INTRODUCTION
Wavelet Transform has been successfully applied to different fields, ranging from pure mathematics to applied science (Nibouche et al., 2000). Numerous studies, carried out on Wavelet Transform, have proven their advantages in image processing and data compression and have made it a basic encoding technique in recent data compression standards. Purely software implementations of the Discrete Wavelet Transform, however, appear to present a problem when real time systems are required in terms of performance. Therefore, hardware acceleration of the DWT has become a topic of recent research.
The basic motivation of our application to radio isotopic images and neutron
tomography projections, (Fig. 1) is to reduce the data volume
and to achieve a low bit rate in the digital representation images without apparent
loss of image quality, in order to reach good 3D reconstruction on a distant
Computer or storage for future retrieval.
Radiological image is a digital 2D image acquired by a certain radiological
procedure which can be: Xrays, neutron radiography, gamma camera image or a
nuclear magnetic reso nance image. For tomography, the gray level value represents
the relative linear attenuation coefficient of the object.

Fig. 1: 
Neutron
imaging system (Kharfi et al., 2005) 
Higher quality images
show finer structural or functional information of body organs and support more
reliable diagnostic outcomes (Wong et al., 1995).
Present approach for image compression helps a lot for 2D image storage and
transmission. It is based on projection images collection, encoding/decoding,
TCP/IP transmission and therefore 3D reconstruction. The design is based on
the MALLAT’s Fast Wavelet Transform Algorithm, which is a fast implementation
of the Discrete Wavelet Transform. The design utilizes various digital techniques
and specific features of the Xilinx VirtexIIV2MB1000 FPGA (Appendix A and B)
to accelerate the transform.
Characteristics of digital radiological images: Digital radiological image is a digital image acquired by a certain radiological procedure which can be Xrays, neutron radiography, gamma camera image or a nuclear magnetic resonance image. It is a twodimensional MxN array of nonnegative integers (gray levels) (Fig. 1).
The performance of radiological image compression depends not only on the compression ratio but also the image quality of reconstructed images. Higher quality images show finer structural or functional information of body organs and support more reliable diagnostic outcomes.
Image compression framework: The general framework for radiological image compression based on wavelet transform (Fig. 2) is similar to other digital compression fields, the framework includes three major stages: Image transformation, quantization (irreversible compression only) and entropy encoding. The relative importance of each stage varies from one technique to another. All reversible compression techniques do not involve the stage of quantization (Wong et al., 1995). It is useful to consider entropy encoding as a twostep process. The first step uses a statistical model to convert coefficients into an intermediate sequence of symbols. The second step converts the symbol to a data stream in which the symbols no longer have externally identifiable boundaries. The entropy coding uses a variablelength code to achieve higher compression rates, such as Huffman coding. The same code tables used to compress an image are needed to decompress it. A promising approach of image compression is progressive transmission that transmits image data in stages and at each stage, an approximation to the original image is reconstructed at the receiver.
A radiological image compression algorithm that provides multiresolution image representation offers many advantages.
Wavelet transform theory: Formally, the wavelet transform is defined by many authors as a mathematical technique in which a particular signal is analyzed (or synthesized) in the time domain by using different versions of a dilated (or contracted) and translated (or shifted) basis function called the wavelet prototype or the mother wavelet.
The CWT as defined by Morlet and Grossman is (Nibouche et al., 2002):

Fig. 2: 
Waveletbased
encoding scheme (Nibouche et al., 2002) 
Where f(t) belongs to the square integrable functions space, L^{2}(R). In the same way, the inverse CWT can be defined as:
Where a and b are, the scaling and the shifting factors and to ease the implementation
of a wavelet system, these factors have been adopted to be a factor of two.
i.e., (Nibouche et al., 2002):
The reconstruction is only possible if the factor C_{ψ} is defined
(<∞), thus ψ(t) (mother wavelet) must have zero mean and finite
energy (admissibility condition). The digital realization of CWT is DWT, given
as the inner product of the signal x(t) being transformed with each of the discrete
orthogonal/orthonormal set of basis wavelets.
Based on Mallat’s work on the relationships between the Quadrature Mirror Filters (QMF), pyramid algorithms and orthonormal wavelet bases, Daubechies succeeded to construct a set of wavelet orthonormal basis functions which have become the cornerstone of many applications such as JPEG2000 compression standard (Mallat, 1989).
Multiresolution, filter banks and DWT implementation: This concept has been introduced first by Mallat (1989) and then Nibouche (2002). It defines clearly the relationships between the QMF, pyramid algorithms and orthonormal wavelet bases. The strength of multiresolution lies in its ability to decompose a signal in finer and finer details. Most importantly, it allows the description of a signal in terms of timefrequency or timescale analysis.
Scaling function: Let the scaling function be defined as:
φ_{k}(t) = φ(tk) K∈Z φ∈L^{2}(R) 
By scaling and translating, a 2dimensional scaling function is generated:
To ease the implementation:
Thus: φ_{j,k}(t) = 2^{j/2} φ(2^{j} tk), a = 2^{j}, b = 2^{j}.k
By spanning the whole space:
Where h(k) are the scaling function coefficients. Equation 3 is called the multiresolution analysis equation
By the same way the wavelet function:
ψ_{k}(t) 
= 
ψ(tk) K∈Z ψ∈L^{2}(R) 
Twodimensional scaled and translated wavelet function is:
ψ_{j,k}(t) = 2^{j/2 }ψ(2^{j}t
– n); j,k∈Ζ such that  ∞<j,k<∞ 
where Z is set of integers.
In terms of the scaling function:
Where g(k) are the wavelet coefficients. Equation 4 leads to a dyadic decomposition.
For practical and computational reasons, discrete time filter banks are required. Such structures decompose a signal into a coarse representation along with added details (Nibouche et al., 2002). We define the relationship between the expansion coefficients at lower and higher scale levels, by using a scaled and shifted version of Eq. 3 and 4:
For finite length signals (which is the case for digital images), the use of
a Finite Impulse Response filter (FIR) is the most appropriate choice. Since
Eq. 5 and 6 compute one output for each
two consecutive inputs, a basic operation required here consists of using a
downsampler or decimator by a factor of two (Fliege, 1994) (Fig.
3).

