diff --git a/EventDetectorUI.py b/EventDetectorUI.py index 5874c5f..5777f11 100644 --- a/EventDetectorUI.py +++ b/EventDetectorUI.py @@ -1,366 +1,390 @@ import sys import os from PyQt5.QtWidgets import * from PyQt5 import QtCore from PyQt5.QtGui import QColor, QIntValidator import matplotlib matplotlib.use('QT5Agg') import numpy as np import matplotlib.pylab as plt from matplotlib.backends.backend_qt5agg import FigureCanvasQTAgg as FigureCanvas from matplotlib.backends.backend_qt5agg import NavigationToolbar2QT as NavigationToolbar font = {'weight' : 'regular', 'size' : 4} matplotlib.rc('font', **font) # pass in the font dict as kwargs matplotlib.rc('lines', linewidth = 0.5) from EventDetection import * from Plotting.EventPlots import * +from LoadData import * class AnalysisUI(QWidget): def __init__(self): super().__init__() self.initUI() def AllVariable(self): # Let's dump all the variable we will need here... self.importedFiles = [] self.selectedFiles = [] self.analyzedFiles = [] self.rawtrace = [] + self.rawtracesamplerate = 0 + self.AllEvents = NC.AllEvents() self.wholetracedrawnforthefirsttime = True + self.dsrawtrace = [] + self.dssamplerate = 0 # Defining the plots self.TurnToolbarsOn = True self.fig_singleEvent = plt.figure(1, figsize=(10, 10)) self.ax_singleEvent = self.fig_singleEvent.add_subplot(111) self.ax_singleEvent.plot(np.arange(1, 10), np.arange(1, 10)) self.fig_wholeTrace = plt.figure(2, figsize=(16, 9)) self.ax_wholeTrace = self.fig_wholeTrace.add_subplot(111) self.line_wholeTrace, = self.ax_wholeTrace.plot(np.arange(1, 10), -np.arange(1, 10)) self.figure_singleEvent = FigureCanvas(figure=self.fig_singleEvent) self.figure_wholeTrace = FigureCanvas(figure=self.fig_wholeTrace) if self.TurnToolbarsOn: self.toolbarWholeTrace = NavigationToolbar(self.figure_wholeTrace, self) self.toolbarSingleEvent = NavigationToolbar(self.figure_singleEvent, self) def initUI(self): # Main Window self.AllVariable() self.setGeometry(0, 0, 1200, 600) self.setWindowTitle('A user interface is like a joke...') # Assemble the different layout pieces self.MakeFileImportLayout() self.MakeAnalysisParametersAndRunLayout() self.MakeEventNavigatorLayout() windowLayout = QHBoxLayout() windowLayout.addWidget(self.FileImportLayout) windowLayout.addWidget(self.AnalysisParameters) windowLayout.addWidget(self.EventNavigatorLayout) self.setLayout(windowLayout) ## All Actions go here: self.button_fileimport.clicked.connect(self.ImportButtonPushed) self.button_clearfilelist.clicked.connect(self.ClearListButtonPushed) self.list_filelist.clicked.connect(self.SelectionInFileListChanged) self.list_filelist.doubleClicked.connect(self.FileDoubleClicked) self.button_startAnalysis.clicked.connect(self.AnalysisButtonPressed) self.text_eventNumber.returnPressed.connect(self.EventNumberChanged) self.button_eventForward.clicked.connect(self.ForwardButtonPushed) self.button_eventBackward.clicked.connect(self.BackwardButtonPushed) + self.setting_LPForWholeTracePlot.editingFinished.connect(self.LPForWholeTracePlotChanged) self.show() def MakeFileImportLayout(self): self.FileImportLayout = QGroupBox("File Import") self.FileImportLayout.setFixedWidth(200) layout = QGridLayout() self.button_fileimport = QPushButton('Open Files') self.button_fileimport.setToolTip('Select the files you want to import into the list. You can do multiple selections.') self.button_clearfilelist = QPushButton('Clear List') self.button_clearfilelist.setToolTip('This deletes your current list of files.') self.list_filelist = QListWidget() self.list_filelist.setSelectionMode(QAbstractItemView.ExtendedSelection) self.list_filelist.setToolTip('This is the list of current files. Entries in green are already analyzed. Your selection defines which files are used in the analysis.') layout.addWidget(self.button_fileimport, 0, 0) layout.addWidget(self.list_filelist, 1, 0) layout.addWidget(self.button_clearfilelist, 2, 0) self.FileImportLayout.setLayout(layout) def MakeAnalysisParametersAndRunLayout(self): self.AnalysisParameters = QGroupBox("Analysis") + self.AnalysisParameters.setSizePolicy(QSizePolicy(QSizePolicy.Minimum, QSizePolicy.Minimum)) self.LowPassSettings = QGroupBox("Low Pass Settings") + self.LowPassSettings.setSizePolicy(QSizePolicy(QSizePolicy.Minimum, QSizePolicy.Minimum)) self.CusumSettings = QGroupBox("Cusum Settings") - self.CusumSettings.setSizePolicy(QSizePolicy(QSizePolicy.Minimum, QSizePolicy.Maximum)) #This allows this part to remain small in vertical position when the UI resizes. + self.CusumSettings.setSizePolicy(QSizePolicy(QSizePolicy.Minimum, QSizePolicy.Minimum)) #This allows this part to remain small in vertical position when the UI resizes. self.OtherSettings = QGroupBox("Other Settings") + self.OtherSettings.setSizePolicy(QSizePolicy(QSizePolicy.Minimum, QSizePolicy.Minimum)) #This allows this part to remain small in vertical position when the UI resizes. self.button_startAnalysis = QPushButton('Start Analysis') self.button_startAnalysis.setToolTip('This button will provoke your computer to sel-destruct. Use at your own risk!!!!!') ## Low Pass Settings self.settings_LP_a = QComboBox() availableoptions = [0.99, 0.999, 0.9999] for i in availableoptions: self.settings_LP_a.addItem(str(i)) self.settings_LP_a.setToolTip('How strong should the filter be? Higher number is stronger filtering, which means less found events...') ## EndValueE self.settings_LP_E = QDoubleSpinBox() self.settings_LP_E.setValue(0) self.settings_LP_E.setSingleStep(0.1) self.settings_LP_E.setToolTip('The current needs to go above mean + E*std to terminate the event') ## StartValueS self.settings_LP_S = QDoubleSpinBox() self.settings_LP_S.setValue(5) self.settings_LP_S.setSingleStep(0.1) self.settings_LP_S.setToolTip('The current needs to go below mean + S*std to initiate an event') layout1 = QGridLayout() - layout1.addWidget(QLabel('Low Pass Filter Coefficient (a)'), 0, 0) + layout1.addWidget(QLabel('Filter Coefficient (a)'), 0, 0) layout1.addWidget(self.settings_LP_a, 0, 1) layout1.addWidget(QLabel('Start Coefficient (S)'), 1, 0) layout1.addWidget(self.settings_LP_S, 1, 1) layout1.addWidget(QLabel('End Coefficient (E)'), 3, 0) layout1.addWidget(self.settings_LP_E, 3, 1) self.LowPassSettings.setLayout(layout1) ##Cusum Settings self.settings_cusum_hBook = QDoubleSpinBox() self.settings_cusum_hBook.setValue(1) self.settings_cusum_hBook.setSingleStep(0.1) self.settings_cusum_hBook.setToolTip('Some random stupid parameter') self.settings_cusum_dI = QDoubleSpinBox() self.settings_cusum_dI.setValue(1) self.settings_cusum_dI.setSuffix(' nA') self.settings_cusum_dI.setSingleStep(0.1) self.settings_cusum_dI.setToolTip('The most probable current drop') layoutCusum = QGridLayout() layoutCusum.addWidget(QLabel('h-Book'), 0, 0) layoutCusum.addWidget(self.settings_cusum_hBook, 0, 1) layoutCusum.addWidget(QLabel('Delta I'), 1, 0) layoutCusum.addWidget(self.settings_cusum_dI, 1, 1) self.CusumSettings.setLayout(layoutCusum) ## Other Settings self.settings_other_ChimeraLP = QDoubleSpinBox() self.settings_other_ChimeraLP.setValue(10) self.settings_other_ChimeraLP.setSingleStep(10) self.settings_other_ChimeraLP.setSuffix(' kHz') self.settings_other_ChimeraLP.setToolTip('Low pass setting for high-bandwidth Chimera Data') self.settings_other_MaxEventLength = QDoubleSpinBox() self.settings_other_MaxEventLength.setValue(0.2) self.settings_other_MaxEventLength.setSingleStep(0.05) self.settings_other_MaxEventLength.setSuffix(' s') self.settings_other_MaxEventLength.setToolTip('maximal event length to be considered an event') self.settings_other_MinEventLength = QDoubleSpinBox() self.settings_other_MinEventLength.setValue(600) self.settings_other_MinEventLength.setSingleStep(100) self.settings_other_MinEventLength.setSuffix(' us') self.settings_other_MinEventLength.setToolTip('maximal event length to be considered an impulse') self.settings_other_fitLength = QDoubleSpinBox() self.settings_other_fitLength.setValue(3) self.settings_other_fitLength.setSingleStep(0.5) self.settings_other_fitLength.setSuffix(' ms') self.settings_other_fitLength.setToolTip('minimal event length to be fitted for impulses') layoutOther = QGridLayout() - layoutOther.addWidget(QLabel('Chimera Low Pass'), 0, 0) + layoutOther.addWidget(QLabel('Chimera LP'), 0, 0) layoutOther.addWidget(self.settings_other_ChimeraLP, 0, 1) layoutOther.addWidget(QLabel('Max event length'), 1, 0) layoutOther.addWidget(self.settings_other_MaxEventLength, 1, 1) layoutOther.addWidget(QLabel('Min event length'), 2, 0) layoutOther.addWidget(self.settings_other_MinEventLength, 2, 1) layoutOther.addWidget(QLabel('Fit length'), 3, 0) layoutOther.addWidget(self.settings_other_fitLength, 3, 1) self.OtherSettings.setLayout(layoutOther) ##Final Assembly layout = QGridLayout() layout.addWidget(self.LowPassSettings, 0, 0) layout.addWidget(self.CusumSettings, 1, 0) layout.addWidget(self.OtherSettings, 2, 0) layout.addWidget(self.button_startAnalysis, 3, 0) self.AnalysisParameters.setLayout(layout) def MakeEventNavigatorLayout(self): self.EventNavigatorLayout = QGroupBox("Event Navigator") self.navigation_buttons =QButtonGroup() self.button_eventForward = QPushButton('Forward') self.button_eventForward.setToolTip('Advance one event forward') self.button_eventBackward = QPushButton('Backward') self.button_eventBackward.setToolTip('Go one event back') + self.setting_LPForWholeTracePlot = QDoubleSpinBox() + self.setting_LPForWholeTracePlot.setValue(10) + self.setting_LPForWholeTracePlot.setSingleStep(0.5) + self.setting_LPForWholeTracePlot.setSuffix(' kHz') self.text_eventNumber = QLineEdit('0') self.text_eventNumber.setValidator(QIntValidator(0, 100000)) self.text_eventNumber.setToolTip('This is the number of the current event. You can edit this field to jump directly to an event.') self.text_eventNumber.setAlignment(QtCore.Qt.AlignCenter) layout = QGridLayout() - layout.addWidget(self.figure_singleEvent, 0, 0, 1, 3) + layout.addWidget(self.figure_singleEvent, 0, 0, 1, 5) self.figure_singleEvent.setFixedHeight(250) - layout.addWidget(self.button_eventForward, 2, 2) - layout.addWidget(self.text_eventNumber, 2, 1) - layout.addWidget(self.button_eventBackward, 2, 0) - layout.addWidget(self.figure_wholeTrace, 3, 0, 1, 3) - self.figure_wholeTrace.setFixedHeight(300) + self.figure_singleEvent.setFixedWidth(600) + layout.addWidget(self.button_eventForward, 2, 4) + layout.addWidget(self.text_eventNumber, 2, 3) + layout.addWidget(self.button_eventBackward, 2, 2) + layout.addWidget(self.setting_LPForWholeTracePlot, 2, 1) + layout.addWidget(QLabel('Downsample plot below: '), 2, 0) + layout.addWidget(self.figure_wholeTrace, 3, 0, 1, 5) + self.figure_wholeTrace.setFixedHeight(250) + self.figure_wholeTrace.setFixedWidth(600) + if self.TurnToolbarsOn: - layout.addWidget(self.toolbarSingleEvent, 1, 0, 1, 3) - layout.addWidget(self.toolbarWholeTrace, 4, 0, 1, 3) + layout.addWidget(self.toolbarSingleEvent, 1, 0, 1, 5) + layout.addWidget(self.toolbarWholeTrace, 4, 0, 1, 5) self.EventNavigatorLayout.setLayout(layout) ## Button Actions def ImportButtonPushed(self, event): files = QFileDialog.getOpenFileNames(self, 'Select Files to Add to list', filter = "Images (*.dat *.log)") ## Add to file list if unique: for i in files[0]: if i not in self.importedFiles: self.importedFiles.append(i) self.UpdateFileList() def ClearListButtonPushed(self, event): self.importedFiles = [] self.UpdateFileList() def UpdateFileList(self): self.list_filelist.clear() for ind, i in enumerate(self.importedFiles): self.list_filelist.addItem(os.path.split(i)[1]) folder = os.path.dirname(i) + os.sep + 'analysisfiles' filename, file_extension = os.path.splitext(os.path.basename(i)) potentialanalysisfile = folder + os.sep + filename + 'data' print(potentialanalysisfile) print(os.path.exists(potentialanalysisfile)) if os.path.isfile(potentialanalysisfile + '.dat') or os.path.isfile(potentialanalysisfile + '.db'): ## If file is present, make the row green. This means analysis was done on it. self.list_filelist.item(ind).setBackground(QColor('green')) def SelectionInFileListChanged(self, event): items = self.list_filelist.selectedItems() x = [] for i in range(len(items)): x.append(str(self.list_filelist.selectedItems()[i].text())) self.selectedFiles = x print('The following files are selected in the list:') print(x) def GetFullFilePath(self, listoffiles): fullfilepaths = [] for i in listoffiles: for j in self.importedFiles: if i in j: fullfilepaths.append(j) return fullfilepaths def AnalysisButtonPressed(self): # Get Full File Paths: fullfilepaths = self.GetFullFilePath(self.selectedFiles) print(fullfilepaths) # Fill coefficient dictionary from the user inputs coefficients = {} coefficients['maxEventLength'] = self.settings_other_MaxEventLength.value() coefficients['minEventLength'] = self.settings_other_MinEventLength.value() coefficients['fitLength'] = self.settings_other_fitLength.value() coefficients['delta'] = self.settings_cusum_dI.value() coefficients['hbook'] = self.settings_cusum_hBook.value() coefficients['a'] = np.float(self.settings_LP_a.currentText()) coefficients['S'] = self.settings_LP_S.value() coefficients['E'] = self.settings_LP_E.value() coefficients['ChimeraLowPass'] = self.settings_other_ChimeraLP.value() coefficients['dt'] = 25 coefficients['deltaRel'] = None for i in fullfilepaths: eventdetectionwithfilegeneration(i, coefficients, forceRun=True, verboseLevel=1) self.UpdateFileList() def BackwardButtonPushed(self): if (int(self.text_eventNumber.text())-1) < 0: self.text_eventNumber.setText('0') else: self.text_eventNumber.setText(str(int(self.text_eventNumber.text())-1)) self.EventNumberChanged() def ForwardButtonPushed(self): if (int(self.text_eventNumber.text())+1) >= len(self.AllEvents.events): self.text_eventNumber.setText(str(len(self.AllEvents.events)-1)) else: self.text_eventNumber.setText(str(int(self.text_eventNumber.text())+1)) self.EventNumberChanged() def EventNumberChanged(self): print(int(self.text_eventNumber.text())) if len(self.AllEvents.events): self.DrawEvent(int(self.text_eventNumber.text())) self.DrawWholeTrace(int(self.text_eventNumber.text())) def DrawEvent(self, number): self.ax_singleEvent.clear() PlotEvent(self.AllEvents.events[number], ax=self.ax_singleEvent, savefile=os.getcwd(), showCUSUM=True, - showCurrent=False, showButtons=False, axisFormatter=False) + showCurrent=False, showButtons=False, axisFormatter=False, plotTitleBool=False) self.figure_singleEvent.draw() def DrawWholeTrace(self, number): - #self.ax_wholeTrace.clear() - ShowEventInTrace_SignalPreloaded(self.rawtrace, self.AllEvents, number, self.ax_wholeTrace, line=self.line_wholeTrace, firstCall=self.wholetracedrawnforthefirsttime) + ShowEventInTrace_SignalPreloaded(self.dsrawtrace, self.AllEvents, number, self.ax_wholeTrace, line=self.line_wholeTrace, firstCall=self.wholetracedrawnforthefirsttime, dscorrection=self.dssamplerate/self.rawtracesamplerate) self.wholetracedrawnforthefirsttime = False #self.ax_wholeTrace.autoscale(enable=True) self.figure_wholeTrace.draw() #self.figure_wholeTrace.flush_events() def FileDoubleClicked(self, event): DoubleclickedFile = self.list_filelist.selectedItems()[0].text() # Get Analysis file fullfilepath = self.GetFullFilePath([DoubleclickedFile])[0] filename, file_extension = os.path.splitext(os.path.basename(fullfilepath)) folder = os.path.dirname(fullfilepath) + os.sep + 'analysisfiles' analysisfilepath = folder + os.sep + filename + 'data' print('Analysis File Path: {}'.format(analysisfilepath)) if os.path.isfile(analysisfilepath + '.dat') or os.path.isfile(analysisfilepath + '.db'): ## Display the analysis print('Lets show the analysis') ## Load the analysis file self.AllEvents = LoadEvents(analysisfilepath) self.text_eventNumber.setValidator(QIntValidator(0, len(self.AllEvents.events)-1)) - self.rawtrace = LoadData.OpenFile(self.AllEvents.events[0].filename, ChimeraLowPass = 10e3, approxImpulseResponse = True, Split = True, verbose = False)['i1'] + inputfile = LoadData.OpenFile(self.AllEvents.events[0].filename, ChimeraLowPass = 10e3, approxImpulseResponse = True, Split = True, verbose = False) + self.rawtrace = inputfile['i1'] + self.rawtracesamplerate = inputfile['samplerate'] + self.LPForWholeTracePlotChanged() self.DrawEvent(int(self.text_eventNumber.text())) self.DrawWholeTrace(int(self.text_eventNumber.text())) print(self.AllEvents) else: error = QMessageBox.critical(self, 'Important Message', 'Please select a file that has been analyzed (turned green).', QMessageBox.Retry) self.wholetracedrawnforthefirsttime = True + def LPForWholeTracePlotChanged(self): + print(self.setting_LPForWholeTracePlot.value()) + self.dsrawtrace, self.dssamplerate = LowPassAndResample(self.rawtrace, self.rawtracesamplerate, self.setting_LPForWholeTracePlot.value() * 1e3, LPtoSR=2) + def closeEvent(self, event): ''' reply = QMessageBox.question(self, 'Message', "Are you sure to quit? Don't make me sad... :( :( :(", QMessageBox.Yes | QMessageBox.No, QMessageBox.No) if reply == QMessageBox.Yes: event.accept() else: event.ignore() ''' ## Execution if __name__ == '__main__': app = QApplication(sys.argv) ex = AnalysisUI() sys.exit(app.exec_()) diff --git a/LoadData.py b/LoadData.py index 177c2ea..e077d20 100644 --- a/LoadData.py +++ b/LoadData.py @@ -1,316 +1,327 @@ import os import platform import shelve from tkinter import filedialog import h5py import numpy as np import pyabf import scipy import scipy.signal as sig from PyQt5 import QtGui from scipy import io from scipy import signal import Functions def ImportABF(datafilename): abf = pyabf.ABF(datafilename) #abf.info() # shows what is available #output={'type': 'Clampfit', 'graphene': 0, 'samplerate': abf.pointsPerSec, 'i1': -20000./65536 * abf.dataY, 'v1': abf.dataC, 'filename': datafilename} output = {'type': 'Clampfit', 'graphene': 0, 'samplerate': abf.dataRate, 'i1': abf.data[0] * 1e-12, 'v1': abf.data[1], 'filename': datafilename} return output def ImportAxopatchData(datafilename): x=np.fromfile(datafilename, np.dtype('>f4')) f=open(datafilename, 'rb') graphene=0 for i in range(0, 10): a=str(f.readline()) #print(a) if 'Acquisition' in a or 'Sample Rate' in a: samplerate=int(''.join(i for i in a if i.isdigit()))/1000 if 'FEMTO preamp Bandwidth' in a: femtoLP=int(''.join(i for i in a if i.isdigit())) if 'I_Graphene' in a: graphene=1 print('This File Has a Graphene Channel!') end = len(x) if graphene: #pore current i1 = x[250:end-3:4] #graphene current i2 = x[251:end-2:4] #pore voltage v1 = x[252:end-1:4] #graphene voltage v2 = x[253:end:4] print('The femto was set to : {} Hz, if this value was correctly entered in the LabView!'.format(str(femtoLP))) output={'FemtoLowPass': femtoLP, 'type': 'Axopatch', 'graphene': 1, 'samplerate': samplerate, 'i1': i1, 'v1': v1, 'i2': i2, 'v2': v2, 'filename': datafilename} else: i1 = np.array(x[250:end-1:2]) v1 = np.array(x[251:end:2]) output={'type': 'Axopatch', 'graphene': 0, 'samplerate': samplerate, 'i1': i1, 'v1': v1, 'filename': datafilename} return output def ImportChimeraRaw(datafilename): matfile=io.loadmat(str(os.path.splitext(datafilename)[0])) #buffersize=matfile['DisplayBuffer'] data = np.fromfile(datafilename, np.dtype('1e-6: if verbose: print('converting to SI units') output['i1']=1e-9*output['i1'] output['v1']=1e-3*output['v1'] elif datafilename[-3::] == 'abf': output = ImportABF(datafilename) if verbose: print('length: ' + str(len(output['i1']))) st = os.stat(datafilename) if platform.system() == 'Darwin': if verbose: print('Platform is ' + platform.system()) output['TimeFileWritten'] = st.st_birthtime output['TimeFileLastModified'] = st.st_mtime output['ExperimentDuration'] = st.st_mtime - st.st_birthtime elif platform.system() == 'Windows': if verbose: print('Platform is Windows') output['TimeFileWritten'] = st.st_ctime output['TimeFileLastModified'] = st.st_mtime output['ExperimentDuration'] = st.st_mtime - st.st_ctime else: if verbose: print('Platform is ' + platform.system() + ', might not get accurate results.') try: output['TimeFileWritten'] = st.st_ctime output['TimeFileLastModified'] = st.st_mtime output['ExperimentDuration'] = st.st_mtime - st.st_ctime except: raise Exception('Platform not detected') return output def SaveToHDF5(inp_file, AnalysisResults, coefficients, outdir): file = str(os.path.split(inp_file['filename'])[1][:-4]) f = h5py.File(outdir + file + '_OriginalDB.hdf5', "w") general = f.create_group("General") general.create_dataset('FileName', data=inp_file['filename']) general.create_dataset('Samplerate', data=inp_file['samplerate']) general.create_dataset('Machine', data=inp_file['type']) general.create_dataset('TransverseRecorded', data=inp_file['graphene']) general.create_dataset('TimeFileWritten', data=inp_file['TimeFileWritten']) general.create_dataset('TimeFileLastModified', data=inp_file['TimeFileLastModified']) general.create_dataset('ExperimentDuration', data=inp_file['ExperimentDuration']) segmentation_LP = f.create_group("LowPassSegmentation") for k,l in AnalysisResults.items(): set1 = segmentation_LP.create_group(k) lpset1 = set1.create_group('LowPassSettings') for o, p in coefficients[k].items(): lpset1.create_dataset(o, data=p) for m, l in AnalysisResults[k].items(): if m is 'AllEvents': eventgroup = set1.create_group(m) for i, val in enumerate(l): eventgroup.create_dataset('{:09d}'.format(i), data=val) elif m is 'Cusum': eventgroup = set1.create_group(m) for i1, val1 in enumerate(AnalysisResults[k]['Cusum']): cusevent = eventgroup.create_group('{:09d}'.format(i1)) cusevent.create_dataset('NumberLevels', data=np.uint64(len(AnalysisResults[k]['Cusum'][i1]['levels']))) if len(AnalysisResults[k]['Cusum'][i1]['levels']): cusevent.create_dataset('up', data=AnalysisResults[k]['Cusum'][i1]['up']) cusevent.create_dataset('down', data=AnalysisResults[k]['Cusum'][i1]['down']) cusevent.create_dataset('both', data=AnalysisResults[k]['Cusum'][i1]['both']) cusevent.create_dataset('fit', data=AnalysisResults[k]['Cusum'][i1]['fit']) # 0: level number, 1: current, 2: length, 3: std cusevent.create_dataset('levels_current', data=AnalysisResults[k]['Cusum'][i1]['levels'][1]) cusevent.create_dataset('levels_length', data=AnalysisResults[k]['Cusum'][i1]['levels'][2]) cusevent.create_dataset('levels_std', data=AnalysisResults[k]['Cusum'][i1]['levels'][3]) else: set1.create_dataset(m, data=l) def SaveVariables(savename, **kwargs): if os.path.isdir(savename): savefile=os.path.join(savename,os.path.basename(savename)+'_Events') else: #cut of .dat extension if savename.lower().endswith('.dat'): savename = savename[0:-4] #Check if file already exists, otherwise popup dialog if os.path.isfile(savename + '.dat'): #root = tkinter.Tk() #root.withdraw() savename = filedialog.asksaveasfile(mode='w', defaultextension=".dat") if savename is None: # asksaveasfile return `None` if dialog closed with "cancel". return savefile=base=os.path.splitext(savename.name)[0] # raise IOError('File ' + savename + '.dat already exists.') else: savefile = savename #Check if directory exists directory = os.path.dirname(savefile) if not os.path.exists(directory): os.makedirs(directory) shelfFile=shelve.open(savefile) for arg_name in kwargs: shelfFile[arg_name]=kwargs[arg_name] shelfFile.close() print('saved as: ' + savefile + '.dat') def LoadVariables(loadname, variableName): if not isinstance(loadname, str): raise Exception('The second argument must be a string') #cut of .dat extension if loadname.lower().endswith('.dat'): loadname = loadname[0:-4] if not os.path.isfile(loadname + '.dat'): raise Exception('File does not exist') shelfFile = shelve.open(loadname) try: Variable = shelfFile[variableName] except KeyError as e: message = 'Key ' + variableName + ' does not exist, available Keys: \n' + "\n".join(list(shelfFile.keys())) print(message) raise else: shelfFile.close() print('Loaded ' + variableName + 'from ' + loadname + '.dat') return Variable + +def LowPassAndResample(inpsignal, samplerate, LP, LPtoSR = 5): + Wn = np.round(2 * LP / samplerate, 4) # [0,1] nyquist frequency + b, a = signal.bessel(4, Wn, btype='low', analog=False) # 4-th order digital filter + z, p, k = signal.tf2zpk(b, a) + eps = 1e-9 + r = np.max(np.abs(p)) + approx_impulse_len = int(np.ceil(np.log(eps) / np.log(r))) + Filt_sig = (signal.filtfilt(b, a, inpsignal, method='gust', irlen=approx_impulse_len)) + ds_factor = np.ceil(samplerate / (2 * LP)) + return (scipy.signal.resample(Filt_sig, int(len(inpsignal) / ds_factor)), samplerate / ds_factor) diff --git a/Plotting/EventPlots.py b/Plotting/EventPlots.py index 726dd9b..c326319 100644 --- a/Plotting/EventPlots.py +++ b/Plotting/EventPlots.py @@ -1,760 +1,773 @@ import numpy as np import matplotlib import matplotlib.pyplot as plt import matplotlib.patches as patches from matplotlib.ticker import EngFormatter from matplotlib.widgets import CheckButtons, Button import argparse import platform import shelve import os import Functions import LoadData import NanoporeClasses as NC from bokeh.models import Legend,LegendItem,Range1d from bokeh.palettes import Spectral4 from bokeh.plotting import figure, output_file, show, output_notebook from bokeh.layouts import row, column, gridplot #needed for plotting in jupyter from bokeh.resources import INLINE import bokeh.io Amp = EngFormatter(unit='A', places=2) Time = EngFormatter(unit='s', places=2) Volt = EngFormatter(unit='V', places=2) Cond = EngFormatter(unit='S', places=2) # Extract y and tau out of events def extractytau(filteredevents, showCurrent): if showCurrent: yVals = [event.currentDrop for event in filteredevents] tau = [event.eventLength for event in filteredevents] else: yVals = [event.currentDrop / event.voltage for event in filteredevents if event.voltage > 0] tau = [event.eventLength for event in filteredevents if event.voltage > 0] return tau, yVals def PlotG_tau(translocationEvents, fig=None, savefile=None, showCurrent=False, normalized=False, showCUSUM=True, showLog=False): """ Function used to produce scatter plots and histograms of the events. The figure produced has 3 subplots: Up: Histogram with the number of events per event length Right: Histogram with the number of events per conductance drop [nS] Center: Scatter-plot of all the events wiht in on x-coordinates the event length [s] and in y-coordinates the conductance drop [nS]. In the plots, the events were distributed into the 3 types of events: CUSUM-fitted in red ('Real' type) Non-fitted in blue ('Rough' type) Impulse in green ('Impulse' type) Parameters ---------- events : AllEvents object All the events to be plotted. savefile : str, optional Full path to file where the plots will be saved. showCurrent : bool, optional False by default. If True, it will change the SI unit in the y-axis from siemens [S] to ampers [A]. So instead of units in conductance drop, it will have current drop. normalized : bool, optional False by default. If True, it will change in the y-axis the unit from siemens [S] to normalized current drop without unit. showCUSUM : bool, optional True by default. """ # categorize events in three types # CUSUM fitted events CUSUMEvents = translocationEvents.GetEventTypes('CUSUM') if len(CUSUMEvents) == 0: CUSUMEvents = translocationEvents.GetEventTypes('Real') # Non-fitted events nonFittedEvents = translocationEvents.GetEventTypes('Rough') # Impulse events impulseEvents = translocationEvents.GetEventTypes('Impulse') catEvents = (CUSUMEvents, nonFittedEvents, impulseEvents) # Save figure def SavePlot(event): # Check if directory exists directory = os.path.dirname(savefile) if showCurrent: fig.savefig(directory + os.sep + 'PlotITau.pdf', transparent=True) else: fig.savefig(directory + os.sep + 'PlotGTau.pdf', transparent=True) def ShowLog(event): name = axScatter.xaxis._scale.name bshowlog.label.set_text('Show ' + name) axScatter.set_xscale('linear') if name == 'log' else axScatter.set_xscale('log') # definitions for the axes left, width = 0.15, 0.55 bottom, height = 0.1, 0.6 left_h = left + width + 0.015 bottom_h = bottom + height + 0.015 rect_scatter = [left, bottom, width, height] rect_histx = [left, bottom_h, width, 0.2] rect_histy = [left_h, bottom, 0.2, height] # start with a rectangular Figure if fig is None: fig = plt.figure(1, figsize=(10, 8)) # define axes and link histogram to scatterplot axScatter = fig.add_axes(rect_scatter) axHistx = fig.add_axes(rect_histx, sharex=axScatter) axHisty = fig.add_axes(rect_histy, sharey=axScatter) # Checkboxes to turn on or off events rax = plt.axes([0.75, 0.73, 0.14, 0.15]) visBool = [True, False, True] labelsCheckBox = ('CUSUM-fitted', 'Not fitted', 'impulse') check = CheckButtons(rax, labelsCheckBox, visBool) # Link button to axes to preserve function rax._check = check bax = plt.axes([0.77, 0.9, 0.1, 0.03]) bnext = Button(bax, 'Save figure') bnext.on_clicked(SavePlot) # Link button to axes to preserve function bax._bnext = bnext baxl = plt.axes([0.77, 0.95, 0.1, 0.03]) txt = 'log' if showLog else 'linear' bshowlog = Button(baxl, 'Show ' + txt) bshowlog.on_clicked(ShowLog) # Link button to axes to preserve function bax._bshowlog = bshowlog # Show labels def setlabels(): plt.setp(axHistx.get_xticklabels(), visible=False) plt.setp(axHisty.get_yticklabels(), visible=False) axScatter.set_xlabel('Event length (s)') axScatter.xaxis.set_major_formatter(Time) if normalized: axScatter.set_ylabel('current drop (normalized)') else: if showCurrent: axScatter.set_ylabel('current drop (A)') axScatter.yaxis.set_major_formatter(Amp) else: axScatter.set_ylabel('Conductance drop (G)') axScatter.yaxis.set_major_formatter(Cond) # Set limits allEvents = [event for event in translocationEvents.events if showCurrent or event.voltage is not 0] bins = 100 extra = 0.1 # 0.1 = 10% allTau, allYVals = extractytau(allEvents, showCurrent) yClean = [y for y in allYVals if str(y) != 'nan'] taurangeHist = np.linspace(min(allTau), max(allTau), num=bins) yValsrangeHist = np.linspace(min(yClean), max(yClean), num=bins) # define colors of the 3 classes colors = ['tomato', 'lightgreen','skyblue'] linecolors = ['red','green', 'blue'] # the scatter plot: scatters = [None]*3 def PlotEvents(visBool): # clear axes axScatter.clear() axHistx.clear() axHisty.clear() nrEvents = 0 for i in range(len(catEvents)): # If checkbox is True, plot events if visBool[i]: # Extract Tau and Y events = catEvents[i] tau, yVals = extractytau(events, showCurrent) scatters[i] = axScatter.scatter(tau, yVals, color=colors[i], marker='o', s=30, linewidths=0.1, edgecolors=linecolors[i], picker=5, visible=visBool[i]) axHistx.hist(tau, bins=taurangeHist, color=colors[i], visible=visBool[i]) axHisty.hist(yVals, bins=yValsrangeHist, orientation='horizontal', color=colors[i], visible=visBool[i]) nrEvents += len(events) if showLog: axScatter.set_xscale('log') # set limits if len(allTau) > 0: tauRange = np.max(allTau) - np.min(allTau) yRange = np.max(yClean) - np.min(yClean) axScatter.set_xlim((np.max([np.min(allTau) - extra * tauRange, 1e-6]), np.max(allTau) + extra * tauRange)) axScatter.set_ylim((np.min(yClean) - extra * yRange, np.max(yClean) + extra * yRange)) setlabels() fig.suptitle('{} number of events'.format(nrEvents)) PlotEvents(visBool) # If click on checkbox, switch Boolean and replot events def func(label): for i in range(len(labelsCheckBox)): if label == labelsCheckBox[i]: visBool[i] = not visBool[i] PlotEvents(visBool) plt.draw() # When clicking on event def onpick(event): for i in range(len(catEvents)): if event.artist == scatters[i]: N = len(event.ind) if not N: return True figi = plt.figure(figsize=(10, 6)) for subplotnum, dataind in enumerate(event.ind): ax = figi.add_subplot(N, 1, subplotnum + 1) PlotEvent(catEvents[i][dataind], ax, savefile, showCUSUM) figi.show() return True fig.canvas.mpl_connect('pick_event', onpick) check.on_clicked(func) if not hasattr(fig, '_AnalysisUI__attached'): plt.show() def PlotGTau(eventClass, xLim = None, yLim = None, showCurrent = False): bokeh.io.output_notebook(INLINE) # categorize events in three types # CUSUM fitted events CUSUMEvents = eventClass.GetEventTypes('CUSUM') # Non-fitted events nonFittedEvents = eventClass.GetEventTypes('Rough') # Impulse events impulseEvents = eventClass.GetEventTypes('Impulse') # define of variables TOOLS = "box_zoom,pan,wheel_zoom,reset" colors = ['tomato', 'lightgreen','skyblue'] linecolors = ['red','green', 'blue'] backgroundColor = "#fafafa" bins = 50 output_notebook() p = figure(plot_width=500, plot_height=500, min_border=10, min_border_left=50, tools=TOOLS, x_axis_location=None, y_axis_location=None,title="Linked Histograms") p.background_fill_color = "#fafafa" # select min max events = [event for event in eventClass.events if showCurrent or event.voltage is not 0] allTau, allYVals = extractytau(events, showCurrent) #Set limits if xLim is None: taurange = np.linspace(min(allTau), max(allTau), num=bins) else: assert(len(xLim) is 2) p.x_range=Range1d(xLim[0], xLim[1]) taurange = np.linspace(xLim[0], xLim[1], num=bins) if yLim is None: yValsrange = np.linspace(min(allYVals), max(allYVals), num=bins) else: assert(len(yLim) is 2) p.y_range=Range1d(yLim[0], yLim[1]) yValsrange = np.linspace(yLim[0], yLim[1], num=bins) # initialize values zerostau = np.zeros(len(taurange) - 1) previoustau = np.zeros(len(taurange) - 1) zerosy = np.zeros(len(yValsrange)-1) previousy = np.zeros(len(yValsrange) - 1) # Define histogram bottom ph = figure(plot_width=p.plot_width, plot_height=100, x_range=p.x_range, y_range=(0, 1), min_border=10, min_border_left=50, y_axis_location="right") ph.background_fill_color = backgroundColor ph.xgrid.grid_line_color = None ph.yaxis.major_label_orientation = np.pi / 4 # Define the second histogram pv = figure(plot_width=300, plot_height=p.plot_height, x_range=(0, 1), y_range=p.y_range, min_border=10, y_axis_location="right") pv.ygrid.grid_line_color = None pv.xaxis.major_label_orientation = np.pi / 4 pv.background_fill_color = backgroundColor legend_it = [] for selectEvents, name, color, linecolor in zip([CUSUMEvents, nonFittedEvents, impulseEvents], ["CUSUM", "Non-fitted", "Impulse"], colors, linecolors): tau, yVals = extractytau(selectEvents,showCurrent) c = p.scatter(tau, yVals, size=3, line_width=2, color=color, alpha=0.8) hhist, hedges = np.histogram(tau, bins=taurange) hh = ph.quad(bottom=zerostau +previoustau, left=hedges[:-1], right=hedges[1:], top=hhist+previoustau, fill_color = color,line_color = linecolor) previoustau +=hhist vhist, vedges = np.histogram(yVals, bins=yValsrange) vh = pv.quad(left=zerosy + previousy, bottom=vedges[:-1], top=vedges[1:], right=vhist+previousy, color=color, line_color= linecolor) previousy +=vhist #Sharing legend is broke #legend_it.append((name, [c, hh, vh])) legend_it.append(LegendItem(label=name, renderers=[vh])) #p.legend.location = "top_right" #legend = Legend(items=legend_it, location=(0, -30)) legend = Legend(items=legend_it) pv.add_layout(legend, 'right') #pv.legend.click_policy = "hide" ph.y_range.end = max(previoustau)*1.1 pv.x_range.end = max(previousy) * 1.1 fig = gridplot([[p, pv],[ph, None]], toolbar_location="left") # show the results show(fig) def PlotGTauVoltage (eventClass, xLim=None, yLim=None, showCurrent=False): bokeh.io.output_notebook(INLINE) #sort the voltages # categorize events in three types # CUSUM fitted events voltageLimits = [0.01, 0.91] voltagesList = eventClass.GetAllVoltages() #define of variables TOOLS = "box_zoom,pan,wheel_zoom,reset" colors = ['tomato', 'lightgreen','skyblue', 'magenta','black'] linecolors = ['red', 'green', 'blue', 'magenta', 'black'] assert(len(voltagesList)<=len(colors)) backgroundColor = "#fafafa" bins = 50 output_notebook() p = figure(plot_width=500, plot_height=500, min_border=10, min_border_left=50, tools=TOOLS, x_axis_location=None, y_axis_location=None,title="Linked Histograms") p.background_fill_color = "#fafafa" #select min max events = [event for event in eventClass.events if showCurrent or event.voltage is not 0] allTau, allYVals = extractytau(events, showCurrent) #Set limits if xLim is None: taurange = np.linspace(min(allTau), max(allTau), num=bins) else: assert(len(xLim) is 2) p.x_range=Range1d(xLim[0], xLim[1]) taurange = np.linspace(xLim[0], xLim[1], num=bins) if yLim is None: yValsrange = np.linspace(min(allYVals), max(allYVals), num=bins) else: assert(len(yLim) is 2) p.y_range=Range1d(yLim[0], yLim[1]) yValsrange = np.linspace(yLim[0], yLim[1], num=bins) #initialize values zerostau = np.zeros(len(taurange) - 1) previoustau = np.zeros(len(taurange) - 1) zerosy = np.zeros(len(yValsrange)-1) previousy = np.zeros(len(yValsrange) - 1) #Define histogram bottom ph = figure(plot_width=p.plot_width, plot_height=100, x_range=p.x_range, y_range=(0, 1), min_border=10, min_border_left=50, y_axis_location="right") ph.background_fill_color = backgroundColor ph.xgrid.grid_line_color = None ph.yaxis.major_label_orientation = np.pi / 4 #Define the second histogram pv = figure(plot_width=300, plot_height=p.plot_height, x_range=(0, 1), y_range=p.y_range, min_border=10, y_axis_location="right") pv.ygrid.grid_line_color = None pv.xaxis.major_label_orientation = np.pi / 4 pv.background_fill_color = backgroundColor legend_it = [] i = 0 #for i in range(len(voltagesList)): for voltage in voltagesList: color = colors[i] linecolor = linecolors[i] i+=1 #for voltage, color, linecolor in zip(voltagesList, colors, linecolors): selectEvents = eventClass.GetEventsforVoltages(voltage) tau, yVals = extractytau(selectEvents, showCurrent) c = p.scatter(tau, yVals, size=3, line_width=2, color=color, alpha=0.8) hhist, hedges = np.histogram(tau, bins=taurange) hh = ph.quad(bottom=zerostau +previoustau, left=hedges[:-1], right=hedges[1:], top=hhist+previoustau, fill_color = color,line_color = linecolor) previoustau +=hhist vhist, vedges = np.histogram(yVals, bins=yValsrange) vh = pv.quad(left=zerosy + previousy, bottom=vedges[:-1], top=vedges[1:], right=vhist+previousy, color=color, line_color= linecolor) previousy +=vhist #Sharing legend is broke #legend_it.append((name, [c, hh, vh])) name = str(Volt.format_data(voltage)) legend_it.append(LegendItem(label=name, renderers=[vh])) #p.legend.location = "top_right" #legend = Legend(items=legend_it, location=(0, -30)) legend = Legend(items=legend_it) pv.add_layout(legend, 'right') #pv.legend.click_policy = "hide" ph.y_range.end = max(previoustau)*1.1 pv.x_range.end = max(previousy) * 1.1 fig = gridplot([[p, pv],[ph, None]], toolbar_location="left") # show the results show(fig) -def PlotEvent(event, ax=None, savefile=os.getcwd(), showCUSUM=True, showCurrent=False, showButtons = True, axisFormatter = True): +def PlotEvent(event, ax=None, savefile=os.getcwd(), showCUSUM=True, showCurrent=False, showButtons = True, axisFormatter = True, plotTitleBool = True): """ Function used to plot a single event passed in argument. The event will be represented in a blue trace and the baseline in a red trace. Parameters ---------- event : TranslocationEvent object Event to be plotted. ax :axes.Axes object or array of Axes object Axes of the plot. savefile : str, optional By default, the current working directory. Full path to the directory to save the plot. showCUSUM : bool, optional False by default. If True, it will plot the CUSUM-fit overlayed in yellow. """ #Link event to axes to keep it around if ax is None: #plt.figure(figsize=(10, 6)) fig, ax = plt.subplots(figsize=(10, 6)) ax._event = event def SavePlot(eventMouse): # Check if directory exists directory = os.path.dirname(savefile) savename=directory + os.sep + 'event.pdf' i = 1 while os.path.exists(directory + os.sep + 'event_{}.pdf'.format(i)): i += 1 savename = directory + os.sep + 'event_{}.pdf'.format(i) eventMouse.inaxes.figure.savefig(savename, transparent=True) def ShowFullTrace(eventMouse): event=eventMouse.inaxes.figure.axes[0]._event ShowEventInTrace(event) if showCUSUM and hasattr(event, 'changeTimes') and len(event.changeTimes)>=2: eventLength = event.eventLengthCUSUM currentDrop = event.currentDropCUSUM else: showCUSUM=False eventLength = event.eventLength currentDrop = event.currentDrop fn=filename_w_ext = os.path.basename(event.filename) if showCurrent: plotTitle = fn + '\n' + 'Event length: {}\nCurrent drop: {} with voltage {}'.format( Time.format_data(eventLength), Cond.format_data(currentDrop), Volt.format_data(event.voltage)) else: plotTitle = fn + '\n' + 'Event length: {}\nConductance drop: {} with voltage {}'.\ format(Time.format_data(eventLength), Cond.format_data(currentDrop/event.voltage),Volt.format_data(event.voltage)) ax.set_xlabel('time (s)') ax.set_ylabel('current (A)') if axisFormatter: ax.xaxis.set_major_formatter(Time) ax.yaxis.set_major_formatter(Amp) - if plotTitle: + if plotTitleBool: plt.title(plotTitle) if showButtons: # Add buttons # Save button bax = plt.axes([0.77, 0.95, 0.15, 0.03]) bsave = Button(bax, 'Save figure') bsave.on_clicked(SavePlot) # Link button to axes to preserve function ax._bsave = bsave # Show original trace button bax2 = plt.axes([0.77, 0.9, 0.15, 0.03]) bfull = Button(bax2, 'Show original Trace') # Link button to axes to preserve function ax._bfull = bfull bfull.on_clicked(ShowFullTrace) #Plotting timeVals1 = np.linspace(0, len(event.before) / event.samplerate, num=len(event.before)) timeVals2 = np.linspace(0 + max(timeVals1), len(event.eventTrace) / event.samplerate + max(timeVals1), num=len(event.eventTrace)) timeVals3 = np.linspace(0 + max(timeVals2), len(event.after) / event.samplerate + max(timeVals2), num=len(event.after)) ax.plot(np.append(timeVals1,timeVals2[0]), np.append(event.before,event.eventTrace[0]), color='tomato') ax.plot(timeVals2, event.eventTrace, color='mediumslateblue') ax.plot(np.append(timeVals2[-1],timeVals3), np.append(event.eventTrace[-1],event.after), color='tomato') if showCUSUM: timeVals = np.linspace(0, len(event.segmentedSignal) / event.samplerate, num=len(event.segmentedSignal)) if hasattr(event,'mcbefore') and hasattr(event,'mcafter') and hasattr(event,'mctrace'): ax.plot(timeVals1, event.mcbefore,'--', color='tomato') x=np.append(np.append(timeVals1[-1],timeVals2),timeVals3[0]) y=np.append(np.append(event.mcbefore[-1],event.mctrace),event.mcafter[0]) ax.plot(x, y, color='yellow') ax.plot(timeVals3, event.mcafter,'--', color='tomato') else: ax.plot(timeVals,event.segmentedSignal, color='yellow') #,timeVals3[0],event.mcafter[0] else: beforeBaseline=np.full(len(event.before), event.baseline) ax.plot(timeVals1,beforeBaseline, '--', color='tomato') afterBaseline = np.full(len(event.after), event.baseline) ax.plot(timeVals3,afterBaseline, '--', color='tomato') meanTrace = np.full(len(event.eventTrace), event.baseline-event.currentDrop) ax.plot(timeVals2,meanTrace, '--', color='mediumslateblue') if 'fig' in locals(): plt.show() -def ShowEventInTrace_SignalPreloaded(FullTrace, AllData, eventnumber, ax, line = None, firstCall = True): +def ShowEventInTrace_SignalPreloaded(FullTrace, AllData, eventnumber, ax, line = None, firstCall = True, dscorrection = None): times = np.linspace(0, len(FullTrace) / AllData.events[eventnumber].samplerate, num=len(FullTrace)) if line: line.set_ydata(FullTrace) line.set_xdata(times) if firstCall: ax.set_xlim(np.min(times), np.max(times)) ax.set_ylim(np.min(FullTrace), np.max(FullTrace)) print('updated lines') else: ax.plot(times, FullTrace, zorder=1) print('Updated plot') ax.set_xlabel('time (s)') ax.set_ylabel('current (A)') # Create a Rectangle patch - if hasattr(AllData.events[eventnumber], 'changeTimes') and len(AllData.events[eventnumber].changeTimes) > 2: - start_i = (AllData.events[eventnumber].beginEventCUSUM - len(AllData.events[eventnumber].before)) / AllData.events[eventnumber].samplerate - end_i = (AllData.events[eventnumber].endEventCUSUM + len(AllData.events[eventnumber].after)) / AllData.events[eventnumber].samplerate + if dscorrection: + if hasattr(AllData.events[eventnumber], 'changeTimes') and len(AllData.events[eventnumber].changeTimes) > 2: + start_i = (AllData.events[eventnumber].beginEventCUSUM - len(AllData.events[eventnumber].before)) / AllData.events[eventnumber].samplerate * dscorrection + end_i = (AllData.events[eventnumber].endEventCUSUM + len(AllData.events[eventnumber].after)) / AllData.events[eventnumber].samplerate * dscorrection + else: + start_i = (AllData.events[eventnumber].beginEvent - len(AllData.events[eventnumber].before)) / AllData.events[eventnumber].samplerate * dscorrection + end_i = (AllData.events[eventnumber].endEvent + len(AllData.events[eventnumber].after)) / AllData.events[eventnumber].samplerate * dscorrection else: - start_i = (AllData.events[eventnumber].beginEvent - len(AllData.events[eventnumber].before)) / AllData.events[eventnumber].samplerate - end_i = (AllData.events[eventnumber].endEvent + len(AllData.events[eventnumber].after)) / AllData.events[eventnumber].samplerate + if hasattr(AllData.events[eventnumber], 'changeTimes') and len(AllData.events[eventnumber].changeTimes) > 2: + start_i = (AllData.events[eventnumber].beginEventCUSUM - len(AllData.events[eventnumber].before)) / \ + AllData.events[eventnumber].samplerate + end_i = (AllData.events[eventnumber].endEventCUSUM + len(AllData.events[eventnumber].after)) / \ + AllData.events[eventnumber].samplerate + else: + start_i = (AllData.events[eventnumber].beginEvent - len(AllData.events[eventnumber].before)) / \ + AllData.events[eventnumber].samplerate + end_i = (AllData.events[eventnumber].endEvent + len(AllData.events[eventnumber].after)) / AllData.events[ + eventnumber].samplerate + minE = np.min(np.append(np.append(AllData.events[eventnumber].eventTrace, AllData.events[eventnumber].before), AllData.events[eventnumber].after)) maxE = np.max(np.append(np.append(AllData.events[eventnumber].eventTrace, AllData.events[eventnumber].before), AllData.events[eventnumber].after)) rect = patches.Rectangle((start_i, minE - 0.1 * (maxE - minE)), end_i - start_i, maxE + 0.2 * (maxE - minE) - minE, linestyle='--', linewidth=1, edgecolor='r', facecolor='none', zorder=10) # Add the patch to the Axes print(ax.patches.clear()) ax.add_patch(rect) def ShowEventInTrace(event): """ Function used to show the event with it's location framed in red in the original full signal trace in blue. Parameters ---------- event : TranslocationEvent object Event to be plotted. """ filename = event.filename loadedData = LoadData.OpenFile(filename, 1e3, True) #, ChimeraLowPass, True, CutTraces) fig, ax = plt.subplots(figsize=(10, 6)) FullTrace = loadedData['i1'] times = np.linspace(0, len(FullTrace) / event.samplerate, num=len(FullTrace)) ax.plot(times, FullTrace, zorder=1) ax.set_xlabel('time (s)') ax.set_ylabel('current (A)') ax.xaxis.set_major_formatter(Time) ax.yaxis.set_major_formatter(Amp) # Create a Rectangle patch if hasattr(event,'changeTimes') and len(event.changeTimes)>2: start_i = (event.beginEventCUSUM - len(event.before))/event.samplerate end_i = (event.endEventCUSUM + len(event.after))/event.samplerate else: start_i = (event.beginEvent - len(event.before))/event.samplerate end_i = (event.endEvent + len(event.after))/event.samplerate minE=np.min(np.append(np.append(event.eventTrace,event.before),event.after)) maxE=np.max(np.append(np.append(event.eventTrace,event.before),event.after)) rect = patches.Rectangle((start_i, minE-0.1*(maxE-minE)), end_i-start_i, maxE+0.2*(maxE-minE)-minE, linestyle='--', linewidth=1, edgecolor='r', facecolor='none',zorder=10) # Add the patch to the Axes ax.add_patch(rect) plt.title(os.path.basename(filename)) plt.show() def PlotCurrentTrace(currentTrace, samplerate): """ Function used in TranslocationEvent class methods in the NanoporeClasses module to plot events. It will plot the TranslocationEvent's currrentTrace passed in argument. Parameters ---------- currentTrace : list of float Data points in current to be plotted. samplerate : float Sampling frequency of the data acquisition. """ timeVals = np.linspace(0, len(currentTrace) / samplerate, num=len(currentTrace)) fig,ax=plt.subplots(figsize=(10, 6)) ax.plot(timeVals, currentTrace) ax.set_xlabel('time (s)') ax.set_ylabel('current (A)') ax.xaxis.set_major_formatter(Time) ax.yaxis.set_major_formatter(Amp) plt.show() def PlotCurrentTraceBaseline(before, currentTrace, after, samplerate, plotTitle=''): """ Function used in TranslocationEvent class methods in the NanoporeClasses module to plot events. It will plot TranslocationEvent's currentTrace and the surrounding baseline (after and before) passed in argument. Parameters ---------- before : list of float Data points in current of the baseline trace before the event. currentTrace : list of float Data points in current of the event trace. after : list of float Data points in current of the baseline trace after the event. samplerate : float Sampling frequency of the data aquisition. plotTitle : str, optional Plot title. """ timeVals1 = np.linspace(0, len(before) / samplerate, num=len(before)) timeVals2 = np.linspace(0 + max(timeVals1), len(currentTrace) / samplerate + max(timeVals1), num=len(currentTrace)) timeVals3 = np.linspace(0 + max(timeVals2), len(after) / samplerate + max(timeVals2), num=len(after)) #plt.figure(figsize=(10, 6)) fig,ax=plt.subplots(figsize=(10, 6)) ax.plot(timeVals1, before, color='red') ax.plot(timeVals2, currentTrace) ax.plot(timeVals3, after, color='red') ax.set_xlabel('time (s)') ax.set_ylabel('current (A)') ax.xaxis.set_major_formatter(Time) ax.yaxis.set_major_formatter(Amp) if plotTitle: plt.title(plotTitle) plt.show() if __name__=='__main__': from tkinter import Tk from tkinter.filedialog import askopenfilenames, askdirectory matplotlib.use('TkAgg') if (platform.system() == 'Darwin'): root = Tk() root.withdraw() if (platform.system() == 'Darwin'): os.system('''/usr/bin/osascript -e 'tell app "Finder" to set frontmost of process "python" to true' ''') parser = argparse.ArgumentParser() parser.add_argument('-i', '--input', help='Input file') args = parser.parse_args() inputData = args.input if inputData == None: inputData = askopenfilenames(filetypes=[('data files', '*.dat')]) #for Mac systems, replace 'Data*.dat' with >> '*.dat' #if inputData: # inputData=os.path.splitext(inputData[0])[0] if (platform.system() == 'Darwin'): root.update() translocationEvents = NC.AllEvents() if inputData: for filename in inputData: shelfFile = shelve.open(os.path.splitext(filename)[0]) translocationEventstemp = shelfFile['TranslocationEvents'] shelfFile.close() translocationEvents.AddEvent(translocationEventstemp) PlotG_tau(translocationEvents, savefile=inputData)