A Finite Bidimensional Wavelet Framework for Computer Graphics: Image Equalization

A lot of material has been written about wavelet theory. Most of these texts provide an elegant framework from the functional and real analysis point of view. The complete infinite dimensional space (the set of all functions such that) is generally used to develop the theory, but this cannot be directly applied to computer software, because the concept of a non denumerable infinite set of vectors or functions is practically useless here. We provide foundations for a finite, linear-algebra based toolkit of wavelets that supply a rich set of tools that can be used to manage image processing, equalization and compression. We test a frequency criterion to design orthonormal wavelet generators and a multirresolution analysis. We show that this criterion can be easily interpreted graphically. Despite our approach only constructs orthonormal wavelet basis; we believe that this approach is general enough to explore possibilities in other computer graphic fields and solution of integer-differential equations on simple domains. We strongly believe that this approach simplifies considerably the wavelet analysis. The frequency criterion expands possibilities of testing two dimensional wavelet bases according to specific graphical needs, and can be applied to different problems that involve regular grid reduction. We propose this criterion as a fundamental basis to design bidimensional ortonormal wavelets for matrix equalization. It gives a wider range of possibilities than restricting only to some well known bases and gives a direct interpretation for computer images.


Introduction
A complete and consistent wavelet theory can be developed from a linear algebra point of view.This approach is characterized for being simpler and more directly applicable to problems such as image compression or realistic computer rendering.For instance, the Frazier Linear Algebra approach [1] can be considered a more adequate approach from the practical point of view, being more applicable in real finite problems.Using the same functional notation as [1], we have that for any , which is not true in infinite dimensional spaces.The notation (as used by [2]) for a N-dimensional vector space on the field of complex numbers is convenient to define the very important operator of convolution.
All translations of the indicating function ( if and ) constitute a basis for .In this context, a wavelet analysis starts from an orthonormal basis in a finite dimensional vector space over the field of complex or real numbers.A vector of this space is called a signal and when it is expanded in terms of a certain basis; some of its coefficients can be small enough to be dropped without affecting the global behavior of the vector.This is the fundamental concept of compression and image equalization.The way of selecting coefficients can be done with digital filters in a method called multirresolution analysis.We will show that there is no need to implement expensive sorting algorithms to make an adequate choice.
This ability to summarize a large amount of data with only a few coefficients is a powerful and useful characteristic of wavelet basis.It makes this concept suitable for image compression and processing, as it will be shown later.This sampling optimization makes wavelets an ideal tool to be used in computer graphics, especially in heavy sampling algorithms such as Ray Tracing or Radiosity [3].
This paper should give a wider vision of what can be done with wavelets, considering them as only orthonormal bases in a finite dimensional vector space.A good approach to computer graphics should consider a finite element method technique combined with finite dimensional wavelet concepts to solve heavy integral equations based on radiosity.Using a more flexible and intuitive way to create basis can help software to accommodate to different graphical needs.Our ideal filter basis provides a simpler way to define basis in the frequency domain, and also it gets the most relevant characteristics of a two dimensional array giving pleasant results.Most of calculations in the spatial domain are left to the software and inverse Fourier transform.
In Section 2, we introduce Fourier transform and frequency domain concepts from a finite dimensional vector space view point, giving as a framework to introduce a bidimensional criterion for wavelet design in Section 3. In Section 4, we give an overview of multirresolution analysis and subband coding in a finite signal context and we describe how such analysis can be applied to bidimensional image processing.
In Section 5 we describe some simple wavelet basis created from very simple geometrical concepts in the frequency domain.We apply these tools to image processing and show that this approach provides a more general environment to design a wavelet toolkit in software, because there are no limits for the wavelet bases that can be designed.This is the essence of image equalization.
In Section 6, we show that a lazy implementation of several convolution operators is enough to produce pleasant and good enough results for two dimensional compressions.This lazy implementation of wavelet design criteria in the frequency domain can also produce new orthonormal wavelet bases which can give excellent compression rates, better than the most popular Haar's basis.
In section 7 we give an overview of how to combine this discrete wavelet analysis with a finite element method technique to solve the radiosity equations and finally, in section 8, we provide some conclusions.

