1 # -*- coding: utf-8 -*-
3 # Copyright (C) 2008-2010 Alberto Gomez-Casado
5 # Massimo Sandal <devicerandom@gmail.com>
6 # W. Trevor King <wking@drexel.edu>
8 # This file is part of Hooke.
10 # Hooke is free software: you can redistribute it and/or
11 # modify it under the terms of the GNU Lesser General Public
12 # License as published by the Free Software Foundation, either
13 # version 3 of the License, or (at your option) any later version.
15 # Hooke is distributed in the hope that it will be useful,
16 # but WITHOUT ANY WARRANTY; without even the implied warranty of
17 # MERCHANTABILITY or FITNESS FOR A PARTICULAR PURPOSE. See the
18 # GNU Lesser General Public License for more details.
20 # You should have received a copy of the GNU Lesser General Public
21 # License along with Hooke. If not, see
22 # <http://www.gnu.org/licenses/>.
24 """The ``convfilt`` module provides :class:`ConvFiltPlugin` and
25 for convolution-based filtering of :mod:`hooke.curve.Curve`\s.
29 :mod:`~hooke.plugin.flatfilt` for a collection of additional filters.
31 :class:`ConfFiltPlugin` was broken out of
32 :class:`~hooke.plugin.flatfilt` to separate its large number of
33 configuration settings.
37 from multiprocessing import Queue
41 from ..command import Command, Argument, Success, Failure
42 from ..config import Setting
43 from ..experiment import VelocityClamp
44 from ..plugin import Plugin, argument_to_setting
45 from ..plugin.curve import CurveArgument
46 from ..plugin.playlist import FilterCommand
47 from ..plugin.vclamp import scale
48 from ..util.fit import PoorFit
49 from ..util.peak import find_peaks, find_peaks_arguments, Peak, _kwargs
52 class ConvFiltPlugin (Plugin):
53 """Convolution-based peak recognition and filtering.
56 super(ConvFiltPlugin, self).__init__(name='convfilt')
57 self._arguments = [ # for Command initialization
58 Argument('convolution', type='float', count='-1',
59 default=[11.0]+[-1.0]*11, help="""
60 Convolution vector roughly matching post-peak cantilever rebound.
61 This should roughly match the shape of the feature you're looking for.
62 """.strip()), # TODO: per-curve convolution vector + interpolation, to
63 # take advantage of the known spring constant.
64 Argument('blind window', type='float', default=20e-9, help="""
65 Meters after the contact point where we do not count peaks to avoid
66 non-specific surface interaction.
68 Argument('min peaks', type='int', default=5, help="""
69 Minimum number of peaks for curve acceptance.
71 ] + copy.deepcopy(find_peaks_arguments)
72 # Set convolution-specific defaults for the fit_peak_arguments.
73 for key,value in [('cut side', 'positive'),
76 ('min deviations', 5.0),
78 ('see double', 10e-9),
80 argument = [a for a in self._arguments if a.name == key][0]
81 argument.default = value
83 Setting(section=self.setting_section, help=self.__doc__)]
84 for argument in self._arguments:
85 self._settings.append(argument_to_setting(
86 self.setting_section, argument))
87 argument.default = None # if argument isn't given, use the config.
88 self._commands = [ConvolutionPeaksCommand(self),
89 ConvolutionFilterCommand(self)]
91 def dependencies(self):
94 def default_settings(self):
98 class ConvolutionPeaksCommand (Command):
99 """Detect peaks in velocity clamp data with a convolution.
103 A simplified version of Kasas' filter [#kasas2000].
104 For any given retracting curve:
106 1) The contact point is found.
107 2) The off-surface data is convolved (using :func:`numpy.convolve`)
108 with a vector that encodes the approximate shape of the target
110 3) Peaks in the convolved curve are extracted with
111 :func:`~hooke.plugins.peak.find_peaks`.
113 The convolution algorithm, with appropriate thresholds, usually
114 recognizes peaks well more than 95% of the time.
116 .. [#kasas2000] S. Kasas, B.M. Riederer, S. Catsicas, B. Cappella,
118 "Fuzzy logic algorithm to extract specific interaction forces
119 from atomic force microscopy data"
120 Rev. Sci. Instrum., 2000.
121 doi: `10.1063/1.1150583 <http://dx.doi.org/10.1063/1.1150583>`_
123 def __init__(self, plugin):
124 config_arguments = [a for a in plugin._arguments
125 if a.name != 'min peaks']
126 # Drop min peaks, since we're not filtering with this
127 # function, just detecting peaks.
128 super(ConvolutionPeaksCommand, self).__init__(
129 name='convolution peaks',
130 arguments=[CurveArgument] + config_arguments,
131 help=self.__doc__, plugin=plugin)
133 def _run(self, hooke, inqueue, outqueue, params):
134 z_data,d_data,params = self._setup(hooke, params)
136 while z_data[start_index] < params['bind window']:
138 conv = numpy.convolve(d_data[start_index:], params['convolution'],
140 peaks = find_peaks(conv, **_kwargs(params, find_peaks_arguments,
141 argument_input_keys=True))
143 peak.name = 'convolution of %s with %s' \
144 % (params['deflection column name'], params['convolution'])
145 peak.index += start_index
148 def _setup(self, hooke, params):
149 """Setup `params` from config and return the z piezo and
152 curve = params['curve']
153 if curve.info['experiment'] != VelocityClamp:
154 raise Failure('%s operates on VelocityClamp experiments, not %s'
155 % (self.name, curve.info['experiment']))
156 for col in ['surface distance (m)', 'deflection (N)']:
157 if col not in curve.data[0].info['columns']:
160 for block in curve.data:
161 if block.info['name'].startswith('retract'):
165 raise Failure('No retraction blocks in %s.' % curve)
166 z_data = data[:,data.info['columns'].index('surface distance (m)')]
167 if 'flattened deflection (N)' in data.info['columns']:
168 params['deflection column name'] = 'flattened deflection (N)'
170 params['deflection column name'] = 'deflection (N)'
171 d_data = data[:,data.info['columns'].index(
172 params['deflection column name'])]
173 for key,value in params.items():
174 if value == None: # Use configured default value.
175 params[key] = self.plugin.config[key]
176 return z_data,d_data,params
179 class ConvolutionFilterCommand (FilterCommand):
180 u"""Create a subset playlist of curves with enough convolution peaks.
184 This filter can reduce a dataset like the one in [#brucale2009] to
185 analyze by hand by about 80-90% (depending on the overall
186 cleanliness of the data set). Thousands of curves can be
187 automatically filtered this way in a few minutes on a standard PC,
188 but the algorithm could still be optimized.
190 .. [#brucale2009] M. Brucale, M. Sandal, S. Di Maio, A. Rampioni,
191 I. Tessari, L. Tosatto, M. Bisaglia, L. Bubacco, B. Samorì.
192 "Pathogenic mutations shift the equilibria of
193 α-Synuclein single molecules towards structured
196 doi: `10.1002/cbic.200800581 <http://dx.doi.org/10.1002/cbic.200800581>`_
200 ConvolutionCommand : Underlying convolution-based peak detection.
202 def __init__(self, plugin):
203 super(ConvolutionFilterCommand, self).__init__(
204 plugin, name='convolution filter playlist')
205 self.arguments.extend(plugin._arguments)
207 def filter(self, curve, hooke, inqueue, outqueue, params):
208 params['curve'] = curve
211 conv_command = [c for c in self.hooke.commands
212 if c.name=='convolution peaks'][0]
213 conv_command.run(hooke, inq, outq, **params)
215 if isinstance(peaks, UncaughtException) \
216 and isinstance(peaks.exception, PoorFit):
218 if not (isinstance(peaks, list) and (len(peaks) == 0
219 or isinstance(peaks[0], Peak))):
220 raise Failure('Expected a list of Peaks, not %s' % peaks)
222 if not isinstance(ret, Success):
224 if params['min peaks'] == None: # Use configured default value.
225 params['min peaks'] = self.plugin.config['min peaks']
226 return len(peaks) >= params['min peaks']