Image Compression:

How Math Led to the JPEG2000 Standard

JPEG2000 Wavelet Transformations

Introduction

The JPEG2000 Image Compression Standard makes use of two biorthogonal wavelet filter pairs. The Cohen-Daubechies-Feaveau filter pair (CDF97)is used for lossy compression and the LeGall filter pair is used to perform lossless image compression. Both filter pairs are odd-length and symmetric. As we will see later on the page, symmetry is important for handling problems that can occur at the edges of an image. The CDF97 pair uses a length 9 lowpass filter and a length 7 highpass filter. The LeGall pair has a length 5 lowpass filter and a length 3 highpass filter. For the inverse transformation, the lengths of the lowpass and highpass filters are reversed (e.g., for the inverse transform using CDF97, the lowpass filter is length 7 and the highpass is length 9). It was important to JPEG committee that each filter pair have similar lengths. Before describing each filter pair, let's see why symmetric filters are important for handling edges in images.

Filter Symmetry and Edge Effects

For the sake of illustration, suppose we have a length 5 symmetric lowpass filter \tilde{{\bf h}} = (\tilde{h}_2, \tilde{h}_1, \tilde{h}_0, \tilde{h}_1, \tilde{h}_2) and a length 5 highpass filter \tilde{{\bf g}} = (\tilde{g}_0, \tilde{g}_1, \tilde{g}_2). Recall, that \tilde{{\bf g}} is constructed from symmetric, length 3 lowpass filter {\bf h} = (h_1, h_0, h_1) with \tilde{g}_k = (-1)^k h_{1-k} so \tilde{{\bf g}} = (\tilde{g}_0, \tilde{g}_1, \tilde{g}_2) = (h_1, -h_0, h_1).

If we were to process a length 8 vector with this transformation, the matrix product would look like

\left[\matrix{ \tilde{h}_0 & \tilde{h}_1 & \tilde{h}_2 & 0 & 0 & 0 & \tilde{h}_2 & \tilde{h}_1 \\ \tilde{h}_2 & \tilde{h}_1 & \tilde{h}_0 & \tilde{h}_1 & \tilde{h}_2 & 0 & 0 & 0 \\ 0 & 0 & \tilde{h}_2 & \tilde{h}_1 & \tilde{h}_0 & \tilde{h}_1 & \tilde{h}_2 & 0 \\ \tilde{h}_2 & 0 & 0 & 0 & \tilde{h}_2 & \tilde{h}_1 & \tilde{h}_0 & \tilde{h}_1 \\ h_1 & -h_0 & h_1 & 0 & 0 & 0 & 0 & 0 \\ 0 & 0 & h_1 & -h_0 & h_1 & 0 & 0 & 0 \\ 0 & 0 & 0 & 0 & h_1 & -h_0 & h_1 & 0 \\ h_1 & 0 & 0 & 0 & 0 & 0 & h_1 & -h_0 }\right] \cdot \left[ \matrix{ v_1 \\ v_2 \\ v_3 \\ v_4 \\ v_5 \\ v_6 \\ v_7 \\ v_8}\right] = \left[\matrix{ \tilde{h}_2 v_7 + \tilde{h}_1 v_8 + \tilde{h}_0 v_1 + \tilde{h}_1 v_2 + \tilde{h}_2 v_3 \\ \tilde{h}_2 v_1 + \tilde{h}_1 v_2 + \tilde{h}_0 v_3 + \tilde{h}_1 v_4 + \tilde{h}_2 v_5 \\ \tilde{h}_2 v_3 + \tilde{h}_1 v_4 + \tilde{h}_0 v_5 + \tilde{h}_1 v_6 + \tilde{h}_2 v_7 \\ \tilde{h}_2 v_5 + \tilde{h}_1 v_6 + \tilde{h}_0 v_7 + \tilde{h}_1 v_8 + \tilde{h}_2 v_1 \\ h_1 v_1 - h_0 v_2 + h_1 v_3 \\ h_1 v_3 - h_0 v_4 + h_1 v_5 \\ h_1 v_5 - h_0 v_6 + h_1 v_7 \\ h_1 v_7 - h_0 v_8 + h_1 v_1 }\right]

