Transcript Signal Processing with Wavelets

Chapter 11
Signal Processing with Wavelets
Objectives
• Define and illustrate the difference between a stationary
and non-stationary signal.
• Describe the relationship between wavelets and subband coding of a signal using quadrature mirror filters
with the property of perfect reconstruction.
• Illustrate the multi-level decomposition of a signal into
approximation and detail components using wavelet
decomposition filters.
• Illustrate the application of wavelet analysis using
MATLAB® to noise suppression, signal compression, and
the identification of transient features in a signal.
Motivation for Wavelet Analysis
• Signals of practical interest are usually nonstationary, meaning that their time-domain and
frequency-domain characteristics vary over short
time intervals (i.e., music, seismic data, etc)
• Classical Fourier analysis (Fourier transforms)
assumes a signal that is either infinite in extent
or stationary within the analysis window.
• Non-stationary analysis requires a different
approach: Wavelet Analysis
• Wavelet analysis also produces better solutions
to important problems such as the transform
compression of images (jpeg versus jpeg2000)
Basic Theory of Wavelets
• Wavelet analysis can be understood as a
form of sub-band coding with quadrature
mirror filters
• The two basic wavelet processes are
decomposition and reconstruction
Wavelet Decomposition
•
•
A single level decomposition puts a signal through 2 complementary low-pass and
high-pass filters
The output of the low-pass filter gives the approximation (A) coefficients, while the
high pass filter gives the detail (D) coefficients
Low-Pass
Filter
DownSample
2X
Approximation (A)
Signal
High-Pass
Filter
Decomposition Filters for Daubechies-8 Wavelets
DownSample
2X
Detail (D)
Wavelet Reconstruction
• The A and D coefficients can be used to
reconstruct the signal perfectly when run
through the mirror reconstruction filters of
the wavelet family
Wavelet Families
• Wavelet families consist of a particular set
of quadrature mirror filters with the
property of perfect reconstruction.
• These families are completely determined
by the impulse responses of the set of 4
filters.
Example:
Filter Set for the Daubechies-5 Wavelet Family
% Set wavelet name.
>> wname = 'db5';
% Compute the four filters associated with wavelet name given
% by the input string wname.
>> [Lo_D,Hi_D,Lo_R,Hi_R] = wfilters(wname);
>> subplot(221); stem(Lo_D);
>> title('Decomposition low-pass filter');
>> subplot(222); stem(Hi_D);
>> title('Decomposition high-pass filter');
>> subplot(223); stem(Lo_R);
>> title('Reconstruction low-pass filter');
>> subplot(224); stem(Hi_R);
>> title('Reconstruction high-pass filter');
>> xlabel('The four filters for db5')
Example:
Filter Set for the Daubechies-5 Wavelet Family
Decomposition low-pass filter
Decomposition high-pass filter
1
1
0.5
0.5
0
0
-0.5
-0.5
0
5
10
-1
0
Reconstruction low-pass filter
5
10
Reconstruction high-pass filter
1
1
0.5
0.5
0
0
-0.5
-0.5
0
5
10
-1
0
5
The four filters for db5
10
Example:
Filter Set for the Daubechies-5 Wavelet Family
>> fvtool(Lo_D,1,Hi_D,1)
>> fvtool(Lo_R,1,Hi_R,1)
Magnitude Response
Magnitude Response
1.5
1
1
Magnitude
Magnitude
1.5
0.5
0.5
0
0
0
0.1
0.2
0.3
0.4
0.5
0.6
0.7
Normalized Frequency ( rad/sample)
Decomposition Filters
0.8
0.9
0
0.1
0.2
0.3
0.4
0.5
0.6
0.7
Normalized Frequency ( rad/sample)
Reconstruction Filters
0.8
0.9
Multi-level Decomposition of a Signal with
Wavelets
Signal
A1
A2
A3
D1
D2
D3
The decomposition tree
can be schematically
described as:
Aj = Aj+ 1 + Dj+ 1
Multi-level Decomposition of a Signal with Wavelets
Frequency Domain (Sub-band Coding)
A2
D1
D2
fs/8
fs/4
fs/2
Example: One-level Decomposition of a Noisy Signal
>> x=analog(100,4,40,10000);
>> xn=x+0.5*randn(size(x));
>> [cA,cD]=dwt(xn,'db8');
% Construct a 100 Hz sinusoid of amplitude 4
% Add Gaussian noise
% Compute the first level decomposition with dwt
% and the Daubechies-8 wavelet
>> subplot(3,1,1),plot(xn),title('Original Signal')
>> subplot(3,1,2),plot(cA),title('One Level Approximation')
>> subplot(3,1,3),plot(cD),title('One Level Detail')
Original Signal
5
0
-5
0
50
100
150
200
250
300
350
400
450
One Level Approximation
10
0
-10
0
50
100
150
200
250
200
250
One Level Detail
2
0
-2
0
50
100
150
Single level discrete
wavelet decomposition
with the Daubechies-8
wavelet family
One-Level Decomposition of a Non-Stationary Signal
>> fs=2500;
>> len=100;
>> [x1,t1]=analog(50,.5,len,fs); % The time vector t1 is
in milliseconds
>> [x2,t2]=analog(100,.25,len,fs);
>> [x3,t3]=analog(200,1,len,fs);
>> y1=cat(2,x1,x2,x3); % Concatenate the signals
>> ty1=[t1,t2+len,t3+2*len]; %Concatenate the time
vectors 1 to len, len to 2*len, etc.
>> [A1,D1]=dwt(y1,'db8');
>> subplot(3,1,1),plot(y1),title('Original Signal')
>> subplot(3,1,2),plot(A1),title('One Level
Approximation')
>> subplot(3,1,3),plot(D1),title('One Level Detail')
Original Signal
1
0
-1
0
100
200
300
400
500
600
700
800
300
350
400
300
350
400
One Level Approximation
2
The detail
coefficients reveal
the transitions in the
non-stationary signal
0
-2
0
50
100
150
200
250
One Level Detail
0.2
0
-0.2
0
50
100
150
200
250
De-Noising a Signal with Multilevel Wavelet Decomposition
>> x=analog(100,4,40,10000);
>> xn=x+0.5*randn(size(x));
>> [C,L] = wavedec(xn,4,'db8'); % Do a multi-level
analysis to four levels with the
% Daubechies-8
wavelet
>> A1 = wrcoef('a',C,L,'db8',1); % Reconstruct the
approximations at various levels
>> A2 = wrcoef('a',C,L,'db8',2);
>> A3 = wrcoef('a',C,L,'db8',3);
>> A4 = wrcoef('a',C,L,'db8',4);
>> subplot(5,1,1),plot(xn),title('Original Signal')
>> subplot(5,1,2),plot(A1),title('Reconstructed
Approximation - Level 1')
>> subplot(5,1,3),plot(A2),title(' Reconstructed
Approximation - Level 2')
>> subplot(5,1,4),plot(A3),title(' Reconstructed
Approximation - Level 3')
>> subplot(5,1,5),plot(A4),title(' Reconstructed
Approximation - Level 4')
Original Signal
10
0
-10
5
0
-5
Significant de-noising occurs
with the level-4 approximation
coefficients (Daubechies-8
wavelets)
5
0
-5
5
0
-5
5
0
-5
0
50
100
150
200
250
300
350
Reconstructed Approximation - Level 1
400
450
0
50
100
150
200
250
300
350
Reconstructed Approximation - Level 2
400
450
0
50
100
150
200
250
300
350
Reconstructed Approximation - Level 3
400
450
0
50
100
150
200
250
300
350
Reconstructed Approximation - Level 4
400
450
0
50
100
400
450
150
200
250
300
350
Finding Signal Discontinuities
>> x=analog(100,4,40,10000);
>> x(302:305)=.25;
>> [A,D]=dwt(x,'db8');
>> subplot(3,1,1),plot(x),title('Original Signal')
>> subplot(3,1,2),plot(A),title('First Level Approximation')
>> subplot(3,1,3),plot(D),title('First Level Detail')
3 sample discontinuity
at sample 302
Original Signal
5
Daubechies-8 wavelets
0
-5
0
50
100
150
200
250
300
350
400
450
First Level Approximation
10
0
-10
0
50
100
150
200
250
200
250
First Level Detail
0.5
0
-0.5
0
50
100
150
Discontinuity response
in the 1st level detail
coefficients (sample
151 because of the 2X
down-sampling)
Simple Signal Compression Using a Wavelet
Approximation
>> load leleccum
>> x=leleccum;
>> w = 'db3';
>> [C,L] = wavedec(x,4,w);
>> A4 = wrcoef('a',C,L,'db3',4);
>> A3 = wrcoef('a',C,L,'db3',3);
>> A2 = wrcoef('a',C,L,'db3',2);
>> A1 = wrcoef('a',C,L,'db3',1);
>> a3 = appcoef(C,L,w,3);
>> subplot(2,1,1),plot(x),axis([0,4000,100,600])
>> title('Original Signal')
>> subplot(2,1,2),plot(A3),axis([0,4000,100,600])
>> title('Approximation Reconstruction at Level 3 Using
the Daubechies-3 Wavelet')
>> (length(a3)/length(x))*100
ans =
12.5926
Original Signal
600
500
400
300
200
100
0
500
1000
1500
2000
2500
3000
3500
4000
Approximation Reconstruction at Level 3 Using the Daubechies-3 Wavelet
600
500
400
300
200
100
0
500
1000
1500
2000
2500
3000
3500
4000
The wavelet approximation at
level-3 contains only 13 % of the
original signal values because of
the wavelet down-sampling, but
still retains the important signal
characteristics.
Compression by Thresholding
>> load leleccum
>> x=leleccum;
>> w = 'db3'; % Specify the Daubechies-4 wavelet
>> [C,L] = wavedec(x,4,w); % Multi-level decomposition to 4 levels.
>> a3 = appcoef(C,L,w,3); % Extract the level 3 approximation coefficients
>> d3 = detcoef(C,L,3);
% Extract the level 3 detail coefficients.
>> subplot(2,1,1), plot(a3),title('Approximation Coefficients at Level 3')
>> subplot(2,1,2), plot(d3),title('Detail Coefficients at Level 3')
Approximation Coefficients at Level 3
2000
1500
1000
500
0
0
100
200
300
400
500
600
500
600
Detail Coefficients at Level 3
50
0
-50
0
100
200
300
400
These are the A3 and D3
coefficients for the signal.
Many of the D3
coefficients could be
“zeroed” without losing
much signal information or
power
Compression by Thresholding
>> load leleccum
>> x=leleccum; % Uncompressed signal
>> w = 'db3';
% Set wavelet family
>> n=3;
% Set decomposition level
>> [C,L] = wavedec(x,n,w); % Find the decomposition
structure of x to level n using w.
>> thr = 10; % Set the threshold value
>> keepapp = 1; %Logical parameter = do not threshold
approximation coefficients
>> sorh='h'; % Use hard thresholding
>> [xd,cxd,lxd, perf0,perfl2]
=wdencmp('gbl',C,L,w,n,thr,sorh,keepapp);
>> subplot(2,1,1), plot(x),title('Original Signal')
>> subplot(2,1,2),plot(xd),title('Compressed Signal
(Detail Thresholding)')
>> perf0 % Percent of coefficients set to zero
>> perfl2 % Percent retained energy in the compressed
signal
perf0 =
83.4064
perfl2 =
99.9943
Original Signal
600
400
In this compression
83% of the coefficients
were set to zero, but
99% of the energy in
the signal was
retained.
200
0
0
500
1000
1500
2000
2500
3000
3500
4000
4500
4000
4500
Compressed Signal (Detail Thresholding)
600
400
200
0
0
500
1000
1500
2000
2500
3000
3500
Compression by Thresholding
>> D1 = wrcoef('d',C,L,w,1);
>> D2 = wrcoef('d',C,L,w,2);
>> D3 = wrcoef('d',C,L,w,3);
>> d1 = wrcoef('d',cxd,lxd,w,1);
>> d2 = wrcoef('d',cxd,lxd,w,2);
>> d3 = wrcoef('d',cxd,lxd,w,3);
Original Detail - Levels 3 to 1
>> subplot(3,2,1),plot(D3),title('Original Detail - Levels 3
to 1')
>> subplot(3,2,2),plot(d3),title('Thresholded Detail Levels 3 to 1')
>> subplot(3,2,3),plot(D2)
>> subplot(3,2,4),plot(d2)
>> subplot(3,2,5),plot(D1)
>> subplot(3,2,6),plot(d1)
40
Thresholded Detail - Levels 3 to 1
40
20
20
0
0
-20
0
2000
4000
6000
-20
40
40
20
20
0
0
-20
0
2000
4000
6000
-20
50
50
0
0
-50
0
2000
4000
6000
-50
0
2000
4000
6000
0
2000
4000
6000
0
2000
4000
6000
Zeroing of coefficients
by thresholding
results in effective
signal compression
Summary
• Wavelet processing is based on the idea of subband decomposition and coding.
• Wavelet “families” are characterized by the lowpass and high-pass filters used for
decomposition and perfect reconstruction of
signals.
• Typical applications of wavelet processing
include elimination of noise, signal compression,
and the identification of transient signal features.