hooke/plugin/generalvclamp.py

   1 # Copyright
   2
   3 """Plugin regarding general velocity clamp measurements
   4 """
   5
   6 from hooke.libhooke import WX_GOOD, ClickedPoint
   7 import wxversion
   8 wxversion.select(WX_GOOD)
   9 from wx import PostEvent
  10 import numpy as np
  11 import scipy as sp
  12 import copy
  13 import os.path
  14 import time
  15
  16 import warnings
  17 warnings.simplefilter('ignore',np.RankWarning)
  18
  19
  20 class generalvclampCommands(object):
  21
  22     def _plug_init(self):
  23         self.basecurrent=None
  24         self.basepoints=None
  25         self.autofile=''
  26
  27     def do_distance(self,args):
  28         '''
  29         DISTANCE
  30         (generalvclamp.py)
  31         Measure the distance (in nm) between two points.
  32         For a standard experiment this is the delta X distance.
  33         For a force clamp experiment this is the delta Y distance (actually becomes
  34         an alias of zpiezo)
  35         -----------------
  36         Syntax: distance
  37         '''
  38         if self.current.curve.experiment == 'clamp':
  39             print 'You wanted to use zpiezo perhaps?'
  40             return
  41         else:
  42             dx,unitx,dy,unity=self._delta(set=1)
  43             print str(dx*(10**9))+' nm'
  44             to_dump='distance '+self.current.path+' '+str(dx*(10**9))+' nm'
  45             self.outlet.push(to_dump)
  46
  47
  48     def do_force(self,args):
  49         '''
  50         FORCE
  51         (generalvclamp.py)
  52         Measure the force difference (in pN) between two points
  53         ---------------
  54         Syntax: force
  55         '''
  56         if self.current.curve.experiment == 'clamp':
  57             print 'This command makes no sense for a force clamp experiment.'
  58             return
  59         dx,unitx,dy,unity=self._delta(set=1)
  60         print str(dy*(10**12))+' pN'
  61         to_dump='force '+self.current.path+' '+str(dy*(10**12))+' pN'
  62         self.outlet.push(to_dump)
  63
  64
  65     def do_forcebase(self,args):
  66         '''
  67         FORCEBASE
  68         (generalvclamp.py)
  69         Measures the difference in force (in pN) between a point and a baseline
  70         took as the average between two points.
  71
  72         The baseline is fixed once for a given curve and different force measurements,
  73         unless the user wants it to be recalculated
  74         ------------
  75         Syntax: forcebase [rebase]
  76                 rebase: Forces forcebase to ask again the baseline
  77                 max: Instead of asking for a point to measure, asks for two points and use
  78                      the maximum peak in between
  79         '''
  80         rebase=False #if true=we select rebase
  81         maxpoint=False #if true=we measure the maximum peak
  82
  83         plot=self._get_displayed_plot()
  84         whatset=1 #fixme: for all sets
  85         if 'rebase' in args or (self.basecurrent != self.current.path):
  86             rebase=True
  87         if 'max' in args:
  88             maxpoint=True
  89
  90         if rebase:
  91             print 'Select baseline'
  92             self.basepoints=self._measure_N_points(N=2, whatset=whatset)
  93             self.basecurrent=self.current.path
  94
  95         if maxpoint:
  96             print 'Select two points'
  97             points=self._measure_N_points(N=2, whatset=whatset)
  98             boundpoints=[points[0].index, points[1].index]
  99             boundpoints.sort()
 100             try:
 101                 y=min(plot.vectors[whatset][1][boundpoints[0]:boundpoints[1]])
 102             except ValueError:
 103                 print 'Chosen interval not valid. Try picking it again. Did you pick the same point as begin and end of interval?'
 104         else:
 105             print 'Select point to measure'
 106             points=self._measure_N_points(N=1, whatset=whatset)
 107             #whatplot=points[0].dest
 108             y=points[0].graph_coords[1]
 109
 110         #fixme: code duplication
 111         boundaries=[self.basepoints[0].index, self.basepoints[1].index]
 112         boundaries.sort()
 113         to_average=plot.vectors[whatset][1][boundaries[0]:boundaries[1]] #y points to average
 114
 115         avg=np.mean(to_average)
 116         forcebase=abs(y-avg)
 117         print str(forcebase*(10**12))+' pN'
 118         to_dump='forcebase '+self.current.path+' '+str(forcebase*(10**12))+' pN'
 119         self.outlet.push(to_dump)
 120
 121     def plotmanip_multiplier(self, plot, current):
 122         '''
 123         Multiplies all the Y values of an SMFS curve by a value stored in the 'force_multiplier'
 124         configuration variable. Useful for calibrations and other stuff.
 125         '''
 126
 127         #not a smfs curve...
 128         if current.curve.experiment != 'smfs':
 129             return plot
 130
 131         #only one set is present...
 132         if len(self.plots[0].vectors) != 2:
 133             return plot
 134
 135         #multiplier is 1...
 136         if (self.config['force_multiplier']==1):
 137             return plot
 138
 139         for i in range(len(plot.vectors[0][1])):
 140             plot.vectors[0][1][i]=plot.vectors[0][1][i]*self.config['force_multiplier']
 141
 142         for i in range(len(plot.vectors[1][1])):
 143             plot.vectors[1][1][i]=plot.vectors[1][1][i]*self.config['force_multiplier']
 144
 145         return plot
 146
 147
 148     def plotmanip_flatten(self, plot, current, customvalue=False):
 149         '''
 150         Subtracts a polynomial fit to the non-contact part of the curve, as to flatten it.
 151         the best polynomial fit is chosen among polynomials of degree 1 to n, where n is
 152         given by the configuration file or by the customvalue.
 153
 154         customvalue= int (>0) --> starts the function even if config says no (default=False)
 155         '''
 156
 157         #not a smfs curve...
 158         if current.curve.experiment != 'smfs':
 159             return plot
 160
 161         #only one set is present...
 162         if len(self.plots[0].vectors) != 2:
 163             return plot
 164
 165         #config is not flatten, and customvalue flag is false too
 166         if (not self.config['flatten']) and (not customvalue):
 167             return plot
 168
 169         max_exponent=12
 170         delta_contact=0
 171
 172         if customvalue:
 173             max_cycles=customvalue
 174         else:
 175             max_cycles=self.config['flatten'] #Using > 1 usually doesn't help and can give artefacts. However, it could be useful too.
 176
 177         contact_index=self.find_contact_point()
 178
 179         valn=[[] for item in range(max_exponent)]
 180         yrn=[0.0 for item in range(max_exponent)]
 181         errn=[0.0 for item in range(max_exponent)]
 182
 183         #Check if we have a proper numerical value
 184         try:
 185             zzz=int(max_cycles)
 186         except:
 187             #Loudly and annoyingly complain if it's not a number, then fallback to zero
 188             print '''Warning: flatten value is not a number!
 189             Use "set flatten" or edit hooke.conf to set it properly
 190             Using zero.'''
 191             max_cycles=0
 192
 193         for i in range(int(max_cycles)):
 194
 195             x_ext=plot.vectors[0][0][contact_index+delta_contact:]
 196             y_ext=plot.vectors[0][1][contact_index+delta_contact:]
 197             x_ret=plot.vectors[1][0][contact_index+delta_contact:]
 198             y_ret=plot.vectors[1][1][contact_index+delta_contact:]
 199             for exponent in range(max_exponent):
 200                 try:
 201                     valn[exponent]=sp.polyfit(x_ext,y_ext,exponent)
 202                     yrn[exponent]=sp.polyval(valn[exponent],x_ret)
 203                     errn[exponent]=sp.sqrt(sum((yrn[exponent]-y_ext)**2)/float(len(y_ext)))
 204                 except Exception,e:
 205                     print 'Cannot flatten!'
 206                     print e
 207                     return plot
 208
 209             best_exponent=errn.index(min(errn))
 210
 211             #extension
 212             ycorr_ext=y_ext-yrn[best_exponent]+y_ext[0] #noncontact part
 213             yjoin_ext=np.array(plot.vectors[0][1][0:contact_index+delta_contact]) #contact part
 214             #retraction
 215             ycorr_ret=y_ret-yrn[best_exponent]+y_ext[0] #noncontact part
 216             yjoin_ret=np.array(plot.vectors[1][1][0:contact_index+delta_contact]) #contact part
 217
 218             ycorr_ext=np.concatenate((yjoin_ext, ycorr_ext))
 219             ycorr_ret=np.concatenate((yjoin_ret, ycorr_ret))
 220
 221             plot.vectors[0][1]=list(ycorr_ext)
 222             plot.vectors[1][1]=list(ycorr_ret)
 223
 224         return plot
 225
 226     #---SLOPE---
 227     def do_slope(self,args):
 228         '''
 229         SLOPE
 230         (generalvclamp.py)
 231         Measures the slope of a delimited chunk on the return trace.
 232         The chunk can be delimited either by two manual clicks, or have
 233         a fixed width, given as an argument.
 234         ---------------
 235         Syntax: slope [width]
 236                 The facultative [width] parameter specifies how many
 237                 points will be considered for the fit. If [width] is
 238                 specified, only one click will be required.
 239         (c) Marco Brucale, Massimo Sandal 2008
 240         '''
 241
 242         # Reads the facultative width argument
 243         try:
 244             fitspan=int(args)
 245         except:
 246             fitspan=0
 247
 248         # Decides between the two forms of user input, as per (args)
 249         if fitspan == 0:
 250             # Gets the Xs of two clicked points as indexes on the current curve vector
 251             print 'Click twice to delimit chunk'
 252             points=self._measure_N_points(N=2,whatset=1)
 253         else:
 254             print 'Click once on the leftmost point of the chunk (i.e.usually the peak)'
 255             points=self._measure_N_points(N=1,whatset=1)
 256
 257         slope=self._slope(points,fitspan)
 258
 259         # Outputs the relevant slope parameter
 260         print 'Slope:'
 261         print str(slope)
 262         to_dump='slope '+self.current.path+' '+str(slope)
 263         self.outlet.push(to_dump)
 264
 265     def _slope(self,points,fitspan):
 266                 # Calls the function linefit_between
 267         parameters=[0,0,[],[]]
 268         try:
 269             clickedpoints=[points[0].index,points[1].index]
 270             clickedpoints.sort()
 271         except:
 272             clickedpoints=[points[0].index-fitspan,points[0].index]
 273
 274         try:
 275             parameters=self.linefit_between(clickedpoints[0],clickedpoints[1])
 276         except:
 277             print 'Cannot fit. Did you click twice the same point?'
 278             return
 279
 280         # Outputs the relevant slope parameter
 281         print 'Slope:'
 282         print str(parameters[0])
 283         to_dump='slope '+self.curve.path+' '+str(parameters[0])
 284         self.outlet.push(to_dump)
 285
 286         # Makes a vector with the fitted parameters and sends it to the GUI
 287         xtoplot=parameters[2]
 288         ytoplot=[]
 289         x=0
 290         for x in xtoplot:
 291             ytoplot.append((x*parameters[0])+parameters[1])
 292
 293         clickvector_x, clickvector_y=[], []
 294         for item in points:
 295             clickvector_x.append(item.graph_coords[0])
 296             clickvector_y.append(item.graph_coords[1])
 297
 298         lineplot=self._get_displayed_plot(0) #get topmost displayed plot
 299
 300         lineplot.add_set(xtoplot,ytoplot)
 301         lineplot.add_set(clickvector_x, clickvector_y)
 302
 303
 304         if lineplot.styles==[]:
 305             lineplot.styles=[None,None,None,'scatter']
 306         else:
 307             lineplot.styles+=[None,'scatter']
 308         if lineplot.colors==[]:
 309             lineplot.colors=[None,None,'black',None]
 310         else:
 311             lineplot.colors+=['black',None]
 312
 313
 314         self._send_plot([lineplot])
 315
 316         return parameters[0]
 317
 318
 319     def linefit_between(self,index1,index2,whatset=1):
 320         '''
 321         Creates two vectors (xtofit,ytofit) slicing out from the
 322         current return trace a portion delimited by the two indexes
 323         given as arguments.
 324         Then does a least squares linear fit on that slice.
 325         Finally returns [0]=the slope, [1]=the intercept of the
 326         fitted 1st grade polynomial, and [2,3]=the actual (x,y) vectors
 327         used for the fit.
 328         (c) Marco Brucale, Massimo Sandal 2008
 329         '''
 330         # Translates the indexes into two vectors containing the x,y data to fit
 331         xtofit=self.plots[0].vectors[whatset][0][index1:index2]
 332         ytofit=self.plots[0].vectors[whatset][1][index1:index2]
 333
 334         # Does the actual linear fitting (simple least squares with numpy.polyfit)
 335         linefit=[]
 336         linefit=np.polyfit(xtofit,ytofit,1)
 337
 338         return (linefit[0],linefit[1],xtofit,ytofit)