Fig. 3: 
One
dimensional discrete wavelet transform 

Fig. 4: 
Twodimensional
Mallat’s analysis and synthesis algorithm 
The filter bank, used to implement Mallat’s Algorithm, is defined as a
combination of a low pass filter and high pass filter, both followed by a factor
of two decimation. To allow further level of decomposition, identical stages
are cascaded leading to a multiresolution analysis. Figure 4
shows analysis/synthesis:
For twodimensional (Image) analysis and synthesis. The filter bank used is the (2,2) CohenDebuchiesFeaveu wavelet filter because of its good compression characteristics (Cohen et al., 1992). It is also the fastest and the simplest one to implement in hardware. Polyphase decomposition is adopted to remedy to the problem of resources lost (Fliege, 1994).
The range of coefficients for each subband is divided into 16 quantization levels. The coefficients are quantized into one of the 16 possible levels and those of small magnitude can be neglected without considerable distortion to the image. Huffman trees are an efficient way of coming up with a variable length encoding for a set of characters, given the relative occurrence probabilities.
FPGA implementation of DWT: Configuring FPGA chips is not a straightforward
operation. Normally, it involves a number of different steps (Fig.
5) in which different tools can be used:
• 
Choosing a design capture tool (VHDL or schematic). 
• 
Converting the file into a netlist format. 
• 
Generating the bit stream necessary for configuring the FPGA, the netlist
file is passed through implementation tools (e.g., map, place and route
tools). 
The input image to the encoder is raw gray scale frames of 512 by 512 pixels
(at maximum). Each pixel is represented by 256 gray scale levels (8 bits).

Fig. 5: 
Design
Flow (Nibouche et al., 2002) 
Input
frames are loaded to the 16 bits word memory by a computer and DWT results are
read back, once the FPGAv has processed it. Dynamic vector quantization and
entropy encoding on the wavelet coefficients are then made.
RESULTS
The VHDL implementation (Sanchez, 2003) of DWT/IDWT is simulated on ModelSim PE 6.0d, for correct functionality, before being synthesized using XILINX ISE8.1. Radiological images are converted to memory bit files and vice versa, initial data to be processed is written inside the memory (8 bits for each pixel) and the results (8 bits for integer part of coefficients and 10 bits for fractional part) are stored back in the same memory. Some processed images are shown in Appendix CF. A comparison based on compression ratio (CR), Root mean square error (RMSE) and compression speed is made between VHDL implementations and Matlab7 (Gonzalez, 2004) programs, the CR and the RMSE are varying with image type and size as in Table 1. Performance results for 1levet DWT are shown in ISE reports in Table 2.
Table 1: 
Compression
ratios and RMS errors for different images 

Table 2: 
ISE
synthesis reports 

CONCLUSION
We have implemented a Wavelet transform based image encoder/decoder on reprogrammable hardware FPGA. The system has multiple configurations which support different compression levels, depending on the user input choices. Lessons learned from this study will help us enhance similar implementations in the future for further improvements using the same reconfigurable IC, especially when the whole tomography system will be achieved.
APPENDIX

Fig. A: 
VirtexIIV2MB1000 System Board 

Fig. C: 
Cameraman image and its level 3 DWT 

Fig. D: 
Neutron radiography image and its level 3 DWT 

Fig. E: 
Radio isotopic image and its level 3 DWT 

Fig. F: 
Computed tomography 3D image (a), its level 3 DWT (b) and
the reconstructed image from DWT coefficients (c) 