Iterative Projection Based Noise Shaping System Using DT-CWT

Iterative Projection Based Noise Shaping System Using DT-CWT

Rupali Gajanan Gawali1, Charudatta V. Kulkarni2 Department of Electronics & Telecommunication, MIT College of Engineering Pune, India

Abstract

In this paper, we propose Dual Tree Complex Wavelet Transform for image coding, so here we verify how the projection properties of the DT-CWT can be exploited to concentrate signal energy in a few wavelet coefficients greater than a given threshold, while maintaining image reconstruction quality. These non-zero coefficients are forced to compensate using an iterative transform projection loop with error feedback. Also, experiment with different thresholding and error gain strategies, here we use IP-NS system (Iterative projection based noise shaping system). The benefits of applying IP-NS to DT-CWT signals to provide preliminary or simulation results & develop a new and unique system for achieving these transform coding aims of coefficient elimination and compensation. The system is based on iterative projection of signals between the image domain and transform domain.

Keywords- Dual Tree Complex Wavelet Transform, DWT Image coding, Iterative Projection, Noise shaping, PSNR.

1. Introduction

   1. Discrete wavelet Transform:

      The DWT is extensively used in its non-redundant form known as standard DWT. The filterbank implementation of standard DWT for images is viewed as 2-D DWT.

      Two Dimensional Discrete Wavelet Transform (2- DDWT)Filterbank structure is the simple implementation of 1-D DWT, whereas image-processing applications requires two-dimensional implementation of wavelet transform. Implementation of 2-D DWT is also referred to as

      multidimensional wavelet transform in literature[15]. The state of the art image coding algorithms such as e.g. the recent JPEG2000 standard make use of the separable dyadic 2-D DWT [16], which is only an extension of 1-D DWT applied separately on rows and columns of an image.

      The single level 2-D DWT structure produces three detailed sub-images (HL, HL, HH) corresponding to three different directional-orientations 00,450, 900 (Horizontal, Diagonal and Vertical) and a lower resolution sub-image LL. The filterbank structure can be iterated in a similar manner on the LL channel to provide multilevel decomposition.

      DWT 2 level decomposition hierarchy of an image is illustrated in figure (1).

      LLL	LLH	Level 1
      LHL	LHH		LH Vertical
      	Approximate
      	HL Horizontal		HH Diagonal

      LLL	LLH	Level 1
      LHL	LHH		LH Vertical
      	Approximate
      	HL Horizontal		HH Diagonal

      Level 2 Decomposition

      LL

      Figure 1. DWT 2- level decomposition

      Each decomposition breaks the parent image into four child images. Each of such sub-images is of one fourth of the size of a parent image.

      1. Limitations of dwt:

         Although the standard DWT is a powerful tool, it has three major disadvantages that undermine its application for certain signal and image processing tasks. These disadvantages [3, 5] are described as below.

         • Lack of shift invariance: This means that small shifts in the input signal can cause major variations in the distribution of energy between DWT coefficients at different scales.

         • Poor directional selectivity: for diagonal features because the wavelet filters are separable & real.

         • Absence of Phase information: DWT implementations use separable filtering with real coefficient filters associated with real wavelets resulting in real-valued approximations and details. Such DWT implementations cannot provide the local phase information. All natural signals are basically real- valued, hence to avail the local phase information,

           complex-valued filtering is required [8]. The difference between real & analytic wavelet is shown in figure (2).

           Figure 2. Presentation of (a) real, and (b) analytic wavelets

   2. Complex wavelets

      Fortunately, there is a simple solution to these DWT shortcomings [3].

      The Fourier transform is based on complex-valued oscillating sinusoids,

      with . The oscillating cosine and sine components (the real and imaginary parts, respectively) form a Hilbert transform pair; i.e., they are out of phase with each other. Together they constitute an analytic signal that is supported on only one-half of the frequency axis ( > 0).

      Inspired by the Fourier representation, imagine a CWT but with a complex-valued scaling function and complex- valued wavelet,

      Here, r(t) is real and even and ji(t ) is imaginary and odd. Moreover, if r(t ) and i(t ) form a Hilbert transform pair (90 out of phase with each other), then c(t ) is an analytic

      signal and supported on only one-half of the frequency axis. The complex scaling function is defined similarly. See Figure(2) for an example of a complex wavelet pair that approximately satisfies these properties.

      Projecting the signal onto as in figure (2), we obtain the complex wavelet coefficient,

      With magnitude

      and phase

      When . As with the Fourier transform, complex wavelets can be used to analyze and represent both real valued signals (resulting in symmetries in the coefficients) and complex-valued signals as shown in figure(2). In either case, the CWT enables new coherent multiscale signal processing algorithms that exploit the complex magnitude and phase.

      Recent research in the development of CWTs can be broadly classified in two groups:

      • RCWT (Redundant CWTs): The RCWT include two almost similar CWTs. They are denoted as DT-DWT (Dual-Tree DWT based CWT) with two almost similar versions namely Kingsburys DT-DWT(K), and Selesnicks DT-DWT(S). These redundant transforms consist of two conventional DWT filterbank trees working in parallel with respective filters of both the trees in approximate quadrature. The filterbank structure of both DT-DWTs is same but the design methods to generate the filter coefficients are different. Both DT- DWTs provide phase information; they are shift- invariant with improved directionality[9,10,11]

      • NRCWT (Non-redundant CWTs) : NRCWT is the DWT of that complex valued projection, The projection (mapping) means converting a real signal to analytic (complex) form through digital filtering [12]

        Here, we introduce the RCWT base Dual Tree Complex Wavelet Transform (DT CWT)

        1. Dual Tree Complex Wavelet Transforms:

           The Dual-Tree Complex Wavelet Transform (DT- CWT) is an overcomplete, perfect reconstruction, separable transform with Gabor-like filters. It uses two trees per dimension, each with short lowpass and highpass filters, to synthesize a single linear-phase complex lowpass/highpass filter pair [1]. The filters in the two trees are just the time- reverse of each other, as are the analysis and reconstruction filters.

        2. Features of the Dual Tree Complex Wavelet Transform:

      • Approximate shift invariance

      • Good directional selectivity in 2-D, 3-D etc.

      • Perfect reconstruction wth short support filters.

      • Limited Redundancy : 2:1 in 1-D,4:1 in 2-D (2m:1 for m-D)

      • Efficient order N computation :Twice the simple DWT for 1-D (2m for m-D)

        • It has the ability to differentiate positive and negative frequencies and produces six directional subbands oriented in ±15, ±45, ±75 [1]

          Tree b

          Figure 4. Analysis filterbank for 1-D DT-DWT

          Level 1 Level 2 Level 3

          DT-CWT 2 level decomposition hierarchy of an image is illustrated in figure (3).

          -750	-750
          	-450	+150
          -450	-150

          	+750	+750
          +150	+450
          +150		+450

          -750	-750
          	-450	+150
          -450	-150

          	+750	+750
          +150	+450
          +150		+450

          Figure 3. DT-CWT 2- level decomposition

          1. Dual Tree Framework:

             One effective approach for implementing an analytic wavelet transform is called the dual-tree CWT. The idea behind the dual-tree approach is quite simple. The dual

             2

             0(n)

             0(n)

             1(n)

             1(n)

             2

             2 0(n)

             2 1(n)

             Imaginary Tree

             Real Tree

             2 0(n)

             2 1(n)

             2 0(n)

             2 1(n)

             Tree b

             Tree a

             2 0(n)

             2 1(n)

             2 0(n)

             2 1(n)

             tree CWT employs two real DWTs; the first DWT gives the real part of the transform while the second DWT gives the imaginary part. The analysis and synthesis FBs used to implement the dual-tree CWT and its inverse is illustrated in Figures (4) and (5). The two real wavelet transforms use two different sets of filters, with each satisfying the PR (Perfect reconstruction) conditions. The two sets of filters are jointly designed so that the overall transform is approximately analytic. Let h0(n), p(n) denote the low-pass/high-pass filter pair for the upper FB, and let g0(n), g1(n) denote the low- pass/high-pass filter pair for the lower FB [3,13].

             The inverse of the dual-tree CWT is as simple as the forward transform. To invert the transform, the real part and the imaginary part are each invertedthe inverse of each of the two real DWTs are used to obtain two real signals.

             Level 1 Level 2 Level 3

             Figure 5. Synthesis filterbank for 1-D DT-DWT

             The 2-D DT-DWT structure has an extension of conjugate filtering in 2-D case. The filterbank structure of 2- D dual-tree is shown in figure (6). 2-D structure needs four trees for analysis as well as for synthesis. The pairs of conjugate filters are applied to two dimensions (x and y), which can be expressed as [14]:

             Analysis Synthesis

             Real Tree

             h0(n) 2

             p(n)

             p(n)

             2

             Tree a

             h0(n)

             h0(n)

             2

             p(n)

             p(n)

             2

             g0(n)

             g0(n)

             2

             2

             h0(n)

             h0(n)

             p(n) 2

             g0(n) 2

             f(t)

             	Tree a(hxhy)		Tree
             	Tree b(gxgy)		Tree
             Real Trees
             	Tree c(hxgy)		Tree
             	Tree d(gxhy)		Tree

             	Tree a(hxhy)		Tree
             	Tree b(gxgy)		Tree
             Real Trees
             	Tree c(hxgy)		Tree
             	Tree d(gxhy)		Tree

             Imaginary Trees

             Figure 6. Filterbank structure for 2-D DT-CWT

             g0(n) 2

             g1(n) 2

             1. Noise Shaping

                The modification of large coefficients to compensate

                g1(n)

                g1(n)

                g1(n)

                g1(n)

                2 Imaginary Tree

                2

                for the loss of small coefficients is the goal of the noise shaping. The noise shaping scheme is used to remove additive, uncorrelated and white noise. For image coding purposes, we would like to find the configurations that concentrate signal energy in as few non-zero wavelet

                coefficients as possible, while still producing a good approximation of the original signal. The scheme presented here attempts to modify large coefficients to compensate for the loss of small coefficients, without substantially changing the original image. The changes to the large coefficients and loss of small coefficients can be modelled as additive complex-valued noise [1].

                For illustrative purposes, we describe first how we shape noise for a simple hard thresholding system. Coefficients whose magnitudes are less than the threshold are discarded (as in simple denoising applications), but there are no constraints on the values of the retained coefficients.

             2. Thresholding

                There are two thresholding methods frequently used

      • Soft Thresholding & Hard Thresholding [13]

        1.4.1 Soft Thresholding

        The soft threshold function also called the shrinkage function is shown in figure(7),

        y

        -T

        x

        subbands only, while keeping the low resolution coefficients unaltered.

        1. Peak Signal to Noise Ratio (PSNR)

        Peak Signal-to-Noise Ratio, often abbreviated PSNR, is an engineering term for the ratio between the maximum possible power of a signal and the power of corrupting noise that affects the fidelity of its representation. Because many signals have a very wide dynamic range, PSNR is usually expressed in terms of the logarithmic decibel scale. PSNR is used to measure difference between two images.

        MSE is mean square error between two images. PSNR is used to measure the quality between the original & reconstructed image. Higher PSNR indicates that the better quality of reconstructed image.

        T

        Figure 7.Soft Thresholding

        Soft thresholding is an extension of hard thresholding, first setting to zero the elements whose absolute values are lower than the threshold, and then shrinking the nonzero coefficients toward zero.

        1.4.2 Hard Thresholding

        The hard thresholding function as shown in figure(8),

        y

