posts/Comparing_velocity_clamp_experiments/compare-unfold.py

   1 #!/usr/bin/env python2
   2 #
   3 # Copyright (C) 2012  W. Trevor King <wking@tremily.us>
   4 #
   5 # This program is free software: you can redistribute it and/or modify
   6 # it under the terms of the GNU General Public License as published by
   7 # the Free Software Foundation, either version 3 of the License, or
   8 # (at your option) any later version.
   9 #
  10 # This program is distributed in the hope that it will be useful,
  11 # but WITHOUT ANY WARRANTY; without even the implied warranty of
  12 # MERCHANTABILITY or FITNESS FOR A PARTICULAR PURPOSE.  See the
  13 # GNU General Public License for more details.
  14 #
  15 # You should have received a copy of the GNU General Public License
  16 # along with this program.  If not, see <http://www.gnu.org/licenses/>.
  17
  18 """Compare unfolding pulls between Prof. Yang's software and my own.
  19 """
  20
  21 import collections as _collections
  22 import os as _os
  23
  24 import crunch as _crunch
  25 import gyang_curve as _gyang_curve
  26 from matplotlib import pylab as _pylab
  27 import numpy as _numpy
  28 from pypiezo import surface as _pypiezo_surface
  29
  30
  31 def load_yang_curves(basedir, single=False):
  32     if _os.path.isfile(basedir):
  33         for curve in _gyang_curve.load_all_traces(basedir):
  34             yield(curve)
  35         return
  36     for series in sorted(_os.listdir(basedir)):
  37         if series.endswith('.ibw'):
  38             path = _os.path.join(basedir, series)
  39             for curve in _gyang_curve.load_all_traces(path):
  40                 yield curve
  41                 if single:
  42                     return
  43
  44 def load_king_curves(basedir, single=False):
  45     if _os.path.isfile(basedir):
  46         yield _crunch.load_pull(filename=basedir)
  47         return
  48     for path in sorted(_os.listdir(basedir)):
  49         if path.endswith('.h5'):
  50             curve = _os.path.join(basedir, path)
  51             yield _crunch.load_pull(filename=curve)
  52             if single:
  53                 return
  54
  55 def plot_comparison(yang_dir, king_dir, single=False, contact_slope=False):
  56     yang = load_yang_curves(yang_dir, single=single)
  57     king = load_king_curves(king_dir, single=single)
  58     figure = _pylab.figure()
  59     axes = figure.add_subplot(1, 1, 1)
  60     axes.set_title('velocity clamp unfolding')
  61     if contact_slope:
  62         axes.set_xlabel('pull index')
  63         axes.set_ylabel('contact slope (bits/bit)')
  64         hist_data = {'b': [], 'r': []}
  65     else:
  66         axes.set_xlabel('z piezo (bits)')
  67         axes.set_ylabel('deflection (bits)')
  68     axes.hold(True)
  69     xyc = []
  70     for curve in yang:
  71         xyc.append(
  72             (curve.piezo_expanded + 2**15, curve.deflection + 2**15, 'b'))
  73     for curve in king:
  74         xyc.append((curve['alpha_v'], curve['alpha_f'], 'r'))
  75
  76     for x, y, c in xyc:
  77         ddict = {
  78             'approach': {  # for _analyze_surface_position_data
  79                 'z': x[::-1],
  80                 'deflection': y[::-1],
  81                 }
  82             }
  83         (non_contact_offset,non_contact_slope,kink,alpha_dp
  84          ) = _pypiezo_surface.analyze_surface_position_data(
  85             ddict, return_all_parameters=True)
  86         y0 = non_contact_offset + non_contact_slope*(kink-x.min())
  87         x0 = x.max() - (y.max()-y0)/alpha_dp
  88         if contact_slope:
  89             hist_data[c].append(alpha_dp)
  90         else:
  91             axes.plot(
  92                 x - x0, y - y0,
  93                 linestyle='none', marker=',', color=c, markeredgecolor=c)
  94     if contact_slope:
  95         for c,data in hist_data.items():
  96             data = [d for d in data if d > 5.5]  # HACK!!
  97             axes.plot(
  98                 data,
  99                 linestyle='none', marker='o', color=c, markeredgecolor=c)
 100             #axes.hist(data, color=c, alpha=0.5)
 101             data = _numpy.array(data)
 102             print('\n'.join(
 103                     '{}: {}'.format(key,value)
 104                     for key,value in {
 105                         'color': c,
 106                         'mean': data.mean(),
 107                         'std': data.std(),
 108                         }.items()))
 109     return figure
 110
 111
 112 if __name__ == '__main__':
 113     import argparse
 114
 115     parser = argparse.ArgumentParser(description=__doc__)
 116     parser.add_argument(
 117         '--single', action='store_const', const=True,
 118         help='compare only a single curve of each type')
 119     parser.add_argument(
 120         '--contact-slope', action='store_const', const=True,
 121         help='compare contact slopes')
 122     parser.add_argument(
 123         '-o', '--output', default=None,
 124         help='save figure instead of showing it')
 125     parser.add_argument(
 126         'yang', metavar='DIR',
 127         help="directory for data taken with Prof. Yang's LabVIEW program")
 128     parser.add_argument(
 129         'king', metavar='DIR',
 130         help='directory for data taken with my pyafm-based unfolder')
 131
 132     args = parser.parse_args()
 133
 134     figure = plot_comparison(
 135         yang_dir=_os.path.expanduser(args.yang),
 136         king_dir=_os.path.expanduser(args.king),
 137         single=args.single,
 138         contact_slope=args.contact_slope)
 139     if args.output:
 140         figure.savefig(args.output)
 141     else:
 142         _pylab.show()