1 # calibcant - tools for thermally calibrating AFM cantilevers
3 # Copyright (C) 2008-2012 W. Trevor King <wking@drexel.edu>
5 # This file is part of calibcant.
7 # calibcant is free software: you can redistribute it and/or modify it under
8 # the terms of the GNU General Public License as published by the Free Software
9 # Foundation, either version 3 of the License, or (at your option) any later
12 # calibcant is distributed in the hope that it will be useful, but WITHOUT ANY
13 # WARRANTY; without even the implied warranty of MERCHANTABILITY or FITNESS FOR
14 # A PARTICULAR PURPOSE. See the GNU General Public License for more details.
16 # You should have received a copy of the GNU General Public License along with
17 # calibcant. If not, see <http://www.gnu.org/licenses/>.
19 """Acquire and analyze cantilever calibration data.
21 The relevent physical quantities are:
23 * Vzp_out Output z-piezo voltage (what we generate)
24 * Vzp Applied z-piezo voltage (after external ZPGAIN)
25 * Zp The z-piezo position
26 * Zcant The cantilever vertical deflection
27 * Vphoto The photodiode vertical deflection voltage (what we measure)
28 * Fcant The force on the cantilever
29 * T The temperature of the cantilever and surrounding solution
30 * (another thing we measure or guess)
31 * k_b Boltzmann's constant
33 Which are related by the parameters:
35 * zpGain Vzp_out / Vzp
36 * zpSensitivity Zp / Vzp
37 * photoSensitivity Vphoto / Zcant
38 * k_cant Fcant / Zcant
40 Cantilever calibration assumes a pre-calibrated z-piezo
41 (zpSensitivity) and a amplifier (zpGain). In our lab, the z-piezo is
42 calibrated by imaging a calibration sample, which has features with
43 well defined sizes, and the gain is set with a knob on the Nanoscope.
45 photoSensitivity is measured by bumping the cantilever against the
46 surface, where Zp = Zcant (see bump_acquire() and the bump_analyze
49 k_cant is measured by watching the cantilever vibrate in free solution
50 (see the vib_acquire() and the vib_analyze submodule). The average
51 energy of the cantilever in the vertical direction is given by the
52 equipartition theorem.
54 .. math:: \frac{1}{2} k_b T = \frac{1}{2} k_cant <Zcant**2>
58 .. math:: k_cant = \frac{k_b T}{Zcant**2}
62 .. math:: Zcant = \frac{Vphoto}{photoSensitivity}
66 .. math:: k_cant = \frac{k_b T * photoSensitivty^2}{<Vphoto**2>}
68 We measured photoSensitivity with the surface bumps. We can either
69 measure T using an external function (see temperature.py), or just
70 estimate it (see T_acquire() and the T_analyze submodule). Guessing
71 room temp ~22 deg C is actually fairly reasonable. Assuming the
72 actual fluid temperature is within +/- 5 deg, the error in the spring
73 constant k_cant is within 5/273.15 ~= 2%. A time series of Vphoto
74 while we're far from the surface and not changing Vzp_out will give us
75 the average variance <Vphoto**2>.
77 We do all these measurements a few times to estimate statistical
80 The functions are layed out in the families::
82 bump_*(), vib_*(), T_*(), and calib_*()
84 For each family, * can be any of:
86 * acquire get real-world data
87 * save store real-world data to disk
88 * load get real-world data from disk
89 * analyze interperate the real-world data.
90 * plot show a nice graphic to convince people we're working :p
92 A family name without any `_*` extension (e.g. `bump()`), runs
93 `*_acquire()`, `*_analyze()`, and `*_save()`. `*_analyze()` will run
94 `*_plot()` if `matplotlib` is set in `calibcant.package_config`.
97 from numpy import zeros as _zeros
98 from numpy import float as _float
99 from time import sleep as _sleep
101 from . import LOG as _LOG
103 from .bump import bump as _bump
104 from .T import T as _T
105 from .vib import vib as _vib
106 from .analyze import calib_analyze as _calib_analyze
107 from .analyze import calib_save as _calib_save
110 def move_far_from_surface(stepper, distance):
111 """Step back approximately `distance` meters.
113 steps = int(distance/stepper.step_size)
114 _LOG.info('step back %d steps (~%g m)' % (steps, distance))
115 stepper.step_relative(-steps)
117 def calib_acquire(afm, calibration_config, filename=None, group='/'):
118 """Acquire data for calibrating a cantilever in one function.
121 afm a pyafm.AFM instance
122 calibration_config a .config._CalibrationConfig instance
124 Outputs (all are arrays of recorded data):
125 bumps measured (V_photodiode / nm_tip) proportionality constant
126 Ts measured temperature (K)
127 vibs measured V_photodiode variance (Volts**2) in free solution
129 The temperatures are collected after moving far from the surface
130 but before and vibrations are measured to give everything time to
131 settle after the big move.
133 assert group.endswith('/'), group
135 bumps = _zeros((calibration_config['num-bumps'],), dtype=_float)
136 for i in range(calibration_config['num-bumps']):
137 _LOG.info('acquire bump %d of %d' % (i, calibration_config['num-bumps']))
138 bumps[i] = _bump(afm=afm, bump_config=calibration_config['bump'],
139 filename=filename, group='%sbump/%d/' % (group, i))
140 _LOG.debug('bumps: %s' % bumps)
142 move_far_from_surface(
143 afm.stepper, distance=calibration_config['vibration-spacing'])
145 Ts = _zeros((calibration_config['num-temperatures'],), dtype=_float)
146 for i in range(calibration_config['num-temperatures']):
147 _LOG.info('acquire T %d of %d'
148 % (i, calibration_config['num-temperatures']))
150 get_T=afm.get_temperature,
151 temperature_config=calibration_config['temperature'],
152 filename=filename, group='%stemperature/%d/' % (group, i))
153 _sleep(calibration_config['temperature-sleep'])
154 _LOG.debug('temperatures: %s' % Ts)
157 vibs = _zeros((calibration_config['num-vibrations'],), dtype=_float)
158 for i in range(calibration_config['num-vibrations']):
160 piezo=afm.piezo, vibration_config=calibration_config['vibration'],
161 filename=filename, group='%svibration/%d/' % (group, i))
162 _LOG.debug('vibrations: %s' % vibs)
164 return (bumps, Ts, vibs)
166 def calib(afm, calibration_config, filename=None, group='/'):
167 """Calibrate a cantilever in one function.
170 (see `calib_acquire()`)
173 k cantilever spring constant (in N/m, or equivalently nN/nm)
174 k_s standard deviation in our estimate of k
177 >>> from pprint import pprint
179 >>> from h5config.storage.hdf5 import pprint_HDF5
180 >>> from pycomedi.device import Device
181 >>> from pycomedi.subdevice import StreamingSubdevice
182 >>> from pycomedi.channel import AnalogChannel, DigitalChannel
183 >>> from pycomedi.constant import AREF, IO_DIRECTION, SUBDEVICE_TYPE, UNIT
184 >>> from pypiezo.afm import AFMPiezo
185 >>> from pypiezo.base import PiezoAxis, InputChannel
186 >>> from pypiezo.config import ChannelConfig, AxisConfig
187 >>> from stepper import Stepper
188 >>> from pyafm.afm import AFM
189 >>> from .config import (CalibrationConfig, BumpConfig,
190 ... TemperatureConfig, VibrationConfig)
191 >>> from .analyze import calib_load_all
193 >>> fd,filename = tempfile.mkstemp(suffix='.h5', prefix='calibcant-')
196 >>> d = Device('/dev/comedi0')
199 Setup an `AFMPiezo` instance.
201 >>> s_in = d.find_subdevice_by_type(SUBDEVICE_TYPE.ai,
202 ... factory=StreamingSubdevice)
203 >>> s_out = d.find_subdevice_by_type(SUBDEVICE_TYPE.ao,
204 ... factory=StreamingSubdevice)
206 >>> axis_channel = s_out.channel(
207 ... 0, factory=AnalogChannel, aref=AREF.ground)
208 >>> input_channel = s_in.channel(0, factory=AnalogChannel, aref=AREF.diff)
209 >>> for chan in [axis_channel, input_channel]:
210 ... chan.range = chan.find_range(unit=UNIT.volt, min=-10, max=10)
212 We set the minimum voltage for the `z` axis to -9 (a volt above
213 the minimum possible voltage) to help with testing
214 `.get_surface_position`. Without this minimum voltage, small
215 calibration errors could lead to a railed -10 V input for the
216 first few surface approaching steps, which could lead to an
217 `EdgeKink` error instead of a `FlatFit` error.
219 >>> axis_config = AxisConfig()
220 >>> axis_config.update(
221 ... {'gain':20, 'sensitivity':8e-9, 'minimum':-9})
222 >>> axis_channel_config = ChannelConfig()
223 >>> axis_channel_config['name'] = 'z'
224 >>> axis_config['channel'] = axis_channel_config
225 >>> input_channel_config = ChannelConfig()
226 >>> input_channel_config['name'] = 'deflection'
228 >>> a = PiezoAxis(config=axis_config, axis_channel=axis_channel)
231 >>> c = InputChannel(config=input_channel_config, channel=input_channel)
234 >>> piezo = AFMPiezo(axes=[a], inputs=[c])
236 Setup a `stepper` instance.
238 >>> s_d = d.find_subdevice_by_type(SUBDEVICE_TYPE.dio)
239 >>> d_channels = [s_d.channel(i, factory=DigitalChannel)
240 ... for i in (0, 1, 2, 3)]
241 >>> for chan in d_channels:
242 ... chan.dio_config(IO_DIRECTION.output)
244 >>> def write(value):
245 ... s_d.dio_bitfield(bits=value, write_mask=2**4-1)
247 >>> stepper = Stepper(write=write)
249 Setup an `AFM` instance.
251 >>> afm = AFM(piezo, stepper)
255 >>> calibration_config = CalibrationConfig()
256 >>> calibration_config['bump'] = BumpConfig()
257 >>> calibration_config['temperature'] = TemperatureConfig()
258 >>> calibration_config['vibration'] = VibrationConfig()
259 >>> calib(afm, calibration_config, filename=filename, group='/')
260 TODO: replace skipped example data with real-world values
261 >>> pprint_HDF5(filename) # doctest: +ELLIPSIS, +REPORT_UDIFF
267 <HDF5 dataset "far-steps": shape (), type "<i4">
270 /bump/0/config/deflection
271 /bump/0/config/deflection/channel
272 <HDF5 dataset "channel": shape (), type "<i4">
276 /bump/0/config/z/axis
277 /bump/0/config/z/axis/channel
278 <HDF5 dataset "channel": shape (), type "<i4">
281 <HDF5 dataset "gain": shape (), type "<i4">
284 <HDF5 dataset "processed": shape (), type "<f8">
287 <HDF5 dataset "deflection": shape (2048,), type "<u2">
289 <HDF5 dataset "z": shape (2048,), type "<u2">
295 /calibration/config/bump
296 <HDF5 dataset "far-steps": shape (), type "<i4">
299 <HDF5 dataset "num-bumps": shape (), type "<i4">
302 /calibration/processed
303 /calibration/processed/spring-constant
304 <HDF5 dataset "data": shape (), type "<f8">
306 <HDF5 dataset "standard-deviation": shape (), type "<f8">
308 <HDF5 dataset "units": shape (), type "|S3">
311 /calibration/raw/photodiode-sensitivity
312 <HDF5 dataset "data": shape (10,), type "<f8">
314 <HDF5 dataset "units": shape (), type "|S3">
316 /calibration/raw/temperature
317 <HDF5 dataset "data": shape (10,), type "<f8">
319 <HDF5 dataset "units": shape (), type "|S1">
321 /calibration/raw/thermal-vibration-variance
322 <HDF5 dataset "data": shape (20,), type "<f8">
324 <HDF5 dataset "units": shape (), type "|S3">
328 /temperature/0/config
329 <HDF5 dataset "default": shape (), type "|b1">
331 <HDF5 dataset "units": shape (), type "|S7">
333 <HDF5 dataset "processed": shape (), type "<f8">
335 <HDF5 dataset "raw": shape (), type "<i4">
342 /vibration/0/config/deflection
343 <HDF5 dataset "channel": shape (), type "<i4">
346 /vibration/0/config/vibration
347 <HDF5 dataset "chunk-size": shape (), type "<i4">
350 <HDF5 dataset "processed": shape (), type "<f8">
353 <HDF5 dataset "deflection": shape (65536,), type "<u2">
360 <HDF5 dataset "deflection": shape (65536,), type "<u2">
362 >>> everything = calib_load_all(filename, '/')
363 >>> pprint(everything) # doctest: +ELLIPSIS, +REPORT_UDIFF
364 {'bump_details': [{'bump_config': <BumpConfig ...>,
365 'deflection_channel_config': <ChannelConfig ...>,
366 'processed_bump': ...,
367 'raw_bump': {'deflection': array([...], dtype=uint16),
368 'z': array([...], dtype=uint16)},
369 'z_axis_config': <AxisConfig ...>},
371 'bumps': array([...]),
372 'calibration_config': <CalibrationConfig ...>,
375 'temperature_details': [{'processed_temperature': ...,
376 'raw_temperature': array(22),
377 'temperature_config': <TemperatureConfig ...>},
379 'temperatures': array([...]),
380 'vibration_details': [{'deflection_channel_config': <ChannelConfig ...>,
381 'processed_vibration': ...,
382 'raw_vibration': array([...], dtype=uint16),
383 'vibration_config': <VibrationConfig ...>},
385 'vibrations': array([...])}
387 Close the Comedi device.
391 Cleanup our temporary config file.
395 bumps, Ts, vibs = calib_acquire(
396 afm, calibration_config, filename=filename, group=group)
397 # TODO: convert vib units?
398 k,k_s = _calib_analyze(bumps, Ts, vibs)
399 _calib_save(filename, group=group+'calibration/', bumps=bumps,
400 temperatures=Ts, vibrations=vibs,
401 calibration_config=calibration_config, k=k, k_s=k_s)