Note that rows 1, 4, and 8 of the output suffer from the "wrapping effect" - the weighted averages and differences are constructed using elements at the top and bottom of the input vector. We exploit the symmetry of the filters by applying the wavelet transformation to a periodized version of the input: (v_1, v_2, v_3, v_4, v_5, v_6, v_7, v_8, v_7, v_6, v_5, v_4, v_3, v_2):

\left[\matrix{ \tilde{h}_0 & \tilde{h}_1 & \tilde{h}_2 & 0 & 0 & 0 & 0 & 0 & 0 & 0 & 0 & 0 & \tilde{h}_2 & \tilde{h}_1 \\ \tilde{h}_2 & \tilde{h}_1 & \tilde{h}_0 & \tilde{h}_1 & \tilde{h}_2 & 0 & 0 & 0 & 0 & 0 & 0 & 0 & 0 & 0 \\ 0 & 0 & \tilde{h}_2 & \tilde{h}_1 & \tilde{h}_0 & \tilde{h}_1 & \tilde{h}_2 & 0 & 0 & 0 & 0 & 0 & 0 & 0 \\ 0 & 0 & 0 & 0 & \tilde{h}_2 & \tilde{h}_1 & \tilde{h}_0 & \tilde{h}_1 & \tilde{h}_2 & 0 & 0 & 0 & 0 & 0 \\ 0 & 0 & 0 & 0 & 0 & 0 & \tilde{h}_2 & \tilde{h}_1 & \tilde{h}_0 & \tilde{h}_1 & \tilde{h}_2 & 0 & 0 & 0 \\ 0 & 0 & 0 & 0 & 0 & 0 & 0 & 0 & \tilde{h}_2 & \tilde{h}_1 & \tilde{h}_0 & \tilde{h}_1 & \tilde{h}_2 & 0 \\ \tilde{h}_2 & 0 & 0 & 0 & 0 & 0 & 0 & 0 & 0 & 0 & \tilde{h}_2 & \tilde{h}_1 & \tilde{h}_0 & \tilde{h}_1 \\ h_1 & -h_0 & h_1 & 0 & 0 & 0 & 0 & 0 & 0 & 0 & 0 & 0 & 0 & 0 \\ 0 & 0 & h_1 & -h_0 & h_1 & 0 & 0 & 0 & 0 & 0 & 0 & 0 & 0 & 0 \\ 0 & 0 & 0 & 0 & h_1 & -h_0 & h_1 & 0 & 0 & 0 & 0 & 0 & 0 & 0 \\ 0 & 0 & 0 & 0 & 0 & 0 & h_1 & -h_0 & h_1 & 0 & 0 & 0 & 0 & 0 \\ 0 & 0 & 0 & 0 & 0 & 0 & 0 & 0 & h_1 & -h_0 & h_1 & 0 & 0 & 0 \\ 0 & 0 & 0 & 0 & 0 & 0 & 0 & 0 & 0 & 0 & h_1 & -h_0 & h_1 & 0 \\ h_1 & 0 & 0 & 0 & 0 & 0 & 0 & 0 & 0 & 0 & 0 & 0 & h_1 & -h_0 }\right] \cdot \left[ \matrix{ v_1 \\ v_2 \\ v_3 \\ v_4 \\ v_5 \\ v_6 \\ v_7 \\ v_8 \\ v_7 \\ v_6 \\ v_5 \\ v_4 \\ v_3 \\ v_2 }\right] =

\left[\matrix{ \tilde{h}_2 v_3 + \tilde{h}_1 v_2 + \tilde{h}_0 v_1 + \tilde{h}_1 v_2 + \tilde{h}_2 v_3 \\ \tilde{h}_2 v_1 + \tilde{h}_1 v_2 + \tilde{h}_0 v_3 + \tilde{h}_1 v_4 + \tilde{h}_2 v_5 \\ \tilde{h}_2 v_3 + \tilde{h}_1 v_4 + \tilde{h}_0 v_5 + \tilde{h}_1 v_6 + \tilde{h}_2 v_7 \\ \tilde{h}_2 v_5 + \tilde{h}_1 v_6 + \tilde{h}_0 v_7 + \tilde{h}_1 v_8 + \tilde{h}_2 v_7 \\ \tilde{h}_2 v_7 + \tilde{h}_1 v_8 + \tilde{h}_0 v_7 + \tilde{h}_1 v_6 + \tilde{h}_2 v_5 \\ \tilde{h}_2 v_7 + \tilde{h}_1 v_6 + \tilde{h}_0 v_5 + \tilde{h}_1 v_4 + \tilde{h}_2 v_3 \\ \tilde{h}_2 v_5 + \tilde{h}_1 v_4 + \tilde{h}_0 v_3 + \tilde{h}_1 v_2 + \tilde{h}_2 v_1 \\ h_1 v_1 - h_0 v_2 + h_1 v_3 \\ h_1 v_3 - h_0 v_4 + h_1 v_5 \\ h_1 v_5 - h_0 v_6 + h_1 v_7 \\ h_1 v_7 - h_0 v_8 + h_1 v_7 \\ h_1 v_7 - h_0 v_6 + h_1 v_5 \\ h_1 v_5 - h_0 v_4 + h_1 v_3 \\ h_1 v_3 - h_0 v_2 + h_1 v_1 }\right]

