AM-FM IMAGE FILTERS

Chuong T. Nguyen and Joseph P. Havlicek

School of Electrical and Computer Engineering, University of Oklahoma, Norman, OK, USA

ABSTRACT

We introduce a multicomponent invertible AM-FM image transform and use it to define new nonlinear AM-FM filters for performing modulation domain image processing. The key elements of the transform are analysis and synthesis filterbanks based on the steerable image pyramid and perfect reconstruction demodulation algorithms based on analytic differentiation of continuous cubic tensor spline models fit to the unwrapped phase samples of a digital image. We demonstrate spatially and spectrally localized orientation and frequency selective filtering, simple image restoration, and image fusion in the modulation domain. These results are also among the first to demonstrate high fidelity image reconstructions from computed multicomponent AM-FM models.

Index Terms— AM-FM image models, AM-FM image filters, modulation domain signal processing, multicomponent models

1. INTRODUCTION

Multicomponent AM-FM image models [1–4] represent an image t: ℝ^{2}→ℝ as a sum

t(x) = ∑_{k=1}^{K} t_{k}(x) = ∑_{k=1}^{K} a_{k}(x) cos[φ_{k}(x)]  (1)

of nonstationary AM-FM functions

t_{k}(x) = a_{k}(x) cos[φ_{k}(x)],  (2)

where x ∈ ℝ^{2}, a_{k}(x) ≥ 0, and φ_{k}(x): ℝ^{2}→ℝ. Given the image t(x), a computed AM-FM model consists of estimates of the K AM functions a_{k}(x) which provide a dense local characterization of the local texture contrast and the K FM functions ∇φ_{k}(x) which provide a dense characterization of the local texture orientation and pattern spacing.

Such models have been used with great success in a variety of image analysis applications, including texture segmentation, 3-D shape from texture, texture-based stereopsis, fingerprint classification, content-based retrieval, and regeneration of occluded and damaged textures [5], as well as for infrared target tracking [6] and in the analysis of (2-D) spectrograms of human speech signals [7]. To date, however, they have been considerably less successful in applications requiring image synthesis in addition to analysis (to the best of our knowledge, reconstruction from a computed AM-FM model has been attempted previously only in [3, 8, 9]). The reason is that some means must be devised for decomposing the image into a sum of components (2) that are locally coherent [3, 4] and for isolating these components from one another on a jointly localized basis in space and spatial frequency prior to demodulation. Because of their Fig. 1. Perfect reconstruction AM-FM model for performing image processing in the modulation domain.

excellent joint localization properties, banks of Gabor or Gabor-like bandpass filters have been used almost universally for this purpose. Unfortunately, the very properties that make these filters attractive, viz., localization and smoothness of the impulse and frequency responses simultaneously, also imply that they cannot provide perfect reconstruction [5].

Our goal is to formulate a new perfect reconstruction (or approximately perfect reconstruction) AM-FM image model and use it to perform signal processing in the modulation domain, as depicted schematically in Fig. 1. Here, the AM-FM model plays a role analogous to that of the Fourier transform in traditional LTI filter theory; in particular, it is a clear requirement that the input and output images be visually indistinguishable when the signal processing block in Fig. 1 is set to the “do nothing” operation.

Previously, we achieved perfect reconstruction in [10] by replacing the typical Gabor filterbank with a nonseparable, orientation selective perfect reconstruction wavelet filterbank. But there were two inherent difficulties with using this approach to formulate practical AM-FM image filters. First, a very large number of filterbank channels was required to generate locally coherent image components t_{k}, on the order of K = 500 for a 256 × 256 image. Second, due to the aliasing that invariably occurs in a critically sampled wavelet filterbank [11], it proved difficult to define signal processing operations on the modulating functions a_{k} and ∇φ_{k} that would lead to predictable, visually meaningful changes in the output image relative to the input.

This work was supported in part by the U.S. Army Research Laboratory and the U.S. Army Research Office under grant W911NF-04-1-0221.

In this paper, we introduce a new invertible AM-FM image transform where the analysis and synthesis filterbanks are based on an appropriate adaptation of the steerable image pyramid decomposition [11, 12]. The steerable pyramid admits several desirable properties that make it attractive for computing AM-FM image models. The bandpass channels are Gabor-like and orientation selective, leading to locally coherent image components that correspond well to visual perception of the significant image structure. Moreover, because the steerable pyramid is essentially alias free and self-inverting, it naturally facilitates the intuitive design of AM-FM filters capable of achieving visually meaningful, perceptually motivated signal processing goals. We briefly describe computation of the invertible AM-FM transform in Section 2, while practical examples of AM-FM image filters are given in Section 3.

2. INVERTIBLE AM-FM TRANSFORM

