hooke/plugin/curve.py

   1 # Copyright (C) 2008-2010 Alberto Gomez-Kasai
   2 #                         Fabiano's Benedetti
   3 #                         Massimo Sandal <devicerandom@gmail.com>
   4 #                         W. Trevor King <wking@drexel.edu>
   5 #
   6 # This file is part of Hooke.
   7 #
   8 # Hooke is free software: you can redistribute it and/or modify it
   9 # under the terms of the GNU Lesser General Public License as
  10 # published by the Free Software Foundation, either version 3 of the
  11 # License, or (at your option) any later version.
  12 #
  13 # Hooke is distributed in the hope that it will be useful, but WITHOUT
  14 # ANY WARRANTY; without even the implied warranty of MERCHANTABILITY
  15 # or FITNESS FOR A PARTICULAR PURPOSE.  See the GNU Lesser General
  16 # Public License for more details.
  17 #
  18 # You should have received a copy of the GNU Lesser General Public
  19 # License along with Hooke.  If not, see
  20 # <http://www.gnu.org/licenses/>.
  21
  22 """The ``curve`` module provides :class:`CurvePlugin` and several
  23 associated :class:`hooke.command.Command`\s for handling
  24 :mod:`hooke.curve` classes.
  25 """
  26
  27 import copy
  28
  29 import numpy
  30
  31 from ..command import Command, Argument, Failure
  32 from ..curve import Data
  33 from ..plugin import Builtin
  34 from ..plugin.playlist import current_playlist_callback
  35 from ..util.calculus import derivative
  36 from ..util.fft import unitary_avg_power_spectrum
  37 from ..util.si import ppSI, join_data_label, split_data_label
  38
  39
  40 # Define common or complicated arguments
  41
  42 def current_curve_callback(hooke, command, argument, value):
  43     if value != None:
  44         return value
  45     playlist = current_playlist_callback(hooke, command, argument, value)
  46     curve = playlist.current()
  47     if curve == None:
  48         raise Failure('No curves in %s' % playlist)
  49     return curve
  50
  51 CurveArgument = Argument(
  52     name='curve', type='curve', callback=current_curve_callback,
  53     help="""
  54 :class:`hooke.curve.Curve` to act on.  Defaults to the current curve
  55 of the current playlist.
  56 """.strip())
  57
  58 def _name_argument(name, default, help):
  59     """TODO
  60     """
  61     return Argument(name=name, type='string', default=default, help=help)
  62
  63 def block_argument(*args, **kwargs):
  64     """TODO
  65     """
  66     return _name_argument(*args, **kwargs)
  67
  68 def column_argument(*args, **kwargs):
  69     """TODO
  70     """
  71     return _name_argument(*args, **kwargs)
  72
  73
  74 # Define useful command subclasses
  75
  76 class CurveCommand (Command):
  77     """A :class:`~hooke.command.Command` operating on a
  78     :class:`~hooke.curve.Curve`.
  79     """
  80     def __init__(self, **kwargs):
  81         if 'arguments' in kwargs:
  82             kwargs['arguments'].insert(0, CurveArgument)
  83         else:
  84             kwargs['arguments'] = [CurveArgument]
  85         super(CurveCommand, self).__init__(**kwargs)
  86
  87     def _curve(self, hooke, params):
  88         # HACK? rely on params['curve'] being bound to the local hooke
  89         # playlist (i.e. not a copy, as you would get by passing a
  90         # curve through the queue).  Ugh.  Stupid queues.  As an
  91         # alternative, we could pass lookup information through the
  92         # queue...
  93         return params['curve']
  94
  95
  96 class BlockCommand (CurveCommand):
  97     """A :class:`CurveCommand` operating on a :class:`~hooke.curve.Data` block.
  98     """
  99     def __init__(self, blocks=None, **kwargs):
 100         if blocks == None:
 101             blocks = [('block', None, 'Name of the data block to act on.')]
 102         block_args = []
 103         for name,default,help in blocks:
 104             block_args.append(block_argument(name, default, help))
 105         self._block_arguments = block_args
 106         if 'arguments' not in kwargs:
 107             kwargs['arguments'] = []
 108         kwargs['arguments'] = block_args + kwargs['arguments']
 109         super(BlockCommand, self).__init__(**kwargs)
 110
 111     def _block_names(self, hooke, params):
 112         curve = self._curve(hooke, params)
 113         return [b.info['name'] for b in curve.data]
 114
 115     def _block_index(self, hooke, params, name=None):
 116         if name == None:
 117             name = self._block_arguments[0].name
 118         block_name = params[name]
 119         if block_name == None:
 120             curve = self._curve(hooke=hooke, params=params)
 121             if len(curve.data) == 0:
 122                 raise Failure('no blocks in %s' % curve)
 123             block_name = curve.data[0].info['name']
 124         names = self._block_names(hooke=hooke, params=params)
 125         try:
 126             return names.index(block_name)
 127         except ValueError, e:
 128             curve = self._curve(hooke, params)
 129             raise Failure('no block named %s in %s (%s): %s'
 130                           % (block_name, curve, names, e))
 131
 132     def _block(self, hooke, params, name=None):
 133         # HACK? rely on params['block'] being bound to the local hooke
 134         # playlist (i.e. not a copy, as you would get by passing a
 135         # block through the queue).  Ugh.  Stupid queues.  As an
 136         # alternative, we could pass lookup information through the
 137         # queue...
 138         curve = self._curve(hooke, params)
 139         index = self._block_index(hooke, params, name)
 140         return curve.data[index]
 141
 142
 143 class ColumnAccessCommand (BlockCommand):
 144     """A :class:`BlockCommand` accessing a :class:`~hooke.curve.Data`
 145     block column.
 146     """
 147     def __init__(self, columns=None, **kwargs):
 148         if columns == None:
 149             columns = [('column', None, 'Name of the data column to act on.')]
 150         column_args = []
 151         for name,default,help in columns:
 152             column_args.append(column_argument(name, default, help))
 153         self._column_arguments = column_args
 154         if 'arguments' not in kwargs:
 155             kwargs['arguments'] = []
 156         kwargs['arguments'] = column_args + kwargs['arguments']
 157         super(ColumnAccessCommand, self).__init__(**kwargs)
 158
 159     def _get_column(self, hooke, params, block_name=None, column_name=None):
 160         if column_name == None:
 161             column_name = self._column_arguments[0].name
 162         column_name = params[column_name]
 163         block = self._block(hooke, params, block_name)
 164         columns = block.info['columns']
 165         try:
 166             column_index = columns.index(column_name)
 167         except ValueError, e:
 168             raise Failure('%s not in %s (%s): %s'
 169                           % (column_name, block.info['name'], columns, e))
 170         return block[:,column_index]
 171
 172
 173 class ColumnAddingCommand (ColumnAccessCommand):
 174     """A :class:`ColumnAccessCommand` that also adds columns.
 175     """
 176     def __init__(self, new_columns=None, **kwargs):
 177         if new_columns == None:
 178             new_columns = []
 179         column_args = []
 180         for name,default,help in new_columns:
 181             column_args.append(column_argument(name, default, help))
 182         self._new_column_arguments = column_args
 183         if 'arguments' not in kwargs:
 184             kwargs['arguments'] = []
 185         kwargs['arguments'] = column_args + kwargs['arguments']
 186         super(ColumnAddingCommand, self).__init__(**kwargs)
 187
 188     def _get_column(self, hooke, params, block_name=None, column_name=None):
 189         if column_name == None and len(self._column_arguments) == 0:
 190             column_name = self._new_column_arguments[0].name
 191         return super(ColumnAddingCommand, self)._get_column(
 192             hooke=hooke, params=params, block_name=block_name,
 193             column_name=column_name)
 194
 195     def _set_column(self, hooke, params, block_name=None, column_name=None,
 196                     values=None):
 197         if column_name == None:
 198             column_name = self._column_arguments[0].name
 199         column_name = params[column_name]
 200         block = self._block(hooke=hooke, params=params, name=block_name)
 201         if column_name not in block.info['columns']:
 202             new = Data((block.shape[0], block.shape[1]+1), dtype=block.dtype)
 203             new.info = copy.deepcopy(block.info)
 204             new[:,:-1] = block
 205             new.info['columns'].append(column_name)
 206             block = new
 207             block_index = self._block_index(hooke, params, name=block_name)
 208             self._curve(hooke, params).data[block_index] = block
 209         column_index = block.info['columns'].index(column_name)
 210         block[:,column_index] = values
 211
 212
 213 # The plugin itself
 214
 215 class CurvePlugin (Builtin):
 216     def __init__(self):
 217         super(CurvePlugin, self).__init__(name='curve')
 218         self._commands = [
 219             GetCommand(self), InfoCommand(self), DeltaCommand(self),
 220             ExportCommand(self), DifferenceCommand(self),
 221             DerivativeCommand(self), PowerSpectrumCommand(self)]
 222
 223
 224 # Define commands
 225
 226 class GetCommand (CurveCommand):
 227     """Return a :class:`hooke.curve.Curve`.
 228     """
 229     def __init__(self, plugin):
 230         super(GetCommand, self).__init__(
 231             name='get curve', help=self.__doc__, plugin=plugin)
 232
 233     def _run(self, hooke, inqueue, outqueue, params):
 234         outqueue.put(self._curve(hooke, params))
 235
 236
 237 class InfoCommand (CurveCommand):
 238     """Get selected information about a :class:`hooke.curve.Curve`.
 239     """
 240     def __init__(self, plugin):
 241         args = [
 242             Argument(name='all', type='bool', default=False, count=1,
 243                      help='Get all curve information.'),
 244             ]
 245         self.fields = ['name', 'path', 'experiment', 'driver', 'filetype', 'note',
 246                        'blocks', 'block sizes']
 247         for field in self.fields:
 248             args.append(Argument(
 249                     name=field, type='bool', default=False, count=1,
 250                     help='Get curve %s' % field))
 251         super(InfoCommand, self).__init__(
 252             name='curve info', arguments=args,
 253             help=self.__doc__, plugin=plugin)
 254
 255     def _run(self, hooke, inqueue, outqueue, params):
 256         curve = self._curve(hooke, params)
 257         fields = {}
 258         for key in self.fields:
 259             fields[key] = params[key]
 260         if reduce(lambda x,y: x and y, fields.values()) == False:
 261             params['all'] = True # No specific fields set, default to 'all'
 262         if params['all'] == True:
 263             for key in self.fields:
 264                 fields[key] = True
 265         lines = []
 266         for key in self.fields:
 267             if fields[key] == True:
 268                 get = getattr(self, '_get_%s' % key.replace(' ', '_'))
 269                 lines.append('%s: %s' % (key, get(curve)))
 270         outqueue.put('\n'.join(lines))
 271
 272     def _get_name(self, curve):
 273         return curve.name
 274
 275     def _get_path(self, curve):
 276         return curve.path
 277
 278     def _get_experiment(self, curve):
 279         return curve.info.get('experiment', None)
 280
 281     def _get_driver(self, curve):
 282         return curve.driver
 283
 284     def _get_filetype(self, curve):
 285         return curve.info.get('filetype', None)
 286
 287     def _get_note(self, curve):
 288         return curve.info.get('note', None)
 289
 290     def _get_blocks(self, curve):
 291         return len(curve.data)
 292
 293     def _get_block_sizes(self, curve):
 294         return [block.shape for block in curve.data]
 295
 296
 297 class DeltaCommand (BlockCommand):
 298     """Get distance information between two points.
 299
 300     With two points A and B, the returned distances are A-B.
 301     """
 302     def __init__(self, plugin):
 303         super(DeltaCommand, self).__init__(
 304             name='delta',
 305             arguments=[
 306                 Argument(name='point', type='point', optional=False, count=2,
 307                          help="""
 308 Indicies of points bounding the selected data.
 309 """.strip()),
 310                 Argument(name='SI', type='bool', default=False,
 311                          help="""
 312 Return distances in SI notation.
 313 """.strip())
 314                 ],
 315             help=self.__doc__, plugin=plugin)
 316
 317     def _run(self, hooke, inqueue, outqueue, params):
 318         data = self._block(hooke, params)
 319         As = data[params['point'][0],:]
 320         Bs = data[params['point'][1],:]
 321         ds = [A-B for A,B in zip(As, Bs)]
 322         if params['SI'] == False:
 323             out = [(name, d) for name,d in zip(data.info['columns'], ds)]
 324         else:
 325             out = []
 326             for name,d in zip(data.info['columns'], ds):
 327                 n,units = split_data_label(name)
 328                 out.append(
 329                   (n, ppSI(value=d, unit=units, decimals=2)))
 330         outqueue.put(out)
 331
 332
 333 class ExportCommand (BlockCommand):
 334     """Export a :class:`hooke.curve.Curve` data block as TAB-delimeted
 335     ASCII text.
 336
 337     A "#" prefixed header will optionally appear at the beginning of
 338     the file naming the columns.
 339     """
 340     def __init__(self, plugin):
 341         super(ExportCommand, self).__init__(
 342             name='export block',
 343             arguments=[
 344                 Argument(name='output', type='file', default='curve.dat',
 345                          help="""
 346 File name for the output data.  Defaults to 'curve.dat'
 347 """.strip()),
 348                 Argument(name='header', type='bool', default=True,
 349                          help="""
 350 True if you want the column-naming header line.
 351 """.strip()),
 352                 ],
 353             help=self.__doc__, plugin=plugin)
 354
 355     def _run(self, hooke, inqueue, outqueue, params):
 356         data = self._block(hooke, params)
 357
 358         with open(params['output'], 'w') as f:
 359             if params['header'] == True:
 360                 f.write('# %s \n' % ('\t'.join(data.info['columns'])))
 361             numpy.savetxt(f, data, delimiter='\t')
 362
 363
 364 class DifferenceCommand (ColumnAddingCommand):
 365     """Calculate the difference between two columns of data.
 366
 367     The difference is added to block A as a new column.
 368
 369     Note that the command will fail if the columns have different
 370     lengths, so be careful when differencing columns from different
 371     blocks.
 372     """
 373     def __init__(self, plugin):
 374         super(DifferenceCommand, self).__init__(
 375             name='difference',
 376             blocks=[
 377                 ('block A', None,
 378                  'Name of block A in A-B.  Defaults to the first block'),
 379                 ('block B', None,
 380                  'Name of block B in A-B.  Defaults to matching `block A`.'),
 381                 ],
 382             columns=[
 383                 ('column A', None,
 384                  """
 385 Column of data from block A to difference.  Defaults to the first column.
 386 """.strip()),
 387                 ('column B', None,
 388                  """
 389 Column of data from block B to difference.  Defaults to matching `column A`.
 390 """.strip()),
 391                 ],
 392             new_columns=[
 393                 ('output column', None,
 394                  """
 395 Name of the new column for storing the difference (without units, defaults to
 396 `difference of <block A> <column A> and <block B> <column B>`).
 397 """.strip()),
 398                 ],
 399             help=self.__doc__, plugin=plugin)
 400
 401     def _run(self, hooke, inqueue, outqueue, params):
 402         params = self.__setup_params(hooke=hooke, params=params)
 403         data_A = self._get_column(hooke=hooke, params=params,
 404                                   block_name='block A',
 405                                   column_name='column A')
 406         data_B = self._get_column(hooke=hooke, params=params,
 407                                   block_name='block B',
 408                                   column_name='column B')
 409         out = data_A - data_B
 410         self._set_column(hooke=hooke, params=params,
 411                          block_name='block A',
 412                          column_name='output column',
 413                          values=out)
 414
 415     def __setup_params(self, hooke, params):
 416         curve = self._curve(hooke, params)
 417         if params['block A'] == None:
 418             params['block A'] = curve.data[0].info['name']
 419         if params['block B'] == None:
 420             params['block B'] = params['block A']
 421         block_A = self._block(hooke, params=params, name='block A')
 422         block_B = self._block(hooke, params=params, name='block B')
 423         if params['column A'] == None:
 424             params['column A'] = block.info['columns'][0]
 425         if params['column B'] == None:
 426             params['column B'] = params['column A']
 427         a_name,a_unit = split_data_label(params['column A'])
 428         b_name,b_unit = split_data_label(params['column B'])
 429         if a_unit != b_unit:
 430             raise Failure('Unit missmatch: %s != %s' % (a_unit, b_unit))
 431         if params['output column'] == None:
 432             params['output column'] = join_data_label(
 433                 'difference of %s %s and %s %s' % (
 434                     block_A.info['name'], a_name,
 435                     block_B.info['name'], b_name),
 436                 a_unit)
 437         else:
 438             params['output column'] = join_data_label(
 439                 params['output column'], a_unit)
 440         return params
 441
 442
 443 class DerivativeCommand (ColumnAddingCommand):
 444     """Calculate the derivative (actually, the discrete differentiation)
 445     of a data column.
 446
 447     See :func:`hooke.util.calculus.derivative` for implementation
 448     details.
 449     """
 450     def __init__(self, plugin):
 451         super(DerivativeCommand, self).__init__(
 452             name='derivative',
 453             columns=[
 454                 ('x column', None,
 455                  'Column of data block to differentiate with respect to.'),
 456                 ('f column', None,
 457                  'Column of data block to differentiate.'),
 458                 ],
 459             new_columns=[
 460                 ('output column', None,
 461                  """
 462 Name of the new column for storing the derivative (without units, defaults to
 463 `derivative of <f column> with respect to <x column>`).
 464 """.strip()),
 465                 ],
 466             arguments=[
 467                 Argument(name='weights', type='dict', default={-1:-0.5, 1:0.5},
 468                          help="""
 469 Weighting scheme dictionary for finite differencing.  Defaults to
 470 central differencing.
 471 """.strip()),
 472                 ],
 473             help=self.__doc__, plugin=plugin)
 474
 475     def _run(self, hooke, inqueue, outqueue, params):
 476         params = self.__setup_params(hooke=hooke, params=params)
 477         x_data = self._get_column(hooke=hooke, params=params,
 478                                   column_name='x column')
 479         f_data = self._get_column(hooke=hooke, params=params,
 480                                   column_name='f column')
 481         d = derivative(
 482             x_data=x_data, f_data=f_data, weights=params['weights'])
 483         self._set_column(hooke=hooke, params=params,
 484                          column_name='output column',
 485                          values=d)
 486
 487     def __setup_params(self, hooke, params):
 488         curve = self._curve(hooke, params)
 489         x_name,x_unit = split_data_label(params['x column'])
 490         f_name,f_unit = split_data_label(params['f column'])
 491         d_unit = '%s/%s' % (f_unit, x_unit)
 492         if params['output column'] == None:
 493             params['output column'] = join_data_label(
 494                 'derivative of %s with respect to %s' % (
 495                     f_name, x_name),
 496                 d_unit)
 497         else:
 498             params['output column'] = join_data_label(
 499                 params['output column'], d_unit)
 500         return params
 501
 502
 503 class PowerSpectrumCommand (ColumnAddingCommand):
 504     """Calculate the power spectrum of a data column.
 505     """
 506     def __init__(self, plugin):
 507         super(PowerSpectrumCommand, self).__init__(
 508             name='power spectrum',
 509             arguments=[
 510                 Argument(name='output block', type='string',
 511                          help="""
 512 Name of the new data block for storing the power spectrum (defaults to
 513 `power spectrum of <source block name> <source column name>`).
 514 """.strip()),
 515                 Argument(name='bounds', type='point', optional=True, count=2,
 516                          help="""
 517 Indicies of points bounding the selected data.
 518 """.strip()),
 519                 Argument(name='freq', type='float', default=1.0,
 520                          help="""
 521 Sampling frequency.
 522 """.strip()),
 523                 Argument(name='freq units', type='string', default='Hz',
 524                          help="""
 525 Units for the sampling frequency.
 526 """.strip()),
 527                 Argument(name='chunk size', type='int', default=2048,
 528                          help="""
 529 Number of samples per chunk.  Use a power of two.
 530 """.strip()),
 531                 Argument(name='overlap', type='bool', default=False,
 532                          help="""
 533 If `True`, each chunk overlaps the previous chunk by half its length.
 534 Otherwise, the chunks are end-to-end, and not overlapping.
 535 """.strip()),
 536                 ],
 537             help=self.__doc__, plugin=plugin)
 538
 539     def _run(self, hooke, inqueue, outqueue, params):
 540         params = self.__setup_params(hooke=hooke, params=params)
 541         data = self._get_column(hooke=hooke, params=params)
 542         bounds = params['bounds']
 543         if bounds != None:
 544             data = data[bounds[0]:bounds[1]]
 545         freq_axis,power = unitary_avg_power_spectrum(
 546             data, freq=params['freq'],
 547             chunk_size=params['chunk size'],
 548             overlap=params['overlap'])
 549         b = Data((len(freq_axis),2), dtype=data.dtype)
 550         b.info['name'] = params['output block']
 551         b.info['columns'] = [
 552             params['output freq column'],
 553             params['output power column'],
 554             ]
 555         self._curve(hooke, params).data.append(b)
 556         self._set_column(hooke, params, block_name='output block',
 557                          column_name='output freq column',
 558                          values=freq_axis)
 559         self._set_column(hooke, params, block_name='output block',
 560                          column_name='output power column',
 561                          values=power)
 562         outqueue.put(b)
 563
 564     def __setup_params(self, hooke, params):
 565         if params['output block'] in self._block_names(hooke, params):
 566             raise Failure('output block %s already exists in %s.'
 567                           % (params['output block'],
 568                              self._curve(hooke, params)))
 569         data = self._get_column(hooke=hooke, params=params)
 570         d_name,d_unit = split_data_label(data.info['name'])
 571         if params['output block'] == None:
 572             params['output block'] = 'power spectrum of %s %s' % (
 573                 data.info['name'], params['column'])
 574         self.params['output freq column'] = join_data_label(
 575             'frequency axis', params['freq units'])
 576         self.params['output power column'] = join_data_label(
 577             'power density', '%s^2/%s' % (data_units, params['freq units']))
 578         return params
 579
 580
 581 class OldCruft (object):
 582
 583     def do_forcebase(self,args):
 584         '''
 585         FORCEBASE
 586         (generalvclamp.py)
 587         Measures the difference in force (in pN) between a point and a baseline
 588         took as the average between two points.
 589
 590         The baseline is fixed once for a given curve and different force measurements,
 591         unless the user wants it to be recalculated
 592         ------------
 593         Syntax: forcebase [rebase]
 594                 rebase: Forces forcebase to ask again the baseline
 595                 max: Instead of asking for a point to measure, asks for two points and use
 596                      the maximum peak in between
 597         '''
 598         rebase=False #if true=we select rebase
 599         maxpoint=False #if true=we measure the maximum peak
 600
 601         plot=self._get_displayed_plot()
 602         whatset=1 #fixme: for all sets
 603         if 'rebase' in args or (self.basecurrent != self.current.path):
 604             rebase=True
 605         if 'max' in args:
 606             maxpoint=True
 607
 608         if rebase:
 609             print 'Select baseline'
 610             self.basepoints=self._measure_N_points(N=2, whatset=whatset)
 611             self.basecurrent=self.current.path
 612
 613         if maxpoint:
 614             print 'Select two points'
 615             points=self._measure_N_points(N=2, whatset=whatset)
 616             boundpoints=[points[0].index, points[1].index]
 617             boundpoints.sort()
 618             try:
 619                 y=min(plot.vectors[whatset][1][boundpoints[0]:boundpoints[1]])
 620             except ValueError:
 621                 print 'Chosen interval not valid. Try picking it again. Did you pick the same point as begin and end of interval?'
 622         else:
 623             print 'Select point to measure'
 624             points=self._measure_N_points(N=1, whatset=whatset)
 625             #whatplot=points[0].dest
 626             y=points[0].graph_coords[1]
 627
 628         #fixme: code duplication
 629         boundaries=[self.basepoints[0].index, self.basepoints[1].index]
 630         boundaries.sort()
 631         to_average=plot.vectors[whatset][1][boundaries[0]:boundaries[1]] #y points to average
 632
 633         avg=np.mean(to_average)
 634         forcebase=abs(y-avg)
 635         print str(forcebase*(10**12))+' pN'
 636         to_dump='forcebase '+self.current.path+' '+str(forcebase*(10**12))+' pN'
 637         self.outlet.push(to_dump)
 638
 639     #---SLOPE---
 640     def do_slope(self,args):
 641         '''
 642         SLOPE
 643         (generalvclamp.py)
 644         Measures the slope of a delimited chunk on the return trace.
 645         The chunk can be delimited either by two manual clicks, or have
 646         a fixed width, given as an argument.
 647         ---------------
 648         Syntax: slope [width]
 649                 The facultative [width] parameter specifies how many
 650                 points will be considered for the fit. If [width] is
 651                 specified, only one click will be required.
 652         (c) Marco Brucale, Massimo Sandal 2008
 653         '''
 654
 655         # Reads the facultative width argument
 656         try:
 657             fitspan=int(args)
 658         except:
 659             fitspan=0
 660
 661         # Decides between the two forms of user input, as per (args)
 662         if fitspan == 0:
 663             # Gets the Xs of two clicked points as indexes on the current curve vector
 664             print 'Click twice to delimit chunk'
 665             points=self._measure_N_points(N=2,whatset=1)
 666         else:
 667             print 'Click once on the leftmost point of the chunk (i.e.usually the peak)'
 668             points=self._measure_N_points(N=1,whatset=1)
 669
 670         slope=self._slope(points,fitspan)
 671
 672         # Outputs the relevant slope parameter
 673         print 'Slope:'
 674         print str(slope)
 675         to_dump='slope '+self.current.path+' '+str(slope)
 676         self.outlet.push(to_dump)
 677
 678     def _slope(self,points,fitspan):
 679         # Calls the function linefit_between
 680         parameters=[0,0,[],[]]
 681         try:
 682             clickedpoints=[points[0].index,points[1].index]
 683             clickedpoints.sort()
 684         except:
 685             clickedpoints=[points[0].index-fitspan,points[0].index]
 686
 687         try:
 688             parameters=self.linefit_between(clickedpoints[0],clickedpoints[1])
 689         except:
 690             print 'Cannot fit. Did you click twice the same point?'
 691             return
 692
 693         # Outputs the relevant slope parameter
 694         print 'Slope:'
 695         print str(parameters[0])
 696         to_dump='slope '+self.curve.path+' '+str(parameters[0])
 697         self.outlet.push(to_dump)
 698
 699         # Makes a vector with the fitted parameters and sends it to the GUI
 700         xtoplot=parameters[2]
 701         ytoplot=[]
 702         x=0
 703         for x in xtoplot:
 704             ytoplot.append((x*parameters[0])+parameters[1])
 705
 706         clickvector_x, clickvector_y=[], []
 707         for item in points:
 708             clickvector_x.append(item.graph_coords[0])
 709             clickvector_y.append(item.graph_coords[1])
 710
 711         lineplot=self._get_displayed_plot(0) #get topmost displayed plot
 712
 713         lineplot.add_set(xtoplot,ytoplot)
 714         lineplot.add_set(clickvector_x, clickvector_y)
 715
 716
 717         if lineplot.styles==[]:
 718             lineplot.styles=[None,None,None,'scatter']
 719         else:
 720             lineplot.styles+=[None,'scatter']
 721         if lineplot.colors==[]:
 722             lineplot.colors=[None,None,'black',None]
 723         else:
 724             lineplot.colors+=['black',None]
 725
 726
 727         self._send_plot([lineplot])
 728
 729         return parameters[0]
 730
 731
 732     def linefit_between(self,index1,index2,whatset=1):
 733         '''
 734         Creates two vectors (xtofit,ytofit) slicing out from the
 735         current return trace a portion delimited by the two indexes
 736         given as arguments.
 737         Then does a least squares linear fit on that slice.
 738         Finally returns [0]=the slope, [1]=the intercept of the
 739         fitted 1st grade polynomial, and [2,3]=the actual (x,y) vectors
 740         used for the fit.
 741         (c) Marco Brucale, Massimo Sandal 2008
 742         '''
 743         # Translates the indexes into two vectors containing the x,y data to fit
 744         xtofit=self.plots[0].vectors[whatset][0][index1:index2]
 745         ytofit=self.plots[0].vectors[whatset][1][index1:index2]
 746
 747         # Does the actual linear fitting (simple least squares with numpy.polyfit)
 748         linefit=[]
 749         linefit=np.polyfit(xtofit,ytofit,1)
 750
 751         return (linefit[0],linefit[1],xtofit,ytofit)