The Discrete Periodic Wavelet Transform in 1D:

A Demonstration of the Routine Wavelet1Dv1.m

by 

James F. Scholl
Mathematical Physics Department 
Rockwell Science Center
1049 Camino Dos Rios
Thousand Oaks, CA 91360

e-mail: jfs@risc.rockwell.com

Calling the Routine, With Some Wavelet Bases to Try

<<"Wavelet1Dv1.m"

the main functions:

?DPWT1D

?IDPWT1D

Some Basis system filters:

haar = {0.5, 0.5};

D4 = 0.5 * {N[(1 + Sqrt[3])/4], N[(3 + Sqrt[3])/4], 
      N[(3 - Sqrt[3])/4], N[(1 - Sqrt[3])/4]};

Coif2 = {-0.051429728471,0.238929728471,
 0.602859456942, 0.272140543058, -0.051429972847,
 -0.011070271529};

Coif4 = {0.011587596739, -0.02932013798, -0.04763959031,
     0.273021046535, 0.574682393857, 0.294867193696,
     -0.054085607092, -0.042026480461, 0.016744410163, 
     0.003967883613, -0.001289203356, -0.000509505399};

Test Data

Random data of length 32 to demonstrate perfect reconstruction:

testdata1 = Table[Random[Real, {-1, 1}], {32}]

A Hermite Function Modulated Gaussian with noise to Demonstrate subpictures and scale:

signal = (1./250.) Table[(HermiteH[6,x] Exp[-x^2/2])^2, 
     {x, -6.35, 6.35, 0.1}];

Length[signal]

noise = Table[Random[Real, {-10, 10}], {128}];

testdata2 = signal + noise;

{Min[testdata2], Max[testdata2]}

Signal plot:

ListPlot[signal, PlotJoined->True, PlotRange->{-20.,120.},
  PlotLabel->"(Hermite[6,x] Exp[-x^2/2])^2"]

Noise Plot:

ListPlot[noise, PlotJoined->True, PlotRange->{-20.,120.},
  PlotLabel->"noise"]

Total Test Signal Plot:

ListPlot[testdata2, PlotJoined->True, PlotRange->{-20.,120.},
  PlotLabel->"signal and noise"]

Test Calculation for Perfect Reconstruction

Random signals make very robust tests for perfect reconstruction:

Wavelet Transform will be done w.r.t the Daubechies D4 system

First I will do the calculation with too many levels to illustrate this error handling condition:

DPWT1D[testdata1, D4, 7]

Next will be the full decomposition (lev=5); this will also demonstrate how the output data is arranged:

wcof1 = DPWT1D[testdata1, D4, 5]

Next I will reconstruct the data:

testdata1R = IDPWT1D[wcof1, D4]

The residual between the reconstructed and original signal; note the perfect reconstruction!!

testdata1 - testdata1R

Test Calculation Illustrating Subsignal Manipulation and Simple-Minded Noise Removal

We can separate signal from noise by looking at the different scale subsignals. We will try to remove the noise.

Wavelet Transform will be done w.r.t the 4th order Coifman basis (Coif4) with 4 levels of decomposition

Wcof2 = DPWT1D[testdata2, Coif4, 4];

These are the component subsignals in the decomposition:

A4 = Wcof2[[1]]; Length[A4]

D4 = Wcof2[[2,1]]; Length[D4]

D3 = Wcof2[[2,2]]; Length[D3]

D2 = Wcof2[[2,3]]; Length[D2]

D1 = Wcof2[[2,4]]; Length[D1]

We will plot the subsignals to illustrate how the output is arranged as well as the multiresolution aspect of wavelets. We will plot by increasing frequency.

ListPlot[A4, PlotJoined->True, PlotRange->All,
  PlotLabel->"A4 subsignal / coefficients of W.T.",
  PlotStyle->RGBColor[1,0,0]]

ListPlot[D4, PlotJoined->True, PlotRange->All,
  PlotLabel->"D4 subsignal / coefficients of W.T.",
  PlotStyle->RGBColor[0,0,1]]

ListPlot[D3, PlotJoined->True, PlotRange->All,
  PlotLabel->"D3 subsignal / coefficients of W.T.",
  PlotStyle->RGBColor[0,0,1]]

ListPlot[D2, PlotJoined->True, PlotRange->All,
  PlotLabel->"D2 subsignal / coefficients of W.T.",
  PlotStyle->RGBColor[0,0,1]]

ListPlot[D1, PlotJoined->True, PlotRange->All,
  PlotLabel->"D1 subsignal / coefficients of W.T.",
  PlotStyle->RGBColor[0,0,1]]

Note the noise manifests itself in the highest frequency subsignal; lets try to remove it by setting these coefficients to zero:

D1prime = Table[0.0, {64}]

The revised coefficients for the inverse wavelet transform: This also illustrates how we can do some signal processing.

Wcof2a = {A4, {D4, D3, D2, D1prime}};

testdata2R = IDPWT1D[Wcof2a, Coif4];

A plot of the reconstructed signal:

recon1 = ListPlot[testdata2R, PlotRange->All,
  PlotJoined->True, PlotStyle->RGBColor[1,0,0],
  PlotLabel->"signal with some noise removed"]

A plot of the original signal without the noise:

original = ListPlot[signal, PlotRange->All,
 PlotJoined->True, PlotStyle->RGBColor[0,0,1],
  PlotLabel->"original signal w/o noise"]

Show[original, recon1, PlotLabel->None]

A Plot of the Residuals in the same scale as the main signals; not all the noise was taken out but you see the idea

ListPlot[signal-testdata2R, PlotRange->{-20., 120.},
 PlotJoined->True, PlotStyle->RGBColor[0,0.5,0],
  PlotLabel->"residual signal"]

Here is how much noise we did remove:

ListPlot[testdata2-testdata2R, PlotRange->{-20., 120.},
 PlotJoined->True, PlotStyle->RGBColor[0,0.5,0],
  PlotLabel->"noise removed"]