Take a close look at the out of the transformation above and compare to the output of the original length 8 vector. Let's start with the lowpass portions. For the periodized vector, there are seven output values and elements 2, 3, and 4 equal elements 7, 6, and 5, respectively. Moreover, elements 2 and 3 are exactly the same as elements 2 and 3 in the original transform AND elements 1 and 4 use consecutive values of the original input vector! So we will take as our new lowpass portion of the length 8 transform the first four values of the output of the length 14 transformation. A similar observation can be made for the highpass portion of the length 14 transformation - elements 1, 2, and 3 of the highpass portion equal elements 8, 7, and 6, respectively and the fourth element uses consecutive values of the original input vector. Thus our length 8 modified wavelet transformation is

\left[ \matrix{ v_1 \\ v_2 \\ v_3 \\ v_4 \\ v_5 \\ v_6 \\ v_7 \\ v_8 }\right] \rightarrow \left[ \matrix{ \tilde{h}_2 v_3 + \tilde{h}_1 v_2 + \tilde{h}_0 v_1 + \tilde{h}_1 v_2 + \tilde{h}_2 v_3 \\ \tilde{h}_2 v_1 + \tilde{h}_1 v_2 + \tilde{h}_0 v_3 + \tilde{h}_1 v_4 + \tilde{h}_2 v_5 \\ \tilde{h}_2 v_3 + \tilde{h}_1 v_4 + \tilde{h}_0 v_5 + \tilde{h}_1 v_6 + \tilde{h}_2 v_7 \\ \tilde{h}_2 v_5 + \tilde{h}_1 v_6 + \tilde{h}_0 v_7 + \tilde{h}_1 v_8 + \tilde{h}_2 v_7 \\ h_1 v_1 - h_0 v_2 + h_1 v_3 \\ h_1 v_3 - h_0 v_4 + h_1 v_5 \\ h_1 v_5 - h_0 v_6 + h_1 v_7 \\ h_1 v_7 - h_0 v_8 + h_1 v_7 }\right]

Certainly, we can write down a matrix for this transformation. Rows 2, 3, 5, 6, and 7 look exactly like the original matrix. Only rows 1, 4, and 8 change:

\left[\matrix{ \tilde{h}_0 & 2\tilde{h}_1 & 2\tilde{h}_2 & 0 & 0 & 0 & 0 & 0 \\ \tilde{h}_2 & \tilde{h}_1 & \tilde{h}_0 & \tilde{h}_1 & \tilde{h}_2 & 0 & 0 & 0 \\ 0 & 0 & \tilde{h}_2 & \tilde{h}_1 & \tilde{h}_0 & \tilde{h}_1 & \tilde{h}_2 & 0 \\ 0 & 0 & 0 & 0 & \tilde{h}_2 & \tilde{h}_1 & \tilde{h}_0 + \tilde{h}_2 & \tilde{h}_1 \\ h_1 & -h_0 & h_1 & 0 & 0 & 0 & 0 & 0 \\ 0 & 0 & h_1 & -h_0 & h_1 & 0 & 0 & 0 \\ 0 & 0 & 0 & 0 & h_1 & -h_0 & h_1 & 0 \\ 0 & 0 & 0 & 0 & 0 & 0 & 2h_1 & -h_0 }\right]

