[python-dtcwt] 185/497: a first stab at moving functionality to a backend class

Ghislain Vaillant ghisvail-guest at moszumanska.debian.org
Tue Jul 21 18:06:02 UTC 2015


This is an automated email from the git hooks/post-receive script.

ghisvail-guest pushed a commit to branch debian/sid
in repository python-dtcwt.

commit 40ee56a31026f0bec8aba3776227e945bd595103
Author: Rich Wareham <rjw57 at cam.ac.uk>
Date:   Mon Nov 11 12:46:23 2013 +0000

    a first stab at moving functionality to a backend class
    
    Move the lowlevel filtering code to a NumPy backend class
---
 dtcwt/backend/numpy/lowlevel.py | 444 +++++++++++++++++-----------------------
 dtcwt/lowlevel.py               |  75 ++++++-
 dtcwt/opencl/lowlevel.py        |   2 +-
 dtcwt/opencl/transform2d.py     |   2 +-
 dtcwt/sampling.py               |   2 +-
 dtcwt/transform1d.py            |   3 +-
 dtcwt/transform2d.py            |   3 +-
 dtcwt/transform3d.py            |   3 +-
 dtcwt/utils.py                  |  69 ++++++-
 tests/testreflect.py            |   2 +-
 10 files changed, 341 insertions(+), 264 deletions(-)

diff --git a/dtcwt/backend/numpy/lowlevel.py b/dtcwt/backend/numpy/lowlevel.py
index 878e767..f5858cd 100644
--- a/dtcwt/backend/numpy/lowlevel.py
+++ b/dtcwt/backend/numpy/lowlevel.py
@@ -1,43 +1,6 @@
 import numpy as np
 from six.moves import xrange
-
-def asfarray(X):
-    """Similar to :py:func:`numpy.asfarray` except that this function tries to
-    preserve the original datatype of X if it is already a floating point type
-    and will pass floating point arrays through directly without copying.
-
-    """
-    X = np.asanyarray(X)
-    return np.asfarray(X, dtype=X.dtype)
-
-def appropriate_complex_type_for(X):
-    """Return an appropriate complex data type depending on the type of X. If X
-    is already complex, return that, if it is floating point return a complex
-    type of the appropriate size and if it is integer, choose an complex
-    floating point type depending on the result of :py:func:`numpy.asfarray`.
-
-    """
-    X = asfarray(X)
-    
-    if np.issubsctype(X.dtype, np.complex64) or np.issubsctype(X.dtype, np.complex128):
-        return X.dtype
-    elif np.issubsctype(X.dtype, np.float32):
-        return np.complex64
-    elif np.issubsctype(X.dtype, np.float64):
-        return np.complex128
-
-    # God knows, err on the side of caution
-    return np.complex128
-
-def as_column_vector(v):
-    """Return *v* as a column vector with shape (N,1).
-    
-    """
-    v = np.atleast_2d(v)
-    if v.shape[0] == 1:
-        return v.T
-    else:
-        return v
+from dtcwt.utils import as_column_vector, asfarray, appropriate_complex_type_for, reflect
 
 def _centered(arr, newsize):
     # Return the center newsize portion of the array.
@@ -77,250 +40,221 @@ def _column_convolve(X, h):
     outShape[0] = abs(X.shape[0] - h_size) + 1
     return _centered(out, outShape)
 
-def reflect(x, minx, maxx):
-    """Reflect the values in matrix *x* about the scalar values *minx* and
-    *maxx*.  Hence a vector *x* containing a long linearly increasing series is
-    converted into a waveform which ramps linearly up and down between *minx* and
-    *maxx*.  If *x* contains integers and *minx* and *maxx* are (integers + 0.5), the
-    ramps will have repeated max and min samples.
-   
-    .. codeauthor:: Rich Wareham <rjw57 at cantab.net>, Aug 2013
-    .. codeauthor:: Nick Kingsbury, Cambridge University, January 1999.
-    
-    """
-
-    # Copy x to avoid in-place modification
-    y = np.array(x, copy=True)
-
-    # Reflect y in maxx.
-    t = y > maxx
-    y[t] = (2*maxx - y[t]).astype(y.dtype)
-
-    while np.any(y < minx):
-        # Reflect y in minx.
-        t = y < minx
-        y[t] = (2*minx - y[t]).astype(y.dtype)
+class LowLevelBackendNumPy(object):
+    def colfilter(self, X, h):
+        """Filter the columns of image *X* using filter vector *h*, without decimation.
+        If len(h) is odd, each output sample is aligned with each input sample
+        and *Y* is the same size as *X*.  If len(h) is even, each output sample is
+        aligned with the mid point of each pair of input samples, and Y.shape =
+        X.shape + [1 0].
+
+        :param X: an image whose columns are to be filtered
+        :param h: the filter coefficients.
+        :returns Y: the filtered image.
+
+        .. codeauthor:: Rich Wareham <rjw57 at cantab.net>, August 2013
+        .. codeauthor:: Cian Shaffrey, Cambridge University, August 2000
+        .. codeauthor:: Nick Kingsbury, Cambridge University, August 2000
+
+        """
+        
+        # Interpret all inputs as arrays
+        X = asfarray(X)
+        h = as_column_vector(h)
+
+        r, c = X.shape
+        m = h.shape[0]
+        m2 = np.fix(m*0.5)
+
+        # Symmetrically extend with repeat of end samples.
+        # Use 'reflect' so r < m2 works OK.
+        xe = reflect(np.arange(-m2, r+m2, dtype=np.int), -0.5, r-0.5)
 
