# 1D Haar Wavelets

This numerical tour explores 1D multiresolution analysis with Haar transform.

## Contents

- Installing toolboxes and setting up the path.
- Forward Haar Transform
- Backward Haar Transform
- Haar Approximation
- The Shape of a Wavelet

## Installing toolboxes and setting up the path.

You need to download the following files: signal toolbox and general toolbox.

You need to unzip these toolboxes in your working directory, so that you have `toolbox_signal` and `toolbox_general` in your directory.

For Scilab user: you must replace the Matlab comment '%' by its Scilab counterpart '//'.

Recommandation: You should create a text file named for instance `numericaltour.sce` (in Scilab) or `numericaltour.m` (in Matlab) to write all the Scilab/Matlab command you want to execute. Then, simply run `exec('numericaltour.sce');` (in Scilab) or `numericaltour;` (in Matlab) to run the commands.

Execute this line only if you are using Matlab.

    getd = @(p)path(p,path); % scilab users must *not* execute this

Then you can add the toolboxes to the path.

    getd('toolbox_signal/');
    getd('toolbox_general/');

## Forward Haar Transform

The Haar transform is the simplest orthogonal wavelet transform. It is computed by iterating difference and averaging between odd and even samples of the signal.

First we load a 1D signal.

    name = 'piece-regular';
    n = 512;
    f = rescale( load_signal(name, n) );

Initialize the transformed coefficients as the signal itself and set the initial scale as the maximum one. `fw` will be iteratively transformated and will contain the Haar coefficients.

    fw = f;
    j = log2(n)-1;

Select the sub-part of the image to transform.

    A = fw(1:2^(j+1));

Compute average and differences to get coarse and details. They are weighted by `1/sqrt(2)` to maintain orthogonality.

    Coarse = ( A(1:2:length(A)) + A(2:2:length(A)) )/sqrt(2);
    Detail = ( A(1:2:length(A)) - A(2:2:length(A)) )/sqrt(2);

Concatenate them to get the result.

    A = [Coarse; Detail];

Display the result of the first step of the transform.

    clf;
    subplot(2,1,1);
    plot(f); axis('tight'); title('Signal');
    subplot(2,1,2);
    plot(A); axis('tight'); title('Transformed');

Display the signal and its coarse coefficients.

    s = 400; t = 40;
    clf;
    subplot(2,1,1);
    plot(f,'.-'); axis([s-t s+t 0 1]); title('Signal (zoom)');
    subplot(2,1,2);
    plot(Coarse,'.-'); axis([(s-t)/2 (s+t)/2 min(A) max(A)]); title('Averages (zoom)');

Exercice 1: (the solution is exo1.m) Implement a full wavelet Haar transform that extract iteratively wavelet coefficients, by repeating these steps. Take care of choosing the correct number of steps.

    exo1;

Check that the transform is orthogonal, which means that the energy of the coefficient is the same as the energy of the signal.

    disp(strcat(['Energy of the signal       = ' num2str(norm(f).^2,3)]));
    disp(strcat(['Energy of the coefficients = ' num2str(norm(fw).^2,3)]));

    Energy of the signal       = 88.6
    Energy of the coefficients = 88.6

We display the whole set of coefficients `fw`, with red vertical separator between the scales. Can you recognize where are the low frequencies and the high frequencies ? You can use the function `plot_wavelet` to help you.

    clf;
    plot_wavelet(fw);
    axis([1 n -2 2]);

## Backward Haar Transform

The backward transform computes a signal `f1` from coefficients `fw`. If `fw` are the coefficients of a signal `f`, then one recovers `f=f1`. It is comuted with itertive up-sampling and filtering.

Initialize the image to recover `f1` as the transformed coefficient, and select the smallezt possible scale.

    f1 = fw;
    j = 0;

Retrieve coarse and detail coefficients in the vertical direction.

    Coarse = f1(1:2^j);
    Detail = f1(2^j+1:2^(j+1));

Undo the transform by computing sums and differences.

    f1(1:2:2^(j+1)) = ( Coarse + Detail )/sqrt(2);
    f1(2:2:2^(j+1)) = ( Coarse - Detail )/sqrt(2);

Exercice 2: (the solution is exo2.m) Write the inverse wavelet transform that computes `f1` from the coefficients `fw`.

    exo2;

Check that we have correctly recovered the signal.

    disp(strcat((['Error |f-f1|/|f| = ' num2str(norm(f-f1)/norm(f))])));

    Error |f-f1|/|f| = 1.0168e-15

## Haar Approximation

Non-linear approximation is obtained by thresholding low amplitude wavelet coefficients.

Set the threshold value.

    T = .5;

Coefficients `fw(i)` smaller in magnitude than `T` are set to zero.

    fwT = fw .* (abs(fw)>T);

Display the coefficients before and after thresholding.

    clf;
    subplot(2,1,1);
    plot_wavelet(fw); axis('tight'); title('Original coefficients');
    subplot(2,1,2);
    plot_wavelet(fwT); axis('tight'); title('Thresholded coefficients');

Exercice 3: (the solution is exo3.m) Find the threshold `T` so that the number of remaining coefficients in `fwT` are a fixzd number `m`. Use this threshold to compute `fwT` and then display the corresponding approximation `f1` of `f`. Try for an increasing number `m` of coeffiients.

    exo3;

## The Shape of a Wavelet

A wavelet coefficient `fw[k]` corresponds to an inner product `f,psi_{j,p}` where `k` depends on `j` and `p`.

The wavelet `f = psi_{j,p}` is computed by applying the inverse wavelet transform to `fw` where `fw[k]=1` and `fw[l]=0` for `l \neq k`.

Exercice 4: (the solution is exo4.m) Compute wavelets at several positions and scales.

    exo4;