A number of good approaches have been proposed for computing continuous and discrete domain multicomponent AM-FM image models, including approximate techniques based on the multidimensional Teager-Kaiser energy operator [2] and on a multidimensional extension of the 1-D analytic signal [3], as well as emerging new techniques based on Max-Gabor analysis over local patches [7] and on the empirical mode decomposition [13]. Here we focus on the analytic signal-based approach because, for a complex extension z_{k}(x) = t_{k}(x) + jℋ[t_{k}(x)] = a_{k}(x) exp[jφ_{k}(x)] of a continuous domain AM-FM image component t_{k} given by (2), where ℋ[·] is the partial Hilbert transform, it leads to an exact demodulation algorithm [3]

∇φ_{k}(x) = Re[ \frac{∇z_{k}(x)}{jz_{k}(x)} ]  (3)
a_{k}(x) = |z_{k}(x)|  (4)

that is free of approximation errors. Moreover, provided that boundary conditions on the phase of each component are saved as part of the AM-FM model, the component t_{k} can be recovered exactly from the modulations (3) and (4) by integrating to recover the phase and then reconstructing according to (2).

To isolate the individual components t_{k}(x), we implement the steerable pyramid as described in [11, 12] to a depth of four levels with a non-orientation selective high-pass filter H_{0} and nonorientation selective low-pass filters L_{n} for 0 ≤ n ≤ 3. For each of L_{0}, L_{1}, and L_{2}, the response is partitioned into eight orientation selective channels using oriented bandpass filters B_{p}, 0 ≤ p ≤ 7. Although this construction completes the usual implementation of the steerable pyramid, it is not optimal for AM-FM modeling because the high-pass and low-pass residuals delivered by the channels H_{0} and L_{3} are not orientation selective and hence not locally coherent. Therefore, we adapt the steerable pyramid by further decomposing the responses of H_{0} and L_{3} into eight orientation selective channels using oriented filters [−j cos(θ − pπ/8)]^{7}, 0 ≤ p ≤ 7, obtained by discarding the magnitude characteristic B(ω) from the right-hand side of equation (5) in [12]. The residual remaining in channels H_{0} and L_{3} after this decomposition tends to be small in magnitude and generally contributes little to visual perception of the image; we therefore discard it. This results in an adapted steerable pyramid decomposition with 40 channels delivering K = 40 AM-FM image components t_{k}.

We apply AM-FM demodulation between the analysis and synthesis filterbanks. For each component, the discrete versions of (3) and (4) given in [3] are first applied to obtain the AM functions a_{k} and estimates of the FM functions ∇φ_{k}. The estimated FM functions are then used to guide a multidimensional phase unwrapping algorithm as described in [14]. The unwrapped phase is fit with cubic tensor product splines that we differentiate analytically to obtain the frequency field ∇φ_{k} [10, 15].

3. EXAMPLES

To define AM-FM image filters, we apply signal processing to the computed modulating functions expressed in polar coordinates according to r_{k} = |∇φ_{k}| and ψ_{k} = arg ∇φ_{k}. The filtered image is then obtained by first reconstructing the filtered image components from the processed modulations according to (2) and then subjecting the filtered components to the synthesis filterbank. By effectively realizing the modulation domain signal processing paradigm depicted in Fig. 1, this approach opens the door to an entirely new class of nonlinear image filters that we have only begun to investigate. We present a few examples in this section, where we denote the processed modulations by â_{k}(x) and ∇φ̂_{k}(x).

A synthetic radial chirp image is shown in Fig. 2(a). By construction, the amplitude is constant and the phase is quadratic along radials emanating from the center of the image. Our signal processing goal is to attenuate nonstationary structure that is oriented at odd multiples of π/4. As a baseline for comparison, we implemented an LTI notch filter with frequency response given in Fig. 2(b). The output of this filter is shown in Fig. 2(c) and exhibits undesirable artifacts as expected: the LTI filter attenuates Fourier components on a spatially global scale, which achieves the desired result but also degrades the subtle constructive and destructive interference between Fourier components that creates the image structure at orientations other than odd multiples of π/4. This is demonstrated by Fig. 2(d), which gives the residual between the original image in Fig. 2(a) and the LTI processed image in Fig. 2(c). We define an AM-FM notch filter by