-        # Reflect y in maxx.
-        t = y > maxx
-        y[t] = (2*maxx - y[t]).astype(y.dtype)
+        # Perform filtering on the columns of the extended matrix X(xe,:), keeping
+        # only the 'valid' output samples, so Y is the same size as X if m is odd.
+        Y = _column_convolve(X[xe,:], h)
 
-    return y
+        return Y
 
-def colfilter(X, h):
-    """Filter the columns of image *X* using filter vector *h*, without decimation.
-    If len(h) is odd, each output sample is aligned with each input sample
-    and *Y* is the same size as *X*.  If len(h) is even, each output sample is
-    aligned with the mid point of each pair of input samples, and Y.shape =
-    X.shape + [1 0].
+    def coldfilt(self, X, ha, hb):
+        """Filter the columns of image X using the two filters ha and hb =
+        reverse(ha).  ha operates on the odd samples of X and hb on the even
+        samples.  Both filters should be even length, and h should be approx linear
+        phase with a quarter sample advance from its mid pt (i.e. :math:`|h(m/2)| >
+        |h(m/2 + 1)|`).
 
-    :param X: an image whose columns are to be filtered
-    :param h: the filter coefficients.
-    :returns Y: the filtered image.
+        .. code-block:: text
 
-    .. codeauthor:: Rich Wareham <rjw57 at cantab.net>, August 2013
-    .. codeauthor:: Cian Shaffrey, Cambridge University, August 2000
-    .. codeauthor:: Nick Kingsbury, Cambridge University, August 2000
+                              ext        top edge                     bottom edge       ext
+            Level 1:        !               |               !               |               !
+            odd filt on .    b   b   b   b   a   a   a   a   a   a   a   a   b   b   b   b   
+            odd filt on .      a   a   a   a   b   b   b   b   b   b   b   b   a   a   a   a
+            Level 2:        !               |               !               |               !
+            +q filt on x      b       b       a       a       a       a       b       b       
+            -q filt on o          a       a       b       b       b       b       a       a
 
-    """
-    
-    # Interpret all inputs as arrays
-    X = asfarray(X)
-    h = as_column_vector(h)
-
-    r, c = X.shape
-    m = h.shape[0]
-    m2 = np.fix(m*0.5)
-
-    # Symmetrically extend with repeat of end samples.
-    # Use 'reflect' so r < m2 works OK.
-    xe = reflect(np.arange(-m2, r+m2, dtype=np.int), -0.5, r-0.5)
-
-    # Perform filtering on the columns of the extended matrix X(xe,:), keeping
-    # only the 'valid' output samples, so Y is the same size as X if m is odd.
-    Y = _column_convolve(X[xe,:], h)
-
-    return Y
-
-def coldfilt(X, ha, hb):
-    """Filter the columns of image X using the two filters ha and hb =
-    reverse(ha).  ha operates on the odd samples of X and hb on the even
-    samples.  Both filters should be even length, and h should be approx linear
-    phase with a quarter sample advance from its mid pt (i.e. :math:`|h(m/2)| >
-    |h(m/2 + 1)|`).
-
-    .. code-block:: text
-
-                          ext        top edge                     bottom edge       ext
-        Level 1:        !               |               !               |               !
-        odd filt on .    b   b   b   b   a   a   a   a   a   a   a   a   b   b   b   b   
-        odd filt on .      a   a   a   a   b   b   b   b   b   b   b   b   a   a   a   a
-        Level 2:        !               |               !               |               !
-        +q filt on x      b       b       a       a       a       a       b       b       
-        -q filt on o          a       a       b       b       b       b       a       a
-
-    The output is decimated by two from the input sample rate and the results
-    from the two filters, Ya and Yb, are interleaved to give Y.  Symmetric
-    extension with repeated end samples is used on the composite X columns
-    before each filter is applied.
-
-    Raises ValueError if the number of rows in X is not a multiple of 4, the
-    length of ha does not match hb or the lengths of ha or hb are non-even.
-
-    .. codeauthor:: Rich Wareham <rjw57 at cantab.net>, August 2013
-    .. codeauthor:: Cian Shaffrey, Cambridge University, August 2000
-    .. codeauthor:: Nick Kingsbury, Cambridge University, August 2000
+        The output is decimated by two from the input sample rate and the results
+        from the two filters, Ya and Yb, are interleaved to give Y.  Symmetric
+        extension with repeated end samples is used on the composite X columns
+        before each filter is applied.
 
