Discrete wavelet transform
From Free net encyclopedia
Template:Unreferenced In numerical analysis and functional analysis, the discrete wavelet transform (DWT) refers to wavelet transforms for which the wavelets are discretely sampled.
The first DWT was invented by Alfréd Haar, a Hungarian mathematician. For an input represented by a list of <math>2^n</math> numbers, the Haar wavelet transform may be considered to simply pair up input values, storing the difference and passing the sum. This process is repeated recursively, pairing up the sums to provide the next scale: finally resulting in <math>2^n-1</math> differences and one final sum.
This simple DWT illustrates the desirable properties of wavelets in general. Firstly, the discrete transform can be performed in <math>O(n)</math> operations. Secondly, the transform captures not only some notion of the frequency content of the input, by examining it at different scales, but also captures temporal content, i.e. the times at which these frequencies occur. Combined, these two properties make the Fast Wavelet Transform (FWT), an alternative to the conventional Fast Fourier Transform.
The most common set of discrete wavelet transforms were formulated by the Belgian mathematician Ingrid Daubechies in 1988. This formulation is based upon the use of recurrence relations to generate progressively finer discrete samplings of an implicit mother wavelet function, each resolution being twice that of the previous scale. In her seminal paper, Daubechies derives a family of wavelets, the first of which is the Haar wavelet. Interest in this field has exploded since then, with the development of many descendants of Daubechies' original family of wavelets.
Other forms of discrete wavelet transform include the non- or undecimated wavelet transform (where downsampling is omitted), the Newland transform (where an orthonormal basis of wavelets is formed from appropriately constructed top-hat filters in frequency space). Wavelet packet transforms are also related to the discrete wavelet transform. Complex wavelet transform is another form.
The discrete wavelet transform has a huge number of applications in science, engineering, mathematics and computer science. Most notably, the discrete wavelet transform is used for signal coding, where the properties of the transform are exploited to represent a discrete signal in a more redundant form, often as a preconditioning for data compression.
Contents |
Definition
One level of the transform
The DWT of a signal <math>x</math> is calculated by passing it through a series of filters. First the samples are passed through a low pass filter with impulse response <math>g</math> resulting in a convolution of the two:
- <math>y[n] = (x * g)[n] = \sum\limits_{k = - \infty }^\infty {x[k] \cdot g[n - k]} </math>
The signal is also decomposed simultaneously using a high-pass filter <math>h</math>. The outputs giving the detail coefficients (from the high-pass filter) and approximation coefficients (from the low-pass). It is important that the two filters are related to each other and they are known as a quadrature mirror filter.
However, since half the frequencies of the signal have now been removed, half the samples can be discarding according to Nyquist’s rule. The filter outputs are then downsampled by 2:
- <math>y_{\mathrm{low}} [n] = \sum\limits_{k = - \infty }^\infty {x[k] \cdot g[2\cdot n - k]} </math>
- <math>y_{\mathrm{high}} [n] = \sum\limits_{k = - \infty }^\infty {x[k] \cdot h[2\cdot n - k]} </math>
This decomposition has halved the time resolution since only half of each filter output characterises the signal. However, each output has half the frequency band of the input so the frequency resolution has been doubled. This is in keeping with the Heisenberg Uncertainty Principle.
With the downsampling operator <math>\downarrow</math>
- <math>(y \downarrow k)[n] = y[k\cdot n] </math>
the above summation can be written more concisely.
- <math>y_{\mathrm{low}} = (x*g)\downarrow 2 </math>
- <math>y_{\mathrm{high}} = (x*h)\downarrow 2 </math>
However computing a complete convolution <math>x*g</math> with subsequent downsampling would waste computation time.
The Lifting scheme is an optimization where these two computations are interleaved.
Cascading and Filter banks
This decomposition is repeated to further increase the frequency resolution and the approximation coefficients decomposed with high and low pass filters and then down-sampled. This is represented as a binary tree with nodes representing a sub-space with a different time-frequency localisation. The tree is known as a filter bank.
Image:Wavelets - Filter Bank.png
At each level in the above diagram the signal is decomposed into low and high frequencies. Due to the decomposition process the input signal must be a multiple of <math>2^n</math> where <math>n</math> is the number of levels.
For example a signal with 16 samples, frequency range 0 to <math>f_n</math> and 3 levels of decomposition, 4 output scales are produced:
| Level | Frequencies | Samples |
|---|---|---|
| 3 | <math>0</math> to <math>Template:F n/8</math> | 4 |
| <math>Template:F n/8</math> to <math>Template:F n/4</math> | 4 | |
| 2 | <math>Template:F n/4</math> to <math>Template:F n/2</math> | 8 |
| 1 | <math>Template:F n/2</math> to <math>f_n</math> | 16 |
Code examples
In its simplest form, the DWT is remarkably easy to compute. The following snippet of the Objective Caml code computes the 1D Haar wavelet transform of an integer-power-of-two-length list of numbers:
# let haar l =
let rec aux l s d = match l, s, d with
[s], [], d -> s :: d
| [], s, d -> aux s [] d
| h1::h2::t, s, d -> aux t (h1 +. h2 :: s) (h1 -. h2 :: d)
| _ -> invalid_arg "haar" in
aux l [] [];;
val haar : float list -> float list = <fun>
For example:
# haar [1.; 2.; 3.; 4.; -4.; -3.; -2.; -1.];; - : float list = [0.; 20.; 4.; 4.; -1.; -1.; -1.; -1.]
The code in an Imperative programming language like Java is more verbose:
public static int[] invoke(int[] input){
//THIS WILL DESTROY THE CONTENTS OF THE input ARRAY
//WE NEGLECT TO VERIFY input.length=2^n, n>1
int[] output=new int[input.length];
for(int length=input.length>>1;;length>>=1){
//length=2^n, WITH DECREASING n
for(int i=0;i<length;i++){
int sum=input[i*2]+input[i*2+1];
int dif=input[i*2]-input[i*2+1];
output[i]=sum;
output[i+length]=dif;
}//for
if (length==1) return output;
//SWAP ARRAYS TO DO NEXT INTERATION
System.arraycopy(output, 0, input, 0, length<<1);
}//for
}//method
WARNING!! These two methods give different answers. The Java program gives {0, 20, -4, -4, -1, -1, -1, -1} with input of {1, 2, 3, 4, -4, -3, -2, -1}. I suspect the Objective Caml program is swapping array elements, making every second pass result in negative differences. (disputed — see [[: talk:Discrete wavelet transform#{{{1|Disputed}}}|talk page]])
See also
fr:Transformée en ondelette discrète ru:Discrete wavelet transform