Go to:
Title Page               Chapter 1           Appendix A   Appendix J
Copyright                Chapter 2           Appendix B   Appendix K
Abstract                 Chapter 3           Appendix C   Appendix L
Acknowledgements         Chapter 4           Appendix D   References
Table of Contents        Chapter 5           Appendix E
List of Figure           Conclusions/        Appendix F
List of Tables           Future Directions   Appendix G
List of Audio Examples                       Appendix H
List of Programs                             Appendix I


Appendix L

Prototype MATLAB code with examples



This appendix contains an informal report on the denoising of audio using lifted interpolated wavelets, submitted as a final project for Math 123 ("Current problems in analysis," a seminar on wavelet analysis taught by Geoff Davis) during the Winter 1996 term at Dartmouth College. This project provided the tools necessary for much of the work in Chapter 3, and provided a foundation for the analyses and experiments in Chapter 4 and Chapter 5.

The aim of this project was first to prototype the fast interpolated wavelet forward and inverse transforms in the MATLAB environment, and then to translate the resulting working algorithms into a faster C++ environment in preparation for attempted real-time audio signal processing. Real-time processing of audio using the fast forward and inverse lifted wavelet transform was not successful; however, the resulting MATLAB framework that was developed during the process proved invaluable in later experiments. In fact, this framework was used to test ideas on small sets of non-audio data before investing the time required to translate the ideas into a faster C++ environment suitable for CD-quality (44.1 kHz) audio. The MATLAB environment was also used to verify the output of the C++ translations of the transforms.

The text, MATLAB code, and C++ code for this project are included in this appendix. The audio cassette tape and software referenced in this report are not included here, as the examples referenced have been improved and accompanying compact disk.

Also note that the page numbers listed in the table of contents of this project refer to the circled page numbers at the bottom-right-hand corner of each page. The uncircled, centered page numbers refer to the overall page numbers in this thesis.


Corey Cheng

Math 123

3/12/96

Final Project: Denoising of Audio files with Interpolating Lifted Wavelets

Summary

The lifting scheme with interpolating wavelets, as discussed by Wim Sweldens, was used to attempt broad-band denoising of stereo, 16-bit, 44.1 kHz AIFF (Audio Interchange File Format) audio files. Soft thresholding (threshold of 100, suggested by Wicherhauser and Coifman) was used with primal and dual wavelets having 3 vanishing moments each, over windows 32768 audio frames long (1 frame = 1 stereo sample = total of 2 samples/frame). Thirteen (13) iterations of the transform were performed for each window.

Prototyping was done in MATLAB, with a final translation done into a quicker C++ environment. Code was written to perform forward and inverse transforms using both syntax's with the help of a polynomial interpolation algorithm (Neville's algorithm), taken from Press' Numerical Recipes handbook (also translated into both MATLAB and C++ syntax's from the original FORTRAN). Denoising was partially successful in that most high frequency noise (tape hiss) was removed; however, residual, low-frequency noise is prevalent in the output (empirically found to be about 200-400 Hz at a mixing console), and is heard as a rather disruptive, low-frequency rumbling.

References used: Wickerhauser and Coifman, Sweldens, Press, Strang/Nguyen.

Theoretical Notes

Lifting and dual lifting are processes that alter a simpler biorthogonal wavelet transform. Dual lifting makes the dual wavelet more complex, which normal lifting makes the primary wavelet more complex. In making the wavelets more complex, it is hoped that the properties of the wavelets can be improved (increase the number of vanishing moments of the more complex dual/primary wavelet, control its shape for feature detection, etc.).

The purpose of lifting and dual lifting is to transform a simple, biorthogonal wavelet transform into a more complex, biorthogonal transform by the intermediate use of a third filter. Dual lifting and lifting can be performed sequentially, so that first either a dual lifting or lifting can be performed to improve the dual/primary wavelet, and then a lifting or a dual lifting can be performed to improve the primary/dual wavelet. In particular, lifting and dual lifting:

1) alter the original biorthogonal transform's high/low forward/inverse filter coefficients

2) alter the shape of the wavelets and scaling functions (and their duals) associated with the original biorthogonal transform

3) provide a predictable way of deriving the more complex fast wavelet transform (associated with the more complex biorthogonal wavelet) from the simpler fast wavelet transform (associated with the less complex biorthogonal wavelet transform) and the third, intermediate filter.

[Please see Sweldens for details.]

The art of lifting and dual lifting comes from the choice of this third intermediate filter, whose coefficients may vary with scale and translation (unlike traditional filters which are constant throughout each scale and translation).

In the case of interpolating scaling functions, one may interpret the entire forward transform as a lifting of the dual lifting of the lazy wavelet transform, which simply separates the input into even and odd samples; similarly, one may interpret the entire inverse transform as an inverse dual lifting of an inverse lifting of the inverse lazy wavelet transform (which simply interleaves two sets of samples).

