import copy
import numpy
+import scipy
from ..command import Command, Argument, Failure, NullQueue
from ..config import Setting
if r2 >= 1:
self._model_data[:r2] = p[0] + p[1] * numpy.arange(r2)
if r2 < len(self._data)-1:
- self._model_data[r2:] = p[0] + p[1]*p[2] + p[3] * numpy.arange(len(self._data)-r2)
+ self._model_data[r2:] = \
+ p[0] + p[1]*p[2] + p[3] * numpy.arange(len(self._data)-r2)
return self._model_data
def set_data(self, data, info=None):
params.append(0.) # restore the non-contact slope parameter
# check that the fit is reasonable, see the :meth:`model` docstring
- # for parameter descriptions.
+ # for parameter descriptions. HACK: hardcoded cutoffs.
if abs(params[3]*10) > abs(params[1]) :
raise PoorFit('Slope in non-contact region, or no slope in contact')
if params[2] < self.info['min position']+0.02*self.info['position range']:
The adjusted data columns `surface distance (m)` and `surface
adjusted deflection (m)` are also added to the block.
+
+ You can select the contact point algorithm with the creatively
+ named `surface contact point algorithm` configuration setting.
+ Currently available options are:
+
+ * fmms (:meth:`find_contact_point_fmms`)
+ * ms (:meth:`find_contact_point_ms`)
+ * wtk (:meth:`find_contact_point_wtk`)
"""
def __init__(self, plugin):
super(SurfaceContactCommand, self).__init__(
['surface distance (m)', 'surface adjusted deflection (m)'])
z_data = data[:,data.info['columns'].index('z piezo (m)')]
d_data = data[:,data.info['columns'].index('deflection (m)')]
- i,deflection_offset,ps = self._find_contact_point(
- z_data, d_data, outqueue)
+ i,deflection_offset,ps = self.find_contact_point(
+ params['curve'], z_data, d_data, outqueue)
surface_offset = z_data[i]
new.info['surface distance offset (m)'] = surface_offset
new.info['surface deflection offset (m)'] = deflection_offset
data = params['curve'].data[int(params['block'])] # HACK, int() should be handled by ui
params['curve'].data[int(params['block'])] = new # HACK, int() should be handled by ui
- def _find_contact_point(self, z_data, d_data, outqueue=None):
- fn = getattr(self, '_find_contact_point_%s'
+ def find_contact_point(self, curve, z_data, d_data, outqueue=None):
+ """Railyard for the `find_contact_point_*` family.
+
+ Uses the `surface contact point algorithm` configuration
+ setting to call the appropriate backend algorithm.
+ """
+ fn = getattr(self, 'find_contact_point_%s'
% self.plugin.config['surface contact point algorithm'])
- return fn(z_data, d_data, outqueue)
+ return fn(curve, z_data, d_data, outqueue)
- def _find_contact_point_fmms(self, z_data, d_data, outqueue=None):
+ def find_contact_point_fmms(self, curve, z_data, d_data, outqueue=None):
"""Algorithm by Francesco Musiani and Massimo Sandal.
Notes
-----
Algorithm:
- * 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.
+ 0) Driver-specific workarounds, e.g. deal with the PicoForce
+ trigger bug by excluding retraction portions with excessive
+ deviation.
+ 1) Select the second half (non-contact side) of the retraction
+ curve.
+ 2) Fit the selection to a line.
+ 3) If the fit is not almost horizontal, halve the selection
+ and retrun to (2).
+ 4) Average the selection and use it as a baseline.
+ 5) Slide in from the start (contact side) of the retraction
+ curve, until you find a point with greater than baseline
+ deflection. That point is the contact point.
"""
- 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_point_ms(self, z_data, d_data, outqueue=None):
+ if curve.info['filetype'] == 'picoforce':
+ # Take care of the picoforce trigger bug (TODO: example
+ # data file demonstrating the bug). We exclude portions
+ # of the curve that have too much standard deviation.
+ # Yes, a lot of magic is here.
+ check_start = len(d_data)-len(d_data)/20
+ monster_start = len(d_data)
+ while True:
+ # look at the non-contact tail
+ non_monster = d_data[check_start:monster_start]
+ if non_monster.std() < 2e-10: # HACK: hardcoded cutoff
+ break
+ else: # move further away from the monster
+ check_start -= len(d_data)/50
+ monster_start -= len(d_data)/50
+ z_data = z_data[:monster_start]
+ d_data = d_data[:monster_start]
+
+ # take half of the thing to start
+ selection_start = len(d_data)/2
+ while True:
+ z_chunk = z_data[selection_start:]
+ d_chunk = d_data[selection_start:]
+ slope,intercept,r,two_tailed_prob,stderr_of_the_estimate = \
+ scipy.stats.linregress(z_chunk, d_chunk)
+ # We stop if we found an almost-horizontal fit or if we're
+ # getting to small a selection. FIXME: 0.1 and 5./6 here
+ # are "magic numbers" (although reasonable)
+ if (abs(slope) < 0.1 # deflection (m) / surface (m)
+ or selection_start > 5./6*len(d_data)):
+ break
+ selection_start += 10
+
+ d_baseline = d_chunk.mean()
+
+ # find the first point above the calculated baseline
+ i = 0
+ while i < len(d_data) and d_data[i] < ymean:
+ i += 1
+ return (i, d_baseline, {})
+
+ def find_contact_point_ms(self, curve, z_data, d_data, outqueue=None):
"""Algorithm by Massimo Sandal.
Notes
* ?
"""
- #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]
else:
return dummy.index
- def _find_contact_point_wtk(self, z_data, d_data, outqueue=None):
+ def find_contact_point_wtk(self, curve, z_data, d_data, outqueue=None):
"""Algorithm by W. Trevor King.
Notes
projects / hooke.git / commitdiff