calibcant/calibrate.py

   1 #!/usr/bin/python
   2 #
   3 # calibcant - tools for thermally calibrating AFM cantilevers
   4 #
   5 # Copyright (C) 2008-2010 W. Trevor King <wking@drexel.edu>
   6 #
   7 # This file is part of CalibCant.
   8 #
   9 # CalibCant is free software: you can redistribute it and/or
  10 # modify it under the terms of the GNU Lesser General Public
  11 # License as published by the Free Software Foundation, either
  12 # version 3 of the License, or (at your option) any later version.
  13 #
  14 # CalibCant is distributed in the hope that it will be useful,
  15 # but WITHOUT ANY WARRANTY; without even the implied warranty of
  16 # MERCHANTABILITY or FITNESS FOR A PARTICULAR PURPOSE.  See the
  17 # GNU Lesser General Public License for more details.
  18 #
  19 # You should have received a copy of the GNU Lesser General Public
  20 # License along with CalibCant.  If not, see
  21 # <http://www.gnu.org/licenses/>.
  22
  23 """
  24 Aquire and analyze cantilever calibration data.
  25
  26 W. Trevor King Dec. 2007-Jan. 2008
  27
  28 GPL BOILERPLATE
  29
  30
  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 or guess)
  40  k_b      Boltzmann's constant
  41
  42 Which are related by the parameters :
  43  zpGain           Vzp_out / Vzp
  44  zpSensitivity    Zp / Vzp
  45  photoSensitivity Vphoto / Zcant
  46  k_cant           Fcant / Zcant
  47
  48 Cantilever calibration assumes a pre-calibrated z-piezo
  49 (zpSensitivity) and a amplifier (zpGain).  In our lab, the z-piezo is
  50 calibrated by imaging a calibration sample, which has features with
  51 well defined sizes, and the gain is set with a knob on the Nanoscope.
  52
  53 photoSensitivity is measured by bumping the cantilever against the
  54 surface, where Zp = Zcant (see bump_aquire() and the bump_analyze
  55 submodule).
  56
  57 k_cant is measured by watching the cantilever vibrate in free solution
  58 (see the vib_aquire() and the vib_analyze submodule).  The average
  59 energy of the cantilever in the vertical direction is given by the
  60 equipartition theorem.
  61     1/2 k_b T   =   1/2 k_cant <Zcant**2>
  62  so     k_cant  = k_b T / Zcant**2
  63  but    Zcant   = Vphoto / photoSensitivity
  64  so     k_cant  = k_b T * photoSensitivty**2 / <Vphoto**2>
  65
  66 We measured photoSensitivity with the surface bumps.  We can either
  67 measure T using an external function (see temperature.py), or just
  68 estimate it (see T_aquire() and the T_analyze submodule).  Guessing
  69 room temp ~22 deg C is actually fairly reasonable.  Assuming the
  70 actual fluid temperature is within +/- 5 deg, the error in the spring
  71 constant k_cant is within 5/273.15 ~= 2%.  A time series of Vphoto
  72 while we're far from the surface and not changing Vzp_out will give us
  73 the average variance <Vphoto**2>.
  74
  75 We do all these measurements a few times to estimate statistical
  76 errors.
  77
  78 The functions are layed out in the families:
  79  bump_*(), vib_*(), T_*(), and calib_*()
  80 where calib_{save|load|analyze}() deal with derived data, not
  81 real-world data.
  82
  83 For each family, * can be any of :
  84  aquire       get real-world data
  85  save         store real-world data to disk
  86  load         get real-world data from disk
  87  analyze      interperate the real-world data.
  88  plot         show a nice graphic to convince people we're working :p
  89  load_analyze_tweaked
  90               read a file with a list of paths to previously saved
  91               real world data load each file using *_load(), analyze
  92               using *_analyze(), and optionally plot using *_plot().
  93               Intended for re-processing old data.
  94 A family name without any _* extension (e.g. bump()), runs *_aquire(),
  95  *_save(), *_analyze(), *_plot().
  96
  97 We also define the two positioning functions:
  98  move_just_onto_surface() and move_far_from_surface()
  99 which make automating the calibration procedure more straightforward.
 100 """
 101
 102 import numpy
 103 import time
 104
 105 import FFT_tools
 106 import piezo.z_piezo_utils as z_piezo_utils
 107 from splittable_kwargs import splittableKwargsFunction, \
 108     make_splittable_kwargs_function
 109
 110 import .common
 111 import .config
 112 import .bump_analyze
 113 import .T_analyze
 114 import .vib_analyze
 115 import .analyze
 116
 117
 118 # bump family
 119
 120 @splittableKwargsFunction()
 121 def bump_aquire(zpiezo, push_depth=200, npoints=1024, push_speed=1000) :
 122     """
 123     Ramps closer push_depth and returns to the original position.
 124     Inputs:
 125      zpiezo     an opened zpiezo.zpiezo instance
 126      push_depth distance to approach, in nm
 127      npoints    number points during the approach and during the retreat
 128      push_speed piezo speed during approach and retreat, in nm/s
 129     Returns the aquired ramp data dictionary, with data in DAC/ADC bits.
 130     """
 131     # generate the bump output
 132     nm_per_step = float(push_depth) / npoints
 133     freq = push_speed / nm_per_step # freq is sample frequency in Hz
 134     start_pos = zpiezo.curPos()
 135     pos_dist = zpiezo.pos_nm2out(push_depth) - zpiezo.pos_nm2out(0)
 136     close_pos = start_pos + pos_dist
 137     appr = numpy.linspace(start_pos, close_pos, npoints)
 138     retr = numpy.linspace(close_pos, start_pos, npoints)
 139     out = numpy.concatenate((appr, retr))
 140     # run the bump, and measure deflection
 141     if config.TEXT_VERBOSE :
 142         print "Bump %g nm at %g nm/s" % (push_depth, push_speed)
 143     data = zpiezo.ramp(out, freq)
 144     return data
 145
 146 @splittableKwargsFunction(bump_aquire,
 147                           (bump_analyze.bump_save, 'data'),
 148                           (bump_analyze.bump_analyze, 'data'))
 149 def bump(**kwargs):
 150     """
 151     Wrapper around bump_aquire(), bump_save(), bump_analyze()
 152     """
 153     bump_aquire_kwargs,bump_save_kwargs,bump_analyze_kwargs = \
 154         bump._splitargs(bump, kwargs)
 155     data = bump_aquire(**bump_aquire_kwargs)
 156     bump_analyze.bump_save(data, **bump_save_kwargs)
 157     photoSensitivity = bump_analyze.bump_analyze(data, **bump_analyze_kwargs)
 158     return photoSensitivity
 159
 160 # T family.
 161 # Fairly stubby, since a one shot Temp measurement is a common thing.
 162 # We just wrap that to provide a consistent interface.
 163
 164 @splittableKwargsFunction()
 165 def T_aquire(get_T=None) :
 166     """
 167     Measure the current temperature of the sample,
 168     or, if get_T == None, fake it by returning config.DEFAULT_TEMP
 169     """
 170     if get_T == None :
 171         if config.TEXT_VERBOSE :
 172             print "Fake temperature %g" % config.DEFAULT_TEMP
 173         return config.DEFAULT_TEMP
 174     else :
 175         if config.TEXT_VERBOSE :
 176             print "Measure temperature"
 177         return get_T()
 178
 179 @splittableKwargsFunction(T_aquire,
 180                           (T_analyze.T_save, 'T'),
 181                           (T_analyze.T_analyze, 'T'))
 182 def T(**kwargs):
 183     """
 184     Wrapper around T_aquire(), T_save(), T_analyze(), T_plot()
 185     """
 186     T_aquire_kwargs,T_save_kwargs,T_analyze_kwargs = \
 187         T._splitargs(T, kwargs)
 188     T_raw = T_aquire(**T_aquire_kwargs)
 189     T_analyze.T_save(T_raw, **T_save_kwargs)
 190     T_ret = T_analyze.T_analyze(T_raw, **T_analyze_kwargs) # returns array
 191     return T_ret[0]
 192
 193 # vib family
 194
 195 @splittableKwargsFunction()
 196 def vib_aquire(zpiezo, time=1, freq=50e3) :
 197     """
 198     Record data for TIME seconds at FREQ Hz from ZPIEZO at it's current position.
 199     """
 200     # round up to the nearest power of two, for efficient FFT-ing
 201     nsamps = FFT_tools.ceil_pow_of_two(time*freq)
 202     time = nsamps / freq
 203     # take some data, keeping the position voltage at it's current value
 204     out = numpy.ones((nsamps,), dtype=numpy.uint16) * zpiezo.curPos()
 205     if config.TEXT_VERBOSE :
 206         print "get %g seconds of data" % time
 207     data = zpiezo.ramp(out, freq)
 208     data['sample frequency Hz'] = numpy.array([freq])
 209     return data
 210
 211 @splittableKwargsFunction(vib_aquire,
 212                           (vib_analyze.vib_save, 'data'),
 213                           (vib_analyze.vib_analyze, 'deflection_bits', 'freq'))
 214 def vib(**kwargs) :
 215     """
 216     Wrapper around vib_aquire(), vib_save(), vib_analyze()
 217     """
 218     vib_aquire_kwargs,vib_save_kwargs,vib_analyze_kwargs = \
 219         vib._splitargs(vib, kwargs)
 220     data = vib_aquire(**vib_aquire_kwargs)
 221     vib_analyze.vib_save(data, **vib_save_kwargs)
 222     freq = data['sample frequency Hz']
 223     deflection_bits = data['Deflection input']
 224     Vphoto_var = vib_analyze.vib_analyze(deflection_bits=deflection_bits,
 225                                          freq=freq, **vib_analyze_kwargs)
 226     return Vphoto_var
 227
 228 # A few positioning functions, so we can run bump_aquire() and vib_aquire()
 229 # with proper spacing relative to the surface.
 230
 231 @splittableKwargsFunction()
 232 def move_just_onto_surface(stepper, zpiezo, Depth_nm=-50, setpoint=2) :
 233     """
 234     Uses z_piezo_utils.getSurfPos() to pinpoint the position of the
 235     surface.  Adjusts the stepper position as required to get within
 236     stepper_tol nm of the surface.  Then set Vzp to place the
 237     cantilever Depth_nm onto the surface.  Negative Depth_nm values
 238     will place the cantilever that many nm _off_ the surface.
 239
 240     If getSurfPos() fails to find the surface, backs off (for safety)
 241     and steps in (without moving the zpiezo) until Vphoto > setpoint.
 242     """
 243     stepper_tol = 250 # nm, generous estimate of the fullstep stepsize
 244
 245     if config.TEXT_VERBOSE :
 246         print "moving just onto surface"
 247     # Zero the piezo
 248     if config.TEXT_VERBOSE :
 249         print "zero the z piezo output"
 250     zpiezo.jumpToPos(zpiezo.pos_nm2out(0))
 251     # See if we're near the surface already
 252     if config.TEXT_VERBOSE :
 253         print "See if we're starting near the surface"
 254     try :
 255         dist = zpiezo.pos_out2nm( \
 256             z_piezo_utils.getSurfPos(zpiezo, zpiezo.def_V2in(setpoint))
 257                                 )
 258     except (z_piezo_utils.tooClose, z_piezo_utils.poorFit), string :
 259         if config.TEXT_VERBOSE :
 260             print "distance failed with: ", string
 261             print "Back off 200 half steps"
 262         # Back away 200 steps
 263         stepper.step_rel(-400)
 264         stepper.step_rel(200)
 265         sp = zpiezo.def_V2in(setpoint) # sp = setpoint in bits
 266         zpiezo.updateInputs()
 267         cd = zpiezo.curDef()           # cd = current deflection in bits
 268         if config.TEXT_VERBOSE :
 269             print "Single stepping approach"
 270         while cd < sp :
 271             if config.TEXT_VERBOSE :
 272                 print "deflection %g < setpoint %g.  step closer" % (cd, sp)
 273             stepper.step_rel(2) # Full step in
 274             zpiezo.updateInputs()
 275             cd = zpiezo.curDef()
 276         # Back off two steps (protecting against backlash)
 277         if config.TEXT_VERBOSE :
 278             print "Step back 4 half steps to get off the setpoint"
 279         stepper.step_rel(-200)
 280         stepper.step_rel(196)
 281         # get the distance to the surface
 282         zpiezo.updateInputs()
 283         if config.TEXT_VERBOSE :
 284             print "get surf pos, with setpoint %g (%d)" % (setpoint, zpiezo.def_V2in(setpoint))
 285         for i in range(20) : # HACK, keep stepping back until we get a distance
 286             try :
 287                 dist = zpiezo.pos_out2nm( \
 288                     z_piezo_utils.getSurfPos(zpiezo,zpiezo.def_V2in(setpoint)))
 289             except (z_piezo_utils.tooClose, z_piezo_utils.poorFit), string :
 290                 stepper.step_rel(-200)
 291                 stepper.step_rel(198)
 292                 continue
 293             break
 294         if i >= 19 :
 295             print "tried %d times, still too close! bailing" % i
 296             print "probably an invalid setpoint."
 297             raise Exception, "weirdness"
 298     if config.TEXT_VERBOSE :
 299         print 'distance to surface ', dist, ' nm'
 300     # fine tune the stepper position
 301     while dist < -stepper_tol : # step back if we need to
 302         stepper.step_rel(-200)
 303         stepper.step_rel(198)
 304         dist = zpiezo.pos_out2nm( \
 305             z_piezo_utils.getSurfPos(zpiezo, zpiezo.def_V2in(setpoint)))
 306         if config.TEXT_VERBOSE :
 307             print 'distance to surface ', dist, ' nm, step back'
 308     while dist > stepper_tol : # and step forward if we need to
 309         stepper.step_rel(2)
 310         dist = zpiezo.pos_out2nm( \
 311             z_piezo_utils.getSurfPos(zpiezo, zpiezo.def_V2in(setpoint)))
 312         if config.TEXT_VERBOSE :
 313             print 'distance to surface ', dist, ' nm, step closer'
 314     # now adjust the zpiezo to place us just onto the surface
 315     target = dist + Depth_nm
 316     zpiezo.jumpToPos(zpiezo.pos_nm2out(target))
 317     # and we're there :)
 318     if config.TEXT_VERBOSE :
 319         print "We're %g nm into the surface" % Depth_nm
 320
 321 @splittableKwargsFunction()
 322 def move_far_from_surface(stepper, um_back=50) :
 323     """
 324     Step back a specified number of microns.
 325     (uses very rough estimate of step distance at the moment)
 326     """
 327     step_nm = 100
 328     steps = int(um_back*1000/step_nm)
 329     print "step back %d steps" % steps
 330     stepper.step_rel(-steps)
 331
 332
 333 # and finally, the calib family
 334
 335 @splittableKwargsFunction((move_just_onto_surface, 'stepper', 'zpiezo'),
 336                           (bump, 'zpiezo', 'log_dir', 'Vphoto_in2V'),
 337                           (move_far_from_surface, 'stepper'),
 338                           (T, 'log_dir'),
 339                           (vib, 'zpiezo', 'log_dir', 'Vphoto_in2V'),
 340                           (analyze.calib_save, 'bumps','Ts','vibs','log_dir'))
 341 def calib_aquire(stepper, zpiezo, num_bumps=10, num_Ts=10, num_vibs=20,
 342                  log_dir=config.LOG_DIR, Vphoto_in2V=config.Vphoto_in2V,
 343                  **kwargs):
 344     """
 345     Aquire data for calibrating a cantilever in one function.
 346     return (bump, T, vib), each of which is an array.
 347     Inputs :
 348      stepper       a stepper.stepper_obj for coarse Z positioning
 349      zpiezo        a z_piezo.z_piezo for fine positioning and deflection readin
 350      num_bumps     number of 'bumps' (see Outputs)
 351      num_temps     number of 'Ts' (see Outputs)
 352      num_vibs      number of 'vib's (see Outputs)
 353      log_dir       directory to log data to.  Default 'None' disables logging.
 354      Vphoto_in2V   function to convert photodiode input bits to Volts
 355
 356      + other kwargs.  Run calib_aquire._kwargs(calib_aquire) to see
 357      all options.  Run calib_aquire._childSplittables to see a list
 358      of kwarg functions that this function calls.
 359
 360     Outputs (all are arrays of recorded data) :
 361      bumps measured (V_photodiode / nm_tip) proportionality constant
 362      Ts    measured temperature (K)
 363      vibs  measured V_photodiode variance in free solution
 364     """
 365     move_just_onto_surface_kwargs,bump_kwargs,move_far_from_surface_kwargs, \
 366         T_kwargs,vib_kwargs,calib_save_kwargs = \
 367         calib_aquire._splitargs(calib_aquire, kwargs)
 368     # get bumps
 369     bumps = numpy.zeros((num_bumps,), dtype=numpy.float)
 370     for i in range(num_bumps) :
 371         move_just_onto_surface(stepper, zpiezo, **move_just_onto_surface_kwargs)
 372         bumps[i] = bump(zpiezo=zpiezo, log_dir=log_dir,
 373                         Vphoto_in2V=Vphoto_in2V, **bump_kwargs)
 374     if config.TEXT_VERBOSE :
 375         print bumps
 376
 377     move_far_from_surface(stepper, **move_far_from_surface_kwargs)
 378
 379     # get Ts
 380     Ts = numpy.zeros((num_Ts,), dtype=numpy.float)
 381     for i in range(num_Ts) :
 382         Ts[i] = T(**T_kwargs)
 383         time.sleep(1) # wait a bit to get an independent temperature measure
 384     print Ts
 385
 386     # get vibs
 387     vibs = numpy.zeros((num_vibs,), dtype=numpy.float)
 388     for i in range(num_vibs) :
 389         vibs[i] = vib(zpiezo=zpiezo, log_dir=log_dir, Vphoto_in2V=Vphoto_in2V,
 390                       **vib_kwargs)
 391     print vibs
 392
 393     analyze.calib_save(bumps, Ts, vibs, log_dir, **calib_save_kwargs)
 394
 395     return (bumps, Ts, vibs)
 396
 397
 398 @splittableKwargsFunction( \
 399     (calib_aquire, 'log_dir'),
 400     (analyze.calib_analyze, 'bumps','Ts','vibs'))
 401 def calib(log_dir=config.LOG_DIR, **kwargs) :
 402     """
 403     Calibrate a cantilever in one function.
 404     The I-don't-care-about-the-details black box version :p.
 405     return (k, k_s)
 406     Inputs:
 407      (see calib_aquire())
 408     Outputs :
 409      k    cantilever spring constant (in N/m, or equivalently nN/nm)
 410      k_s  standard deviation in our estimate of k
 411     Notes :
 412      See get_calibration_data() for the data aquisition code
 413      See analyze_calibration_data() for the analysis code
 414     """
 415     calib_aquire_kwargs,calib_analyze_kwargs = \
 416         calib._splitargs(calib, kwargs)
 417     a, T, vib = calib_aquire(**calib_aquire_kwargs)
 418     k,k_s,ps2_m, ps2_s,T_m,T_s,one_o_Vp2_m,one_o_Vp2_s = \
 419         analyze.calib_analyze(a, T, vib, **calib_analyze_kwargs)
 420     analyze.calib_save_analysis(k, k_s, ps2_m, ps2_s, T_m, T_s,
 421                                 one_o_Vp2_m, one_o_Vp2_s, log_dir)
 422     return (k, k_s)
 423
 424
 425