hooke/driver/picoforce.py

   1 # Copyright (C) 2006-2010 Alberto Gomez-Casado
   2 #                         Massimo Sandal <devicerandom@gmail.com>
   3 #                         W. Trevor King <wking@drexel.edu>
   4 #
   5 # This file is part of Hooke.
   6 #
   7 # Hooke is free software: you can redistribute it and/or
   8 # modify it under the terms of the GNU Lesser General Public
   9 # License as published by the Free Software Foundation, either
  10 # version 3 of the License, or (at your option) any later version.
  11 #
  12 # Hooke is distributed in the hope that it will be useful,
  13 # but WITHOUT ANY WARRANTY; without even the implied warranty of
  14 # MERCHANTABILITY or FITNESS FOR A PARTICULAR PURPOSE.  See the
  15 # GNU Lesser General Public License for more details.
  16 #
  17 # You should have received a copy of the GNU Lesser General Public
  18 # License along with Hooke.  If not, see
  19 # <http://www.gnu.org/licenses/>.
  20
  21 """Driver for Veeco PicoForce force spectroscopy files.
  22 """
  23
  24 import pprint
  25 import re
  26 import time
  27
  28 import numpy
  29
  30 from .. import curve as curve # this module defines data containers.
  31 from .. import experiment as experiment # this module defines expt. types
  32 from ..config import Setting # configurable setting class
  33 from . import Driver as Driver # this is the Driver base class
  34
  35
  36 __version__='0.0.0.20100516'
  37
  38 class PicoForceDriver (Driver):
  39     """Handle Veeco Picoforce force spectroscopy files.
  40     """
  41     def __init__(self):
  42         super(PicoForceDriver, self).__init__(name='picoforce')
  43
  44     def is_me(self, path):
  45         f = file(path, 'r')
  46         header = f.read(30)
  47         f.close()
  48
  49         return header[2:17] == 'Force file list'
  50
  51     def read(self, path):
  52         info = self._read_header_path(path)
  53         self._check_version(info)
  54         data = self._read_data_path(path, info)
  55         info['filetype'] = self.name
  56         info['experiment'] = experiment.VelocityClamp
  57         return (data, info)
  58
  59     def _read_header_path(self, path):
  60         """Read curve information from the PicoForce file at `path`.
  61
  62         See :meth:`._read_header_file`.
  63         """
  64         return self._read_header_file(file(path, 'rb'))
  65
  66     def _read_header_file(self, file):
  67         r"""Read curve information from a PicoForce file.
  68
  69         Return a dict of dicts representing the information.  If a
  70         field is repeated multiple times, it's value is replaced by a
  71         list of the values for each occurence.
  72
  73         Examples
  74         --------
  75
  76         >>> import pprint
  77         >>> import StringIO
  78         >>> f = StringIO.StringIO('\r\n'.join([
  79         ...             '\*Force file list',
  80         ...             '\Version: 0x06120002',
  81         ...             '\Date: 04:42:34 PM Tue Sep 11 2007',
  82         ...             '\Start context: FOL2',
  83         ...             '\Data length: 40960',
  84         ...             '\Text: ',
  85         ...             '\*Equipment list',
  86         ...             '\Description: Extended PicoForce',
  87         ...             '\Controller: IIIA',
  88         ...             '\*Ciao force image list',
  89         ...             '\Data offset: 40960',
  90         ...             '\Data length: 8192',
  91         ...             '\*Ciao force image list',
  92         ...             '\Data offset: 49152',
  93         ...             '\Data length: 8192',
  94         ...             '\*Ciao force image list',
  95         ...             '\Data offset: 57344',
  96         ...             '\Data length: 8192',
  97         ...             ]))
  98         >>> p = PicoForceDriver()
  99         >>> d = p._read_header_file(f)
 100         >>> pprint.pprint(d, width=60)
 101         {'Ciao force image list': [{'Data length': '8192',
 102                                     'Data offset': '40960'},
 103                                    {'Data length': '8192',
 104                                     'Data offset': '49152'},
 105                                    {'Data length': '8192',
 106                                     'Data offset': '57344'}],
 107          'Equipment list': {'Controller': 'IIIA',
 108                             'Description': 'Extended PicoForce'},
 109          'Force file list': {'Data length': '40960',
 110                              'Date': '04:42:34 PM Tue Sep 11 2007',
 111                              'Start context': 'FOL2',
 112                              'Text:': None,
 113                              'Version': '0x06120002'}}
 114         """
 115         info = {}
 116         header_field = None
 117         for line in file:
 118             line = line.strip()
 119             if line.startswith('\*File list end'):
 120                 break
 121             if line.startswith(r'\*'):
 122                 header_field = line[len(r'\*'):]
 123                 if header_field in info:
 124                     if isinstance(info[header_field], list):
 125                         info[header_field].append({}) # >=3rd appearance
 126                     else: # Second appearance
 127                         info[header_field] = [info[header_field], {}]
 128                 else: # First appearance
 129                     info[header_field] = {}
 130             else:
 131                 assert line.startswith('\\'), line
 132                 fields = line[len('\\'):].split(': ', 1)
 133                 key = fields[0]
 134                 if len(fields) == 1: # fields = [key]
 135                     value = None
 136                 else: # fields = [key, value]
 137                     value = fields[1]
 138                 if isinstance(info[header_field], list): # >=2nd header_field
 139                     target_dict = info[header_field][-1]
 140                 else: # first appearance of header_field
 141                     target_dict = info[header_field]
 142                 if key in target_dict and target_dict[key] != value:
 143                     raise NotImplementedError(
 144                         'Overwriting %s: %s -> %s'
 145                         % (key, target_dict[key], value))
 146                 target_dict[key] = value
 147         return (info)
 148
 149     def _check_version(self, info):
 150         """Ensure the input file is a version we understand.
 151
 152         Otherwise, raise `ValueError`.
 153         """
 154         version = info['Force file list'].get('Version', None)
 155         if version not in ['0x06120002', '0x06130001', '0x07200000']:
 156             raise NotImplementedError(
 157                 '%s file version %s not supported (yet!)\n%s'
 158                 % (self.name, version,
 159                    pprint.pformat(info['Force file list'])))
 160
 161     def _read_data_path(self, path, info):
 162         """Read curve data from the PicoForce file at `path`.
 163
 164         See :meth:`._read_data_file`.
 165         """
 166         f = file(path, 'rb')
 167         data = self._read_data_file(f, info)
 168         f.close()
 169         return data
 170
 171     def _read_data_file(self, file, info):
 172         file.seek(0)
 173         traces = self._extract_traces(buffer(file.read()), info)
 174         for k,t in traces.items():
 175             t.tofile('t-%s' % k, sep='\n')
 176         self._validate_traces(
 177             traces['Z sensor'], traces['Deflection'])
 178         L = len(traces['Deflection'])
 179         approach = self._extract_block(
 180             info, traces['Z sensor'], traces['Deflection'], 0, L/2, 'approach')
 181         retract = self._extract_block(
 182             info, traces['Z sensor'], traces['Deflection'], L/2, L, 'retract')
 183         data = [approach, retract]
 184         return data
 185
 186     def _extract_traces(self, buffer, info):
 187         """Extract each of the three vector blocks in a PicoForce file.
 188
 189         The blocks are (in variable order):
 190
 191         * Z piezo sensor input
 192         * Deflection input
 193         * Deflection again?
 194
 195         And their headers are marked with 'Ciao force image list'.
 196         """
 197         traces = {}
 198         version = info['Force file list']['Version']
 199         type_re = re.compile('S \[(\w*)\] "([\w\s.]*)"')
 200         for image in info['Ciao force image list']:
 201             offset = int(image['Data offset'])
 202             length = int(image['Data length'])
 203             sample_size = int(image['Bytes/pixel'])
 204             if sample_size != 2:
 205                 raise NotImplementedError('Size: %s' % sample_size)
 206             rows = length / sample_size
 207             d = curve.Data(
 208                 shape=(rows),
 209                 dtype=numpy.int16,
 210                 buffer=buffer,
 211                 offset=offset,
 212                 info=image,
 213                 )
 214             if version in ['0x06120002', '0x06130001']:
 215                 match = type_re.match(image['@4:Image Data'])
 216                 assert match != None, 'Bad regexp for %s, %s' \
 217                     % ('@4:Image Data', image['@4:Image Data'])
 218                 assert match.group(1).lower() == match.group(2).replace(' ','').lower(), \
 219                     'Name missmatch: "%s", "%s"' % (match.group(1), match.group(2))
 220                 tname = match.group(2)
 221             else:
 222                 assert version == '0x07200000', version
 223                 match = type_re.match(image['@4:Image Data'])
 224                 assert match != None, 'Bad regexp for %s, %s' \
 225                     % ('@4:Image Data', image['@4:Image Data'])
 226                 if match.group(1) == 'PulseFreq1':
 227                     assert match.group(2) == 'Freq. 1', match.group(2)
 228                 else:
 229                     assert match.group(1).lower() == match.group(2).replace(' ','').lower(), \
 230                         'Name missmatch: "%s", "%s"' % (match.group(1), match.group(2))
 231                 tname = match.group(2)
 232                 if tname == 'Freq. 1':  # Normalize trace names between versions
 233                     tname = 'Z sensor'
 234                 elif tname == 'Deflection Error':
 235                     tname = 'Deflection'
 236             if tname in traces:
 237                 continue
 238                 msg = None
 239                 if (traces[tname].info != image):
 240                     ik = set(image.keys())
 241                     ok = set(traces[tname].info.keys())
 242                     if ik != ok:
 243                         mmsg = 'extra keys: %s, missing keys %s' % (ik-ok, ok-ik)
 244                     else:
 245                         mmsg = []
 246                         for key in image.keys():
 247                             if image[key] != traces[tname].info[key]:
 248                                 mmsg.append(
 249                                     '%s (%s != %s)'
 250                                     % (key, image[key], traces[tname].info[key]))
 251                         mmsg = ', '.join(mmsg)
 252                     msg = 'info difference: %s' % mmsg
 253                 elif not (traces[tname] == d).all():
 254                     msg = 'data difference'
 255                 if msg != None:
 256                     raise NotImplementedError(
 257                         'Missmatched duplicate traces for %s: %s'
 258                         % (tname, msg))
 259             traces[tname] = d
 260         return traces
 261
 262     def _validate_traces(self, z_piezo, deflection):
 263         if len(z_piezo) != len(deflection):
 264             raise ValueError('Trace length missmatch: %d != %d'
 265                              % (len(z_piezo), len(deflection)))
 266
 267     def _extract_block(self, info, z_piezo, deflection, start, stop, name):
 268         block = curve.Data(
 269             shape=(stop-start, 2),
 270             dtype=numpy.float)
 271         block[:,0] = z_piezo[start:stop]
 272         block[:,1] = deflection[start:stop]
 273         block.info = self._translate_block_info(
 274             info, z_piezo.info, deflection.info, name)
 275         block.info['columns'] = ['z piezo (m)', 'deflection (m)']
 276         block = self._scale_block(block)
 277         return block
 278
 279     def _translate_block_info(self, info, z_piezo_info, deflection_info, name):
 280         version = info['Force file list']['Version']
 281         ret = {
 282             'name': name,
 283             'raw info': info,
 284             'raw z piezo info': z_piezo_info,
 285             'raw deflection info': deflection_info,
 286             'spring constant (N/m)': float(z_piezo_info['Spring Constant']),
 287             }
 288
 289         t = info['Force file list']['Date'] # 04:42:34 PM Tue Sep 11 2007
 290         ret['time'] = time.strptime(t, '%I:%M:%S %p %a %b %d %Y')
 291
 292         volt_re = re.compile(
 293             'V \[Sens. ([\w\s.]*)\] \(([.0-9]*) V/LSB\) (-?[.0-9]*) V')
 294         hz_re = re.compile(
 295             'V \[Sens. ([\w\s.]*)\] \(([.0-9]*) kHz/LSB\) (-?[.0-9]*) kHz')
 296         if version in ['0x06120002', '0x06130001']:
 297             match = volt_re.match(z_piezo_info['@4:Z scale'])
 298             assert match != None, 'Bad regexp for %s, %s' \
 299                 % ('@4:Z scale', z_piezo_info['@4:Z scale'])
 300             assert match.group(1) == 'ZSensorSens', z_piezo_info['@4:Z scale']
 301         else:
 302             assert version == '0x07200000', version
 303             match = hz_re.match(z_piezo_info['@4:Z scale'])
 304             assert match != None, 'Bad regexp for %s, %s' \
 305                 % ('@4:Z scale', z_piezo_info['@4:Z scale'])
 306             assert match.group(1) == 'Freq. 1', z_piezo_info['@4:Z scale']
 307         ret['z piezo sensitivity (V/bit)'] = float(match.group(2))
 308         ret['z piezo range (V)'] = float(match.group(3))
 309         ret['z piezo offset (V)'] = 0.0
 310         # offset assumed if raw data is signed...
 311
 312         match = volt_re.match(deflection_info['@4:Z scale'])
 313         assert match != None, 'Bad regexp for %s, %s' \
 314             % ('@4:Z scale', deflection_info['@4:Z scale'])
 315         assert match.group(1) == 'DeflSens', z_piezo_info['@4:Z scale']
 316         ret['deflection sensitivity (V/bit)'] = float(match.group(2))
 317         ret['deflection range (V)'] = float(match.group(3))
 318         ret['deflection offset (V)'] = 0.0
 319         # offset assumed if raw data is signed...
 320
 321         nm_sens_re = re.compile('V ([.0-9]*) nm/V')
 322         match = nm_sens_re.match(info['Scanner list']['@Sens. Zsens'])
 323         assert match != None, 'Bad regexp for %s/%s, %s' \
 324             % ('Scanner list', '@Sens. Zsens', info['Scanner list']['@4:Z scale'])
 325         ret['z piezo sensitivity (m/V)'] = float(match.group(1))*1e-9
 326
 327         match = nm_sens_re.match(info['Ciao scan list']['@Sens. DeflSens'])
 328         assert match != None, 'Bad regexp for %s/%s, %s' \
 329             % ('Ciao scan list', '@Sens. DeflSens', info['Ciao scan list']['@Sens. DeflSens'])
 330         ret['deflection sensitivity (m/V)'] = float(match.group(1))*1e-9
 331
 332         match = volt_re.match(info['Ciao force list']['@Z scan start'])
 333         assert match != None, 'Bad regexp for %s/%s, %s' \
 334             % ('Ciao force list', '@Z scan start', info['Ciao force list']['@Z scan start'])
 335         ret['z piezo scan (V/bit)'] = float(match.group(2))
 336         ret['z piezo scan start (V)'] = float(match.group(3))
 337
 338         match = volt_re.match(info['Ciao force list']['@Z scan size'])
 339         assert match != None, 'Bad regexp for %s/%s, %s' \
 340             % ('Ciao force list', '@Z scan size', info['Ciao force list']['@Z scan size'])
 341         ret['z piezo scan size (V)'] = float(match.group(3))
 342
 343         const_re = re.compile('C \[([:\w\s]*)\] ([.0-9]*)')
 344         match = const_re.match(z_piezo_info['@Z magnify'])
 345         assert match != None, 'Bad regexp for %s, %s' \
 346             % ('@Z magnify', info['@Z magnify'])
 347         assert match.group(1) == '4:Z scale', match.group(1)
 348         ret['z piezo gain'] = float(match.group(2))
 349
 350         if version in ['0x06120002', '0x06130001']:
 351             match = volt_re.match(z_piezo_info['@4:Z scale'])
 352             assert match != None, 'Bad regexp for %s, %s' \
 353                 % ('@4:Z scale', info['@4:Z scale'])
 354             assert match.group(1) == 'ZSensorSens', match.group(1)
 355             ret['z piezo sensitivity (V/bit)'] = float(match.group(2))
 356             ret['z piezo range (V)'] = float(match.group(3))
 357         else:
 358             assert version == '0x07200000', version
 359             pass
 360
 361         match = volt_re.match(z_piezo_info['@4:Ramp size'])
 362         assert match != None, 'Bad regexp for %s, %s' \
 363             % ('@4:Ramp size', info['@4:Ramp size'])
 364         assert match.group(1) == 'Zsens', match.group(1)
 365         ret['z piezo ramp size (V/bit)'] = float(match.group(2))
 366         ret['z piezo ramp size (V)'] = float(match.group(3))
 367
 368         match = volt_re.match(z_piezo_info['@4:Ramp offset'])
 369         assert match != None, 'Bad regexp for %s, %s' \
 370             % ('@4:Ramp offset', info['@4:Ramp offset'])
 371         assert match.group(1) == 'Zsens', match.group(1)
 372         ret['z piezo ramp offset (V/bit)'] = float(match.group(2))
 373         ret['z piezo ramp offset (V)'] = float(match.group(3))
 374
 375         # Unaccounted for:
 376         #   Samps*
 377
 378         return ret
 379
 380     def _scale_block(self, data):
 381         """Convert the block from its native format to a `numpy.float`
 382         array in SI units.
 383         """
 384         ret = curve.Data(
 385             shape=data.shape,
 386             dtype=numpy.float,
 387             )
 388         info = data.info
 389         ret.info = info
 390         ret.info['raw data'] = data # store the raw data
 391         data.info = {} # break circular reference info <-> data
 392
 393         z_col = info['columns'].index('z piezo (m)')
 394         d_col = info['columns'].index('deflection (m)')
 395
 396         # Leading '-' because Veeco's z increases towards the surface
 397         # (positive indentation), but it makes more sense to me to
 398         # have it increase away from the surface (positive
 399         # separation).
 400         ret[:,z_col] = -(
 401             (data[:,z_col].astype(ret.dtype)
 402              * info['z piezo sensitivity (V/bit)']
 403              - info['z piezo offset (V)'])
 404             * info['z piezo gain']
 405             * info['z piezo sensitivity (m/V)']
 406             )
 407
 408         # Leading '-' because deflection voltage increases as the tip
 409         # moves away from the surface, but it makes more sense to me
 410         # to have it increase as it moves toward the surface (positive
 411         # tension on the protein chain).
 412         ret[:,d_col] = -(
 413             (data[:,d_col]
 414              * info['deflection sensitivity (V/bit)']
 415              - info['deflection offset (V)'])
 416             * info['deflection sensitivity (m/V)']
 417             )
 418
 419         return ret