==============================================================================
C-Scene Issue #2
Implementation of a lossless wavelet transform in C
Paul Brannan, a.k.a "cout"
==============================================================================

Introduction

Simply put, a wavelet can be thought of as a function that oscillates, similar to a wave, but as it tends toward plus or minus infinity, it approaches zero. The functions are used at different scales, because in functional analysis, different artifacts show up at different scales (see the picture of the Daubechies wavelet at the right). There are cert ainly more characteristics that mathematicians use for defining which wavelets are the most useful, but these two characteristics are what make wavelets in general a very exciting area to explore. The study of wavelets generally takes a strong mathematical background, but a simple Haar-based lossless compressor is a good place to start for those who haven't yet taken Calculus XVII and Super Advanced Linear Algebra.

Wavelets are used for many reasons. Wavelets allow signals to be decomposed into various parts, some of which can simply be thrown away. This is called quantization, and with DCT (the algorithm used in JPEG compression) can often result in "blockiness" in an image. With wavelet-based compression, such blockiness is minimized, thus allowing more of the image to be thrown away without suffering consequences in terms of image quality.

However, wavelets have another key characteristic. They have the ability to store an image at many resolutions, while minimizing overhead. Using wavelets, a smaller (lower resolution) version of the original image can be transmitted first, providing a rough sketch of what the image looks like. Next, instead of transmitting the entire version of a more detailed image, only the information necessary to derive a higher resolution image from the low resolution image is translated. As more of the image is transmitted, the decoded image increases in quality. This kind of transmission is known as progressive transmission.

There are a number of lossless methods for progressively transmitting images already existing today. One of the more common methods is used in the GIF format, a graphic file format common on the Internet. While to the user this kind of progressive transmission may seem to be very similar to wavelet-based compression, it is actually a completely different method. An exploration into the heart of the Haar transform will reveal the differences.

Background Information

History of Wavelets

In 1847, a man by the name of Joseph Fourier discovered a new way to transform functions, by representing the functions as a composite of sine waves. Fourier analysis is good for describing signals in a very concise manner. Sine waves, used in Fourier analysis, have the characteristic of repeating periodically forever. It is this property of sine waves that allows Fourier analysis to represent long, regular signals very well. However, Fourier analysis does not deal very well with signals that last a short amount of time (Von Baeyer 71).

In 1909, a man by the name of Haar discovered an alternative. He began with a function h(x) which was equal to 1 on [0,¹/₂), -1 on [¹/₂,1) and 0 outside the interval [0,1). Shown in Figure 1, this was the basis function for various "atoms" h_n(x) which were used to construct f(x) (Meyer 15-16).

Transforming Images with the Haar system

A 1-D signal can easily be transformed using the Haar system. The Haar "step function" can be amplified and shifted vertically to form the function shown in Figure 2. The values which represent the transformations (the amount the function is amplified/shifted) are actually equal to be the averages and differences of the two pixels. The average represents the amount to shift h(x) up, and the difference represents how much h(x) is to be amplified. For the pixels in Figure 2, h(x) would be shifted up 96 units, and amplified by -32 units (which is the difference divided by 2).

Images are often continuous functions, except near "edges," or areas of sharp contrast. Thus, in order to be able to use the Haar construction with images, it is first necessary to pass the image through a low-pass filter to soften the edges; in the Haar system, this is done using averages (Meyer 46-47). Next, the image is also passed through a high-pass filter, computed using the differences. The results of the two filters can actually be calculated in a single pass, as was done in Figure 3b. The results of the low-pass filter are in fact a scaled-down version of the picture, so they can be sent through the process again, as shown in Figure 3c (Gilbert).

In images, there is usually a lot more low frequency (changing slowly over the image) information than high frequency (quickly changing) information. Because of this, most of the values resulting from the high-pass filter are close to 0. The more of these values that are close to 0, the more efficiently the image can be compressed (Gilbert).

Generally, the difficult part of wavelet-based compression is defining the low-pass and high-pass filters. The filters should be orthogonal, meaning the output of each filter represents a different part of the signal, and they should have narrow support, meaning that they only affect a small area around the current pixel being analyzed (Meyer 59 and Gilbert). The Haar system is the only orthogonal system that meets these criteria. The complaint about the Haar system is that it produces harsh "edge effects" after the image has been quantized (Meyer 59). The generally accepted solution to this problem is a type of wavelet called a bi-orthogonal wavelet, but the Haar system should work well when the data is not quantized.

Comparison of Wavelets and Interlacing

