12345678910111213141516171819202122232425262728293031323334353637383940414243444546474849505152535455565758596061626364656667686970717273747576777879808182838485868788899091929394959697989910010110210310410510610710810911011111211311411511611711811912012112212312412512612712812913013113213313413513613713813914014114214314414514614714814915015115215315415515615715815916016116216316416516616716816917017117217317417517617717817918018118218318418518618718818919019119219319419519619719819920020120220320420520620720820921021121221321421521621721821922022122222322422522622722822923023123223323423523623723823924024124224324424524624724824925025125225325425525625725825926026126226326426526626726826927027127227327427527627727827928028128228328428528628728828929029129229329429529629729829930030130230330430530630730830931031131231331431531631731831932032132232332432532632732832933033133233333433533633733833934034134234334434534634734834935035135235335435535635735835936036136236336436536636736836937037137237337437537637737837938038138238338438538638738838939039139239339439539639739839940040140240340440540640740840941041141241341441541641741841942042142242342442542642742842943043143243343443543643743843944044144244344444544644744844945045145245345445545645745845946046146246346446546646746846947047147247347447547647747847948048148248348448548648748848949049149249349449549649749849950050150250350450550650750850951051151251351451551651751851952052152252352452552652752852953053153253353453553653753853954054154254354454554654754854955055155255355455555655755855956056156256356456556656756856957057157257357457557657757857958058158258358458558658758858959059159259359459559659759859960060160260360460560660760860961061161261361461561661761861962062162262362462562662762862963063163263363463563663763863964064164264364464564664764864965065165265365465565665765865966066166266366466566666766866967067167267367467567667767867968068168268368468568668768868969069169269369469569669769869970070170270370470570670770870971071171271371471571671771871972072172272372472572672772872973073173273373473573673773873974074174274374474574674774874975075175275375475575675775875976076176276376476576676776876977077177277377477577677777877978078178278378478578678778878979079179279379479579679779879980080180280380480580680780880981081181281381481581681781881982082182282382482582682782882983083183283383483583683783883984084184284384484584684784884985085185285385485585685785885986086186286386486586686786886987087187287387487587687787887988088188288388488588688788888989089189289389489589689789889990090190290390490590690790890991091191291391491591691791891992092192292392492592692792892993093193293393493593693793893994094194294394494594694794894995095195295395495595695795895996096196296396496596696796896997097197297397497597697797897998098198298398498598698798898999099199299399499599699799899910001001100210031004100510061007100810091010101110121013101410151016101710181019102010211022102310241025102610271028102910301031103210331034103510361037103810391040104110421043104410451046104710481049105010511052105310541055105610571058105910601061106210631064106510661067106810691070107110721073107410751076107710781079108010811082108310841085108610871088108910901091109210931094109510961097109810991100110111021103110411051106110711081109111011111112111311141115111611171118111911201121112211231124112511261127112811291130113111321133113411351136113711381139114011411142114311441145114611471148114911501151115211531154115511561157115811591160116111621163116411651166116711681169117011711172117311741175117611771178117911801181118211831184118511861187118811891190119111921193119411951196119711981199120012011202120312041205120612071208120912101211121212131214121512161217121812191220122112221223122412251226122712281229123012311232123312341235123612371238123912401241124212431244124512461247124812491250125112521253125412551256125712581259126012611262126312641265126612671268126912701271127212731274127512761277127812791280128112821283128412851286128712881289129012911292129312941295129612971298129913001301130213031304130513061307130813091310131113121313131413151316131713181319132013211322132313241325132613271328132913301331133213331334133513361337133813391340134113421343134413451346134713481349135013511352135313541355135613571358135913601361136213631364136513661367136813691370137113721373137413751376137713781379138013811382138313841385138613871388138913901391139213931394139513961397139813991400140114021403140414051406140714081409141014111412141314141415141614171418141914201421142214231424142514261427142814291430143114321433143414351436143714381439144014411442144314441445144614471448144914501451145214531454145514561457145814591460146114621463146414651466146714681469147014711472147314741475147614771478147914801481148214831484148514861487148814891490149114921493149414951496149714981499150015011502150315041505150615071508150915101511151215131514151515161517151815191520152115221523152415251526152715281529153015311532153315341535153615371538153915401541154215431544154515461547154815491550155115521553155415551556155715581559156015611562156315641565156615671568156915701571157215731574157515761577157815791580158115821583158415851586158715881589159015911592159315941595159615971598159916001601160216031604160516061607160816091610161116121613161416151616161716181619162016211622162316241625162616271628162916301631163216331634163516361637163816391640164116421643164416451646164716481649165016511652165316541655165616571658165916601661166216631664166516661667166816691670167116721673167416751676167716781679168016811682168316841685168616871688168916901691169216931694169516961697169816991700170117021703170417051706170717081709171017111712171317141715171617171718171917201721172217231724172517261727172817291730173117321733173417351736173717381739174017411742174317441745174617471748174917501751175217531754175517561757175817591760176117621763176417651766176717681769177017711772177317741775177617771778177917801781178217831784178517861787178817891790179117921793179417951796179717981799180018011802180318041805180618071808180918101811181218131814181518161817181818191820182118221823182418251826182718281829183018311832183318341835183618371838183918401841184218431844184518461847184818491850185118521853185418551856185718581859186018611862186318641865186618671868186918701871187218731874187518761877187818791880188118821883188418851886188718881889189018911892189318941895189618971898189919001901190219031904190519061907190819091910191119121913191419151916191719181919192019211922192319241925192619271928192919301931193219331934193519361937193819391940194119421943194419451946194719481949195019511952195319541955195619571958195919601961196219631964196519661967196819691970197119721973197419751976197719781979198019811982198319841985198619871988198919901991199219931994199519961997199819992000200120022003200420052006200720082009201020112012201320142015201620172018201920202021202220232024202520262027202820292030203120322033203420352036203720382039204020412042204320442045204620472048204920502051205220532054205520562057205820592060206120622063206420652066206720682069207020712072207320742075207620772078207920802081208220832084208520862087208820892090209120922093209420952096209720982099210021012102210321042105210621072108210921102111211221132114211521162117211821192120212121222123212421252126212721282129213021312132213321342135213621372138213921402141214221432144214521462147214821492150215121522153215421552156215721582159216021612162216321642165216621672168216921702171217221732174217521762177217821792180218121822183218421852186218721882189219021912192219321942195219621972198219922002201220222032204220522062207220822092210221122122213221422152216221722182219222022212222222322242225222622272228222922302231223222332234223522362237223822392240224122422243224422452246224722482249225022512252225322542255225622572258225922602261226222632264226522662267226822692270227122722273227422752276227722782279228022812282228322842285228622872288228922902291229222932294229522962297229822992300230123022303230423052306230723082309231023112312231323142315231623172318231923202321232223232324232523262327232823292330233123322333233423352336233723382339234023412342234323442345234623472348234923502351235223532354235523562357235823592360236123622363236423652366236723682369237023712372237323742375237623772378237923802381238223832384238523862387238823892390239123922393239423952396239723982399240024012402240324042405240624072408240924102411241224132414241524162417241824192420242124222423242424252426242724282429243024312432243324342435243624372438243924402441244224432444244524462447244824492450245124522453 |
- from PyQt5.QtCore import *
- import numpy as np
- import scipy.signal as signal
- import pylab
- import sys
- import scipy
- import copy
- import struct
- from scipy.io.matlab import mio
- from numpy import pi
- from math import floor
- import matplotlib as mpl
- from matplotlib.ticker import FuncFormatter
- import matplotlib.font_manager as fm
- import matplotlib.pyplot as plt
- import matplotlib.ticker
- from matplotlib.ticker import MaxNLocator
-
- import multiprocessing
- import itertools
-
- import akvo.tressel.adapt as adapt
-
- import akvo.tressel.decay as decay
- import akvo.tressel.pca as pca
- import akvo.tressel.rotate as rotate
- import akvo.tressel.cmaps as cmaps
-
- import cmocean
- plt.register_cmap(name='viridis', cmap=cmaps.viridis)
- plt.register_cmap(name='inferno', cmap=cmaps.inferno)
- plt.register_cmap(name='inferno_r', cmap=cmaps.inferno_r)
-
- plt.register_cmap(name='magma', cmap=cmaps.magma)
- plt.register_cmap(name='magma_r', cmap=cmaps.magma_r)
-
-
- def loadGMRBinaryFID( rawfname, istack, info ):
- """ Reads a single binary GMR file and fills into DATADICT
- """
-
-
-
-
- nps = (int)((info["prePulseDelay"])*info["samp"])
- npul = (int)(self.pulseLength[0]*self.samp)
-
-
- nds = nps+npul+(int)((self.deadTime)*self.samp);
- nd1 = (int)(1.*self.samp)
-
- invGain = 1./self.RxGain
- invCGain = self.CurrentGain
-
- pulse = "Pulse 1"
- chan = self.DATADICT[pulse]["chan"]
- rchan = self.DATADICT[pulse]["rchan"]
-
- rawFile = open( rawfname, 'rb')
-
- T = N_samp * self.dt
- TIMES = np.arange(0, T, self.dt) - .0002
-
- for ipm in range(self.nPulseMoments):
- buf1 = rawFile.read(4)
- buf2 = rawFile.read(4)
-
- N_chan = struct.unpack('>i', buf1 )[0]
- N_samp = struct.unpack('>i', buf2 )[0]
-
- DATA = np.zeros([N_samp, N_chan+1])
- for ichan in range(N_chan):
- DATADUMP = rawFile.read(4*N_samp)
- for irec in range(N_samp):
- DATA[irec,ichan] = struct.unpack('>f', DATADUMP[irec*4:irec*4+4])[0]
-
- return DATA, TIMES
-
- class SNMRDataProcessor(QObject):
- """ Revised class for preprocessing sNMR Data.
- Derived types can read GMR files
- """
- def __init__(self):
- QObject.__init__(self)
- self.numberOfMoments = 0
- self.numberOfPulsesPerMoment = 0
- self.pulseType = "NONE"
- self.transFreq = 0
- self.pulseLength = np.zeros(1)
- self.nPulseMoments = 0
- self.dt = 0
-
- def mfreqz(self, b,a=1):
- """ Plots the frequency response of a filter specified with a and b weights
- """
- import scipy.signal as signal
- pylab.figure(124)
- w,h = signal.freqz(b,a)
- w /= max(w)
- w *= .5/self.dt
- h_dB = 20 * pylab.log10 (abs(h))
- pylab.subplot(211)
-
- pylab.plot(w,h_dB)
- pylab.ylim(-150, 5)
- pylab.ylabel('Magnitude (dB)')
-
- pylab.xlabel(r'Hz')
- pylab.title(r'Frequency response')
- pylab.subplot(212)
- h_Phase = pylab.unwrap(pylab.arctan2(pylab.imag(h), pylab.real(h)))
-
- pylab.plot(w,h_Phase)
- pylab.ylabel('Phase (radians)')
- pylab.xlabel(r'Hz')
-
- pylab.title(r'Phase response')
- pylab.subplots_adjust(hspace=0.5)
-
- def mfreqz2(self, b, a, canvas):
- "for analysing filt-filt"
- import scipy.signal as signal
- canvas.reAx2(False,False)
-
- canvas.ax1.tick_params(axis='both', which='major', labelsize=8)
- canvas.ax2.tick_params(axis='both', which='major', labelsize=8)
-
-
-
- w,h = signal.freqz(b,a)
- w /= max(w)
- w *= .5/self.dt
- h_dB = 20 * pylab.log10(abs(h*h) + 1e-16)
-
-
- canvas.ax1.plot(w,h_dB)
- canvas.ax1.set_ylim(-150, 5)
- canvas.ax1.set_ylabel('Magnitude [db]', fontsize=8)
-
- canvas.ax1.set_xlabel(r'[Hz]', fontsize=8)
- canvas.ax1.set_title(r'Frequency response', fontsize=8)
- canvas.ax1.grid(True)
-
- tt = np.arange(0, .02, self.dt)
- impulse = signal.dimpulse((self.filt_z, self.filt_p, self.filt_k, self.dt), t=tt)
-
-
-
- impulse_dB = 20.*np.log10(np.abs(np.array(impulse[1][0])))
-
- canvas.ax2.plot(np.array(impulse[0]), impulse[1][0])
-
-
- canvas.ax2.set_ylabel('response [%]', fontsize=8)
- canvas.ax2.set_xlabel(r'time [s]', fontsize=8)
- canvas.ax2.set_title(r'impulse response', fontsize=8)
-
- canvas.draw()
-
- return impulse
-
-
- class GMRDataProcessor(SNMRDataProcessor):
-
-
- progressTrigger = pyqtSignal("int")
- doneTrigger = pyqtSignal()
- enableDSPTrigger = pyqtSignal()
- updateProcTrigger = pyqtSignal()
-
- def __init__(self):
-
- SNMRDataProcessor.__init__(self)
- self.maxBusV = 0.
- self.samp = 50000.
- self.dt = 2e-5
- self.deadTime = .0055
- self.prePulseDelay = 0.05
- self.windead = 0.
- self.pulseType = -1
- self.transFreq = -1
- self.maxBusV = -1
- self.pulseLength = -1
- self.interpulseDelay = -1
- self.repetitionDelay = -1
- self.nPulseMoments = -1
- self.TuneCapacitance = -1
- self.nTransVersion = -1
- self.nDAQVersion = -1
- self.nInterleaves = -1
-
-
- self.RotatedAmplitude = False
-
-
-
- def Print(self):
- print ("pulse type", self.pulseType)
- print ("maxBusV", self.maxBusV)
- print ("inner pulse delay", self.interpulseDelay)
- print ("tuning capacitance", self.TuneCapacitance)
- print ("sampling rate", self.samp)
- print ("dt", self.dt)
- print ("dead time", self.deadTime)
- print ("pre pulse delay", self.prePulseDelay)
- print ("number of pulse moments", self.nPulseMoments)
- print ("pulse Length", self.pulseLength)
- print ("trans freq", self.transFreq)
-
- def readHeaderFile(self, FileName):
-
- HEADER = np.loadtxt(FileName)
-
- pulseTypeDict = {
- 1 : lambda: "FID",
- 2 : lambda: "T1",
- 3 : lambda: "SPINECHO",
- 4 : lambda: "4PhaseT1"
- }
-
- pulseLengthDict = {
- 1 : lambda x: np.ones(1) * x,
- 2 : lambda x: np.ones(2) * x,
- 3 : lambda x: np.array([x, 2.*x]),
- 4 : lambda x: np.ones(2) * x
- }
-
- self.pulseType = pulseTypeDict.get((int)(HEADER[0]))()
- self.transFreq = HEADER[1]
- self.maxBusV = HEADER[2]
- self.pulseLength = pulseLengthDict.get((int)(HEADER[0]))(1e-3*HEADER[3])
- self.interpulseDelay = 1e-3*HEADER[4]
- self.repetitionDelay = HEADER[5]
- self.nPulseMoments = (int)(HEADER[6])
- self.TuneCapacitance = HEADER[7]
- self.nTransVersion = HEADER[8]
- self.nDAQVersion = HEADER[9]
- self.nInterleaves = HEADER[10]
-
- self.gain()
-
-
- self.samp = 50000.
- self.dt = 2e-5
-
-
- if self.nDAQVersion >= 2:
-
-
-
- self.samp = HEADER[14]
- self.dt = 1./self.samp
- self.deadTime = .0055
- self.prePulseDelay = 0.05
-
-
- def gain(self):
-
-
-
-
- w = 2*np.pi*self.transFreq
-
- L_coil = 1e6/(self.TuneCapacitance*(w**2))
- R_coil = 1.
- Z1_in = .5 + 1j*.5*w
- Z2_in = 1./(1j*w*.000001616)
- Z_eq_inv = (1./Z1_in) + (1./Z2_in)
- Zeq = 1./Z_eq_inv
- Zsource = R_coil + 1j*w*L_coil
- voltage_in = Zeq / (Zsource + Zeq)
- self.circuitGain = np.abs(voltage_in)
- self.circuitPhase_deg = (180/np.pi)+np.angle(voltage_in)
- circuitImpedance_ohms = np.abs(Zsource + Zeq)
-
-
-
- if self.nTransVersion == 4:
- self.PreAmpGain = 1000.
- elif self.nTransVersion == 1 or self.nTransVersion == 2 or self.nTransVersion == 3 or self.nTransVersion == 6:
- self.PreAmpGain = 500.
- else:
- print ("unsupported transmitter version")
- exit(1)
-
-
- self.RxGain = self.circuitGain * self.PreAmpGain
-
-
-
- if floor(self.nDAQVersion) == 1:
- self.CurrentGain = 150.
- elif floor(self.nDAQVersion) == 2:
- self.CurrentGain = 180.
-
- def updateProgress(self):
- pass
-
- def TDSmartStack(self, outlierTest, MADcutoff, canvas):
-
- Stack = {}
-
- for pulse in self.DATADICT["PULSES"]:
- stack = np.zeros(( len(self.DATADICT[pulse]["chan"]), self.DATADICT["nPulseMoments"],\
- len(self.DATADICT["stacks"]), len(self.DATADICT[pulse]["TIMES"]) ))
- for ipm in range(self.DATADICT["nPulseMoments"]):
- istack = 0
- for sstack in self.DATADICT["stacks"]:
- if self.pulseType == "FID" or pulse == "Pulse 2":
- mod = (-1)**(ipm%2) * (-1)**(sstack%2)
- elif self.pulseType == "4PhaseT1":
- mod = (-1)**(ipm%2) * (-1)**(((sstack-1)/2)%2)
- ichan = 0
- for chan in self.DATADICT[pulse]["chan"]:
- stack[ichan,ipm,istack,:] += mod*self.DATADICT[pulse][chan][ipm][sstack]
- ichan += 1
- istack += 1
- Stack[pulse] = stack
-
-
-
-
- canvas.reAxH2(np.shape(stack)[0], False, False)
- axes = canvas.fig.axes
- SimpleStack = {}
- VarStack = {}
- for pulse in self.DATADICT["PULSES"]:
- SimpleStack[pulse] = {}
- VarStack[pulse] = {}
- ichan = 0
- for chan in self.DATADICT[pulse]["chan"]:
- SimpleStack[pulse][chan] = 1e9*np.average( Stack[pulse][ichan], 1 )
- VarStack[pulse][chan] = 1e9*np.std( Stack[pulse][ichan], 1 )
- ax1 = axes[ 2*ichan ]
-
-
- y_formatter = matplotlib.ticker.ScalarFormatter(useOffset=False)
- ax1.yaxis.set_major_formatter(y_formatter)
-
- ax1.plot( 1e3*self.DATADICT[pulse]["TIMES"], np.average( SimpleStack[pulse][chan], 0 ))
- ax1.set_title("Ch." + str(chan) + ": avg FID", fontsize=8)
- ax1.set_xlabel(r"time (ms)", fontsize=8)
-
- if ichan == 0:
- ax1.set_ylabel(r"signal [nV]", fontsize=8)
- else:
- plt.setp(ax1.get_yticklabels(), visible=False)
- plt.setp(ax1.get_yaxis().get_offset_text(), visible=False)
-
-
-
-
-
-
-
-
-
-
-
-
- ichan += 1
-
-
-
-
- if outlierTest == "MAD":
- MADStack = {}
- VarStack = {}
-
- madstack = np.zeros(( len(self.DATADICT[pulse]["chan"]),\
- self.DATADICT["nPulseMoments"], len(self.DATADICT[pulse]["TIMES"]) ))
- varstack = np.zeros(( len(self.DATADICT[pulse]["chan"]),\
- self.DATADICT["nPulseMoments"], len(self.DATADICT[pulse]["TIMES"]) ))
- for pulse in self.DATADICT["PULSES"]:
- MADStack[pulse] = {}
- VarStack[pulse] = {}
- ichan = 0
- for chan in self.DATADICT[pulse]["chan"]:
- ax1 = axes[ 2*ichan ]
- for ipm in range(self.DATADICT["nPulseMoments"]):
-
-
-
-
-
-
-
-
-
-
-
-
-
-
-
-
- x = 1e9*copy.deepcopy(Stack[pulse][ichan][ipm,:,:])
- tile_med = np.tile( np.median(x, axis=0), (np.shape(x)[0],1))
- MAD = MADcutoff * np.median(np.abs(x - tile_med), axis=0)
- tile_MAD = np.tile( MAD, (np.shape(x)[0],1))
- good = np.abs(x-tile_med)/tile_MAD < 2.
-
- madstack[ichan][ipm] = copy.deepcopy( np.ma.masked_array(x, good != True).mean(axis=0) )
- varstack[ichan][ipm] = copy.deepcopy( np.ma.masked_array(x, good != True).std(axis=0) )
-
-
- percent = int(1e2* (float)((ipm)+ichan*self.DATADICT["nPulseMoments"]) /
- (float)(self.DATADICT["nPulseMoments"] * len(self.DATADICT[pulse]["chan"])))
- self.progressTrigger.emit(percent)
-
- ax1.plot( 1e3*self.DATADICT[pulse]["TIMES"], np.average( madstack[ichan], 0 ))
-
- MADStack[pulse][chan] = madstack[ichan]
- VarStack[pulse][chan] = varstack[ichan]
- ichan += 1
-
- self.DATADICT["stack"] = MADStack
-
- else:
- self.DATADICT["stack"] = SimpleStack
-
-
-
-
-
-
-
-
- im2 = []
- im1 = []
- for pulse in self.DATADICT["PULSES"]:
- ichan = 0
- axes = canvas.fig.axes
- vvmin = 1e10
- vvmax = 0
- for chan in self.DATADICT[pulse]["chan"]:
- ax1 = axes[2*ichan ]
- ax2 = axes[2*ichan+1]
- if outlierTest == "MAD":
- X = np.fft.rfft( MADStack[pulse][chan][0,:] )
- nu = np.fft.fftfreq(len( MADStack[pulse][chan][0,:]), d=self.dt)
- else:
- X = np.fft.rfft( SimpleStack[pulse][chan][0,:] )
- nu = np.fft.fftfreq(len( SimpleStack[pulse][chan][0,:]), d=self.dt)
-
- nu = nu[0:len(X)]
- nu[-1] = np.abs(nu[-1])
- df = nu[1] - nu[0]
- of = 0
-
- istart = int((self.transFreq-50.)/df)
- iend = int((self.transFreq+50.)/df)
- of = nu[istart]
-
- def freqlabel(xxx, pos):
- return '%1.0f' %(of + xxx*df)
- formatter = FuncFormatter(freqlabel)
-
- SFFT = np.zeros( (self.DATADICT["nPulseMoments"], len(X)), dtype=np.complex64 )
- SFFT[0,:] = X
- for ipm in range(1, self.DATADICT["nPulseMoments"]):
- if outlierTest == "MAD":
- SFFT[ipm,:] = np.fft.rfft( MADStack[pulse][chan][ipm,:] )
- else:
- SFFT[ipm,:] = np.fft.rfft( SimpleStack[pulse][chan][ipm,:] )
-
-
-
- db = (np.abs(SFFT[:,istart:iend]))
-
-
-
-
- vvmin = min(vvmin, np.min (db))
- vvmax = max(vvmax, np.max (db))
- im2.append(ax2.matshow( db, aspect='auto', cmap=cmocean.cm.ice, vmin=vvmin, vmax=vvmax))
-
-
-
- if ichan == 0:
- ax2.set_ylabel(r"$q$ (A $\cdot$ s)", fontsize=8)
- else:
-
- plt.setp(ax2.get_yticklabels(), visible=False)
-
- ax2.xaxis.set_major_formatter(formatter)
- ax2.xaxis.set_ticks_position('bottom')
- ax2.xaxis.set_major_locator(MaxNLocator(3))
-
- y_formatter = matplotlib.ticker.ScalarFormatter(useOffset=False)
- ax2.yaxis.set_major_formatter(y_formatter)
-
-
-
-
-
-
-
-
-
-
- ichan += 1
-
-
- canvas.fig.subplots_adjust(hspace=.1, wspace=.05, left=.075, right=.95 )
-
-
-
-
-
- cb2 = canvas.fig.colorbar(im2[-1], ax=axes[1::2], format='%1.0e', orientation='horizontal', shrink=.35, aspect=30)
- cb2.ax.tick_params(axis='both', which='major', labelsize=8)
- cb2.set_label("$\mathcal{V}_N$ (nV)", fontsize=8)
-
-
- canvas.draw()
- self.doneTrigger.emit()
-
- def FDSmartStack(self, outlierTest, MADcutoff, canvas):
-
- print("FFT stuff")
- self.dataCubeFFT()
-
- Stack = {}
-
- for pulse in self.DATADICT["PULSES"]:
- stack = np.zeros(( len(self.DATADICT[pulse]["chan"]), \
- self.DATADICT["nPulseMoments"],\
- len(self.DATADICT["stacks"]),\
- len(self.DATADICT[pulse][self.DATADICT[pulse]["chan"][0] ]["FFT"]["nu"])//2 + 1),\
- dtype=np.complex )
- for ipm in range(self.DATADICT["nPulseMoments"]):
- istack = 0
- for sstack in self.DATADICT["stacks"]:
- if self.pulseType == "FID" or pulse == "Pulse 2":
- mod = (-1)**(ipm%2) * (-1)**(sstack%2)
- elif self.pulseType == "4PhaseT1":
- mod = (-1)**(ipm%2) * (-1)**(((sstack-1)/2)%2)
- ichan = 0
- for chan in self.DATADICT[pulse]["chan"]:
-
- stack[ichan,ipm,istack,:] += mod*self.DATADICT[pulse][chan]["FFT"][sstack][ipm,:]
- ichan += 1
- istack += 1
- Stack[pulse] = stack
-
-
-
-
- canvas.reAxH2(np.shape(stack)[0], False, False)
- axes = canvas.fig.axes
- SimpleStack = {}
- VarStack = {}
- for pulse in self.DATADICT["PULSES"]:
- SimpleStack[pulse] = {}
- VarStack[pulse] = {}
- ichan = 0
- for chan in self.DATADICT[pulse]["chan"]:
- SimpleStack[pulse][chan] = 1e9*np.average( Stack[pulse][ichan], 1 )
- VarStack[pulse][chan] = 1e9*np.std( Stack[pulse][ichan], 1 )
- ax1 = axes[ 2*ichan ]
-
-
- y_formatter = matplotlib.ticker.ScalarFormatter(useOffset=False)
- ax1.yaxis.set_major_formatter(y_formatter)
-
-
-
- ax1.matshow( np.real(SimpleStack[pulse][chan]), aspect='auto')
- ax1.set_title("Ch." + str(chan) + ": avg FID", fontsize=8)
- ax1.set_xlabel(r"time (ms)", fontsize=8)
-
- if ichan == 0:
- ax1.set_ylabel(r"signal [nV]", fontsize=8)
- else:
- plt.setp(ax1.get_yticklabels(), visible=False)
- plt.setp(ax1.get_yaxis().get_offset_text(), visible=False)
- ichan += 1
-
-
-
-
- if outlierTest == "MAD":
- MADStack = {}
- VarStack = {}
-
- madstack = np.zeros(( len(self.DATADICT[pulse]["chan"]),\
- self.DATADICT["nPulseMoments"],\
- len(self.DATADICT[pulse][self.DATADICT[pulse]["chan"][0] ]["FFT"]["nu"])//2 + 1))
- varstack = np.zeros(( len(self.DATADICT[pulse]["chan"]),\
- self.DATADICT["nPulseMoments"],\
- len(self.DATADICT[pulse][self.DATADICT[pulse]["chan"][0] ]["FFT"]["nu"])//2 + 1))
- for pulse in self.DATADICT["PULSES"]:
- MADStack[pulse] = {}
- VarStack[pulse] = {}
- ichan = 0
- for chan in self.DATADICT[pulse]["chan"]:
- ax1 = axes[ 2*ichan ]
- for ipm in range(self.DATADICT["nPulseMoments"]):
-
-
-
-
-
-
-
-
-
-
-
-
-
-
-
-
- x = 1e9*copy.deepcopy(Stack[pulse][ichan][ipm,:,:])
- tile_med = np.tile( np.median(x, axis=0), (np.shape(x)[0],1))
- MAD = MADcutoff * np.median(np.abs(x - tile_med), axis=0)
- tile_MAD = np.tile( MAD, (np.shape(x)[0],1))
- good = np.abs(x-tile_med)/tile_MAD < 2.
-
- madstack[ichan][ipm] = copy.deepcopy( np.ma.masked_array(x, good != True).mean(axis=0) )
- varstack[ichan][ipm] = copy.deepcopy( np.ma.masked_array(x, good != True).std(axis=0) )
-
-
- percent = int(1e2* (float)((ipm)+ichan*self.DATADICT["nPulseMoments"]) /
- (float)(self.DATADICT["nPulseMoments"] * len(self.DATADICT[pulse]["chan"])))
- self.progressTrigger.emit(percent)
-
- ax2 = axes[2*ichan+1]
-
- MADStack[pulse][chan] = madstack[ichan]
- VarStack[pulse][chan] = varstack[ichan]
- ax2.matshow( np.real(MADStack[pulse][chan]), aspect='auto')
- ichan += 1
-
- self.DATADICT["stack"] = MADStack
-
- else:
- self.DATADICT["stack"] = SimpleStack
-
-
-
-
-
-
-
-
-
-
-
-
-
-
-
-
-
-
-
-
-
-
-
-
-
-
-
-
-
-
-
-
-
-
-
-
-
-
-
-
-
-
-
-
-
-
-
-
-
-
-
-
-
-
-
-
-
-
-
-
-
-
-
-
-
-
-
-
-
-
-
-
-
-
-
-
-
-
-
-
-
-
-
-
-
-
-
-
-
-
-
-
-
-
-
- canvas.draw()
- self.doneTrigger.emit()
-
-
- def sumData(self, canvas, fred):
- chans = copy.deepcopy(self.DATADICT[self.DATADICT["PULSES"][0]]["chan"])
- nchan = len(chans)
-
- for ich in range(nchan-1):
- for ch in chans[ich+1:]:
- chsum = chans[ich] + "+" + ch
- for pulse in self.DATADICT["PULSES"]:
- self.DATADICT[pulse][chsum] = {}
- for ipm in range(self.DATADICT["nPulseMoments"]):
- self.DATADICT[pulse][chsum][ipm] = {}
- for istack in self.DATADICT["stacks"]:
- self.DATADICT[pulse][chsum][ipm][istack] = self.DATADICT[pulse][chans[ich]][ipm][istack] - self.DATADICT[pulse][ch][ipm][istack]
- if chsum == "1+2":
-
-
- self.DATADICT[pulse]["chan"].append(chsum)
-
-
- sumall = False
- if sumall:
- chsum = ""
- for ch in chans:
- chsum += ch + "+"
- chsum = chsum[0:-1]
-
- for pulse in self.DATADICT["PULSES"]:
- self.DATADICT[pulse][chsum] = {}
- for ipm in range(self.DATADICT["nPulseMoments"]):
- self.DATADICT[pulse][chsum][ipm] = {}
- for istack in self.DATADICT["stacks"]:
- self.DATADICT[pulse][chsum][ipm][istack] = copy.deepcopy(self.DATADICT[pulse][chans[0]][ipm][istack])
- for ch in chans[1:]:
- self.DATADICT[pulse][chsum][ipm][istack] += self.DATADICT[pulse][ch][ipm][istack]
- self.DATADICT[pulse]["chan"].append(chsum)
-
-
-
-
-
-
-
-
-
- self.doneTrigger.emit()
-
- def quadDet(self, clip, phase, canvas):
-
- from scipy import signal
- self.RotatedAmplitude = True
-
- wL = self.transFreq * 2*np.pi
- vL = self.transFreq
-
- dt = self.dt
-
-
- CA = {}
- IP = {}
- NR = {}
- RE = {}
- IM = {}
-
-
- CAmax = {}
- NRmax = {}
- REmax = {}
- IMmax = {}
- E0,phi,df,T2 = 100.,0,0,.2
- first = False
- self.sigma = {}
- for pulse in self.DATADICT["PULSES"]:
- CA[pulse] = {}
- IP[pulse] = {}
- NR[pulse] = {}
- RE[pulse] = {}
- IM[pulse] = {}
- CAmax[pulse] = 0
- NRmax[pulse] = 0
- REmax[pulse] = 0
- IMmax[pulse] = 0
- ichan = 0
- self.sigma[pulse] = {}
- for chan in self.DATADICT[pulse]["chan"]:
- CA[pulse][chan] = np.zeros( (self.DATADICT["nPulseMoments"], len(self.DATADICT[pulse]["TIMES"])-clip ) )
- IP[pulse][chan] = np.zeros( (self.DATADICT["nPulseMoments"], len(self.DATADICT[pulse]["TIMES"])-clip ) )
- NR[pulse][chan] = np.zeros( (self.DATADICT["nPulseMoments"], len(self.DATADICT[pulse]["TIMES"])-clip ) )
- RE[pulse][chan] = np.zeros( (self.DATADICT["nPulseMoments"], len(self.DATADICT[pulse]["TIMES"])-clip ) )
- IM[pulse][chan] = np.zeros( (self.DATADICT["nPulseMoments"], len(self.DATADICT[pulse]["TIMES"])-clip ) )
- for ipm in range(0, self.DATADICT["nPulseMoments"]):
-
- xn = self.DATADICT["stack"][pulse][chan][ipm,:]
- ht = signal.hilbert(xn)*np.exp(-1j*wL*self.DATADICT[pulse]["TIMES"])
-
-
- RE[pulse][chan][ipm,:] = np.real(ht[clip::])
- IM[pulse][chan][ipm,:] = np.imag(ht[clip::])
- REmax[pulse] = max(REmax[pulse], np.max(np.real(ht[clip::])))
- IMmax[pulse] = max(IMmax[pulse], np.max(np.imag(ht[clip::])))
-
-
- IP[pulse][chan][ipm,:] = np.angle(ht)[clip::]
-
-
-
-
-
-
- [success, E0, df, phi, T2] = decay.quadratureDetect2( ht.real, ht.imag, self.DATADICT[pulse]["TIMES"])
-
-
-
-
- D = self.RotateAmplitude( ht.real, ht.imag, phi, df, self.DATADICT[pulse]["TIMES"] )
- CA[pulse][chan][ipm,:] = D.imag[clip::]
- NR[pulse][chan][ipm,:] = D.real[clip::]
- CAmax[pulse] = max(CAmax[pulse], np.max(D.imag[clip::]) )
- NRmax[pulse] = max(NRmax[pulse], np.max(D.real[clip::]) )
- self.sigma[pulse][chan] = np.std(NR[pulse][chan])
-
- percent = int(1e2* (float)((ipm)+ichan*self.DATADICT["nPulseMoments"]) /
- (float)(self.DATADICT["nPulseMoments"] * len(self.DATADICT[pulse]["chan"])))
- self.progressTrigger.emit(percent)
- ichan += 1
-
- self.DATADICT["CA"] = CA
- self.DATADICT["IP"] = IP
- self.DATADICT["NR"] = NR
- self.DATADICT["RE"] = RE
- self.DATADICT["IM"] = IM
-
- self.DATADICT["CAmax"] = CAmax
- self.DATADICT["NRmax"] = NRmax
- self.DATADICT["REmax"] = REmax
- self.DATADICT["IMmax"] = IMmax
-
- self.doneTrigger.emit()
-
- def plotQuadDet(self, clip, phase, canvas):
-
- canvas.reAxH2( len(self.DATADICT[ self.DATADICT["PULSES"][0] ]["chan"] ), False, False)
-
-
-
-
- dcmap = cmocean.cm.curl_r
- canvas.reAxH2( len(self.DATADICT[ self.DATADICT["PULSES"][0] ]["chan"] ), False, False)
- for pulse in self.DATADICT["PULSES"]:
- ichan = 0
- axes = canvas.fig.axes
- mmaxr = 0.
- mmaxi = 0.
- if clip > 0:
- time_sp = 1e3 * (self.DATADICT[pulse]["TIMES"][clip-1::] - self.DATADICT[pulse]["PULSE_TIMES"][-1] )
- else:
- time_sp = 1e3 * (self.DATADICT[pulse]["TIMES"] - self.DATADICT[pulse]["PULSE_TIMES"][-1] )
-
- QQ = np.average(self.DATADICT[pulse]["Q"], axis=1 )
-
- for chan in self.DATADICT[pulse]["chan"]:
- ax1 = axes[2*ichan ]
- ax2 = axes[2*ichan+1]
- if phase == 0:
- im1 = ax1.pcolormesh( time_sp, QQ, self.DATADICT["RE"][pulse][chan], cmap=dcmap, rasterized=True,\
- vmin=-self.DATADICT["REmax"][pulse] , vmax=self.DATADICT["REmax"][pulse] )
- im2 = ax2.pcolormesh( time_sp, QQ, self.DATADICT["IM"][pulse][chan], cmap=dcmap, rasterized=True,\
- vmin=-self.DATADICT["IMmax"][pulse], vmax=self.DATADICT["IMmax"][pulse] )
- if phase == 1:
- im1 = ax1.pcolormesh( time_sp, QQ, self.DATADICT["CA"][pulse][chan], cmap=dcmap, rasterized=True,
- vmin=-self.DATADICT["CAmax"][pulse] , vmax=self.DATADICT["CAmax"][pulse] )
- im2 = ax2.pcolormesh( time_sp, QQ, self.DATADICT["IP"][pulse][chan], cmap=cmocean.cm.phase, rasterized=True,\
- vmin=-np.pi, vmax=np.pi)
- if phase == 2:
- im1 = ax1.pcolormesh( time_sp, QQ, self.DATADICT["CA"][pulse][chan], cmap=dcmap, rasterized=True,\
- vmin=-self.DATADICT["CAmax"][pulse] , vmax=self.DATADICT["CAmax"][pulse] )
- im2 = ax2.pcolormesh( time_sp, QQ, self.DATADICT["NR"][pulse][chan], cmap=dcmap, rasterized=True,\
- vmin=-self.DATADICT["NRmax"][pulse] , vmax=self.DATADICT["NRmax"][pulse] )
-
-
-
-
-
-
-
-
-
-
-
-
-
-
-
-
-
-
-
-
- percent = int(1e2* (float)(ichan)/len(self.DATADICT[pulse]["chan"]))
- self.progressTrigger.emit(percent)
-
- if ichan == 0:
- ax1.set_ylabel(r"$q$ ( $\mathrm{A}\cdot\mathrm{s}$)", fontsize=8)
- ax2.set_ylabel(r"$q$ ( $\mathrm{A}\cdot\mathrm{s}$)", fontsize=8)
- else:
-
-
- plt.setp(ax1.get_yticklabels(), visible=False)
- plt.setp(ax2.get_yticklabels(), visible=False)
- ichan += 1
-
- ax1.set_yscale('log')
- ax2.set_yscale('log')
-
- plt.setp(ax1.get_xticklabels(), visible=False)
-
- ax1.set_ylim( np.min(QQ), np.max(QQ) )
- ax2.set_ylim( np.min(QQ), np.max(QQ) )
-
- ax1.set_xlim( np.min(time_sp), np.max(time_sp) )
- ax2.set_xlim( np.min(time_sp), np.max(time_sp) )
-
-
- ax2.set_xlabel(r"Time (ms)", fontsize=8)
-
- canvas.fig.subplots_adjust(hspace=.15, wspace=.05, left=.075, right=.95, bottom=.1, top=.95 )
- tick_locator = MaxNLocator(nbins=3)
-
- cb1 = canvas.fig.colorbar(im1, ax=axes[0::2], format='%1.0e', orientation='horizontal', shrink=.35, aspect=30)
- cb1.ax.tick_params(axis='both', which='major', labelsize=8)
- cb1.set_label("$\mathcal{V}_N$ (nV)", fontsize=8)
- cb1.locator = tick_locator
- cb1.update_ticks()
-
- tick_locator2 = MaxNLocator(nbins=3)
- cb2 = canvas.fig.colorbar(im2, ax=axes[1::2], format='%1.0e', orientation='horizontal', shrink=.35, aspect=30, pad=.2)
- cb2.ax.tick_params(axis='both', which='major', labelsize=8)
- cb2.set_label("$\mathcal{V}_N$ (nV)", fontsize=8)
-
-
- cb2.locator = tick_locator2
- cb2.update_ticks()
-
- canvas.draw()
- self.doneTrigger.emit()
-
- def RotateAmplitude(self, X, Y, zeta, df, t):
- V = X + 1j*Y
- return np.abs(V) * np.exp( 1j * ( np.angle(V) - zeta - 2.*np.pi*df*t ) )
-
- def gateIntegrate( self, gpd, clip, canvas ):
- """ Gate integrate the real, imaginary, phased, and noise residual channels
- """
-
- self.gated = True
- self.GATED = {}
-
- for pulse in self.DATADICT["PULSES"]:
- QQ = np.average(self.DATADICT[pulse]["Q"], axis=1 )
- ichan = 0
- for chan in self.DATADICT[pulse]["chan"]:
- self.GATED[chan] = {}
- for ipm in range(0, self.DATADICT["nPulseMoments"]):
-
-
-
-
-
- time_sp = 1e3 * (self.DATADICT[pulse]["TIMES"] - self.DATADICT[pulse]["PULSE_TIMES"][-1] )
-
-
-
-
- GT, GCA, GTT, sig_stack, isum = rotate.gateIntegrate( self.DATADICT["CA"][pulse][chan][ipm], time_sp, gpd, self.sigma[pulse][chan], 1.5 )
- GT, GNR, GTT, sig_stack, isum = rotate.gateIntegrate( self.DATADICT["NR"][pulse][chan][ipm], time_sp, gpd, self.sigma[pulse][chan], 1.5 )
- GT, GRE, GTT, sig_stack, isum = rotate.gateIntegrate( self.DATADICT["RE"][pulse][chan][ipm], time_sp, gpd, self.sigma[pulse][chan], 1.5 )
- GT, GIM, GTT, sig_stack, isum = rotate.gateIntegrate( self.DATADICT["IM"][pulse][chan][ipm], time_sp, gpd, self.sigma[pulse][chan], 1.5 )
- GT, GIP, GTT, sig_stack, isum = rotate.gateIntegrate( self.DATADICT["IP"][pulse][chan][ipm], time_sp, gpd, self.sigma[pulse][chan], 1.5 )
-
- if ipm == 0:
-
-
-
- self.GATED[chan]["CA"] = np.zeros( ( self.DATADICT["nPulseMoments"], len(GT)-clip) )
- self.GATED[chan]["NR"] = np.zeros( ( self.DATADICT["nPulseMoments"], len(GT)-clip) )
- self.GATED[chan]["RE"] = np.zeros( ( self.DATADICT["nPulseMoments"], len(GT)-clip) )
- self.GATED[chan]["IM"] = np.zeros( ( self.DATADICT["nPulseMoments"], len(GT)-clip) )
- self.GATED[chan]["IP"] = np.zeros( ( self.DATADICT["nPulseMoments"], len(GT)-clip) )
- self.GATED[chan]["isum"] = isum
-
-
- self.GATEDABSCISSA = GT[clip:]
- self.GATEDWINDOW = GTT[clip:]
-
-
-
- self.GATED[chan]["CA"][ipm] = GCA.real[clip:]
- self.GATED[chan]["NR"][ipm] = GNR.real[clip:]
- self.GATED[chan]["RE"][ipm] = GRE.real[clip:]
- self.GATED[chan]["IM"][ipm] = GIM.real[clip:]
- self.GATED[chan]["IP"][ipm] = GIP.real[clip:]
-
- percent = int(1e2* (float)((ipm)+ichan*self.DATADICT["nPulseMoments"]) /
- (float)(self.DATADICT["nPulseMoments"] * len(self.DATADICT[pulse]["chan"])))
- self.progressTrigger.emit(percent)
-
- self.GATED[chan]["GTT"] = GTT[clip:]
- self.GATED[chan]["GT"] = GT[clip:]
- self.GATED[chan]["QQ"] = QQ
- ichan += 1
- self.doneTrigger.emit()
-
- def bootstrap_resample(self, X, n=None):
-
- """ Bootstrap resample an array_like
- Parameters
- ----------
- X : array_like
- data to resample
- n : int, optional
- length of resampled array, equal to len(X) if n==None
- Results
- -------
- returns X_resamples
- """
- if n == None:
- n = len(X)
-
- resample_i = np.floor(np.random.rand(n)*len(X)).astype(int)
- return X[resample_i]
-
- def bootstrap_sigma(self, pulse, chan):
-
-
- nt = len(self.GATED[chan]["GT"])
- nb = 5000
- XS = np.zeros( (nb, nt) )
- for ii in range(nb):
- for it in range(nt):
- if self.GATED[chan]["isum"][it] < 8:
- XS[ii, it] = self.sigma[pulse][chan]
- else:
-
- if it == 0:
- X = self.bootstrap_resample( np.concatenate( (self.GATED[chan]["NR"][:,it], self.GATED[chan]["NR"][:,it+1], \
- self.GATED[chan]["NR"][:,it+2], self.GATED[chan]["NR"][:,it+3] ) ), n=nt )
- elif it == 1:
- X = self.bootstrap_resample( np.concatenate( (self.GATED[chan]["NR"][:,it-1], self.GATED[chan]["NR"][:,it], \
- self.GATED[chan]["NR"][:,it+1], self.GATED[chan]["NR"][:,it+2] ) ), n=nt )
- elif it == nt-2:
- X = self.bootstrap_resample( np.concatenate( (self.GATED[chan]["NR"][:,it+1], self.GATED[chan]["NR"][:,it], \
- self.GATED[chan]["NR"][:,it-1], self.GATED[chan]["NR"][:,it-2] ) ), n=nt )
- elif it == nt-1:
- X = self.bootstrap_resample( np.concatenate( (self.GATED[chan]["NR"][:,it], self.GATED[chan]["NR"][:,it-1], \
- self.GATED[chan]["NR"][:,it-2], self.GATED[chan]["NR"][:,it-3] ) ), n=nt )
- else:
- X = self.bootstrap_resample( np.concatenate( (self.GATED[chan]["NR"][:,it-2] , self.GATED[chan]["NR"][:,it-1], \
- self.GATED[chan]["NR"][:,it], self.GATED[chan]["NR"][:,it+1], self.GATED[chan]["NR"][:,it+2] )), n=nt )
- XS[ii,it] = np.std(X)
- return XS
-
- def plotGateIntegrate( self, gpd, clip, phase, canvas ):
- """ Plot the gate integration
- """
-
- canvas.reAxH2( len(self.DATADICT[ self.DATADICT["PULSES"][0] ]["chan"] ), False, False)
- axes = canvas.fig.axes
- cmap = cmocean.cm.balance_r
-
-
- vmax1 = 0
- vmax2 = 0
- for pulse in self.DATADICT["PULSES"]:
- for chan in self.DATADICT[pulse]["chan"]:
- if phase == 0:
- vmax1 = max(vmax1, np.max(np.abs(self.GATED[chan]["RE"])))
- vmax2 = max(vmax2, np.max(np.abs(self.GATED[chan]["IM"])))
- elif phase == 1:
- vmax1 = max(vmax1, np.max(np.abs(self.GATED[chan]["CA"])))
- vmax2 = np.pi
- elif phase == 2:
- vmax1 = max(vmax1, np.max(np.abs(self.GATED[chan]["CA"])))
- vmax2 = max(vmax2, np.max(np.abs(self.GATED[chan]["NR"])))
-
-
- for pulse in self.DATADICT["PULSES"]:
- ichan = 0
- for chan in self.DATADICT[pulse]["chan"]:
-
- ax1 = axes[2*ichan ]
- ax2 = axes[2*ichan+1]
-
- if phase == 0:
- im1 = ax1.pcolormesh(self.GATED[chan]["GTT"], self.GATED[chan]["QQ"], self.GATED[chan]["RE"], cmap=cmap, vmin=-vmax1, vmax=vmax1)
- im2 = ax2.pcolormesh(self.GATED[chan]["GTT"], self.GATED[chan]["QQ"], self.GATED[chan]["IM"], cmap=cmap, vmin=-vmax2, vmax=vmax2)
- elif phase == 1:
- im1 = ax1.pcolormesh(self.GATED[chan]["GTT"], self.GATED[chan]["QQ"], self.GATED[chan]["CA"], cmap=cmap, vmin=-vmax1, vmax=vmax1)
- im2 = ax2.pcolormesh(self.GATED[chan]["GTT"], self.GATED[chan]["QQ"], self.GATED[chan]["IP"], cmap=cmocean.cm.phase, vmin=-vmax2, vmax=vmax2)
- elif phase == 2:
- im1 = ax1.pcolormesh(self.GATED[chan]["GTT"], self.GATED[chan]["QQ"], self.GATED[chan]["CA"], cmap=cmap, vmin=-vmax1, vmax=vmax1)
- XS = self.bootstrap_sigma(pulse, chan)
-
-
-
-
-
-
-
-
-
-
-
-
-
-
-
-
-
-
-
-
-
-
-
-
-
-
-
-
-
-
-
- ax2.plot( self.GATED[chan]["GT"], np.std(self.GATED[chan]["NR"], axis=0), color='darkgrey', linewidth=2, label="gate std" )
- ax2.plot( np.tile(self.GATED[chan]["GT"], (5000,1) ), XS, '.', color='lightgrey', linewidth=1, alpha=.5 )
- ax2.plot( self.GATED[chan]["GT"], np.average(XS, axis=0), color='black', linewidth=2, label="bootstrap avg." )
- ax2.plot( self.GATED[chan]["GT"], np.median(XS, axis=0), color='black', linewidth=2, label="bootstrap med." )
- ax2.legend()
-
- im1.set_edgecolor('face')
- if phase != 2:
- im2.set_edgecolor('face')
-
-
- plt.setp(ax1.get_xticklabels(), visible=False)
-
- ax1.set_ylim( np.min(self.GATED[chan]["QQ"]), np.max(self.GATED[chan]["QQ"]) )
-
- if phase != 2:
- ax2.set_ylim( np.min(self.GATED[chan]["QQ"]), np.max(self.GATED[chan]["QQ"]) )
-
- ax1.set_xlim( np.min(self.GATED[chan]["GTT"]), np.max(self.GATED[chan]["GTT"]) )
- ax2.set_xlim( np.min(self.GATED[chan]["GTT"]), np.max(self.GATED[chan]["GTT"]) )
-
- ax1.set_yscale('log')
- ax2.set_yscale('log')
-
- qlabs = np.append(np.concatenate( ( self.GATED[chan]["QQ"][0:1], self.GATED[chan]["QQ"][9::10] )), self.GATED[chan]["QQ"][-1] )
- ax1.yaxis.set_ticks( qlabs )
- if phase != 2:
- ax2.yaxis.set_ticks( qlabs )
-
- formatter = matplotlib.ticker.FuncFormatter(lambda x, pos: str((round(x,1))))
-
- ax1.yaxis.set_major_formatter(formatter)
- ax2.yaxis.set_major_formatter(formatter)
-
- ax1.xaxis.set_major_formatter(formatter)
- ax2.xaxis.set_major_formatter(formatter)
-
- ax1.set_xscale('log')
- ax2.set_xscale('log')
-
- if ichan == 0:
- ax1.set_ylabel(r"$q$ ( $\mathrm{A}\cdot\mathrm{s}$)", fontsize=8)
- if phase == 2:
- ax2.set_ylabel(r"noise est. (nV)", fontsize=8)
- else:
- ax2.set_ylabel(r"$q$ ( $\mathrm{A}\cdot\mathrm{s}$)", fontsize=8)
- else:
- plt.setp(ax1.get_yticklabels(), visible=False)
- plt.setp(ax2.get_yticklabels(), visible=False)
-
- ax2.set_xlabel(r"$t-\tau_p$ (ms)", fontsize=8)
- ichan += 1
-
-
-
-
- canvas.fig.subplots_adjust(hspace=.15, wspace=.05, left=.075, right=.95, bottom=.1, top=.95 )
- tick_locator = MaxNLocator(nbins=5)
-
- cb1 = canvas.fig.colorbar(im1, ax=axes[0::2], format='%1.0e', orientation='horizontal', shrink=.35, aspect=30)
- cb1.ax.tick_params(axis='both', which='major', labelsize=8)
- cb1.set_label("$\mathcal{V}_N$ (nV)", fontsize=8)
- cb1.locator = tick_locator
- cb1.update_ticks()
-
- if phase != 2:
- cb2 = canvas.fig.colorbar(im2, ax=axes[1::2], format='%1.0e', orientation='horizontal', shrink=.35, aspect=30, pad=.2)
- cb2.ax.tick_params(axis='both', which='major', labelsize=8)
- cb2.set_label("$\mathcal{V}_N$ (nV)", fontsize=8)
-
- cb2.locator = tick_locator
- cb2.update_ticks()
-
- canvas.draw()
- self.doneTrigger.emit()
-
-
- def FDSmartStack(self, cv, canvas):
- from matplotlib.colors import LogNorm
- from matplotlib.ticker import MaxNLocator
- """
- Currently this stacks 4-phase second pulse data only, we need to generalise
- """
-
- try:
- canvas.fig.clear()
- except:
- pass
-
-
- self.dataCubeFFT( )
-
-
- canvas.ax1 = canvas.fig.add_axes([.1, .1, .2, .8])
- canvas.ax2 = canvas.fig.add_axes([.325, .1, .2, .8])
- canvas.ax3 = canvas.fig.add_axes([.55, .1, .2, .8])
- canvas.ax4 = canvas.fig.add_axes([.815, .1, .05, .8])
- canvas.ax1.tick_params(axis='both', which='major', labelsize=8)
- canvas.ax2.tick_params(axis='both', which='major', labelsize=8)
- canvas.ax3.tick_params(axis='both', which='major', labelsize=8)
- canvas.ax4.tick_params(axis='both', which='major', labelsize=8)
- canvas.ax1.set_ylabel("pulse index", fontsize=8)
- canvas.ax1.set_xlabel(r"$\omega$ bin", fontsize=8)
- canvas.ax2.set_xlabel(r"$\omega$ bin", fontsize=8)
- canvas.ax3.set_xlabel(r"$\omega$ bin", fontsize=8)
- canvas.ax2.yaxis.set_ticklabels([])
- canvas.ax3.yaxis.set_ticklabels([])
-
-
-
-
-
-
-
-
-
-
-
-
-
-
- for pulse in self.DATADICT["PULSES"]:
- for ichan in self.DATADICT[pulse]["chan"]:
-
- CONTAINER = {}
- CONTAINER["Cycle 1"] = []
- CONTAINER["Cycle 2"] = []
- for istack in self.DATADICT["stacks"]:
-
- ipm = 8
-
-
-
- if not istack%4%4:
-
- CONTAINER["Cycle 1"].append(-self.DATADICT[pulse][ichan]["FFT"][istack])
-
- elif not istack%4%3:
-
- CONTAINER["Cycle 2"].append(-self.DATADICT[pulse][ichan]["FFT"][istack])
-
- elif not istack%4%2:
-
- CONTAINER["Cycle 2"].append( self.DATADICT[pulse][ichan]["FFT"][istack])
-
- else:
-
- CONTAINER["Cycle 1"].append( self.DATADICT[pulse][ichan]["FFT"][istack])
-
-
-
-
-
-
-
-
-
-
-
-
-
-
-
-
- """ All this phase cycling wreaks havoc on a normal calculation of std. and variance. Instead, we resort to calculating
- a pooled variance. In this assumption is that the precision of the measurment is constant. This is a poor choice for
- any type of moving sensor.
- """
-
-
-
-
-
- CONTAINER = .5*(np.array(CONTAINER["Cycle 2"]) - np.array(CONTAINER["Cycle 1"]))
- print ("container shape", np.shape( CONTAINER), self.iWindowStart+1, self.iWindowEnd-1)
- dmin = np.min(np.abs(np.average(np.array(CONTAINER)[:,:,self.iWindowStart+1:self.iWindowEnd-1], axis=0)))
- dmax = np.max(np.abs(np.average(np.array(CONTAINER)[:,:,self.iWindowStart+1:self.iWindowEnd-1], axis=0)))
- mn = canvas.ax1.matshow( 20.*np.log10(np.abs(np.average(np.array(CONTAINER)[:,:, self.iWindowStart+1:self.iWindowEnd-1], axis=0))), aspect='auto', vmin=-120, vmax=-40)
-
- canvas.ax2.matshow( 20*np.log10(np.std(np.real(np.array(CONTAINER)[:,:,self.iWindowStart+1:self.iWindowEnd-1]), axis=0)), aspect='auto', vmin=-120, vmax=-40)
- canvas.ax3.matshow( 20*np.log10(np.std(np.imag(np.array(CONTAINER)[:,:,self.iWindowStart+1:self.iWindowEnd-1]), axis=0)), aspect='auto', vmin=-120, vmax=-40)
-
-
- cb1 = mpl.colorbar.Colorbar(canvas.ax4, mn)
- cb1.ax.tick_params(labelsize=8)
- cb1.set_label("power [dB]", fontsize=8)
- canvas.ax1.xaxis.set_major_locator(MaxNLocator(4))
- canvas.ax2.xaxis.set_major_locator(MaxNLocator(4))
- canvas.ax3.xaxis.set_major_locator(MaxNLocator(4))
- canvas.draw()
-
- self.doneTrigger.emit()
-
- def effectivePulseMoment(self, cv, canvas):
-
- canvas.reAxH(2)
- nstack = len(self.DATADICT["stacks"])
- canvas.ax1.set_yscale('log')
-
- for pulse in self.DATADICT["PULSES"]:
- self.DATADICT[pulse]["qeff"] = {}
- self.DATADICT[pulse]["q_nu"] = {}
- for ipm in range(self.DATADICT["nPulseMoments"]):
- self.DATADICT[pulse]["qeff"][ipm] = {}
- self.DATADICT[pulse]["q_nu"][ipm] = {}
-
-
-
-
-
-
-
-
- scolours = cmocean.cm.ice(np.linspace(0,1,1.5*len(self.DATADICT["stacks"])))
- iistack = 0
- for istack in self.DATADICT["stacks"]:
-
- x = self.DATADICT[pulse]["CURRENT"][ipm][istack]
- X = np.fft.rfft(x)
- v = np.fft.fftfreq(len(x), self.dt)
- v = v[0:len(X)]
- v[-1] = np.abs(v[-1])
-
- I0 = np.abs(X)/len(X)
- qeff = I0 * (self.DATADICT[pulse]["PULSE_TIMES"][-1]-self.DATADICT[pulse]["PULSE_TIMES"][0])
- canvas.ax1.set_title(r"pulse moment index " +str(ipm), fontsize=10)
- canvas.ax1.set_xlabel(r"$\nu$ [Hz]", fontsize=8)
- canvas.ax1.set_ylabel(r"$q_{eff}$ [A$\cdot$sec]", fontsize=8)
- canvas.ax1.plot(v, qeff, color=scolours[iistack] )
- self.DATADICT[pulse]["qeff"][ipm][istack] = qeff
- self.DATADICT[pulse]["q_nu"][ipm][istack] = v
- iistack += 1
- canvas.draw()
-
- percent = int(1e2* (float)((istack)+ipm*self.DATADICT["nPulseMoments"]) /
- (float)(len(self.DATADICT["PULSES"])*self.DATADICT["nPulseMoments"]*nstack))
- self.progressTrigger.emit(percent)
-
- canvas.draw()
-
-
- self.plotQeffNu(cv, canvas.ax2)
- canvas.draw()
- self.doneTrigger.emit()
-
- def plotQeffNu(self, cv, ax):
-
-
-
- nstack = len(self.DATADICT["stacks"])
- iFID = 0
- for pulse in self.DATADICT["PULSES"]:
- self.DATADICT[pulse]["Q"] = np.zeros( (self.DATADICT["nPulseMoments"], len(self.DATADICT["stacks"])) )
- ilabel = True
- for ipm in range(self.DATADICT["nPulseMoments"]):
-
- scolours = cmocean.cm.ice(np.linspace(0,1,1.5*len(self.DATADICT["stacks"])))
-
-
- istack = 0
- for stack in self.DATADICT["stacks"]:
-
- icv = int (round(cv / self.DATADICT[pulse]["q_nu"][ipm][stack][1]))
- self.DATADICT[pulse]["Q"][ipm,istack] = self.DATADICT[pulse]["qeff"][ipm][stack][icv]
- if ilabel:
- ax.scatter(ipm, self.DATADICT[pulse]["qeff"][ipm][stack][icv], facecolors='none', edgecolors=scolours[istack], label=(str(pulse)))
- ilabel = False
- else:
- ax.scatter(ipm, self.DATADICT[pulse]["qeff"][ipm][stack][icv], facecolors='none', edgecolors=scolours[istack])
-
-
- percent = int(1e2* (float)((istack)+ipm*self.DATADICT["nPulseMoments"]) /
- (float)(len(self.DATADICT["PULSES"])*self.DATADICT["nPulseMoments"]*nstack))
- self.progressTrigger.emit(percent)
- istack += 1
- iFID += 1
- ax.set_xlabel(r"pulse moment index", fontsize=8)
- ax.set_ylabel(r"$q_{eff}$ [A$\cdot$sec]", fontsize=8)
- ax.set_yscale('log')
- ax.set_xlim(0, ax.get_xlim()[1])
- ax.legend(loc='upper right', scatterpoints = 1, prop={'size':6})
-
- def enableDSP(self):
- self.enableDSPTrigger.emit()
-
- def adaptiveFilter(self, M, flambda, truncate, mu, PCA, canvas):
-
- canvas.reAx2(shx=False, shy=False)
-
- canvas.ax2.tick_params(axis='both', which='major', labelsize=8)
- canvas.ax2.ticklabel_format(style='sci', scilimits=(0,0), axis='y')
-
- canvas.ax2.tick_params(axis='both', which='major', labelsize=8)
- canvas.ax2.ticklabel_format(style='sci', scilimits=(0,0), axis='y')
-
- if truncate:
- itrunc =(int) ( round( 1e-3*truncate*self.samp ) )
-
- print( "adaptive filter size", 1e3*self.dt*M, " [ms]" )
-
- Filt = adapt.AdaptiveFilter(flambda)
- H = {}
- for pulse in self.DATADICT["PULSES"]:
- H[pulse] = {}
- for ichan in self.DATADICT[pulse]["chan"]:
- print("setting H", pulse, ichan)
- H[pulse][ichan] = np.zeros(M)
-
- iFID = 0
-
-
-
-
-
- for istack in self.DATADICT["stacks"]:
- for ipm in range(self.DATADICT["nPulseMoments"]):
- for pulse in self.DATADICT["PULSES"]:
- canvas.ax1.clear()
- canvas.ax2.clear()
- for ichan in self.DATADICT[pulse]["chan"]:
-
- canvas.ax1.plot( self.DATADICT[pulse]["TIMES"], 1e9* self.DATADICT[pulse][ichan][ipm][istack],\
- label = "noisy")
- RX = []
- for irchan in self.DATADICT[pulse]["rchan"]:
- RX.append(self.DATADICT[pulse][irchan][ipm][istack][::-1])
- if all(H[pulse][ichan]) == 0:
-
- [e,H[pulse][ichan]] = Filt.adapt_filt_Ref( self.DATADICT[pulse][ichan][ipm][istack][::-1],\
- RX,\
- M, mu, PCA, flambda, H[pulse][ichan])
- [e,H[pulse][ichan]] = Filt.adapt_filt_Ref( self.DATADICT[pulse][ichan][ipm][istack][::-1],\
- RX,\
- M, mu, PCA, flambda, H[pulse][ichan])
- else:
- [e,H[pulse][ichan]] = Filt.adapt_filt_Ref( self.DATADICT[pulse][ichan][ipm][istack][::-1],\
- RX,\
- M, mu, PCA, flambda, H[pulse][ichan])
-
- if truncate:
- canvas.ax1.plot( self.DATADICT[pulse]["TIMES"][0:itrunc], 1e9* e[::-1][0:itrunc],\
- label = pulse + " ipm=" + str(ipm) + " istack=" + str(istack) + " ichan=" + str(ichan))
- self.DATADICT[pulse][ichan][ipm][istack] = e[::-1][0:itrunc]
- else:
- canvas.ax1.plot( self.DATADICT[pulse]["TIMES"], 1e9* e[::-1],\
- label = pulse + " ipm=" + str(ipm) + " istack=" + str(istack) + " ichan=" + str(ichan))
- self.DATADICT[pulse][ichan][ipm][istack] = e[::-1]
-
- canvas.ax2.plot( H[pulse][ichan] , label="taps")
- canvas.ax1.legend(prop={'size':6})
- canvas.ax2.legend(prop={'size':6})
-
- canvas.ax1.set_xlabel(r"time [s]", fontsize=8)
- canvas.ax1.set_ylabel(r"signal [nV]", fontsize=8)
-
- canvas.ax2.set_xlabel(r"filter index", fontsize=8)
- canvas.ax2.set_ylabel(r"scale factor", fontsize=8)
-
- canvas.draw()
-
-
-
- for rchan in self.DATADICT[pulse]["rchan"]:
- if truncate:
- self.DATADICT[pulse][rchan][ipm][istack] = self.DATADICT[pulse][rchan][ipm][istack][0:itrunc]
-
-
- percent = (int)(1e2*((float)(istack*self.DATADICT["nPulseMoments"]+(ipm))/( len(self.DATADICT["PULSES"])*self.nPulseMoments*(len(self.DATADICT["stacks"])+1) )))
- self.progressTrigger.emit(percent)
-
-
-
-
-
-
-
-
-
- if truncate:
- self.DATADICT[pulse]["TIMES"] = self.DATADICT[pulse]["TIMES"][0:itrunc]
-
- self.doneTrigger.emit()
- self.updateProcTrigger.emit()
-
-
-
- def plotFT(self, canvas, istart=0, iend=0):
-
- try:
- canvas.fig.clear()
- except:
- pass
-
- canvas.ax1 = canvas.fig.add_axes([.1, .1, .65, .8])
- canvas.ax1c = canvas.fig.add_axes([.8, .1, .05, .8])
- canvas.ax1.tick_params(axis='both', which='major', labelsize=8)
-
- for pulse in self.DATADICT["PULSES"]:
- for istack in self.DATADICT["stacks"]:
- for ichan in self.DATADICT[pulse]["chan"]:
-
- XA = np.zeros((self.DATADICT["nPulseMoments"] , len(self.DATADICT[pulse][ichan][0][istack])/2+1))
- nu = np.fft.fftfreq(self.DATADICT[pulse][ichan][0][istack].size, d=self.dt)
- nu[-1] *= -1
- df = nu[1]
- of = 0
- if istart:
- of = nu[istart]
- def freqlabel(x, pos):
- return '%1.0f' %(of + x*df)
- formatter = FuncFormatter(freqlabel)
- canvas.ax1.clear()
- for ipm in range(self.DATADICT["nPulseMoments"]):
- X = np.fft.rfft(self.DATADICT[pulse][ichan][ipm][istack])
- XA[ipm,:] = np.abs(X)
- if istart:
- mn = canvas.ax1.matshow(20.*np.log10(XA[:,istart:iend+1]), aspect='auto', vmax=-40, vmin=-120)
- else:
- mn = canvas.ax1.matshow(20.*np.log10(XA), aspect='auto', vmax=-40, vmin=-120)
- smin = np.min(20.*np.log10(XA))
- smax = np.max(20.*np.log10(XA))
- canvas.ax1.xaxis.set_major_formatter(formatter)
- cb1 = mpl.colorbar.Colorbar(canvas.ax1c, mn)
- cb1.ax.tick_params(labelsize=8)
- cb1.set_label("signal [dB]", fontsize=8)
- canvas.ax1.set_xlabel(r"$\nu$ [Hz]", fontsize=10)
- canvas.ax1.set_ylabel(r"$q_{index}$", fontsize=10)
- canvas.draw()
-
- def plotFT(self, canvas, istart=0, iend=0):
-
- try:
- canvas.fig.clear()
- except:
- pass
-
- canvas.ax1 = canvas.fig.add_axes([.1, .1, .65, .8])
- canvas.ax1c = canvas.fig.add_axes([.8, .1, .05, .8])
- canvas.ax1.tick_params(axis='both', which='major', labelsize=8)
-
- for pulse in self.DATADICT["PULSES"]:
- for istack in self.DATADICT["stacks"]:
- for ichan in self.DATADICT[pulse]["chan"]:
-
- XA = np.zeros((self.DATADICT["nPulseMoments"] , len(self.DATADICT[pulse][ichan][0][istack])//2+1))
- nu = np.fft.fftfreq(self.DATADICT[pulse][ichan][0][istack].size, d=self.dt)
- nu[-1] *= -1
- df = nu[1]
- of = 0
- if istart:
- of = nu[istart]
- def freqlabel(x, pos):
- return '%1.0f' %(of + x*df)
- formatter = FuncFormatter(freqlabel)
- canvas.ax1.clear()
- for ipm in range(self.DATADICT["nPulseMoments"]):
- X = np.fft.rfft(self.DATADICT[pulse][ichan][ipm][istack])
- XA[ipm,:] = np.abs(X)
- if istart:
- mn = canvas.ax1.matshow(20.*np.log10(XA[:,istart:iend+1]), aspect='auto', vmax=-40, vmin=-120, cmap='viridis')
- else:
- mn = canvas.ax1.matshow(20.*np.log10(XA), aspect='auto', vmax=-40, vmin=-120, cmap='viridis')
- canvas.ax1.xaxis.set_major_formatter(formatter)
- cb1 = mpl.colorbar.Colorbar(canvas.ax1c, mn)
- cb1.ax.tick_params(labelsize=8)
- cb1.set_label("signal [dB]", fontsize=8)
- canvas.ax1.set_xlabel(r"$\nu$ [Hz]", fontsize=10)
- canvas.ax1.set_ylabel(r"$q_{index}$", fontsize=10)
- canvas.draw()
-
-
- def dataCubeFFT(self):
- """
- Performs FFT on entire cube of DATA, and REFERENCE channels, but not pulse currents,
- Results are saved to a new field in the data structure
-
- The GMR varies phase as a function of pulse moment index, so that the first pusle moment is zero phase,
- the second is pi/2 the third is zero. This method corrects for this, so that all pulse moments are in phase.
-
- Technically we may not want to do this, if there is some system response that this cycles away, and we lose track of
- how many of each cycle we have, could this be problomatic? I think it will come out in the wash as we keep track of the
- rest of the phase cycles. Holy phase cycling batman.
- """
- for pulse in self.DATADICT["PULSES"]:
- for ichan in np.append(self.DATADICT[pulse]["chan"], self.DATADICT[pulse]["rchan"]):
-
- self.DATADICT[pulse][ichan]["FFT"] = {}
- self.DATADICT[pulse][ichan]["FFT"]["nu"] = np.fft.fftfreq(self.DATADICT[pulse][ichan][0][self.DATADICT["stacks"][0]].size, d=self.dt)
- self.DATADICT[pulse][ichan]["FFT"]["nu"][-1] *= -1
- for istack in self.DATADICT["stacks"]:
- self.DATADICT[pulse][ichan]["FFT"][istack] = np.zeros((self.DATADICT["nPulseMoments"] , len(self.DATADICT[pulse][ichan][0][istack])//2+1), dtype=complex)
- for ipm in range(self.DATADICT["nPulseMoments"]):
-
- mod = (-1)**(ipm%2) * (-1)**(istack%2)
- self.DATADICT[pulse][ichan]["FFT"][istack][ipm,:] = np.fft.rfft( self.DATADICT[pulse][ichan][ipm][istack] )
-
-
-
-
-
-
-
- def adaptiveFilterFD(self, ftype, band, centre, canvas):
-
- try:
- canvas.fig.clear()
- except:
- pass
-
- canvas.ax1 = canvas.fig.add_axes([.1, .5, .7, .4])
- canvas.ax1c = canvas.fig.add_axes([.85, .5, .05, .4])
- canvas.ax1.tick_params(axis='both', which='major', labelsize=8)
-
-
- canvas.ax2 = canvas.fig.add_axes([.1, .05, .7, .4])
- canvas.ax2c = canvas.fig.add_axes([.85, .05, .05, .4])
- canvas.ax2.tick_params(axis='both', which='major', labelsize=8)
-
-
- self.dataCubeFFT()
-
- Filt = adapt.AdaptiveFilter(0.)
- for pulse in self.DATADICT["PULSES"]:
-
- [WINDOW, nd, wstart, wend, dead] = self.computeWindow(pulse, band, centre, ftype)
- for istack in self.DATADICT["stacks"]:
- for ichan in self.DATADICT[pulse]["chan"]:
-
- nd = len(self.DATADICT[pulse][ichan][0][istack])
- XX = np.zeros((self.DATADICT["nPulseMoments"] , len(self.DATADICT[pulse][ichan][0][istack])//2+1), dtype=complex)
- nu = np.fft.fftfreq(self.DATADICT[pulse][ichan][0][istack].size, d=self.dt)
- nu[-1] *= -1
-
- def freqlabel(x, pos):
- return '%1.0f' %((wstart)*nu[1] + x*nu[1])
- formatter = FuncFormatter(freqlabel)
- canvas.ax1.clear()
- for ipm in range(self.DATADICT["nPulseMoments"]):
- X = np.fft.rfft(self.DATADICT[pulse][ichan][ipm][istack])
- XX[ipm,:] = X
-
- XX = XX*WINDOW
- XX = XX[:,wstart:wend]
- smin = np.min(20.*np.log10(np.abs(XX)))
- smax = np.max(20.*np.log10(np.abs(XX)))
-
- smax = -40
- smin = -120
- mn = canvas.ax1.matshow(20.*np.log10(np.abs(XX)), aspect='auto', vmin=smin, vmax=smax)
-
- canvas.ax1.xaxis.set_major_formatter(formatter)
- cb1 = mpl.colorbar.Colorbar(canvas.ax1c, mn)
-
- RX = []
- for ichan in self.DATADICT[pulse]["rchan"]:
- R = np.zeros((self.DATADICT["nPulseMoments"] , len(self.DATADICT[pulse][ichan][0][istack])//2+1), dtype=complex)
- for ipm in range(self.DATADICT["nPulseMoments"]):
- R[ipm,:] = np.fft.rfft(self.DATADICT[pulse][ichan][ipm][istack])
- RX.append(R[:,wstart:wend])
- XC = Filt.transferFunctionFFT(XX, RX)
-
-
-
-
-
- mc = canvas.ax2.matshow(20.*np.log10(np.abs(XC)), aspect='auto', vmin=smin, vmax=smax)
- cb2 = mpl.colorbar.Colorbar(canvas.ax2c, mc)
- cmin = np.min(20.*np.log10(np.abs(XC)))
- cmax = np.max(20.*np.log10(np.abs(XC)))
- canvas.ax2.xaxis.set_major_formatter(formatter)
-
- canvas.draw()
-
-
-
-
-
- self.doneTrigger.emit()
-
-
- def findSpikes(self, x, width, threshold, rollOn):
- import scipy.ndimage as im
- spikes = np.zeros( len(x) )
- med = im.median_filter(x, width,mode='nearest')
- std = np.std(x)
- spikes = (np.abs(x-med) > threshold * std)
- return np.array(np.where(spikes[rollOn::])) + rollOn
-
-
-
-
-
-
-
-
-
-
-
-
-
-
-
-
-
-
-
-
-
-
-
-
-
-
-
-
-
-
-
-
-
-
-
-
-
-
-
-
-
- def designFilter(self, cf, PB, SB, gpass, gstop, ftype, canvas):
- ''' cf is central frequency
- pb is pass band
- sb is stop band
- '''
- TS = (cf) / (.5/self.dt)
- PB = PB / (.5/self.dt)
- SB = SB / (.5/self.dt)
-
-
- if ftype=="Butterworth":
- [bord, wn] = signal.buttord([TS-PB,TS+PB], [TS-SB,TS+SB], gpass, gstop)
- [self.filt_b, self.filt_a] = signal.butter(bord, wn, btype='bandpass', output='ba')
- [self.filt_z, self.filt_p, self.filt_k] = signal.butter(bord, wn, btype='band', output='zpk')
- elif ftype == "Chebychev Type II":
- [bord, wn] = signal.cheb2ord([TS-PB,TS+PB], [TS-SB,TS+SB], gpass, gstop)
- [self.filt_b, self.filt_a] = signal.cheby2(bord, gstop, wn, btype='bandpass', output='ba')
- [self.filt_z, self.filt_p, self.filt_k] = signal.cheby2(bord, gstop, wn, btype='band', output='zpk')
- elif ftype == "Elliptic":
- [bord, wn] = signal.ellipord([TS-PB,TS+PB], [TS-SB,TS+SB], gpass, gstop)
- [self.filt_b, self.filt_a] = signal.ellip(bord, gpass, gstop, wn, btype='bandpass', output='ba')
- [self.filt_z, self.filt_p, self.filt_k] = signal.ellip(bord, gpass, gstop, wn, btype='band', output='zpk')
-
-
- impulse = self.mfreqz2(self.filt_b, self.filt_a, canvas)
- self.fe = -5
- for it in range(len(impulse[0])):
- if abs(impulse[1][0][it][0]) >= .1 * gpass:
- self.fe = impulse[0][it]
-
- return [bord, self.fe]
-
- def downsample(self, truncate, dec, canvas):
- """ Downsamples and truncates the raw signal.
- """
- canvas.reAx2()
-
- canvas.ax1.tick_params(axis='both', which='major', labelsize=8)
- canvas.ax1.ticklabel_format(style='sci', scilimits=(0,0), axis='y')
- canvas.ax2.tick_params(axis='both', which='major', labelsize=8)
- canvas.ax2.ticklabel_format(style='sci', scilimits=(0,0), axis='y')
-
- self.samp /= dec
- self.dt = 1./self.samp
-
- if truncate:
- itrunc = (int)( 1e-3*truncate*self.samp )
- iFID = 0
- for pulse in self.DATADICT["PULSES"]:
- for ipm in range(self.DATADICT["nPulseMoments"]):
- for istack in self.DATADICT["stacks"]:
- canvas.ax1.clear()
- canvas.ax2.clear()
- for ichan in np.append(self.DATADICT[pulse]["chan"], self.DATADICT[pulse]["rchan"]):
-
- ndi = np.shape(self.DATADICT[pulse][ichan][ipm][istack])[0]%dec
- if ndi:
- [self.DATADICT[pulse][ichan][ipm][istack], RSTIMES] = signal.resample(self.DATADICT[pulse][ichan][ipm][istack][0:-ndi],\
- len(self.DATADICT[pulse][ichan][ipm][istack][0:-ndi])//dec,\
- self.DATADICT[pulse]["TIMES"][0:-ndi], window='hamm')
- else:
- [self.DATADICT[pulse][ichan][ipm][istack], RSTIMES] = signal.resample(self.DATADICT[pulse][ichan][ipm][istack],\
- len(self.DATADICT[pulse][ichan][ipm][istack])//dec,\
- self.DATADICT[pulse]["TIMES"], window='hamm')
- if truncate:
- self.DATADICT[pulse][ichan][ipm][istack] = self.DATADICT[pulse][ichan][ipm][istack][0:itrunc]
- RSTIMES = RSTIMES[0:itrunc]
-
-
- for ichan in self.DATADICT[pulse]["chan"]:
- canvas.ax2.plot( RSTIMES, 1e9*self.DATADICT[pulse][ichan][ipm][istack], \
- label = pulse + " ipm=" + str(ipm) + " istack=" + str(istack) + " ichan=" + str(ichan))
- for ichan in self.DATADICT[pulse]["rchan"]:
- canvas.ax1.plot( RSTIMES, 1e9*self.DATADICT[pulse][ichan][ipm][istack], \
- label = pulse + " ipm=" + str(ipm) + " istack=" + str(istack) + " ichan=" + str(ichan))
-
- canvas.ax1.set_xlabel(r"time [s]", fontsize=8)
- canvas.ax1.set_ylabel(r"signal [nV]", fontsize=8)
- canvas.ax2.set_xlabel(r"time [s]", fontsize=8)
- canvas.ax2.set_ylabel(r"signal [nV]", fontsize=8)
-
- canvas.ax1.legend(prop={'size':6})
- canvas.ax2.legend(prop={'size':6})
- canvas.ax1.ticklabel_format(style='sci', scilimits=(0,0), axis='y')
- canvas.ax2.ticklabel_format(style='sci', scilimits=(0,0), axis='y')
- canvas.draw()
- percent = (int)(1e2*((float)(iFID*self.DATADICT["nPulseMoments"]+(ipm))/( len(self.DATADICT["PULSES"])*self.nPulseMoments)))
- self.progressTrigger.emit(percent)
- iFID += 1
- self.DATADICT[pulse]["TIMES"] = RSTIMES
-
-
-
- for pulse in self.DATADICT["PULSES"]:
- for ipm in range(self.DATADICT["nPulseMoments"]):
- for istack in self.DATADICT["stacks"]:
- ndi = np.shape(self.DATADICT[pulse]["CURRENT"][ipm][istack])[0]%dec
- if ndi:
- [self.DATADICT[pulse]["CURRENT"][ipm][istack], RSPTIMES] = signal.resample(self.DATADICT[pulse]["CURRENT"][ipm][istack][0:-ndi],\
- len(self.DATADICT[pulse]["CURRENT"][ipm][istack][0:-ndi])//dec,\
- self.DATADICT[pulse]["PULSE_TIMES"][0:-ndi], window='hamm')
- else:
- [self.DATADICT[pulse]["CURRENT"][ipm][istack], RSPTIMES] = signal.resample(self.DATADICT[pulse]["CURRENT"][ipm][istack],\
- len(self.DATADICT[pulse]["CURRENT"][ipm][istack])//dec,\
- self.DATADICT[pulse]["PULSE_TIMES"], window='hamm')
- self.DATADICT[pulse]["PULSE_TIMES"] = RSPTIMES
- self.doneTrigger.emit()
- self.updateProcTrigger.emit()
-
- def computeWindow(self, pulse, band, centre, ftype, canvas=None):
-
-
- nd = len(self.DATADICT[pulse][self.DATADICT[pulse]["chan"][0]][0][self.DATADICT["stacks"][0]])
- fft1 = np.fft.rfft(self.DATADICT[pulse][self.DATADICT[pulse]["chan"][0]][0][self.DATADICT["stacks"][0]])
- freqs = np.fft.fftfreq(nd, self.dt)
- df = freqs[1] - freqs[0]
- N = int((round)(band/df))
-
- if ftype == "Hamming":
- window = np.hamming(N)
- elif ftype == "Hanning":
- window = np.hanning(N)
- elif ftype == "Rectangular":
- window = np.ones(N)
- elif ftype == "Flat top":
- window = signal.flattop(N)
- else:
- print ("in windowFilter, window type undefined")
-
- WINDOW = np.zeros(len(fft1))
- ifreq = int(round(centre/df))
- istart = ifreq-len(window)//2
- iend = 0
- if N%2:
- WINDOW[ifreq-N//2:ifreq+N//2+1] = window
- iend = ifreq+N//2+1
- else:
- WINDOW[ifreq-N//2:ifreq+N//2] = window
- iend = ifreq+N//2
-
- self.WINDOW = WINDOW
- self.iWindowStart = istart
- self.iWindowEnd = iend
- self.FFTtimes = nd
-
- fft1 = np.fft.irfft(WINDOW)
-
- self.windead = 0.
- for ift in np.arange(100,0,-1):
-
- if fft1[ift] > .001:
-
- dead = 1e3*self.DATADICT[pulse]["TIMES"][ift] - 1e3*self.DATADICT[pulse]["TIMES"][0]
- self.windead = self.DATADICT[pulse]["TIMES"][ift] - self.DATADICT[pulse]["TIMES"][0]
- break
-
- if canvas != None:
- canvas.fig.clear()
- canvas.ax1 = canvas.fig.add_axes([.1, .6, .75, .35])
- canvas.ax2 = canvas.fig.add_axes([.1, .1, .75, .35])
- canvas.ax1.plot(WINDOW)
- canvas.ax2.plot( 1e3* self.DATADICT[pulse]["TIMES"][0:100] - 1e3*self.DATADICT[pulse]["TIMES"][0], fft1[0:100] )
- canvas.ax2.set_xlabel("time (ms)")
- canvas.ax2.set_title("IFFT")
- canvas.draw()
-
- return [WINDOW, nd, istart, iend, dead]
-
- def windowFilter(self, ftype, band, centre, canvas):
-
-
-
-
- iFID = 0
- for pulse in self.DATADICT["PULSES"]:
- [WINDOW, nd, istart, iend,dead] = self.computeWindow(pulse, band, centre, ftype)
- for istack in self.DATADICT["stacks"]:
- for ipm in range(self.DATADICT["nPulseMoments"]):
- for ichan in np.append(self.DATADICT[pulse]["chan"], self.DATADICT[pulse]["rchan"]):
- fft = np.fft.rfft( self.DATADICT[pulse][ichan][ipm][istack] )
- fft *= WINDOW
- self.DATADICT[pulse][ichan][ipm][istack] = np.fft.irfft(fft , nd)
-
- percent = (int)(1e2*((float)(iFID*self.DATADICT["nPulseMoments"]+(ipm))/(len(self.DATADICT["PULSES"])*self.nPulseMoments)))
- self.progressTrigger.emit(percent)
- iFID += 1
-
- self.plotFT(canvas, istart, iend)
- self.doneTrigger.emit()
-
- def bandpassFilter(self, canvas, blank, plot=True):
-
- if plot:
- canvas.reAx2()
- canvas.ax1.tick_params(axis='both', which='major', labelsize=8)
- canvas.ax1.ticklabel_format(style='sci', scilimits=(0,0), axis='y')
- canvas.ax2.tick_params(axis='both', which='major', labelsize=8)
- canvas.ax2.ticklabel_format(style='sci', scilimits=(0,0), axis='y')
-
- ife = (int)( max(self.fe, self.windead) * self.samp )
-
- iFID = 0
- for pulse in self.DATADICT["PULSES"]:
- self.DATADICT[pulse]["TIMES"] = self.DATADICT[pulse]["TIMES"][ife:-ife]
- for ipm in range(self.DATADICT["nPulseMoments"]):
- for istack in self.DATADICT["stacks"]:
- canvas.ax1.clear()
- canvas.ax2.clear()
-
- for ichan in self.DATADICT[pulse]["rchan"]:
-
-
-
-
-
-
-
- self.DATADICT[pulse][ichan][ipm][istack] = \
- signal.filtfilt(self.filt_b, self.filt_a, self.DATADICT[pulse][ichan][ipm][istack])[ife:-ife]
-
-
- if plot:
- canvas.ax1.plot( self.DATADICT[pulse]["TIMES"], 1e9*self.DATADICT[pulse][ichan][ipm][istack], \
- label = pulse + " ipm=" + str(ipm) + " istack=" + str(istack) + " rchan=" + str(ichan))
-
- for ichan in self.DATADICT[pulse]["chan"]:
-
-
-
-
-
-
-
- self.DATADICT[pulse][ichan][ipm][istack] = \
- scipy.signal.filtfilt(self.filt_b, self.filt_a, self.DATADICT[pulse][ichan][ipm][istack])[ife:-ife]
-
-
- if plot:
- canvas.ax2.plot( self.DATADICT[pulse]["TIMES"], 1e9*self.DATADICT[pulse][ichan][ipm][istack], \
- label = "data " + pulse + " ipm=" + str(ipm) + " istack=" + str(istack) + " chan=" + str(ichan))
-
- if plot:
- canvas.ax1.set_xlabel(r"time [s]", fontsize=8)
- canvas.ax1.set_ylabel(r"signal [nV]", fontsize=8)
- canvas.ax2.set_xlabel(r"time [s]", fontsize=8)
- canvas.ax2.set_ylabel(r"signal [nV]", fontsize=8)
- canvas.ax1.legend(prop={'size':6})
- canvas.ax2.legend(prop={'size':6})
- canvas.draw()
-
- percent = (int)(1e2*((float)(iFID*self.DATADICT["nPulseMoments"]+(ipm))/(len(self.DATADICT["PULSES"])*self.nPulseMoments)))
- self.progressTrigger.emit(percent)
- iFID += 1
- self.doneTrigger.emit()
- self.updateProcTrigger.emit()
-
- def loadGMRBinaryFID( self, rawfname, istack ):
- """ Reads a single binary GMR file and fills into DATADICT
- """
-
-
-
-
- nps = (int)((self.prePulseDelay)*self.samp)
- npul = (int)(self.pulseLength[0]*self.samp)
-
-
- nds = nps+npul+(int)((self.deadTime)*self.samp);
- nd1 = (int)(1.*self.samp)
-
- invGain = 1./self.RxGain
- invCGain = self.CurrentGain
-
- pulse = "Pulse 1"
- chan = self.DATADICT[pulse]["chan"]
- rchan = self.DATADICT[pulse]["rchan"]
-
- rawFile = open( rawfname, 'rb')
-
- for ipm in range(self.nPulseMoments):
- buf1 = rawFile.read(4)
- buf2 = rawFile.read(4)
-
- N_chan = struct.unpack('>i', buf1 )[0]
- N_samp = struct.unpack('>i', buf2 )[0]
-
- T = N_samp * self.dt
- TIMES = np.arange(0, T, self.dt) - .0002
-
- DATA = np.zeros([N_samp, N_chan+1])
- for ichan in range(N_chan):
- DATADUMP = rawFile.read(4*N_samp)
- for irec in range(N_samp):
- DATA[irec,ichan] = struct.unpack('>f', DATADUMP[irec*4:irec*4+4])[0]
-
-
- for ichan in chan:
- self.DATADICT["Pulse 1"][ichan][ipm][istack] = DATA[:,eval(ichan)+3][nds:nds+nd1] * invGain
- self.DATADICT["Pulse 1"]["TIMES"] = TIMES[nds:nds+nd1]
- self.DATADICT["Pulse 1"]["CURRENT"][ipm][istack] = DATA[:,1][nps:nps+npul] * invCGain
- self.DATADICT["Pulse 1"]["PULSE_TIMES"] = TIMES[nps:nps+npul]
-
-
- for ichan in rchan:
- self.DATADICT["Pulse 1"][ichan][ipm][istack] = DATA[:,eval(ichan)+3][nds:nds+nd1] * invGain
- self.DATADICT["Pulse 1"]["TIMES"] = TIMES[nds:nds+nd1]
-
- def loadGMRASCIIFID( self, rawfname, istack ):
- """Based on the geoMRI instrument manufactured by VistaClara. Imports
- a suite of raw .lvm files with the following format (on one line)
-
- time(s) DC_Bus/100(V) Current+/75(A) Curr-/75(A) Voltage+/200(V) \
- Ch1(V) Ch2(V) Ch3(V) Ch4(V)
-
- Sampling rate is assumed at 50 kHz
- """
- import pandas as pd
-
-
-
- nps = (int)((self.prePulseDelay)*self.samp)
- npul = (int)(self.pulseLength[0]*self.samp)
-
-
- nds = nps+npul+(int)((self.deadTime)*self.samp);
- nd1 = (int)(1.*self.samp) - nds
- ndr = (int)(1.*self.samp)
-
- invGain = 1./self.RxGain
- invCGain = self.CurrentGain
-
- pulse = "Pulse 1"
- chan = self.DATADICT[pulse]["chan"]
- rchan = self.DATADICT[pulse]["rchan"]
-
- T = 1.5
- TIMES = np.arange(0, T, self.dt) - .0002
-
- self.DATADICT["Pulse 1"]["TIMES"] = TIMES[nds:nds+nd1]
- self.DATADICT["Pulse 1"]["PULSE_TIMES"] = TIMES[nps:nps+npul]
-
-
-
- DATA = pd.read_csv(rawfname, header=None, sep="\t").values
- for ipm in range(self.nPulseMoments):
- for ichan in np.append(chan,rchan):
- self.DATADICT["Pulse 1"][ichan][ipm][istack] = DATA[:, eval(ichan)+4][nds:(nds+nd1)] * invGain
- self.DATADICT["Pulse 1"]["CURRENT"][ipm][istack] = DATA[:,2][nps:nps+npul] * invCGain
- nds += ndr
- nps += ndr
-
-
- def loadFIDData(self, base, procStacks, chanin, rchanin, FIDProc, canvas, deadTime, plot):
- '''
- Loads a GMR FID dataset, reads binary and ASCII format files
- '''
-
- canvas.reAx3(True,False)
-
- chan = []
- for ch in chanin:
- chan.append(str(ch))
-
- rchan = []
- for ch in rchanin:
- rchan.append(str(ch))
-
-
-
- self.deadTime = deadTime
- self.samp = 50000.
- self.dt = 1./self.samp
-
-
-
-
- PULSES = [FIDProc]
- PULSES = ["Pulse 1"]
-
- self.DATADICT = {}
- self.DATADICT["nPulseMoments"] = self.nPulseMoments
- self.DATADICT["stacks"] = procStacks
- self.DATADICT["PULSES"] = PULSES
- for pulse in PULSES:
- self.DATADICT[pulse] = {}
- self.DATADICT[pulse]["chan"] = chan
- self.DATADICT[pulse]["rchan"] = rchan
- self.DATADICT[pulse]["CURRENT"] = {}
- for ichan in np.append(chan,rchan):
- self.DATADICT[pulse][ichan] = {}
- for ipm in range(self.nPulseMoments):
- self.DATADICT[pulse][ichan][ipm] = {}
- self.DATADICT[pulse]["CURRENT"][ipm] = {}
- for istack in procStacks:
- self.DATADICT[pulse][ichan][ipm][istack] = np.zeros(3)
- self.DATADICT[pulse]["CURRENT"][ipm][istack] = np.zeros(3)
-
-
-
- iistack = 0
- fnames = []
- for istack in procStacks:
- if self.nDAQVersion < 2.3:
-
- self.loadGMRASCIIFID( base + "_" + str(istack), istack )
- else:
- self.loadGMRBinaryFID( base + "_" + str(istack) + ".lvm", istack )
-
-
- percent = (int) (1e2*((float)((iistack*self.nPulseMoments+ipm+1)) / (len(procStacks)*self.nPulseMoments)))
- self.progressTrigger.emit(percent)
- iistack += 1
-
-
-
-
-
-
-
-
- if plot:
- iistack = 0
- for istack in procStacks:
- for ipm in range(self.nPulseMoments):
- canvas.ax1.clear()
- canvas.ax2.clear()
- canvas.ax3.clear()
-
- for ichan in chan:
- canvas.ax1.plot(self.DATADICT["Pulse 1"]["PULSE_TIMES"], self.DATADICT["Pulse 1"]["CURRENT"][ipm][istack] , color='black')
- canvas.ax3.plot(self.DATADICT["Pulse 1"]["TIMES"], self.DATADICT["Pulse 1"][ichan][ipm][istack], label="Pulse 1 FID data ch. "+str(ichan))
-
- for ichan in rchan:
- canvas.ax2.plot(self.DATADICT["Pulse 1"]["TIMES"], self.DATADICT["Pulse 1"][ichan][ipm][istack], label="Pulse 1 FID ref ch. "+str(ichan))
-
- canvas.ax3.legend(prop={'size':6})
- canvas.ax2.legend(prop={'size':6})
-
- canvas.ax1.set_title("stack "+str(istack)+" pulse index " + str(ipm), fontsize=8)
- canvas.ax1.set_xlabel("time [s]", fontsize=8)
- canvas.ax1.set_ylabel("Current [A]", fontsize=8)
- canvas.ax1.ticklabel_format(style='sci', scilimits=(0,0), axis='y')
-
- canvas.ax2.set_ylabel("RAW signal [V]", fontsize=8)
- canvas.ax2.tick_params(axis='both', which='major', labelsize=8)
- canvas.ax2.tick_params(axis='both', which='minor', labelsize=6)
- canvas.ax2.set_xlabel("time [s]", fontsize=8)
- canvas.ax2.ticklabel_format(style='sci', scilimits=(0,0), axis='y')
- canvas.ax3.ticklabel_format(style='sci', scilimits=(0,0), axis='y')
- canvas.draw()
-
- percent = (int) (1e2*((float)((iistack*self.nPulseMoments+ipm+1)) / (len(procStacks)*self.nPulseMoments)))
- self.progressTrigger.emit(percent)
- iistack += 1
-
- self.enableDSP()
- self.doneTrigger.emit()
-
- def load4PhaseT1Data(self, base, procStacks, chan, rchan, FIDProc, canvas, deadTime, plot):
-
- """
- Designed to load GMR 4-phase data which use the following convention for phase cycles
- P1 P2
- Stack 1 -> 0 0 <-- <--
- Stack 2 -> 0 pi/2 | <-- <--
- Stack 3 -> pi/2 0 <-- | <--
- Stack 4 -> pi/2 pi/2 <-- <--
- The cycle is determined by stack indice. Walbrecker proposes for pulse2 data (Stack2 - Stack1) / 2
- equivalently (Stack 4 - Stack3) will yield the same voltage response wrt. the second pulse.
- Alternatively Stack 4 can be converted to be aligned with Stack 1 by negating, and Stack 3 Can be aligned with Stack 2 by negating
- Then there are just the two phase cycles that can be stacked like normal.
- Unfortunately, we need to stack each cycle first, then perform corrections for phase cycling. The reason for this is that otherwise,
- the entire point is lost, as the signal that is desired to be cancelled out may not be balanced evenly across the stacks. That is to say,
- if there is an uneven number of a certain phase cycle.
-
- We could, I suppose impose this condition, but I think I would rather not?
- + more samples for std. deviation calculation
- + single spikes will have less residual effect
- - can no longer do normality tests etc. and remove data that are suspect.
- - requires a dumb stack, and may also require removal of entire stacks of data
-
- Additonally, the GMR varies phase as a function of pulse moment index, so that the first pusle moment is zero phase, the second is pi/2 the third is zero ...
- This however, is altered by the above convention. It gets a little complicated...
- """
-
- import struct
- canvas.reAx2()
-
-
- self.prePulseDelay = 0.01
- self.deadTime = deadTime
- self.samp = 50000.
- self.dt = 1./self.samp
- invGain = 1./self.RxGain
- invCGain = self.CurrentGain
-
-
-
-
- nps = (int)((self.prePulseDelay)*self.samp)
- nps2 = (int)((self.prePulseDelay+self.interpulseDelay)*self.samp)
- npul = (int)(self.pulseLength[0]*self.samp)
- np2 = (int)(self.pulseLength[1]*self.samp)
-
-
- nds = nps+npul+(int)((self.deadTime)*self.samp);
- nd1 = (int)((self.interpulseDelay)*self.samp)
-
- nd2s = nps+npul+nd1+(int)((self.deadTime)*self.samp);
- nd2 = (int)((1.)*self.samp)
- nd1 -= (int)((.028)*self.samp) + nps
-
-
-
- PULSES = [FIDProc]
- if FIDProc == "Both":
- PULSES = ["Pulse 1","Pulse 2"]
-
- self.DATADICT = {}
- self.DATADICT["nPulseMoments"] = self.nPulseMoments
- self.DATADICT["stacks"] = procStacks
- self.DATADICT["PULSES"] = PULSES
- for pulse in PULSES:
- self.DATADICT[pulse] = {}
- self.DATADICT[pulse]["chan"] = chan
- self.DATADICT[pulse]["rchan"] = rchan
- self.DATADICT[pulse]["CURRENT"] = {}
- for ichan in np.append(chan,rchan):
- self.DATADICT[pulse][ichan] = {}
- for ipm in range(self.nPulseMoments):
- self.DATADICT[pulse][ichan][ipm] = {}
- self.DATADICT[pulse]["CURRENT"][ipm] = {}
- for istack in procStacks:
- self.DATADICT[pulse][ichan][ipm][istack] = np.zeros(3)
- self.DATADICT[pulse]["CURRENT"][ipm][istack] = np.zeros(3)
-
-
-
- iistack = 0
- for istack in procStacks:
- rawFile = open(base + "_" + str(istack) + ".lvm", 'rb')
- for ipm in range(self.nPulseMoments):
-
- N_chan = struct.unpack('>i', rawFile.read(4))[0]
- N_samp = struct.unpack('>i', rawFile.read(4))[0]
-
- T = N_samp * self.dt
- TIMES = np.arange(0, T, self.dt) - .0002
-
- DATA = np.zeros([N_samp, N_chan+1])
- for ichan in range(N_chan):
- DATADUMP = rawFile.read(4*N_samp)
- for irec in range(N_samp):
- DATA[irec,ichan] = struct.unpack('>f', DATADUMP[irec*4:irec*4+4])[0]
- if plot:
- canvas.ax1.clear()
- canvas.ax2.clear()
- li = np.shape( DATA[:,4][nd2s:nd2s+nd2] )[0]
-
-
-
-
- for ichan in chan:
- if FIDProc == "Pulse 1":
- self.DATADICT["Pulse 1"][ichan][ipm][istack] = DATA[:,ichan+3][nds:nds+nd1] * invGain
- self.DATADICT["Pulse 1"]["TIMES"] = TIMES[nds:nds+nd1]
- self.DATADICT["Pulse 1"]["CURRENT"][ipm][istack] = DATA[:,1][nps:nps+npul] * invCGain
- self.DATADICT["Pulse 1"]["PULSE_TIMES"] = TIMES[nps:nps+npul]
- if plot:
- canvas.ax1.plot(self.DATADICT["Pulse 1"]["TIMES"], self.DATADICT["Pulse 1"][ichan][ipm][istack], label="Pulse 1 FID data ch. "+str(ichan))
- canvas.ax2.plot(self.DATADICT["Pulse 1"]["PULSE_TIMES"], self.DATADICT["Pulse 1"]["CURRENT"][ipm][istack] , color='black')
- elif FIDProc == "Pulse 2":
- self.DATADICT["Pulse 2"][ichan][ipm][istack] = DATA[:,ichan+3][nd2s:nd2s+nd2] *invGain
- self.DATADICT["Pulse 2"]["TIMES"] = TIMES[nd2s:nd2s+nd2]
- self.DATADICT["Pulse 2"]["CURRENT"][ipm][istack] = DATA[:,1][nps2:nps2+np2] * invCGain
- self.DATADICT["Pulse 2"]["PULSE_TIMES"] = TIMES[nps2:nps2+np2]
- if plot:
- canvas.ax1.plot(self.DATADICT["Pulse 2"]["TIMES"], self.DATADICT["Pulse 2"][ichan][ipm][istack], label="Pulse 2 FID data ch. "+str(ichan))
- canvas.ax2.plot( self.DATADICT["Pulse 2"]["PULSE_TIMES"], self.DATADICT["Pulse 2"]["CURRENT"][ipm][istack], color='black' )
- else:
- self.DATADICT["Pulse 1"][ichan][ipm][istack] = DATA[:,ichan+3][nds:nds+nd1] * invGain
- self.DATADICT["Pulse 2"][ichan][ipm][istack] = DATA[:,ichan+3][nd2s:nd2s+nd2] * invGain
- self.DATADICT["Pulse 1"]["TIMES"] = TIMES[nds:nds+nd1]
- self.DATADICT["Pulse 2"]["TIMES"] = TIMES[nd2s:nd2s+nd2]
- self.DATADICT["Pulse 1"]["CURRENT"][ipm][istack] = DATA[:,1][nps:nps+npul] * invCGain
- self.DATADICT["Pulse 1"]["PULSE_TIMES"] = TIMES[nps:nps+npul]
- self.DATADICT["Pulse 2"]["CURRENT"][ipm][istack] = DATA[:,1][nps2:nps2+np2] * invCGain
- self.DATADICT["Pulse 2"]["PULSE_TIMES"] = TIMES[nps2:nps2+np2]
- if plot:
- canvas.ax1.plot(self.DATADICT["Pulse 1"]["TIMES"], self.DATADICT["Pulse 1"][ichan][ipm][istack], label="Pulse 1 FID data ch. "+str(ichan))
- canvas.ax1.plot(self.DATADICT["Pulse 2"]["TIMES"], self.DATADICT["Pulse 2"][ichan][ipm][istack], label="Pulse 2 FID data ch. "+str(ichan))
- canvas.ax2.plot( self.DATADICT["Pulse 1"]["PULSE_TIMES"], self.DATADICT["Pulse 1"]["CURRENT"][ipm][istack] , color='black' )
- canvas.ax2.plot( self.DATADICT["Pulse 2"]["PULSE_TIMES"], self.DATADICT["Pulse 2"]["CURRENT"][ipm][istack] , color='black')
-
- for ichan in rchan:
- if FIDProc == "Pulse 1":
- self.DATADICT["Pulse 1"][ichan][ipm][istack] = DATA[:,ichan+3][nds:nds+nd1] * invGain
- self.DATADICT["Pulse 1"]["TIMES"] = TIMES[nds:nds+nd1]
- if plot:
- canvas.ax1.plot(self.DATADICT["Pulse 1"]["TIMES"], self.DATADICT["Pulse 1"][ichan][ipm][istack], label="Pulse 1 FID ref ch. "+str(ichan))
- elif FIDProc == "Pulse 2":
- self.DATADICT["Pulse 2"][ichan][ipm][istack] = DATA[:,ichan+3][nd2s:nd2s+nd2] * invGain
- self.DATADICT["Pulse 2"]["TIMES"] = TIMES[nd2s:nd2s+nd2]
- if plot:
- canvas.ax1.plot(self.DATADICT["Pulse 2"]["TIMES"], self.DATADICT["Pulse 2"][ichan][ipm][istack], label="Pulse 2 FID ref ch. "+str(ichan))
- else:
- self.DATADICT["Pulse 1"][ichan][ipm][istack] = DATA[:,ichan+3][nds:nds+nd1] * invGain
- self.DATADICT["Pulse 2"][ichan][ipm][istack] = DATA[:,ichan+3][nd2s:nd2s+nd2] * invGain
- self.DATADICT["Pulse 1"]["TIMES"] = TIMES[nds:nds+nd1]
- self.DATADICT["Pulse 2"]["TIMES"] = TIMES[nd2s:nd2s+nd2]
- if plot:
- canvas.ax1.plot(self.DATADICT["Pulse 1"]["TIMES"], self.DATADICT["Pulse 1"][ichan][ipm][istack], label="Pulse 1 FID ref ch. "+str(ichan))
- canvas.ax1.plot(self.DATADICT["Pulse 2"]["TIMES"], self.DATADICT["Pulse 2"][ichan][ipm][istack], label="Pulse 2 FID ref ch. "+str(ichan))
-
- if plot:
- canvas.ax1.legend(prop={'size':6})
- canvas.ax1.set_title("stack "+str(istack)+" pulse index " + str(ipm), fontsize=8)
- canvas.ax1.set_xlabel("time [s]", fontsize=8)
- canvas.ax1.set_ylabel("RAW signal [V]", fontsize=8)
- canvas.ax2.set_ylabel("Current [A]", fontsize=8)
- canvas.ax2.tick_params(axis='both', which='major', labelsize=8)
- canvas.ax2.tick_params(axis='both', which='minor', labelsize=6)
- canvas.ax2.ticklabel_format(style='sci', scilimits=(0,0), axis='y')
- canvas.ax1.ticklabel_format(style='sci', scilimits=(0,0), axis='y')
-
- canvas.draw()
-
-
- percent = (int) (1e2*((float)((iistack*self.nPulseMoments+ipm+1)) / (len(procStacks)*self.nPulseMoments)))
- self.progressTrigger.emit(percent)
-
- iistack += 1
-
- self.enableDSP()
- self.doneTrigger.emit()
-
-
- if __name__ == "__main__":
-
- if len(sys.argv) < 4:
- print( "mrsurvey path/to/header <stack1> <stackN> ")
- exit()
-
- GMR = GMRDataProcessor()
- GMR.readHeaderFile(sys.argv[1])
- GMR.Print()
-
- if GMR.pulseType == "FID":
- GMR.loadFIDData(sys.argv[1], sys.argv[2], sys.argv[3], 5)
-
- if GMR.pulseType == "4PhaseT1":
- GMR.load4PhaseT1Data(sys.argv[1], sys.argv[2], sys.argv[3], 5)
-
- pylab.show()
|