Lazy Wavelet Transform Files

Lazy (linear prediction) Wavelet Transform with Lifting
Author: Arkadi Kagan

This project consists from the follow parts:
1. lazy_wavelet.c, lazy_wavelet_2d.c, lazy_wavelet_byte.c and lazy_wavelet_2d_byte.c
   implementation files. Each implementation is entairly stand-along library. There is no
   inter-dependencies of any kind or common utilities. The only dependency allowed is on
   standard libraries.
2. test_1d to check 1-D algorithm validity
3. test_2d to check 2-D algorithm validity
4. test_video and all its dependencies

The "test_video" sample is intended to demonstrate wavelet transform on a video stream.
This sample code is doing the follow:
* Read Microsoft provided AVI "SAMPLE.AVI" frame by frame
* Each frame is converted from RGB to gray levels 0-255,
  the transformed video is stored in "sample.gray.avi"
* Each gray sample is transformed into array of wavelet coeffecients, stored in
  "sample.coeff.avi" for visual review
* The computed coeffecients are downsampled to fit into a byte
* The sample cut-off 1 / 16 of all coefficients and store in "sample.fwt", set zero all the rest
* Restore the original image with inverse wavelet transform and save into "sample.out.avi"
Note that "sample.fwt" contain all information that was used to restore the video into
"sample.out.avi". However, this project do not have an aim to create real compression codec. Instead,
we demonstrate that wavelet transform by itself can perform as a competitative image compression
and analysis tool.

The wavelet transform code strictly follows recipe from:
Ingrid Doubaches and Wim Sweldens, "Factoring Wavelet Transform into Lifting Steps",
Journal of Fourier Analysis and Applications, Volume 4, Number 3, May 1998, pages 247-269,
DOI: 10.1007/BF02476026

The simplified process of Fast Discrete Wavelet Transform can be reformulated as follow:
* Split the source sequence into odd and even members.
* Predict odd members by corresponding values of the surrounding even members.
  In our case, we compute mean x[i - 1] and x[i + 1].
  We store delta = x[i] - Predicted(x[i]) as a coefficient.
* Update the even members to restore the mean value to be equal to the overall mean value.
  In our case we add (delta[i - 1] + delta[i + 1]) / 4 to each even member.
After this process is completted, we have a set of updated original values (even members)
and a separate set of coefficients to be stored. The updated original values are farther
processed in recursion.

The defined above method computes Lazy (linear prediction) Wavelet Transform. In this project
we do not attempt to create general code that fit all possible Wavelet Transforms. Instead,
our aim is to create a single efficient implementation that can be used as an auxiliary
tool for other tasks.

As a natural extension to the 1-D wavelet transform, we define here a 2-D transform. Our choice
of decomposition algorithm is piramidal decomposition. The pyramid decomposition is defined in
David Salomon, "Data Compression: The Complete Reference, Fourth Edition", section 5.6.1, p. 20
Here is a quotation:
"The latter method computes the wavelet transform of the image by alternating between rows and columns. The first step is
 to calculate averages and differences for all the rows (just one iteration, not the entire wavelet transform). This creates
 averages in the left half of the image and differences in the right half. The second step is to calculate averages and
 differences (just one iteration) for all the columns, which results in averages in the top-left quadrant of the image and
 differences elsewhere."

In our implementation great care was taken to allow image size that is not limited to power of two. To clarify
the problem solved, lets take a 1-D sample of length 5:
 x = {0, 1, 2, 3, 4}
If we split this sequence into odd/even as suggested by Ingrid Doubaches and Wim Sweldens, we get
 {0, 2, 4, 1, 3}
As we can see, there are 3 even numbers and only 2 odd numbers. When we are trying to restore
the original, how do we split into {0, 2, 4}v{1, 3}, and not into {0, 2}v{4, 1, 3}?
In one-step process the answer is obvious, but lets take one step farther: split {0, 2, 4} into
{0, 4} and {2}. In this case we have {0, 4, 2, 1, 3} and we want to know that {0, 4}v{2} is the first
block of data to process. Solution to the above problem is to calculate all the block lengths ahead
before any transform is performed. This implies that we need a log_2(length) buffer to store the
block lengths.

Integer arithmetic transform.

A rather important modification of floating point transform can be found in lazy_wavelet_byte.c and
lazy_wavelet_2d_byte.c. Instead of doing all the calculation with float or double, we do all the same
algorithm with signed char and intermediate results we keep in integer variables.

This transform is less generic than floating point variant, because values are limited into a single
signed char, that is ranging from -127 to +127. However, this transform is sufficient for video coding,
as has been shown in test_video.c sample.

Use of integer arithmetic only opens a door for optimizations like Intel MMX or GPU API's for parallel
processing.