IEEE SIGNAL PROCESSING LETTERS, VOL. 6, NO. 7, JULY 1999 165Generalized Perona–MalikEquation for Image RestorationGuo W. Wei, Member, IEEEAbstract— This letter introduces generalizations of the Per-ona–Malik equation. An edge enhancing functional is proposedfor direct edge enhancement. A number of super diffusion op-erators is introduced for fast and effective smoothing. Statisti-cal information is utilized for robust edge-stopping. Numericalintegration is conducted by using a recently developed quasiin-terpolating wavelet method. Computer experiments indicate thatthe present algorithm is very efficient for edge-detecting andnoise-removing.I. INTRODUCTIONTHE PERONA–MALIK equation [1], proposed in 1990,has recently stimulated a great deal of interest in imageprocessing community [2]–[12]. It is commonly believed thatthe Perona–Malik equation provides a potential algorithmfor image segmentation, noise removing, edge detection, andimage enhancement. The basic idea behind the Perona–Malikalgorithm is to evolve an original image,, under an edge-controlled diffusion operator [1](1)Here,is generalized diffusion coefficient which isso designed that its values are very small at the edge ofan image. Perona and Malik argued that their anisotropicdiffusion equation has no additional maxima (minima) whichdo not belong to the initial image data. This point has beenchallenged recently [2], [10]. It is well-known [2], [3], [5],[6] that this anisotropic diffusion algorithm may break downwhen the gradient generated by noise is comparable to imageedges and features. However, this can be alleviated by usinga regularization procedure.The purpose of this letter is to introduce generalizations ofthe Perona–Malik equation. We introduce an edge-controlledenhancing functional for image enhancement. Super diffusionoperators are introduced for fast and effective smoothing.Statistical information is used for robust edge-stopping. Arecently developed quasi interpolating wavelet algorithm [13]is utilized for the numerical integration of generalized Perona-Malik equations. The efficiency and robustness of the presentManuscript received December 17, 1998. This work was supported in partby the National University of Singapore. The associate editor coordinatingthe review of this manuscript and approving it for publication was Prof. S.Reeves.The author is with the Department of Computational Science, NationalUniversity of Singapore, Singapore 119260, Republic of Singapore (e-mail:[email protected]).Publisher Item Identifier S 1070-9908(99)04962-7.approach is illustrated by the restoration and enhancement ofa noisy Lena image.II. THEORYIn many cases, edges of an image may be worn down forvarious reasons. In other cases, the resolution of an imagecan be very low. Hence, a direct feature enhancement is oftenrequired. We propose a real-valued, bounded edge enhancingfunctional(2)Hereis appropriately chosen so that it isedge sensitive and the contrast of image edges is enhanced.This leads to a generalized Perona–Malik equation(3)Obviously, for constantand a simple form of , (3) reducesto the inhomogeneous diffusion equation.The heat equation can be derived from Fourier’s law for heatflux,, with being a constant. This,from the point of view of kinetic theory, is an approximationto a quasi-homogeneous system that is near equilibrium. Abetter approximation can be expressed as a super flux(4)whereare constants and high order terms describethe influence of inhomogeneity in temperature field and flux-flux correlations to the heat flux. In particular,can beregarded as an energy flux operator. Energy conservation leadsto(5)whereis a source term which can be a nonlinear function.Equation (5) is a generalized reaction-diffusion equation whichincludes not only the usual diffusion and production terms, butalso super diffusion terms. The case of second order super flux,, has been used forthe decription of a number of physical phenomena, such aspattern formation in alloys, glasses, polymer, combustion, andbiological systems.1070–9908/99$10.00  1999 IEEE166 IEEE SIGNAL PROCESSING LETTERS, VOL. 6, NO. 7, JULY 1999(a) (b)(c) (d)Fig. 1. Restoration of a noisy Lena image. (a) Noisy Lena image. (b) Restored with the Perona–Malik diffusion operator. (c) Restored with thePerona–Malik diffusion operator and the edge enhancing functional. (d) Restored with the Perona–Malik diffusion operator, the super-diffusion operator,and the edge enhancing functional.In image systems, the distribution of image pixels can behighly inhomogeneous. Hence the generalized Perona–Malikequation (3) can be made more efficient for image segmen-tation and noise removing by incorporating an edge sensitivesuper diffusion operatorHere are edge sensitive diffusion functions. Aspecial form used in this Letter for numerical experiments is(7)The diffusion functions can be appropriately chosen in manydifferent ways. In this letter, we choose the Gaussian for bothand(8)whereis automatically determined by the noise level, .Here bothand are chosen as a local statistical varianceofand(9)wheredenotes the local average of centered atpoint. The area of local average used in this work is 17 17pixel. In this work, we simply choose the edge enhancingfunctional as(10)withvarying as a function of and .WEI: GENERALIZED PERONA–MALIK EQUATION 167III. NUMERICAL ASPECTS AND RESULTSThe Perona–Malik equation (1) and generalized Per-ona–Malik equations (3) and (7) are spatially discretizedusing a quasi-interpolating wavelet algorithm [13](11)whereis the th derivative of the quasiinterpolating wavelet scaling function [13]. Here we chooseand for remov-ing noise and accurately detecting edges. Since our quasi-interpolating wavelet algorithm has a Schwartz class kernel,no additional regularization procedure is needed. Reflectingboundary condition is used in the time integration. The fourthorder Runge–Kutta method is employed for the time discretiza-tion of the Perona–Malik equation and its generalizations. Atotal of 50 iterations is used in all integrations.We use a 512512 Lena image for our demonstration. TheLena image is degraded with Gaussian noise to obtain a peaksignal-to-noise ratio (PSNR) of 25 dB [Fig. 1(a)]. The PSNRused here is calculated as(12)whereand are the original image and noisyimage samples, respectively, andand are the numberof pixels horizontally and vertically, respectively. A thresholdvalue of 38 is used. Fig. 1(b) shows the