Convolutions and frequency domain
One of the most useful operators in digital signal processing is convolution.It has been proven that this operator can represent physical light phenomena such as reflection [4], or local diffuse shading [5].
Since we are dealing with finite dimensional signals, a circular convolution is defined as a binary operator * which transforms two vectors into a third vector by the rule [6]: A Fourier basis, defined by the component wise vectors , is a set of eigenvectors of any convolution operator that therefore diagonalizes it.A vector of all eigenvalues of is called a Fourier transform of and is denoted by .These bases can be easily extended to multidimensional arrays, using tensorial products [7].A dimensional convolution is described by the relation:

Where
. A P-dimensional Fourier basis is given by: Fourier bases have a strong physical meaning.Low frequencies (corresponding to values of with near 0) are associated with smooth transitions and high frequencies are associated with abrupt transitions.In a P-dimensional environment, there can be frequency combinations of low and high frequencies in different dimensional projections.In P-dimensional environments, signal equalization becomes a richer concept, where this frequency interpretation satisfies different user interests.A wavelet basis can be designed in this context.We have avoided the use of integrals for defining Fourier transforms and general convolutions.This finite sums context is enough to achieve useful results in Computer Graphics as it will be shown later.
Fourier transforms can also be defined as sums [6] but we have focused on their fundamental property.They diagonalize any convolution operator.Fourier basis has this nice characteristic on the finite dimensional domain.Other orthonormal basis can be built using softer characterizations.

Bidimensional wavelets for image compression
A wavelet basis in the bidimensional space is generated by 4 vectors whose translations are a basis for the entire original vector space [2].Creating such a basis can be facilitated by the next characterization in the frequency domain, in which we extend Frazier results [1] to the bidimensional case:

Theorem 1:
The set of all translations of vectors is an orthonormal set of vectors if and only if the matrix: Is unitary (its columns are a set of orthonormal vectors) for each .
A direct proof of this theorem requires some calculations which we do not provide here.It can be found in [8][9].This theorem is applicable only when and are even numbers.It permits to design an entire wavelet basis in the frequency domain (a filter) by just constructing a matrix whose column vectors are orthonormal.The reason why we need 4 vectors to create an orthonormal basis is that this is the minimum quantity of vectors necessary to achieve certain frequency localization [9].The relation between spatial and frequency localization will always be limited by Heisenberg uncertainty principle [10].

Bidimensional subband coding and multirresolution analysis
We define a downsampling operator as a function such that .Thus, a downsampling operator takes only even components of a vector.Similarly, an upsampling operator enlarges a vector by inserting zeros between its components [11].Let us define for any the vector .
With this notation, we have the reconstruction equation for a filter bank with 4 branches, using a wavelet basis : A wavelet basis generator can be defined in each subspace .We define, new vectors and inductively as and , being . These vectors are generated from smaller versions of themselves.All translations of generate a subspace .We say: We have .This set of subspaces is called a multirresolution analysis.If we define: Then we have .This gives us a representation of a signal in different level of detail, according to a perfect reconstruction bidimensional equation: Applying equation ( 3), we have an iterative filter bank which first decomposes (downsamples) the signal, reducing four times the dimension of the original vector space.It can always be reconstructed with a convolution and an upsampling operator, none of which depends on the original signal.
Let us analyze equation ( 3) according to the subspace representation given above.The convolution is nothing but a projection of vector on space for it represents a linear combination of translations of vector .Adding vectors (linear combinations of orthonormal bases for ) we get a vector which is the projection of on subspace .Iterating this process, we get a decreasing sequence of subspaces at end of which the complete space is obtained.This is the representation of a signal at different levels of detail.Projections on the spaces are called detail branches and presumably they can be suppressed without affecting significantly the global behavior of the signal.This is a conceptual issue that should not be related to the size of the coefficients but the kind of details we consider unnecessary.The heaviest parts of these calculations are convolutions.This is the reason we prefer to use a wavelet basis that differs of 0 in only in a few entries (a finite interpretation of compact support).
In figure 1, we show a basic diagram for a two dimensional finite resolution analysis.Down arrows are downsampling operators and up arrows are upsampling operators.Color boxes represent convolutions.Blue boxes are internal convolutions for the downsampling and orange boxes represent convolutions to recover a downsampled signal.Only one of the four branches should be selected to project on the increasing sequence of vector subspaces of the multirresolution analysis.