*** The use of Neville's algorithm for both the dual and inverse lifting in both the forward and inverse composite transforms arises from the linear proportionality of the primary intermediate filter to dual intermediate filter [Sweldens, "The Lifting Scheme: A custom-design construction of biorthogonal wavelets, Theorem 12, p 18]. In essence, Neville's algorithm chooses the intermediate filter for us, and evaluates the result of the filtering without the user's needing to know the s coefficients exactly. At the boundary, we still use Neville's algorithm to evaluate the wavelet and scaling coefficients, even though interpolation is done in a different way at the boundaries. This can be interpreted as follows: even though the intermediate filter is different at the boundaries, Neville's algorithm chooses and evaluates this filter for us; Sweldens assures us that if the dual and primary lifting is done the same way (i.e. using Neville's algorithm in the same way at the boundaries for both lifting and dual lifting), then we can still achieve perfect reconstruction.

Please see the extensive code comments (especially in the file denoise.cc) for deeper explanations of specific, implementation-oriented issues.

Discussion and Conclusions

The lifting scheme with interpolated wavelets is a powerful technique for wavelet analysis: it allows the custom design of biorthogonal, symmetric wavelets with an arbitrary (2n -1) number of vanishing moments that preserve these vanishing moment properties even at the boundaries of an analyzed window. However, despite these advantages, when the interpolation at the boundaries does not produce a wavelet coefficient with zero value, the resulting boundary wavelet coefficient tends to have a very large magnitude due to the extreme sensitivity of the polynomial interpolation model outside and near the edge of the given set of points around which the interpolation is based.

The wavelet coefficients produced at the boundary, in this manner, tend to have little use, and therefore not only do not solve the problem of wavelet analysis at the boundary, but even make the problem worse! In lay terms, when the polynomial approximation in the lifting scheme works exactly, it works beautifully; however, when it the approximation is off as is, by far, most often the case, it tends to be off by large amounts (a few orders of magnitude is an empirical best guess!). Artifacts can be heard as audible clicks at the boundary.

Luckily, the broad-band denoising mostly alters wavelet coefficients with small values, so that these large values are ignored in the denoising process. The perfect reconstruction property of the lifted wavelet transform just transforms these large spikes back into the original audio, so the listener is unaware of boundary effects due to extremely large wavelet coefficients. However, these large, intermediate wavelet coefficients at the boundaries can sometimes be heard as clicks when the user requests both a primary/dual wavelet with a high number of vanishing moments and a fairly high number of iterations.

In essence, this implies that this style of wavelet analysis and processing will only be useful for applications which deal directly with small-valued wavelet coefficients, since it is almost guaranteed that there will be (artificially) large-valued coefficients at all samples near window boundaries.

Since the standard wavelet decomposition scheme was used in this noise reduction (recursive splitting of the low-pass output), it is no wonder that there are low-frequency noise artifacts left in the denoised output. Soft thresholding only eliminates wavelet coefficients below a certain (threshold) value and adjust others to minimize discontinuity in the final reconstruction; however, since the wavelet decomposition must stop after only a few iterations to prevent the surfacing of boundary artifacts, not all low-frequency noise can be removed. A best basis algorithm (as proposed and implemented by Coifman and Wickerhauser), or a non-standard division of the frequency space would do better in eliminating these low-frequency artifacts.

Enclosures

1) Printed code examples in C++ and MATLAB, printed diagrams
8"test2.m - takes input signal specified in setup.m and analyzes and resynthesizes the signal, setting various levels of wavelet coefficients to zero249"
category description page
1diagram boundary artifacts using circular convolution 3
2" elimination of boundary artifacts using interpolated lifted wavelet transforms 4
3" "Giant Steps" excerpt reconstructions 5
4" "Giant Steps" excerpt wavelet coeffs. 6
5C++ code C++ denoising program for the SGI Indy, includes forward and reverse lifted wavelet transform procedures, neville's algorithm, and driver programs for the command-line "denoise" compiled command. 7
6MATLAB code setup.m - provides easy access for filter selection and test signal selection for commands test1 and test2 20
7" test1.m - test of perfect reconstruction of test signal in setup.m using circular convolution of filters selected in setup.m 23

2) Audio examples

In order on the enclosed audio cassette:
status name source
1original risset, Inharmonique, excerpt LP
2denoised " -
3original Blooming Flowers by the River in Moonlight Night, excerpt Cass.
4denoised " -
5original Killing Joke, Good Samaritan CD
6denoised " -

3) Executables and on-line code

I have included a compressed, tar file of my MATLAB working directory (via Blitzmail). I have also included my C++ code and a working makefile, along with all the necessary files to compile the denoise command under an SGI Irix 5.3 system with development libraries libaudio.a and libaudiofile.a under GNU G++ 2.6.2 (careful with the Indys in Sudikoff 004 - they run GNU G++ 2.7.0, which is stricter about the ANSI C++ standard than 2.6.2). I have also included a pre-compiled binary of the C++ command "denoise." MATLAB Files of interest are listed in the above table, but most information can be found in setup.m, test2.m and flt_test2.m


Go to:
Title Page               Chapter 1           Appendix A   Appendix J
Copyright                Chapter 2           Appendix B   Appendix K
Abstract                 Chapter 3           Appendix C   Appendix L
Acknowledgements         Chapter 4           Appendix D   References
Table of Contents        Chapter 5           Appendix E
List of Figure           Conclusions/        Appendix F
List of Tables           Future Directions   Appendix G
List of Audio Examples                       Appendix H
List of Programs                             Appendix I