1 # -*- coding: utf-8 -*-
3 from hooke.libhooke import WX_GOOD, ClickedPoint
6 wxversion.select(WX_GOOD)
7 from wx import PostEvent
15 warnings.simplefilter('ignore',np.RankWarning)
18 class autopeakCommands(object):
20 def do_autopeak(self,args):
24 Automatically performs a number of analyses on the peaks of the given curve.
25 Currently it automatically:
26 - fits peaks with WLC or FJC function (depending on how the fit_function variable is set)
27 - measures peak maximum forces with a baseline
28 - measures slope in proximity of peak maximum
29 Requires flatten plotmanipulator , fit.py plugin , flatfilts.py plugin with convfilt
32 autopeak [rebase] [pl=value] [t=value] [noauto] [reclick]
34 rebase : Re-asks baseline interval
36 pl=[value] : Use a fixed persistent length for the fit. If pl is not given,
37 the fit will be a 2-variable
38 fit. DO NOT put spaces between 'pl', '=' and the value.
39 The value must be in meters.
40 Scientific notation like 0.35e-9 is fine.
42 t=[value] : Use a user-defined temperature. The value must be in
43 kelvins; by default it is 293 K.
44 DO NOT put spaces between 't', '=' and the value.
46 noauto : allows for clicking the contact point by
47 hand (otherwise it is automatically estimated) the first time.
48 If subsequent measurements are made, the same contact point
49 clicked the first time is used
51 reclick : redefines by hand the contact point, if noauto has been used before
52 but the user is unsatisfied of the previously choosen contact point.
54 usepoints : fit interval by number of points instead than by nanometers
56 noflatten : does not use the "flatten" plot manipulator
58 When you first issue the command, it will ask for the filename. If you are giving the filename
59 of an existing file, autopeak will resume it and append measurements to it. If you are giving
60 a new filename, it will create the file and append to it until you close Hooke.
63 Useful variables (to set with SET command):
65 fit_function = type of function to use for elasticity. If "wlc" worm-like chain is used, if "fjc" freely jointed
68 temperature= temperature of the system for wlc/fjc fit (in K)
70 auto_slope_span = number of points on which measure the slope, for slope
72 auto_fit_nm = number of nm to fit before the peak maximum, for WLC/FJC (if usepoints false)
73 auto_fit_points = number of points to fit before the peak maximum, for WLC/FJC (if usepoints true)
75 baseline_clicks = -1: no baseline, f=0 at the contact point (whether hand-picked or automatically found)
77 1: decide baseline with a single click and length defined in auto_left_baseline
78 2: let user click points of baseline
79 auto_left_baseline = length in nm to use as baseline from the right point (if baseline_clicks=0 , 1)
80 auto_right_baseline = distance in nm of peak-most baseline point from last peak (if baseline_clicks = 0)
82 auto_min_p ; auto_max_p = Minimum and maximum persistence length (if using WLC) or Kuhn length (if using FJC)
83 outside of which the peak is automatically discarded (in nm)
86 #default fit etc. variables
88 T=self.config['temperature']
90 slope_span=int(self.config['auto_slope_span'])
92 rebase=False #if true=we select rebase
93 noflatten=False #if true=we avoid flattening
95 #initialize output data vectors
104 displayed_plot=self._get_displayed_plot(0)
106 #COMMAND LINE PARSING
107 #--Using points instead of nm interval
108 if 'usepoints' in args.split():
109 fit_points=int(self.config['auto_fit_points'])
114 #--Recalculate baseline
115 if 'rebase' in args or (self.basecurrent != self.current.path):
118 if 'noflatten' in args:
121 #--Custom persistent length / custom temperature
122 for arg in args.split():
123 #look for a persistent length argument.
125 pl_expression=arg.split('=')
126 pl_value=float(pl_expression[1]) #actual value
127 #look for a T argument. FIXME: spaces are not allowed between 'pl' and value
128 if ('t=' in arg[0:2]) or ('T=' in arg[0:2]):
129 t_expression=arg.split('=')
130 T=float(t_expression[1])
131 #--Contact point arguments
132 if 'reclick' in args.split():
133 print 'Click contact point'
134 contact_point, contact_point_index = self.pickup_contact_point()
135 elif 'noauto' in args.split():
136 if self.wlccontact_index==None or self.wlccurrent != self.current.path:
137 print 'Click contact point'
138 contact_point , contact_point_index = self.pickup_contact_point()
140 contact_point=self.wlccontact_point
141 contact_point_index=self.wlccontact_index
143 #Automatically find contact point
144 cindex=self.find_contact_point()
145 contact_point=self._clickize(displayed_plot.vectors[1][0], displayed_plot.vectors[1][1], cindex)
146 #--END COMMAND LINE PARSING--
149 peak_location, peak_size = self.find_current_peaks(noflatten)
151 if len(peak_location) == 0:
152 print 'No peaks to fit.'
155 fitplot=copy.deepcopy(displayed_plot)
157 #Pick up force baseline
159 self.basepoints=self.baseline_points(peak_location, displayed_plot)
161 boundaries=[self.basepoints[0].index, self.basepoints[1].index]
163 to_average=displayed_plot.vectors[1][1][boundaries[0]:boundaries[1]] #y points to average
164 avg=np.mean(to_average)
166 clicks=self.config['baseline_clicks']
169 avg=displayed_plot.vectors[1][1][contact_point_index]
171 avg=displayed_plot.vectors[1][1][cindex]
173 for peak in peak_location:
177 fit_points=self.fit_interval_nm(peak, displayed_plot, self.config['auto_fit_nm'], True)
178 peak_point=self._clickize(displayed_plot.vectors[1][0],displayed_plot.vectors[1][1],peak)
179 other_fit_point=self._clickize(displayed_plot.vectors[1][0],displayed_plot.vectors[1][1],peak-fit_points)
182 points=[contact_point, peak_point, other_fit_point]
184 if abs(peak_point.index-other_fit_point.index) < 2:
187 if self.config['fit_function']=='wlc':
189 params, yfit, xfit, fit_errors = self.wlc_fit(points, displayed_plot.vectors[1][0], displayed_plot.vectors[1][1], pl_value, T, return_errors=True)
190 elif self.config['fit_function']=='fjc':
191 params, yfit, xfit, fit_errors = self.fjc_fit(points, displayed_plot.vectors[1][0], displayed_plot.vectors[1][1], pl_value, T, return_errors=True)
193 print 'Unknown fit function'
194 print 'Please set fit_function as wlc or fjc'
199 delta_to_measure=displayed_plot.vectors[1][1][peak-delta_force:peak+delta_force]
200 y=min(delta_to_measure)
201 #save force values (pN)
203 slope=self.linefit_between(peak-slope_span,peak)[0]
206 #check fitted data and, if right, add peak to the measurement
207 #FIXME: code duplication
208 if len(params)==1: #if we did choose 1-value fit
209 p_lengths.append(pl_value)
210 c_lengths.append(params[0]*(1.0e+9))
211 sigma_p_lengths.append(0)
212 sigma_c_lengths.append(fit_errors[0]*(1.0e+9))
213 forces.append(abs(y-avg)*(1.0e+12))
215 #Add WLC fit lines to plot
216 fitplot.add_set(xfit,yfit)
217 if len(fitplot.styles)==0:
221 fitplot.styles.append(None)
222 fitplot.colors.append(None)
224 p_leng=params[1]*(1.0e+9)
225 #check if persistent length makes sense. otherwise, discard peak.
226 if p_leng>self.config['auto_min_p'] and p_leng<self.config['auto_max_p']:
227 p_lengths.append(p_leng)
228 c_lengths.append(params[0]*(1.0e+9))
229 sigma_c_lengths.append(fit_errors[0]*(1.0e+9))
230 sigma_p_lengths.append(fit_errors[1]*(1.0e+9))
231 forces.append(abs(y-avg)*(1.0e+12))
234 #Add WLC fit lines to plot
235 fitplot.add_set(xfit,yfit)
236 if len(fitplot.styles)==0:
240 fitplot.styles.append(None)
241 fitplot.colors.append(None)
246 #add basepoints to fitplot
247 fitplot.add_set([self.basepoints[0].graph_coords[0],self.basepoints[1].graph_coords[0]],[self.basepoints[0].graph_coords[1],self.basepoints[1].graph_coords[1]])
248 fitplot.styles.append('scatter')
249 fitplot.colors.append(None)
251 #Show wlc fits and peak locations
252 self._send_plot([fitplot])
254 print 'Using fit function: ',self.config['fit_function']
255 print 'Measurements for all peaks detected:'
256 print 'contour (nm)', c_lengths
257 print 'sigma contour (nm)',sigma_c_lengths
258 print 'p (nm)',p_lengths
259 print 'sigma p (nm)',sigma_p_lengths
260 print 'forces (pN)',forces
261 print 'slopes (N/m)',slopes
264 while controller==False:
265 #Ask the user what peaks to ignore from analysis.
266 print 'Peaks to ignore (0,1...n from contact point,return to take all)'
267 print 'N to discard measurement'
268 exclude_raw=raw_input('Input:')
274 if not exclude_raw=='':
275 exclude=exclude_raw.split(',')
277 exclude=[int(item) for item in exclude]
283 sigma_c_lengths[i]=None
284 sigma_p_lengths[i]=None
290 #Clean data vectors from ignored peaks
291 #FIXME:code duplication
292 c_lengths=[item for item in c_lengths if item != None]
293 p_lengths=[item for item in p_lengths if item != None]
294 forces=[item for item in forces if item != None]
295 slopes=[item for item in slopes if item != None]
296 sigma_c_lengths=[item for item in sigma_c_lengths if item != None]
297 sigma_p_lengths=[item for item in sigma_p_lengths if item != None]
299 print 'Measurements for chosen peaks:'
300 print 'contour (nm)',c_lengths
301 print 'sigma contour (nm)',sigma_c_lengths
302 print 'p (nm)',p_lengths
303 print 'sigma p (nm)',sigma_p_lengths
304 print 'forces (pN)',forces
305 print 'slopes (N/m)',slopes
308 if self.autofile=='':
309 self.autofile=raw_input('Autopeak filename? (return to ignore) ')
310 if self.autofile=='':
314 if not os.path.exists(self.autofile):
315 f=open(self.autofile,'w+')
316 f.write('Analysis started '+time.asctime()+'\n')
317 f.write('----------------------------------------\n')
318 f.write('; Contour length (nm) ; Persistence length (nm) ; Max.Force (pN) ; Slope (N/m) ; Sigma contour (nm) ; Sigma persistence (nm)\n')
322 f=open(self.autofile,'a+')
324 f.write(self.current.path+'\n')
325 for i in range(len(c_lengths)):
326 f.write(' ; '+str(c_lengths[i])+' ; '+str(p_lengths[i])+' ; '+str(forces[i])+' ; '+str(slopes[i])+' ; '+str(sigma_c_lengths[i])+' ; '+str(sigma_p_lengths[i])+'\n')
329 self.do_note('autopeak')