Some frequency filters for wavelet basis
We can create a simple wavelet basis by defining a bidimensional matrix whose entries are 1 in some octants of the bidimensional array and 0 everywhere else.The four generators of a wavelet basis are defined as complements of this simple array.Figures 3,4 and 5 give the results of this wavelet basis applied to an image.Its Fourier transform is defined trivially according to the following equations: All other wavelet generators are defined as translations of this simple configuration.At the end, the entire finite dimensional box is covered as shown in Figure 3.This basis should not be interpreted to have compact support in an infinite dimensional space.But this is not an issue in finite dimensional vector spaces as Figure 3 shows.This filter takes all possible combinations of Low, High-Low, High-High frequencies and therefore considers smooth variations along the image and also borders of the objects, as shown in figure 3, right.This equalizing multirresolutive basis provides a representation that constitutes a very accurate summary of the original signal (Figure 4).This is a projection of the original signal on the space .Only 1.5% of total complex coefficients are used to get this.Inverse Fourier transform of these generators are obtained via software (Fast Fourier Transform) and we are not interested in defining this wavelet in spatial domain.To get an idea of how simple is to create a wavelet basis using the Theorem 1 of section 3, we compare these results with the well known Haar's basis.
Haar basis is one of the most widely used in computer graphics because its simplicity [12].It is considered one of the best bases to introduce the topic [13].But according to Heisenberg uncertainty principle, its small spatial localization will produce a smooth signal in frequency domain and compression will be sparse.Figure 5 shows the Haar wavelet basis and its Fourier transform.It takes a sample of all frequency components in the signal.This two dimensional Haar basis is defined as a tensorial product of one dimensional Haar basis: So, , , and .A more traditional approach of low-pass and high-pass can also produce interesting results as shown in Figure 7.

A lazy wavelet implementation
The essence of wavelet analysis for computer graphics is the ability to manipulate the level of detail of a signal and use it in different levels.Here we have used an easy and intuitive method to create wavelets base on their frequency behavior.This leads to a multirresolutive image equalizer.This is a lazy application.We "blindly" apply the perfect reconstruction equation.This means that we are not concerned, during the compression process, with the magnitude of coefficients in the wavelet basis expansion or any other trick to choose the best basis elements.The effectiveness of the method depends on a fast implementation of the convolution operator.A fast wavelet transform is based on the quantity of elements in the vector array of the basis that are different of 0 [10] .Another approach should use Fast Fourier Transform to calculate convolutions [14].
A two dimensional equalization can be used on an image to obtain even more effects, and several wavelet bases can be combined in different levels of the multirresolution analysis as shown in figure 8.These images combine a global approximation in low frequencies and then add some details coefficients to get borders.Theorem (1) should be seen as a design criterion.Despite we are restricting to orthonormal wavelet basis, leaving off a considerable range of possibilities, the concept of creating it from four branches should give us enough freedom and abstraction to combine several basis through a multirresolution analysis and adapt them to a particular application.This theorem could be a basic criterion for orthonormal wavelet design.With this approach it is easier to adapt the theory to algorithmic needs using applied Linear Algebra.
Using elementary operations (such as finite dimensional projections, dot products, change of basis), it is also easy to readapt the theory to other branches of numerical analysis, such as finite element method.This will be studied in the next section.

