From: W. Trevor King Date: Mon, 17 May 2010 18:40:02 +0000 (-0400) Subject: Merged hooke.plugin.plotmanip into hooke.plugin.curve. X-Git-Url: http://git.tremily.us/?p=hooke.git;a=commitdiff_plain;h=8273b2acd0162fd19e79cf52ab3822454d5b2c50;hp=82d0396f036eedbcb7e04ac16a8c8cdf73d564f5 Merged hooke.plugin.plotmanip into hooke.plugin.curve. Also: * Fixed a number of typos in hooke.plugin.curve. * Rewrote derivative code for Numpy in hooke.util.calculus.derivative. * Removed FFT code (I'm replacing it with my FFT_tools module shortly). * Moved v-clamp-specific code into generalvclamp. * Moved plotmanips to hooke.plugin.plotmanip holding area. --- diff --git a/hooke/plugin/curve.py b/hooke/plugin/curve.py index 13ee496..bf21918 100644 --- a/hooke/plugin/curve.py +++ b/hooke/plugin/curve.py @@ -1,5 +1,7 @@ -# Copyright (C) 2010 Fibrin's Benedetti -# W. Trevor King +# Copyright (C) 2008-2010 Alberto Gomez-Casado +# Fabrizio Benedetti +# Massimo Sandal +# W. Trevor King # # This file is part of Hooke. # @@ -23,8 +25,10 @@ associated :class:`hooke.command.Command`\s for handling """ from ..command import Command, Argument, Failure +from ..curve import Data from ..plugin import Builtin from ..plugin.playlist import current_playlist_callback +from ..util.calculus import derivative class CurvePlugin (Builtin): @@ -32,7 +36,7 @@ class CurvePlugin (Builtin): super(CurvePlugin, self).__init__(name='curve') def commands(self): - return [InfoCommand(), ] + return [InfoCommand(), ExportCommand()] # Define common or complicated arguments @@ -114,18 +118,17 @@ class InfoCommand (Command): def _get_block_sizes(self, curve): return [block.shape for block in curve.data] - class ExportCommand (Command): """Export a :class:`hooke.curve.Curve` data block as TAB-delimeted ASCII text. """ def __init__(self): - super(InfoCommand, self).__init__( - name='curve info', + super(ExportCommand, self).__init__( + name='export block', arguments=[ CurveArgument, Argument(name='block', aliases=['set'], type='int', default=0, - help=""" + help=""" Data block to save. For an approach/retract force curve, `0` selects the approacing curve and `1` selects the retracting curve. """.strip()), @@ -137,7 +140,83 @@ File name for the output data. Defaults to 'curve.dat' help=self.__doc__) def _run(self, hooke, inqueue, outqueue, params): - data = params['curve'].data[params['index']] + data = params['curve'].data[params['block']] f = open(params['output'], 'w') data.tofile(f, sep='\t') f.close() + +class DifferenceCommand (Command): + """Calculate the derivative (actually, the discrete differentiation) + of a curve data block. + + See :func:`hooke.util.calculus.derivative` for implementation + details. + """ + def __init__(self): + super(DifferenceCommand, self).__init__( + name='block difference', + arguments=[ + CurveArgument, + Argument(name='block one', aliases=['set one'], type='int', + default=1, + help=""" +Block A in A-B. For an approach/retract force curve, `0` selects the +approacing curve and `1` selects the retracting curve. +""".strip()), + Argument(name='block two', aliases=['set two'], type='int', + default=0, + help='Block B in A-B.'), + Argument(name='x column', type='int', default=0, + help=""" +Column of data block to differentiate with respect to. +""".strip()), + Argument(name='y column', type='int', default=1, + help=""" +Column of data block to differentiate. +""".strip()), + ], + help=self.__doc__) + + def _run(self, hooke, inqueue, outqueue, params): + a = params['curve'].data[params['block one']] + b = params['curve'].data[params['block two']] + assert a[:,params['x column']] == b[:,params['x column']]: + out = Data((a.shape[0],2), dtype=a.dtype) + out[:,0] = a[:,params['x column']] + out[:,1] = a[:,params['y column']] - b[:,params['y column']]: + outqueue.put(out) + +class DerivativeCommand (Command): + """Calculate the difference between two blocks of data. + """ + def __init__(self): + super(DerivativeCommand, self).__init__( + name='block derivative', + arguments=[ + CurveArgument, + Argument(name='block', aliases=['set'], type='int', default=0, + help=""" +Data block to differentiate. For an approach/retract force curve, `0` +selects the approacing curve and `1` selects the retracting curve. +""".strip()), + Argument(name='x column', type='int', default=0, + help=""" +Column of data block to differentiate with respect to. +""".strip()), + Argument(name='y column', type='int', default=1, + help=""" +Column of data block to differentiate. +""".strip()), + Argument(name='weights', type='dict', default={-1:-0.5, 1:0.5}, + help=""" +Weighting scheme dictionary for finite differencing. Defaults to +central differencing. +""".strip()), + ], + help=self.__doc__) + + def _run(self, hooke, inqueue, outqueue, params): + data = params['curve'].data[params['block']] + outqueue.put(derivative( + block, x_col=params['x column'], y_col=params['y column'], + weights=params['weights'])) diff --git a/hooke/plugin/generalvclamp.py b/hooke/plugin/generalvclamp.py index a672143..d152608 100644 --- a/hooke/plugin/generalvclamp.py +++ b/hooke/plugin/generalvclamp.py @@ -39,6 +39,25 @@ warnings.simplefilter('ignore',np.RankWarning) class generalvclampCommands(object): + def do_subtplot(self, args): + ''' + SUBTPLOT + (procplots.py plugin) + Plots the difference between ret and ext current curve + ------- + Syntax: subtplot + ''' + #FIXME: sub_filter and sub_order must be args + + if len(self.plots[0].vectors) != 2: + print 'This command only works on a curve with two different plots.' + pass + + outplot=self.subtract_curves(sub_order=1) + + plot_graph=self.list_of_events['plot_graph'] + wx.PostEvent(self.frame,plot_graph(plots=[outplot])) + def _plug_init(self): self.basecurrent=None self.basepoints=None diff --git a/hooke/plugin/plotmanip.py b/hooke/plugin/plotmanip.py new file mode 100644 index 0000000..0455c7a --- /dev/null +++ b/hooke/plugin/plotmanip.py @@ -0,0 +1,92 @@ + +class Plotmanip (object): +#-----PLOT MANIPULATORS + def plotmanip_median(self, plot, current, customvalue=None): + ''' + does the median of the y values of a plot + ''' + if customvalue: + median_filter=customvalue + else: + median_filter=self.config['medfilt'] + + if median_filter==0: + return plot + + if float(median_filter)/2 == int(median_filter)/2: + median_filter+=1 + + nplots=len(plot.vectors) + c=0 + while c -# W. Trevor King -# -# This file is part of Hooke. -# -# Hooke is free software: you can redistribute it and/or -# modify it under the terms of the GNU Lesser General Public -# License as published by the Free Software Foundation, either -# version 3 of the License, or (at your option) any later version. -# -# Hooke is distributed in the hope that it will be useful, -# but WITHOUT ANY WARRANTY; without even the implied warranty of -# MERCHANTABILITY or FITNESS FOR A PARTICULAR PURPOSE. See the -# GNU Lesser General Public License for more details. -# -# You should have received a copy of the GNU Lesser General Public -# License along with Hooke. If not, see -# . - -"""Processed plots plugin for force curves. -""" - -from ..libhooke import WX_GOOD -import wxversion -wxversion.select(WX_GOOD) - -import wx -import numpy as np -import scipy as sp -import scipy.signal -import copy - -from .. import curve as lhc - - -class procplotsCommands(object): - - def _plug_init(self): - pass - - def do_derivplot(self,args): - ''' - DERIVPLOT - (procplots.py plugin) - Plots the derivate (actually, the discrete differentiation) of the current force curve - --------- - Syntax: derivplot - ''' - dplot=self.derivplot_curves() - plot_graph=self.list_of_events['plot_graph'] - wx.PostEvent(self.frame,plot_graph(plots=[dplot])) - - def derivplot_curves(self): - ''' - do_derivplot helper function - - create derivate plot curves for force curves. - ''' - dplot=lhc.PlotObject() - dplot.vectors=[] - - for vector in self.plots[0].vectors: - dplot.vectors.append([]) - dplot.vectors[-1].append(vector[0][:-1]) - dplot.vectors[-1].append(np.diff(vector[1])) - - dplot.destination=1 - dplot.units=self.plots[0].units - - return dplot - - def do_subtplot(self, args): - ''' - SUBTPLOT - (procplots.py plugin) - Plots the difference between ret and ext current curve - ------- - Syntax: subtplot - ''' - #FIXME: sub_filter and sub_order must be args - - if len(self.plots[0].vectors) != 2: - print 'This command only works on a curve with two different plots.' - pass - - outplot=self.subtract_curves(sub_order=1) - - plot_graph=self.list_of_events['plot_graph'] - wx.PostEvent(self.frame,plot_graph(plots=[outplot])) - - def subtract_curves(self, sub_order): - ''' - subtracts the extension from the retraction - ''' - xext=self.plots[0].vectors[0][0] - yext=self.plots[0].vectors[0][1] - xret=self.plots[0].vectors[1][0] - yret=self.plots[0].vectors[1][1] - - #we want the same number of points - maxpoints_tot=min(len(xext),len(xret)) - xext=xext[0:maxpoints_tot] - yext=yext[0:maxpoints_tot] - xret=xret[0:maxpoints_tot] - yret=yret[0:maxpoints_tot] - - if sub_order: - ydiff=[yretval-yextval for yretval,yextval in zip(yret,yext)] - else: #reverse subtraction (not sure it's useful, but...) - ydiff=[yextval-yretval for yextval,yretval in zip(yext,yret)] - - outplot=copy.deepcopy(self.plots[0]) - outplot.vectors[0][0], outplot.vectors[1][0] = xext,xret #FIXME: if I use xret, it is not correct! - outplot.vectors[1][1]=ydiff - outplot.vectors[0][1]=[0 for item in outplot.vectors[1][0]] - - return outplot - - -#-----PLOT MANIPULATORS - def plotmanip_median(self, plot, current, customvalue=None): - ''' - does the median of the y values of a plot - ''' - if customvalue: - median_filter=customvalue - else: - median_filter=self.config['medfilt'] - - if median_filter==0: - return plot - - if float(median_filter)/2 == int(median_filter)/2: - median_filter+=1 - - nplots=len(plot.vectors) - c=0 - while cstartvalue: - cutindex=index - else: - break - - plot.vectors[0][0]=plot.vectors[0][0][:cutindex] - plot.vectors[0][1]=plot.vectors[0][1][:cutindex] - - return plot - ''' - - - -#FFT--------------------------- - def fft_plot(self, vector): - ''' - calculates the fast Fourier transform for the selected vector in the plot - ''' - fftplot=lhc.PlotObject() - fftplot.vectors=[[]] - - fftlen=len(vector)/2 #need just 1/2 of length - fftplot.vectors[-1].append(np.arange(1,fftlen).tolist()) - - try: - fftplot.vectors[-1].append(abs(np.fft(vector)[1:fftlen]).tolist()) - except TypeError: #we take care of newer NumPy (1.0.x) - fftplot.vectors[-1].append(abs(np.fft.fft(vector)[1:fftlen]).tolist()) - - - fftplot.destination=1 - - - return fftplot - - - def do_fft(self,args): - ''' - FFT - (procplots.py plugin) - Plots the fast Fourier transform of the selected plot - --- - Syntax: fft [top,bottom] [select] [0,1...] - - By default, fft performs the Fourier transform on all the 0-th data set on the - top plot. - - [top,bottom]: which plot is the data set to fft (default: top) - [select]: you pick up two points on the plot and fft only the segment between - [0,1,...]: which data set on the selected plot you want to fft (default: 0) - ''' - - #args parsing - #whatplot = plot to fft - #whatset = set to fft in the plot - select=('select' in args) - if 'top' in args: - whatplot=0 - elif 'bottom' in args: - whatplot=1 - else: - whatplot=0 - whatset=0 - for arg in args: - try: - whatset=int(arg) - except ValueError: - pass - - if select: - points=self._measure_N_points(N=2, whatset=whatset) - boundaries=[points[0].index, points[1].index] - boundaries.sort() - y_to_fft=self.plots[whatplot].vectors[whatset][1][boundaries[0]:boundaries[1]] #y - else: - y_to_fft=self.plots[whatplot].vectors[whatset][1] #y - - fftplot=self.fft_plot(y_to_fft) - fftplot.units=['frequency', 'power'] - plot_graph=self.list_of_events['plot_graph'] - wx.PostEvent(self.frame,plot_graph(plots=[fftplot])) diff --git a/hooke/util/calculus.py b/hooke/util/calculus.py new file mode 100644 index 0000000..504840a --- /dev/null +++ b/hooke/util/calculus.py @@ -0,0 +1,106 @@ +# Copyright + +"""The `calculus` module provides functions for calculating +derivatives and integrals of discrete data. +""" + +import copy + +import numpy + +from ..curve import Data + + +def derivative(data, x_col=0, f_col=1, weights={-1:-0.5, 1:0.5}): + """Calculate the discrete derivative (finite difference) of + data[:,f_col] with respect to data[:,x_col]. + + Examples + -------- + + >>> import pprint + >>> d = Data((5,2), dtype=numpy.float, + ... info={'columns':['x', 'x**2']}) + >>> for i in range(5): + ... d[i,0] = i + ... d[i,1] = i**2 + >>> d + Data([[ 0., 0.], + [ 1., 1.], + [ 2., 4.], + [ 3., 9.], + [ 4., 16.]]) + >>> dd = derivative(d) + >>> dd + Data([[ 0., 1.], + [ 1., 2.], + [ 2., 4.], + [ 3., 6.], + [ 4., 7.]]) + >>> pprint.pprint(dd.info) + {'columns': ['x', 'deriv x**2 with respect to x']} + + Notes + ----- + + Weights + ~~~~~~~ + + The returned :class:`Data` block shares its x vector with the + input data. The ith df/dx value in the returned data is + caclulated with:: + + (df/dx)[i] = (SUM_j w[j] f[i+j]) / h + + where ``h = x[i+1]-x[i]`` is the x coordinate spacing (assumed + constant) and ``j`` ranges over the keys of `weights`. + + There standard schemes translate as follows: + + ======== ====================== =================== + scheme formula weights + ======== ====================== =================== + forward ``(f[i+1]-f[i])/h`` ``{0:-1,1:1}`` + backward ``(f[i]-f[i-1])/h`` ``{0:1,-1:-1}`` + central ``(f[i+1]-f[i-1])/2h`` ``{-1:-0.5,1:0.5}`` + ======== ====================== =================== + + The default scheme is central differencing. + + Boundary conditions + ~~~~~~~~~~~~~~~~~~~ + + These could be configurable in principle. The current scheme just + extrapolates virtual points out to negative `i` following:: + + f[i<0] = 2*f[0] - f[-i] + + With analogous treatment for `i > data.shape[0]`. This ensures that + `f[i]-f[0]` is odd about `i=0`, which keeps derivatives smooth.:: + + f[i] - f[0] = f[0] - f[-i] == -(f[-i] - f[0]) + """ + output = Data((data.shape[0],2), dtype=data.dtype) + output.info = copy.copy(data.info) + output.info['columns'] = [ + data.info['columns'][x_col], + 'deriv %s with respect to %s' \ + % (data.info['columns'][f_col], data.info['columns'][x_col]), + ] + h = data[1,x_col] - data[0,x_col] + chunks = [] + for i,w in weights.items(): + chunk = numpy.roll(w*data[:,f_col], -i) + if i > 0: # chunk shifted down, replace the high `i`s + zero = len(chunk) - 1 - i + for j in range(1,i+1): + chunk[zero+j] = 2*chunk[zero] - chunk[zero-j] + elif i < 0: # chunk shifted up, replace the low `i`s + zero = -i + for j in range(1,zero+1): + chunk[zero-j] = 2*chunk[zero] - chunk[zero+j] + chunks.append(chunk) + print chunks + output[:,0] = data[:,x_col] + output[:,1] = sum(chunks) + return output