It is natural to ask whether or not the above matrix is invertible. The answer is yes! To see the effect of the modified transformation, consider the following vector below. The values are v_k = k, k=1,\ldots,32. The second plot is the regular biorthogonal wavelet transformation using the 5-3 filter pair described in the General Wavelet Transformations subsection. You can see the wrapping effects in the first and last element of the lowpass portion and the last element of the highpass portion. The last plot is the modified transformation. Note that the wrapping effect is much less apparent.


		Original input vector.


		Biorthogonal wavelet transformation.			Modified biorthogonal wavelet transformation.

The same periodization trick works for odd-length filter pairs and for any even length input vectors. A similar periodization exists for even-length filter pairs. Please see Discrete Wavelet Transformations: An Elementary Approach with Applications for more details.

The Cohen-Daubechies-Feauveau Filter

A biorthogonal filter pair of lengths 9 and 7 can be constructed in the same manner described in the General Wavelet Transformations subsection - for the length 7 filter, we start with row 6 of Pascal's triangle, divide each term by 64 and then multiply each term by \sqrt{2}. The result is

\tilde{{\bf h}} = (\sqrt{2}/64, 3\sqrt{2}/32, 15\sqrt{2}/64, \sqrt{2}/3, 15\sqrt{2}/64, 3\sqrt{2}/32, 1\sqrt{2}/64)

The length 9 filter is

{\bf h} = (-5\sqrt{2}/64, 15\sqrt{2}/64, -7\sqrt{2}/8, -7\sqrt{2}/32, 77\sqrt{2}/32, -7\sqrt{2}/32, -7\sqrt{2}/8, 15\sqrt{2}/64, -5\sqrt{2}/64 )

Let's plot the absolute value of the Fourier series associated with both filters:


		\|\, \tilde{H}(\omega) \, \|.			\|\, H(\omega) \, \|.

From the plots of the absolute values of the Fourier series, you can see the problem with this filter pair. The filter \tilde{{\bf h}} is extremely flat at \omega = \pi while the plot of |\, H(\omega)\, | is not as flat at \omega = \pi but all has a maximum near the middle of the interval. To serve as a lowpass filter, we desire a graph more like that for |\, \tilde{H}(\omega)\, |. For this reason, this 9/7 pair is not desirable.

Daubechies, Cohen, and Feauveau came up with a method to balance the number of zeros in the Fourier series at \omega = \pi for the filters. We do not discuss their approach here - please see Discrete Wavelet Transformations: An Elementary Approach with Applications for more details. The filters are

\tilde{{\bf h}} = ( 0.0378285,-0.0238495,-0.110624,0.377403,0.852699,0.377403,-0.110624,-0.0238495,0.0378285)

and

{\bf h } = ( -0.0645389,-0.0406894,0.418092,0.788486,0.418092,-0.0406894,-0.0645389)

The plots of the absolute values of their Fourier series appear below:


		\|\, \tilde{H}(\omega) \, \|.			\|\, H(\omega) \, \|.

The LeGall Filter

The LeGall filter pair is actually a modified version of the biorthogonal filter pair developed in the General Wavelet Transformations subsection. If you recall, we had \tilde{{\bf h}} = (\tilde{h}_1,\tilde{h}_0,\tilde{h}_1) = (\sqrt{2}/4,\sqrt{2}/2,\sqrt{2}/4) and {\bf h} = (h_2, h_1, h_0, h_1, h_2) = (-\sqrt{2}/8,\sqrt{2}/4,3\sqrt{2}/4,\sqrt{2}/4,-\sqrt{2}/8). The LeGall filter pair is obtained by multiplying \tilde{{\bf h}} by \sqrt{2}, dividing the second filter by \sqrt{2} and then switching the filters so that the longer filter is the lowpass filter. We have

\tilde{{\bf h}} = (h_2,h_1,h_0,h_1,h_2) = (-1/8, 1/4, 3/4, 1/4, -1/8)

{\bf h} = (h_1, h_0, h_1) = (1/2,1,1/2)

We construct \tilde{{\bf g}} from {\bf h} using the formula g_k = (-1)^{1-k} h_{1-k} so that \tilde{{\bf g}} = (-1/2, 1, -1/2).

The filters lengths here are shorter - this will reduce computation time of the transformation. Here is the modified wavelet transformation (see above) for a length 8 vector:

\left[ \matrix{ 3/4 & 1/2 & -1/4 & 0 & 0 & 0 & 0 & 0 \\ -1/8 & 1/4 & 3/4 & 1/4 & -1/8 & 0 & 0 & 0 \\ 0 & 0 & -1/8 & 1/4 & 3/4 & 1/4 & -1/8 & 0 \\ 0 & 0 & 0 & 0 & -1/8 & 1/4 & 5/8 & 1/4 \\ -1/2 & 1 & -1/2 & 0 & 0 & 0 & 0 & 0 \\ 0 & -1/2 & 1 & -1/2 & 0 & 0 & 0 & 0 \\ 0 & 0 & 0 & -1/2 & 1 & -1/2 & 0 & 0 \\ 0 & 0 & 0 & 0 & 0 & 0 & -1 & 1 } \right] \cdot \left[\matrix{ v_1 \\ v_2 \\ v_3 \\ v_4 \\ v_5 \\ v_6 \\ v_7 \\ v_8 }\right] = \left[ \matrix{ -v_3/8 + v_2/4 + 3v_1/4 + v_2/4 - v_3/8 \\ -v_1/8 + v_2/4 + 3v_3/4 + v_4/4 - v_5/8 \\ -v_3/8 + v_4/4 + 3v_5/4 + v_6/4 - v_7/8 \\ -v_5/8 + v_6/4 + 3v_7/4 + v_8/4 - v_7/8 \\ -v_1/2 + v_2 - v_3/2 \\ -v_3/2 + v_4 - v_5/2 \\ -v_5/2 + v_6 - v_7/2 \\ -v_7/2 + v_8 - v_7/2 }\right]

Wim Sweldens developed a method for quickly computing the above transformation. His method is called lifting and it involves computing the highpass portion of the transformation and then lifting the lowpass portion from the highpass. The method is also useful since we can modify it to map integers to integers AND invert the process!

First, let's create two new vectors {\bf e}, {\bf o} from our original input. We have

{\bf o} = (o_1, o_2, o_3, o_4) = (v_1, v_3, v_5, v_7)

{\bf e} = (e_1, e_2, e_3, e_4) = (v_2, v_4, v_6, v_8)

If we call the elements in the highpass portion d_k, then d_k = e_k - (o_k + o_{k+1})/2 where k=1,2,3,4 and o_5 = o_4. We now lift the lowpass elements s_k from the highpass elements.

Look what happens when we add d_1 and d_1:

\begin{eqnarray} d_1 + d_2 &=& e_1 - (o_1 + o_2)/2 + e_2 - (o_2 + o_3)/2 \\ &=& v_2 - (v_1 + v_3)/2 + v_4 - (v_3 + v_5)/2 \\ &=& -v_1/2 + v_2 - v_3 + v_4 - v_5/2 \end{eqnarray}

This sum is very close to s_2 = -v_1/8 + v_2/4 + 3v_3/4 + v_4/4 - v_5/8. Indeed if we divide by 4 and rearrange some terms we obtain

\begin{eqnarray} (d_1 + d_2)/4 &=& -v_1/8 + v_2/4 - v_3/4 + v_4/4 - v_5/8 \\ &=& \left( -v_1/8 + v_2/4 + 3v_3/4 + v_4/4 - v_5/8\right) - v_3 \\ &=& s_2 - v_3 \\ &=& s_2 - o_2 \end{eqnarray}

We can solve for s_2 and obtain s_2 = o_2 + (d_1 + d_2)/4. In general, we have s_k = o_k + (d_{k-1}+d_k)/4, k=1,2,3,4 where d_0 = d_1. The method for computing the lowpass values is faster than the inner product.

We can also answer the inverse question. Given d_k, s_k, we can compute o_k = s_k - (d_{k-1}+d_k)/4 and once we have the o_k, we can obtain e_k via e_k = d_k + (o_k + o_{k+1})/2.

JPEG2000 maps integers to integers using a slight modification of this transformation:

\begin{eqnarray} d_k &=& e_k - \lfloor\, (o_k+o_{k+1})/2 \, \rfloor \\ s_k &=& o_k + \lfloor\, (d_{k-1}+d_k)/4 + 1/2 \, \rfloor \end{eqnarray}

The 1/2 that appears in the second formula prevents bias in the flooring process. It is not entirely clear, but the method is completely invertible! This seems counterintuitive given we have transformed data and then floored it (a nonlinear process). But the inverse is unique - please see Discrete Wavelet Transformations: An Elementary Approach with Applications for more details.