One-dimensional example of the DWT with the Haar Wavelet

Setup

Define Quadrature Mirror Filters for the Haar Wavelet

(Subscript[h, 0])  =  (Subscript[h, 1])  =  1/(Sqrt[2]) and (Subscript[g, 0]) = -(Subscript[g, 1]) =  1/(Sqrt[2])

\(h[x_Integer]\  := \ 0; \ h[0]\  = \ 1/Sqrt[2]; \ h[1]\  = \ 1/Sqrt[2]; 
    
g[x_Integer]\  := \ 0; \ g[0]\  = \ 1/Sqrt[2]; \ 
    g[1]\  = \ \(-1\)/Sqrt[2]; \)

Define the "Operator Form" of the DWT

\(H[a_List, n_Integer]\  := \ 
      Table[h[k\  - \ 2  n], {k, 0, Length[a] - 1}] . a // N; 

    G[a_List, n_Integer]\  := \ 
      Table[g[k - 2  n], {k, 0, Length[a] - 1}] . a // N; \)

Define the "Operator Form" of the IDWT

\(Hstar[a_List, k_Integer]\  := \ 
      Table[h[k\  - \ 2  n], {n, 0, Length[a] - 1}] . a // N; 

    Gstar[a_List, k_Integer]\  := \ 
      Table[g[k\  - \ 2\ n], {n, 0, Length[a]\  - \ 1}] . a // N; \)

H[a_List]  :=  Table[Simplify[H[a, n]], {n, 0, Length[a]/2 - 1}]G[a_List]  :=  Table[Simplify[G[a, n]], {n, 0, Length[a]/2 - 1}]

\(Hstar[a_List]\  := \ 
      Table[Simplify[Hstar[a, n]], {n, 0, Length[a]*2 - 1}]; 

    Gstar[a_List]\  := \ 
      Table[Simplify[Gstar[a, n]], {n, 0, Length[a]*2 - 1}]; \)

Define some auxiliary functions that will compose other functions

This is a somewhat inefficient way of doing things, but it works for the 
demonstration.

compose[f_, g_]  :=  f[g[#]]&;

compose[f_, 0]  =  #&;

compose[f_, n_Integer]  :=  
      compose[f, n]  =  compose[f, compose[f, n - 1]];

Define the functions that will be called

DWT[a_List]  :=  
      Flatten[Append[
          Table[compose[G, compose[H, i]][a], {i, 0, 
              Log[2, Length[a]] - 1}], compose[H, Log[2, Length[a]]][a]]]

\(IDWT[a_List]\  := \ 
      Module[{l\  = \ Log[2, Length[a]]}, 
		
        Sum[\(compose[Hstar, j]\)[
              Gstar[a[\([
                    Range[2^l - 2^\((l - j)\) + 1, 
                      2^l\  - \ 2^\((l - j - 1)\)]]\)]]], {j, 0, l - 1}] + 
          \(compose[Hstar, l]\)[{a[\([Length[a]]\)]}]]\)

Define some auxiliary functions for denoising

thresh[a_List, x_]  :=  If[# >= x, #, 0]&/@a;

thresh2[a_List, x_]  :=  If[Abs[#]  >=  x, #, 0]&/@a;

threshmax[a_List, x_]  :=  thresh[a, Max[a]*x];

threshmax2[a_List, x_]  :=  thresh2[a, Max[Abs/@a]*x]

Plot2Lists[list1_List, list2_List, opts___]  :=  
      Show[ListPlot[list1, PlotJoined -> True, 
          PlotStyle -> {RGBColor[1, 0, 0]}, DisplayFunction -> Identity], 
        ListPlot[list2, PlotStyle -> {RGBColor[0, 0, 1]}, PlotJoined -> True, 
          DisplayFunction -> Identity], DisplayFunction -> $DisplayFunction, 
        opts]

Demonstration

The underlying curve is just a sine curve

udata  =  Table[Sin[x*2 Pi/512], {x, 0, 511}];

ListPlot[udata, PlotJoined -> True]

Data from that curve

data  =  
      Table[Sin[x*2 Pi/512], {x, 0, 511}]  +  
        (Table[Random[], {x, 0, 511}]* .5  -   .25);

ListPlot[data, PlotJoined -> True]

The first smooth is grossly oversmoothed

Note the sparsity of the wavelet coefficient vector.

wdata  =  DWT[data]

swdata  =  thresh[wdata, 10]

sdata  =  IDWT[swdata];

Plot2Lists[udata, sdata]

wdata  =  DWT[data];

Max[wdata]Min[wdata]

This is a better smoothing of the data

swdata  =  thresh2[wdata, 1]

sdata  =  IDWT[swdata];

Plot2Lists[sdata, udata]

A worse (or better) smooth?

Though we keep more coefficients, the smooth here seems to be worse. However, 
note that it is actually better in the middle. We happened to keep the 
coefficients of those wavelets that are approximating the middle of the plot 
range. Perhaps simple thresholding is not the best way to denoise.

swdata2  =  thresh[wdata,  .5]

sdata2  =  IDWT[swdata2];

Plot2Lists[sdata2, udata]

Undersmoothed approximation

This approximation below shows gross undersmoothing. The noise is "obvious" 
(especially since we are comparing the smooth with the underlying curve).

swdata3  =  thresh[wdata,  .2]

sdata3  =  IDWT[swdata3];

Plot2Lists[sdata3, udata]