Thursday, April 11, 2019

discrete signals - Daubechies wavelet transform


i have N samples obtained by sampling a signal with lot of frequency contents. How will i apply daubechies wavelet transform to obtain the frequency and its location? i need to write a program which will process the signal and gives the frequency and location as the result.



Answer



Looks like you need a general explanation of the discrete wavelet transform (DWT). DWT breaks a signal down into subbands distributed evenly in a logarithmic frequency scale, each subband sampled at a rate proportional to the frequencies in that band. The traditional Fourier transformation has no time domain resolution at all, or when done using many short windows on a longer data, equal resolution at all frequencies. The distribution of samples in the time and frequency domain by DWT is of form:


log f
|XXXXXXXXXXXXXXXX X = a sample
|X X X X X X X X f = frequency
|X X X X t = time
|X X
|X

|X
----------------t

Single subband decomposition and reconstruction:


    -> high -> decimate -------------> dilute -> high
| pass by 2 high subband by 2 pass \
in | + out
| / =in
-> low -> decimate -------------> dilute -> low
pass by 2 low subband by 2 pass


This creates two subbands from the input signal, both sampled at half the original frequency. The filters approximate halfband finite impulse response (FIR) filters and are determined by the choice of wavelet. Using Daubechies wavelets (and most others), the data can be reconstructed to the exact original even when the halfband filters are not perfect. Note that in the above scheme, the total amount of information (samples) stays the same throughout.


Decimation by 2: ABCDEFGHIJKLMNOPQR -> ACEGIKMOQ
Dilution by 2: ACEGIKMOQ -> A0C0E0G0I0K0M0O0Q0

To get the logarithmic resolution in frequency, the low subband is re-transformed, and again, the low subband from this transformation gets the same treatment etc.


Decomposition:


    -> high -> decimate --------------------------------> subband0
| pass by 2
in | -> high -> decimate ---------------> subband1

| | pass by 2
-> low -> decim | -> high -> decim -> subband2
pass by 2 | | pass by 2
-> low -> decim |
pass by 2 | . down to what suffices
-> . or if periodic data,
. until short of data

Reconstruction:


subband0 -----------------------------------> dilute -> high

by 2 pass \
subband1 ------------------> dilute -> high + out
by 2 pass \ / =in
subband2 -> dilute -> high + dilute -> low
by 2 pass \ / by 2 pass
+ dilute -> low
Start . / by 2 pass
here! . -> dilute -> low
. by 2 pass


In a real-time application, the filters introduce delays, so you need to compensate them by adding additional delays to less-delayed higher bands, to get the summation work as intended.


For periodic signals or windowed operation, this problem doesn't exist - a single subband transformation is a matrix multiplication, with wrapping implemented in the matrix:


Decomposition:


|L0|   |C0  C1  C2  C3                | |I0|     L = lowpass output
|H0| |C3 -C2 C1 -C0 | |I1| H = highpass output
|L1| | C0 C1 C2 C3 | |I2| I = input
|H1| = | C3 -C2 C1 -C0 | |I3| C = coefficients
|L2| | C0 C1 C2 C3| |I4|
|H2| | C3 -C2 C1 -C0| |I5|
|L3| |C2 C3 C0 C1| |I6|

|H3| |C1 -C0 C3 -C2| |I7| Daubechies 4-coef:

1+sqrt(3) 3+sqrt(3) 3-sqrt(3) 1-sqrt(3)
C0 = --------- C1 = --------- C2 = --------- C3 = ---------
4 sqrt(2) 4 sqrt(2) 4 sqrt(2) 4 sqrt(2)

Reconstruction:


|I0|   |C0  C3                  C2  C1| |L0|
|I1| |C1 -C2 C3 -C0| |H0|
|I2| |C2 C1 C0 C3 | |L1|

|I3| = |C3 -C0 C1 -C2 | |H1|
|I4| | C2 C1 C0 C3 | |L2|
|I5| | C3 -C0 C1 -C2 | |H2|
|I6| | C2 C1 C0 C3| |L3|
|I7| | C3 -C0 C1 -C2| |H3|

C0, C1, C2, C3 are the "db2" lowpass FIR filter coefficients. Highpass coefficients you get by reversing tap order and multiplying by sequence 1,-1, 1,-1, ... Because these are orthogonal wavelets, the analysis and reconstruction coefficients are the same.


A coefficient set convolved by its reverse is an ideal halfband lowpass filter multiplied by a symmetric windowing function. This creates the kind of symmetry in the frequency domain that enables aliasing-free reconstruction. Daubechies wavelets are the minimum-phase, minimum number of taps solutions for a number of vanishing moments (seven in "db7" etc.), which determines their frequency selectivity.


I was asked to show the matrices for 6 coefficients, so here they are, made a bit larger for clarity but could be the same size as before too. Decomposition:


|L0|   |C0  C1  C2  C3  C4  C5                | |I0|     

|H0| |C5 -C4 C3 -C2 C1 -C0 | |I1|
|L1| | C0 C1 C2 C3 C4 C5 | |I2|
|H1| | C5 -C4 C3 -C2 C1 -C0 | |I3|
|L2| = | C0 C1 C2 C3 C4 C5| |I4|
|H2| | C5 -C4 C3 -C2 C1 -C0| |I5|
|L3| |C4 C5 C0 C1 C2 C3| |I6|
|H3| |C1 -C0 C5 -C4 C3 -C2| |I7|
|L4| |C2 C3 C4 C5 C0 C1| |I8|
|H4| |C3 -C2 C1 -C0 C5 -C4| |I9|


Reconstruction:


|I0|   |C0  C5                  C4  C1  C2  C3| |L0|
|I1| |C1 -C4 C5 -C0 C3 -C2| |H0|
|I2| |C2 C3 C0 C5 C4 C1| |L1|
|I3| |C3 -C2 C1 -C4 C5 -C0| |H1|
|I4| = |C4 C1 C2 C3 C0 C5 | |L2|
|I5| |C5 -C0 C3 -C2 C1 -C4 | |H2|
|I6| | C4 C1 C2 C3 C0 C5 | |L3|
|I7| | C5 -C0 C3 -C2 C1 -C4 | |H3|
|I8| | C4 C1 C2 C3 C0 C5| |L4|

|I9| | C5 -C0 C3 -C2 C1 -C4| |H4|

With:


C0 = 3.326705529500826159985115891390056300129233992450683597084705e-01
C1 = 8.068915093110925764944936040887134905192973949948236181650920e-01
C2 = 4.598775021184915700951519421476167208081101774314923066433867e-01
C3 = -1.350110200102545886963899066993744805622198452237811919756862e-01
C4 = -8.544127388202666169281916918177331153619763898808662976351748e-02
C5 = 3.522629188570953660274066471551002932775838791743161039893406e-02


More coefficient sets can be found here.


No comments:

Post a Comment

digital communications - Understanding the Matched Filter

I have a question about matched filtering. Does the matched filter maximise the SNR at the moment of decision only? As far as I understand, ...