hooke/plugin/fclamp.py

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