From: W. Trevor King Date: Thu, 17 Dec 2009 19:32:21 +0000 (-0500) Subject: Removed sha-bang from non-executable python files + whitespace cleanups. X-Git-Url: http://git.tremily.us/?p=hooke.git;a=commitdiff_plain;h=c0bd161ab41e18f33f4497b8e3206f0a12bf9616 Removed sha-bang from non-executable python files + whitespace cleanups. Whitespace cleanups mostly involved removing trailing whitespace. --- diff --git a/hooke/driver/csvdriver.py b/hooke/driver/csvdriver.py index 4222255..fb465ee 100644 --- a/hooke/driver/csvdriver.py +++ b/hooke/driver/csvdriver.py @@ -1,12 +1,10 @@ -#!/usr/bin/env python - ''' csvdriver.py Simple driver to read general comma-separated values in Hooke Columns are read this way: - + X1 , Y1 , X2 , Y2 , X3 , Y3 ... If the number of columns is odd, the last column is ignored. @@ -19,57 +17,55 @@ import libhooke as lh import csv class csvdriverDriver(lhc.Driver): - + def __init__(self, filename): - + self.filedata = open(filename,'r') - self.data = list(self.filedata) + self.data = list(self.filedata) self.filedata.close() - + self.filetype = 'generic' self.experiment = '' - + self.filename=filename - + def is_me(self): myfile=file(self.filename) headerline=myfile.readlines()[0] myfile.close() - + #using a custom header makes things much easier... #(looking for raw CSV data is at strong risk of confusion) if headerline[:-1]=='Hooke data': return True else: return False - + def close_all(self): self.filedata.close() - + def default_plots(self): rrows=csv.reader(self.data) rows=list(rrows) #transform the csv.reader iterator in a normal list columns=lh.transposed2(rows[1:]) - + main_plot=lhc.PlotObject() main_plot.vectors=[] - + for index in range(0,len(columns),2): main_plot.vectors.append([]) temp_x=columns[index] temp_y=columns[index+1] - + #convert to float (the csv gives strings) temp_x=[float(item) for item in temp_x] temp_y=[float(item) for item in temp_y] - + main_plot.vectors[-1].append(temp_x) main_plot.vectors[-1].append(temp_y) - + main_plot.units=['x','y'] main_plot.title=self.filename main_plot.destination=0 - + return [main_plot] - - \ No newline at end of file diff --git a/hooke/driver/hemingclamp.py b/hooke/driver/hemingclamp.py index e2f2e1c..41e721a 100755 --- a/hooke/driver/hemingclamp.py +++ b/hooke/driver/hemingclamp.py @@ -1,11 +1,9 @@ -#!/usr/bin/env python - ''' libhemingclamp.py Library for interpreting Hemingway force spectroscopy files. -Copyright (C) 2008 Massimo Sandal, Marco Brucale (University of Bologna, Italy) +Copyright (C) 2008 Massimo Sandal, Marco Brucale (University of Bologna, Italy) This program is released under the GNU General Public License version 2. ''' @@ -16,35 +14,35 @@ __changelog__=''' 2007_02_07: Initial implementation ''' import string -import libhookecurve as lhc +import libhookecurve as lhc class DataChunk(list): '''Dummy class to provide ext and ret methods to the data list. In this case ext and self can be equal. ''' - + def ext(self): return self - + def ret(self): return self class hemingclampDriver(lhc.Driver): - + def __init__(self, filename): - + self.filedata = open(filename,'r') self.data = self.filedata.readlines()[6:] self.filedata.close() - + self.filetype = 'hemingclamp' self.experiment = 'clamp' - + self.filename=filename - + def __del__(self): - self.filedata.close() - + self.filedata.close() + def is_me(self): ''' we define our magic heuristic for HemingClamp files @@ -56,7 +54,7 @@ class hemingclampDriver(lhc.Driver): return True else: return False - + def _getdata_all(self): time = [] phase = [] @@ -65,7 +63,7 @@ class hemingclampDriver(lhc.Driver): imposed = [] trim_indexes = [] trim_counter = 0.0 - + for i in self.data: temp = string.split(i) #time.append(float(temp[0])*(1.0e-3)) # This is managed differently now, since each data point = 1ms: see below @@ -78,24 +76,24 @@ class hemingclampDriver(lhc.Driver): if phase[x] != trim_counter: trim_indexes.append(x) trim_counter = phase[x] - + #we rebuild the time counter assuming 1 point = 1 millisecond c=0.0 for z in zpiezo: time.append(c) - c+=(1.0e-3) - + c+=(1.0e-3) + return time,phase,zpiezo,defl,imposed,trim_indexes - + def time(self): return DataChunk(self._getdata_all()[0]) def phase(self): return DataChunk(self._getdata_all()[1]) - + def zpiezo(self): return DataChunk(self._getdata_all()[2]) - + def deflection(self): return DataChunk(self._getdata_all()[3]) @@ -104,31 +102,31 @@ class hemingclampDriver(lhc.Driver): def trimindexes(self): return DataChunk(self._getdata_all()[5]) - + def close_all(self): ''' Explicitly closes all files ''' self.filedata.close() - + def default_plots(self): main_plot=lhc.PlotObject() defl_plot=lhc.PlotObject() - + time=self.time() phase=self.phase() zpiezo=self.zpiezo() deflection=self.deflection() imposed=self.imposed() - + main_plot.vectors=[[time,zpiezo],[time,phase]] main_plot.units=['seconds','meters'] main_plot.destination=0 main_plot.title=self.filename - + defl_plot.vectors=[[time,deflection],[time,imposed]] defl_plot.units=['seconds','Newtons'] defl_plot.destination=1 - + return [main_plot, defl_plot] - \ No newline at end of file + diff --git a/hooke/driver/jpk.py b/hooke/driver/jpk.py index 4f64331..5fd5f5d 100644 --- a/hooke/driver/jpk.py +++ b/hooke/driver/jpk.py @@ -1,15 +1,13 @@ -#!/usr/bin/env python - import string -import libhookecurve as lhc +import libhookecurve as lhc class DataChunk(list): #Dummy class to provide ext and ret methods to the data list. - + def ext(self): halflen=(len(self)/2) return self[0:halflen] - + def ret(self): halflen=(len(self)/2) return self[halflen:] @@ -20,26 +18,26 @@ class jpkDriver(lhc.Driver): self.filename=filename #self.filename can always be useful, and should be defined self.filedata = open(filename,'r') #We open the file self.filelines=self.filedata.readlines() - self.filedata.close() + self.filedata.close() '''These are two strings that can be used by Hooke commands/plugins to understand what they are looking at. They have no other meaning. They have to be somehow defined however - commands often look for those variables. - + self.filetype should contain the name of the exact filetype defined by the driver (so that filetype-specific commands can know if they're dealing with the correct filetype) self.experiment should contain instead the type of data involved (for example, various drivers can be used for force-clamp experiments, - but hooke commands could like to know if we're looking at force clamp data, regardless of their origin, and not other + but hooke commands could like to know if we're looking at force clamp data, regardless of their origin, and not other kinds of data) - + Of course, all other variables you like can be defined in the class. ''' self.filetype = 'jpk' self.experiment = 'smfs' - - + + def __del__(self): - self.filedata.close() - + self.filedata.close() + def is_me(self): ''' we define our magic heuristic for jpk files @@ -51,24 +49,24 @@ class jpkDriver(lhc.Driver): return True else: return False - + def close_all(self): self.filedata.close() - + def _read_data_segment(self): #routine that actually reads the data - + height_ms=[] height_m=[] height=[] v_deflection=[] h_deflection=[] - + self.springconstant=0 #if we don't meet any spring constant, use deflection... - + for line in self.filelines: #we meet the segment defining the order of data columns - + if line[0:9]=='# columns': splitline=line.split()[2:] height_ms_index=splitline.index('smoothedStrainGaugeHeight') @@ -76,10 +74,10 @@ class jpkDriver(lhc.Driver): height_index=splitline.index('height') v_deflection_index=splitline.index('vDeflection') #h_deflection=splitline.index('hDeflection') - + if line[0:16]=='# springConstant': self.springconstant=float(line.split()[2]) - + if line[0] != '#' and len(line.split())>1: dataline=line.split() height_ms.append(float(dataline[height_ms_index])) @@ -87,31 +85,31 @@ class jpkDriver(lhc.Driver): height.append(float(dataline[height_index])) v_deflection.append(float(dataline[v_deflection_index])) #h_deflection.append(float(dataline[h_deflection_index])) - + if self.springconstant != 0: force=[item*self.springconstant for item in v_deflection] else: #we have measured no spring constant :( force=v_deflection - + height_ms=DataChunk([item*-1 for item in height_ms]) height_m=DataChunk([item*-1 for item in height_m]) height=DataChunk([item*-1 for item in height]) deflection=DataChunk(v_deflection) force=DataChunk(force) - + return height_ms,height_m,height,deflection,force - + def deflection(self): height_ms,height_m,height,deflection,force=self._read_data_segment() deflection_ext=deflection.ext() deflection_ret=deflection.ret() deflection_ret.reverse() return deflection_ext,deflection_ret - + def default_plots(self): - + height_ms,height_m,height,deflection,force=self._read_data_segment() - + height_ms_ext=height_ms.ext() height_ms_ret=height_ms.ret() force_ext=force.ext() @@ -119,21 +117,21 @@ class jpkDriver(lhc.Driver): #reverse the return data, to make it coherent with hooke standard height_ms_ret.reverse() force_ret.reverse() - - main_plot=lhc.PlotObject() + + main_plot=lhc.PlotObject() main_plot.add_set(height_ms_ext,force_ext) main_plot.add_set(height_ms_ret,force_ret) - - - + + + if self.springconstant != 0: main_plot.units=['meters','force'] else: main_plot.units=['meters','meters'] - + main_plot.normalize_vectors() - + main_plot.destination=0 main_plot.title=self.filename - - return [main_plot] \ No newline at end of file + + return [main_plot] diff --git a/hooke/driver/mcs.py b/hooke/driver/mcs.py index 711e8c3..48836f2 100644 --- a/hooke/driver/mcs.py +++ b/hooke/driver/mcs.py @@ -1,5 +1,3 @@ -#!/usr/bin/env python - ''' mcs.py @@ -13,7 +11,7 @@ import libhooke as lh import struct class mcsDriver(lhc.Driver): - + def __init__(self, filename): ''' Open the RED (A) ones; the BLUE (D) mirror ones will be automatically opened @@ -25,37 +23,37 @@ class mcsDriver(lhc.Driver): oth[-8]='d' othername=''.join(oth) self.filename=filename - self.othername=othername - + self.othername=othername + #print self.filename, self.othername - + self.filedata=open(filename,'rb') self.reddata=self.filedata.read() self.filedata.close() - + self.filebluedata=open(othername,'rb') #open also the blue ones self.bluedata=self.filebluedata.read() self.filebluedata.close() - + self.filetype = 'mcs' self.experiment = 'smfluo' - + def is_me(self): if self.filename[-3:].lower()=='mcs': return True else: return False - + def close_all(self): self.filedata.close() self.filebluedata.close() - - + + def default_plots(self): red_data=self.read_file(self.reddata) blue_data=self.read_file(self.bluedata) blue_data=[-1*float(item) for item in blue_data] #visualize blue as "mirror" of red - + main_plot=lhc.PlotObject() main_plot.add_set(range(len(red_data)),red_data) main_plot.add_set(range(len(blue_data)),blue_data) @@ -64,17 +62,17 @@ class mcsDriver(lhc.Driver): main_plot.destination=0 main_plot.title=self.filename main_plot.colors=['red','blue'] - + return [main_plot] - - def read_file(self, raw_data): + + def read_file(self, raw_data): real_data=[] intervalsperfile=struct.unpack('h', raw_data[10:12])[0] #read in number of intervals in this file #this data is contained in bit offset 10-12 in mcs file #see http://docs.python.org/library/struct.html#module-struct for additional explanation - + numbytes=len(raw_data) #data is stored in 4-byte chunks, starting with pos 256 for j in range(0,intervalsperfile): #read in all intervals in file temp=raw_data[256+j*4:256+j*4+4] #data starts at byte offset 256 real_data.append(struct.unpack('i', temp)[0]) #[0] because it returns a 1-element tuple - return real_data \ No newline at end of file + return real_data diff --git a/hooke/driver/mfp1dexport.py b/hooke/driver/mfp1dexport.py index 57e75d9..5cc18e5 100644 --- a/hooke/driver/mfp1dexport.py +++ b/hooke/driver/mfp1dexport.py @@ -1,5 +1,3 @@ -#!/usr/bin/env python - ''' mfp1dexport.py @@ -12,49 +10,49 @@ import libhookecurve as lhc import libhooke as lh class mfp1dexportDriver(lhc.Driver): - + def __init__(self, filename): - + self.filename=filename self.filedata=open(filename,'rU') self.lines=list(self.filedata.readlines()) self.filedata.close() - + self.filetype='mfp1dexport' self.experiment='smfs' - + def close_all(self): self.filedata.close() - + def is_me(self): try: self.raw_header=self.lines[0:38] except: #Not enough lines for a header; not a good file return False - + #FIXME: We want a more reasonable header recognition if self.raw_header[0][0:4]=='Wave': return True else: return False - + def _read_columns(self): - + self.raw_columns=self.lines[39:] - + kline=None for line in self.lines: if line[:7]=='SpringC': kline=line break - + kline=kline.split(':') - + #self.k=float(self.raw_header[23][8:]) self.k=float(kline[1]) - - + + xext=[] xret=[] yext=[] @@ -65,15 +63,15 @@ class mfp1dexportDriver(lhc.Driver): yext.append(float(spline[1])) xret.append(float(spline[2])) yret.append(float(spline[3])) - + return [[xext,yext],[xret,yret]] - + def deflection(self): self.data=self._read_columns() return self.data[0][1],self.data[1][1] - - - def default_plots(self): + + + def default_plots(self): main_plot=lhc.PlotObject() defl_ext,defl_ret=self.deflection() yextforce=[i*self.k for i in defl_ext] diff --git a/hooke/driver/picoforce.py b/hooke/driver/picoforce.py index d331767..e7df889 100755 --- a/hooke/driver/picoforce.py +++ b/hooke/driver/picoforce.py @@ -1,5 +1,3 @@ -#!/usr/bin/env python - ''' libpicoforce.py @@ -20,11 +18,11 @@ __version__='0.0.0.20080404' class DataChunk(list): #Dummy class to provide ext and ret methods to the data list. - + def ext(self): halflen=(len(self)/2) return self[0:halflen] - + def ret(self): halflen=(len(self)/2) return self[halflen:] @@ -32,39 +30,39 @@ class DataChunk(list): class picoforceDriver(lhc.Driver): #Construction and other special methods - + def __init__(self,filename): ''' constructor method ''' - + self.textfile=file(filename) self.binfile=file(filename,'rb') - + #The 0,1,2 data chunks are: #0: D (vs T) #1: Z (vs T) #2: D (vs Z) - - + + self.filepath=filename self.debug=False - + self.filetype='picoforce' self.experiment='smfs' - - + + #Hidden methods. These are meant to be used only by API functions. If needed, however, #they can be called just like API methods. - + def _get_samples_line(self): ''' Gets the samples per line parameters in the file, to understand trigger behaviour. ''' self.textfile.seek(0) - + samps_expr=re.compile(".*Samps") - + samps_values=[] for line in self.textfile.readlines(): if samps_expr.match(line): @@ -73,25 +71,25 @@ class picoforceDriver(lhc.Driver): samps_values.append(samps) except: pass - + #We raise a flag for the fact we meet an offset, otherwise we would take spurious data length arguments. - + return int(samps_values[0]) - + def _get_chunk_coordinates(self): ''' This method gets the coordinates (offset and length) of a data chunk in our Picoforce file. - - It returns a list containing two tuples: - the first element of each tuple is the data_offset, the second is the corresponding + + It returns a list containing two tuples: + the first element of each tuple is the data_offset, the second is the corresponding data size. - - In near future probably each chunk will get its own data structure, with + + In near future probably each chunk will get its own data structure, with offset, size, type, etc. ''' self.textfile.seek(0) - + offset_expr=re.compile(".*Data offset") length_expr=re.compile(".*Data length") @@ -105,42 +103,42 @@ class picoforceDriver(lhc.Driver): offset=int(line.split()[2]) #the third word splitted is the offset (in bytes) data_offsets.append(offset) #We raise a flag for the fact we meet an offset, otherwise we would take spurious data length arguments. - flag_offset=1 - + flag_offset=1 + #same for the data length - if length_expr.match(line) and flag_offset: + if length_expr.match(line) and flag_offset: size=int(line.split()[2]) data_sizes.append(size) #Put down the offset flag until the next offset is met. flag_offset=0 return zip(data_offsets,data_sizes) - + def _get_data_chunk(self,whichchunk): ''' reads a data chunk and converts it in 16bit signed int. ''' offset,size=self._get_chunk_coordinates()[whichchunk] - - + + self.binfile.seek(offset) raw_chunk=self.binfile.read(size) - + my_chunk=[] for data_position in range(0,len(raw_chunk),2): data_unit_bytes=raw_chunk[data_position:data_position+2] #The unpack function converts 2-bytes in a signed int ('h'). #we use output[0] because unpack returns a 1-value tuple, and we want the number only data_unit=struct.unpack('h',data_unit_bytes)[0] - my_chunk.append(data_unit) - + my_chunk.append(data_unit) + return DataChunk(my_chunk) - + def _get_Zscan_info(self,index): ''' gets the Z scan informations needed to interpret the data chunk. These info come from the general section, BEFORE individual chunk headers. - + By itself, the function will parse for three parameters. (index) that tells the function what to return when called by exposed API methods. @@ -149,20 +147,20 @@ class picoforceDriver(lhc.Driver): index=2 : returns Zscan_V_size ''' self.textfile.seek(0) - + ciaoforcelist_expr=re.compile(".*Ciao force") zscanstart_expr=re.compile(".*@Z scan start") zscansize_expr=re.compile(".*@Z scan size") - + ciaoforce_flag=0 theline=0 for line in self.textfile.readlines(): if ciaoforcelist_expr.match(line): ciaoforce_flag=1 #raise a flag: zscanstart and zscansize params to read are later - + if ciaoforce_flag and zscanstart_expr.match(line): raw_Zscanstart_line=line.split() - + if ciaoforce_flag and zscansize_expr.match(line): raw_Zscansize_line=line.split() @@ -171,35 +169,35 @@ class picoforceDriver(lhc.Driver): for itemscanstart,itemscansize in zip(raw_Zscanstart_line,raw_Zscansize_line): Zscanstart_line.append(itemscanstart.strip('[]()')) Zscansize_line.append(itemscansize.strip('[]()')) - + Zscan_V_LSB=float(Zscanstart_line[6]) Zscan_V_start=float(Zscanstart_line[8]) Zscan_V_size=float(Zscansize_line[8]) - + return (Zscan_V_LSB,Zscan_V_start,Zscan_V_size)[index] - + def _get_Z_magnify_scale(self,whichchunk): ''' gets Z scale and Z magnify Here we get Z scale/magnify from the 'whichchunk' only. whichchunk=1,2,3 TODO: make it coherent with data_chunks syntaxis (0,1,2) - + In future, should we divide the *file* itself into chunk descriptions and gain true chunk data structures? ''' self.textfile.seek(0) - + z_scale_expr=re.compile(".*@4:Z scale") z_magnify_expr=re.compile(".*@Z magnify") - + ramp_size_expr=re.compile(".*@4:Ramp size") ramp_offset_expr=re.compile(".*@4:Ramp offset") - + occurrences=0 found_right=0 - - + + for line in self.textfile.readlines(): if z_magnify_expr.match(line): occurrences+=1 @@ -208,141 +206,141 @@ class picoforceDriver(lhc.Driver): raw_z_magnify_expression=line.split() else: found_right=0 - + if found_right and z_scale_expr.match(line): raw_z_scale_expression=line.split() if found_right and ramp_size_expr.match(line): raw_ramp_size_expression=line.split() if found_right and ramp_offset_expr.match(line): raw_ramp_offset_expression=line.split() - + return float(raw_z_magnify_expression[5]),float(raw_z_scale_expression[7]), float(raw_ramp_size_expression[7]), float(raw_ramp_offset_expression[7]), float(raw_z_scale_expression[5][1:]) - - + + #Exposed APIs. - #These are the methods that are meant to be called from external apps. - + #These are the methods that are meant to be called from external apps. + def LSB_to_volt(self,chunknum,voltrange=20): ''' Converts the LSB data of a given chunk (chunknum=0,1,2) in volts. First step to get the deflection and the force. - + SYNTAXIS: item.LSB_to_volt(chunknum, [voltrange]) - + The voltrange is by default set to 20 V. ''' return DataChunk([((float(lsb)/65535)*voltrange) for lsb in self.data_chunks[chunknum]]) - + def LSB_to_deflection(self,chunknum,deflsensitivity=None,voltrange=20): ''' Converts the LSB data in deflection (meters). - + SYNTAXIS: item.LSB_to_deflection(chunknum, [deflection sensitivity], [voltrange]) - + chunknum is the chunk you want to parse (0,1,2) - - The deflection sensitivity by default is the one parsed from the file. + + The deflection sensitivity by default is the one parsed from the file. The voltrange is by default set to 20 V. ''' if deflsensitivity is None: deflsensitivity=self.get_deflection_sensitivity() - - lsbvolt=self.LSB_to_volt(chunknum) + + lsbvolt=self.LSB_to_volt(chunknum) return DataChunk([volt*deflsensitivity for volt in lsbvolt]) - + def deflection(self): ''' Get the actual force curve deflection. ''' deflchunk= self.LSB_to_deflection(2) return deflchunk.ext(),deflchunk.ret() - + def LSB_to_force(self,chunknum=2,Kspring=None,voltrange=20): ''' Converts the LSB data (of deflection) in force (newtons). - + SYNTAXIS: item.LSB_to_force([chunknum], [spring constant], [voltrange]) - + chunknum is the chunk you want to parse (0,1,2). The chunk used is by default 2. The spring constant by default is the one parsed from the file. The voltrange is by default set to 20 V. ''' if Kspring is None: Kspring=self.get_spring_constant() - - lsbdefl=self.LSB_to_deflection(chunknum) + + lsbdefl=self.LSB_to_deflection(chunknum) return DataChunk([(meter*Kspring) for meter in lsbdefl]) - + def get_Zscan_V_start(self): return self._get_Zscan_info(1) - + def get_Zscan_V_size(self): return self._get_Zscan_info(2) - + def get_Z_scan_sensitivity(self): ''' gets Z sensitivity ''' self.textfile.seek(0) - + z_sensitivity_expr=re.compile(".*@Sens. Zsens") - + for line in self.textfile.readlines(): if z_sensitivity_expr.match(line): z_sensitivity=float(line.split()[3]) #return it in SI units (that is: m/V, not nm/V) return z_sensitivity*(10**(-9)) - + def get_Z_magnify(self,whichchunk): ''' Gets the Z magnify factor. Normally it is 1, unknown exact use as of 2006-01-13 ''' return self._get_Z_magnify_scale(whichchunk)[0] - + def get_Z_scale(self,whichchunk): ''' Gets the Z scale. ''' return self._get_Z_magnify_scale(whichchunk)[1] - + def get_ramp_size(self,whichchunk): ''' Gets the -user defined- ramp size ''' return self._get_Z_magnify_scale(whichchunk)[2] - + def get_ramp_offset(self,whichchunk): ''' Gets the ramp offset ''' return self._get_Z_magnify_scale(whichchunk)[3] - + def get_Z_scale_LSB(self,whichchunk): ''' Gets the LSB-to-volt conversion factor of the Z data. (so called hard-scale in the Nanoscope documentation) - + ''' return self._get_Z_magnify_scale(whichchunk)[4] - + def get_deflection_sensitivity(self): ''' gets deflection sensitivity - ''' + ''' self.textfile.seek(0) - + def_sensitivity_expr=re.compile(".*@Sens. DeflSens") - + for line in self.textfile.readlines(): if def_sensitivity_expr.match(line): def_sensitivity=float(line.split()[3]) break #return it in SI units (that is: m/V, not nm/V) return def_sensitivity*(10**(-9)) - + def get_spring_constant(self): ''' gets spring constant. @@ -350,35 +348,35 @@ class picoforceDriver(lhc.Driver): They are normally all equal, but we retain all three for future... ''' self.textfile.seek(0) - + springconstant_expr=re.compile(".*Spring Constant") - + constants=[] - + for line in self.textfile.readlines(): if springconstant_expr.match(line): constants.append(float(line.split()[2])) - + return constants[0] - + def get_Zsensorsens(self): ''' gets Zsensorsens for Z data. - + This is the sensitivity needed to convert the LSB data in nanometers for the Z-vs-T data chunk. - ''' + ''' self.textfile.seek(0) - + zsensorsens_expr=re.compile(".*Sens. ZSensorSens") - + for line in self.textfile.readlines(): if zsensorsens_expr.match(line): zsensorsens_raw_expression=line.split() #we must take only first occurrence, so we exit from the cycle immediately break - + return (float(zsensorsens_raw_expression[3]))*(10**(-9)) - + def Z_data(self): ''' returns converted ext and ret Z curves. @@ -390,18 +388,18 @@ class picoforceDriver(lhc.Driver): #rampsize_zt=self.get_ramp_size(2) #rampoffset_zt=self.get_ramp_offset(2) zsensorsens=self.get_Zsensorsens() - + ''' The magic formula that converts the Z data is: - + meters = LSB * V_lsb_conversion_factor * ZSensorSens ''' - + #z_curves=[item*Zlsb_zt*zsensorsens for item in self.data_chunks[1].pair['ext']],[item*Zlsb_zt*zsensorsens for item in self.data_chunks[1].pair['ret']] - z_curves=[item*Zlsb_zt*zsensorsens for item in self.data_chunks[1].ext()],[item*Zlsb_zt*zsensorsens for item in self.data_chunks[1].ret()] + z_curves=[item*Zlsb_zt*zsensorsens for item in self.data_chunks[1].ext()],[item*Zlsb_zt*zsensorsens for item in self.data_chunks[1].ret()] z_curves=[DataChunk(item) for item in z_curves] return z_curves - + def Z_extremes(self): ''' returns the extremes of the Z values @@ -410,9 +408,9 @@ class picoforceDriver(lhc.Driver): z_extremes={} z_extremes['ext']=zcurves[0][0],zcurves[0][-1] z_extremes['ret']=zcurves[1][0],zcurves[1][-1] - + return z_extremes - + def Z_step(self): ''' returns the calculated step between the Z values @@ -421,10 +419,10 @@ class picoforceDriver(lhc.Driver): zpoints={} z_extremes=self.Z_extremes() - + zrange['ext']=abs(z_extremes['ext'][0]-z_extremes['ext'][1]) zrange['ret']=abs(z_extremes['ret'][0]-z_extremes['ret'][1]) - + #We must take 1 from the calculated zpoints, or when I use the arange function gives me a point more #with the step. That is, if I have 1000 points, and I use arange(start,stop,step), I have 1001 points... #For cleanness, solution should really be when using arange, but oh well... @@ -432,17 +430,17 @@ class picoforceDriver(lhc.Driver): zpoints['ret']=len(self.Z_data()[1])-1 #this syntax must become coherent!! return (zrange['ext']/zpoints['ext']),(zrange['ret']/zpoints['ret']) - + def Z_domains(self): ''' returns the Z domains on which to plot the force data. - + The Z domains are returned as a single long DataChunk() extended list. The extension and retraction part - can be extracted using ext() and ret() methods. - ''' + can be extracted using ext() and ret() methods. + ''' x1step=self.Z_step()[0] - x2step=self.Z_step()[1] - + x2step=self.Z_step()[1] + try: xext=arange(self.Z_extremes()['ext'][0],self.Z_extremes()['ext'][1],-x1step) xret=arange(self.Z_extremes()['ret'][0],self.Z_extremes()['ret'][1],-x2step) @@ -450,7 +448,7 @@ class picoforceDriver(lhc.Driver): xext=arange(0,1) xret=arange(0,1) print 'picoforce.py: Warning. xext, xret domains cannot be extracted.' - + if not (len(xext)==len(xret)): if self.debug: #print warning @@ -460,36 +458,36 @@ class picoforceDriver(lhc.Driver): print "You cannot trust the resulting curve." print "Until a solution is found, I substitute the ext domain with the ret domain. Sorry." xext=xret - + return DataChunk(xext.tolist()+xret.tolist()) - + def Z_scan_size(self): return self.get_Zscan_V_size()*self.get_Z_scan_sensitivity() - + def Z_start(self): return self.get_Zscan_V_start()*self.get_Z_scan_sensitivity() - + def ramp_size(self,whichchunk): ''' to be implemented if needed ''' raise "Not implemented yet." - - + + def ramp_offset(self,whichchunk): ''' to be implemented if needed ''' raise "Not implemented yet." - + def detriggerize(self, forcext): ''' Cuts away the trigger-induced s**t on the extension curve. DEPRECATED cutindex=2 startvalue=forcext[0] - - for index in range(len(forcext)-1,2,-2): + + for index in range(len(forcext)-1,2,-2): if forcext[index]>startvalue: cutindex=index else: @@ -498,7 +496,7 @@ class picoforceDriver(lhc.Driver): return cutindex ''' return 0 - + def is_me(self): ''' self-identification of file type magic @@ -506,40 +504,40 @@ class picoforceDriver(lhc.Driver): curve_file=file(self.filepath) header=curve_file.read(30) curve_file.close() - + if header[2:17] == 'Force file list': #header of a picoforce file self.data_chunks=[self._get_data_chunk(num) for num in [0,1,2]] return True else: return False - + def close_all(self): ''' Explicitly closes all files ''' self.textfile.close() self.binfile.close() - + def default_plots(self): ''' creates the default PlotObject ''' - - + + force=self.LSB_to_force() zdomain=self.Z_domains() - + samples=self._get_samples_line() #cutindex=0 #cutindex=self.detriggerize(force.ext()) - + main_plot=lhc.PlotObject() - + main_plot.vectors=[[zdomain.ext()[0:samples], force.ext()[0:samples]],[zdomain.ret()[0:samples], force.ret()[0:samples]]] main_plot.normalize_vectors() main_plot.units=['meters','newton'] main_plot.destination=0 main_plot.title=self.filepath - - + + return [main_plot] diff --git a/hooke/driver/picoforcealt.py b/hooke/driver/picoforcealt.py index 078430b..ec69118 100644 --- a/hooke/driver/picoforcealt.py +++ b/hooke/driver/picoforcealt.py @@ -1,5 +1,3 @@ -#!/usr/bin/env python - ''' libpicoforce.py @@ -22,11 +20,11 @@ __version__='0.0.0.20081706' class DataChunk(list): #Dummy class to provide ext and ret methods to the data list. - + def ext(self): halflen=(len(self)/2) return self[0:halflen] - + def ret(self): halflen=(len(self)/2) return self[halflen:] @@ -34,40 +32,40 @@ class DataChunk(list): class picoforcealtDriver(lhc.Driver): #Construction and other special methods - + def __init__(self,filename): ''' constructor method ''' - + self.textfile=file(filename) self.binfile=file(filename,'rb') - + #The 0,1,2 data chunks are: #0: D (vs T) #1: Z (vs T) #2: D (vs Z) - + self.forcechunk=0 self.distancechunk=1 #TODO eliminate the need to set chunk numbers - + self.filepath=filename self.debug=True - + self.filetype='picoforce' self.experiment='smfs' - - - + + + def _get_samples_line(self): ''' Gets the samples per line parameters in the file, to understand trigger behaviour. ''' self.textfile.seek(0) - + samps_expr=re.compile(".*Samps") - + samps_values=[] for line in self.textfile.readlines(): if samps_expr.match(line): @@ -76,25 +74,25 @@ class picoforcealtDriver(lhc.Driver): samps_values.append(samps) except: pass - + #We raise a flag for the fact we meet an offset, otherwise we would take spurious data length arguments. - + return int(samps_values[0]) - + def _get_chunk_coordinates(self): ''' This method gets the coordinates (offset and length) of a data chunk in our Picoforce file. - - It returns a list containing two tuples: - the first element of each tuple is the data_offset, the second is the corresponding + + It returns a list containing two tuples: + the first element of each tuple is the data_offset, the second is the corresponding data size. - - In near future probably each chunk will get its own data structure, with + + In near future probably each chunk will get its own data structure, with offset, size, type, etc. ''' self.textfile.seek(0) - + offset_expr=re.compile(".*Data offset") length_expr=re.compile(".*Data length") @@ -108,35 +106,35 @@ class picoforcealtDriver(lhc.Driver): offset=int(line.split()[2]) #the third word splitted is the offset (in bytes) data_offsets.append(offset) #We raise a flag for the fact we meet an offset, otherwise we would take spurious data length arguments. - flag_offset=1 - + flag_offset=1 + #same for the data length - if length_expr.match(line) and flag_offset: + if length_expr.match(line) and flag_offset: size=int(line.split()[2]) data_sizes.append(size) #Put down the offset flag until the next offset is met. flag_offset=0 return zip(data_offsets,data_sizes) - + def _get_data_chunk(self,whichchunk): ''' reads a data chunk and converts it in 16bit signed int. ''' offset,size=self._get_chunk_coordinates()[whichchunk] - - + + self.binfile.seek(offset) raw_chunk=self.binfile.read(size) - + my_chunk=[] for data_position in range(0,len(raw_chunk),2): data_unit_bytes=raw_chunk[data_position:data_position+2] #The unpack function converts 2-bytes in a signed int ('h'). #we use output[0] because unpack returns a 1-value tuple, and we want the number only data_unit=struct.unpack('h',data_unit_bytes)[0] - my_chunk.append(data_unit) - + my_chunk.append(data_unit) + return DataChunk(my_chunk) def _force(self): @@ -149,72 +147,72 @@ class picoforcealtDriver(lhc.Driver): voltrange=1 z_scale=self._get_Z_scale() deflsensitivity=self.get_deflection_sensitivity() - volts=[((float(lsb))*voltrange*z_scale) for lsb in self.data_chunks[self.forcechunk]] + volts=[((float(lsb))*voltrange*z_scale) for lsb in self.data_chunks[self.forcechunk]] deflect=[volt*deflsensitivity for volt in volts] - + return deflect - - + + def _Z(self): - #returns distance vector (calculated instead than from data chunk) + #returns distance vector (calculated instead than from data chunk) rampsize=self._get_rampsize() sampsline=self._get_samples_line() senszscan=self._get_Z_scan_sens() - + xstep=senszscan*rampsize/sampsline*10**(-9) xext=arange(sampsline*xstep,0,-xstep) xret=arange(sampsline*xstep,0,-xstep) - + return DataChunk(xext.tolist()+xret.tolist()) - + def _get_Z_scale(self): self.textfile.seek(0) expr=re.compile(".*@4:Z scale") - + for line in self.textfile.readlines(): if expr.match(line): zscale=float((line.split()[5]).strip("() []")) break return zscale - + def _get_rampsize(self): self.textfile.seek(0) expr=re.compile(".*@4:Ramp size:") - + for line in self.textfile.readlines(): if expr.match(line): zsens=float((line.split()[7]).strip("() []")) break return zsens - + def _get_Z_scan_sens(self): self.textfile.seek(0) expr=re.compile(".*@Sens. Zsens") - + for line in self.textfile.readlines(): if expr.match(line): zsens=float((line.split()[3]).strip("() []")) break return zsens - - - + + + def get_deflection_sensitivity(self): ''' gets deflection sensitivity - ''' + ''' self.textfile.seek(0) - + def_sensitivity_expr=re.compile(".*@Sens. DeflSens") - + for line in self.textfile.readlines(): if def_sensitivity_expr.match(line): def_sensitivity=float(line.split()[3]) break #return it in SI units (that is: m/V, not nm/V) return def_sensitivity*(10**(-9)) - + def get_spring_constant(self): ''' gets spring constant. @@ -222,17 +220,17 @@ class picoforcealtDriver(lhc.Driver): They are normally all equal, but we retain all three for future... ''' self.textfile.seek(0) - + springconstant_expr=re.compile(".*Spring Constant") - + constants=[] - + for line in self.textfile.readlines(): if springconstant_expr.match(line): constants.append(float(line.split()[2])) - + return constants[0] - + def is_me(self): ''' self-identification of file type magic @@ -240,21 +238,21 @@ class picoforcealtDriver(lhc.Driver): curve_file=file(self.filepath) header=curve_file.read(30) curve_file.close() - + if header[2:17] == 'Force file list': #header of a picoforce file #here DONT translate chunk self.data_chunks=[self._get_data_chunk(num) for num in [0,1,2]] return True else: return False - + def close_all(self): ''' Explicitly closes all files ''' self.textfile.close() self.binfile.close() - + def default_plots(self): ''' creates the default PlotObject @@ -268,8 +266,8 @@ class picoforcealtDriver(lhc.Driver): main_plot.units=['meters','newton'] main_plot.destination=0 main_plot.title=self.filepath - - + + return [main_plot] def deflection(self): diff --git a/hooke/driver/tutorialdriver.py b/hooke/driver/tutorialdriver.py index b16512a..9d38c47 100644 --- a/hooke/driver/tutorialdriver.py +++ b/hooke/driver/tutorialdriver.py @@ -1,5 +1,3 @@ -#!/usr/bin/env python - ''' tutorialdriver.py @@ -53,7 +51,7 @@ class tutorialdriverDriver(lhc.Driver): ''' This is a *generic* driver, not a specific force spectroscopy driver. See the written documentation to see what a force spectroscopy driver must be defined to take into account Hooke facilities. - + Our driver is a class with the name convention nameofthedriverDriver, where "nameofthedriver" is the filename. The driver must inherit from the parent class lhc.Driver, so the syntax is class nameofthedriverDriver(lhc.Driver) @@ -68,40 +66,40 @@ class tutorialdriverDriver(lhc.Driver): ''' In this case, we have a data format that is just a list of ASCII values, so we can just divide that in rows, and generate a list with each item being a row. - Of course if your data files are binary, or follow a different approach, do whatever you need. :) + Of course if your data files are binary, or follow a different approach, do whatever you need. :) ''' self.data = list(self.filedata) self.filedata.close() #remember to close the file - + '''These are two strings that can be used by Hooke commands/plugins to understand what they are looking at. They have no other meaning. They have to be somehow defined however - commands often look for those variables. - + self.filetype should contain the name of the exact filetype defined by the driver (so that filetype-specific commands can know if they're dealing with the correct filetype) self.experiment should contain instead the type of data involved (for example, various drivers can be used for force-clamp experiments, - but hooke commands could like to know if we're looking at force clamp data, regardless of their origin, and not other + but hooke commands could like to know if we're looking at force clamp data, regardless of their origin, and not other kinds of data) - + Of course, all other variables you like can be defined in the class. ''' self.filetype = 'tutorial' self.experiment = 'generic' - + def is_me(self): ''' THIS METHOD MUST BE DEFINED. RETURNS: Boolean (True or False) - This method must be an heuristic that looks at the file content and decides if the file can be opened by the driver itself. - It returns True if the file opened can be interpreted by the current driver, False otherwise. + This method must be an heuristic that looks at the file content and decides if the file can be opened by the driver itself. + It returns True if the file opened can be interpreted by the current driver, False otherwise. Defining this method allows Hooke to understand what kind of files we're looking at automatically. - + We have to re-open/re-close the file here. ''' - + myfile=open(self.filename, 'r') headerline=myfile.readlines()[0] #we take the first line myfile.close() - + ''' Here, our "magic fingerprint" is the TUTORIAL_FILE header. Of course, depending on the data file, you can have interesting headers, or patterns, etc. that you can use to guess the data format. What matters is successful recognizing, and returning @@ -111,12 +109,12 @@ class tutorialdriverDriver(lhc.Driver): return True else: return False - + def _generate_vectors(self): ''' Here we parse the data and generate the raw vectors. This method has only to do with the peculiar file format here, so it's of - no big interest (I just wrote it to present a functional driver). - + no big interest (I just wrote it to present a functional driver). + Only thing to remember, it can be nice to call methods that are used only "internally" by the driver (or by plugins) with a "_" prefix, so to have a visual remark. But it's just an advice. ''' @@ -134,10 +132,10 @@ class tutorialdriverDriver(lhc.Driver): elif item[:-1]=='PLOT2': whatplot=item[:-1] elif item[0]=='X' or item[0]=='Y': - pos=positions[item[:-1]] + pos=positions[item[:-1]] else: pass - + return vectors def close_all(self): @@ -147,62 +145,62 @@ class tutorialdriverDriver(lhc.Driver): In this method, all file objects defined in the driver must be closed. ''' self.filename.close() - - + + def default_plots(self): ''' THIS METHOD MUST BE DEFINED. RETURNS: [ lhc.PlotObject ] or [ lhc.PlotObject, lhc.PlotObject] - + This is the method that returns the plots to Hooke. It must return a list with *one* or *two* PlotObjects. - + See the libhookecurve.py source code to see how PlotObjects are defined and work in detail. ''' gen_vectors=self._generate_vectors() - + #Here we create the main plot PlotObject and initialize its vectors. main_plot=lhc.PlotObject() main_plot.vectors=[] #The same for the other plot. other_plot=lhc.PlotObject() other_plot.vectors=[] - + ''' Now we fill the plot vectors with our data. set 1 set 2 The "correct" shape of the vector is [ [[x1,x2,x3...],[y1,y2,y3...]] , [[x1,x2,x3...],[y1,y2,y3...]] ], so we have to put stuff in this way into it. - + The add_set() method takes care of this , just use plot.add_set(x,y). ''' main_plot.add_set(gen_vectors['PLOT1'][0],gen_vectors['PLOT1'][1]) main_plot.add_set(gen_vectors['PLOT1'][2],gen_vectors['PLOT1'][3]) - + other_plot.add_set(gen_vectors['PLOT2'][0],gen_vectors['PLOT2'][1]) other_plot.add_set(gen_vectors['PLOT2'][2],gen_vectors['PLOT2'][3]) - + ''' normalize_vectors() trims the vectors, so that if two x/y couples are of different lengths, the latest points are trimmed (otherwise we have a python error). Always a good idea to run it, to avoid crashes. ''' main_plot.normalize_vectors() other_plot.normalize_vectors() - + ''' Here we define: - units: [string, string], define the measure units of X and Y axes - destination: 0/1 , defines where to plot the plot (0=top, 1=bottom), default=0 - title: string , the plot title. - + for each plot. Again, see libhookecurve.py comments for details. ''' main_plot.units=['unit of x','unit of y'] main_plot.destination=0 main_plot.title=self.filename+' main' - + other_plot.units=['unit of x','unit of y'] other_plot.destination=1 other_plot.title=self.filename+' other' - - return [main_plot, other_plot] \ No newline at end of file + + return [main_plot, other_plot] diff --git a/hooke/hooke.py b/hooke/hooke.py index 6046847..55f23b1 100755 --- a/hooke/hooke.py +++ b/hooke/hooke.py @@ -811,4 +811,5 @@ def main(): app.MainLoop() -main() +if __name__ == '__main__': + main() diff --git a/hooke/hooke_cli.py b/hooke/hooke_cli.py index d3c9693..fc0cee3 100755 --- a/hooke/hooke_cli.py +++ b/hooke/hooke_cli.py @@ -10,7 +10,6 @@ Copyright (C) 2006 Massimo Sandal (University of Bologna, Italy). This program is released under the GNU General Public License version 2. ''' - from libhooke import * #FIXME import libhookecurve as lhc @@ -50,39 +49,39 @@ import platform class HookeCli(cmd.Cmd): - + def __init__(self,frame,list_of_events,events_from_gui,config,drivers): cmd.Cmd.__init__(self) - + self.prompt = 'hooke: ' - - + + self.current_list=[] #the playlist we're using - - self.current=None #the current curve under analysis. + + self.current=None #the current curve under analysis. self.plots=None ''' The actual hierarchy of the "current curve" is a bit complex: - + self.current = the lhc.HookeCurve container object of the current curve self.current.curve = the current "real" curve object as defined in the filetype driver class self.current.curve.default_plots() = the default plots of the filetype driver. - - The plot objects obtained by mean of self.current.curve.default_plots() + + The plot objects obtained by mean of self.current.curve.default_plots() then undergoes modifications by the plotmanip - modifier functions. The modified plot is saved in self.plots and used if needed by other functions. + modifier functions. The modified plot is saved in self.plots and used if needed by other functions. ''' - - + + self.pointer=0 #a pointer to navigate the current list - + #Things that come from outside self.frame=frame #the wx frame we refer to self.list_of_events=list_of_events #a list of wx events we use to interact with the GUI self.events_from_gui=events_from_gui #the Queue object we use to have messages from the GUI self.config=config #the configuration dictionary self.drivers=drivers #the file format drivers - + #get plot manipulation functions plotmanip_functions=[] for object_name in dir(self): @@ -97,8 +96,8 @@ class HookeCli(cmd.Cmd): self.plotmanip[nameindex] = item else: pass - - + + self.playlist_saved=0 #Did we save the playlist? self.playlist_name='' #Name of playlist self.notes_saved=1 #Did we save the notes? @@ -106,10 +105,10 @@ class HookeCli(cmd.Cmd): #create outlet self.outlet=lout.Outlet() - + #Data that must be saved in the playlist, related to the whole playlist (not individual curves) - self.playlist_generics={} - + self.playlist_generics={} + #make sure we execute _plug_init() for every command line plugin we import for plugin_name in self.config['plugins']: try: @@ -123,7 +122,7 @@ class HookeCli(cmd.Cmd): #load default list, if possible self.do_loadlist(self.config['defaultlist']) - + #HELPER FUNCTIONS #Everything sending an event should be here def _measure_N_points(self, N, whatset=1): @@ -138,7 +137,7 @@ class HookeCli(cmd.Cmd): except Empty: pass return points - + def _get_displayed_plot(self,dest=0): ''' returns the currently displayed plot. @@ -152,32 +151,32 @@ class HookeCli(cmd.Cmd): if displayed_plot: break return displayed_plot - + def _send_plot(self,plots): ''' sends a plot to the GUI ''' wx.PostEvent(self.frame, self.list_of_events['plot_graph'](plots=plots)) return - + def _find_plotmanip(self, name): ''' returns a plot manipulator function from its name ''' return self.plotmanip[self.config['plotmanips'].index(name)] - + def _clickize(self, xvector, yvector, index): ''' - returns a ClickedPoint() object from an index and vectors of x, y coordinates + returns a ClickedPoint() object from an index and vectors of x, y coordinates ''' point=ClickedPoint() point.index=index point.absolute_coords=xvector[index],yvector[index] point.find_graph_coords(xvector,yvector) return point - + #HERE COMMANDS BEGIN - + def help_set(self): print ''' SET @@ -209,13 +208,13 @@ Syntax: set [variable] [value] value=None else: value=args[1] - + self.config[key]=value self.do_plot(0) - + #PLAYLIST MANAGEMENT AND NAVIGATION #------------------------------------ - + def help_loadlist(self): print ''' LOADLIST @@ -227,41 +226,41 @@ Syntax: loadlist [playlist file] #checking for args: if nothing is given as input, we warn and exit. while len(args)==0: args=linp.safeinput('File to load?') - + arglist=args.split() play_to_load=arglist[0] - + #We assume a Hooke playlist has the extension .hkp if play_to_load[-4:] != '.hkp': play_to_load+='.hkp' - - try: + + try: playxml=PlaylistXML() self.current_list, self.playlist_generics=playxml.load(play_to_load) self.current_playxml=playxml except IOError: print 'File not found.' return - + print 'Loaded %s curves' %len(self.current_list) - + if 'pointer' in self.playlist_generics.keys(): self.pointer=int(self.playlist_generics['pointer']) else: #if no pointer is found, set the current curve as the first curve of the loaded playlist self.pointer=0 print 'Starting at curve ',self.pointer - + self.current=self.current_list[self.pointer] - + #resets saved/notes saved state self.playlist_saved=0 self.playlist_name='' - self.notes_saved=0 - + self.notes_saved=0 + self.do_plot(0) - - + + def help_genlist(self): print ''' GENLIST @@ -277,29 +276,29 @@ genlist dir/*.* are all equivalent syntax. ------------ Syntax: genlist [input files] - + ''' def do_genlist(self,args): #args list is: input path, output name if len(args)==0: args=linp.safeinput('Input files?') - - arglist=args.split() + + arglist=args.split() list_path=arglist[0] - + #if it's a directory, is like /directory/*.* #FIXME: probably a bit kludgy. - if os.path.isdir(list_path): + if os.path.isdir(list_path): if platform.system == 'Windows': SLASH="\\" else: SLASH="/" if list_path[-1] == SLASH: list_path=list_path+'*.*' - else: + else: list_path=list_path+SLASH+'*.*' - - #expanding correctly the input list with the glob module :) + + #expanding correctly the input list with the glob module :) list_files=glob.glob(list_path) list_files.sort() @@ -307,25 +306,25 @@ Syntax: genlist [input files] for item in list_files: try: if os.path.isfile(item): - self.current_list.append(lhc.HookeCurve(os.path.abspath(item))) + self.current_list.append(lhc.HookeCurve(os.path.abspath(item))) except: pass - - self.pointer=0 + + self.pointer=0 if len(self.current_list)>0: self.current=self.current_list[self.pointer] else: print 'Empty list!' return - + #resets saved/notes saved state self.playlist_saved=0 self.playlist_name='' - self.notes_saved=0 - + self.notes_saved=0 + self.do_plot(0) - - + + def do_savelist(self,args): ''' SAVELIST @@ -335,22 +334,22 @@ Syntax: genlist [input files] ''' while len(args)==0: args=linp.safeinput('Output file?',['savedlist.txt']) - + output_filename=args - + self.playlist_generics['pointer']=self.pointer - + #autocomplete filename if not specified if output_filename[-4:] != '.hkp': output_filename+='.hkp' - + playxml=PlaylistXML() playxml.export(self.current_list, self.playlist_generics) - playxml.save(output_filename) - + playxml.save(output_filename) + #remembers we have saved playlist self.playlist_saved=1 - + def help_addtolist(self): print ''' ADDTOLIST @@ -364,14 +363,14 @@ Syntax: addtolist [filename] print 'You must give the input filename you want to add' self.help_addtolist() return - + filenames=glob.glob(args) - + for filename in filenames: self.current_list.append(lhc.HookeCurve(os.path.abspath(filename))) #we need to save playlist self.playlist_saved=0 - + def help_printlist(self): print ''' PRINTLIST @@ -382,8 +381,8 @@ Syntax: printlist def do_printlist(self,args): for item in self.current_list: print item.path - - + + def help_jump(self): print ''' JUMP @@ -392,31 +391,31 @@ Jumps to a given curve. Syntax: jump {$curve} If the curve is not in the current playlist, it politely asks if we want to add it. - ''' + ''' def do_jump(self,filename): ''' jumps to the curve with the given filename. if the filename is not in the playlist, it asks if we must add it or not. ''' - + if filename=='': filename=linp.safeinput('Jump to?') - + filepath=os.path.abspath(filename) print filepath - + c=0 item_not_found=1 while item_not_found: try: - + if self.current_list[c].path == filepath: self.pointer=c self.current=self.current_list[self.pointer] item_not_found=0 self.do_plot(0) else: - c+=1 + c+=1 except IndexError: #We've found the end of the list. answer=linp.safeinput('Curve not found in playlist. Add it to list?',['y']) @@ -429,10 +428,10 @@ If the curve is not in the current playlist, it politely asks if we want to add self.current=self.current_list[-1] self.pointer=(len(current_list)-1) self.do_plot(0) - + item_not_found=0 - - + + def do_index(self,args): ''' INDEX @@ -440,9 +439,9 @@ If the curve is not in the current playlist, it politely asks if we want to add ----- Syntax: index ''' - print self.pointer+1, 'of', len(self.current_list) - - + print self.pointer+1, 'of', len(self.current_list) + + def help_next(self): print ''' NEXT @@ -458,22 +457,22 @@ Syntax: next, n print 'No curve file loaded, currently!' print 'This should not happen, report to http://code.google.com/p/hooke' return - + if self.pointer == (len(self.current_list)-1): self.pointer=0 print 'Playlist finished; back to first curve.' else: self.pointer+=1 - + self.current=self.current_list[self.pointer] self.do_plot(0) - - + + def help_n(self): self.help_next() def do_n(self,args): self.do_next(args) - + def help_previous(self,args): print ''' PREVIOUS @@ -491,22 +490,22 @@ Syntax: previous, p return if self.pointer == 0: self.pointer=(len(self.current_list)-1) - print 'Start of playlist; jump to last curve.' + print 'Start of playlist; jump to last curve.' else: self.pointer-=1 - + self.current=self.current_list[self.pointer] self.do_plot(args) - - + + def help_p(self): self.help_previous() def do_p(self,args): self.do_previous(args) - + #PLOT INTERACTION COMMANDS -#------------------------------- +#------------------------------- def help_plot(self): print ''' PLOT @@ -515,7 +514,7 @@ Plots the current force curve Syntax: plot ''' def do_plot(self,args): - + self.current.identify(self.drivers) self.plots=self.current.curve.default_plots() try: @@ -524,21 +523,21 @@ Syntax: plot print 'Unexpected error occurred in do_plot().' print e return - + #apply the plotmanip functions eventually present nplots=len(self.plots) c=0 while c 1: - dest=int(args[1]) - + dest=int(args[1]) + export_image=self.list_of_events['export_image'] wx.PostEvent(self.frame, export_image(name=name, dest=dest)) - - + + def help_txt(self): print ''' TXT @@ -660,7 +659,7 @@ X1 , Y1 , X2 , Y2 , X3 , Y3 ... Syntax: txt [filename] {plot to export} ''' def do_txt(self,args): - + def transposed2(lists, defval=0): ''' transposes a list of lists, i.e. from [[a,b,c],[x,y,z]] to [[a,x],[b,y],[c,z]] without losing @@ -669,7 +668,7 @@ Syntax: txt [filename] {plot to export} ''' if not lists: return [] return map(lambda *row: [elem or defval for elem in row], *lists) - + whichplot=0 args=args.split() if len(args)==0: @@ -680,39 +679,39 @@ Syntax: txt [filename] {plot to export} whichplot=int(args[1]) except: pass - - columns=[] + + columns=[] for dataset in self.plots[whichplot].vectors: - for i in range(0,len(dataset)): + for i in range(0,len(dataset)): columns.append([]) for value in dataset[i]: - columns[-1].append(str(value)) - + columns[-1].append(str(value)) + rows=transposed2(columns, 'nan') rows=[' , '.join(item) for item in rows] text='\n'.join(rows) - + txtfile=open(filename,'w+') #Save units of measure in header txtfile.write('X:'+self.plots[whichplot].units[0]+'\n') txtfile.write('Y:'+self.plots[whichplot].units[1]+'\n') txtfile.write(text) txtfile.close() - - + + #LOGGING, REPORTING, NOTETAKING - + def do_note_old(self,args): ''' NOTE_OLD **deprecated**: Use note instead. Will be removed in 0.9 - + Writes or displays a note about the current curve. If [anything] is empty, it displays the note, otherwise it adds a note. The note is then saved in the playlist if you issue a savelist command --------------- - Syntax: note_old [anything] + Syntax: note_old [anything] ''' if args=='': @@ -725,20 +724,20 @@ Syntax: txt [filename] {plot to export} args=args.decode('ascii','ignore') if len(args)==0: args='?' - + self.current_list[self.pointer].notes=args self.notes_saved=0 - - + + def do_note(self,args): ''' NOTE - + Writes or displays a note about the current curve. If [anything] is empty, it displays the note, otherwise it adds a note. The note is then saved in the playlist if you issue a savelist command. --------------- - Syntax: note_old [anything] + Syntax: note_old [anything] ''' if args=='': @@ -750,8 +749,8 @@ Syntax: txt [filename] {plot to export} f=open(self.notes_filename,'w') f.write(title_line) f.close() - - #bypass UnicodeDecodeError troubles + + #bypass UnicodeDecodeError troubles try: args=args.decode('ascii') except: @@ -759,25 +758,25 @@ Syntax: txt [filename] {plot to export} if len(args)==0: args='?' self.current_list[self.pointer].notes=args - + f=open(self.notes_filename,'a+') note_string=(self.current.path+' | '+self.current.notes+'\n') f.write(note_string) f.close() - + def help_notelog(self): print ''' NOTELOG Writes a log of the notes taken during the session for the current playlist --------------- +-------------- Syntax notelog [filename] -''' +''' def do_notelog(self,args): - + if len(args)==0: args=linp.safeinput('Notelog filename?',['notelog.txt']) - + note_lines='Notes taken at '+time.asctime()+'\n' for item in self.current_list: if len(item.notes)>0: @@ -785,7 +784,7 @@ Syntax notelog [filename] #FIXME: file path should be truncated... note_string=(item.path+' | '+item.notes+'\n') note_lines+=note_string - + try: f=open(args,'a+') f.write(note_lines) @@ -793,7 +792,7 @@ Syntax notelog [filename] except IOError, (ErrorNumber, ErrorMessage): print 'Error: notes cannot be saved. Catched exception:' print ErrorMessage - + self.notes_saved=1 def help_copylog(self): @@ -804,15 +803,15 @@ Moves the annotated curves to another directory Syntax copylog [directory] ''' def do_copylog(self,args): - + if len(args)==0: args=linp.safeinput('Destination directory?') #TODO default - + mydir=os.path.abspath(args) if not os.path.isdir(mydir): print 'Destination is not a directory.' return - + for item in self.current_list: if len(item.notes)>0: try: @@ -848,7 +847,7 @@ Syntax copylog [directory] self.outlet.delete(args) #OS INTERACTION COMMANDS -#----------------- +#----------------- def help_dir(self): print ''' DIR, LS @@ -858,16 +857,16 @@ Syntax: dir [path] ls [path] ''' def do_dir(self,args): - + if len(args)==0: args='*' print glob.glob(args) - + def help_ls(self): self.help_dir(self) def do_ls(self,args): self.do_dir(args) - + def help_pwd(self): print ''' PWD @@ -876,8 +875,8 @@ Gives the current working directory. Syntax: pwd ''' def do_pwd(self,args): - print os.getcwd() - + print os.getcwd() + def help_cd(self): print ''' CD @@ -891,8 +890,8 @@ Syntax: cd os.chdir(mypath) except OSError: print 'I cannot access that directory.' - - + + def help_system(self): print ''' SYSTEM @@ -902,15 +901,15 @@ Syntax system [command line] ''' pass def do_system(self,args): - waste=os.system(args) - + waste=os.system(args) + def do_debug(self,args): ''' this is a dummy command where I put debugging things ''' print self.config['plotmanips'] pass - + def help_current(self): print ''' CURRENT @@ -920,7 +919,7 @@ Syntax: current ''' def do_current(self,args): print self.current.path - + def do_info(self,args): ''' INFO @@ -934,13 +933,13 @@ Syntax: current for set in plot.vectors: lengths=[len(item) for item in set] print 'Data set size: ',lengths - + def do_version(self,args): ''' VERSION ------ Prints the current version and codename, plus library version. Useful for debugging. - ''' + ''' print 'Hooke '+__version__+' ('+__codename__+')' print 'Released on: '+__releasedate__ print '---' @@ -954,7 +953,7 @@ Syntax: current print 'Platform: '+str(platform.uname()) print '---' print 'Loaded plugins:',self.config['loaded_plugins'] - + def help_exit(self): print ''' EXIT, QUIT @@ -962,30 +961,27 @@ Exits the program cleanly. ------ Syntax: exit Syntax: quit -''' +''' def do_exit(self,args): we_exit='N' - + if (not self.playlist_saved) or (not self.notes_saved): we_exit=linp.safeinput('You did not save your playlist and/or notes. Exit?',['n']) else: we_exit=linp.safeinput('Exit?',['y']) - + if we_exit[0].upper()=='Y': wx.CallAfter(self.frame.Close) sys.exit(0) else: return - + def help_quit(self): self.help_exit() def do_quit(self,args): self.do_exit(args) - - - if __name__ == '__main__': mycli=HookeCli(0) mycli.cmdloop() diff --git a/hooke/libhooke.py b/hooke/libhooke.py index f668635..69101c8 100755 --- a/hooke/libhooke.py +++ b/hooke/libhooke.py @@ -1,5 +1,4 @@ #!/usr/bin/env python - ''' libhooke.py @@ -11,8 +10,6 @@ With algorithms contributed by Francesco Musiani (University of Bologna, Italy) This program is released under the GNU General Public License version 2. ''' - - import libhookecurve as lhc import scipy @@ -26,20 +23,20 @@ import string import csv HOOKE_VERSION=['0.8.3_devel', 'Seinei', '2008-04-16'] -WX_GOOD=['2.6','2.8'] - +WX_GOOD=['2.6','2.8'] + class PlaylistXML: ''' This module allows for import/export of an XML playlist into/out of a list of HookeCurve objects ''' - + def __init__(self): - + self.playlist=None #the DOM object representing the playlist data structure self.playpath=None #the path of the playlist XML file self.plaything=None self.hidden_attributes=['curve'] #This list contains hidden attributes that we don't want to go into the playlist. - + def export(self, list_of_hooke_curves, generics): ''' Creates an initial playlist from a list of files. @@ -48,21 +45,21 @@ class PlaylistXML: - ''' - + ''' + #create the output playlist, a simple XML document impl=xml.dom.minidom.getDOMImplementation() #create the document DOM object and the root element newdoc=impl.createDocument(None, "playlist",None) top_element=newdoc.documentElement - + #save generics variables playlist_generics=newdoc.createElement("generics") top_element.appendChild(playlist_generics) for key in generics.keys(): newdoc.createAttribute(key) playlist_generics.setAttribute(key,str(generics[key])) - + #save curves and their attributes for item in list_of_hooke_curves: #playlist_element=newdoc.createElement("curve") @@ -71,17 +68,17 @@ class PlaylistXML: for key in item.__dict__: if not (key in self.hidden_attributes): newdoc.createAttribute(key) - playlist_element.setAttribute(key,str(item.__dict__[key])) - + playlist_element.setAttribute(key,str(item.__dict__[key])) + self.playlist=newdoc - + def load(self,filename): ''' loads a playlist file ''' myplay=file(filename) self.playpath=filename - + #the following 3 lines are needed to strip newlines. otherwise, since newlines #are XML elements too (why?), the parser would read them (and re-save them, multiplying #newlines...) @@ -89,24 +86,24 @@ class PlaylistXML: the_file=myplay.read() the_file_lines=the_file.split('\n') the_file=''.join(the_file_lines) - - self.playlist=xml.dom.minidom.parseString(the_file) - + + self.playlist=xml.dom.minidom.parseString(the_file) + #inner parsing functions def handlePlaylist(playlist): list_of_files=playlist.getElementsByTagName("element") generics=playlist.getElementsByTagName("generics") return handleFiles(list_of_files), handleGenerics(generics) - + def handleGenerics(generics): generics_dict={} if len(generics)==0: return generics_dict - + for attribute in generics[0].attributes.keys(): generics_dict[attribute]=generics[0].getAttribute(attribute) return generics_dict - + def handleFiles(list_of_files): new_playlist=[] for myfile in list_of_files: @@ -115,11 +112,11 @@ class PlaylistXML: for attribute in myfile.attributes.keys(): #extract attributes for the single curve the_curve.__dict__[attribute]=myfile.getAttribute(attribute) new_playlist.append(the_curve) - + return new_playlist #this is the true thing returned at the end of this function...(FIXME: clarity) - + return handlePlaylist(self.playlist) - + def save(self,output_filename): ''' @@ -130,7 +127,7 @@ class PlaylistXML: except IOError: print 'libhooke.py : Cannot save playlist. Wrong path or filename' return - + self.playlist.writexml(outfile,indent='\n') outfile.close() @@ -138,22 +135,22 @@ class PlaylistXML: class HookeConfig: ''' Handling of Hooke configuration file - + Mostly based on the simple-yet-useful examples of the Python Library Reference about xml.dom.minidom - + FIXME: starting to look a mess, should require refactoring ''' - + def __init__(self): self.config={} self.config['plugins']=[] self.config['drivers']=[] self.config['plotmanips']=[] - + def load_config(self, filename): myconfig=file(filename) - + #the following 3 lines are needed to strip newlines. otherwise, since newlines #are XML elements too, the parser would read them (and re-save them, multiplying #newlines...) @@ -161,9 +158,9 @@ class HookeConfig: the_file=myconfig.read() the_file_lines=the_file.split('\n') the_file=''.join(the_file_lines) - - self.config_tree=xml.dom.minidom.parseString(the_file) - + + self.config_tree=xml.dom.minidom.parseString(the_file) + def getText(nodelist): #take the text from a nodelist #from Python Library Reference 13.7.2 @@ -172,7 +169,7 @@ class HookeConfig: if node.nodeType == node.TEXT_NODE: rc += node.data return rc - + def handleConfig(config): display_elements=config.getElementsByTagName("display") plugins_elements=config.getElementsByTagName("plugins") @@ -186,12 +183,12 @@ class HookeConfig: handleWorkdir(workdir_elements) handleDefaultlist(defaultlist_elements) handlePlotmanip(plotmanip_elements) - + def handleDisplay(display_elements): for element in display_elements: for attribute in element.attributes.keys(): self.config[attribute]=element.getAttribute(attribute) - + def handlePlugins(plugins): for plugin in plugins[0].childNodes: try: @@ -205,28 +202,28 @@ class HookeConfig: self.config['drivers'].append(str(driver.tagName)) except: #if we allow fancy formatting of xml, there is a text node, so tagName fails for it... pass - + def handlePlotmanip(plotmanips): for plotmanip in plotmanips[0].childNodes: try: self.config['plotmanips'].append(str(plotmanip.tagName)) except: #if we allow fancy formatting of xml, there is a text node, so tagName fails for it... pass - + def handleWorkdir(workdir): ''' default working directory ''' wdir=getText(workdir[0].childNodes) self.config['workdir']=wdir.strip() - + def handleDefaultlist(defaultlist): ''' default playlist ''' dflist=getText(defaultlist[0].childNodes) self.config['defaultlist']=dflist.strip() - + handleConfig(self.config_tree) #making items in the dictionary more machine-readable for item in self.config.keys(): @@ -242,13 +239,13 @@ class HookeConfig: self.config[item]=None else: pass - + return self.config - - + + def save_config(self, config_filename): print 'Not Implemented.' - pass + pass class ClickedPoint: @@ -256,38 +253,38 @@ class ClickedPoint: this class defines what a clicked point on the curve plot is ''' def __init__(self): - + self.is_marker=None #boolean ; decides if it is a marker self.is_line_edge=None #boolean ; decides if it is the edge of a line (unused) self.absolute_coords=(None,None) #(float,float) ; the absolute coordinates of the clicked point on the graph self.graph_coords=(None,None) #(float,float) ; the coordinates of the plot that are nearest in X to the clicked point self.index=None #integer ; the index of the clicked point with respect to the vector selected self.dest=None #0 or 1 ; 0=top plot 1=bottom plot - - + + def find_graph_coords_old(self, xvector, yvector): ''' Given a clicked point on the plot, finds the nearest point in the dataset (in X) that corresponds to the clicked point. OLD & DEPRECATED - to be removed ''' - + #FIXME: a general algorithm using min() is needed! print '---DEPRECATED FIND_GRAPH_COORDS_OLD---' best_index=0 best_dist=10**9 #should be more than enough given the scale - + for index in scipy.arange(1,len(xvector),1): dist=((self.absolute_coords[0]-xvector[index])**2)+(100*((self.absolute_coords[1]-yvector[index])))**2 #TODO, generalize? y coordinate is multiplied by 100 due to scale differences in the plot if dist0: #if valid values are string we use alphainput, if it is only one we take as default if type(valid[0]) is StringType: @@ -30,15 +27,14 @@ def safeinput (message, valid=[]): return alphainput(message, valid[0], 0,[]) else: return alphainput(message,'', 1,valid) - + #if valid values are numbers we use numinput if type(valid[0]) is IntType: if len(valid)==1: return numinput(message,valid[0],1,[]) else: return numinput(message,'',1,valid) - - + def alphainput (message, default, repeat, valid): ''' @@ -46,12 +42,12 @@ def alphainput (message, default, repeat, valid): default: return value if user input was not correct (and repeat=0) repeat: keeps asking user till it gets a valid input valid: list of allowed answers, empty list for "anything" - ''' + ''' if default and not repeat: print 'Press [enter] for default: ('+str(default)+')' reply=raw_input(message) if len(valid)>0: - if reply in valid: + if reply in valid: return reply else: if repeat==1: @@ -72,12 +68,11 @@ def alphainput (message, default, repeat, valid): reply=raw_input(message) return reply - def checkalphainput (test, default, valid): #useful when input was taken form command args if len(valid)>0: - if test in valid: + if test in valid: return test else: return default @@ -95,17 +90,17 @@ def numinput(message, default, repeat, limits): default: return value if user input was not correct (and repeat=0) repeat: keeps asking user till it gets a valid input limits: pair of values, input is checked to be between them, empty list for "any number" - ''' + ''' if default and not repeat: print 'Press [enter] for default: '+str(default) - + reply=raw_input(message) - + try: intreply=int(reply) except: intreply=None - + if len(limits)==2: high=int(limits.pop()) low=int(limits.pop()) @@ -138,6 +133,7 @@ def numinput(message, default, repeat, limits): intreply=None return intreply + def checknuminput(test,default,limits): #useful when input was taken from command args if len(limits)==2: @@ -152,4 +148,3 @@ def checknuminput(test,default,limits): return int(test) else: return default - diff --git a/hooke/liboutlet.py b/hooke/liboutlet.py index 06e0940..58edf32 100644 --- a/hooke/liboutlet.py +++ b/hooke/liboutlet.py @@ -1,5 +1,3 @@ -#!/usr/bin/env python - ''' Basic outlet object @@ -8,17 +6,16 @@ Copyright (C) 2008 Alberto Gomez-Casado (University of Twente). This program is released under the GNU General Public License version 2. ''' - import re class Outlet(object): - + def __init__(self): self.buffer=[] #relations is still unused self.relations=[] - + def push(self, args): #queue new entry self.buffer.append(args) @@ -40,15 +37,15 @@ class Outlet(object): def empty(self): self.buffer=[] - + def read_last(self): return self.buffer[len(self.buffer)-1] - + def read_first(self): return self.buffer[0] def read_type(self,dtype): - #returns entries matching a given type (force, distance, point...) + #returns entries matching a given type (force, distance, point...) aux=[] index=0 if dtype=='all': @@ -57,5 +54,3 @@ class Outlet(object): if re.match(dtype+'*',i): aux.append(i) return aux - - diff --git a/hooke/libpeakspot.py b/hooke/libpeakspot.py index 398fe84..a81112d 100644 --- a/hooke/libpeakspot.py +++ b/hooke/libpeakspot.py @@ -1,5 +1,3 @@ -#!/usr/bin/env python - ''' a library of helping functions for spotting force spectroscopy peaks @@ -14,43 +12,43 @@ def conv_dx(data,vect): dim=len(data) window=len(vect) temparr=[0.0]*dim - + end=dim-window for j in range(end): for k in range(window): temparr[j]+=data[j+k]*vect[k] - return temparr - + return temparr + def absdev(arr): ''' Calculates the absolute deviation of a vector ''' med=0.0 absD=0.0 - + med=np.mean(arr) for j in arr: absD+=abs(j-med) return absD/len(arr) - + def noise_absdev(data,positive=False,maxcut=0.2,stable=0.005): ''' Returns the standard deviation of the noise. The logic is: we cut the most negative (or positive) data points until the absolute deviation becomes stable (it doesn't vary more than 0.005) or we have cut more than maxcut*len(data) points. Then calculate the absolute deviation. - + If positive=True we cut the most positive data points, False=we cut the negative ones. ''' out=[item for item in data] #we copy just to be sure... out.sort() if positive: out.reverse() - + temp_absdev=absdev(out) - + for index in range(len(out)): cutindex=(index+1)*5 cut_absdev=absdev(out[cutindex:]) #we jump five points after five points... @@ -58,18 +56,18 @@ def noise_absdev(data,positive=False,maxcut=0.2,stable=0.005): temp_absdev=cut_absdev else: break - + return cut_absdev - + def abovenoise(convoluted,noise_level,cut_index=0,abs_devs=4): ''' Generates a vector which is 0 where the vector is less than abs_devs*noise_level ; 1 if not (spike). ''' #calculate absolute noise deviation #noise_level=noise_absdev(convoluted[cut_index:]) - + above=[] - + for index in range(len(convoluted)): if indexself.config['auto_min_p'] and p_leng CLI COMMUNICATION class fitCommands: - + def _plug_init(self): self.wlccurrent=None self.wlccontact_point=None self.wlccontact_index=None - + def wlc_fit(self,clicked_points,xvector,yvector, pl_value, T=293, return_errors=False): ''' Worm-like chain model fitting. The function is the simple polynomial worm-like chain as proposed by C.Bustamante, J.F.Marko, E.D.Siggia and S.Smith (Science. 1994 Sep 9;265(5178):1599-600.) ''' - + ''' clicked_points[0] = contact point (calculated or hand-clicked) clicked_points[1] and [2] are edges of chunk ''' - + #STEP 1: Prepare the vectors to apply the fit. - - + + if pl_value is not None: pl_value=pl_value/(10**9) - + #indexes of the selected chunk first_index=min(clicked_points[1].index, clicked_points[2].index) last_index=max(clicked_points[1].index, clicked_points[2].index) - + #getting the chunk and reverting it xchunk,ychunk=xvector[first_index:last_index],yvector[first_index:last_index] xchunk.reverse() - ychunk.reverse() + ychunk.reverse() #put contact point at zero and flip around the contact point (the fit wants a positive growth for extension and force) xchunk_corr_up=[-(x-clicked_points[0].graph_coords[0]) for x in xchunk] ychunk_corr_up=[-(y-clicked_points[0].graph_coords[1]) for y in ychunk] - + #make them arrays xchunk_corr_up=scipy.array(xchunk_corr_up) ychunk_corr_up=scipy.array(ychunk_corr_up) - - + + #STEP 2: actually do the fit - + #Find furthest point of chunk and add it a bit; the fit must converge #from an excess! xchunk_high=max(xchunk_corr_up) xchunk_high+=(xchunk_high/10) - + #Here are the linearized start parameters for the WLC. #[lambd=1/Lo , pii=1/P] - + p0=[(1/xchunk_high),(1/(3.5e-10))] p0_plfix=[(1/xchunk_high)] ''' ODR STUFF fixme: remove these comments after testing ''' - - + + def f_wlc(params,x,T=T): ''' wlc function for ODR fitting @@ -98,7 +96,7 @@ class fitCommands: therm=Kb*T y=(therm*pii/4.0) * (((1-(x*lambd))**-2) - 1 + (4*x*lambd)) return y - + def f_wlc_plfix(params,x,pl_value=pl_value,T=T): ''' wlc function for ODR fitting @@ -109,7 +107,7 @@ class fitCommands: therm=Kb*T y=(therm*pii/4.0) * (((1-(x*lambd))**-2) - 1 + (4*x*lambd)) return y - + #make the ODR fit realdata=scipy.odr.RealData(xchunk_corr_up,ychunk_corr_up) if pl_value: @@ -118,11 +116,11 @@ class fitCommands: else: model=scipy.odr.Model(f_wlc) o = scipy.odr.ODR(realdata, model, p0) - + o.set_job(fit_type=2) out=o.run() fit_out=[(1/i) for i in out.beta] - + #Calculate fit errors from output standard deviations. #We must propagate the error because we fit the *inverse* parameters! #The error = (error of the inverse)*(value**2) @@ -130,8 +128,8 @@ class fitCommands: for sd,value in zip(out.sd_beta, fit_out): err_real=sd*(value**2) fit_errors.append(err_real) - - def wlc_eval(x,params,pl_value,T): + + def wlc_eval(x,params,pl_value,T): ''' Evaluates the WLC function ''' @@ -139,17 +137,17 @@ class fitCommands: lambd, pii = params else: lambd = params - + if pl_value: pii=1/pl_value - + Kb=(1.38065e-23) #boltzmann constant therm=Kb*T #so we have thermal energy - + return ( (therm*pii/4.0) * (((1-(x*lambd))**-2.0) - 1 + (4.0*x*lambd)) ) - + #STEP 3: plotting the fit - + #obtain domain to plot the fit - from contact point to last_index plus 20 points thule_index=last_index+10 if thule_index > len(xvector): #for rare cases in which we fit something at the END of whole curve. @@ -159,17 +157,17 @@ class fitCommands: xfit_chunk.reverse() xfit_chunk_corr_up=[-(x-clicked_points[0].graph_coords[0]) for x in xfit_chunk] xfit_chunk_corr_up=scipy.array(xfit_chunk_corr_up) - + #the fitted curve: reflip, re-uncorrect yfit=wlc_eval(xfit_chunk_corr_up, out.beta, pl_value,T) yfit_down=[-y for y in yfit] yfit_corr_down=[y+clicked_points[0].graph_coords[1] for y in yfit_down] - + if return_errors: return fit_out, yfit_corr_down, xfit_chunk, fit_errors else: return fit_out, yfit_corr_down, xfit_chunk, None - + def fjc_fit(self,clicked_points,xvector,yvector, pl_value, T=293, return_errors=False): ''' Freely-jointed chain function @@ -179,39 +177,39 @@ class fitCommands: clicked_points[0] is the contact point (calculated or hand-clicked) clicked_points[1] and [2] are edges of chunk ''' - + #STEP 1: Prepare the vectors to apply the fit. if pl_value is not None: pl_value=pl_value/(10**9) - + #indexes of the selected chunk first_index=min(clicked_points[1].index, clicked_points[2].index) last_index=max(clicked_points[1].index, clicked_points[2].index) - + #getting the chunk and reverting it xchunk,ychunk=xvector[first_index:last_index],yvector[first_index:last_index] xchunk.reverse() - ychunk.reverse() + ychunk.reverse() #put contact point at zero and flip around the contact point (the fit wants a positive growth for extension and force) xchunk_corr_up=[-(x-clicked_points[0].graph_coords[0]) for x in xchunk] ychunk_corr_up=[-(y-clicked_points[0].graph_coords[1]) for y in ychunk] - - + + #make them arrays xchunk_corr_up=scipy.array(xchunk_corr_up) ychunk_corr_up=scipy.array(ychunk_corr_up) - - + + #STEP 2: actually do the fit - + #Find furthest point of chunk and add it a bit; the fit must converge #from an excess! xchunk_high=max(xchunk_corr_up) xchunk_high+=(xchunk_high/10) - + #Here are the linearized start parameters for the WLC. #[lambd=1/Lo , pii=1/P] - + p0=[(1/xchunk_high),(1/(3.5e-10))] p0_plfix=[(1/xchunk_high)] ''' @@ -223,7 +221,7 @@ class fitCommands: hyperbolic cotangent ''' return (np.exp(2*z)+1)/(np.exp(2*z)-1) - + def x_fjc(params,f,T=T): ''' fjc function for ODR fitting @@ -231,11 +229,11 @@ class fitCommands: lambd,pii=params Kb=(1.38065e-23) therm=Kb*T - + #x=(therm*pii/4.0) * (((1-(x*lambd))**-2) - 1 + (4*x*lambd)) x=(1/lambd)*(coth(f*(1/pii)/therm) - (therm*pii)/f) return x - + def x_fjc_plfix(params,f,pl_value=pl_value,T=T): ''' fjc function for ODR fitting @@ -247,7 +245,7 @@ class fitCommands: #y=(therm*pii/4.0) * (((1-(x*lambd))**-2) - 1 + (4*x*lambd)) x=(1/lambd)*(coth(f*(1/pii)/therm) - (therm*pii)/f) return x - + #make the ODR fit realdata=scipy.odr.RealData(ychunk_corr_up,xchunk_corr_up) if pl_value: @@ -256,11 +254,11 @@ class fitCommands: else: model=scipy.odr.Model(x_fjc) o = scipy.odr.ODR(realdata, model, p0) - + o.set_job(fit_type=2) out=o.run() fit_out=[(1/i) for i in out.beta] - + #Calculate fit errors from output standard deviations. #We must propagate the error because we fit the *inverse* parameters! #The error = (error of the inverse)*(value**2) @@ -268,8 +266,8 @@ class fitCommands: for sd,value in zip(out.sd_beta, fit_out): err_real=sd*(value**2) fit_errors.append(err_real) - - def fjc_eval(y,params,pl_value,T): + + def fjc_eval(y,params,pl_value,T): ''' Evaluates the WLC function ''' @@ -277,16 +275,16 @@ class fitCommands: lambd, pii = params else: lambd = params - + if pl_value: pii=1/pl_value - + Kb=(1.38065e-23) #boltzmann constant therm=Kb*T #so we have thermal energy #return ( (therm*pii/4.0) * (((1-(x*lambd))**-2.0) - 1 + (4.0*x*lambd)) ) return (1/lambd)*(coth(y*(1/pii)/therm) - (therm*pii)/y) - - + + #STEP 3: plotting the fit #obtain domain to plot the fit - from contact point to last_index plus 20 points thule_index=last_index+10 @@ -302,54 +300,54 @@ class fitCommands: #or other buggy situations. Kludge to live with it now... ychunk=yvector[:thule_index] y_evalchunk=np.linspace(min(ychunk),max(ychunk),100) - + yfit_down=[-y for y in y_evalchunk] yfit_corr_down=[y+clicked_points[0].graph_coords[1] for y in yfit_down] yfit_corr_down=scipy.array(yfit_corr_down) - + #the fitted curve: reflip, re-uncorrect xfit=fjc_eval(yfit_corr_down, out.beta, pl_value,T) xfit=list(xfit) xfit.reverse() xfit_chunk_corr_up=[-(x-clicked_points[0].graph_coords[0]) for x in xfit] - + #xfit_chunk_corr_up=scipy.array(xfit_chunk_corr_up) #deltay=yfit_down[0]-yvector[clicked_points[0].index] - + #This is a terrible, terrible kludge to find the point where it should normalize (and from where it should plot) xxxdists=[] for index in scipy.arange(1,len(xfit_chunk_corr_up),1): - xxxdists.append((clicked_points[0].graph_coords[0]-xfit_chunk_corr_up[index])**2) + xxxdists.append((clicked_points[0].graph_coords[0]-xfit_chunk_corr_up[index])**2) normalize_index=xxxdists.index(min(xxxdists)) #End of kludge - + deltay=yfit_down[normalize_index]-clicked_points[0].graph_coords[1] yfit_corr_down=[y-deltay for y in yfit_down] - + if return_errors: return fit_out, yfit_corr_down[normalize_index+1:], xfit_chunk_corr_up[normalize_index+1:], fit_errors else: return fit_out, yfit_corr_down[normalize_index+1:], xfit_chunk_corr_up[normalize_index+1:], None - - + + def do_wlc(self,args): ''' WLC (fit.py plugin) - + See the fit command ''' self.do_fit(args) - + def do_fjc(self,args): ''' FJC (fit.py plugin) - + See the fit command ''' self.do_fit(args) - + def do_fit(self,args): ''' FIT @@ -358,32 +356,32 @@ class fitCommands: First you have to click a contact point. Then you have to click the two edges of the data you want to fit. - + The fit function depends on the fit_function variable. You can set it with the command "set fit_function wlc" or "set fit_function fjc" depending on the function you prefer. - - For WLC, the function is the simple polynomial worm-like chain as proposed by - C.Bustamante, J.F.Marko, E.D.Siggia and S.Smith (Science. 1994 + + For WLC, the function is the simple polynomial worm-like chain as proposed by + C.Bustamante, J.F.Marko, E.D.Siggia and S.Smith (Science. 1994 Sep 9;265(5178):1599-600.) - - For FJC, ref: + + For FJC, ref: C.Ray and B.B. Akhremitchev; http://www.chem.duke.edu/~boris/research/force_spectroscopy/fit_efjc.pdf Arguments: - pl=[value] : Use a fixed persistent length (WLC) or Kuhn length (FJC) for the fit. If pl is not given, - the fit will be a 2-variable + pl=[value] : Use a fixed persistent length (WLC) or Kuhn length (FJC) for the fit. If pl is not given, + the fit will be a 2-variable fit. DO NOT put spaces between 'pl', '=' and the value. - The value must be in nanometers. - + The value must be in nanometers. + t=[value] : Use a user-defined temperature. The value must be in kelvins; by default it is 293 K. DO NOT put spaces between 't', '=' and the value. - - noauto : allows for clicking the contact point by + + noauto : allows for clicking the contact point by hand (otherwise it is automatically estimated) the first time. If subsequent measurements are made, the same contact point clicked is used - + reclick : redefines by hand the contact point, if noauto has been used before but the user is unsatisfied of the previously choosen contact point. --------- @@ -400,10 +398,10 @@ class fitCommands: if ('t=' in arg[0:2]) or ('T=' in arg[0:2]): t_expression=arg.split('=') T=float(t_expression[1]) - + #use the currently displayed plot for the fit displayed_plot=self._get_displayed_plot() - + #handle contact point arguments correctly if 'reclick' in args.split(): print 'Click contact point' @@ -429,7 +427,7 @@ class fitCommands: contact_point.absolute_coords=displayed_plot.vectors[1][0][cindex], displayed_plot.vectors[1][1][cindex] contact_point.find_graph_coords(displayed_plot.vectors[1][0], displayed_plot.vectors[1][1]) contact_point.is_marker=True - + print 'Click edges of chunk' points=self._measure_N_points(N=2, whatset=1) points=[contact_point]+points @@ -444,11 +442,11 @@ class fitCommands: print 'No recognized fit function defined!' print 'Set your fit function to wlc or fjc.' return - + except: print 'Fit not possible. Probably wrong interval -did you click two *different* points?' return - + #FIXME: print "Kuhn length" for FJC print 'Fit function:',self.config['fit_function'] print 'Contour length: ',params[0]*(1.0e+9),' nm' @@ -458,64 +456,64 @@ class fitCommands: print name_of_charlength+': ',params[1]*(1.0e+9),' nm' to_dump='persistent '+self.current.path+' '+str(params[1]*(1.0e+9))+' nm' self.outlet.push(to_dump) - + if fit_errors: fit_nm=[i*(10**9) for i in fit_errors] print 'Standard deviation (contour length)', fit_nm[0] if len(fit_nm)>1: print 'Standard deviation ('+name_of_charlength+')', fit_nm[1] - - + + #add the clicked points in the final PlotObject clickvector_x, clickvector_y=[], [] for item in points: clickvector_x.append(item.graph_coords[0]) clickvector_y.append(item.graph_coords[1]) - + #create a custom PlotObject to gracefully plot the fit along the curves - + fitplot=copy.deepcopy(displayed_plot) fitplot.add_set(xfit,yfit) fitplot.add_set(clickvector_x,clickvector_y) - + #FIXME: this colour/styles stuff must be solved at the root! if fitplot.styles==[]: fitplot.styles=[None,None,None,'scatter'] else: fitplot.styles+=[None,'scatter'] - + if fitplot.colors==[]: fitplot.colors=[None,None,None,None] else: fitplot.colors+=[None,None] - + self._send_plot([fitplot]) - - + + #---------- - - + + def find_contact_point(self,plot=False): ''' Finds the contact point on the curve. - + The current algorithm (thanks to Francesco Musiani, francesco.musiani@unibo.it and Massimo Sandal) is: - take care of the PicoForce trigger bug - exclude retraction portions with too high standard deviation - fit the second half of the retraction curve to a line - if the fit is not almost horizontal, take a smaller chunk and repeat - otherwise, we have something horizontal - so take the average of horizontal points and use it as a baseline - + Then, start from the rise of the retraction curve and look at the first point below the baseline. - + FIXME: should be moved, probably to generalvclamp.py ''' - + if not plot: plot=self.plots[0] - + outplot=self.subtract_curves(1) xret=outplot.vectors[1][0] ydiff=outplot.vectors[1][1] @@ -524,7 +522,7 @@ class fitCommands: yext=plot.vectors[0][1] xret2=plot.vectors[1][0] yret=plot.vectors[1][1] - + #taking care of the picoforce trigger bug: we exclude portions of the curve that have too much #standard deviation. yes, a lot of magic is here. monster=True @@ -537,13 +535,13 @@ class fitCommands: else: #move away from the monster monlength-=int(len(xret)/50) finalength-=int(len(xret)/50) - - + + #take half of the thing endlength=int(len(xret)/2) - + ok=False - + while not ok: xchunk=yext[endlength:monlength] ychunk=yext[endlength:monlength] @@ -553,36 +551,36 @@ class fitCommands: if (abs(regr[1]) > 0.1) and ( endlength < len(xret)-int(len(xret)/6) ) : endlength+=10 else: - ok=True - - + ok=True + + ymean=np.mean(ychunk) #baseline - + index=0 point = ymean+1 - + #find the first point below the calculated baseline while point > ymean: try: point=yret[index] - index+=1 + index+=1 except IndexError: #The algorithm didn't find anything below the baseline! It should NEVER happen - index=0 + index=0 return index - + return index - - - + + + def find_contact_point2(self, debug=False): ''' TO BE DEVELOPED IN THE FUTURE Finds the contact point on the curve. - + FIXME: should be moved, probably to generalvclamp.py ''' - + #raw_plot=self.current.curve.default_plots()[0] raw_plot=self.plots[0] '''xext=self.plots[0].vectors[0][0] @@ -594,43 +592,43 @@ class fitCommands: yext=raw_plot.vectors[0][1] xret2=raw_plot.vectors[1][0] yret=raw_plot.vectors[1][1] - + first_point=[xext[0], yext[0]] last_point=[xext[-1], yext[-1]] - + #regr=scipy.polyfit(first_point, last_point,1)[0:2] diffx=abs(first_point[0]-last_point[0]) diffy=abs(first_point[1]-last_point[1]) - + #using polyfit results in numerical errors. good old algebra. a=diffy/diffx b=first_point[1]-(a*first_point[0]) baseline=scipy.polyval((a,b), xext) - + ysub=[item-basitem for item,basitem in zip(yext,baseline)] - + contact=ysub.index(min(ysub)) - + return xext,ysub,contact - + #now, exploit a ClickedPoint instance to calculate index... dummy=ClickedPoint() dummy.absolute_coords=(x_intercept,y_intercept) dummy.find_graph_coords(xret2,yret) - + if debug: return dummy.index, regr, regr_contact else: return dummy.index - - + + def x_do_contact(self,args): ''' DEBUG COMMAND to be activated in the future ''' xext,ysub,contact=self.find_contact_point2(debug=True) - + contact_plot=self.plots[0] contact_plot.add_set(xext,ysub) contact_plot.add_set([xext[contact]],[self.plots[0].vectors[0][1][contact]]) @@ -639,8 +637,8 @@ class fitCommands: contact_plot.styles=[None,None,None,'scatter'] self._send_plot([contact_plot]) return - - + + index,regr,regr_contact=self.find_contact_point2(debug=True) print regr print regr_contact @@ -649,8 +647,8 @@ class fitCommands: #nc_line=[(item*regr[0])+regr[1] for item in x_nc] nc_line=scipy.polyval(regr,xret) c_line=scipy.polyval(regr_contact,xret) - - + + contact_plot=self.current.curve.default_plots()[0] contact_plot.add_set(xret, nc_line) contact_plot.add_set(xret, c_line) @@ -658,4 +656,3 @@ class fitCommands: #contact_plot.styles.append(None) contact_plot.destination=1 self._send_plot([contact_plot]) - \ No newline at end of file diff --git a/hooke/plugin/flatfilts.py b/hooke/plugin/flatfilts.py index 3e9f9c8..a2cc569 100755 --- a/hooke/plugin/flatfilts.py +++ b/hooke/plugin/flatfilts.py @@ -1,5 +1,3 @@ -#!/usr/bin/env python - ''' FLATFILTS @@ -27,19 +25,19 @@ import libhookecurve as lhc class flatfiltsCommands: - + def _plug_init(self): #configurate convfilt variables convfilt_configurator=ConvfiltConfig() - + #different OSes have different path conventions if self.config['hookedir'][0]=='/': slash='/' #a Unix or Unix-like system else: slash='\\' #it's a drive letter, we assume it's Windows - + self.convfilt_config=convfilt_configurator.load_config(self.config['hookedir']+slash+'convfilt.conf') - + def do_flatfilt(self,args): ''' FLATFILT @@ -63,33 +61,33 @@ class flatfiltsCommands: median_filter=7 min_npks=4 min_deviation=9 - + args=args.split(' ') if len(args) == 2: min_npks=int(args[0]) min_deviation=int(args[1]) else: pass - + print 'Processing playlist...' notflat_list=[] - + c=0 for item in self.current_list: c+=1 - + try: notflat=self.has_features(item, median_filter, min_npks, min_deviation) print 'Curve',item.path, 'is',c,'of',len(self.current_list),': features are ',notflat except: notflat=False print 'Curve',item.path, 'is',c,'of',len(self.current_list),': cannot be filtered. Probably unable to retrieve force data from corrupt file.' - + if notflat: item.features=notflat item.curve=None #empty the item object, to further avoid memory leak notflat_list.append(item) - + if len(notflat_list)==0: print 'Found nothing interesting. Check your playlist, could be a bug or criteria could be too much stringent' return @@ -100,28 +98,28 @@ class flatfiltsCommands: self.current_list=notflat_list self.current=self.current_list[self.pointer] self.do_plot(0) - + def has_features(self,item,median_filter,min_npks,min_deviation): ''' decides if a curve is flat enough to be rejected from analysis: it sees if there are at least min_npks points that are higher than min_deviation times the absolute value of noise. - + Algorithm original idea by Francesco Musiani, with my tweaks and corrections. ''' retvalue=False - - item.identify(self.drivers) + + item.identify(self.drivers) #we assume the first is the plot with the force curve #do the median to better resolve features from noise flat_plot=self.plotmanip_median(item.curve.default_plots()[0], item, customvalue=median_filter) - flat_vects=flat_plot.vectors + flat_vects=flat_plot.vectors item.curve.close_all() #needed to avoid *big* memory leaks! del item.curve del item - - #absolute value of derivate + + #absolute value of derivate yretdiff=diff(flat_vects[1][1]) yretdiff=[abs(value) for value in yretdiff] #average of derivate values @@ -134,19 +132,19 @@ class flatfiltsCommands: c_pks+=1 else: break - + if c_pks>=min_npks: retvalue = c_pks - + del flat_plot, flat_vects, yretdiff - + return retvalue ################################################################ #-----CONVFILT------------------------------------------------- #-----Convolution-based peak recognition and filtering. #Requires the libpeakspot.py library - + def has_peaks(self, plot, abs_devs=None): ''' Finds peak position in a force curve. @@ -154,13 +152,13 @@ class flatfiltsCommands: ''' if abs_devs==None: abs_devs=self.convfilt_config['mindeviation'] - - + + xret=plot.vectors[1][0] yret=plot.vectors[1][1] #Calculate convolution. convoluted=lps.conv_dx(yret, self.convfilt_config['convolution']) - + #surely cut everything before the contact point cut_index=self.find_contact_point(plot) #cut even more, before the blind window @@ -172,34 +170,34 @@ class flatfiltsCommands: blind_index+=1 cut_index+=blind_index #do the dirty convolution-peak finding stuff - noise_level=lps.noise_absdev(convoluted[cut_index:], self.convfilt_config['positive'], self.convfilt_config['maxcut'], self.convfilt_config['stable']) - above=lps.abovenoise(convoluted,noise_level,cut_index,abs_devs) + noise_level=lps.noise_absdev(convoluted[cut_index:], self.convfilt_config['positive'], self.convfilt_config['maxcut'], self.convfilt_config['stable']) + above=lps.abovenoise(convoluted,noise_level,cut_index,abs_devs) peak_location,peak_size=lps.find_peaks(above,seedouble=self.convfilt_config['seedouble']) - + #take the maximum for i in range(len(peak_location)): peak=peak_location[i] maxpk=min(yret[peak-10:peak+10]) index_maxpk=yret[peak-10:peak+10].index(maxpk)+(peak-10) peak_location[i]=index_maxpk - + return peak_location,peak_size - - + + def exec_has_peaks(self,item,abs_devs): ''' encapsulates has_peaks for the purpose of correctly treating the curve objects in the convfilt loop, to avoid memory leaks ''' - item.identify(self.drivers) + item.identify(self.drivers) #we assume the first is the plot with the force curve plot=item.curve.default_plots()[0] - + if 'flatten' in self.config['plotmanips']: #If flatten is present, use it for better recognition of peaks... flatten=self._find_plotmanip('flatten') #extract flatten plot manipulator plot=flatten(plot, item, customvalue=1) - + peak_location,peak_size=self.has_peaks(plot,abs_devs) #close all open files item.curve.close_all() @@ -207,10 +205,10 @@ class flatfiltsCommands: del item.curve del item return peak_location, peak_size - + #------------------------ #------commands---------- - #------------------------ + #------------------------ def do_peaks(self,args): ''' PEAKS @@ -223,37 +221,37 @@ class flatfiltsCommands: ''' if len(args)==0: args=self.convfilt_config['mindeviation'] - + try: abs_devs=float(args) except: print 'Wrong argument, using config value' abs_devs=float(self.convfilt_config['mindeviation']) - + defplots=self.current.curve.default_plots()[0] #we need the raw, uncorrected plots - + if 'flatten' in self.config['plotmanips']: flatten=self._find_plotmanip('flatten') #extract flatten plot manipulator defplots=flatten(defplots, self.current) else: print 'You have the flatten plot manipulator not loaded. Enabling it could give you better results.' - + peak_location,peak_size=self.has_peaks(defplots,abs_devs) print 'Found '+str(len(peak_location))+' peaks.' to_dump='peaks '+self.current.path+' '+str(len(peak_location)) self.outlet.push(to_dump) #print peak_location - + #if no peaks, we have nothing to plot. exit. if len(peak_location)==0: return - + #otherwise, we plot the peak locations. xplotted_ret=self.plots[0].vectors[1][0] yplotted_ret=self.plots[0].vectors[1][1] xgood=[xplotted_ret[index] for index in peak_location] ygood=[yplotted_ret[index] for index in peak_location] - + recplot=self._get_displayed_plot() recplot.vectors.append([xgood,ygood]) if recplot.styles==[]: @@ -262,9 +260,9 @@ class flatfiltsCommands: else: recplot.styles+=['scatter'] recplot.colors+=[None] - + self._send_plot([recplot]) - + def do_convfilt(self,args): ''' CONVFILT @@ -284,26 +282,26 @@ class flatfiltsCommands: If called without arguments, it uses default values. ''' - + min_npks=self.convfilt_config['minpeaks'] min_deviation=self.convfilt_config['mindeviation'] - + args=args.split(' ') if len(args) == 2: min_npks=int(args[0]) min_deviation=int(args[1]) else: pass - + print 'Processing playlist...' print '(Please wait)' notflat_list=[] - + c=0 for item in self.current_list: c+=1 - - try: + + try: peak_location,peak_size=self.exec_has_peaks(item,min_deviation) if len(peak_location)>=min_npks: isok='+' @@ -313,18 +311,18 @@ class flatfiltsCommands: except: peak_location,peak_size=[],[] print 'Curve',item.path, 'is',c,'of',len(self.current_list),': cannot be filtered. Probably unable to retrieve force data from corrupt file.' - + if len(peak_location)>=min_npks: item.peak_location=peak_location item.peak_size=peak_size item.curve=None #empty the item object, to further avoid memory leak notflat_list.append(item) - + #Warn that no flattening had been done. if not ('flatten' in self.config['plotmanips']): print 'Flatten manipulator was not found. Processing was done without flattening.' print 'Try to enable it in your configuration file for better results.' - + if len(notflat_list)==0: print 'Found nothing interesting. Check your playlist, could be a bug or criteria could be too much stringent' return @@ -335,8 +333,8 @@ class flatfiltsCommands: self.current_list=notflat_list self.current=self.current_list[self.pointer] self.do_plot(0) - - + + def do_setconv(self,args): ''' SETCONV @@ -354,7 +352,7 @@ class flatfiltsCommands: print 'This is not an internal convfilt variable!' print 'Run "setconv" without arguments to see a list of defined variables.' return - + if len(args)==1: print self.convfilt_config[args[0]] elif len(args)>1: @@ -369,19 +367,19 @@ class flatfiltsCommands: class ConvfiltConfig: ''' Handling of convfilt configuration file - + Mostly based on the simple-yet-useful examples of the Python Library Reference about xml.dom.minidom - + FIXME: starting to look a mess, should require refactoring ''' - + def __init__(self): self.config={} - - + + def load_config(self, filename): - myconfig=file(filename) + myconfig=file(filename) #the following 3 lines are needed to strip newlines. otherwise, since newlines #are XML elements too, the parser would read them (and re-save them, multiplying #newlines...) @@ -389,9 +387,9 @@ class ConvfiltConfig: the_file=myconfig.read() the_file_lines=the_file.split('\n') the_file=''.join(the_file_lines) - - self.config_tree=xml.dom.minidom.parseString(the_file) - + + self.config_tree=xml.dom.minidom.parseString(the_file) + def getText(nodelist): #take the text from a nodelist #from Python Library Reference 13.7.2 @@ -400,23 +398,23 @@ class ConvfiltConfig: if node.nodeType == node.TEXT_NODE: rc += node.data return rc - + def handleConfig(config): noiseabsdev_elements=config.getElementsByTagName("noise_absdev") convfilt_elements=config.getElementsByTagName("convfilt") handleAbsdev(noiseabsdev_elements) handleConvfilt(convfilt_elements) - + def handleAbsdev(noiseabsdev_elements): for element in noiseabsdev_elements: for attribute in element.attributes.keys(): self.config[attribute]=element.getAttribute(attribute) - + def handleConvfilt(convfilt_elements): for element in convfilt_elements: for attribute in element.attributes.keys(): self.config[attribute]=element.getAttribute(attribute) - + handleConfig(self.config_tree) #making items in the dictionary machine-readable for item in self.config.keys(): @@ -424,5 +422,5 @@ class ConvfiltConfig: self.config[item]=eval(self.config[item]) except NameError: #if it's an unreadable string, keep it as a string pass - - return self.config \ No newline at end of file + + return self.config diff --git a/hooke/plugin/generalclamp.py b/hooke/plugin/generalclamp.py index 4ac30e7..327e585 100644 --- a/hooke/plugin/generalclamp.py +++ b/hooke/plugin/generalclamp.py @@ -1,73 +1,71 @@ -#!/usr/bin/env python - ''' GENERALCLAMP.py Plugin regarding general force clamp measurements ''' from libhooke import WX_GOOD, ClickedPoint -import wxversion +import wxversion import libhookecurve as lhc wxversion.select(WX_GOOD) -from wx import PostEvent +from wx import PostEvent + +class generalclampCommands: -class generalclampCommands: - def plotmanip_clamp(self, plot, current, customvalue=False): ''' - Handles some viewing options for the "force clamp" data format, depending on the state of these configuration variables: - (1) If self.config['fc_showphase'] != 0, the 'phase' data column (i.e. the 2nd) is shown in the 0th graph (else it isn't) - (2) If self.config['fc_showimposed'] != 0, the 'imposed deflection' data column (i.e. the 5th) is shown in the 1st graph (else it isn't) + Handles some viewing options for the "force clamp" data format, depending on the state of these configuration variables: + (1) If self.config['fc_showphase'] != 0, the 'phase' data column (i.e. the 2nd) is shown in the 0th graph (else it isn't) + (2) If self.config['fc_showimposed'] != 0, the 'imposed deflection' data column (i.e. the 5th) is shown in the 1st graph (else it isn't) (3) If self.config['fc_interesting'] == 0, the entire curve is shown in the graphs; if it has a non-zero value N, only phase N is shown. - - NOTE - my implementation of point(3) feels quite awkward - someone smarter than me plz polish that! - + + NOTE - my implementation of point(3) feels quite awkward - someone smarter than me plz polish that! + ''' - + #not a fclamp curve... if current.curve.experiment != 'clamp': - return plot - - if self.config['fc_interesting'] != 0 and plot.destination==0: - lower = int((self.config['fc_interesting'])-1) - upper = int((self.config['fc_interesting'])+1) - trim = current.curve.trimindexes()[lower:upper] - newtime = [] - newzpiezo = [] - newphase = [] - for x in range(trim[0],trim[1]): - newtime.append(self.plots[0].vectors[0][0][x]) - newzpiezo.append(self.plots[0].vectors[0][1][x]) - newphase.append(self.plots[0].vectors[1][1][x]) - self.plots[0].vectors[0][0] = newtime - self.plots[0].vectors[0][1] = newzpiezo - self.plots[0].vectors[1][0] = newtime - self.plots[0].vectors[1][1] = newphase - - if self.config['fc_interesting'] != 0 and plot.destination==1: - lower = int((self.config['fc_interesting'])-1) - upper = int((self.config['fc_interesting'])+1) - trim = current.curve.trimindexes()[lower:upper] - newtime = [] - newdefl = [] - newimposed = [] - for x in range(trim[0],trim[1]): - newtime.append(self.plots[1].vectors[0][0][x]) - newdefl.append(self.plots[1].vectors[0][1][x]) - newimposed.append(self.plots[1].vectors[1][1][x]) - self.plots[1].vectors[0][0] = newtime - self.plots[1].vectors[0][1] = newdefl - self.plots[1].vectors[1][0] = newtime - self.plots[1].vectors[1][1] = newimposed - - if self.config['fc_showphase'] == 0 and plot.destination==0: - self.plots[0].remove_set(1) - - if self.config['fc_showimposed'] == 0 and plot.destination==1: - self.plots[1].remove_set(1) - + return plot + + if self.config['fc_interesting'] != 0 and plot.destination==0: + lower = int((self.config['fc_interesting'])-1) + upper = int((self.config['fc_interesting'])+1) + trim = current.curve.trimindexes()[lower:upper] + newtime = [] + newzpiezo = [] + newphase = [] + for x in range(trim[0],trim[1]): + newtime.append(self.plots[0].vectors[0][0][x]) + newzpiezo.append(self.plots[0].vectors[0][1][x]) + newphase.append(self.plots[0].vectors[1][1][x]) + self.plots[0].vectors[0][0] = newtime + self.plots[0].vectors[0][1] = newzpiezo + self.plots[0].vectors[1][0] = newtime + self.plots[0].vectors[1][1] = newphase + + if self.config['fc_interesting'] != 0 and plot.destination==1: + lower = int((self.config['fc_interesting'])-1) + upper = int((self.config['fc_interesting'])+1) + trim = current.curve.trimindexes()[lower:upper] + newtime = [] + newdefl = [] + newimposed = [] + for x in range(trim[0],trim[1]): + newtime.append(self.plots[1].vectors[0][0][x]) + newdefl.append(self.plots[1].vectors[0][1][x]) + newimposed.append(self.plots[1].vectors[1][1][x]) + self.plots[1].vectors[0][0] = newtime + self.plots[1].vectors[0][1] = newdefl + self.plots[1].vectors[1][0] = newtime + self.plots[1].vectors[1][1] = newimposed + + if self.config['fc_showphase'] == 0 and plot.destination==0: + self.plots[0].remove_set(1) + + if self.config['fc_showimposed'] == 0 and plot.destination==1: + self.plots[1].remove_set(1) + return plot - + def do_time(self,args): ''' Measures the time difference (in seconds) between two points @@ -75,12 +73,12 @@ class generalclampCommands: ---- Syntax: time ''' - if self.current.curve.experiment == 'clamp': - time=self._delta(set=0)[0] + if self.current.curve.experiment == 'clamp': + time=self._delta(set=0)[0] print str(time*1000)+' ms' else: print 'This command makes no sense for a non-force clamp experiment.' - + def do_zpiezo(self,args): ''' Measures the zpiezo difference (in nm) between two points @@ -89,11 +87,11 @@ class generalclampCommands: Syntax: zpiezo ''' if self.current.curve.experiment == 'clamp': - zpiezo=self._delta(set=0)[2] + zpiezo=self._delta(set=0)[2] print str(zpiezo*(10**9))+' nm' else: print 'This command makes no sense for a non-force clamp experiment.' - + def do_defl(self,args): ''' Measures the deflection difference (in nm) between two points @@ -108,7 +106,7 @@ class generalclampCommands: print str(defl*(10**12))+' pN' else: print 'This command makes no sense for a non-force clamp experiment.' - + def do_step(self,args): ''' Measures the length and time duration of a time-Z step @@ -126,126 +124,126 @@ class generalclampCommands: dt=abs(points[1].graph_coords[0]-points[0].graph_coords[0]) print 'dZ: ',dz,' nm' print 'dT: ',dt,' s' - + + else: + print 'This command makes no sense for a non-force clamp experiment.' + + def do_fcfilt(self,args): + ''' + Filters out featureless force clamp curves of the current playlist. + It's very similar to 'flatfilt' for velocity clamp curves. + Creates a new playlist only containing non-empty curves. + + WARNING - Only works if you set an appropriate fc_interesting config variable! + WARNING - arguments are NOT optional at the moment! + + Syntax: fcfilt maxretraction(nm) mindeviation (pN) + + Suggested values for an (i27)8 experiment with our setup are 200nm and 10-15 pN + ''' + + if self.config['fc_interesting'] == 0: + print 'You must specify the phase of interest (using set fc_interesing X) prior to running fcfilt!' + return + + maxretraction=0 + threshold=0 + args=args.split(' ') + if len(args)==2: + maxretraction=int(args[0]) + threshold=int(args[1]) + else: + print 'Arguments are not optional for fcfilt. You should pass two numbers:' + print '(1) the maximum plausible piezo retraction in NANOMETERS (e.g. the length of the protein)' + print "(2) the threshold, in PICONEWTONS. If signal deviates from imposed more than this, it's an event" + return + + + print 'Processing playlist... go get yourself a cup of coffee.' + notflat_list=[] + + c=0 + + for item in self.current_list: + c+=1 + try: + notflat=self.has_stuff(item,maxretraction,threshold) + print 'Curve',item.path,'is',c,'of',len(self.current_list),'--->Has Stuff =',notflat + except: + notflat=False + print 'Curve',item.path,'is',c,'of',len(self.current_list),'--->could not be processed' + if notflat: + item.features=notflat + item.curve=None + notflat_list.append(item) + + if len(notflat_list)==0: + print 'Nothing interesting here. Reconsider either your filtering criteria or your experimental data' + return + else: + print 'Found',len(notflat_list),'potentially interesting curves.' + print 'Regenerating Playlist...' + self.pointer=0 + self.current_list=notflat_list + self.current=self.current_list[self.pointer] + self.do_plot(0) + + def has_stuff(self,item,maxretraction,threshold): + ''' + Decides whether a curve has some features in the interesting phase. + Algorithm: + - clip the interesting phase portion of the curve. + - discard the first 20 milliseconds (this is due to a quirk of our hardware). + - look at the zpiezo plot and note down when (if) retratcs more than [maxretraction] nm away from the first point. + - clip off any data after this point, with an excess of 100 points (again, an hardware quirk) + - if the remainder is less than 100 points, ditch the curve. + - now look at the deflection plot and check if there are points more than [threshold] pN over the 'flat zone'. + - if you find such points, bingo! + ''' + + item.identify(self.drivers) + + lower = int((self.config['fc_interesting'])-1) + upper = int((self.config['fc_interesting'])+1) + trim_idxs = item.curve.trimindexes()[lower:upper] + lo=trim_idxs[0]+20 #clipping the first 20 points off... + hi=trim_idxs[1] + trimmed_zpiezo=item.curve.default_plots()[0].vectors[0][1][lo:hi] + trimmed_defl=item.curve.default_plots()[1].vectors[0][1][lo:hi] + trimmed_imposed=item.curve.default_plots()[1].vectors[1][1][lo:hi] + imposed=trimmed_imposed[21] #just to match the 20-pts clipping... + + item.curve.close_all() + del item.curve + del item + + starting_z=trimmed_zpiezo[0] + plausible=starting_z-(maxretraction*1e-9) + det_trim=0 + while trimmed_zpiezo[det_trim]>plausible: + det_trim+=1 + if det_trim >= len(trimmed_zpiezo): #breaking cycles makes me shiver... + det_trim=len(trimmed_zpiezo) #but I cannot think of anything better now. + break + further_trim=det_trim-100 + if further_trim<100: + return False + trimmed_defl=trimmed_defl[:further_trim] + + trimmed_defl.sort() + ninetypercent=int(0.9*len(trimmed_defl)) + j=0 + sum=0 + for j in trimmed_defl[:ninetypercent]: + sum+=j + avg=float(sum/ninetypercent) + sweetspot=float(avg+(threshold*1e-12)) + if trimmed_defl[-1]>sweetspot: + flag=True else: - print 'This command makes no sense for a non-force clamp experiment.' - - def do_fcfilt(self,args): - ''' - Filters out featureless force clamp curves of the current playlist. - It's very similar to 'flatfilt' for velocity clamp curves. - Creates a new playlist only containing non-empty curves. - - WARNING - Only works if you set an appropriate fc_interesting config variable! - WARNING - arguments are NOT optional at the moment! - - Syntax: fcfilt maxretraction(nm) mindeviation (pN) - - Suggested values for an (i27)8 experiment with our setup are 200nm and 10-15 pN - ''' - - if self.config['fc_interesting'] == 0: - print 'You must specify the phase of interest (using set fc_interesing X) prior to running fcfilt!' - return - - maxretraction=0 - threshold=0 - args=args.split(' ') - if len(args)==2: - maxretraction=int(args[0]) - threshold=int(args[1]) - else: - print 'Arguments are not optional for fcfilt. You should pass two numbers:' - print '(1) the maximum plausible piezo retraction in NANOMETERS (e.g. the length of the protein)' - print "(2) the threshold, in PICONEWTONS. If signal deviates from imposed more than this, it's an event" - return - - - print 'Processing playlist... go get yourself a cup of coffee.' - notflat_list=[] - - c=0 - - for item in self.current_list: - c+=1 - try: - notflat=self.has_stuff(item,maxretraction,threshold) - print 'Curve',item.path,'is',c,'of',len(self.current_list),'--->Has Stuff =',notflat - except: - notflat=False - print 'Curve',item.path,'is',c,'of',len(self.current_list),'--->could not be processed' - if notflat: - item.features=notflat - item.curve=None - notflat_list.append(item) - - if len(notflat_list)==0: - print 'Nothing interesting here. Reconsider either your filtering criteria or your experimental data' - return - else: - print 'Found',len(notflat_list),'potentially interesting curves.' - print 'Regenerating Playlist...' - self.pointer=0 - self.current_list=notflat_list - self.current=self.current_list[self.pointer] - self.do_plot(0) - - def has_stuff(self,item,maxretraction,threshold): - ''' - Decides whether a curve has some features in the interesting phase. - Algorithm: - - clip the interesting phase portion of the curve. - - discard the first 20 milliseconds (this is due to a quirk of our hardware). - - look at the zpiezo plot and note down when (if) retratcs more than [maxretraction] nm away from the first point. - - clip off any data after this point, with an excess of 100 points (again, an hardware quirk) - - if the remainder is less than 100 points, ditch the curve. - - now look at the deflection plot and check if there are points more than [threshold] pN over the 'flat zone'. - - if you find such points, bingo! - ''' - - item.identify(self.drivers) - - lower = int((self.config['fc_interesting'])-1) - upper = int((self.config['fc_interesting'])+1) - trim_idxs = item.curve.trimindexes()[lower:upper] - lo=trim_idxs[0]+20 #clipping the first 20 points off... - hi=trim_idxs[1] - trimmed_zpiezo=item.curve.default_plots()[0].vectors[0][1][lo:hi] - trimmed_defl=item.curve.default_plots()[1].vectors[0][1][lo:hi] - trimmed_imposed=item.curve.default_plots()[1].vectors[1][1][lo:hi] - imposed=trimmed_imposed[21] #just to match the 20-pts clipping... - - item.curve.close_all() - del item.curve - del item - - starting_z=trimmed_zpiezo[0] - plausible=starting_z-(maxretraction*1e-9) - det_trim=0 - while trimmed_zpiezo[det_trim]>plausible: - det_trim+=1 - if det_trim >= len(trimmed_zpiezo): #breaking cycles makes me shiver... - det_trim=len(trimmed_zpiezo) #but I cannot think of anything better now. - break - further_trim=det_trim-100 - if further_trim<100: - return False - trimmed_defl=trimmed_defl[:further_trim] - - trimmed_defl.sort() - ninetypercent=int(0.9*len(trimmed_defl)) - j=0 - sum=0 - for j in trimmed_defl[:ninetypercent]: - sum+=j - avg=float(sum/ninetypercent) - sweetspot=float(avg+(threshold*1e-12)) - if trimmed_defl[-1]>sweetspot: - flag=True - else: - flag=False - - del trimmed_defl,trimmed_zpiezo,trimmed_imposed - - return flag - \ No newline at end of file + flag=False + + del trimmed_defl,trimmed_zpiezo,trimmed_imposed + + return flag + diff --git a/hooke/plugin/generaltccd.py b/hooke/plugin/generaltccd.py index 70ad513..1b95049 100644 --- a/hooke/plugin/generaltccd.py +++ b/hooke/plugin/generaltccd.py @@ -1,5 +1,3 @@ -#!/usr/bin/env python - ''' generaltccd.py @@ -7,24 +5,24 @@ General utilities for TCCD stuff ''' class generaltccdCommands: - + def plotmanip_threshold(self, plot, current, customvalue=False): ''' Cuts from the plot everything below the threshold. Set the threshold with "set tccd_threshold" ''' - + if current.curve.experiment != 'smfluo': return plot - + if not self.config['tccd_threshold'] and (not customvalue): return plot - + if customvalue: thresh=customvalue else: thresh=self.config['tccd_threshold'] - + for set in plot.vectors: newy=[] for value in set[1]: @@ -32,11 +30,11 @@ class generaltccdCommands: newy.append(0) else: newy.append(value) - + set[1]=newy - + return plot - + def plotmanip_coincident(self,plot,current, customvalue=False): ''' @@ -44,10 +42,10 @@ class generaltccdCommands: ''' if current.curve.experiment != 'smfluo': return plot - + if not self.config['tccd_coincident'] and (not customvalue): return plot - + newred=[] newblue=[] for index in range(len(plot.vectors[0][1])): @@ -57,8 +55,8 @@ class generaltccdCommands: else: newred.append(0) newblue.append(0) - + plot.vectors[0][1]=newred plot.vectors[1][1]=newblue - - return plot \ No newline at end of file + + return plot diff --git a/hooke/plugin/generalvclamp.py b/hooke/plugin/generalvclamp.py index 1c71de0..32beaf0 100644 --- a/hooke/plugin/generalvclamp.py +++ b/hooke/plugin/generalvclamp.py @@ -1,5 +1,3 @@ -#!/usr/bin/env python - ''' generalvclamp.py @@ -21,12 +19,12 @@ warnings.simplefilter('ignore',np.RankWarning) class generalvclampCommands: - + def _plug_init(self): self.basecurrent=None self.basepoints=None self.autofile='' - + def do_distance(self,args): ''' DISTANCE @@ -55,7 +53,7 @@ class generalvclampCommands: Measure the force difference (in pN) between two points --------------- Syntax: force - ''' + ''' if self.current.curve.experiment == 'clamp': print 'This command makes no sense for a force clamp experiment.' return @@ -63,15 +61,15 @@ class generalvclampCommands: print str(dy*(10**12))+' pN' to_dump='force '+self.current.path+' '+str(dy*(10**12))+' pN' self.outlet.push(to_dump) - - + + def do_forcebase(self,args): ''' FORCEBASE (generalvclamp.py) Measures the difference in force (in pN) between a point and a baseline took as the average between two points. - + The baseline is fixed once for a given curve and different force measurements, unless the user wants it to be recalculated ------------ @@ -82,19 +80,19 @@ class generalvclampCommands: ''' rebase=False #if true=we select rebase maxpoint=False #if true=we measure the maximum peak - + plot=self._get_displayed_plot() whatset=1 #fixme: for all sets if 'rebase' in args or (self.basecurrent != self.current.path): rebase=True if 'max' in args: maxpoint=True - + if rebase: print 'Select baseline' self.basepoints=self._measure_N_points(N=2, whatset=whatset) self.basecurrent=self.current.path - + if maxpoint: print 'Select two points' points=self._measure_N_points(N=2, whatset=whatset) @@ -109,12 +107,12 @@ class generalvclampCommands: points=self._measure_N_points(N=1, whatset=whatset) #whatplot=points[0].dest y=points[0].graph_coords[1] - + #fixme: code duplication boundaries=[self.basepoints[0].index, self.basepoints[1].index] boundaries.sort() to_average=plot.vectors[whatset][1][boundaries[0]:boundaries[1]] #y points to average - + avg=np.mean(to_average) forcebase=abs(y-avg) print str(forcebase*(10**12))+' pN' @@ -126,15 +124,15 @@ class generalvclampCommands: Multiplies all the Y values of an SMFS curve by a value stored in the 'force_multiplier' configuration variable. Useful for calibrations and other stuff. ''' - + #not a smfs curve... if current.curve.experiment != 'smfs': return plot - + #only one set is present... if len(self.plots[0].vectors) != 2: return plot - + #multiplier is 1... if (self.config['force_multiplier']==1): return plot @@ -146,45 +144,45 @@ class generalvclampCommands: plot.vectors[1][1][i]=plot.vectors[1][1][i]*self.config['force_multiplier'] return plot - - + + def plotmanip_flatten(self, plot, current, customvalue=False): ''' Subtracts a polynomial fit to the non-contact part of the curve, as to flatten it. - the best polynomial fit is chosen among polynomials of degree 1 to n, where n is + the best polynomial fit is chosen among polynomials of degree 1 to n, where n is given by the configuration file or by the customvalue. - + customvalue= int (>0) --> starts the function even if config says no (default=False) ''' - + #not a smfs curve... if current.curve.experiment != 'smfs': return plot - + #only one set is present... if len(self.plots[0].vectors) != 2: return plot - + #config is not flatten, and customvalue flag is false too if (not self.config['flatten']) and (not customvalue): return plot - + max_exponent=12 delta_contact=0 - + if customvalue: max_cycles=customvalue else: max_cycles=self.config['flatten'] #Using > 1 usually doesn't help and can give artefacts. However, it could be useful too. - + contact_index=self.find_contact_point() - + valn=[[] for item in range(max_exponent)] yrn=[0.0 for item in range(max_exponent)] errn=[0.0 for item in range(max_exponent)] - + for i in range(int(max_cycles)): - + x_ext=plot.vectors[0][0][contact_index+delta_contact:] y_ext=plot.vectors[0][1][contact_index+delta_contact:] x_ret=plot.vectors[1][0][contact_index+delta_contact:] @@ -200,22 +198,22 @@ class generalvclampCommands: return plot best_exponent=errn.index(min(errn)) - + #extension ycorr_ext=y_ext-yrn[best_exponent]+y_ext[0] #noncontact part - yjoin_ext=np.array(plot.vectors[0][1][0:contact_index+delta_contact]) #contact part + yjoin_ext=np.array(plot.vectors[0][1][0:contact_index+delta_contact]) #contact part #retraction ycorr_ret=y_ret-yrn[best_exponent]+y_ext[0] #noncontact part yjoin_ret=np.array(plot.vectors[1][1][0:contact_index+delta_contact]) #contact part - + ycorr_ext=np.concatenate((yjoin_ext, ycorr_ext)) ycorr_ret=np.concatenate((yjoin_ret, ycorr_ret)) - + plot.vectors[0][1]=list(ycorr_ext) plot.vectors[1][1]=list(ycorr_ret) - + return plot - + #---SLOPE--- def do_slope(self,args): ''' @@ -251,7 +249,7 @@ class generalvclampCommands: clickedpoints=[] points=self._measure_N_points(N=1,whatset=1) clickedpoints=[points[0].index-fitspan,points[0].index] - + # Calls the function linefit_between parameters=[0,0,[],[]] try: @@ -259,31 +257,31 @@ class generalvclampCommands: except: print 'Cannot fit. Did you click twice the same point?' return - + # Outputs the relevant slope parameter print 'Slope:' print str(parameters[0]) to_dump='slope '+self.current.path+' '+str(parameters[0]) self.outlet.push(to_dump) - + # Makes a vector with the fitted parameters and sends it to the GUI xtoplot=parameters[2] ytoplot=[] x=0 for x in xtoplot: ytoplot.append((x*parameters[0])+parameters[1]) - + clickvector_x, clickvector_y=[], [] for item in points: clickvector_x.append(item.graph_coords[0]) clickvector_y.append(item.graph_coords[1]) lineplot=self._get_displayed_plot(0) #get topmost displayed plot - + lineplot.add_set(xtoplot,ytoplot) lineplot.add_set(clickvector_x, clickvector_y) - - + + if lineplot.styles==[]: lineplot.styles=[None,None,None,'scatter'] else: @@ -292,8 +290,8 @@ class generalvclampCommands: lineplot.styles=[None,None,None,None] else: lineplot.colors+=[None,None] - - + + self._send_plot([lineplot]) def linefit_between(self,index1,index2,whatset=1): @@ -310,12 +308,9 @@ class generalvclampCommands: # Translates the indexes into two vectors containing the x,y data to fit xtofit=self.plots[0].vectors[whatset][0][index1:index2] ytofit=self.plots[0].vectors[whatset][1][index1:index2] - + # Does the actual linear fitting (simple least squares with numpy.polyfit) linefit=[] linefit=np.polyfit(xtofit,ytofit,1) return (linefit[0],linefit[1],xtofit,ytofit) - - - diff --git a/hooke/plugin/macro.py b/hooke/plugin/macro.py index bdeafab..800d311 100644 --- a/hooke/plugin/macro.py +++ b/hooke/plugin/macro.py @@ -1,5 +1,3 @@ -#!/usr/bin/env python - ''' COMMAND MACRO PLUGIN FOR HOOKE @@ -89,7 +87,7 @@ class macroCommands: if args=='stop': self.pause=0 self.prompt=self.auxprompt - if len(self.currentmacro) != 0: + if len(self.currentmacro) != 0: answer=linput.safeinput('Do you want to save this macro? ',['y']) if answer[0].lower() == 'y': self.do_savemacro('') @@ -105,7 +103,7 @@ class macroCommands: self.pause=0 self.collect() else: - if len(self.currentmacro) != 0: + if len(self.currentmacro) != 0: answer=linput.safeinput('Another macro is already beign recorded\nDo you want to save it?',['y']) if answer[0].lower() == 'y': self.do_savemacro('') @@ -122,14 +120,14 @@ class macroCommands: Saves previously recorded macro into a script file for future use ------- Syntax: savemacro [macroname] - If no macroname is supplied one will be interactively asked + If no macroname is supplied one will be interactively asked ''' saved_ok=0 if self.currentmacro==None: print 'No macro is being recorded!' return 0 - if len(macroname)==0: + if len(macroname)==0: macroname=linput.safeinput('Enter new macro name: ') if len(macroname) == 0: print 'Invalid name' @@ -156,10 +154,10 @@ class macroCommands: macroname.hkm should be present at [hooke]/macros directory By default the macro will be executed over current curve - passing 'playlist' word as second argument executes macroname + passing 'playlist' word as second argument executes macroname over all curves - By default curve(s) will be processed silently, passing 'v' - as second/third argument will print each command that is + By default curve(s) will be processed silently, passing 'v' + as second/third argument will print each command that is executed Note that macros applied to playlists should end by export @@ -195,11 +193,11 @@ class macroCommands: print 'Could not find a macro named '+macropath return 0 txtfile=open(macropath) - if cycle ==1: + if cycle ==1: #print self.current_list for item in self.current_list: self.current=item - self.do_plot(0) + self.do_plot(0) for command in txtfile: @@ -225,9 +223,3 @@ class macroCommands: if verbose==1: print 'Executing command '+' '.join(testcmd) self.onecmd(' '.join(testcmd)) - - - - - - diff --git a/hooke/plugin/massanalysis.py b/hooke/plugin/massanalysis.py index e3d3830..6fcb03d 100644 --- a/hooke/plugin/massanalysis.py +++ b/hooke/plugin/massanalysis.py @@ -1,5 +1,3 @@ -#!/usr/bin/env python - ''' massanalysis.py @@ -21,10 +19,10 @@ import csv class massanalysisCommands: def _plug_init(self): - self.mass_variables={} + self.mass_variables={} self.interesting_variables=['curve','firstpeak_distance','lastpeak_distance','Npeaks','median_distance','mean_distance'] self._clean_data() - + def _clean_data(self): for variable in self.interesting_variables: self.mass_variables[variable]=[] @@ -33,18 +31,18 @@ class massanalysisCommands: ''' calculates X distance of a peak from the contact point ''' - item.identify(self.drivers) - + item.identify(self.drivers) + real_positions=[] cut_index=self.find_contact_point() - + #we assume the first is the plot with the force curve plot=item.curve.default_plots()[0] xret=plot.vectors[1][0] - + start_x=xret[cut_index] - - real_positions=[abs((xret[index])-(start_x)) for index in locations] + + real_positions=[abs((xret[index])-(start_x)) for index in locations] #close all open files item.curve.close_all() #needed to avoid *big* memory leaks! @@ -61,39 +59,39 @@ class massanalysisCommands: ''' self._clean_data() #if we recall it, clean previous data! min_deviation=self.convfilt_config['mindeviation'] - - + + c=0 for item in self.current_list: try: peak_location,peak_size=self.exec_has_peaks(item, min_deviation) real_positions=self.peak_position_from_contact(item, peak_location) - + self.mass_variables['Npeaks'].append(len(peak_location)) - + if len(peak_location) > 1: self.mass_variables['firstpeak_distance'].append(min(real_positions)) self.mass_variables['lastpeak_distance'].append(max(real_positions)) - + distancepeaks=[] for index in range(len(real_positions)-1): distancepeaks.append(real_positions[index+1]-real_positions[index]) else: self.mass_variables['firstpeak_distance'].append(0) - self.mass_variables['lastpeak_distance'].append(0) - + self.mass_variables['lastpeak_distance'].append(0) + if len(peak_location) > 2: self.mass_variables['median_distance'].append(np.median(distancepeaks)) - self.mass_variables['mean_distance'].append(np.mean(distancepeaks)) + self.mass_variables['mean_distance'].append(np.mean(distancepeaks)) else: self.mass_variables['median_distance'].append(0) - self.mass_variables['mean_distance'].append(0) - + self.mass_variables['mean_distance'].append(0) + print 'curve',c except SyntaxError: print 'curve',c,'not mapped' pass - + c+=1 def do_plotmap(self,args): @@ -107,37 +105,35 @@ class massanalysisCommands: print 'Give me two arguments between those:' print self.interesting_variables return - + scattermap=lhc.PlotObject() scattermap.vectors=[[]] scattermap.vectors[0].append(x) scattermap.vectors[0].append(y) - + scattermap.units=[args[0],args[1]] scattermap.styles=['scatter'] scattermap.destination=1 - + self._send_plot([scattermap]) - + def do_savemaps(self,args): ''' args=filename ''' - + ''' def csv_write_cols(data, f): - + #from Bruno Desthuillers on comp.lang.python - + writer = csv.writer(f) keys = data.keys() writer.writerow(dict(zip(keys,keys))) for row in zip(*data.values()): - writer.writerow(dict(zip(keys, row))) + writer.writerow(dict(zip(keys, row))) ''' - + f=open(args,'wb') lh.csv_write_dictionary(f,self.mass_variables) f.close() - - \ No newline at end of file diff --git a/hooke/plugin/multidistance.py b/hooke/plugin/multidistance.py index c93900f..adda783 100644 --- a/hooke/plugin/multidistance.py +++ b/hooke/plugin/multidistance.py @@ -14,7 +14,7 @@ warnings.simplefilter('ignore',np.RankWarning) class multidistanceCommands: - + def do_multidistance(self,args): ''' MULTIDISTANCE @@ -25,7 +25,7 @@ class multidistanceCommands: of an existing file, autopeak will resume it and append measurements to it. If you are giving a new filename, it will create the file and append to it until you close Hooke. You can also define a minimun deviation of the peaks. - + Syntax: multidistance [deviation] deviation = number of times the convolution signal is above the noise absolute deviation. @@ -48,20 +48,20 @@ class multidistanceCommands: pk_loc,peak_size=self.has_peaks(defplot, abs_devs) return pk_loc, peak_size - + noflatten=False peaks_location, peak_size=find_current_peaks(noflatten, args) - + #if no peaks, we have nothing to plot. exit. if len(peaks_location)==0: return - + #otherwise, we plot the peak locations. xplotted_ret=self.plots[0].vectors[1][0] yplotted_ret=self.plots[0].vectors[1][1] xgood=[xplotted_ret[index] for index in peaks_location] ygood=[yplotted_ret[index] for index in peaks_location] - + recplot=self._get_displayed_plot() recplot.vectors.append([xgood,ygood]) if recplot.styles==[]: @@ -79,7 +79,7 @@ class multidistanceCommands: if exclude_raw=='N': print 'Discarded.' return - + if not exclude_raw=='': exclude=exclude_raw.split(',') #we convert in numbers the input @@ -100,13 +100,13 @@ class multidistanceCommands: peaks_location= peaks_location[0:new_a]+peaks_location[new_a+1:] peak_size= peak_size[0:new_a]+peak_size[new_a+1:] count+=1 - + #we calculate the distance vector dist=[] for i in range(len(peaks_location)-1): dist.append(xplotted_ret[peaks_location[i]]-xplotted_ret[peaks_location[i+1]]) - - + + @@ -116,7 +116,7 @@ class multidistanceCommands: if self.autofile=='': print 'Not saved.' return - + if not os.path.exists(self.autofile): f=open(self.autofile,'w+') f.write('Analysis started '+time.asctime()+'\n') @@ -128,14 +128,14 @@ class multidistanceCommands: f.write(str(o)) f.write("\n") f.close() - + print 'Saving...' f=open(self.autofile,'a+') - + f.write(self.current.path+'\n') for i in dist: f.write(";") f.write(str(i)) - f.write("\n") + f.write("\n") f.close() diff --git a/hooke/plugin/pcluster.py b/hooke/plugin/pcluster.py index 5084b83..5f51939 100644 --- a/hooke/plugin/pcluster.py +++ b/hooke/plugin/pcluster.py @@ -1,5 +1,3 @@ -#!/usr/bin/env python - from mdp import pca from libhooke import WX_GOOD, ClickedPoint import wxversion @@ -18,7 +16,7 @@ warnings.simplefilter('ignore',np.RankWarning) class pclusterCommands: - + def _plug_init(self): self.clustplot1=None self.clustplot2=None @@ -39,7 +37,7 @@ class pclusterCommands: self.my_work_dir = os.getcwd()+slash+pclus_dir+slash self.my_curr_dir = os.path.basename(os.getcwd()) os.mkdir(self.my_work_dir) - + #--Custom persistent length pl_value=None for arg in args.split(): @@ -49,33 +47,33 @@ class pclusterCommands: pl_value=float(pl_expression[1]) #actual value else: pl_value=None - + #configuration variables min_npks = self.convfilt_config['minpeaks'] min_deviation = self.convfilt_config['mindeviation'] - + pclust_filename = "automeasure_"+self.my_curr_dir+"_blind"+blindw+".txt" #raw_input('Automeasure filename? ') realclust_filename = "coordinate_"+self.my_curr_dir+"_blind"+blindw+".txt" #raw_input('Coordinates filename? ') peackforce_filename = "peakforce_"+self.my_curr_dir+"_blind"+blindw+".txt" #raw_input('Peacks and Forces filename? ') - + f=open(self.my_work_dir+pclust_filename,'w+') f.write('Analysis started '+time.asctime()+'\n') f.write('----------------------------------------\n') f.write('; Contour length (nm) ; Persistence length (nm) ; Max.Force (pN) ; Slope (N/m) ; Sigma contour (nm) ; Sigma persistence (nm)\n') f.close() - + f=open(self.my_work_dir+realclust_filename,'w+') f.write('Analysis started '+time.asctime()+'\n') f.write('----------------------------------------\n') f.write('; Peak number ; Mean delta (nm) ; Median delta (nm) ; Mean force (pN) ; Median force (pN) ; First peak length (nm) ; Last peak length (nm) ; Max force (pN) ; Min force (pN) ; Max delta (nm) ; Min delta (nm) ; Peaks Diff\n') f.close() - + f=open(self.my_work_dir+peackforce_filename,'w+') f.write('Analysis started '+time.asctime()+'\n') f.write('----------------------------------------\n') f.write('; Peak number ; 1 peak Length (nm) ; 1 peak Force (pN) ; 2 peak Length (nm) ; 2 peak Force (pN) ; 3 peak Length (nm) ; 3 peak Force (pN) ; 4 peak Length (nm) ; 4 peak Force (pN) ; 5 peak Length (nm) ; 5 peak Force (pN) ; 6 peak Length (nm) ; 6 peak Force (pN) ; 7 peak Length (nm) ; 7 peak Force (pN) ; 8 peak Length (nm) ; 8 peak Force (pN)\n') f.close() - + # ------ FUNCTION ------ def fit_interval_nm(start_index,plot,nm,backwards): ''' @@ -86,7 +84,7 @@ class pclusterCommands: ''' whatset=1 #FIXME: should be decidable x_vect=plot.vectors[1][0] - + c=0 i=start_index start=x_vect[start_index] @@ -113,7 +111,7 @@ class pclusterCommands: fit_errors [6.5817195369767644e-010, 2.4415923138871498e-011] ''' fit_points=int(self.config['auto_fit_points']) # number of points to fit before the peak maximum <50> - + T=self.config['temperature'] #temperature of the system in kelvins. By default it is 293 K. <301.0> cindex=self.find_contact_point(itplot[0]) #Automatically find contact point <158, libhooke.ClickedPoint> contact_point=self._clickize(itplot[0].vectors[1][0], itplot[0].vectors[1][1], cindex) @@ -131,7 +129,7 @@ class pclusterCommands: def features_peaks(itplot, peak, fit_points, contact_point, pl_value, T, cindex, avg): ''' - calculate informations for each peak and add they in + calculate informations for each peak and add they in c_lengths, p_lengths, sigma_c_lengths, sigma_p_lengths, forces, slopes ''' c_leng=None @@ -140,17 +138,17 @@ class pclusterCommands: sigma_p_leng=None force=None slope=None - + delta_force=10 slope_span=int(self.config['auto_slope_span']) - + peak_point=self._clickize(itplot[0].vectors[1][0],itplot[0].vectors[1][1],peak) other_fit_point=self._clickize(itplot[0].vectors[1][0],itplot[0].vectors[1][1],peak-fit_points) - + points=[contact_point, peak_point, other_fit_point] - + params, yfit, xfit, fit_errors = self.wlc_fit(points, itplot[0].vectors[1][0], itplot[0].vectors[1][1], pl_value, T, return_errors=True) - + #Measure forces delta_to_measure=itplot[0].vectors[1][1][peak-delta_force:peak+delta_force] y=min(delta_to_measure) @@ -168,12 +166,12 @@ class pclusterCommands: #check if persistent length makes sense. otherwise, discard peak. if p_leng>self.config['auto_min_p'] and p_leng "pca gaeta_coor_blind50.txt 1,3,6" Automatically measures pca from coordinates filename and shows two interactives plots - With the second argument (arbitrary) you can select the columns and the multiplier factor + With the second argument (arbitrary) you can select the columns and the multiplier factor to use for the pca (for es "1,3*50,6,8x10,9"). Dont use spaces. "*" or "x" are the same thing. Without second argument it reads pca_config.txt file (c)Paolo Pancaldi, Massimo Sandal 2009 ''' - + # reads the columns of pca if self.config['hookedir'][0]=='/': slash='/' #a Unix or Unix-like system @@ -355,7 +353,7 @@ class pclusterCommands: conf=open(self.my_hooke_dir+"pca_config.txt") config = conf.readlines() conf.close() - + self.plot_myCoord = [] # tiene le coordinate prese direttamente dal file creato con pCluster self.plot_origCoord = [] # tiene le coordinate solo delle colonne scelte e moltiplicate per i valori scelti self.plot_pcaCoord = [] # tiene le due colonne della pca @@ -366,8 +364,8 @@ class pclusterCommands: self.plot_NewPcaCoord = [] # tiene le due colonne della pca filtrate per densita self.plot_NewPcaCoordTr=[] # tiene le due colonne della pca trasposta filtrate per densita plot_path_temp = "" - - # prende in inpunt un arg (nome del file) + + # prende in inpunt un arg (nome del file) # e il secondo le colonne su cui lavorare (e' arbitrario, riceve x es "1,2,3") arg = args.split(" ") if arg[0]==args: @@ -375,7 +373,7 @@ class pclusterCommands: else: file_name=arg[0] config[0] = arg[1] - + # creo l'array "plot_myCoord" con tutte le coordinate dei plots # e l'array plot_paths con tutti i percorsi dei plots nPlotTot = -3 #tolgo le prime 3 righe iniziali del file @@ -388,7 +386,7 @@ class pclusterCommands: if row[0]==" " and row.find('nan')==-1 and row.find("-1.#IND")==-1: row = row[row.index(";",2)+2:].split(" ; ") # non considero la prima colonna col #picchi row = [float(i) for i in row] - + #0:Mean delta, 1:Median delta, 2:Mean force, 3:Median force, 4:First peak length, 5:Last peak length #6:Max delta 7:Min delta 8:Max force 9:Min force 10:Std delta 11:Std force if (row[0]<500 and row[1]<500 and row[2]<500 and row[3]<500 and row[4]<500 and row[5]<500 and row[6]<500 and row[7]<500 and row[8]<500 and row[9]<500 and row[10]<500 and row[11]<500): @@ -397,8 +395,8 @@ class pclusterCommands: self.plot_myCoord.append(row) self.plot_paths.append(plot_path_temp) f.close() - - # creo l'array con alcune colonne e pure moltiplicate + + # creo l'array con alcune colonne e pure moltiplicate for row in self.plot_myCoord: res=[] for cols in config[0].split(","): @@ -413,7 +411,7 @@ class pclusterCommands: molt = 1 res.append(row[col]*molt) self.plot_origCoord.append(res) - + # array convert, calculate PCA, transpose self.plot_origCoord = np.array(self.plot_origCoord,dtype='float') #print self.plot_origCoord.shape @@ -421,11 +419,11 @@ class pclusterCommands: self.plot_pcaCoordTr = np.transpose(self.plot_pcaCoord) pca_X=np.array(self.plot_pcaCoordTr[0],dtype='float') pca_Y=np.array(self.plot_pcaCoordTr[1],dtype='float') - + ''' # Start section of testing with good plots # 4 TESTING! Xsyn_1=[] - Ysyn_1=[] + Ysyn_1=[] Xgb1_1=[] Ygb1_1=[] Xbad_1=[] @@ -434,7 +432,7 @@ class pclusterCommands: goodnames=goodnamefile.readlines() nPlotGood = len(goodnames)-2 #tolgo prima e ultima riga goodnames=[i.split()[0] for i in goodnames[1:]] - + for index in range(len(self.plot_paths)): if self.plot_paths[index][:-1] in goodnames: Xsyn_1.append(pca_X[index]) @@ -444,7 +442,7 @@ class pclusterCommands: Ybad_1.append(pca_Y[index]) # Stop section of testing with good plots # 4 TESTING! ''' - + # print first plot clustplot1=lhc.PlotObject() clustplot1.add_set(pca_X,pca_Y) @@ -456,7 +454,7 @@ class pclusterCommands: clustplot1.destination=0 self._send_plot([clustplot1]) self.clustplot1=clustplot1 - + # density and filer estimation kernel = sp.stats.kde.gaussian_kde(sp.c_[pca_X,pca_Y].T) tallest = 0 @@ -479,7 +477,7 @@ class pclusterCommands: axis_X = np.linspace(xmin,xmax,num=100) axis_Y = np.linspace(ymin,ymax,num=100) ''' - + # density filtering: # tramite "kernel.evaluate" trovo lo score (altezza) di ogni coordinata e decido se mantenerla o no filtered_pca_X = [] @@ -494,16 +492,16 @@ class pclusterCommands: filtered_PcaCoordTr.append(filtered_pca_X) filtered_PcaCoordTr.append(filtered_pca_Y) filtered_PcaCoord = np.transpose(filtered_PcaCoordTr) - + # creo i due array "plot_FiltOrigCoord" e "plot_FiltPaths" contenenti solo i dati filtrati con alta densita for index in range(len(self.plot_pcaCoord)): if self.plot_pcaCoord[index] in filtered_PcaCoord: self.plot_FiltOrigCoord.append(self.plot_myCoord[index]) self.plot_FiltPaths.append(self.plot_paths[index]) - + ''' # START PCA#2: USELESS!!! - + # creo l array con alcune colonne e pure moltiplicate temp_coord = [] for row in self.plot_FiltOrigCoord: @@ -521,7 +519,7 @@ class pclusterCommands: res.append(row[col]*molt) temp_coord.append(res) self.plot_FiltOrigCoord = temp_coord - + # ricalcolo la PCA: array convert, calculate PCA, transpose self.plot_FiltOrigCoord = np.array(self.plot_FiltOrigCoord,dtype='float') #print self.plot_FiltOrigCoord.shape @@ -529,7 +527,7 @@ class pclusterCommands: self.plot_NewPcaCoordTr = np.transpose(self.plot_NewPcaCoord) pca_X2=np.array(self.plot_NewPcaCoordTr[0],dtype='float') pca_Y2=np.array(self.plot_NewPcaCoordTr[1],dtype='float') - + # Start section of testing with good plots # 4 TESTING! Xsyn_2=[] Ysyn_2=[] @@ -542,7 +540,7 @@ class pclusterCommands: else: Xbad_2.append(pca_X2[index]) Ybad_2.append(pca_Y2[index]) - + # print second plot clustplot2=lhc.PlotObject() #clustplot2.add_set(pca_X2,pca_Y2) @@ -555,7 +553,7 @@ class pclusterCommands: self._send_plot([clustplot2]) self.clustplot2=clustplot2 ''' - + # PRINT density plot clustplot2=lhc.PlotObject() clustplot2.add_set(filtered_pca_X,filtered_pca_Y) @@ -565,7 +563,7 @@ class pclusterCommands: clustplot2.destination=1 self._send_plot([clustplot2]) self.clustplot2=clustplot2 - + # printing results config_pca1 = config[0].replace("*", "x").rstrip("\n") config_pca2 = config[2].replace("*", "x").rstrip("\n") @@ -582,7 +580,7 @@ class pclusterCommands: print "Piu alta: ", tallest #print "- FILTRO 3: 2'PCA:"+config_pca2 print "" - + # -- exporting coordinates and plot of PCA in debug mode! -- if config[3].find("true")!=-1: #1' PCA: save plot and build coordinate s file @@ -610,7 +608,7 @@ class pclusterCommands: f.write (" ; " + str(cel)) f.write ("\n") f.close() - + # pCLUSTER SAVING!!! import shutil pcl_name = file_name.replace("coordinate_", "goodplots_").replace('.txt','_'+config_pca1+'_'+str(my_filter).replace(".",",")) @@ -622,12 +620,12 @@ class pclusterCommands: f.write (myfile+"\n") shutil.copy2(myfile, pcl_name) f.close() - - + + def do_multipca(self,args): ''' MULTIPCA -> "multipca gaeta_coor_blind50.txt 3" - Automatically multiply the column suggest in second argument for value between 1-100 (step of 2), + Automatically multiply the column suggest in second argument for value between 1-100 (step of 2), measures pca from coordinates filename and save the png plots. (c)Paolo Pancaldi, Massimo Sandal 2009 ''' @@ -655,7 +653,7 @@ class pclusterCommands: for j in range(1, 13): if i!=j: self.do_pca(file_name + " " + str(i) + "," + str(j)) - + def do_triplepca(self,args): ''' TRIPLEPCA -> "triplepca gaeta_coor_blind50.txt" @@ -675,7 +673,7 @@ class pclusterCommands: ''' It returns id, coordinates and file name of a clicked dot on a PCA graphic ''' - + self._send_plot([self.clustplot1]) #quick workaround for BAD problems in the GUI print 'Click point' point = self._measure_N_points(N=1, whatset=0) @@ -687,7 +685,7 @@ class pclusterCommands: print "coord: " + str(dot_coord) self.do_genlist(str(plot_file)) #self.do_jump(str(plot_file)) - + # indea iniziata e messa da parte... def do_peakforce(self, args): ''' @@ -695,14 +693,14 @@ class pclusterCommands: Automatically measures peack and force plots (c)Paolo Pancaldi, Massimo Sandal 2009 ''' - + # prende in inpunt un arg (nome del file) file_name=args f=open(file_name) - + # scrivo un file temp g = open('_prove.txt','w') - + plot_file = ''; rows = f.readlines() for row in rows: @@ -731,5 +729,3 @@ class pclusterCommands: g.write(writeme) g.close() f.close() - - \ No newline at end of file diff --git a/hooke/plugin/procplots.py b/hooke/plugin/procplots.py index f1cc0a3..68e7f8f 100755 --- a/hooke/plugin/procplots.py +++ b/hooke/plugin/procplots.py @@ -1,5 +1,3 @@ -#!/usr/bin/env python - ''' PROCPLOTS Processed plots plugin for force curves. @@ -18,10 +16,10 @@ import scipy.signal import copy class procplotsCommands: - + def _plug_init(self): pass - + def do_derivplot(self,args): ''' DERIVPLOT @@ -29,30 +27,30 @@ class procplotsCommands: Plots the derivate (actually, the discrete differentiation) of the current force curve --------- Syntax: derivplot - ''' - dplot=self.derivplot_curves() - plot_graph=self.list_of_events['plot_graph'] + ''' + dplot=self.derivplot_curves() + plot_graph=self.list_of_events['plot_graph'] wx.PostEvent(self.frame,plot_graph(plots=[dplot])) - + def derivplot_curves(self): ''' do_derivplot helper function - + create derivate plot curves for force curves. - ''' + ''' dplot=lhc.PlotObject() - dplot.vectors=[] - + dplot.vectors=[] + for vector in self.plots[0].vectors: dplot.vectors.append([]) dplot.vectors[-1].append(vector[0][:-1]) dplot.vectors[-1].append(np.diff(vector[1])) - + dplot.destination=1 dplot.units=self.plots[0].units - - return dplot - + + return dplot + def do_subtplot(self,args): ''' SUBTPLOT @@ -62,42 +60,42 @@ class procplotsCommands: Syntax: subtplot ''' #FIXME: sub_filter and sub_order must be args - + if len(self.plots[0].vectors) != 2: print 'This command only works on a curve with two different plots.' pass - + outplot=self.subtract_curves(sub_order=1) - - plot_graph=self.list_of_events['plot_graph'] + + plot_graph=self.list_of_events['plot_graph'] wx.PostEvent(self.frame,plot_graph(plots=[outplot])) - + def subtract_curves(self, sub_order): ''' subtracts the extension from the retraction - ''' + ''' xext=self.plots[0].vectors[0][0] yext=self.plots[0].vectors[0][1] xret=self.plots[0].vectors[1][0] yret=self.plots[0].vectors[1][1] - + #we want the same number of points maxpoints_tot=min(len(xext),len(xret)) xext=xext[0:maxpoints_tot] yext=yext[0:maxpoints_tot] xret=xret[0:maxpoints_tot] yret=yret[0:maxpoints_tot] - + if sub_order: ydiff=[yretval-yextval for yretval,yextval in zip(yret,yext)] else: #reverse subtraction (not sure it's useful, but...) - ydiff=[yextval-yretval for yextval,yretval in zip(yext,yret)] - + ydiff=[yextval-yretval for yextval,yretval in zip(yext,yret)] + outplot=copy.deepcopy(self.plots[0]) outplot.vectors[0][0], outplot.vectors[1][0] = xext,xret #FIXME: if I use xret, it is not correct! outplot.vectors[1][1]=ydiff outplot.vectors[0][1]=[0 for item in outplot.vectors[1][0]] - + return outplot @@ -110,21 +108,21 @@ class procplotsCommands: median_filter=customvalue else: median_filter=self.config['medfilt'] - + if median_filter==0: return plot - + if float(median_filter)/2 == int(median_filter)/2: median_filter+=1 - + nplots=len(plot.vectors) c=0 while cstartvalue: cutindex=index else: @@ -174,12 +172,12 @@ class procplotsCommands: plot.vectors[0][0]=plot.vectors[0][0][:cutindex] plot.vectors[0][1]=plot.vectors[0][1][:cutindex] - + return plot - ''' - - - + ''' + + + #FFT--------------------------- def fft_plot(self, vector): ''' @@ -187,22 +185,22 @@ class procplotsCommands: ''' fftplot=lhc.PlotObject() fftplot.vectors=[[]] - + fftlen=len(vector)/2 #need just 1/2 of length - fftplot.vectors[-1].append(np.arange(1,fftlen).tolist()) - + fftplot.vectors[-1].append(np.arange(1,fftlen).tolist()) + try: fftplot.vectors[-1].append(abs(np.fft(vector)[1:fftlen]).tolist()) except TypeError: #we take care of newer NumPy (1.0.x) fftplot.vectors[-1].append(abs(np.fft.fft(vector)[1:fftlen]).tolist()) - - + + fftplot.destination=1 - - + + return fftplot - - + + def do_fft(self,args): ''' FFT @@ -210,15 +208,15 @@ class procplotsCommands: Plots the fast Fourier transform of the selected plot --- Syntax: fft [top,bottom] [select] [0,1...] - + By default, fft performs the Fourier transform on all the 0-th data set on the top plot. - + [top,bottom]: which plot is the data set to fft (default: top) [select]: you pick up two points on the plot and fft only the segment between [0,1,...]: which data set on the selected plot you want to fft (default: 0) ''' - + #args parsing #whatplot = plot to fft #whatset = set to fft in the plot @@ -235,7 +233,7 @@ class procplotsCommands: whatset=int(arg) except ValueError: pass - + if select: points=self._measure_N_points(N=2, whatset=whatset) boundaries=[points[0].index, points[1].index] @@ -243,9 +241,8 @@ class procplotsCommands: y_to_fft=self.plots[whatplot].vectors[whatset][1][boundaries[0]:boundaries[1]] #y else: y_to_fft=self.plots[whatplot].vectors[whatset][1] #y - - fftplot=self.fft_plot(y_to_fft) + + fftplot=self.fft_plot(y_to_fft) fftplot.units=['frequency', 'power'] - plot_graph=self.list_of_events['plot_graph'] + plot_graph=self.list_of_events['plot_graph'] wx.PostEvent(self.frame,plot_graph(plots=[fftplot])) - \ No newline at end of file diff --git a/hooke/plugin/superimpose.py b/hooke/plugin/superimpose.py index fdd0a33..ce10e99 100644 --- a/hooke/plugin/superimpose.py +++ b/hooke/plugin/superimpose.py @@ -1,5 +1,3 @@ -#!/usr/bin/env python - from libhooke import WX_GOOD import wxversion wxversion.select(WX_GOOD) @@ -9,38 +7,38 @@ import libhookecurve as lhc from numpy import arange, mean class superimposeCommands: - + def _plug_init(self): self.imposed=[] - + def do_selimpose(self,args): ''' SELIMPOSE (superimpose.py plugin) Hand-selects the curve portion to superimpose ''' #fixme: set to superimpose should be in args - + if args=='clear': self.imposed=[] return - + current_set=1 - + points=self._measure_two_points() boundaries=[points[0].index, points[1].index] boundaries.sort() - + theplot=self.plots[0] #append the selected section self.imposed.append([]) self.imposed[-1].append(theplot.vectors[1][0][boundaries[0]:boundaries[1]]) #x self.imposed[-1].append(theplot.vectors[1][1][boundaries[0]:boundaries[1]]) #y - + #align X first point self.imposed[-1][0] = [item-self.imposed[-1][0][0] for item in self.imposed[-1][0]] #align Y first point self.imposed[-1][1] = [item-self.imposed[-1][1][0] for item in self.imposed[-1][1]] - + def do_plotimpose(self,args): ''' PLOTIMPOSE (sumperimpose.py plugin) @@ -49,16 +47,16 @@ class superimposeCommands: imposed_object=lhc.PlotObject() imposed_object.vectors=self.imposed print 'Plotting',len(imposed_object.vectors),'imposed curves' - + imposed_object.normalize_vectors() - + imposed_object.units=self.plots[0].units imposed_object.title='Imposed curves' imposed_object.destination=1 - - plot_graph=self.list_of_events['plot_graph'] + + plot_graph=self.list_of_events['plot_graph'] PostEvent(self.frame,plot_graph(plots=[imposed_object])) - + def do_plotavgimpose(self,args): ''' PLOTAVGIMPOSE (superimpose.py plugin) @@ -69,10 +67,10 @@ class superimposeCommands: min_x=[] for curve in self.imposed: min_x.append(min(curve[0])) - + #find minimum extension min_ext_limit=max(min_x) - + x_avg=arange(step,min_ext_limit,step) y_avg=[] for value in x_avg: @@ -80,19 +78,19 @@ class superimposeCommands: for curve in self.imposed: for xvalue, yvalue in zip(curve[0],curve[1]): if xvalue >= (value+step) and xvalue <= (value-step): - to_avg.append(yvalue) + to_avg.append(yvalue) y_avg.append(mean(to_avg)) - + print 'len x len y' print len(x_avg), len(y_avg) print y_avg - + avg_object=lhc.PlotObject() avg_object.vectors=[[x_avg, y_avg]] avg_object.normalize_vectors() avg_object.units=self.plots[0].units avg_object.title="Average curve" avg_object.destination=1 - - plot_graph=self.list_of_events['plot_graph'] - PostEvent(self.frame,plot_graph(plots=[avg_object])) \ No newline at end of file + + plot_graph=self.list_of_events['plot_graph'] + PostEvent(self.frame,plot_graph(plots=[avg_object])) diff --git a/hooke/plugin/tutorial.py b/hooke/plugin/tutorial.py index f9f2f33..ad7a6fe 100755 --- a/hooke/plugin/tutorial.py +++ b/hooke/plugin/tutorial.py @@ -1,5 +1,3 @@ -#!/usr/bin/env python - ''' TUTORIAL PLUGIN FOR HOOKE @@ -25,14 +23,14 @@ class tutorialCommands: ''' Here we define the class containing all the Hooke commands we want to define in the plugin. - + Notice the class name!! The syntax is filenameCommands. That is, if your plugin is pluggy.py, your class name is pluggyCommands. - + Otherwise, the class will be ignored by Hooke. - ''' - + ''' + def _plug_init(self): ''' This is the plugin initialization. @@ -40,57 +38,57 @@ class tutorialCommands: If there is something you need to do when Hooke starts, code it in this function. ''' print 'I am the Tutorial plugin initialization!' - + #Here we initialize a local configuration variable; see plotmanip_absvalue() for explanation. - self.config['tutorial_absvalue']=0 + self.config['tutorial_absvalue']=0 pass - + def do_nothing(self,args): ''' This is a boring but working example of an actual Hooke command. A Hooke command is a function of the xxxxCommands class, which is ALWAYS defined this way: - + def do_nameofcommand(self,args) - + *do_ is needed to make Hooke understand this function is a command *nameofcommand is how the command will be called in the Hooke command line. *self is, well, self *args is ALWAYS needed (otherwise Hooke will crash executing the command). We will see later what args is. - + Note that if you now start Hooke with this plugin activated and you type in the Hooke command line "help nothing" you will see this very text as output. So the help of a command is a string comment below the function definition, like this one. - + Commands usually return None. ''' print 'I am a Hooke command. I do nothing.' - + def do_printargs(self,args): ''' This command prints the args you give to it. args is always a string, that contains everything you write after the command. So if you issue "mycommand blah blah 12345" args is "blah blah 12345". - + Again, args is needed in the definition even if your command does not use it. ''' print 'You gave me those args: '+args - + def help_tutorial(self): ''' - This is a help function. + This is a help function. If you want a help function for something that is not a command, you can write a help function like this. Calling "help tutorial" will execute this function. ''' print 'You called help_tutorial()' - + def do_environment(self,args): ''' This plugin contains a panoramic of the Hooke command line environment variables, and prints their current value. ''' - + '''self.current_list TYPE: [ libhookecurve.HookeCurve ], len=variable contains the actual playlist of Hooke curve objects. @@ -99,38 +97,38 @@ class tutorialCommands: ''' print 'current_list length:',len(self.current_list) print 'current_list 0th:',self.current_list[0] - + '''self.pointer TYPE: int contains the index of the current curve in the playlist ''' print 'pointer: ',self.pointer - + '''self.current TYPE: libhookecurve.HookeCurve contains the current curve displayed. We will see later how it works. ''' print 'current:',self.current - + '''self.plots TYPE: [ libhookecurve.PlotObject ], len=1,2 contains the current default plots. - Each PlotObject contains all info needed to display + Each PlotObject contains all info needed to display the plot: apart from the data vectors, the title, destination etc. Usually self.plots[0] is the default topmost plot, self.plots[1] is the accessory bottom plot. ''' print 'plots:',self.plots - + '''self.config TYPE: { string:anything } contains the current Hooke configuration variables, in form of a dictionary. ''' print 'config:',self.config - + '''self.plotmanip TYPE: [ function ] Contains the ordered plot manipulation functions. @@ -140,21 +138,21 @@ class tutorialCommands: *YOU SHOULD NEVER MODIFY THAT* ''' print 'plotmanip: ',self.plotmanip - + '''self.drivers TYPE: [ class ] Contains the plot reading drivers. *YOU SHOULD NEVER MODIFY THAT* ''' print 'drivers: ',self.drivers - + '''self.frame TYPE: wx.Frame Contains the wx Frame of the GUI. ***NEVER, EVER TOUCH THAT.*** ''' print 'frame: ',self.frame - + '''self.list_of_events TYPE: { string:wx.Event } Contains the wx.Events to communicate with the GUI. @@ -162,7 +160,7 @@ class tutorialCommands: to create a GUI plugin. ''' print 'list of events:',self.list_of_events - + '''self.events_from_gui TYPE: Queue.Queue Contains the Queue where data from the GUI is put. @@ -170,50 +168,50 @@ class tutorialCommands: to create a GUI plugin. ''' print 'events from gui:',self.events_from_gui - + '''self.playlist_saved TYPE: Int (0/1) ; Boolean Flag that tells if the playlist has been saved or not. ''' print 'playlist saved:',self.playlist_saved - + '''self.playlist_name TYPE: string Name of current playlist ''' print 'playlist name:',self.playlist_name - + '''self.notes_saved TYPE: Int (0/1) ; Boolean Flag that tells if the playlist has been saved or not. ''' print 'notes saved:',self.notes_saved - + def do_myfirstplot(self,args): ''' In this function, we see how to create a PlotObject and send it to the screen. ***Read the code of PlotObject in libhookecurve.py before!***. ''' - + #We generate some interesting data to plot for this example. xdata1=np.arange(-5,5,0.1) xdata2=np.arange(-5,5,0.1) ydata1=[item**2 for item in xdata1] ydata2=[item**3 for item in xdata2] - + #Create the object. #The PlotObject class lives in the libhookecurve library. myplot=lhc.PlotObject() ''' - The *data* of the plot live in the .vectors list. - + The *data* of the plot live in the .vectors list. + plot.vectors is a multidimensional array: plot.vectors[0]=set1 plot.vectors[1]=set2 plot.vectors[2]=sett3 etc. - + 2 curves in a x,y plot are: [[[x1],[y1]],[[x2],[y2]]] for example: @@ -224,28 +222,28 @@ class tutorialCommands: x2 = self.vectors[1][0] y2 = self.vectors[1][1] ''' - #Pour 0-th dataset into myplot: + #Pour 0-th dataset into myplot: myplot.add_set(xdata1,ydata1) - - #Pour 1-st dataset into myplot: + + #Pour 1-st dataset into myplot: myplot.add_set(xdata2,ydata2) - + #Add units to x and y axes #units=[string, string] myplot.units=['x axis','y axis'] - + #Where do we want the plot? 0=top, 1=bottom myplot.destination=1 - + '''Send it to the GUI. Note that you *have* to encapsulate it into a list, so you have to send [myplot], not simply myplot. - + You can also send more two plots at once self.send_plot([plot1,plot2]) ''' self._send_plot([myplot]) - + def do_myfirstscatter(self,args): ''' @@ -256,40 +254,40 @@ class tutorialCommands: xdata2=np.arange(-5,5,1) ydata1=[item**2 for item in xdata1] ydata2=[item**3 for item in xdata2] - + myplot=lhc.PlotObject() myplot.add_set(xdata1,ydata1) myplot.add_set(xdata2,ydata2) - - + + #Add units to x and y axes myplot.units=['x axis','y axis'] - + #Where do we want the plot? 0=top, 1=bottom myplot.destination=1 - + '''None=standard line plot 'scatter'=scatter plot By default, the styles attribute is an empty list. If you want to define a scatter plot, you must define all other plots as None or 'scatter', depending on what you want. - + Here we define the second set to be plotted as scatter, and the first to be plotted as line. - + Here we define also the colors to be the default Matplotlib colors ''' myplot.styles=[None,'scatter'] myplot.colors=[None,None] self._send_plot([myplot]) - + def do_clickaround(self,args): ''' Here we click two points on the curve and take some parameters from the points we have clicked. ''' - + ''' points = self._measure_N_points(N=Int, whatset=Int) *N = number of points to measure(1...n) @@ -298,17 +296,17 @@ class tutorialCommands: ''' points=self._measure_N_points(N=2,whatset=1) print 'You clicked the following points.' - + ''' These are the absolute coordinates of the - point clicked. + point clicked. [float, float] = x,y ''' print 'Absolute coordinates:' print points[0].absolute_coords print points[1].absolute_coords print - + ''' These are the coordinates of the points clicked, remapped on the graph. @@ -322,7 +320,7 @@ class tutorialCommands: print points[0].graph_coords print points[1].graph_coords print - + ''' These are the indexes of the clicked points on the dataset vector. @@ -330,39 +328,39 @@ class tutorialCommands: print 'Index of points on the graph:' print points[0].index print points[1].index - - + + def help_thedifferentplots(self): ''' The *three* different default plots you should be familiar with in Hooke. - + Each plot contains of course the respective data in their vectors attribute, so here you learn also which data access for each situation. ''' print ''' 1. THE RAW, CURRENT PLOTS - + self.current --- Contains the current libhookecurve.HookeCurve container object. A HookeCurve object defines only two default attributes: - + * self.current.path = string The path of the current displayed curve - + * self.current.curve = libhookecurve.Driver The curve object. This is not only generated by the driver, this IS a driver instance in itself. This means that in self.current.curve you can access the specific driver APIs, if you know them. - + And defines only one method: * self.current.identify() Fills in the self.current.curve object. See in the cycling tutorial. - + ***** The REAL curve data actually lives in: --- @@ -370,75 +368,75 @@ class tutorialCommands: Contains the raw PlotObject-s, as "spitted out" by the driver, without any intervention. This is as close to the raw data as Hooke gets. - + One or two plots can be spit out; they are always enclosed in a list. ***** - + Methods of self.current.curve are: --- - + * self.current.curve.is_me() (Used by identify() only.) - + * self.current.curve.close_all() Closes all driver open files; see the cycling tutorial. ''' - + print ''' 2. THE PROCESSED, DEFAULT PLOT - + The plot that is spitted out by the driver is *not* the usual default plot that is displayed by calling "plot" at the Hooke prompt. - + This is because the raw, driver-generated plot is usually *processed* by so called *plot processing* functions. We will see in the tutorial how to define - them. - + them. + For example, in force spectroscopy force curves, raw data are automatically corrected for deflection. Other data can be, say, filtered by default. - - The default plots are accessible in + + The default plots are accessible in self.plots = [ libhooke.PlotObject ] - + self.plots[0] is usually the topmost plot self.plots[1] is usually the bottom plot (if present) ''' - + print ''' 3. THE PLOT DISPLAYED RIGHT NOW. - + Sometimes the plots you are displaying *right now* is different from the previous two. You may have a fit trace, you may have issued some command that spits out - a custom plot and you want to rework that, whatever. - + a custom plot and you want to rework that, whatever. + You can obtain in any moment the plot currently displayed by Hooke by issuing - + PlotObject = self._get_displayed_plot(dest) * dest = Int (0/1) dest=0 : top plot dest=1 : bottom plot ''' - - + + def do_cycling(self,args): ''' Here we cycle through our playlist and print some info on the curves we find. Cycling through the playlist needs a bit of care to avoid memory leaks and dangling open files... - + Look at the source code for more information. ''' - + def things_when_cycling(item): ''' We encapsulate here everything has to open the actual curve file. By doing it all here, we avoid to do acrobacies when deleting objects etc. in the main loop: we do the dirty stuff here. ''' - + ''' identify() - + This method looks for the correct driver in self.drivers to use; and puts the curve content in the .curve attribute. Basically, until identify() is called, the HookeCurve object @@ -446,17 +444,17 @@ class tutorialCommands: the Hooke plot routine), the HookeCurve object is "filled" with the actual curve. ''' - + item.identify(self.drivers) - + ''' After the identify(), item.curve contains the curve, and item.curve.default_plots() behaves exactly like self.current.curve.default_plots() -but for the given item. ''' itplot=item.curve.default_plots() - + print 'length of X1 vector:',len(itplot[0].vectors[0][0]) #just to show something - + ''' The following three lines are a magic spell you HAVE to do before closing the function. @@ -465,30 +463,30 @@ class tutorialCommands: item.curve.close_all() #Avoid open files dangling del item.curve #Avoid memory leaks del item #Just be paranoid. Don't ask. - + return - - + + c=0 for item in self.current_list: print 'Looking at curve ',c,'of',len(self.current_list) things_when_cycling(item) c+=1 - + return - - - + + + def plotmanip_absvalue(self, plot, current, customvalue=None): ''' This function defines a PLOT MANIPULATOR. A plot manipulator is a function that takes a plot in input, does something to the plot and returns the modified plot in output. The function, once plugged, gets automatically called everytime self.plots is updated - + For example, in force spectroscopy force curves, raw data are automatically corrected for deflection. Other data can be, say, filtered by default. - + To create and activate a plot manipulator you have to: * Write a function (like this) which name starts with "plotmanip_" (just like commands start with "do_") @@ -500,20 +498,20 @@ class tutorialCommands: * The function must return a plot object. * Add an entry in hooke.conf: if your function is "plotmanip_something" you will have to add in the plotmanips section: example - + - + - + Important: Plot manipulators are *in pipe*: each plot manipulator output becomes the input of the next one. The order in hooke.conf *is the order* in which plot manipulators are connected, so in the example above we have: self.current.curve.default_plots() --> detriggerize --> correct --> median --> something --> self.plots ''' - + ''' Here we see what is in a configuration variable to enable/disable the plot manipulator as user will using the Hooke "set" command. @@ -521,16 +519,16 @@ class tutorialCommands: ''' if not self.config['tutorial_absvalue']: return plot - + #We do something to the plot, for demonstration's sake #If we needed variables, we would have used customvalue. plot.vectors[0][1]=[abs(i) for i in plot.vectors[0][1]] plot.vectors[1][1]=[abs(i) for i in plot.vectors[1][1]] - + #Return the plot object. return plot - - + + #TODO IN TUTORIAL: #how to add lines to an existing plot!! #peaks diff --git a/hooke/plugin/viewer.py b/hooke/plugin/viewer.py index 42f4b6d..961b6a2 100644 --- a/hooke/plugin/viewer.py +++ b/hooke/plugin/viewer.py @@ -1,5 +1,3 @@ -#!/usr/bin/env python - ''' Viewer test case @@ -23,7 +21,7 @@ class viewerCommands: def do_vwnew(self,args): #creates a new viewer self.viewerlist.append(lview.Ascii(self.outlet)) - dt=linput.safeinput('What type of data will this viewer handle? (force/distance/all)',['force', 'distance', 'all']) + dt=linput.safeinput('What type of data will this viewer handle? (force/distance/all)',['force', 'distance', 'all']) #TODO update types, make a list somewhere? print dt self.viewerlist[-1].setdtype(dt) @@ -36,4 +34,4 @@ class viewerCommands: if len(args)==0: args=0 - self.viewerlist[int(args)].action() \ No newline at end of file + self.viewerlist[int(args)].action()