-    """
-    # Make sure all inputs are arrays
-    X = asfarray(X)
-    ha = asfarray(ha)
-    hb = asfarray(hb)
-
-    r, c = X.shape
-    if r % 4 != 0:
-        raise ValueError('No. of rows in X must be a multiple of 4')
-
-    if ha.shape != hb.shape:
-        raise ValueError('Shapes of ha and hb must be the same')
-
-    if ha.shape[0] % 2 != 0:
-        raise ValueError('Lengths of ha and hb must be even')
-
-    m = ha.shape[0]
-    m2 = np.fix(m*0.5)
-
-    # Set up vector for symmetric extension of X with repeated end samples.
-    xe = reflect(np.arange(-m, r+m), -0.5, r-0.5)
-
-    # Select odd and even samples from ha and hb. Note that due to 0-indexing
-    # 'odd' and 'even' are not perhaps what you might expect them to be.
-    hao = as_column_vector(ha[0:m:2])
-    hae = as_column_vector(ha[1:m:2])
-    hbo = as_column_vector(hb[0:m:2])
-    hbe = as_column_vector(hb[1:m:2])
-    t = np.arange(5, r+2*m-2, 4)
-    r2 = r/2;
-    Y = np.zeros((r2,c), dtype=X.dtype)
-
-    if np.sum(ha*hb) > 0:
-       s1 = slice(0, r2, 2)
-       s2 = slice(1, r2, 2)
-    else:
-       s2 = slice(0, r2, 2)
-       s1 = slice(1, r2, 2)
-    
-    # Perform filtering on columns of extended matrix X(xe,:) in 4 ways. 
-    Y[s1,:] = _column_convolve(X[xe[t-1],:],hao) + _column_convolve(X[xe[t-3],:],hae)
-    Y[s2,:] = _column_convolve(X[xe[t],:],hbo) + _column_convolve(X[xe[t-2],:],hbe)
-
-    return Y
-
-def colifilt(X, ha, hb):
-    """ Filter the columns of image X using the two filters ha and hb =
-    reverse(ha).  ha operates on the odd samples of X and hb on the even
-    samples.  Both filters should be even length, and h should be approx linear
-    phase with a quarter sample advance from its mid pt (i.e `:math:`|h(m/2)| >
-    |h(m/2 + 1)|`).
-    
-    .. code-block:: text
-
-                          ext       left edge                      right edge       ext
-        Level 2:        !               |               !               |               !
-        +q filt on x      b       b       a       a       a       a       b       b       
-        -q filt on o          a       a       b       b       b       b       a       a
-        Level 1:        !               |               !               |               !
-        odd filt on .    b   b   b   b   a   a   a   a   a   a   a   a   b   b   b   b   
-        odd filt on .      a   a   a   a   b   b   b   b   b   b   b   b   a   a   a   a
-   
-    The output is interpolated by two from the input sample rate and the
-    results from the two filters, Ya and Yb, are interleaved to give Y.
-    Symmetric extension with repeated end samples is used on the composite X
-    columns before each filter is applied.
-   
-    .. codeauthor:: Rich Wareham <rjw57 at cantab.net>, August 2013
-    .. codeauthor:: Cian Shaffrey, Cambridge University, August 2000
-    .. codeauthor:: Nick Kingsbury, Cambridge University, August 2000
+        Raises ValueError if the number of rows in X is not a multiple of 4, the
+        length of ha does not match hb or the lengths of ha or hb are non-even.
 
-    """
-    # Make sure all inputs are arrays
-    X = asfarray(X)
-    ha = asfarray(ha)
-    hb = asfarray(hb)
+        .. codeauthor:: Rich Wareham <rjw57 at cantab.net>, August 2013
+        .. codeauthor:: Cian Shaffrey, Cambridge University, August 2000
+        .. codeauthor:: Nick Kingsbury, Cambridge University, August 2000
 
-    r, c = X.shape
-    if r % 2 != 0:
-        raise ValueError('No. of rows in X must be a multiple of 2')
+        """
+        # Make sure all inputs are arrays
+        X = asfarray(X)
+        ha = asfarray(ha)
+        hb = asfarray(hb)
 
-    if ha.shape != hb.shape:
-        raise ValueError('Shapes of ha and hb must be the same')
+        r, c = X.shape
+        if r % 4 != 0:
+            raise ValueError('No. of rows in X must be a multiple of 4')
 
-    if ha.shape[0] % 2 != 0:
-        raise ValueError('Lengths of ha and hb must be even')
+        if ha.shape != hb.shape:
+            raise ValueError('Shapes of ha and hb must be the same')
 
-    m = ha.shape[0]
-    m2 = np.fix(m*0.5)
+        if ha.shape[0] % 2 != 0:
+            raise ValueError('Lengths of ha and hb must be even')
 
-    Y = np.zeros((r*2,c), dtype=X.dtype)
-    if not np.any(np.nonzero(X[:])[0]):
-        return Y
+        m = ha.shape[0]
+        m2 = np.fix(m*0.5)
 
-    if m2 % 2 == 0:
-        # m/2 is even, so set up t to start on d samples.
         # Set up vector for symmetric extension of X with repeated end samples.
-        # Use 'reflect' so r < m2 works OK.
-        xe = reflect(np.arange(-m2, r+m2, dtype=np.int), -0.5, r-0.5)
-       
-        t = np.arange(3, r+m, 2)
-        if np.sum(ha*hb) > 0:
-            ta = t
-            tb = t - 1
-        else:
-            ta = t - 1
-            tb = t
-       
+        xe = reflect(np.arange(-m, r+m), -0.5, r-0.5)
+
         # Select odd and even samples from ha and hb. Note that due to 0-indexing
         # 'odd' and 'even' are not perhaps what you might expect them to be.
         hao = as_column_vector(ha[0:m:2])
         hae = as_column_vector(ha[1:m:2])
         hbo = as_column_vector(hb[0:m:2])
         hbe = as_column_vector(hb[1:m:2])
-       
-        s = np.arange(0,r*2,4)
-       
-        Y[s,:]   = _column_convolve(X[xe[tb-2],:],hae)
-        Y[s+1,:] = _column_convolve(X[xe[ta-2],:],hbe)
-        Y[s+2,:] = _column_convolve(X[xe[tb  ],:],hao)
-        Y[s+3,:] = _column_convolve(X[xe[ta  ],:],hbo)
-    else:
-        # m/2 is odd, so set up t to start on b samples.
-        # Set up vector for symmetric extension of X with repeated end samples.
-        # Use 'reflect' so r < m2 works OK.
-        xe = reflect(np.arange(-m2, r+m2, dtype=np.int), -0.5, r-0.5)
+        t = np.arange(5, r+2*m-2, 4)
+        r2 = r/2;
+        Y = np.zeros((r2,c), dtype=X.dtype)
 
-        t = np.arange(2, r+m-1, 2)
         if np.sum(ha*hb) > 0:
-            ta = t
-            tb = t - 1
+           s1 = slice(0, r2, 2)
+           s2 = slice(1, r2, 2)
         else:
-            ta = t - 1
-            tb = t
-       
-        # Select odd and even samples from ha and hb. Note that due to 0-indexing
-        # 'odd' and 'even' are not perhaps what you might expect them to be.
-        hao = as_column_vector(ha[0:m:2])
-        hae = as_column_vector(ha[1:m:2])
-        hbo = as_column_vector(hb[0:m:2])
-        hbe = as_column_vector(hb[1:m:2])
+           s2 = slice(0, r2, 2)
+           s1 = slice(1, r2, 2)
+        
+        # Perform filtering on columns of extended matrix X(xe,:) in 4 ways. 
+        Y[s1,:] = _column_convolve(X[xe[t-1],:],hao) + _column_convolve(X[xe[t-3],:],hae)
+        Y[s2,:] = _column_convolve(X[xe[t],:],hbo) + _column_convolve(X[xe[t-2],:],hbe)
+
+        return Y
+
+    def colifilt(self, X, ha, hb):
+        """ Filter the columns of image X using the two filters ha and hb =
+        reverse(ha).  ha operates on the odd samples of X and hb on the even
+        samples.  Both filters should be even length, and h should be approx linear
+        phase with a quarter sample advance from its mid pt (i.e `:math:`|h(m/2)| >
+        |h(m/2 + 1)|`).
+        
+        .. code-block:: text
+
+                              ext       left edge                      right edge       ext
+            Level 2:        !               |               !               |               !
+            +q filt on x      b       b       a       a       a       a       b       b       
+            -q filt on o          a       a       b       b       b       b       a       a
+            Level 1:        !               |               !               |               !
+            odd filt on .    b   b   b   b   a   a   a   a   a   a   a   a   b   b   b   b   
+            odd filt on .      a   a   a   a   b   b   b   b   b   b   b   b   a   a   a   a
        