2. Iterative – Projection System

   The DT-CWT is an energy-preserving transform; the signals energy is reduced through thresholding, and as a result of the projection of the noisy signal back into the image domain. The iterative projection forms the basis of noise shaping system [1]. Figure (9) shows the system block diagram.

   DT-CWT i IDT-CWT

   -T x

   T

   Original Input Image

   (X)

   X Xh_cap Xi_cap

   Analysis

   (A)

   Analysis

   (A)

   Hard Thresholding

   Hard Thresholding

   Synthesis

   (R)

   Synthesis

   (R)

   ei =

   =X-Xi_cap

   Figure 8. Hard Thresholding

   which keeps the input if it is larger than the threshold, otherwise it is set to zero.

   The wavelet thresholding procedure removes noise by thresholding the wavelet coefficients of the detailed

   +

   DT-CWT

   k

   k

   Analysis

   Wi (A)

   Figure 9. Block diagram of Iterative Projection System

   As shown in block diagram figure (9), A is Analysis Operator

   R is Synthesis or reconstruction operator k is Energy Gain Factor

   X is Original image

   Xi_cap is Threshold reconstructed image i is Hard threshold

   1. Basic Algorithm:

      To achieve sparsity following steps of Basic algorithm given below:

      • Set iteration i=1 & take the DT CWT of the input image.

      • Set to zero all wavelet coefficients with magnitude smaller than a threshold.

      • Take the inverse DT CWT & measure the error due to loss of smaller coefficients.

      • Take the DT CWT of the error image i.e. error image is projected back into DTCWT domain & added to the non zero wavelet coefficients(Xh_cap) from step 2 to reduce the error & to form the wavelet coefficients Xi+1 at the start of iteration i+1.

      • Increment i, reduce threshold a little & repeat steps 2 to 4

      • When there are sufficient non zero coefficients to give the required rate distortion tradeoff, keep threshold constant & iterate a few more times until converged.

