pypiezo/z_piezo_utils.py

   1 # Copyright (C) 2008-2011 W. Trevor King <wking@drexel.edu>
   2 #
   3 # This file is part of pypiezo.
   4 #
   5 # pypiezo is free software; you can redistribute it and/or modify it
   6 # under the terms of the GNU General Public License as published by the
   7 # Free Software Foundation, either version 3 of the License, or (at your
   8 # option) any later version.
   9 #
  10 # pypiezo is distributed in the hope that it will be useful, but
  11 # WITHOUT ANY WARRANTY; without even the implied warranty of
  12 # MERCHANTABILITY or FITNESS FOR A PARTICULAR PURPOSE.  See the GNU
  13 # General Public License for more details.
  14 #
  15 # You should have received a copy of the GNU General Public License
  16 # along with pypiezo.  If not, see <http://www.gnu.org/licenses/>.
  17
  18 TEXT_VERBOSE = False
  19 PYLAB_INTERACTIVE = True
  20 PYLAB_VERBOSE = False
  21 BASE_FIG_NUM = 60
  22 LOG_DATA = True
  23 LOG_DIR = '${DEFAULT}/z_piezo_utils'
  24
  25 NULL_STEPS_BEFORE_SINGLE_PT_SWEEP = False
  26 SLEEP_DURING_SURF_POS_AQUISITION = False # doesn't help
  27
  28 from numpy import array, arange, ones, zeros, linspace, uint16, float, sin, cos, pi
  29 from scipy.stats import linregress
  30 import curses_check_for_keypress
  31 import scipy.optimize
  32 import data_logger
  33 import time
  34
  35 _flush_plot = None
  36 _final_flush_plot = None
  37 _pylab = None
  38 def _dummy_fn(): pass
  39
  40 def _import_pylab():
  41     """Import pylab plotting functions for when we need to plot.  This
  42     function can be called multiple times, and ensures that the pylab
  43     setup is only imported once.  It defines the functions
  44       _flush_plot() to be called after incremental changes,
  45       _final_flush_plot(), to be called at the end of any graphical processing, and
  46       _pylab, a reference to the pylab module."""
  47     global _pylab
  48     global _flush_plot
  49     global _final_flush_plot
  50     if _pylab == None :
  51         import pylab as _pylab
  52     if _flush_plot == None :
  53         if PYLAB_INTERACTIVE :
  54             _flush_plot = _pylab.draw
  55         else :
  56             _flush_plot = _pylab.show
  57     if _final_flush_plot == None :
  58         if PYLAB_INTERACTIVE :
  59             _final_flush_plot = _pylab.draw
  60         else :
  61             _final_flush_plot = _dummy_fn
  62
  63 class error (Exception) :
  64     pass
  65 class poorFit (error) :
  66     pass
  67 class tooClose (error) :
  68     pass
  69
  70 def moveToPosOrDef(zpiezo, pos, defl, step, return_data=False) :
  71     if return_data or LOG_DATA or PYLAB_VERBOSE :
  72         aquire_data = True
  73     else :
  74         aquire_data = False
  75     zpiezo._check_range(pos)
  76     if step == 0 :
  77         raise error, "Must have non-zero step size"
  78     elif step < 0 and pos > zpiezo.curPos() :
  79         step = -step
  80     elif step > 0 and pos < zpiezo.curPos() :
  81         step = -step
  82
  83     zpiezo._jump.reserve()
  84     zpiezo.updateInputs()
  85     if TEXT_VERBOSE :
  86         print "current z_piezo position %6d, current deflection %6d <? %6d" % (zpiezo.curPos(), zpiezo.curDef(), defl)
  87     if aquire_data :
  88         def_array=[zpiezo.curDef()]
  89         pos_array=[zpiezo.curPos()]
  90     if NULL_STEPS_BEFORE_SINGLE_PT_SWEEP == True :
  91         # take a few null steps to let the input deflection stabalize
  92         for i in range(100) :
  93             zpiezo.jumpToPos(zpiezo.curPos())
  94             if aquire_data :
  95                 def_array.append(zpiezo.curDef())
  96                 pos_array.append(zpiezo.curPos())
  97     # step in until we hit our target position or exceed our target deflection
  98     while zpiezo.curPos() != pos and zpiezo.curDef() < defl :
  99         distTo = pos - zpiezo.curPos()
 100         if abs(distTo) < abs(step) :
 101             jumpTo = pos
 102         else :
 103             jumpTo = zpiezo.curPos() + step
 104         zpiezo.jumpToPos(jumpTo)
 105         if TEXT_VERBOSE :
 106             print "current z_piezo position %6d, current deflection %6d >=< %6d" % (zpiezo.curPos(), zpiezo.curDef(), defl)
 107         if aquire_data :
 108             def_array.append(zpiezo.curDef())
 109             pos_array.append(zpiezo.curPos())
 110     zpiezo._jump.release()
 111
 112     if zpiezo.setExternalZPiezo != None :
 113         zpiezo.setExternalZPiezo(zpiezo.pos_out2nm(zpiezo.curPos()))
 114     if PYLAB_VERBOSE :
 115         _import_pylab()
 116         _pylab.figure(BASE_FIG_NUM)
 117         timestamp = time.strftime('%H%M%S')
 118         _pylab.plot(pos_array, def_array, '.', label=timestamp)
 119         _pylab.title('step approach')
 120         _pylab.legend(loc='best')
 121         _flush_plot()
 122     if LOG_DATA :
 123         log = data_logger.data_log(LOG_DIR, noclobber_logsubdir=False,
 124                                    log_name="moveToPosOrDef")
 125         data = {"Z piezo output":array(pos_array),
 126                 "Deflection input":array(def_array)}
 127         log.write_dict_of_arrays(data)
 128     if return_data :
 129         data = {"Z piezo output":array(pos_array),
 130                 "Deflection input":array(def_array)}
 131         return (zpiezo.curVals(), data)
 132     else :
 133         return zpiezo.curVals()
 134
 135
 136 def wiggleForInterferenceMin(zpiezo, wig_amp=20000, wig_freq=100,
 137                              plotVerbose=True) :
 138     "Output sin to measure interference"
 139     if PYLAB_VERBOSE or plotVerbose:
 140         _import_pylab()
 141     npoints = 1000
 142     scanfreq = wig_freq*npoints
 143     out = zeros((npoints,), dtype=uint16)
 144     for i in range(npoints) :
 145         out[i] = int(sin(4*pi*i/npoints)*wig_amp+32000)
 146         # 4 for 2 periods
 147     if PYLAB_VERBOSE or plotVerbose:
 148         _pylab.figure(BASE_FIG_NUM+1)
 149         a = _pylab.axes()
 150         _pylab.title('wiggle for interference')
 151         _pylab.hold(False)
 152         p = _pylab.plot(out, out, 'b.-')
 153     if LOG_DATA :
 154         log = data_logger.data_log(LOG_DIR, noclobber_logsubdir=False,
 155                                    log_name="wiggle")
 156     c = curses_check_for_keypress.check_for_keypress(test_mode=False)
 157     while c.input() == None :
 158         data = zpiezo.ramp(out, scanfreq)
 159         dmin = data["Deflection input"].min()
 160         dmax = data["Deflection input"].max()
 161         if PYLAB_VERBOSE or plotVerbose:
 162             p[0].set_ydata(data["Deflection input"])
 163             a.set_ylim([dmin,dmax])
 164             _flush_plot()
 165         if LOG_DATA :
 166             log.write_dict_of_arrays(data)
 167     if PYLAB_VERBOSE or plotVerbose:
 168         _pylab.hold(True)
 169
 170 def getSurfPosData(zpiezo, maxDefl, textVerbose=False, plotVerbose=False) :
 171     "Measure the distance to the surface"
 172     if plotVerbose and not PYLAB_VERBOSE :
 173         _import_pylab()
 174     origPos = zpiezo.curPos()
 175     # fully retract the piezo
 176     if TEXT_VERBOSE or textVerbose :
 177         print "\tRetract the piezo"
 178     out = linspace(zpiezo.curPos(), zpiezo.zpMin, 2000)
 179     ret1 = zpiezo.ramp(out, 10e3)
 180     del(out) # start with a clean output vector, seems to help reduce bounce?
 181     zpiezo.updateInputs() # opening output seems to cause bounce?
 182     # locate high deflection position
 183     if TEXT_VERBOSE or textVerbose :
 184         print "\tApproach until there is dangerous deflection (> %d)" % maxDefl
 185     if SLEEP_DURING_SURF_POS_AQUISITION == True :
 186         time.sleep(.2) # sleeping briefly seems to reduce bounce?
 187     cp, mtpod = moveToPosOrDef(zpiezo, zpiezo.zpMax, maxDefl,
 188                                (zpiezo.zpMax-zpiezo.zpMin)/2000,
 189                                return_data=True)
 190     highContactPos = zpiezo.curPos()
 191     if highContactPos <= origPos :   # started out too close, we need some room to move in...
 192         if TEXT_VERBOSE or textVerbose :
 193             print "\tToo close, return to the original position: ", origPos
 194         out = linspace(zpiezo.curPos(), origPos, 2000)
 195         ret3=zpiezo.ramp(out, 10e3)
 196         del(out)
 197         zpiezo.updateInputs()
 198         if PYLAB_VERBOSE or plotVerbose :
 199             _pylab.figure(BASE_FIG_NUM+2)
 200             def plot_dict(d, label) :
 201                 print d.keys()
 202                 _pylab.plot(d["Z piezo output"], d["Deflection input"], '.',label=label)
 203             plot_dict(ret1, 'First retract')
 204             plot_dict(mtpod, 'Single step in')
 205             plot_dict(ret3, 'Return to start')
 206             _pylab.legend(loc='best')
 207             _pylab.title('Too close')
 208         raise tooClose, "Way too close!\n(highContactPos %d <= origPos %d)" % (highContactPos, origPos)
 209     zpiezo.updateInputs()
 210     # fully retract the piezo again
 211     if TEXT_VERBOSE or textVerbose :
 212         print "\tRetract the piezo again"
 213     if SLEEP_DURING_SURF_POS_AQUISITION == True :
 214         time.sleep(.2)
 215     out = linspace(highContactPos, zpiezo.zpMin, 2000)
 216     ret2 = zpiezo.ramp(out, 10e3)
 217     del(out)
 218     zpiezo.updateInputs()
 219     if zpiezo.curPos() != zpiezo.zpMin :
 220         raise Exception, "bad code"
 221     # scan to the high contact position
 222     if TEXT_VERBOSE or textVerbose :
 223         print "\tRamp in to the deflected position: ", highContactPos
 224     if SLEEP_DURING_SURF_POS_AQUISITION == True :
 225         time.sleep(.2)
 226     out = linspace(zpiezo.curPos(), highContactPos, 4000)
 227     data = zpiezo.ramp(out, 10e3)
 228     del(out)
 229     zpiezo.updateInputs()
 230     if LOG_DATA :
 231         log = data_logger.data_log(LOG_DIR, noclobber_logsubdir=False,
 232                                    log_name="getSurfPos")
 233         log.write_dict_of_arrays(data)
 234     if SLEEP_DURING_SURF_POS_AQUISITION == True :
 235         time.sleep(.2)
 236     # return to the original position
 237     if TEXT_VERBOSE or textVerbose :
 238         print "\tReturn to the original position: ", origPos
 239     out = linspace(highContactPos, origPos, 2000)
 240     ret3=zpiezo.ramp(out, 10e3)
 241     del(out)
 242     zpiezo.updateInputs()
 243     return {"ret1":ret1, "mtpod":mtpod, "ret2":ret2,
 244             "approach":data, "ret3":ret3}
 245
 246 def bilinear(x, params):
 247     """bilinear fit for surface bumps.  Model has two linear regimes
 248     which meet at x=kink_position and have independent slopes.
 249     """
 250     assert type(x) == type(array([0,1])), "Invalid x array: %s" % type(x)
 251     left_offset,left_slope,kink_position,right_slope = params
 252     left_mask = x < kink_position
 253     right_mask = x >= kink_position # = not left_mask
 254     left_y = left_offset + left_slope*x
 255     right_y = left_offset + left_slope*kink_position \
 256         + right_slope*(x-kink_position)
 257     return left_mask * left_y + right_mask * right_y
 258
 259 def analyzeSurfPosData(ddict, textVerbose=False, plotVerbose=False, retAllParams=False) :
 260     # ususes ddict["approach"] for analysis
 261     # the others are just along to be plotted
 262
 263     data = ddict["approach"]
 264     # analyze data, using bilinear model
 265     # y = p0 + p1 x                for x <= p2
 266     #   = p0 + p1 p2 + p3 (x-p2)   for x >= p2
 267     if TEXT_VERBOSE or textVerbose :
 268         print "\tAnalyze the data"
 269     dump_before_index = 0 # 25 # HACK!!
 270     # Generate a reasonable guess...
 271     startPos = int(data["Z piezo output"].min())
 272     finalPos = int(data["Z piezo output"].max())
 273     startDef = int(data["Deflection input"].min())
 274     finalDef = int(data["Deflection input"].max())
 275     # startDef and startPos are probably for 2 different points
 276
 277     left_offset   = startDef
 278     left_slope    = 0
 279     kink_position = (finalPos+startPos)/2.0
 280     right_slope   = 2.0*(finalDef-startDef)/(finalPos-startPos)
 281     pstart = [left_offset, left_slope, kink_position, right_slope]
 282
 283     offset_scale = (finalPos - startPos)/100
 284     left_slope_scale = right_slope/10
 285     kink_scale = (finalPos-startPos)/100
 286     right_slope_scale = right_slope
 287     scale = [offset_scale, left_slope_scale, kink_scale, right_slope_scale]
 288     if TEXT_VERBOSE or textVerbose :
 289         print "start Def %d, final def %d, start pos %d, final pos %d" % (startDef, finalDef, startPos, finalPos)
 290         print "Guessed params ", pstart
 291     def residual(p, y, x):
 292         Y = bilinear(x, p)
 293         return Y-y
 294     leastsq = scipy.optimize.leastsq
 295     params,cov,info,mesg,ier = leastsq(
 296         residual, pstart,
 297         args=(data["Deflection input"][dump_before_index:],
 298               data["Z piezo output"][dump_before_index:]),
 299         full_output=True, maxfev=10000)
 300     if TEXT_VERBOSE or textVerbose :
 301         print "Best fit parameters: ", params
 302
 303     if PYLAB_VERBOSE or plotVerbose :
 304         _import_pylab()
 305         _pylab.figure(BASE_FIG_NUM+3)
 306         def plot_dict(d, label) :
 307             _pylab.plot(d["Z piezo output"], d["Deflection input"], '.',label=label)
 308         try :
 309             plot_dict(ddict['ret1'], 'First retract')
 310             plot_dict(ddict['mtpod'], 'Single step in')
 311             plot_dict(ddict['ret2'], 'Second retract')
 312         except KeyError :
 313             pass # if we don't have the surrounding data, don't worry
 314         plot_dict(data, 'Main approach')
 315         try :
 316             plot_dict(ddict['ret3'], 'Return to start')
 317         except KeyError :
 318             pass # if we don't have the surrounding data, don't worry
 319         def fit_fn(x, params) :
 320             if x <= params[2] :
 321                 return params[0] + params[1]*x
 322             else :
 323                 return params[0] + params[1]*params[2] + params[3]*(x-params[2])
 324         _pylab.plot([startPos, params[2], finalPos],
 325                     [fit_fn(startPos, params), fit_fn(params[2], params), fit_fn(finalPos, params)],
 326                     '-',label='fit')
 327         _pylab.title('surf_pos %5g %5g %5g %5g' % (params[0], params[1], params[2], params[3]))
 328         _pylab.legend(loc='best')
 329         _flush_plot()
 330
 331     # check that the fit is reasonable
 332     # params[1] is slope in non-contact region
 333     # params[3] is slope in contact region
 334     if abs(params[1]*10) > abs(params[3]) :
 335         raise poorFit, "Slope in non-contact region, or no slope in contact"
 336     if params[2] < startPos+0.02*(finalPos-startPos) : # params[2] is kink position
 337         raise poorFit, "No kink (kink %g less than %g, need more space to left)" % \
 338             (params[2], startPos+0.02*(finalPos-startPos))
 339     if params[2] > finalPos-0.02*(finalPos-startPos) :
 340         raise poorFit, "No kink (kink %g more than %g, need more space to right)" % \
 341             (params[2], finalPos-0.02*(finalPos-startPos))
 342     if params[3] < 0.5 * right_slope : # params[3] is slope in contact region
 343         raise poorFit, "Too far"
 344     if TEXT_VERBOSE or textVerbose :
 345         print "\tSurface position: ", params[2]
 346     if retAllParams :
 347         return params
 348     return params[2]
 349
 350 def getSurfPos(zpiezo, maxDefl, textVerbose=False, plotVerbose=False) :
 351     ddict = getSurfPosData(zpiezo, maxDefl, textVerbose, plotVerbose)
 352     return analyzeSurfPosData(ddict, textVerbose, plotVerbose)
 353
 354 def test_smoothness(zp, plotVerbose=True) :
 355     posA = 20000
 356     posB = 50000
 357     setpoint = zp.def_V2in(3)
 358     steps = 200
 359     outfreq = 1e5
 360     outarray = linspace(posB, posA, 1000)
 361     indata=[]
 362     outdata=[]
 363     curVals = zp.jumpToPos(posA)
 364     zp.pCurVals(curVals)
 365     time.sleep(1) # let jitters die down
 366     for i in range(10) :
 367         print "ramp %d to %d" % (zp.curPos(), posB)
 368         curVals, data = moveToPosOrDef(zp, posB, setpoint, step=steps,
 369                                        return_data = True)
 370         indata.append(data)
 371         out = zp.ramp(outarray, outfreq)
 372         outdata.append(out)
 373     if plotVerbose :
 374         from pylab import figure, plot, title, legend, hold, subplot
 375     if PYLAB_VERBOSE or plotVerbose :
 376         _import_pylab()
 377         _pylab.figure(BASE_FIG_NUM+4)
 378         for i in range(10) :
 379             _pylab.plot(indata[i]['Z piezo output'],
 380                         indata[i]['Deflection input'], '+--', label='in')
 381             _pylab.plot(outdata[i]['Z piezo output'],
 382                         outdata[i]['Deflection input'], '.-', label='out')
 383         _pylab.title('test smoothness (step in, ramp out)')
 384         #_pylab.legend(loc='best')
 385
 386 def test() :
 387     import z_piezo
 388     zp = z_piezo.z_piezo()
 389     curVals = zp.moveToPosOrDef(zp.pos_nm2out(600), defl=zp.curDef()+6000, step=(zp.pos_nm2out(10)-zp.pos_nm2out(0)))
 390     if TEXT_VERBOSE :
 391         zp.pCurVals(curVals)
 392     pos = zp.getSurfPos(maxDefl=zp.curDef()+6000)
 393     if TEXT_VERBOSE :
 394         print "Surface at %g nm", pos
 395     print "success"
 396     if PYLAB_VERBOSE and _final_flush_plot != None:
 397         _final_flush_plot()