3 # calibcant - tools for thermally calibrating AFM cantilevers
5 # Copyright (C) 2007,2008, William Trevor King
7 # This program is free software; you can redistribute it and/or
8 # modify it under the terms of the GNU General Public License as
9 # published by the Free Software Foundation; either version 3 of the
10 # License, or (at your option) any later version.
12 # This program is distributed in the hope that it will be useful, but
13 # WITHOUT ANY WARRANTY; without even the implied warranty of
14 # MERCHANTABILITY or FITNESS FOR A PARTICULAR PURPOSE.
15 # See the GNU General Public License for more details.
17 # You should have received a copy of the GNU General Public License
18 # along with this program; if not, write to the Free Software
19 # Foundation, Inc., 59 Temple Place - Suite 330, Boston, MA
22 # The author may be contacted at <wking@drexel.edu> on the Internet, or
23 # write to Trevor King, Drexel University, Physics Dept., 3141 Chestnut St.,
24 # Philadelphia PA 19104, USA.
27 Separate the more general calib_analyze() from the other calib_*()
28 functions in calibcant. Also provide a command line interface
29 for analyzing data acquired through other workflows.
31 The relevent physical quantities are :
32 Vzp_out Output z-piezo voltage (what we generate)
33 Vzp Applied z-piezo voltage (after external ZPGAIN)
34 Zp The z-piezo position
35 Zcant The cantilever vertical deflection
36 Vphoto The photodiode vertical deflection voltage (what we measure)
37 Fcant The force on the cantilever
38 T The temperature of the cantilever and surrounding solution
39 (another thing we measure)
40 k_b Boltzmann's constant
42 Which are related by the parameters :
44 zpSensitivity Zp / Vzp
45 photoSensitivity Vphoto / Zcant
52 from splittable_kwargs import splittableKwargsFunction, \
53 make_splittable_kwargs_function
60 kb = 1.3806504e-23 # Boltzmann's constant in J/K
62 #@splittableKwargsFunction((calib_plot, 'bumps', 'Ts', 'vibs'))
63 # Some of the child functions aren't yet defined, so postpone
64 # make-splittable until later in the module.
65 def calib_analyze(bumps, Ts, vibs, **kwargs) :
67 Analyze data from get_calibration_data()
68 return (k, k_s, !!!a2_r, T_r, one_o_Vp2_r)
69 Inputs (all are arrays of recorded data) :
70 bumps measured (V_photodiode / nm_tip) proportionality constant
71 Ts measured temperature (K)
72 vibs measured V_photodiode variance in free solution (V**2)
74 k cantilever spring constant (in N/m, or equivalently nN/nm)
75 k_s standard deviation in our estimate of k
76 !!!a2_r relative error in a**2
77 !!!T_r relative error in T
78 !!!one_o_Vp2_r relative error in 1/Vphotodiode_variance
80 We're assuming vib is mostly from thermal cantilever vibrations
81 (and then only from vibrations in the single vertical degree of freedom),
82 and not from other noise sources.
83 The various relative errors are returned to help you gauge the
84 largest source of random error in your measurement of k.
85 If one of them is small, don't bother repeating that measurment too often.
86 If one is large, try repeating that measurement more.
87 Remember that you need enough samples to have a valid error estimate in
88 the first place, and that none of this addresses any systematic errors.
90 calib_plot_kwargs, = calib_analyze._splitargs(calib_analyze, kwargs)
91 photoSensitivity2 = bumps**2
92 one_o_Vphoto2 = 1/vibs
94 photoSensitivity2_m = photoSensitivity2.mean()
96 one_o_Vphoto2_m = one_o_Vphoto2.mean()
97 # Vphoto / photoSensitivity = x
98 # k = kb T / <x**2> = kb T photoSensitiviy**2 * (1e9nm/m)**2 /
101 # units, photoSensitivity = Vphoto/(Zcant in nm),
102 # so Vphoto/photoSensitivity = Zcant in nm
103 # so k = J/K * K / nm^2 * (1e9nm/m)**2 = N/m
104 k = kb * T_m * photoSensitivity2_m * one_o_Vphoto2_m * 1e18
106 # propogation of errors !!!
107 # first, get standard deviations
108 photoSensitivity2_s = photoSensitivity2.std()
110 one_o_Vphoto2_s = one_o_Vphoto2.std()
111 # !!!!now, get relative errors
112 photoSensitivity2_r = photoSensitivity2_s / photoSensitivity2_m
114 one_o_Vphoto2_r = one_o_Vphoto2_s / one_o_Vphoto2_m
116 k_s = k*(photoSensitivity2_r**2 + T_r**2 + one_o_Vphoto2_r**2)**0.5
118 calib_plot(bumps, Ts, vibs, **calib_plot_kwargs)
121 photoSensitivity2_m, photoSensitivity2_s,
122 T_m, T_s, one_o_Vphoto2_m, one_o_Vphoto2_s)
124 @splittableKwargsFunction()
125 def string_errors(k_m, k_s,
126 photoSensitivity2_m, photoSensitivity2_s,
128 one_o_Vphoto2_m, one_o_Vphoto2_s) :
130 photoSensitivity2_r = photoSensitivity2_s / photoSensitivity2_m
132 one_o_Vphoto2_r = one_o_Vphoto2_s / one_o_Vphoto2_m
133 string = "Variable (units) : mean +/- std. dev. (relative error)\n"
134 string += "Cantilever k (N/m) : %g +/- %g (%g)\n" \
136 string += "photoSensitivity**2 (V/nm)**2 : %g +/- %g (%g)\n" \
137 % (photoSensitivity2_m, photoSensitivity2_s, photoSensitivity2_r)
138 string += "T (K) : %g +/- %g (%g)\n" \
140 string += "1/Vphoto**2 (1/V)**2 : %g +/- %g (%g)" \
141 % (one_o_Vphoto2_m, one_o_Vphoto2_s, one_o_Vphoto2_r)
144 @splittableKwargsFunction()
145 def calib_save(bumps, Ts, vibs, log_dir=None) :
147 Save a dictonary with the bump, T, and vib data.
150 data = {'bump':bumps, 'T':Ts, 'vib':vibs}
151 log = data_logger.data_log(log_dir, noclobber_logsubdir=False,
153 log.write_dict_of_arrays(data)
155 def calib_load(datafile) :
156 "Load the dictionary data, using data_logger.date_load()"
157 dl = data_logger.data_load()
158 data = dl.read_dict_of_arrays(path)
159 return (data['bump'], data['T'], data['vib'])
161 def calib_save_analysis(k, k_s,
162 photoSensitivity2_m, photoSensitivity2_s,
163 T_m, T_s, one_o_Vphoto2_m, one_o_Vphoto2_s,
166 log = data_logger.data_log(log_dir, noclobber_logsubdir=False,
167 log_name="calib_analysis_text")
168 log.write_binary(string_errors(k, k_s,
169 photoSensitivity2_m, photoSensitivity2_s,
171 one_o_Vphoto2_m, one_o_Vphoto2_s)
174 @splittableKwargsFunction()
175 def calib_plot(bumps, Ts, vibs, plotVerbose=False) :
176 if plotVerbose or config.PYLAB_VERBOSE :
177 common._import_pylab()
178 common._pylab.figure(config.BASE_FIGNUM+4)
179 common._pylab.subplot(311)
180 common._pylab.plot(bumps, 'g.-')
181 common._pylab.title('Photodiode sensitivity (V/nm)')
182 common._pylab.subplot(312)
183 common._pylab.plot(Ts, 'r.-')
184 common._pylab.title('Temperature (K)')
185 common._pylab.subplot(313)
186 common._pylab.plot(vibs, 'b.-')
187 common._pylab.title('Thermal deflection variance (Volts**2)')
190 make_splittable_kwargs_function(calib_analyze,
191 (calib_plot, 'bumps', 'Ts', 'vibs'))
193 @splittableKwargsFunction((calib_analyze, 'bumps', 'Ts', 'vibs'))
194 def calib_load_analyze_tweaks(bump_tweaks, vib_tweaks, T_tweaks=None) :
195 raise NotImplementedError
196 a = read_tweaked_bumps(bump_tweaks)
197 vib = V_photo_variance_from_file(vib_tweaks)
200 return analyze_calibration_data(a, T, vib, log_dir=log_dir)
202 # commandline interface functions
205 def array_from_string(string):
207 for num in string.split(',') :
208 ret.append(float(num))
210 return numpy.array(ret)
212 def read_data(ifile):
213 "ifile can be a filename string or open (seekable) file object"
214 unlabeled_data=scipy.io.read_array(file)
215 return unlabeled_data
217 def get_array(string, filename, name) :
218 "get an array from supplied command line options"
220 array = array_from_string(string)
221 elif filename != None :
222 array = read_data(filename)
224 raise Exception, "no %s information given" % (name)
227 if __name__ == '__main__' :
228 # command line interface
229 from optparse import OptionParser
231 usage_string = ('%prog <bumps> <temps> <vibs>\n'
232 '2008, W. Trevor King.\n'
234 'Takes arrays of Vphotodiode sensitivity (V/nm), Temperature (K), \n'
235 'and Vibration variance (V**2) as comma seperated lists.\n'
236 'Returns the cantilever spring constant (pN/nm).\n'
238 ' $ %prog -b 0.02,0.03,0.025 -t 298.2,300.1 -v 6e-9,5.5e-9\n'
240 parser = OptionParser(usage=usage_string, version='%prog '+common.VERSION)
241 parser.add_option('-b', '--bump-string', dest='bump_string',
242 help='comma seperated photodiode sensitivities (V/nm)',
243 type='string', metavar='BUMPS')
244 parser.add_option('-t', '--temp-string', dest='temp_string',
245 help='comma seperated temperatures (K)',
246 type='string', metavar='TEMPS')
247 parser.add_option('-v', '--vib-string', dest='vib_string',
248 help='comma seperated vibration variances (V**2)',
249 type='string', metavar='VIBS')
250 parser.add_option('-B', '--bump-file', dest='bump_file',
251 help='comma seperated photodiode sensitivities (V/nm)',
252 type='string', metavar='BUMPFILE')
253 parser.add_option('-T', '--temp-file', dest='temp_file',
254 help='comma seperated temperatures (K)',
255 type='string', metavar='TEMPFILE')
256 parser.add_option('-V', '--vib-file', dest='vib_file',
257 help='comma seperated vibration variances (V**2)',
258 type='string', metavar='VIBFILE')
259 parser.add_option('-C', '--celsius', dest='celsius',
260 help='Use Celsius input temperatures instead of Kelvin (defaul %default)\n',
261 action='store_true', default=False)
262 parser.add_option('-o', '--output-file', dest='ofilename',
263 help='write output to FILE (default stdout)',
264 type='string', metavar='FILE')
265 parser.add_option('-p', '--plot-inputs', dest='plot',
266 help='plot the input arrays to check their distribution',
267 action='store_true', default=False)
268 parser.add_option('', '--verbose', dest='verbose', action='store_true',
269 help='print lots of debugging information',
272 options,args = parser.parse_args()
275 config.TEXT_VERBOSE = options.verbose
276 config.PYLAB_INTERACTIVE = False
277 config.PYLAB_VERBOSE = options.plot
279 bumps = get_array(options.bump_string, options.bump_file, "bump")
280 temps = get_array(options.temp_string, options.temp_file, "temp")
281 vibs = get_array(options.vib_string, options.vib_file, "vib")
284 # calib_plot(bumps, temps, vibs)
287 for i in range(len(temps)) :
288 temps[i] = T_analyze.C_to_K(temps[i])
290 km,ks,ps2m,ps2s,Tm,Ts,ooVp2m,ooVp2s = \
291 calib_analyze(bumps, temps, vibs, plotVerbose=options.plot)
294 print >> sys.stderr, \
295 string_errors(km,ks,ps2m,ps2s,Tm,Ts,ooVp2m,ooVp2s)
297 if options.ofilename != None :
298 print >> file(options.ofilename, 'w'), km