hooke/plugin/generalclamp.py

   1 # Copyright (C) 2008-2010 Alberto Gomez-Casado
   2 #                         Marco Brucale
   3 #                         Massimo Sandal <devicerandom@gmail.com>
   4 #                         W. Trevor King <wking@drexel.edu>
   5 #
   6 # This file is part of Hooke.
   7 #
   8 # Hooke is free software: you can redistribute it and/or
   9 # modify it under the terms of the GNU Lesser General Public
  10 # License as published by the Free Software Foundation, either
  11 # version 3 of the License, or (at your option) any later version.
  12 #
  13 # Hooke is distributed in the hope that it will be useful,
  14 # but WITHOUT ANY WARRANTY; without even the implied warranty of
  15 # MERCHANTABILITY or FITNESS FOR A PARTICULAR PURPOSE.  See the
  16 # GNU Lesser General Public License for more details.
  17 #
  18 # You should have received a copy of the GNU Lesser General Public
  19 # License along with Hooke.  If not, see
  20 # <http://www.gnu.org/licenses/>.
  21
  22 """Plugin regarding general force clamp measurements
  23 """
  24
  25 from hooke.libhooke import WX_GOOD, ClickedPoint
  26
  27 import wxversion
  28 wxversion.select(WX_GOOD)
  29 from wx import PostEvent
  30
  31 from .. import curve as lhc
  32
  33
  34 class generalclampCommands(object):
  35
  36     def plotmanip_clamp(self, plot, current, customvalue=False):
  37         '''
  38         Handles some viewing options for the "force clamp" data format, depending on the state of these configuration variables:
  39         (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)
  40         (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)
  41         (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.
  42
  43         NOTE - my implementation of point(3) feels quite awkward - someone smarter than me plz polish that!
  44
  45         '''
  46
  47         #not a fclamp curve...
  48         if current.curve.experiment != 'clamp':
  49             return plot
  50
  51         if self.config['fc_interesting'] != 0 and plot.destination==0:
  52             lower = int((self.config['fc_interesting'])-1)
  53             upper = int((self.config['fc_interesting'])+1)
  54             trim = current.curve.trimindexes()[lower:upper]
  55             newtime = []
  56             newzpiezo = []
  57             newphase = []
  58             for x in range(trim[0],trim[1]):
  59                 newtime.append(self.plots[0].vectors[0][0][x])
  60                 newzpiezo.append(self.plots[0].vectors[0][1][x])
  61                 newphase.append(self.plots[0].vectors[1][1][x])
  62             self.plots[0].vectors[0][0] = newtime
  63             self.plots[0].vectors[0][1] = newzpiezo
  64             self.plots[0].vectors[1][0] = newtime
  65             self.plots[0].vectors[1][1] = newphase
  66
  67         if self.config['fc_interesting'] != 0 and plot.destination==1:
  68             lower = int((self.config['fc_interesting'])-1)
  69             upper = int((self.config['fc_interesting'])+1)
  70             trim = current.curve.trimindexes()[lower:upper]
  71             newtime = []
  72             newdefl = []
  73             newimposed = []
  74             for x in range(trim[0],trim[1]):
  75                 newtime.append(self.plots[1].vectors[0][0][x])
  76                 newdefl.append(self.plots[1].vectors[0][1][x])
  77                 newimposed.append(self.plots[1].vectors[1][1][x])
  78             self.plots[1].vectors[0][0] = newtime
  79             self.plots[1].vectors[0][1] = newdefl
  80             self.plots[1].vectors[1][0] = newtime
  81             self.plots[1].vectors[1][1] = newimposed
  82
  83         if self.config['fc_showphase'] == 0 and plot.destination==0:
  84             self.plots[0].remove_set(1)
  85
  86         if self.config['fc_showimposed'] == 0 and plot.destination==1:
  87             self.plots[1].remove_set(1)
  88
  89         return plot
  90
  91     def do_time(self,args):
  92         '''
  93         Measures the time difference (in seconds) between two points
  94         Implemented only for force clamp
  95         ----
  96         Syntax: time
  97         '''
  98         if self.current.curve.experiment == 'clamp':
  99             time=self._delta(set=0)[0]
 100             print str(time*1000)+' ms'
 101         else:
 102             print 'This command makes no sense for a non-force clamp experiment.'
 103
 104     def do_zpiezo(self,args):
 105         '''
 106         Measures the zpiezo difference (in nm) between two points
 107         Implemented only for force clamp
 108         ----
 109         Syntax: zpiezo
 110         '''
 111         if self.current.curve.experiment == 'clamp':
 112             zpiezo=self._delta(set=0)[2]
 113             print str(zpiezo*(10**9))+' nm'
 114         else:
 115             print 'This command makes no sense for a non-force clamp experiment.'
 116
 117     def do_defl(self,args):
 118         '''
 119         Measures the deflection difference (in nm) between two points
 120         Implemented only for force clamp
 121         NOTE: It makes sense only on the time VS defl plot; it is still not masked for the other plot...
 122         -----
 123         Syntax: defl
 124         '''
 125         if self.current.curve.experiment == 'clamp':
 126             print "Warning - don't use on the zpiezo plot!"
 127             defl=self._delta(set=1)[2]
 128             print str(defl*(10**12))+' pN'
 129         else:
 130             print 'This command makes no sense for a non-force clamp experiment.'
 131
 132     def do_step(self,args):
 133         '''
 134         Measures the length and time duration of a time-Z step
 135         -----
 136         Syntax: step
 137         '''
 138         if self.current.curve.experiment == 'clamp':
 139             print 'Click three points in this fashion:'
 140             print ' (0)-------(1)'
 141             print '           |'
 142             print '           |'
 143             print '           (2)----------'
 144             points=self._measure_N_points(N=3,whatset=0)
 145             dz=abs(points[2].graph_coords[1]-points[1].graph_coords[1])*(10e+8)
 146             dt=abs(points[1].graph_coords[0]-points[0].graph_coords[0])
 147             print 'dZ: ',dz,' nm'
 148             print 'dT: ',dt,' s'
 149
 150         else:
 151             print 'This command makes no sense for a non-force clamp experiment.'
 152
 153     def do_fcfilt(self,args):
 154         '''
 155         Filters out featureless force clamp curves of the current playlist.
 156         It's very similar to 'flatfilt' for velocity clamp curves.
 157         Creates a new playlist only containing non-empty curves.
 158
 159         WARNING - Only works if you set an appropriate fc_interesting config variable!
 160         WARNING - arguments are NOT optional at the moment!
 161
 162         Syntax: fcfilt maxretraction(nm) mindeviation (pN)
 163
 164         Suggested values for an (i27)8 experiment with our setup are 200nm and 10-15 pN
 165         '''
 166
 167         if self.config['fc_interesting'] == 0:
 168             print 'You must specify the phase of interest (using set fc_interesing X) prior to running fcfilt!'
 169             return
 170
 171         maxretraction=0
 172         threshold=0
 173         args=args.split(' ')
 174         if len(args)==2:
 175             maxretraction=int(args[0])
 176             threshold=int(args[1])
 177         else:
 178             print 'Arguments are not optional for fcfilt. You should pass two numbers:'
 179             print '(1) the maximum plausible piezo retraction in NANOMETERS (e.g. the length of the protein)'
 180             print "(2) the threshold, in PICONEWTONS. If signal deviates from imposed more than this, it's an event"
 181             return
 182
 183
 184         print 'Processing playlist... go get yourself a cup of coffee.'
 185         notflat_list=[]
 186
 187         c=0
 188
 189         for item in self.current_list:
 190             c+=1
 191             try:
 192                 notflat=self.has_stuff(item,maxretraction,threshold)
 193                 print 'Curve',item.path,'is',c,'of',len(self.current_list),'--->Has Stuff =',notflat
 194             except:
 195                 notflat=False
 196                 print 'Curve',item.path,'is',c,'of',len(self.current_list),'--->could not be processed'
 197             if notflat:
 198                 item.features=notflat
 199                 item.curve=None
 200                 notflat_list.append(item)
 201
 202         if len(notflat_list)==0:
 203             print 'Nothing interesting here. Reconsider either your filtering criteria or your experimental data'
 204             return
 205         else:
 206             print 'Found',len(notflat_list),'potentially interesting curves.'
 207             print 'Regenerating playlist...'
 208             self.pointer=0
 209             self.current_list=notflat_list
 210             self.current=self.current_list[self.pointer]
 211             self.do_plot(0)
 212
 213     def has_stuff(self,item,maxretraction,threshold):
 214         '''
 215         Decides whether a curve has some features in the interesting phase.
 216         Algorithm:
 217             - clip the interesting phase portion of the curve.
 218             - discard the first 20 milliseconds (this is due to a quirk of our hardware).
 219             - look at the zpiezo plot and note down when (if) retratcs more than [maxretraction] nm away from the first point.
 220             - clip off any data after this point, with an excess of 100 points (again, an hardware quirk)
 221             - if the remainder is less than 100 points, ditch the curve.
 222             - now look at the deflection plot and check if there are points more than [threshold] pN over the 'flat zone'.
 223             - if you find such points, bingo!
 224         '''
 225
 226         item.identify(self.drivers)
 227
 228         lower = int((self.config['fc_interesting'])-1)
 229         upper = int((self.config['fc_interesting'])+1)
 230         trim_idxs = item.curve.trimindexes()[lower:upper]
 231         lo=trim_idxs[0]+20                                                  #clipping the first 20 points off...
 232         hi=trim_idxs[1]
 233         trimmed_zpiezo=item.curve.default_plots()[0].vectors[0][1][lo:hi]
 234         trimmed_defl=item.curve.default_plots()[1].vectors[0][1][lo:hi]
 235         trimmed_imposed=item.curve.default_plots()[1].vectors[1][1][lo:hi]
 236         imposed=trimmed_imposed[21]                                         #just to match the 20-pts clipping...
 237
 238         item.curve.close_all()
 239         del item.curve
 240         del item
 241
 242         starting_z=trimmed_zpiezo[0]
 243         plausible=starting_z-(maxretraction*1e-9)
 244         det_trim=0
 245         while trimmed_zpiezo[det_trim]>plausible:
 246             det_trim+=1
 247             if det_trim >= len(trimmed_zpiezo):                              #breaking cycles makes me shiver...
 248                 det_trim=len(trimmed_zpiezo)                                 #but I cannot think of anything better now.
 249                 break
 250         further_trim=det_trim-100
 251         if further_trim<100:
 252             return False
 253         trimmed_defl=trimmed_defl[:further_trim]
 254
 255         trimmed_defl.sort()
 256         ninetypercent=int(0.9*len(trimmed_defl))
 257         j=0
 258         sum=0
 259         for j in trimmed_defl[:ninetypercent]:
 260             sum+=j
 261         avg=float(sum/ninetypercent)
 262         sweetspot=float(avg+(threshold*1e-12))
 263         if trimmed_defl[-1]>sweetspot:
 264             flag=True
 265         else:
 266             flag=False
 267
 268         del trimmed_defl,trimmed_zpiezo,trimmed_imposed
 269
 270         return flag
 271