3. Implementation of Iterative Projection System

   • Our iterative projection based noise shaping system applied to 256×256 pixel Lena Image.

   • Let X be the data from an N x N (256×256) real-valued image is as shown in figure(10), rearranged into an

     N 2 x 1 column matrix.

   • The DT-CWT analysis operator matrix A has dimensions N M x N 2, with M = 2N.

   • Level one DT-CWT of image is taken, then the thresholding operation is done. Let the DT-CWT of image is represented by Xi

   • Applying the hard threshold to Xi eliminates insignificant coefficients, leaving the significant coefficients unaffected and producing Xh_cap,where matrix element Xh_cap(m,n) is given by,

     Then is reconstructed, let we call .

   • For ease of result display the value of hard threshold is chosen from 64 to 60.then thresholded & reconstructed image as shown in figure(10). Since we are discarding some of the coefficients, it will obviously introduce some kind of loss of information i.e. error in reconstructed image after thresholding ( .

   • So the error between the original and re-constructed image is ei = X Xi_cap. Image reconstructed from error coefficient ei is shown in figure (10).

   • On the ith iteration, some non-linear operation, such as thresholding or quantising, is applied to wavelet coefficients Xi. This operation can be represented as the addition of noise ni. When the system starts the input to the non-linear operator is Xi, the DT-CWT of image X. After this (i.e. for i>0), the input to the non-linear operator switches to the feedback loop.

   • The noisy coefficients Xh_cap are then projected back into image domain, where the error ei between the original image X and the most recent reconstruction Xi_cap calculated. The error is projected into the DT- CWT domain and added to the noisy coefficients Xh_cap to form the wavelet coefficients Xi+1 at the start of iteration i+1.

   • Although the DT-CWT is an energy-preserving transform, the signals energy is reduced through thresholding, and as a result of the projection of the noisy signal back into the image domain. We obviously need to introduce energy gain.

   • A gain factor k>=1 is applied to the error signal to replace energy lost in the system.

   • The choice of k affects which noise components remain in the system. With unit gain (k=1), the noise components of the DT-CWT analysis operator disappear during the loop.

   • The effect of k=1 is that the number of nonzero coefficients remaining after thresholding decreases with each iteration, while the PSNR of the reconstructed image increases, until the signal converges.

   • When k>1, the energy of the final reconstructed image is closer to the energy of the original. This improves the reconstructed images PSNR, but at the expense of more non-zero coefficients. We consider the number of non-zero coefficients to be a measure of efficiency because zero coefficients are very efficient to code, and so bit rate tends to be proportional to the number of non-zero coefficients

   • We can partially circumvent this trade-off between the reconstruction quality and number of significant coefficients by starting with a large threshold and decreasing it each iteration until our target threshold is reached. A large starting threshold causes many coefficients to be eliminated initially. With the right balance between the energy gain k and the amount by which the threshold decreases each iteration, most insignificant coefficients remain insignificant, while the signals energy is maintained and image reconstruction improves.

   • Our iterative projection-based noise shaping system was applied to the 256 x 256 8 -bit Lena image. A

     threshold value = 32 was applied to the wavelet and scaling coefficients. The number of non-zero coefficients clearly decreases with each iteration, while image quality improves

     • The control system, where the initial threshold 64 is reduced by 1 each iteration until theta = 32, performs best as theta decrease the PSNR increase as shown in table(I & II) below:

       Table I: As Theta decreases PSNR increases of threshold reconstructed image

       Theta	PSNR in DB	Theta	PSNR in DB	Theta	PSNR in DB	Theta	PSNR in DB
       64	24.85	56	27.53	48	29.69	40	32.16
       63	29.71	55	32.25	47	34.13	39	36.31
       62	25.53	54	28.10	46	30.22	38	33.13
       61	30.41	53	32.72	45	34.59	37	37.20
       60	26.23	52	28.66	44	30.78	36	34.54
       59	31.09	51	33.21	43	35.10	35	38.42
       58	26.92	50	29.16	42	31.40	34	36.41
       57	31.71	49	33.66	41	35.65	33	40.00
       						32	38.49

       Table II: As Theta decreases PSNR increases of final reconstructed image

       Theta	PSNR in DB	Theta	PSNR in DB	Theta	PSNR in DB	Theta	PSNR in DB
       64	20.54	56	23.17		48	25.28	40	27.73
       63	25.52	55	27.88	47	29.72	39	31.88
       62	21.25	54	23.73	46	25.81	38	28.69
       61	26.13	53	28.33	45	30.18	37	32.76
       60	21.91	52	24.27	44	26.36	36	30.10
       59	26.76	51	28.82	43	30.68	35	33.98
       58	22.58	50	24.76	42	26.98	34	31.96
       57	27.36	49	29.26	41	31.22	33	35.56
       						32	34.04

     • Image showing output of implementation of iterative projection system:

   • At initial threshold = 64

     Figure 10. Output Image of IP-NS at = 64

   • At final Threshold = 32

     Figure 11. Output Image of IP-NS at = 32

   • Graph between Theta vs. PSNR of threshold reconstructed Lena.bmp image & final reconstructed Lena.bmp image using DT CWT

     Figure 12. Theta Vs PSNR using DT CWT

   • Reduce threshold by 1 at each iteration, while the PSNR of the Threshold reconstructed image& final reconstructed image (New image) increases as shown in above figure (12).

   • As PSNR increases, while image reconstruction quality improves.

4. Conclusion

   Iterative projection based noise shaping (IP-NS) system produces gain over conventional single-stage thresholding. Also, the IP-NS system produces a version of an images DT-CWT as shown in figure (10, 11) that results in the same PSNR reconstruction with 60 to 70% fewer non-zero coefficients than single-stage thresholding, making it more amenable to efficient coding.

   The initial threshold 64 is reduced by 1 each iteration until theta = 32, performs best because slower Threshold reduction rate i.e. 1/iteration produces a DT-CWT signal with fewer non zero coefficients but same reconstruction quality, here as the theta decreases PSNR increases shown in table (II) & image reconstruction improves shown in figure (10,11) with different threshold ().

5. References

1. T.H. Reeves and N.G. Kingsbury, Overcomplete Image Coding Using Iterative Projection-Based Noise Shaping, Signal Processing Group, Department of Engineering, University of Cambridge, U.K., IEEE 2002

2. H. B¨olcskei, F. Hlawatsch, Oversampled filter banks: Optimal noise shaping, design freedom, and noise analysis, Proc. IEEE ICASSP 97, Vol. 3, pp. 24532456, April 1997.

3. Ivan W. Selesnick, Richard G. Baraniuk, and Nick G. Kingsbury, The Dual Tree Complex Wavelet Transform, IEEE signal processing magazine ,november 2005

4. N.G. Kingsbury, A dual-tree complex wavelet transform with improved orthogonality and symmetry properties, Proc. IEEE ICIP 2000, paper 1429, Sept. 2000.

5. N.G. Kingsbury, Complex wavelets for shift invariant analysis and filtering of signals, Journal of Applied Computational Harmonic Analysis, Vol. 10, pp. 234253, May 2001.

6. S.G.Mallat,Z.Zhang,Matching pursuits with time frequency dictionaries,IEEE Trans. Signal Processing, Vol. 41, pp. 33973415, Dec. 1993.

7. A. Abbas and T. Tran, Fast approximations of the orthogonal dual- tree wavelet bases, in Proc. IEEE Int. Conf. Acoust., Speech, Signal Processing (ICASSP), Philadelphia, Mar. 2005, vol. 4, pp. 605608.

8. J M Lina, Image Processing with Complex Daubechies Wavelets, J. of Math. Imaging and Vision, 7, 211-223, 1997

9. N G Kingsbury, Shift Invariant Properties of the Dual-Tree Complex Wavelet Transform, Proc. IEEE Conf. on Acoustics, Speech and Signal Processing, Phoenix, AZ, paper SPTM 3.6, March 16-19, 1999

10. I W Selesnick, Hilbert Transform Pairs of Wavelet Bases, IEEE Sig. Proc. Letters, 8,6,170-173, 2001

11. I W Selesnick, The Design of Hilbert Transform Pairs of Wavelet Bases via the Flat Delay Filter, Proc. of ICASSP, Salt lake City, UT, 2001

12. R Spaendonck, F Fernandes, M Coates, and C Burrus, Non- redundant, Directionally Selective, Complex Wavelets, IEEE Proc. of Int. Conf. Image Proc., 379-382, 2000

13. Dr. Salih Husain Ali & Aymen Dawood Salman, Image Compression Based on 2D Dual Tree Complex Wavelet Transform (2D DT-CWT), Eng. & Tech. Journal,Vol. 28,No. 7,2010

14. D. Srinivasulu Reddy, S. Varadarajan, M. N. GiriPrasad, 2D Dual- Tree Complex Wavelet Transform Based Image Analysis ,

    Contemporary Engineering Sciences, Vol. 5, 2012, no. 3, 127 – 136

15. S G Mallat, A Theory for Multiresolution Signal Decomposition: The Wavelet Representation , IEEE Trans. Pattern Reco. and Machine Int., 11(7), 674-693, 1989

16. ISO/IEC JTC 1/SC 29/WG1 N943, JPEG2000 Requirements and Profiles, International Organization of Standardization, 1999

17. Robi Polikar, The Wavelet Tutorial, Available: http://users.rowan.edu/~polikar/WAVELETS/WTtutorial.html