MATLABsignal processingwavelets
Part of my assignment about signal processing says the following:
- Compute the Discrete Wavelet Transform for the input signals
- Group the wavelet coefficients in trees growing across scales
- Classify the trees in high- and low-energy trees
I've managed to implement my own DWT, calculate the coefficients and group them in trees.
Now I'm left with calculating the energy of the wavelet coefficients.
Question: How to calculate the energy of wavelet coefficients?
I am doing this in Matlab.
Well, since the wavelet transform is not precisely a time-frequency analysis but more a time-scale analysis one can not really display the wavelet coefficients in the time frequency plane. (However, even for the STFT the time-frequency plane is only a crutch for illustration purposes.)
Maybe it helps to think of $N=2^n$. Then you can analyze your signal into $n$ scales with the wavelet transform. See, e.g., this image with $n=3$:
Then the wavelet-transform stores the $c_{31}$ and all the $d_{ij}$'s. If you want to put the coefficients in the time-frequency plane you would have to stretch the tiles as you indicated in your figures. There would be $n$ vertical colums (in $\omega$), the rightmost having $2^{n/2}$ boxes and the number of boxes halving in each column (the last two containing just one box). This adds up to $2^n$.
However, I suggest to read Mallat's "A wavelet tour of signal processing - the sparse way" Chapters 1.3, 4.3 and 8.1. It's quite hard to produce illustrating figures for the wavelet transform here...
Considering the efficiency: To be a bit more precise, the discrete wavelet transform (DWT) of a signal of length $2^n$ and a filter with length $k$ takes $O(k2^n)$, while FFT is $O(n2^n)$ with comparable constants. Hence, DWT is faster it the filter is not too long.
The problem lies with the convolution. Standard Matlab conv function will use full convolution. So the length will always be dependent on the longest of the 2 arguments. This can be resolved by using the same parameter for convolution.
'same': Central part of the convolution of the same size as u.
So the correct code is:
function R=myDWT(sig, count)
[Lo_D, Hi_D] = wfilter_bior44();
input = sig;
while(count ~= 0) % While count not equal to 0
% Pass through filters by using convolution
Ca = conv(input, Lo_D, 'same');
Cd = conv(input, Hi_D, 'same');
% Downsample by 2
Ca = downsample(Ca, 2);
Cd = downsample(Cd, 2);
% TODO: Save Ca and Cd somewhere
count = count - 1;
input = Ca;
end
R = input;
end
Best Answer
If such coefficients are calculated from an orthonormal basis, rather than from a redundant wavelet transform like the discrete case of a CWT, then the energy is the sum of the square of the absoute values of those coefficients, because orthonormal wavelet bases converge in norm, and the energy of a signal is defined as the integral of its absolute value. See https://en.wikipedia.org/wiki/Energy_(signal_processing)