-#!/usr/bin/env python
-# -*- coding: utf-8 -*-
-
-'''
-FIT
-
-Force spectroscopy curves basic fitting plugin.
-Licensed under the GNU GPL version 2
+# Copyright (C) 2008-2010 Alberto Gomez-Casado
+# Fabrizio Benedetti
+# Massimo Sandal <devicerandom@gmail.com>
+# W. Trevor King <wking@drexel.edu>
+#
+# This file is part of Hooke.
+#
+# Hooke is free software: you can redistribute it and/or
+# modify it under the terms of the GNU Lesser General Public
+# License as published by the Free Software Foundation, either
+# version 3 of the License, or (at your option) any later version.
+#
+# Hooke is distributed in the hope that it will be useful,
+# but WITHOUT ANY WARRANTY; without even the implied warranty of
+# MERCHANTABILITY or FITNESS FOR A PARTICULAR PURPOSE. See the
+# GNU Lesser General Public License for more details.
+#
+# You should have received a copy of the GNU Lesser General Public
+# License along with Hooke. If not, see
+# <http://www.gnu.org/licenses/>.
+
+"""Force spectroscopy curves basic fitting plugin.
Non-standard Dependencies:
procplots.py (plot processing plugin)
-'''
+"""
+
from ..libhooke import WX_GOOD, ClickedPoint
import wxversion
fitplot.colors+=[None,None]
self._send_plot([fitplot])
-
-
-
- #----------
-
-
- def find_contact_point(self,plot=False):
- '''
- Finds the contact point on the curve.
-
- The current algorithm (thanks to Francesco Musiani, francesco.musiani@unibo.it and Massimo Sandal) is:
- - take care of the PicoForce trigger bug - exclude retraction portions with too high standard deviation
- - fit the second half of the retraction curve to a line
- - if the fit is not almost horizontal, take a smaller chunk and repeat
- - otherwise, we have something horizontal
- - so take the average of horizontal points and use it as a baseline
-
- Then, start from the rise of the retraction curve and look at the first point below the
- baseline.
-
- FIXME: should be moved, probably to generalvclamp.py
- '''
-
- if not plot:
- plot=self.plots[0]
-
- outplot=self.subtract_curves(1)
- xret=outplot.vectors[1][0]
- ydiff=outplot.vectors[1][1]
-
- xext=plot.vectors[0][0]
- yext=plot.vectors[0][1]
- xret2=plot.vectors[1][0]
- yret=plot.vectors[1][1]
-
- #taking care of the picoforce trigger bug: we exclude portions of the curve that have too much
- #standard deviation. yes, a lot of magic is here.
- monster=True
- monlength=len(xret)-int(len(xret)/20)
- finalength=len(xret)
- while monster:
- monchunk=scipy.array(ydiff[monlength:finalength])
- if abs(np.std(monchunk)) < 2e-10:
- monster=False
- else: #move away from the monster
- monlength-=int(len(xret)/50)
- finalength-=int(len(xret)/50)
-
-
- #take half of the thing
- endlength=int(len(xret)/2)
-
- ok=False
-
- while not ok:
- xchunk=yext[endlength:monlength]
- ychunk=yext[endlength:monlength]
- regr=scipy.stats.linregress(xchunk,ychunk)[0:2]
- #we stop if we found an almost-horizontal fit or if we're going too short...
- #FIXME: 0.1 and 6 here are "magic numbers" (although reasonable)
- if (abs(regr[1]) > 0.1) and ( endlength < len(xret)-int(len(xret)/6) ) :
- endlength+=10
- else:
- ok=True
-
-
- ymean=np.mean(ychunk) #baseline
-
- index=0
- point = ymean+1
-
- #find the first point below the calculated baseline
- while point > ymean:
- try:
- point=yret[index]
- index+=1
- except IndexError:
- #The algorithm didn't find anything below the baseline! It should NEVER happen
- index=0
- return index
-
- return index
-
-
-
- def find_contact_point2(self, debug=False):
- '''
- TO BE DEVELOPED IN THE FUTURE
- Finds the contact point on the curve.
-
- FIXME: should be moved, probably to generalvclamp.py
- '''
-
- #raw_plot=self.current.curve.default_plots()[0]
- raw_plot=self.plots[0]
- '''xext=self.plots[0].vectors[0][0]
- yext=self.plots[0].vectors[0][1]
- xret2=self.plots[0].vectors[1][0]
- yret=self.plots[0].vectors[1][1]
- '''
- xext=raw_plot.vectors[0][0]
- yext=raw_plot.vectors[0][1]
- xret2=raw_plot.vectors[1][0]
- yret=raw_plot.vectors[1][1]
-
- first_point=[xext[0], yext[0]]
- last_point=[xext[-1], yext[-1]]
-
- #regr=scipy.polyfit(first_point, last_point,1)[0:2]
- diffx=abs(first_point[0]-last_point[0])
- diffy=abs(first_point[1]-last_point[1])
-
- #using polyfit results in numerical errors. good old algebra.
- a=diffy/diffx
- b=first_point[1]-(a*first_point[0])
- baseline=scipy.polyval((a,b), xext)
-
- ysub=[item-basitem for item,basitem in zip(yext,baseline)]
-
- contact=ysub.index(min(ysub))
-
- return xext,ysub,contact
-
- #now, exploit a ClickedPoint instance to calculate index...
- dummy=ClickedPoint()
- dummy.absolute_coords=(x_intercept,y_intercept)
- dummy.find_graph_coords(xret2,yret)
-
- if debug:
- return dummy.index, regr, regr_contact
- else:
- return dummy.index
-
-
-
- def x_do_contact(self,args):
- '''
- DEBUG COMMAND to be activated in the future
- '''
- xext,ysub,contact=self.find_contact_point2(debug=True)
-
- contact_plot=self.plots[0]
- contact_plot.add_set(xext,ysub)
- contact_plot.add_set([xext[contact]],[self.plots[0].vectors[0][1][contact]])
- #contact_plot.add_set([first_point[0]],[first_point[1]])
- #contact_plot.add_set([last_point[0]],[last_point[1]])
- contact_plot.styles=[None,None,None,'scatter']
- self._send_plot([contact_plot])
- return
-
-
- index,regr,regr_contact=self.find_contact_point2(debug=True)
- print regr
- print regr_contact
- raw_plot=self.current.curve.default_plots()[0]
- xret=raw_plot.vectors[0][0]
- #nc_line=[(item*regr[0])+regr[1] for item in x_nc]
- nc_line=scipy.polyval(regr,xret)
- c_line=scipy.polyval(regr_contact,xret)
-
-
- contact_plot=self.current.curve.default_plots()[0]
- contact_plot.add_set(xret, nc_line)
- contact_plot.add_set(xret, c_line)
- contact_plot.styles=[None,None,None,None]
- #contact_plot.styles.append(None)
- contact_plot.destination=1
- self._send_plot([contact_plot])
projects / hooke.git / blobdiff