auto_min_p ; auto_max_p = Minimum and maximum persistence length (if using WLC) or Kuhn length (if using FJC)
outside of which the peak is automatically discarded (in nm)
'''
-
- #MACROS.
- #FIXME: to move outside function
- def fit_interval_nm(start_index,plot,nm,backwards):
- '''
- Calculates the number of points to fit, given a fit interval in nm
- start_index: index of point
- plot: plot to use
- backwards: if true, finds a point backwards.
- '''
- whatset=1 #FIXME: should be decidable
- x_vect=plot.vectors[1][0]
- c=0
- i=start_index
- start=x_vect[start_index]
- maxlen=len(x_vect)
- while abs(x_vect[i]-x_vect[start_index])*(10**9) < nm:
- if i==0 or i==maxlen-1: #we reached boundaries of vector!
- return c
-
- if backwards:
- i-=1
- else:
- i+=1
- c+=1
- return c
-
- def pickup_contact_point():
- '''macro to pick up the contact point by clicking'''
- contact_point=self._measure_N_points(N=1, whatset=1)[0]
- contact_point_index=contact_point.index
- self.wlccontact_point=contact_point
- self.wlccontact_index=contact_point.index
- self.wlccurrent=self.current.path
- return contact_point, contact_point_index
-
- def find_current_peaks(noflatten):
- #Find peaks.
- defplot=self.current.curve.default_plots()[0]
- if not noflatten:
- flatten=self._find_plotmanip('flatten') #Extract flatten plotmanip
- defplot=flatten(defplot, self.current, customvalue=1) #Flatten curve before feeding it to has_peaks
- peak_location,peak_size=self.has_peaks(defplot, self.convfilt_config['mindeviation'])
- return peak_location, peak_size
-
#default fit etc. variables
pl_value=None
T=self.config['temperature']
#--END COMMAND LINE PARSING--
- peak_location, peak_size = find_current_peaks(noflatten)
+ peak_location, peak_size = self.find_current_peaks(noflatten)
if len(peak_location) == 0:
print 'No peaks to fit.'
#Pick up force baseline
if rebase:
- clicks=self.config['baseline_clicks']
- if clicks==0:
- self.basepoints=[]
- base_index_0=peak_location[-1]+fit_interval_nm(peak_location[-1], displayed_plot, self.config['auto_right_baseline'],False)
- self.basepoints.append(self._clickize(displayed_plot.vectors[1][0],displayed_plot.vectors[1][1],base_index_0))
- base_index_1=self.basepoints[0].index+fit_interval_nm(self.basepoints[0].index, displayed_plot, self.config['auto_left_baseline'],False)
- self.basepoints.append(self._clickize(displayed_plot.vectors[1][0],displayed_plot.vectors[1][1],base_index_1))
- elif clicks>0:
- print 'Select baseline'
- if clicks==1:
- self.basepoints=self._measure_N_points(N=1, whatset=whatset)
- base_index_1=self.basepoints[0].index+fit_interval_nm(self.basepoints[0].index, displayed_plot, self.config['auto_left_baseline'], False)
- self.basepoints.append(self._clickize(displayed_plot.vectors[1][0],displayed_plot.vectors[1][1],base_index_1))
- else:
- self.basepoints=self._measure_N_points(N=2, whatset=whatset)
-
- self.basecurrent=self.current.path
+ self.basepoints=self.baseline_points(peak_location, displayed_plot)
boundaries=[self.basepoints[0].index, self.basepoints[1].index]
boundaries.sort()
#WLC FITTING
#define fit interval
if not usepoints:
- fit_points=fit_interval_nm(peak, displayed_plot, self.config['auto_fit_nm'], True)
+ fit_points=self.fit_interval_nm(peak, displayed_plot, self.config['auto_fit_nm'], True)
peak_point=self._clickize(displayed_plot.vectors[1][0],displayed_plot.vectors[1][1],peak)
other_fit_point=self._clickize(displayed_plot.vectors[1][0],displayed_plot.vectors[1][1],peak-fit_points)
#Show wlc fits and peak locations
self._send_plot([fitplot])
- #self.do_peaks('')
print 'Using fit function: ',self.config['fit_function']
print 'Measurements for all peaks detected:'
--- /dev/null
+# -*- coding: utf-8 -*-
+from libhooke import WX_GOOD, ClickedPoint
+import wxversion
+wxversion.select(WX_GOOD)
+from wx import PostEvent
+import numpy as np
+import scipy as sp
+import copy
+import os.path
+import time
+
+
+class curvetoolsCommands:
+
+ def fit_interval_nm(self,start_index,plot,nm,backwards):
+ '''
+ Calculates the number of points to fit, given a fit interval in nm
+ start_index: index of point
+ plot: plot to use
+ backwards: if true, finds a point backwards.
+ '''
+ whatset=1 #FIXME: should be decidable
+ x_vect=plot.vectors[1][0]
+
+ c=0
+ i=start_index
+ start=x_vect[start_index]
+ maxlen=len(x_vect)
+ while abs(x_vect[i]-x_vect[start_index])*(10**9) < nm:
+ if i==0 or i==maxlen-1: #we reached boundaries of vector!
+ return c
+
+ if backwards:
+ i-=1
+ else:
+ i+=1
+ c+=1
+ return c
+
+
+
+ def find_current_peaks(self,noflatten, a=True, maxpeak=True):
+ #Find peaks.
+ if a==True:
+ a=self.convfilt_config['mindeviation']
+ try:
+ abs_devs=float(a)
+ except:
+ print "Bad input, using default."
+ abs_devs=self.convfilt_config['mindeviation']
+
+ defplot=self.current.curve.default_plots()[0]
+ if not noflatten:
+ flatten=self._find_plotmanip('flatten') #Extract flatten plotmanip
+ defplot=flatten(defplot, self.current, customvalue=1) #Flatten curve before feeding it to has_peaks
+ pk_location,peak_size=self.has_peaks(defplot, abs_devs, maxpeak)
+ return pk_location, peak_size
+
+
+ def pickup_contact_point(self,N=1,whatset=1):
+ '''macro to pick up the contact point by clicking'''
+ contact_point=self._measure_N_points(N=1, whatset=1)[0]
+ contact_point_index=contact_point.index
+ self.wlccontact_point=contact_point
+ self.wlccontact_index=contact_point.index
+ self.wlccurrent=self.current.path
+ return contact_point, contact_point_index
+
+
+
+ def baseline_points(self,peak_location, displayed_plot):
+ clicks=self.config['baseline_clicks']
+ if clicks==0:
+ self.basepoints=[]
+ base_index_0=peak_location[-1]+self.fit_interval_nm(peak_location[-1], displayed_plot, self.config['auto_right_baseline'],False)
+ self.basepoints.append(self._clickize(displayed_plot.vectors[1][0],displayed_plot.vectors[1][1],base_index_0))
+ base_index_1=self.basepoints[0].index+self.fit_interval_nm(self.basepoints[0].index, displayed_plot, self.config['auto_left_baseline'],False)
+ self.basepoints.append(self._clickize(displayed_plot.vectors[1][0],displayed_plot.vectors[1][1],base_index_1))
+ elif clicks>0:
+ print 'Select baseline'
+ if clicks==1:
+ self.basepoints=self._measure_N_points(N=1, whatset=1)
+ base_index_1=self.basepoints[0].index+self.fit_interval_nm(self.basepoints[0].index, displayed_plot, self.config['auto_left_baseline'], False)
+ self.basepoints.append(self._clickize(displayed_plot.vectors[1][0],displayed_plot.vectors[1][1],base_index_1))
+ else:
+ self.basepoints=self._measure_N_points(N=2, whatset=1)
+
+ self.basecurrent=self.current.path
+ return self.basepoints
+
+
+
#!/usr/bin/env python
+# -*- coding: utf-8 -*-
'''
FLATFILTS
#-----Convolution-based peak recognition and filtering.
#Requires the libpeakspot.py library
- def has_peaks(self, plot, abs_devs=None):
+ def has_peaks(self, plot, abs_devs=None, maxpeak=True, window=10):
'''
Finds peak position in a force curve.
FIXME: should be moved in libpeakspot.py
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
+ #take the minimum or the maximum of a peak
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
+ valpk=min(yret[peak-window:peak+window]) #maximum in force (near the unfolding point)
+ index_pk=yret[peak-window:peak+window].index(valpk)+(peak-window)
+
+ if maxpeak==False:
+ valpk=max(yret[peak:peak+window]) #minimum in force, near the baseline
+ index_pk=yret[peak:peak+window].index(valpk)+(peak)
+
+# Let's explain that for the minimum. Immaging that we know that there is a peak at position/region 100 and you have found its y-value,
+# Now you look in the array, from 100-10 to 100+10 (if the window is 10).
+# This "100-10 to 100+10" is substancially a new array with its index. In this array you have 20
+# elements, so the index of your y-value will be 10.
+# Now to find the index in the TOTAL array you have to add the "position" of the "region" (that in this case
+# correspond to 100) and also substract the window size ---> (+100-10)
+
+ peak_location[i]=index_pk
return peak_location,peak_size
item.peak_size=peak_size
item.curve=None #empty the item object, to further avoid memory leak
notflat_list.append(item)
-
+
+ for i in range(1000):
+ k=0
+
#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.'
-->\r
<plugins>\r
<fit/>\r
+ <curvetools/>\r
<procplots/>\r
<flatfilts/>\r
<generalclamp/>\r
<pcluster/>\r
<generaltccd/>\r
<multidistance/>\r
+ <jumpstat/>\r
</plugins>\r
\r
<!--\r
--- /dev/null
+# -*- coding: utf-8 -*-
+from libhooke import WX_GOOD, ClickedPoint
+import wxversion
+wxversion.select(WX_GOOD)
+from wx import PostEvent
+import numpy as np
+import scipy as sp
+import copy
+import os.path
+import time
+import sys
+import warnings
+warnings.simplefilter('ignore',np.RankWarning)
+
+
+class jumpstatCommands():
+
+ def do_jumpstat(self,args):
+ '''
+ JUMPSTAT
+ jumpstat.py
+ Based on the convolution recognition automatically give:
+ - the delta distance between the peaks,
+ - the delta-force from the top of the peaks and subsequent relaxation,
+ - the delta-force from the top of the peaks and the baseline
+ The command allow also to remove the unwanted peaks that can be due to interference.
+ When you first issue the command, it will ask for the filename. If you are giving the filename
+ 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:
+ jumpstat [deviation]
+ deviation = number of times the convolution signal is above the noise absolute deviation.
+ '''
+
+
+ #finding the max and the minimum positions for all the peaks
+ noflatten=False
+ #we use if else only to avoid a "bad input" message from find_current_peaks
+ if (len(args)==0):
+ max_peaks_location, peak_size=self.find_current_peaks(noflatten)
+ min_peaks_location, pks2=self.find_current_peaks(noflatten, True, False)
+ else:
+ max_peaks_location, peak_size=self.find_current_peaks(noflatten, args)
+ min_peaks_location, pks2=self.find_current_peaks(noflatten, args, False)
+
+
+ #print "max_peaks_location: "+str(len(max_peaks_location))
+ #print "min_peaks_location: "+str(len(min_peaks_location))
+
+ #if no peaks, we have nothing to plot. exit.
+ if len(max_peaks_location)==0:
+ print "No peaks on this curve."
+ return
+
+ if len(max_peaks_location)!=len(min_peaks_location):
+ print "Something went wrong in peaks recognition, number of minima is different from number of maxima. Exiting."
+ 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 max_peaks_location]
+ ygood=[yplotted_ret[index] for index in max_peaks_location]
+
+ xafter=[xplotted_ret[index] for index in min_peaks_location]
+ yafter=[yplotted_ret[index] for index in min_peaks_location]
+
+ recplot=self._get_displayed_plot()
+ recplot2=self._get_displayed_plot()
+ recplot.vectors.append([xgood,ygood])
+ recplot2.vectors.append([xafter,yafter])
+
+ if recplot.styles==[]:
+ recplot.styles=[None,None,'scatter']
+ recplot.colors=[None,None,None]
+ else:
+ recplot.styles+=['scatter']
+ recplot.colors+=[None]
+
+ if recplot2.styles==[]:
+ recplot2.styles=[None,None,None]
+ recplot2.colors=[None,'1.0',None]
+ else:
+ recplot2.styles+=['scatter']
+ recplot2.colors+=['0.5']
+
+ self._send_plot([recplot])
+ self._send_plot([recplot2])
+
+
+ #finding the baseline
+ self.basepoints=self.baseline_points(max_peaks_location, recplot)
+ boundaries=[self.basepoints[0].index, self.basepoints[1].index]
+ boundaries.sort()
+ to_average=recplot.vectors[1][1][boundaries[0]:boundaries[1]] #y points to average
+ avg=np.mean(to_average)
+
+
+ dist=[]
+ jumpforce=[]
+ force=[]
+
+ #we calculate the distance vector
+ for g in range(len(max_peaks_location)-1):
+ dist.append((10**9)*(xplotted_ret[max_peaks_location[g]]-xplotted_ret[max_peaks_location[g+1]]))
+ print "Distance values for the peaks in nm:"
+ print dist
+
+ #the jump-force vector
+ for g in range(len(max_peaks_location)):
+ jumpforce.append((10**12) *(yplotted_ret[min_peaks_location[g]] -yplotted_ret[max_peaks_location[g]]) )
+ print "Force values for the jumps of the peaks in pN:"
+ print jumpforce
+
+ #the force from baseline vector
+ for g in range(len(max_peaks_location)):
+ force.append((10**12)*(avg-yplotted_ret[max_peaks_location[g]]))
+ print "Force values for the peaks in pN:"
+ print force
+
+
+
+ #Now ask for the peaks that we don't want
+ print 'Peaks to ignore (0,1...n from contact point,return to take all)'
+ print 'N to discard measurement'
+ exclude_raw=raw_input('Input:')
+ if exclude_raw=='N':
+ print 'Discarded.'
+ return
+
+ if not exclude_raw=='':
+ exclude=exclude_raw.split(',')
+ #we convert in numbers the input
+ try:
+ exclude=[int(item) for item in exclude]
+ except:
+ print 'Bad input, taking nothing.'
+ return
+
+# we remove the peaks that we don't want from the list, we need a counter beacause if we remove
+# a peaks the other peaks in the list are shifted by one at each step
+ count=0
+ for a in exclude:
+ if (a==0):
+ max_peaks_location=max_peaks_location[1:]
+ min_peaks_location=min_peaks_location[1:]
+ else:
+ new_a=a-count
+ max_peaks_location= max_peaks_location[0:new_a]+max_peaks_location[new_a+1:]
+ min_peaks_location= min_peaks_location[0:new_a]+min_peaks_location[new_a+1:]
+ peak_size= peak_size[0:new_a]+peak_size[new_a+1:]
+ count+=1
+
+
+ #print "max_peaks_location: "+str(len(max_peaks_location))
+ #print "min_peaks_location: "+str(len(min_peaks_location))
+
+
+ dist=[]
+ jumpforce=[]
+ force=[]
+ #we recalculate the distances and the forces after the removing of the unwanted peaks
+ for g in range(len(max_peaks_location)-1):
+ dist.append(xplotted_ret[max_peaks_location[g]]-xplotted_ret[max_peaks_location[g+1]])
+ for g in range(len(max_peaks_location)):
+ jumpforce.append( yplotted_ret[min_peaks_location[g]] - yplotted_ret[max_peaks_location[g]] )
+ for g in range(len(max_peaks_location)):
+ force.append(avg - yplotted_ret[max_peaks_location[g]])
+
+
+
+
+
+ #Save file info
+ if self.autofile=='':
+ self.autofile=raw_input('Jumpstat filename? (return to ignore) ')
+ 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')
+ f.write('----------------------------------------\n')
+ f.write('; Delta Distance length (m); Jump Force pN; Standard Force pN\n')
+ f.write(self.current.path+'\n')
+ for k in range(len(dist)):
+ f.write(";")
+ f.write(str(dist[k])+";"+str(jumpforce[k])+";"+str(force[k])+"\n" )
+ f.write("\n")
+ f.close()
+
+ else:
+ f=open(self.autofile,'a+')
+ f.write(self.current.path+'\n')
+ for k in range(len(dist)):
+ f.write(";")
+ f.write(str(dist[k])+";"+str(jumpforce[k])+";"+str(force[k])+"\n" )
+ f.write("\n")
+ f.close()
+
+ print 'Saving...'
\ No newline at end of file
deviation = number of times the convolution signal is above the noise absolute deviation.
'''
- def find_current_peaks(noflatten, a):
- #Find peaks.
- if len(a)==0:
- a=self.convfilt_config['mindeviation']
- try:
- abs_devs=float(a)
- except:
- print "Bad input, using default."
- abs_devs=self.convfilt_config['mindeviation']
-
- defplot=self.current.curve.default_plots()[0]
- if not noflatten:
- flatten=self._find_plotmanip('flatten') #Extract flatten plotmanip
- defplot=flatten(defplot, self.current, customvalue=1) #Flatten curve before feeding it to has_peaks
- 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)
+ peaks_location, peak_size=self.find_current_peaks(noflatten)
#if no peaks, we have nothing to plot. exit.
if len(peaks_location)==0:
projects / hooke.git / commitdiff