The type of progressive transmission used today is known as interlacing and is used in the GIF format. Transforming the image this way produces four scaled-down images. The first image is derived from the first set of every fourth line, the second image from the second set of every fourth line, and so on. This results in an image similar to the one shown in Figure 4.

The Haar Transform

The general idea behind a wavelet-based transform is to use a high-pass and a low-pass filter. The high-pass (difference) filter is defined as:

D = a - b

where a and b are the two values to transform, and D is the difference between a and b. The low-pass (average) filter is defined as:

A = (a + b) / 2

where a and b are again the two values to transform, and A is the average between a and b. By taking values of every pair of vertically adjacent pixels in the image and transforming them into their average and difference, the entire image can be transformed. In the example program, the transform is only applied vertically, as a 1-D transform, so the results can be compared with GIF-style interlaced images.

However, in practice this algorithm does not work "as is" if an 8-bit char is used for the high-pass filter. One reason for this is simple: an 8-bit char cannot hold a signed 9-bit value. Invariably, problems will arise when one attempts to do this. It has been suggested that the high-pass filter should actually be stored as a 9-bit value; doing so can actually result in very good compression. There is, however, an alternative to this.

Implementation of the Haar transform

Three types of errors can occur in the transform, keeping the algorithm from being successfully reversed. On a computer, when the average of two values is calculated using integer arithmetic, some precision errors may occur. This is due to the decimal portion of the number being truncated. In fact, precision errors will always occur if the sum or difference of the two numbers is odd. This is a very simple problem to deal with. When the Haar transform is applied to an image, two numbers, the average and difference of the original two numbers, are stored for every pair of numbers in an image. The precision errors occur in the average. The difference can be used to determine if a precision error occurred; if this value is odd an error did in fact occur.

Another type of error also occurs, due in part to the way in which rounding errors are dealt with. The range of numbers is between -128 and 127, for a given value in the images used in this paper. If the sum of the two numbers is odd, then a precision error will occur, as already discussed. However, if the sum of the two numbers is also negative, then it is important to realize how the numbers will be truncated. With computers, numbers are always truncated by simply chopping off the decimal portion of the number. In effect, all numbers are rounded toward zero using this method. This means that negative numbers are always rounded up, instead of down, like positive numbers. In order to use the fix described previously for dealing with precision errors, a statement must be added to the program to allow negative numbers to be rounded down instead of up. A "floor" function (greatest integer less than or equal to) must be implemented instead of truncation.

Finally, if the difference between two numbers is larger than the set limits of the averaging function, then the average of the two numbers will be off in the eighth bit. This was an extremely tough problem to diagnose. However, once the problem was recognized, the fix was easy. If the difference between the two numbers is less than -128 or greater than 127, then the eighth bit of the average must be set to represent the overflow.

The question might be asked "if these errors are occurring in the image, then the transform really isn't lossless, is it?" Well, of course it is lossless. If the errors were left uncorrected, then the transform would be lossy. Try the program haartest.c to see that the transform really is lossless.

Pseudocode for a two-pixel transform:

Encoding:

