[Scipy-svn] r3227 - in trunk/Lib/signal: . tests

scipy-svn@scip... scipy-svn@scip...
Fri Aug 10 19:57:02 CDT 2007


Author: stefan
Date: 2007-08-10 19:56:45 -0500 (Fri, 10 Aug 2007)
New Revision: 3227

Modified:
   trunk/Lib/signal/tests/test_wavelets.py
   trunk/Lib/signal/wavelets.py
Log:
Add complex Morlet wavelet (based on contribution by iCy-fLaME).


Modified: trunk/Lib/signal/tests/test_wavelets.py
===================================================================
--- trunk/Lib/signal/tests/test_wavelets.py	2007-08-11 00:08:07 UTC (rev 3226)
+++ trunk/Lib/signal/tests/test_wavelets.py	2007-08-11 00:56:45 UTC (rev 3227)
@@ -22,5 +22,14 @@
                 assert len(x) == len(phi) == len(psi)
                 assert_equal(len(x),(k-1)*2**J)
 
+    def check_morlet(self):
+        x = wavelets.morlet(50,4.1,complete=True)
+        y = wavelets.morlet(50,4.1,complete=False)
+        assert_equal(len(x),len(y))
+
+        x = wavelets.morlet(10,50,complete=False)
+        y = wavelets.morlet(10,50,complete=True)
+        assert_equal(x,y)
+
 if __name__ == "__main__":
     NumpyTest().run()

Modified: trunk/Lib/signal/wavelets.py
===================================================================
--- trunk/Lib/signal/wavelets.py	2007-08-11 00:08:07 UTC (rev 3226)
+++ trunk/Lib/signal/wavelets.py	2007-08-11 00:56:45 UTC (rev 3227)
@@ -1,8 +1,9 @@
-__all__ = ['daub','qmf','cascade']
+__all__ = ['daub','qmf','cascade','morlet']
 
 import numpy as sb
 from numpy.dual import eig
 from scipy.misc import comb
+from scipy import linspace, pi, exp, zeros
 
 def daub(p):
     """The coefficients for the FIR low-pass filter producing Daubechies wavelets.
@@ -167,3 +168,50 @@
         prevkeys = newkeys
 
     return x, phi, psi
+
+def morlet(M, w=5.0, s=1.0, complete=True):
+    """Complex Morlet wavelet.
+
+    :Parameters:
+        M : int
+            Length of the wavelet.
+        w : float
+            Omega0
+        s : float
+            Scaling factor, windowed from -s*2*pi to +s*2*pi.
+        complete : bool
+            Whether to use the complete or the standard version.
+
+    Notes:
+    ------
+
+    The standard version:
+        pi**-0.25 * exp(1j*w*x) * exp(-0.5*(x**2))
+
+        This commonly used wavelet is often referred to simply as the
+        Morlet wavelet.  Note that, this simplified version can cause
+        admissibility problems at low values of w.
+
+    The complete version:
+        pi**-0.25 * (exp(1j*w*x) - exp(-0.5*(w**2))) * exp(-0.5*(x**2))
+
+        The complete version of the Morlet wavelet, with a correction
+        term to improve admissibility. For w greater than 5, the
+        correction term is negligible.
+
+    Note that the energy of the return wavelet is not normalised
+    according to s.
+
+    The fundamental frequency of this wavelet in Hz is given
+    by f = 2*s*w*r / M where r is the sampling rate.
+
+    """
+    x = linspace(-s*2*pi,s*2*pi,M)
+    output = exp(1j*w*x)
+
+    if complete:
+        x -= exp(-0.5*(w**2))
+
+    output *= exp(-0.5*(x**2)) * pi**(-0.25)
+
+    return output



More information about the Scipy-svn mailing list