Discrete Wavelet Transform¶
Wraps GSL's gsl_wavelet_transform_forward()
,
and gsl_wavelet_transform_inverse()
and creates plot methods.
AUTHOR:
Josh Kantor (2006-10-07) - initial version
David Joyner (2006-10-09) - minor changes to docstrings and examples.
-
sage.calculus.transforms.dwt.
DWT
(n, wavelet_type, wavelet_k)¶ This function initializes an GSLDoubleArray of length n which can perform a discrete wavelet transform.
INPUT:
n
– a power of 2T
– the data in the GSLDoubleArray must be realwavelet_type
– the name of the type of wavelet, valid choices are:'daubechies'
'daubechies_centered'
'haar'
'haar_centered'
'bspline'
'bspline_centered'
For daubechies wavelets,
wavelet_k
specifies a daubechie wavelet with \(k/2\) vanishing moments. \(k = 4,6,...,20\) for \(k\) even are the only ones implemented.For Haar wavelets,
wavelet_k
must be 2.For bspline wavelets,
wavelet_k
of \(103,105,202,204,206,208,301, 305,307,309\) will give biorthogonal B-spline wavelets of order \((i,j)\) wherewavelet_k
is \(100*i+j\). The wavelet transform uses \(J = \log_2(n)\) levels.OUTPUT:
An array of the form \((s_{-1,0}, d_{0,0}, d_{1,0}, d_{1,1}, d_{2,0}, \ldots, d_{J-1,2^{J-1}-1})\) for \(d_{j,k}\) the detail coefficients of level \(j\). The centered forms align the coefficients of the sub-bands on edges.
EXAMPLES:
sage: a = WaveletTransform(128,'daubechies',4) sage: for i in range(1, 11): ....: a[i] = 1 ....: a[128-i] = 1 sage: a.plot().show(ymin=0) sage: a.forward_transform() sage: a.plot().show() sage: a = WaveletTransform(128,'haar',2) sage: for i in range(1, 11): a[i] = 1; a[128-i] = 1 sage: a.forward_transform() sage: a.plot().show(ymin=0) sage: a = WaveletTransform(128,'bspline_centered',103) sage: for i in range(1, 11): a[i] = 1; a[100+i] = 1 sage: a.forward_transform() sage: a.plot().show(ymin=0)
This example gives a simple example of wavelet compression:
sage: a = DWT(2048,'daubechies',6) sage: for i in range(2048): a[i]=float(sin((i*5/2048)**2)) sage: a.plot().show() # long time (7s on sage.math, 2011) sage: a.forward_transform() sage: for i in range(1800): a[2048-i-1] = 0 sage: a.backward_transform() sage: a.plot().show() # long time (7s on sage.math, 2011)
-
class
sage.calculus.transforms.dwt.
DiscreteWaveletTransform
¶ Bases:
sage.libs.gsl.array.GSLDoubleArray
Discrete wavelet transform class.
-
backward_transform
()¶
-
forward_transform
()¶
-
plot
(xmin=None, xmax=None, **args)¶
-
-
sage.calculus.transforms.dwt.
WaveletTransform
(n, wavelet_type, wavelet_k)¶ This function initializes an GSLDoubleArray of length n which can perform a discrete wavelet transform.
INPUT:
n
– a power of 2T
– the data in the GSLDoubleArray must be realwavelet_type
– the name of the type of wavelet, valid choices are:'daubechies'
'daubechies_centered'
'haar'
'haar_centered'
'bspline'
'bspline_centered'
For daubechies wavelets,
wavelet_k
specifies a daubechie wavelet with \(k/2\) vanishing moments. \(k = 4,6,...,20\) for \(k\) even are the only ones implemented.For Haar wavelets,
wavelet_k
must be 2.For bspline wavelets,
wavelet_k
of \(103,105,202,204,206,208,301, 305,307,309\) will give biorthogonal B-spline wavelets of order \((i,j)\) wherewavelet_k
is \(100*i+j\). The wavelet transform uses \(J = \log_2(n)\) levels.OUTPUT:
An array of the form \((s_{-1,0}, d_{0,0}, d_{1,0}, d_{1,1}, d_{2,0}, \ldots, d_{J-1,2^{J-1}-1})\) for \(d_{j,k}\) the detail coefficients of level \(j\). The centered forms align the coefficients of the sub-bands on edges.
EXAMPLES:
sage: a = WaveletTransform(128,'daubechies',4) sage: for i in range(1, 11): ....: a[i] = 1 ....: a[128-i] = 1 sage: a.plot().show(ymin=0) sage: a.forward_transform() sage: a.plot().show() sage: a = WaveletTransform(128,'haar',2) sage: for i in range(1, 11): a[i] = 1; a[128-i] = 1 sage: a.forward_transform() sage: a.plot().show(ymin=0) sage: a = WaveletTransform(128,'bspline_centered',103) sage: for i in range(1, 11): a[i] = 1; a[100+i] = 1 sage: a.forward_transform() sage: a.plot().show(ymin=0)
This example gives a simple example of wavelet compression:
sage: a = DWT(2048,'daubechies',6) sage: for i in range(2048): a[i]=float(sin((i*5/2048)**2)) sage: a.plot().show() # long time (7s on sage.math, 2011) sage: a.forward_transform() sage: for i in range(1800): a[2048-i-1] = 0 sage: a.backward_transform() sage: a.plot().show() # long time (7s on sage.math, 2011)
-
sage.calculus.transforms.dwt.
is2pow
(n)¶