A = (a + b) / 2;		/* Low-pass (average) filter */
D = a - b;			/* High-pass (difference) filter */
if((a + b < 0) and (D & 1)	/* Use floor function if A was */
A = A - 1;   			/* Rounded off */
if((a - b < -128) or (a - b > 127)
A = A + 128;   			/* See if the difference is out of */
				/* range */

Decoding:

a = A + (D / 2);   		/* Calculate the first value */
b = A - (D / 2);   		/* Calculate the second value */
if (D & 1)    			/* Test to see if D is odd */
if (D < 0) b++ else a++; 	/* Increment the larger value */

The Wavelet Transform Programs

The programs tga2har.c and tga2int.c are designed to transform a TGA image into a Haar-transformed image and an interlaced image, respectively. They were not originally written for readability, but they do employ some interesting techniques. They have not been tested on images other than 320x200x256 uncompressed TGA's, so don't expect them to work on anything else (though they may, but the results could be unpredictable). Try running the programs using one of these images as an input file. You will see that the programs write directly to a buffer in screen memory. Speed-wise this is very inefficient, but it allows the program to compile properly on a 16-bit real-mode compiler.

Note that neither of these programs perform any compression. All that is done is that the image is tranformed into another form, in hopes that better compression can be achieved by an external program. In reality, these programs do not generally crea te more compressible images, but they do work as good examples. See http://www.musc.edu/~brannanp/wavelets/ for more information on how the algorithms perform in the real world.

Oddly enough, people don't like to transform images into forms that they cannot later undo. Thus, it is impractical to have a program that performs a wavelet transform on an image, but not have a program that can display the transformed image. That i s what the program loadhar.c is for. Again, the entire tranform is performed inside of screen memory, so a compiler thatcan directly write to segment 0xA000 is needed in order to run this program without modification.

These programs are, however, fairly easy to port to other platforms. By simply rewriting the functions text(), graphics(), and putline(), the programs should run on any platform that supports graphics display.

Discussion

Hopefully, this article has given presented a small bit of what wavelets are. There are a number of wavelet applications that are waiting to be discovered. Wavelets, since the discovery Daubechies made of an orthogonal wavelet system with compact support have been finding their way into image compression, audio compression, function analysis, artificia l intelligence and sight, fractal geometry studies, musical decomposition, and all sorts of other areas. Some sources of information have been provided here for those who may be interested in pursuing wavelet research.

The wavelet system used here is hardly state-of-the-art. It was originally designed as a scientific comparison between GIF interlacing and Haar wavelets, and is now useful for little more than a simple example. Much is left for the reader to complete as an exercise :). For example, a wavelet transform is not generally performed in a single dimension in images; it is only performed in one dimension at a time, but can be alternated between x and y to achieve a two-dimensional transform (see the CREW l ink for a better explanation of this). Also, the transform may create even more compressible data if the palette is sorted before the transform is performed; see Jakulin's PROPOSAL for more information on this.

Bibliographic Information

Apiki, Steve. "Lossless Data Compression." Byte Mar. 1991: 309-387.

This article discusses two common methods of data compression, Huffman and LZW, as well as common variations of these methods. Both of these methods are lossless, meaning that the decompressed version of the file is the same as the original. Sample source code is also provided.

E.R.I.C.: An Algorithm for Efficient Reversible Image Compression. URL: http://www-ias.jpl.nasa.gov/HPCC/eric/eric.html.

A proposal for a lossless wavelet-based image coder.

Title: Compression Pointers. URL: http://www.internz.com/compression-pointers.html.

A very good list of comrpession links. Includes links to FBI Fingerprint research (which uses wavelets), along with a wealth of other links.

Gilbert, David Alan. "Re: Math Retard Needs Help." Usenet newsgroup comp.compression.research 5 Jun 1995.

This is an extremely well written explanation of compression algorithms using wavelets. The algorithms described in this document are the algorithms upon which the programs in this project were based.

Jakulin, Alex. "PROPOSAL for lossless image compression." Usenet newsgroup comp.compression 13 Jan 1995.

Using the delta transform on an image prior to compressing it with the deflate algorithm can decrease the quality of compression. A new method is presented which integrates RLE with the deflate algorithm, by introducing a new symbol that indicates an RLE code follows. Also recommended is sorting the palette by intensities and then delta-coding it; the palette would then be more compressible, as would an image which is delta-coded.

Meyer, Yves. Wavelets: Algorithms and Applications. Translated and revised by Robert D. Ryan. Philadelphia: Society for Industrial and Applied Mathematics, 1993.

This book is an extensive collection of chapters on various topics relating to wavelets, including background on wavelets, audio and video compression applications, and visual systems of robots. The book is easy for a novice to comprehend, provided he/she has some sort of a background in basic calculus.

Title: RICOH CRC: CREW technology. URL: http://yellow.crc.ricoh.com/CREW/.

This is a very good explanation of how an efficient lossless wavelet-based image coder should work. It claims better compression than JBIG.

Rimmer, Steve. Supercharged Bitmapped Graphics. USA: Windcrest Books, 1992.

This book explains the details of a number of graphics file formats that are commonly used. The GIF format (which uses interlacing, a method of progressive transmission) and the TGA format (a format that is easy to implement, and was used in this project) were both discussed in detail.

Stix, Gary. "Encoding the 'Neatness' of Ones and Zeroes." Scientific American Sept. 1991: 54-58.

In 1951, David Huffman developed, for a term paper, a compression methods that was the most efficient method at the time. His professor, Fano, and Claude F. Shannon had struggled with the problem for many years before Huffman's discovery. This article discusses basic coding trees, used in Huffman compression, though the main purpose of the article is to provide a brief biography of David Huffman.

Von Baeyer, Hans Christian. "Wave of the future," Discover May 1995: 69-74.

This source describes the history of Daubechies's wavelets. In addition, and explanation of wavelets that is easy to understand is also included. This source is almost a necessity in order to understand the other sources, which require high levels of mathematics in order to comprehend. You can download example source for lossless wavelets here.