-        s = np.arange(0,r*2,4)
+        The output is interpolated by two from the input sample rate and the
+        results from the two filters, Ya and Yb, are interleaved to give Y.
+        Symmetric extension with repeated end samples is used on the composite X
+        columns before each filter is applied.
        
-        Y[s,:]   = _column_convolve(X[xe[tb],:],hao)
-        Y[s+1,:] = _column_convolve(X[xe[ta],:],hbo)
-        Y[s+2,:] = _column_convolve(X[xe[tb],:],hae)
-        Y[s+3,:] = _column_convolve(X[xe[ta],:],hbe)
+        .. codeauthor:: Rich Wareham <rjw57 at cantab.net>, August 2013
+        .. codeauthor:: Cian Shaffrey, Cambridge University, August 2000
+        .. codeauthor:: Nick Kingsbury, Cambridge University, August 2000
+
+        """
+        # Make sure all inputs are arrays
+        X = asfarray(X)
+        ha = asfarray(ha)
+        hb = asfarray(hb)
+
+        r, c = X.shape
+        if r % 2 != 0:
+            raise ValueError('No. of rows in X must be a multiple of 2')
+
+        if ha.shape != hb.shape:
+            raise ValueError('Shapes of ha and hb must be the same')
+
+        if ha.shape[0] % 2 != 0:
+            raise ValueError('Lengths of ha and hb must be even')
+
+        m = ha.shape[0]
+        m2 = np.fix(m*0.5)
+
+        Y = np.zeros((r*2,c), dtype=X.dtype)
+        if not np.any(np.nonzero(X[:])[0]):
+            return Y
+
+        if m2 % 2 == 0:
+            # m/2 is even, so set up t to start on d samples.
+            # Set up vector for symmetric extension of X with repeated end samples.
+            # Use 'reflect' so r < m2 works OK.
+            xe = reflect(np.arange(-m2, r+m2, dtype=np.int), -0.5, r-0.5)
+           
+            t = np.arange(3, r+m, 2)
+            if np.sum(ha*hb) > 0:
+                ta = t
+                tb = t - 1
+            else:
+                ta = t - 1
+                tb = t
+           
+            # Select odd and even samples from ha and hb. Note that due to 0-indexing
+            # 'odd' and 'even' are not perhaps what you might expect them to be.
+            hao = as_column_vector(ha[0:m:2])
+            hae = as_column_vector(ha[1:m:2])
+            hbo = as_column_vector(hb[0:m:2])
+            hbe = as_column_vector(hb[1:m:2])
+           
+            s = np.arange(0,r*2,4)
+           
+            Y[s,:]   = _column_convolve(X[xe[tb-2],:],hae)
+            Y[s+1,:] = _column_convolve(X[xe[ta-2],:],hbe)
+            Y[s+2,:] = _column_convolve(X[xe[tb  ],:],hao)
+            Y[s+3,:] = _column_convolve(X[xe[ta  ],:],hbo)
+        else:
+            # m/2 is odd, so set up t to start on b samples.
+            # Set up vector for symmetric extension of X with repeated end samples.
+            # Use 'reflect' so r < m2 works OK.
+            xe = reflect(np.arange(-m2, r+m2, dtype=np.int), -0.5, r-0.5)
+
+            t = np.arange(2, r+m-1, 2)
+            if np.sum(ha*hb) > 0:
+                ta = t
+                tb = t - 1
+            else:
+                ta = t - 1
+                tb = t
+           
+            # Select odd and even samples from ha and hb. Note that due to 0-indexing
+            # 'odd' and 'even' are not perhaps what you might expect them to be.
+            hao = as_column_vector(ha[0:m:2])
+            hae = as_column_vector(ha[1:m:2])
+            hbo = as_column_vector(hb[0:m:2])
+            hbe = as_column_vector(hb[1:m:2])
+           
+            s = np.arange(0,r*2,4)
+           
+            Y[s,:]   = _column_convolve(X[xe[tb],:],hao)
+            Y[s+1,:] = _column_convolve(X[xe[ta],:],hbo)
+            Y[s+2,:] = _column_convolve(X[xe[tb],:],hae)
+            Y[s+3,:] = _column_convolve(X[xe[ta],:],hbe)
 
-    return Y
+        return Y
 
 # vim:sw=4:sts=4:et
 
diff --git a/dtcwt/lowlevel.py b/dtcwt/lowlevel.py
index c440bf9..77307e5 100644
--- a/dtcwt/lowlevel.py
+++ b/dtcwt/lowlevel.py
@@ -1 +1,74 @@
-from dtcwt.backend.numpy.lowlevel import *
+from dtcwt.backend.numpy.lowlevel import LowLevelBackendNumPy
+
+_BACKEND = LowLevelBackendNumPy()
+
+def colfilter(X, h):
+    """Filter the columns of image *X* using filter vector *h*, without decimation.
+    If len(h) is odd, each output sample is aligned with each input sample
+    and *Y* is the same size as *X*.  If len(h) is even, each output sample is
+    aligned with the mid point of each pair of input samples, and Y.shape =
+    X.shape + [1 0].
+
+    :param X: an image whose columns are to be filtered
+    :param h: the filter coefficients.
+    :returns Y: the filtered image.
+
+    """
+    return _BACKEND.colfilter(X, h)
+
+def coldfilt(X, ha, hb):
+    """Filter the columns of image X using the two filters ha and hb =
+    reverse(ha).  ha operates on the odd samples of X and hb on the even
+    samples.  Both filters should be even length, and h should be approx linear
+    phase with a quarter sample advance from its mid pt (i.e. :math:`|h(m/2)| >
+    |h(m/2 + 1)|`).
+
+    .. code-block:: text
+
+                          ext        top edge                     bottom edge       ext
+        Level 1:        !               |               !               |               !
+        odd filt on .    b   b   b   b   a   a   a   a   a   a   a   a   b   b   b   b   
+        odd filt on .      a   a   a   a   b   b   b   b   b   b   b   b   a   a   a   a
+        Level 2:        !               |               !               |               !
+        +q filt on x      b       b       a       a       a       a       b       b       
+        -q filt on o          a       a       b       b       b       b       a       a
+
+    The output is decimated by two from the input sample rate and the results
+    from the two filters, Ya and Yb, are interleaved to give Y.  Symmetric
+    extension with repeated end samples is used on the composite X columns
+    before each filter is applied.
+
+    Raises ValueError if the number of rows in X is not a multiple of 4, the
+    length of ha does not match hb or the lengths of ha or hb are non-even.
+
+    """
+    return _BACKEND.coldfilt(X, ha, hb)
+
+def colifilt(X, ha, hb):
+    """ Filter the columns of image X using the two filters ha and hb =
+    reverse(ha).  ha operates on the odd samples of X and hb on the even
+    samples.  Both filters should be even length, and h should be approx linear
+    phase with a quarter sample advance from its mid pt (i.e `:math:`|h(m/2)| >
+    |h(m/2 + 1)|`).
+    
+    .. code-block:: text
+
+                          ext       left edge                      right edge       ext
+        Level 2:        !               |               !               |               !
+        +q filt on x      b       b       a       a       a       a       b       b       
+        -q filt on o          a       a       b       b       b       b       a       a
+        Level 1:        !               |               !               |               !
+        odd filt on .    b   b   b   b   a   a   a   a   a   a   a   a   b   b   b   b   
+        odd filt on .      a   a   a   a   b   b   b   b   b   b   b   b   a   a   a   a
+   
+    The output is interpolated by two from the input sample rate and the
+    results from the two filters, Ya and Yb, are interleaved to give Y.
+    Symmetric extension with repeated end samples is used on the composite X
+    columns before each filter is applied.
+   
+    .. codeauthor:: Rich Wareham <rjw57 at cantab.net>, August 2013
+    .. codeauthor:: Cian Shaffrey, Cambridge University, August 2000
+    .. codeauthor:: Nick Kingsbury, Cambridge University, August 2000
+
+    """
+    return _BACKEND.colifilt(X, ha, hb)
diff --git a/dtcwt/opencl/lowlevel.py b/dtcwt/opencl/lowlevel.py
index 8de1501..3a350ad 100644
--- a/dtcwt/opencl/lowlevel.py
+++ b/dtcwt/opencl/lowlevel.py
@@ -17,7 +17,7 @@ from six.moves import xrange
 import struct
 import functools
 
-from dtcwt.lowlevel import asfarray, as_column_vector, reflect
+from dtcwt.utils import asfarray, as_column_vector
 
 # note that this decorator ignores **kwargs
 # From https://wiki.python.org/moin/PythonDecoratorLibrary#Alternate_memoize_as_nested_functions
diff --git a/dtcwt/opencl/transform2d.py b/dtcwt/opencl/transform2d.py
index 32feb51..81f665a 100644
--- a/dtcwt/opencl/transform2d.py
+++ b/dtcwt/opencl/transform2d.py
@@ -6,7 +6,7 @@ from six.moves import xrange
 
 from dtcwt import biort as _biort, qshift as _qshift
 from dtcwt.defaults import DEFAULT_BIORT, DEFAULT_QSHIFT
-from dtcwt.lowlevel import appropriate_complex_type_for, asfarray
+from dtcwt.utils import appropriate_complex_type_for, asfarray
 from dtcwt.opencl.lowlevel import colfilter, coldfilt, colifilt
 from dtcwt.opencl.lowlevel import axis_convolve, axis_convolve_dfilter, q2c
 from dtcwt.opencl.lowlevel import to_device, to_queue, to_array, empty
diff --git a/dtcwt/sampling.py b/dtcwt/sampling.py
index 16fd8f8..8816a94 100644
--- a/dtcwt/sampling.py
+++ b/dtcwt/sampling.py
@@ -8,7 +8,7 @@ __all__ = (
     'upsample', 'upsample_highpass',
 )
 
-from dtcwt.lowlevel import reflect, asfarray
+from dtcwt.utils import reflect, asfarray
 
 import numpy as np
 
diff --git a/dtcwt/transform1d.py b/dtcwt/transform1d.py
index b6c54e0..2b34ef0 100644
--- a/dtcwt/transform1d.py
+++ b/dtcwt/transform1d.py
@@ -5,7 +5,8 @@ from six.moves import xrange
 
 from dtcwt import biort as _biort, qshift as _qshift
 from dtcwt.defaults import DEFAULT_BIORT, DEFAULT_QSHIFT
-from dtcwt.lowlevel import colfilter, coldfilt, colifilt, as_column_vector, asfarray
+from dtcwt.lowlevel import colfilter, coldfilt, colifilt
+from dtcwt.utils import as_column_vector, asfarray
 
 def dtwavexfm(X, nlevels=3, biort=DEFAULT_BIORT, qshift=DEFAULT_QSHIFT, include_scale=False):
     """Perform a *n*-level DTCWT decompostion on a 1D column vector *X* (or on
diff --git a/dtcwt/transform2d.py b/dtcwt/transform2d.py
index 05ee079..3915098 100644
--- a/dtcwt/transform2d.py
+++ b/dtcwt/transform2d.py
@@ -5,7 +5,8 @@ from six.moves import xrange
 
 from dtcwt import biort as _biort, qshift as _qshift
 from dtcwt.defaults import DEFAULT_BIORT, DEFAULT_QSHIFT
-from dtcwt.lowlevel import colfilter, coldfilt, colifilt, appropriate_complex_type_for, asfarray
+from dtcwt.lowlevel import colfilter, coldfilt, colifilt
+from dtcwt.utils import appropriate_complex_type_for, asfarray
 
 def dtwavexfm2(X, nlevels=3, biort=DEFAULT_BIORT, qshift=DEFAULT_QSHIFT, include_scale=False):
     """Perform a *n*-level DTCWT-2D decompostion on a 2D matrix *X*.
diff --git a/dtcwt/transform3d.py b/dtcwt/transform3d.py
index 7c0c0fc..d013638 100644
--- a/dtcwt/transform3d.py
+++ b/dtcwt/transform3d.py
@@ -5,7 +5,8 @@ from six.moves import xrange
 
 from dtcwt import biort as _biort, qshift as _qshift
 from dtcwt.defaults import DEFAULT_BIORT, DEFAULT_QSHIFT
-from dtcwt.lowlevel import colfilter, coldfilt, colifilt, asfarray
+from dtcwt.lowlevel import colfilter, coldfilt, colifilt
+from dtcwt.utils import asfarray
 
 def dtwavexfm3(X, nlevels=3, biort=DEFAULT_BIORT, qshift=DEFAULT_QSHIFT, ext_mode=4, discard_level_1=False):
     """Perform a *n*-level DTCWT-3D decompostion on a 3D matrix *X*.
diff --git a/dtcwt/utils.py b/dtcwt/utils.py
index 8b716f1..185f40c 100644
--- a/dtcwt/utils.py
+++ b/dtcwt/utils.py
@@ -4,7 +4,6 @@ __all__ = ( 'drawedge', 'drawcirc',)
 
 import numpy as np
 
-
 def drawedge(theta,r,w,N):
     """Generate an image of size N * N pels, of an edge going from 0 to 1
     in height at theta degrees to the horizontal (top of image = 1 if angle = 0).
@@ -57,3 +56,71 @@ def drawcirc(r,w,du,dv,N):
     # Final circle image plane
     p = 0.5 + 0.5 * np.sin(np.minimum(np.maximum((np.exp(np.array([-0.5]) * (x**2 + y**2)).T - np.exp((-0.5))) * (r * 3 / w), np.pi/(-2)), np.pi/2))
     return p
+
+def asfarray(X):
+    """Similar to :py:func:`numpy.asfarray` except that this function tries to
+    preserve the original datatype of X if it is already a floating point type
+    and will pass floating point arrays through directly without copying.
+
+    """
+    X = np.asanyarray(X)
+    return np.asfarray(X, dtype=X.dtype)
+
+def appropriate_complex_type_for(X):
+    """Return an appropriate complex data type depending on the type of X. If X
+    is already complex, return that, if it is floating point return a complex
+    type of the appropriate size and if it is integer, choose an complex
+    floating point type depending on the result of :py:func:`numpy.asfarray`.
+
+    """
+    X = asfarray(X)
+    
+    if np.issubsctype(X.dtype, np.complex64) or np.issubsctype(X.dtype, np.complex128):
+        return X.dtype
+    elif np.issubsctype(X.dtype, np.float32):
+        return np.complex64
+    elif np.issubsctype(X.dtype, np.float64):
+        return np.complex128
+
+    # God knows, err on the side of caution
+    return np.complex128
+
+def as_column_vector(v):
+    """Return *v* as a column vector with shape (N,1).
+    
+    """
+    v = np.atleast_2d(v)
+    if v.shape[0] == 1:
+        return v.T
+    else:
+        return v
+
+def reflect(x, minx, maxx):
+    """Reflect the values in matrix *x* about the scalar values *minx* and
+    *maxx*.  Hence a vector *x* containing a long linearly increasing series is
+    converted into a waveform which ramps linearly up and down between *minx* and
+    *maxx*.  If *x* contains integers and *minx* and *maxx* are (integers + 0.5), the
+    ramps will have repeated max and min samples.
+   
+    .. codeauthor:: Rich Wareham <rjw57 at cantab.net>, Aug 2013
+    .. codeauthor:: Nick Kingsbury, Cambridge University, January 1999.
+    
+    """
+
+    # Copy x to avoid in-place modification
+    y = np.array(x, copy=True)
+
+    # Reflect y in maxx.
+    t = y > maxx
+    y[t] = (2*maxx - y[t]).astype(y.dtype)
+
+    while np.any(y < minx):
+        # Reflect y in minx.
+        t = y < minx
+        y[t] = (2*minx - y[t]).astype(y.dtype)
+
+        # Reflect y in maxx.
+        t = y > maxx
+        y[t] = (2*maxx - y[t]).astype(y.dtype)
+
+    return y
diff --git a/tests/testreflect.py b/tests/testreflect.py
index 3e55258..1467042 100644
--- a/tests/testreflect.py
+++ b/tests/testreflect.py
@@ -1,6 +1,6 @@
 import numpy as np
 
-from dtcwt.lowlevel import reflect
+from dtcwt.utils import reflect
 
 def setup():
     global ramp, reflected

-- 
Alioth's /usr/local/bin/git-commit-notice on /srv/git.debian.org/git/debian-science/packages/python-dtcwt.git



More information about the debian-science-commits mailing list