Eliminate Signal Noise With Discrete Wavelet Transformation
The wavelet transform is a mathematical tool that's becoming quite useful for analyzing many types of signals. It has been proven especially useful in data compression, as well as in adaptive equalizer and transmultiplexer applications.
A wavelet is a small, localized wave of a particular shape and finite duration. Several families, or collections of similar types of wavelets, are in use today. A few go by the names of Haar, Daubechies, and Biorthogonal. Wavelets within each of these families share common properties. For instance, the Biorthogonal wavelet family exhibits linear phase, which is an important characteristic for signal and image reconstruction.
Wavelet analysis is simply the process of decomposing a signal into shifted and scaled versions of a particular wavelet. An important property of wavelet analysis is perfect reconstruction, which is the process of reassembling a decomposed signal or image into its original form without loss of information. By examining wavelet theory as it applies to three specific applications, we find that it works so well because these examples rely on perfect reconstruction for their fundamental operation.
There are no set rules for the choice of the mother wavelet used in wavelet analysis. The choice depends on the properties of the mother wavelet, the properties of the signal to be examined, and the requirements of the analysis. For this reason, it's convenient to have tools that let you easily explore and experiment with many different wavelets and input signals. The following examples use MATLAB, the Wavelet Toolbox, and Simulink to make exploration of wavelet concepts convenient.
In this article, the wavelet we use as an example (called the "mother" wavelet) is the Daubechies wavelet, db4. The 4 in the name represents the order of the filter, which corresponds to eight coefficients.
The Discrete Wavelet Transform (DWT) is commonly employed using dyadic multirate filter banks, which are sets of filters that divide a signal frequency band into subbands. These filter banks are comprised of low-pass, high-pass, or bandpass filters. If the filter banks are wavelet filter banks that consist of special low-pass and high-pass wavelet filters, then the outputs of the low-pass filter are the approximation coefficients. Also, the outputs of the high-pass filter are the detail coefficients.
The process of obtaining the approximation and detail coefficients is called decomposition. Termed multilevel decomposition, this process can be repeated, with successive approximations (the output of the low-pass filter in the first bank) being decomposed in turn, so that one signal is broken down into a number of components.
A two-level decomposition is shown in Figure 1. In this illustration, a2 represents the approximation coefficients, while d2 and d1 represent the detail coefficients resulting from the two-level decomposition. After each decomposition, we employ decimation by two to remove every other sample and, therefore, reduce the amount of data present.
The Inverse Discrete Wavelet Transform (IDWT) reconstructs a signal from the approximation and detail coefficients derived from decomposition. The IDWT differs from the DWT in that it requires upsampling and filtering, in that order. Upsampling, also known as interpolating, means the insertion of zeros between samples in a signal. The right side of the figure shows an example of reconstruction.
Another way to interpret the figure is that the analysis filter bank on the left reduces the rate of an input signal and produces multiple output signals with varying rates. The analysis filter bank performs the DWT represented by the decomposition. The synthesis filter bank on the right increases the rates of multiple input signals while combining them into a single output signal. It performs the IDWT represented by the reconstruction.
The Filters Are The Key Now one might ask, what's unique about wavelet filter banks? The magic is in the filters themselves. By choosing filters that are intimately related for both decomposition and reconstruction processes, the effects of aliasing, which can be introduced by the decimation, are removed.When the signal is reconstructed, it doesn't exhibit any aliasing or distortion (right side of Fig. 1). As a result, the output is said to be a perfect reconstruction.
Wavelet filters have finite length. They aren't truncated versions of infinitely long filter re-sponses. Because of this property, wavelet filter banks can perform local analysis, or the examination of a localized area of a larger signal. Local analysis is an important consideration when dealing with signals that have discontinuities. Wavelet transforms can be applied to these kinds of signals with excellent results. This is due to their ability to locate short-time (local) high-frequency features of a signal and resolve low-frequency behavior at the same time.
As stated earlier, perfect reconstruction is an important property of wavelet filter banks. When the analysis filter bank output is connected to the synthesis filter bank input and the proper delays for alignment are used, as in Figure 1, then the output of the entire system is identical to the input. If a threshold operation is applied to the output of the DWT and wavelet coefficients that are below a specified value are removed, then the system will perform a "de-noising" function.
Two different threshold operations can be viewed in Figure 2. In the first, hard thresholding, coefficients whose absolute values are lower than the threshold are set to zero. Hard thresholding is extended by the second technique, soft thresholding, by shrinking the remaining nonzero coefficients toward zero.
The de-noising process consists of decomposing the original signal, thresholding the detail coefficients, and reconstructing the signal. The decomposition portion of our de-noising example is accomplished via the DWT. The Wavelet Toolbox provides various parameters from which one must pick in order to decompose the signal. These parameters include loading the original signal, selecting the wavelet family, and specifying the level of decomposition.
We have picked Daubechies 4 (db4) as our analysis wavelet, a three-level decomposition. We could have elected to perform more levels of decomposition, as the more levels we chose to decompose our signal, the more detail coefficients we get. But for de-noising our signal in this example, a three-level decomposition provides sufficient noise reduction.
By employing the Wavelet 1-D Discrete Wavelet Analysis Tool from the Wavelet Toolbox, one can calculate the DWT by clicking on the Analyze button. The results of the decomposition are illustrated by Figure 3. The original noisy block signal is shown in the s trace. The a3 trace represents the third-level approximation coefficients, which are the high-scale, low-frequency components. Note how the approximation a3 is similar to the original signal s. The other three waveforms (d3, d2, and d1) are the detail coefficients, which are low-scale, high-frequency components.
After the decomposition was complete, we opened the Wavelet 1-D De-Noising Tool so that we could define our de-noising parameters (Fig. 4) We chose to implement the defaults of the Fixed Form Soft Threshold and Unscaled White Noise for our soft thresholding method. The threshold values were fixed at exactly 4.000 for each decomposition level, which gave us good noise suppression. Both the de-noised and original signals are shown in Figure 5.
Simulation Of Real-Time De-Noising The previous example explored the theory of wavelet decomposition and thresholding in a data-flow environment. While this theoretical exploration is an important part of the de-noising implementation process, the implementation of a real-time system has other aspects that must be considered as well. This is best performed with a simulation tool like Simulink, which allows the use of individual subsytem components in a multirate time-flow environment.The first and most familiar technique for implementing multirate DSP systems is to propagate scalar data samples with varying rates within a system model. A Simulink model containing analysis and synthesis filter banks propagating scalar data samples is shown in Figure 6. This example uses the same input signal s from Figure 3.
For simulation in a data-flow environment, such as MATLAB, processing signals of differing lengths due to changes in sample rates needs to be addressed. The data alignment would be carried out using the appropriate indexing schemes. In a multirate time-flow environment, like Simulink, the data is conceptually infinite in length, and indexing cannot perform the data alignment. Instead, delay elements are introduced after the analysis filter bank to achieve data alignment. This is accomplished by the "Delay Alignment" subsystem (Fig. 6, again).
Furthermore, to compare the output signal with the input, additional delays are introduced into the input signal path. Data alignment is a significant aspect of a practical, real-time implementation. The input, output, and residual signals shown in Figure 6 can be viewed in the scope display in Figure 7.
The wavelet transmultiplexer (WTM) provides an interesting example of the perfect reconstruction property of the DWT. The transmultiplexer combines two source signals for transmission over a single link, then separates the two signals at the receiving end of the channel (Fig. 8). The inputs are assumed to be baseband signals.
The ability of wavelets to provide perfect reconstruction of independent signals, transmitted over a single communications link, is demonstrated in Figure 9. Channels 1 and 2 are perfectly recreated, as indicated by the combined error plot. The error trace is plotted with an expanded vertical scale to demonstrate the absence of signal corruption.
The model shown in Figure 8 demonstrates a two-channel transmultiplexer. But the method can be extended to an arbitrary number of channels. Note that the total data rate is still limited by the Nyquist rate of the high-speed data link.
Similarities With FDM Operation The operation of a WTM is analogous to a frequency-domain multiplexer (FDM) in several respects. In an FDM, baseband input signals are filtered and modulated into adjacent frequency bands, summed together, and then transmitted over a single link. On the receiving end, the transmitted signal is filtered to separate the two adjacent frequency channels. The signals are then demodulated back to baseband.The filters need to pass the desired signal through the filter passband with as little distortion as possible. In addition, the filters must strongly attenuate the adjacent signal to provide a sharp transition from the filter passband to its stopband. This process limits the amount of crosstalk, or signal leakage, from one frequency band to the next. These constraints generally require longer and more expensive filters.
Often, FDM employs an unused frequency band, known as a guard band, between the two modulated frequency bands to relax the requirements on the FDM filters. This decreases spectral efficiency, thereby reducing the usable bandwidth for each input signal.
In a WTM, the filtering performed by the synthesis and analysis wavelet filters is analogous to the filtering steps in the FDM. Plus, the interpolation in the IDWT is equivalent to frequency modulation. From a frequency-domain perspective, the wavelet filters are fairly poor spectral filters, exhibiting slow transitions from passband to stopband, and providing significant distortion in their response.
What makes the WTM special, though, is that the analysis and synthesis filters together completely cancel the filter distortions and signal aliasing. That produces perfect reconstruction of the input signals and, thus, perfect extraction of the multiplexed inputs.
Ideal spectral efficiency can be achieved with the WTM, because no guard band is required. Practical limitations of implementing the channel filter create out-of-band leakage and distortion. In the conventional FDM approach, every channel within the same communications system requires its own filter and is susceptible to crosstalk from neighboring channels. But the WTM method only requires a single bandpass filter for the entire communications channel, and the channel-to-channel interference is eliminated.
Keep in mind that a noisy link can cause imperfect reconstruction of the input signals. Furthermore, the effects of channel noise and other impairments on the recovered signals can differ in FDM- and WTM-based systems.
Image compression is becoming increasingly important as the efficient use of available transmission bandwidth becomes more complex. As complexity increases, system resources must be optimized to use minimal bandwidth and memory. One way to optimize these resources is to employ image compression. The method and amount of compression needs to be such that it's still possible to achieve a reasonable reconstruction of the image. Wavelet transforms have this capability.
The compression procedure is similar to that of de-noising used in an earlier example. The only difference lies in the thresholding applied to the detail coefficients. Two ap-proaches are available in the Wavelet Toolbox for thresholding detail coefficients when compressing two-dimensional data. These are global thresholding and level thresholding.
With the global-thresholding ap-proach, we define a global-threshold method, a compression performance factor, and a relative square norm recovery performance factor. The Wavelet Toolbox derives a global threshold from an equal balance between the percentages of retained energy and number of zero coefficients. With the level-thresholding approach, one would need to visually determine level-dependent thresholds.
In this example, we allow the Wavelet Toolbox to derive a global threshold for our example image. The image shown in Figure 10 was decomposed using the two-dimensional discrete wavelet analysis tool (similar to the one-dimensional tool found in Figure 3). For this example, we decided to perform a two-level decomposition using the biorthogonal spline wavelet bior3.7, which specifies a third-order reconstruction filter and a seventh-order decomposition filter.
The compression tools available in the Wavelet Toolbox perform only the thresholding portion of the compression process. Its performance is measured by the percentage of remaining nonzero elements in the wavelet decomposition. When implementing a real-world compression scheme, one would need to further consider quantization and bit-allocation factors.
The two-dimensional wavelet compression tool automatically generates a threshold based on the thresholding method selected (Fig. 10, again). We picked "Remove near 0," which sets this global threshold to 4. When we click on the Compress button, all coefficients whose values are less than 4 (in this case, 49.81%) are forced to zero. In spite of this case, 98.98% of the original image energy is retained. See the Wavelet Toolbox User's Guide for more information on how these percentages were calculated.
Wavelet analysis is a new and promising tool which complements traditional signal processing techniques. It can offer significant advantages for real-time systems, and it opens the door to new and exciting communications applications.
For Further Information:
- C. Taswell, "The What, How, and Why of Wavelet Shrinkage Denoising," Computing In Science And Engineering, vol. 2, no. 3, May/June 2000, p. 12-19.
- Department of Applied Science Wavelet Group, Lawrence Livermore National Laboratory, www.llnl.gov/das/wavelet/wavelet.html.
- G. Cherubini; J. Cioffi; E. Eleftheriou; S. Olcer; "Filter Bank Modulation Techniques for Very High-Speed Digital Subscriber Lines," IEEE Communications Magazine, vol. 38, no. 5, May 2000, p. 98-104.
- G. Strang and T. Nguyen, Wavelets And Filter Banks, Wellesley-Cambridge Press, 1997.
- M. Misiti; Y. Misiti; G. Oppenheim; J.M. Poggi, Wavelet Toolbox User's Guide, The MathWorks Inc., 1996.
Nickolas Adler is a technical wordforger in the documentation department at The MathWorks Inc., Natick, Mass. He possesses a BS in Physics from Tufts University in Medford, Mass.
Paul Costa is a DSP engineer at The MathWorks Inc. He possesses a BS in electrical engineering, which he received in 1994 from the University of Massachusetts Dartmouth, North Dartmouth, Mass. Currently, he's pursuing an MS in electrical engineering at the University of Massachusetts Dartmouth. Costa may be reached by phone at (508) 647-7349, or via e-mail at [email protected].
Andrew V. Dowd is senior DSP engineer at The MathWorks Inc. He holds a BS in electrical engineering from the University of Virginia in Charlottesville, as well as an MS in electrical engineering from the University of Arizona at Tucson. Dowd may be reached by phone at (508) 647-7916, or via e-mail at [email protected].
Anne Mascarin is the DSP market segment manager at The MathWorks Inc. She holds a BA from Boston University and an MS in electrical engineering from Northeastern University, Boston, Mass. Mascarin may be contacted by phone at (508) 647-7598, or through e-mail at [email protected].
Don Orofino is manager of DSP development at The MathWorks Inc. He earned a PhD in electrical engineering in 1992 from Worcester Polytechnic Institute, Worcester, Mass., as well as a BS/MS in electrical engineering in 1988 from Union College, Schenectady, N.Y. Orofino can be reached by phone at (508) 647-7000, or via e-mail at [email protected].