Other applications: Wavelet finite element method for radiosity
In modeling realistic phenomena, many complicated differential and integral equations emerge.The best way to discretize such kind of problems is to suppose that the approximate solution of the integerdifferential equation is a linear combination of basic functions .That is: Where the are scalars.Basically, these numbers are chosen so that is the orthonormal projection of the real solution over the space .This is the method of weighted residuals [15].Wavelets are considered for the functions .A finite dimensional approach should be considered for this case as well.For example, a method for global illumination called radiosity [16] is particularly representative for the usefulness of this theory [3].As stated by [12], the equilibrium distribution of radiosity satisfies the following equation: The function of Radiosity B(y) is defined as the emitted part of the surface, and is the diffuse reflectance factor.The geometrical component is the relation between the normal vectors of the two surfaces X and Y, and it is defined as follows: To solve this, several methods for approximate a solution have been used.Previous work proposed an unified method to describe this as a rendering equation and several numerical algorithms are described intending to approximate a solution [17].The first step is to discretize the problem.All surfaces are supposed to be divided in patches and there is a bijection between every patch and a rectangular array.This leads us to the cases described above.As more patches are considered, the best approximation finite element methods will give.We can suppose that the correct solution is a linear combination of piecewise constant functions, and this improves according to how many patches are considered.The shapes of these sections of surface correspond to a finite element method issue [18] and we suppose that we are given a two dimensional array of patches that interact with each other according to the discretized equation: Where is a linear operator that is computationally heavy to calculate (it is associated with the radiosity equation kernel [3]).To avoid calculations for so many patches (a good approximation with finite piecewise linear functions depend on it) we can expand the unknown vector as a linear combination of a wavelet basis which produces a multirresolution analysis in the original bidimensional patch space .So, we suppose that: Where is an orthonormal basis for the multirresolution analysis .Substituting this expansion in the discrete equation for radiosity we get: If we make an inner product of this equation with each one of the elements of the basis we can express each one of the coefficients as a linear combination of themselves in which the coefficients are given by internal products .Solving this linear system gives an approximation of .This approximation has been projected, a priori, on the multirresolution analysis subspace .According to the perfect reconstruction equation, basis elements are given by which are an orthogonal basis for .This process can be repeated for the generators of spaces and according to the magnitude of the obtained coefficients, it can be valued if it is necessary to continue refining the process or if there is already a good approximation for the solution .There can be cases in which it may be necessary to refine the solution to the entire space, and sometimes, it will be enough to calculate the approximation in the space with .A good choice of a wavelet basis is important in this issue and some different filters should be tested in this particular problem.So any integer-differential equation can be solved at different levels of detail using multirresolution.Through equalization we can choose what kind of details we want to manipulate.

Conclusions and future work
A strictly finite dimensional approach to wavelet theory and multirresolution analysis methods can be directly applied in a programming language or hardware design [11].Computer graphics algorithms should benefit largely from an optimized fast wavelet transform that implement a finite multirresolution analysis in two or more dimensions.A wavelet toolkit should not be restricted to one kind of wavelet basis such as Haar, and apply theorem of section 3 as a criterion to design different orthonormal wavelet basis.This will expand possibilities to create a real adaptive wavelet filter.We have shown in our experiments that other basis besides Haar are suitable for image compression; despite they are not so lightly supported in the spatial domain.This compression is useful as a resume of a two dimensional array of numbers and not necessarily an optimization of well known compression algorithms.A wavelet filter should adapt to the needs of the application.
A finite dimensional approach can be easily implemented and many optimizations are possible.The possibility of combining this richness of wavelet bases with finite element method can be exploited in a future work.A finite element method, based on multirresolution analysis, provide many exciting possibilities for solving, hierarchically, some differential equations that have a physical interpretation.This not only can help in the compromise between velocity and precision in a Ray Tracing algorithm, but in other physical phenomena as thermodynamics, physical relativity and Maxwell Equations for electromagnetism.The flexibility of this approach depend on the freedom of selecting a wavelet basis, and in the ability to interpret each one of these problems in a discrete form.Finding an adequate basis and a fast implementation for each phase of ray tracing is a possibility that will be explored in the future.
Ray tracing is computationally heavy algorithm that simulates light displacement to the observer eye [19] .A lazy implementation of the method samples every pixel of the screen and displays its corresponding color.In the last years, new hardware development are making feasible to implement this process in real time [20] and a physically based global illumination model emerges as an important tool to add realism [16].This problem can be solved with probabilistic tools [21] or more classical finite element methods.A ray tracing entirely based on finite wavelet signal processing is a very interesting challenge.
Our drawback is that this framework for multirresolutive equalization can be applied only when it is feasible a matrix representation.Second generation wavelets have been defined on more general domains [22] but defining a frequency criterion to construct them is still (as far as we know) an open problem.A very interesting challenge will be to extend theorem 1's concepts and apply them on general domains.

Figure 8 .
Figure 8. Different wavelet bases at different levels of multirresolution.Right figure shows a combination of low-pass high-pass filters and left figure shows a combination of Haar and a band pass filter