1 # Copyright (C) 2011-2012 W. Trevor King <wking@drexel.edu>
3 # This file is part of pypiezo.
5 # pypiezo is free software: you can redistribute it and/or modify it under the
6 # terms of the GNU General Public License as published by the Free Software
7 # Foundation, either version 3 of the License, or (at your option) any later
10 # pypiezo is distributed in the hope that it will be useful, but WITHOUT ANY
11 # WARRANTY; without even the implied warranty of MERCHANTABILITY or FITNESS FOR
12 # A PARTICULAR PURPOSE. See the GNU General Public License for more details.
14 # You should have received a copy of the GNU General Public License along with
15 # pypiezo. If not, see <http://www.gnu.org/licenses/>.
17 "Basic piezo control."
20 from time import sleep as _sleep
22 import numpy as _numpy
23 from scipy.stats import linregress as _linregress
26 import matplotlib as _matplotlib
27 import matplotlib.pyplot as _matplotlib_pyplot
28 import time as _time # for timestamping lines on plots
29 except (ImportError, RuntimeError), e:
31 _matplotlib_import_error = e
33 from pycomedi.constant import AREF, TRIG_SRC, SDF
34 from pycomedi.utility import inttrig_insn, Reader, Writer
36 from . import LOG as _LOG
37 from . import config as _config
38 from . import package_config as _package_config
41 def convert_bits_to_volts(config, data):
42 """Convert bit-valued data to volts.
44 >>> config = _config.ChannelConfig()
45 >>> config['conversion-coefficients'] = [1, 2, 3]
46 >>> config['conversion-origin'] = -1
47 >>> convert_bits_to_volts(config, -1)
49 >>> convert_bits_to_volts(
50 ... config, _numpy.array([-1, 0, 1, 2], dtype=_numpy.float))
51 array([ 1., 6., 17., 34.])
53 coefficients = config['conversion-coefficients']
54 origin = config['conversion-origin']
55 return _numpy.polyval(list(reversed(coefficients)), data-origin)
57 def convert_volts_to_bits(config, data):
58 """Convert bit-valued data to volts.
60 >>> config = _config.ChannelConfig()
61 >>> config['inverse-conversion-coefficients'] = [1, 2, 3]
62 >>> config['inverse-conversion-origin'] = -1
63 >>> convert_volts_to_bits(config, -1)
65 >>> convert_volts_to_bits(
66 ... config, _numpy.array([-1, 0, 1, 2], dtype=_numpy.float))
67 array([ 1., 6., 17., 34.])
69 Note that the inverse coeffiecient and offset are difficult to
70 derive from the forward coefficient and offset. The current
71 Comedilib conversion functions, `comedi_to_physical()` and
72 `comedi_from_physical()` take `comedi_polynomial_t` conversion
73 arguments. `comedi_polynomial_t` is defined in `comedilib.h`_,
74 and holds a polynomial of length
75 `COMEDI_MAX_NUM_POLYNOMIAL_COEFFICIENTS`, currently set to 4. The
76 inverse of this cubic polynomial might be another polynomial, but
77 it might also have a more complicated form. A general inversion
78 solution is considered too complicated, so when you're setting up
79 your configuration, you should use Comedilib to save both the
80 forward and inverse coefficients and offsets.
82 .. _comedilib.h: http://www.comedi.org/git?p=comedi/comedilib.git;a=blob;f=include/comedilib.h;hb=HEAD
84 origin = config['inverse-conversion-origin']
85 inverse_coefficients = config['inverse-conversion-coefficients']
86 if len(inverse_coefficients) == 0:
87 raise NotImplementedError('cubic polynomial inversion')
88 return _numpy.polyval(list(reversed(inverse_coefficients)), data-origin)
90 def convert_volts_to_meters(config, data):
91 """Convert volt-valued data to meters.
93 >>> config = _config.AxisConfig()
94 >>> config['gain'] = 20.0
95 >>> config['sensitivity'] = 8e-9
96 >>> convert_volts_to_meters(config, 1)
98 >>> convert_volts_to_meters(
99 ... config, _numpy.array([1, 6, 17, 34], dtype=_numpy.float))
100 ... # doctest: +ELLIPSIS
101 array([ 1.6...e-07, 9.6...e-07, 2.7...e-06,
104 return data * config['gain'] * config['sensitivity']
106 def convert_meters_to_volts(config, data):
107 """Convert bit-valued data to volts.
109 >>> config = _config.AxisConfig()
110 >>> config['gain'] = 20.0
111 >>> config['sensitivity'] = 8e-9
112 >>> convert_meters_to_volts(config, 1.6e-7)
114 >>> convert_meters_to_volts(
115 ... config, _numpy.array([1.6e-7, 9.6e-7, 2.72e-6, 5.44e-6],
116 ... dtype=_numpy.float))
117 array([ 1., 6., 17., 34.])
119 return data / (config['gain'] * config['sensitivity'])
121 def convert_bits_to_meters(axis_config, data):
122 """Convert bit-valued data to meters.
124 >>> channel_config = _config.ChannelConfig()
125 >>> channel_config['conversion-coefficients'] = [1, 2, 3]
126 >>> channel_config['conversion-origin'] = -1
127 >>> axis_config = _config.AxisConfig()
128 >>> axis_config['gain'] = 20.0
129 >>> axis_config['sensitivity'] = 8e-9
130 >>> axis_config['channel'] = channel_config
131 >>> convert_bits_to_meters(axis_config, 1)
132 ... # doctest: +ELLIPSIS
134 >>> convert_bits_to_meters(
135 ... axis_config, _numpy.array([-1, 0, 1, 2], dtype=_numpy.float))
136 ... # doctest: +ELLIPSIS
137 array([ 1.6...e-07, 9.6...e-07, 2.7...e-06,
140 data = convert_bits_to_volts(axis_config['channel'], data)
141 return convert_volts_to_meters(axis_config, data)
143 def convert_meters_to_bits(axis_config, data):
144 """Convert meter-valued data to volts.
146 >>> channel_config = _config.ChannelConfig()
147 >>> channel_config['inverse-conversion-coefficients'] = [1, 2, 3]
148 >>> channel_config['inverse-conversion-origin'] = -1
149 >>> axis_config = _config.AxisConfig()
150 >>> axis_config['gain'] = 20.0
151 >>> axis_config['sensitivity'] = 8e-9
152 >>> axis_config['channel'] = channel_config
153 >>> convert_meters_to_bits(axis_config, 1.6e-7)
155 >>> convert_meters_to_bits(
157 ... _numpy.array([1.6e-7, 9.6e-7, 2.72e-6, 5.44e-6],
158 ... dtype=_numpy.float))
159 array([ 17., 162., 1009., 3746.])
161 data = convert_meters_to_volts(axis_config, data)
162 return convert_volts_to_bits(axis_config['channel'], data)
164 def get_axis_name(axis_config):
165 """Return the name of an axis from the `AxisConfig`
167 This is useful, for example, with
168 `Config.select_config(get_attribute=get_axis_name)`.
170 channel_config = axis_config['channel']
171 return channel_config['name']
173 def _setup_channel_config(config, channel):
174 """Initialize the `ChannelConfig` `config` using the
175 `AnalogChannel` `channel`.
177 config['device'] = channel.subdevice.device.filename
178 config['subdevice'] = channel.subdevice.index
179 config['channel'] = channel.index
180 config['maxdata'] = channel.get_maxdata()
181 config['range'] = channel.range.value
182 config['analog-reference'] = AREF.index_by_value(channel.aref.value)
183 converter = channel.get_converter()
184 config['conversion-origin'
185 ] = converter.get_to_physical_expansion_origin()
186 config['conversion-coefficients'
187 ] = converter.get_to_physical_coefficients()
188 config['inverse-conversion-origin'
189 ] = converter.get_from_physical_expansion_origin()
190 config['inverse-conversion-coefficients'
191 ] = converter.get_from_physical_coefficients()
194 class PiezoAxis (object):
195 """A one-dimensional piezoelectric actuator.
197 If used, the montoring channel must (as of now) be on the same
198 device as the controlling channel.
200 >>> from pycomedi.device import Device
201 >>> from pycomedi.subdevice import StreamingSubdevice
202 >>> from pycomedi.channel import AnalogChannel
203 >>> from pycomedi.constant import (AREF, SUBDEVICE_TYPE, UNIT)
205 >>> d = Device('/dev/comedi0')
208 >>> s_in = d.find_subdevice_by_type(SUBDEVICE_TYPE.ai,
209 ... factory=StreamingSubdevice)
210 >>> s_out = d.find_subdevice_by_type(SUBDEVICE_TYPE.ao,
211 ... factory=StreamingSubdevice)
213 >>> axis_channel = s_out.channel(
214 ... 0, factory=AnalogChannel, aref=AREF.ground)
215 >>> monitor_channel = s_in.channel(
216 ... 0, factory=AnalogChannel, aref=AREF.diff)
217 >>> for chan in [axis_channel, monitor_channel]:
218 ... chan.range = chan.find_range(unit=UNIT.volt, min=-10, max=10)
220 >>> config = _config.AxisConfig()
221 >>> config.update({'gain':20, 'sensitivity':8e-9})
222 >>> config['channel'] = _config.OutputChannelConfig()
223 >>> config['monitor'] = _config.InputChannelConfig()
224 >>> config['monitor']['device'] = '/dev/comediX'
226 >>> p = PiezoAxis(config=config)
227 ... # doctest: +NORMALIZE_WHITESPACE
228 Traceback (most recent call last):
230 NotImplementedError: piezo axis control and monitor on different devices
231 (/dev/comedi0 and /dev/comediX)
233 >>> config['monitor']['device'] = config['channel']['device']
234 >>> p = PiezoAxis(config=config,
235 ... axis_channel=axis_channel, monitor_channel=monitor_channel)
238 >>> print(config['channel'].dump())
245 analog-reference: ground
246 conversion-coefficients: -10.0, 0.000305180437934
247 conversion-origin: 0.0
248 inverse-conversion-coefficients: 0.0, 3276.75
249 inverse-conversion-origin: -10.0
250 >>> print(config['monitor'].dump())
257 analog-reference: diff
258 conversion-coefficients: -10.0, 0.000305180437934
259 conversion-origin: 0.0
260 inverse-conversion-coefficients: 0.0, 3276.75
261 inverse-conversion-origin: -10.0
263 >>> convert_bits_to_meters(p.config, 0)
264 ... # doctest: +ELLIPSIS
269 def __init__(self, config, axis_channel=None, monitor_channel=None):
271 if (config['monitor'] and
272 config['channel']['device'] != config['monitor']['device']):
273 raise NotImplementedError(
274 ('piezo axis control and monitor on different devices '
276 config['channel']['device'],
277 config['monitor']['device']))
279 raise NotImplementedError(
280 'pypiezo not yet capable of opening its own axis channel')
281 #axis_channel = pycomedi...
282 self.axis_channel = axis_channel
283 if config['monitor'] and not monitor_channel:
284 raise NotImplementedError(
285 'pypiezo not yet capable of opening its own monitor channel')
286 #monitor_channel = pycomedi...
287 self.monitor_channel = monitor_channel
288 self.name = config['channel']['name']
290 def setup_config(self):
291 "Initialize the axis (and monitor) configs."
292 _setup_channel_config(self.config['channel'], self.axis_channel)
293 if self.monitor_channel:
294 _setup_channel_config(
295 self.config['monitor'], self.monitor_channel)
296 if self.config['minimum'] is None:
297 self.config['minimum'] = convert_bits_to_volts(
298 self.config['channel'], 0)
299 if self.config['maximum'] is None:
300 self.config['maximum'] = convert_bits_to_volts(
301 self.config['channel'], self.axis_channel.get_maxdata())
304 class InputChannel(object):
305 """An input channel monitoring some interesting parameter.
307 >>> from pycomedi.device import Device
308 >>> from pycomedi.subdevice import StreamingSubdevice
309 >>> from pycomedi.channel import AnalogChannel
310 >>> from pycomedi.constant import AREF, SUBDEVICE_TYPE, UNIT
312 >>> d = Device('/dev/comedi0')
315 >>> s = d.find_subdevice_by_type(SUBDEVICE_TYPE.ai,
316 ... factory=StreamingSubdevice)
318 >>> channel = s.channel(0, factory=AnalogChannel, aref=AREF.diff)
319 >>> channel.range = channel.find_range(unit=UNIT.volt, min=-10, max=10)
321 >>> channel_config = _config.InputChannelConfig()
323 >>> c = InputChannel(config=channel_config, channel=channel)
325 >>> print(channel_config.dump())
332 analog-reference: diff
333 conversion-coefficients: -10.0, 0.000305180437934
334 conversion-origin: 0.0
335 inverse-conversion-coefficients: 0.0, 3276.75
336 inverse-conversion-origin: -10.0
338 >>> convert_bits_to_volts(c.config, 0)
343 def __init__(self, config, channel=None):
346 raise NotImplementedError(
347 'pypiezo not yet capable of opening its own channel')
348 #channel = pycomedi...
349 self.channel = channel
350 self.name = config['name']
352 def setup_config(self):
353 _setup_channel_config(self.config, self.channel)
356 class Piezo (object):
357 """A piezo actuator-controlled experiment.
359 >>> from pprint import pprint
360 >>> from pycomedi.device import Device
361 >>> from pycomedi.subdevice import StreamingSubdevice
362 >>> from pycomedi.channel import AnalogChannel
363 >>> from pycomedi.constant import AREF, SUBDEVICE_TYPE, UNIT
365 >>> d = Device('/dev/comedi0')
368 >>> s_in = d.find_subdevice_by_type(SUBDEVICE_TYPE.ai,
369 ... factory=StreamingSubdevice)
370 >>> s_out = d.find_subdevice_by_type(SUBDEVICE_TYPE.ao,
371 ... factory=StreamingSubdevice)
373 >>> axis_channel = s_out.channel(
374 ... 0, factory=AnalogChannel, aref=AREF.ground)
375 >>> monitor_channel = s_in.channel(
376 ... 0, factory=AnalogChannel, aref=AREF.diff)
377 >>> input_channel = s_in.channel(1, factory=AnalogChannel, aref=AREF.diff)
378 >>> for chan in [axis_channel, monitor_channel, input_channel]:
379 ... chan.range = chan.find_range(unit=UNIT.volt, min=-10, max=10)
381 >>> axis_config = _config.AxisConfig()
382 >>> axis_config.update({'gain':20, 'sensitivity':8e-9})
383 >>> axis_config['channel'] = _config.OutputChannelConfig()
384 >>> axis_config['channel']['name'] = 'z'
385 >>> axis_config['monitor'] = _config.InputChannelConfig()
386 >>> input_config = _config.InputChannelConfig()
387 >>> input_config['name'] = 'some-input'
389 >>> a = PiezoAxis(config=axis_config, axis_channel=axis_channel,
390 ... monitor_channel=monitor_channel)
393 >>> c = InputChannel(config=input_config, channel=input_channel)
396 >>> p = Piezo(axes=[a], inputs=[c], name='Charlie')
397 >>> inputs = p.read_inputs()
398 >>> pprint(inputs) # doctest: +SKIP
399 {'some-input': 34494L, 'z-monitor': 32669L}
401 >>> pos = convert_volts_to_bits(p.config['axes'][0]['channel'], 0)
405 >>> p.last_output == {'z': int(pos)}
408 :meth:`ramp` raises an error if passed an invalid `data` `dtype`.
410 >>> output_data = _numpy.arange(0, int(pos), step=int(pos/10),
411 ... dtype=_numpy.float)
412 >>> output_data = output_data.reshape((len(output_data), 1))
413 >>> input_data = p.ramp(data=output_data, frequency=10,
414 ... output_names=['z'], input_names=['z-monitor', 'some-input'])
415 Traceback (most recent call last):
417 ValueError: output dtype float64 does not match expected <type 'numpy.uint16'>
418 >>> output_data = _numpy.arange(0, int(pos), step=int(pos/10),
419 ... dtype=p.channel_dtype('z', direction='output'))
420 >>> output_data = output_data.reshape((len(output_data), 1))
421 >>> input_data = p.ramp(data=output_data, frequency=10,
422 ... output_names=['z'], input_names=['z-monitor', 'some-input'])
423 >>> input_data # doctest: +SKIP
434 [32646, 22686]], dtype=uint16)
436 >>> p.last_output == {'z': output_data[-1]}
439 >>> data = p.named_ramp(data=output_data, frequency=10,
440 ... output_names=['z'], input_names=['z-monitor', 'some-input'])
441 >>> pprint(data) # doctest: +ELLIPSIS, +SKIP
442 {'some-input': array([21666, 20566, ..., 22395], dtype=uint16),
443 'z': array([ 0, 3276, ..., 32760], dtype=uint16),
444 'z-monitor': array([ 3102, 6384, ..., 32647], dtype=uint16)}
448 def __init__(self, axes, inputs, name=None):
451 self.config = _config.PiezoConfig()
453 self.config['name'] = name
454 self.config['axes'] = [x.config for x in axes]
455 self.config['inputs'] = [x.config for x in inputs]
456 self.last_output = {}
458 def axis_by_name(self, name):
459 "Get an axis by its name."
460 for axis in self.axes:
461 if axis.name == name:
463 raise ValueError(name)
465 def input_channel_by_name(self, name):
466 "Get an input channel by its name."
467 for input_channel in self.inputs:
468 if input_channel.name == name:
470 raise ValueError(name)
472 def channels(self, direction=None):
473 """Iterate through all `(name, channel)` tuples.
475 =========== ===================
476 `direction` Returned channels
477 =========== ===================
478 'input' all input channels
479 'output' all output channels
481 =========== ===================
483 if direction not in ('input', 'output', None):
484 raise ValueError(direction)
486 if direction != 'input':
487 yield (a.name, a.axis_channel)
488 if a.monitor_channel and direction != 'output':
489 yield ('%s-monitor' % a.name, a.monitor_channel)
490 if direction != 'output':
491 for c in self.inputs:
492 yield (c.name, c.channel)
494 def channel_by_name(self, name, direction=None):
495 """Get a channel by its name.
497 Setting `direction` (see :meth:`channels`) may allow a more
500 for n,channel in self.channels(direction=direction):
503 raise ValueError(name)
505 def channel_dtype(self, channel_name, direction=None):
506 """Get a channel's data type by name.
508 Setting `direction` (see :meth:`channels`) may allow a more
511 channel = self.channel_by_name(name=channel_name, direction=direction)
512 return channel.subdevice.get_dtype()
514 def read_inputs(self):
515 "Read all inputs and return a `name`->`value` dictionary."
516 # There is no multi-channel read instruction, so preform reads
518 ret = dict([(n, c.data_read())
519 for n,c in self.channels(direction='input')])
520 _LOG.debug('current position: %s' % ret)
523 def jump(self, axis_name, position, steps=1, sleep=None):
524 "Move the output named `axis_name` to `position`."
525 _LOG.debug('jump %s to %s in %d steps' % (axis_name, position, steps))
528 orig_pos = self.last_output[axis_name]
531 ("cannot make a soft jump to {} because we don't have a "
532 'last-output position for {}').format(
533 position, axis_name))
536 for pos in _numpy.linspace(orig_pos, position, steps+1)[1:]:
537 self.jump(axis_name=axis_name, position=pos)
541 position = int(position)
542 channel = self.channel_by_name(name=axis_name)
543 channel.data_write(position)
544 self.last_output[axis_name] = position
546 def ramp(self, data, frequency, output_names, input_names=()):
547 """Synchronized IO ramp writing `data` and reading `in_names`.
551 data : numpy array-like
552 Row for each cycle, column for each output channel.
554 Target cycle frequency in Hz.
555 output_names : list of strings
556 Names of output channels in the same order as the columns
558 input_names : list of strings
559 Names of input channels to monitor in the same order as
560 the columns of the returned array.
562 if len(data.shape) != 2:
564 'ramp data must be two dimensional, not %d' % len(data.shape))
565 if data.shape[1] != len(output_names):
567 'ramp data should have on column for each input, '
568 'but has %d columns for %d inputs'
569 % (data.shape[1], len(output_names)))
570 n_samps = data.shape[0]
571 log_string = 'ramp %d samples at %g Hz. out: %s, in: %s' % (
572 n_samps, frequency, output_names, input_names)
573 _LOG.debug(log_string) # _LOG on one line for easy commenting-out
575 output_channels = [self.channel_by_name(name=n, direction='output')
576 for n in output_names]
577 inputs = [self.channel_by_name(name=n, direction='input')
578 for n in input_names]
580 ao_subdevice = output_channels[0].subdevice
581 ai_subdevice = inputs[0].subdevice
582 device = ao_subdevice.device
584 output_dtype = ao_subdevice.get_dtype()
585 if data.dtype != output_dtype:
586 raise ValueError('output dtype %s does not match expected %s'
587 % (data.dtype, output_dtype))
588 input_data = _numpy.ndarray(
589 (n_samps, len(inputs)), dtype=ai_subdevice.get_dtype())
591 _LOG.debug('setup ramp commands')
592 scan_period_ns = int(1e9 / frequency)
593 ai_cmd = ai_subdevice.get_cmd_generic_timed(
594 len(inputs), scan_period_ns)
595 ao_cmd = ao_subdevice.get_cmd_generic_timed(
596 len(output_channels), scan_period_ns)
598 ai_cmd.start_src = TRIG_SRC.int
600 ai_cmd.stop_src = TRIG_SRC.count
601 ai_cmd.stop_arg = n_samps
602 ai_cmd.chanlist = inputs
603 #ao_cmd.start_src = TRIG_SRC.ext
604 #ao_cmd.start_arg = 18 # NI card AI_START1 internal AI start signal
605 ao_cmd.start_src = TRIG_SRC.int
607 ao_cmd.stop_src = TRIG_SRC.count
608 ao_cmd.stop_arg = n_samps-1
609 ao_cmd.chanlist = output_channels
611 ai_subdevice.cmd = ai_cmd
612 ao_subdevice.cmd = ao_cmd
614 rc = ai_subdevice.command_test()
616 _LOG.debug('analog input test %d: %s' % (i, rc))
618 rc = ao_subdevice.command_test()
620 _LOG.debug('analog output test %d: %s' % (i, rc))
622 _LOG.debug('lock subdevices for ramp')
627 _LOG.debug('load ramp commands')
628 ao_subdevice.command()
629 ai_subdevice.command()
631 writer = Writer(ao_subdevice, data)
633 reader = Reader(ai_subdevice, input_data)
636 _LOG.debug('arm analog output')
637 device.do_insn(inttrig_insn(ao_subdevice))
638 _LOG.debug('trigger ramp (via analog input)')
639 device.do_insn(inttrig_insn(ai_subdevice))
640 _LOG.debug('ramp running')
644 _LOG.debug('ramp complete')
646 #_LOG.debug('AI flags: %s' % ai_subdevice.get_flags())
647 ai_subdevice.cancel()
648 ai_subdevice.unlock()
650 # release busy flag, which seems to not be cleared
652 # http://groups.google.com/group/comedi_list/browse_thread/thread/4c7040989197abad/
653 #_LOG.debug('AO flags: %s' % ao_subdevice.get_flags())
654 ao_subdevice.cancel()
655 ao_subdevice.unlock()
656 _LOG.debug('unlocked subdevices after ramp')
658 for i,name in enumerate(output_names):
659 self.last_output[name] = data[-1,i]
661 if _package_config['matplotlib']:
663 raise _matplotlib_import_error
664 figure = _matplotlib_pyplot.figure()
665 axes = figure.add_subplot(1, 1, 1)
667 timestamp = _time.strftime('%H%M%S')
668 axes.set_title('piezo ramp %s' % timestamp)
669 for d,names in [(data, output_names),
670 (input_data, input_names)]:
671 for i,name in enumerate(names):
672 axes.plot(d[:,i], label=name)
676 def named_ramp(self, data, frequency, output_names, input_names=()):
677 input_data = self.ramp(
678 data=data, frequency=frequency, output_names=output_names,
679 input_names=input_names)
681 for i,name in enumerate(output_names):
682 ret[name] = data[:,i]
683 for i,name in enumerate(input_names):
684 ret[name] = input_data[:,i]