â_{k}(x) = \{ \begin{matrix}16δ_{k}(x)a_{k}(x)/π, & δ_{k}(x) < π/16, \\ a_{k}(x), & otherwise \end{matrix}  (5)

and ∇φ̂_{k}(x) = ∇φ_{k}(x), where δ_{k}(x) = |(|ψ_{k}(x)| − π/4)| is an amplitude scaling factor equal to the radian angular distance between ψ_{k}(x) and ±π/4 (because the partial Hilbert transform has action in the horizontal direction, ψ_{k}(x) is restricted to quadrants I and IV of the frequency plane). The resulting processed image is given in Fig. 2(e), where it may be seen that the signal processing goal has been achieved. Moreover, because the AM-FM filter is capable of attenuating oriented structure on a spatially local basis, it delivers a result that is free of the undesirable artifacts seen in the LTI filter output.

In addition to orientation selective processing, AM-FM filters can also be used to perform spatially local amplification and attenuation based on magnitude frequency. To amplify a band of radial frequencies, we implemented a frequency selective filter according to

â_{k}(x) = \{ \begin{matrix}2a_{k}(x), & 0.2 < r_{k}(x) < 0.35, \\ a_{k}(x), & otherwise \end{matrix}  (6)

and ∇φ̂_{k}(x) = ∇φ_{k}(x), where r_{k}(x) is given in units of cycles per pixel. The AM-FM filtering result is given in Fig. 2(f).

Simplistic AM-FM image restoration is illustrated in Fig. 2(g)- (j). The original image is given in Fig. 2(g), while the image in Fig. 2(h) was obtained by applying a low-pass linear blur and adding Gaussian white noise. As a baseline comparison, the result of a

(a)  (b)  (c)  (d)

(e)  (f)  (g)  (h)

(i)  (j)  (k)  (l)

(m)  (n)  (o)  (p)

Fig. 2. Examples. (a) radial chirp image. (b) LTI notch filter frequency response. (c) Output of LTI notch filter. (d) Difference image between (a) and (c) showing lack of spatial localization in the LTI filtering operation. (e) Orientation selective attenuation performed by AM-FM notch filter. (f) Frequency selective enhancement performed by AM-FM bandpass filter. (g) Original Lena image. (h) Corrupted by linear blur and additive noise. (i) Result of naïve LTI high-pass restoration filter. (j) Result of elementary AM-FM enhancement/restoration filter. (k),(l) Clocks input image pair. (m) AM-FM image fusion result. (n),(o) Registered CT and MR images (imagefusion.org). (p) AM-FM image fusion result. naïve high-pass LTI filter approximating the pseudo-inverse is shown in Fig. 2(i). The design concept for the AM-FM restoration filter combines elements similar to both wavelet shrinkage and unsharp masking. The noise power is distributed widely throughout the steerable pyramid channels resulting in a relatively small contribution to the individual AM functions a_{k}, whereas the coherent image structure tends to be jointly localized resulting in strong contributions to the AM functions, particularly in the vicinity of edges. We applied a simple threshold to the amplitude modulations computed from Fig. 2(h) to attenuate the noise, reconstructed, and then applied the same high-pass filter that was used in Fig. 2(i) to generate a high-pass mask. This mask was added back to the degraded image in Fig. 2(h) to obtain the enhanced/restored result shown in Fig. 2(j).

Finally, elementary examples of AM-FM image fusion based on local contrast are given Fig. 2(k)-(m) and Fig. 2(n)-(p). The main idea is that sharp intensity edges are associated with local AM values a_{k}(x) that are relatively large, whereas defocused surfaces are associated with AM values that are relatively smaller than those of their in-focus counterparts. Thus, for a problem like the well-known pair of clock images shown in Fig. 2(k) and (l), a fused image showing both clocks in focus can be obtained as follows. We first compute a multicomponent AM-FM model for each of the two input images shown in Fig. 2(k) and (l). For each pixel t̂_{k}(x) of component k in the fused image, we take â_{k}(x) and ∇φ̂_{k}(x) directly from the input image for which a_{k}(x) is larger on a pixel-by-pixel basis. The fused image result is given in Fig. 2(m). An identical AM-FM algorithm was used to obtain the result shown in Fig. 2(p) by fusing the CT image of Fig. 2(n) and the MR image of Fig. 2(o) (CT and MR images

courtesy of imagefusion.org).

4. CONCLUSION

In this paper, we introduced an invertible AM-FM image transform based on an adapted steerable image pyramid in combination with powerful perfect reconstruction joint demodulation algorithms. We used this transform to develop new nonlinear AM-FM filters for performing image processing in the modulation domain. Compared to LTI filtering, the advantages of AM-FM filters are that they are capable of performing spatially and spectrally localized processing directly in terms of the visually important nonstationary image structure. The results given here represent an entirely new class of nonlinear image filters that we have only begun to investigate. We demonstrated several simple AM-FM filters, including notch and bandpass filters for amplifying and attenuating spatially local structure based on orientation and granularity (i.e., magnitude frequency) and an elementary AM-FM image restoration filter. We also applied AM-FM filtering to perform image fusion based on local contrast in two typical scenarios. Compared to previous AM-FM image analysis techniques, these results are significant in that they are among the first to demonstrate high fidelity AM-FM image reconstructions and intuitive, practical modulation domain filter design.

Most of the filters considered here relied heavily on AM-based processing. FM-based processing is generally more difficult due to the line singularity associated with the partial Hilbert transform, which tends to perturb the computed frequencies moreso than the amplitudes. We are currently developing strategies for rotating the Hilbert transform on a component-by-component basis to steer the singularity away from the instantaneous frequency vectors ∇φ_{k}. It will also be interesting to investigate a wide variety of AM-FM filtering techniques, including the application of LTI and nonlinear spatial filters directly to the computed modulations.

5. REFERENCES

[1] M.S. Pattichis and A.C. Bovik, “Analyzing image structure by multidimensional frequency modulation,” IEEE Trans. Pattern Anal., Machine Intel., vol. 29, no. 5, pp. 753–766, May 2007.

[2] P. Maragos and A. C. Bovik, “Image demodulation using multidimensional energy separation,” J. Opt. Soc. Amer. A, vol. 12, no. 9, pp. 1867–1876, Sep. 1995.

[3] J. P. Havlicek, D. S. Harding, and A. C. Bovik, “Multidimensional quasi-eigenfunction approximations and multicomponent AM-FM models,” IEEE Trans. Image Proc., vol. 9, no. 2, pp. 227–242, Feb. 2000.

[4] A. C. Bovik, N. Gopal, T. Emmoth, and A. Restrepo, “Localized measurement of emergent image frequencies by Gabor wavelets,” IEEE Trans. Info. Theory, vol. 38, no. 2, pp. 691– 712, Mar. 1992.

[5] J. P. Havlicek, P. C. Tay, and A. C. Bovik, “AM-FM image models: Fundamental techniques and emerging trends,” in Handbook of Image and Video Processing, A. C. Bovik, Ed., pp. 377–395. Elsevier Academic Press, Burlington, MA, 2nd edition, 2005.

[6] C. T. Nguyen, J. P. Havlicek, and M. Yeary, “Modulation domain template tracking,” in Proc. IEEE Int’l. Conf. Comput. Vision, Pattern Recog., Minneapolis, MN, Jun. 17-22, 2007, 8 pp.

[7] T. Ezzat, J. Bouvrie, and T. Poggio, “AM-FM demodulation of spectrograms using localized 2D max-Gabor analysis,” in Proc. IEEE Int’l. Conf. Acoust., Speech, Signal Process., Honolulu, HI, Apr. 15-20, 2007, vol. IV, pp. 1061–1064.

[8] S. Lu and P. C. Doerschuk, “Nonlinear modeling and processing of speech based on sums of AM-FM formant models,” IEEE Trans. Signal Proc., vol. 44, no. 4, pp. 773–782, Apr. 1996.

[9] J. P. Havlicek, D. S. Harding, and A. C. Bovik, “The multicomponent AM-FM image representation,” IEEE Trans. Image Proc., vol. 5, no. 6, pp. 1094–1100, Jun. 1996.

[10] R. A. Sivley and J. P. Havlicek, “Perfect reconstruction AM-FM image models,” in Proc. IEEE Int’l. Conf. Image Proc., Atlanta, GA, Oct. 8-11, 2006, pp. 2125–2128.

[11] E. P. Simoncelli, W. T. Freeman, E. H. Adelson, and D. J. Heeger, “Shiftable multi-scale transform,” IEEE Trans. Info. Theory, vol. 38, no. 2, pp. 587–607, March. 1992.

[12] E. P. Simoncelli and W. T. Freeman, “The steerable pyramid: a flexible architecture for multi-scale derivative computation,” in Proc. IEEE Int’l. Conf. Image Proc., Washington, DC., Oct. 23-26, 1995, pp. 444–447.

[13] N. E. Huang, Z. Shen, S. R. Long, M. C. Wu, H. H. Shih, Q. Zheng, N. -C. Yen, C. C. Tung, and H. H. Liu, “The empirical mode decomposition and the Hilbert spectrum for nonlinear and non-stationary time series analysis,” Proc. R. Soc. Lond. A, vol. 454, pp. 903–995, 1998.

[14] R. A. Sivley and J. P. Havlicek, “Multidimensional phase unwrapping for consistent APF estimation,” in Proc. IEEE Int’l. Conf. Image Proc., Genoa, Italy, Sep. 11-14, 2005, vol. II, pp. 458–461.

[15] R. A. Sivley and J. P. Havlicek, “A spline-based framework for perfrect reconstruction AM-FM models,” in Proc. IEEE Southwest Symp. Image Anal., Interp., Denver, CO, Mar. 26-28 2006, pp. 198–202.