From 3c5cce3b5d67d72c6db719eee07251690b517606 Mon Sep 17 00:00:00 2001
From: "W. Trevor King" <wking@drexel.edu>
Date: Sat, 7 Aug 2010 18:11:37 -0400
Subject: [PATCH] Create extra block holding PowerSpectrumCommand output data

---
 hooke/plugin/curve.py | 45 ++++++++++++++++++++++++++++++++++++-------
 hooke/util/fft.py     |  7 +------
 2 files changed, 39 insertions(+), 13 deletions(-)

diff --git a/hooke/plugin/curve.py b/hooke/plugin/curve.py
index 3208afa..5ba69d7 100644
--- a/hooke/plugin/curve.py
+++ b/hooke/plugin/curve.py
@@ -32,7 +32,8 @@ from ..plugin import Builtin
 from ..plugin.playlist import current_playlist_callback
 from ..util.calculus import derivative
 from ..util.fft import unitary_avg_power_spectrum
-from ..util.si import ppSI, split_data_label
+from ..util.si import ppSI, join_data_label, split_data_label
+
 
 
 class CurvePlugin (Builtin):
@@ -305,13 +306,21 @@ class PowerSpectrumCommand (Command):
 Data block to act on.  For an approach/retract force curve, `0`
 selects the approaching curve and `1` selects the retracting curve.
 """.strip()),
-                Argument(name='column', type='int', default=1,
+                Argument(name='column', type='string', optional=False,
                          help="""
-Column of data block containing time-series data.
+Name of the data block column containing to-be-transformed data.
+""".strip()),
+                Argument(name='bounds', type='point', optional=True, count=2,
+                         help="""
+Indicies of points bounding the selected data.
 """.strip()),
                 Argument(name='freq', type='float', default=1.0,
                          help="""
 Sampling frequency.
+""".strip()),
+                Argument(name='freq units', type='string', default='Hz',
+                         help="""
+Units for the sampling frequency.
 """.strip()),
                 Argument(name='chunk size', type='int', default=2048,
                          help="""
@@ -321,13 +330,35 @@ Number of samples per chunk.  Use a power of two.
                          help="""
 If `True`, each chunk overlaps the previous chunk by half its length.
 Otherwise, the chunks are end-to-end, and not overlapping.
+""".strip()),
+                Argument(name='output block name', type='string',
+                         default='power spectrum',
+                         help="""
+Name of the new data block (without `of <source block name> <source column name>`) for storing the power spectrum.
 """.strip()),
                 ],
             help=self.__doc__, plugin=plugin)
 
     def _run(self, hooke, inqueue, outqueue, params):
         data = params['curve'].data[params['block']]
-        outqueue.put(unitary_avg_power_spectrum(
-                data[:,params['column']], freq=params['freq'],
-                chunk_size=params['chunk size'],
-                overlap=params['overlap']))
+        col = data.info['columns'].index(params['column'])
+        d = data[:,col]
+        if bounds != None:
+            d = d[params['bounds'][0]:params['bounds'][1]]
+        freq_axis,power = unitary_avg_power_spectrum(
+            d, freq=params['freq'],
+            chunk_size=params['chunk size'],
+            overlap=params['overlap'])
+        name,data_units = split_data_label(params['column'])
+        b = Data((len(freq_axis),2), dtype=data.dtype)
+        b.info['name'] = '%s of %s %s' % (
+            params['output block name'], data.info['name'], params['column'])
+        b.info['columns'] = [
+            join_data_label('frequency axis', params['freq units']),
+            join_data_label('power density',
+                            '%s^2/%s' % (data_units, params['freq units'])),
+            ]
+        b[:,0] = freq_axis
+        b[:,1] = power
+        params['curve'].data.append(b)
+        outqueue.put(b)
diff --git a/hooke/util/fft.py b/hooke/util/fft.py
index 4b13b31..0f7f7b0 100644
--- a/hooke/util/fft.py
+++ b/hooke/util/fft.py
@@ -114,7 +114,6 @@ def unitary_rfft(data, freq=1.0):
 
     Notes
     -----
-    
     If the units on your data are Volts,
     and your sampling frequency is in Hz,
     then `freq_axis` will be in Hz,
@@ -307,7 +306,7 @@ def avg_power_spectrum(data, freq=1.0, chunk_size=2048,
 
 def unitary_avg_power_spectrum(data, freq=1.0, chunk_size=2048,
                                overlap=True, window=window_hann):
-    """Compute the unitary avgerage power spectrum of `data`.
+    """Compute the unitary average power spectrum of `data`.
 
     See Also
     --------
@@ -334,7 +333,6 @@ class TestRFFT (unittest.TestCase):
 
     Notes
     -----
-
     The expected return values are [#dft]_:
 
     .. math:: X_k = \sum_{m=0}^{n-1} x_m \exp^{-2\pi imk/n}
@@ -372,7 +370,6 @@ class TestUnitaryRFFT (unittest.TestCase):
 
         Notes
         -----
-
         Which is:
 
         .. math:: \sum_{m=0}^{n-1} |x_m|^2 dt = \sum_{k=0}^{n-1} |X_k|^2 df
@@ -421,7 +418,6 @@ class TestUnitaryRFFT (unittest.TestCase):
 
         Notes
         -----
-
         Analytic result:
 
         .. math:: \rfft(\rect(at)) = 1/|a|\cdot\sinc(f/a)
@@ -486,7 +482,6 @@ class TestUnitaryRFFT (unittest.TestCase):
 
         Notes
         -----
-
         Analytic result:
 
         .. math::
-- 
2.26.2