Abstract
The traditional integer-order partial differential equations and gradient regularization based image denoising techniques often suffer from staircase effect, speckle artifacts, and the loss of image contrast and texture details. To address these issues, in this paper, a difference curvature driven fractional anisotropic diffusion for image noise removal is presented, which uses two new techniques, fractional calculus and difference curvature, to describe the intensity variations in images. The fractional-order derivatives information of an image can deal well with the textures of the image and achieve a good tradeoff between eliminating speckle artifacts and restraining staircase effect. The difference curvature constructed by the second order derivatives along the direction of gradient of an image and perpendicular to the gradient can effectively distinguish between ramps and edges. Fourier transform technique is also proposed to compute the fractional-order derivative. Experimental results demonstrate that the proposed denoising model can avoid speckle artifacts and staircase effect and preserve important features such as curvy edges, straight edges, ramps, corners, and textures. They are obviously superior to those of traditional integral based methods. The experimental results also reveal that our proposed model yields a good visual effect and better values of MSSIM and PSNR.
1. Introduction
A challenging task in designing image denoising model is to retain significant features, such as edges and texture details, while eliminating noise. Among a variety of integer-order partial differential equations based image denoising techniques [1–3], energy minimization based variational technique in two well-known forms of the nonlinear diffusion with anisotropic diffusion proposed by Perona and Malik [4] (Perona-Malik or PM) and total variation (TV) regularization proposed by Rudin et al. [5] (Rudin-Osher-Fatemi or ROF) are the ones that have shown a good edge preservation capability. Most of these nonlinear diffusion models are developed from the following model: by changing the form of , which decides the degree of denoising and preservation of singularities (e.g., noise, edges, and texture). The isotropic diffusion expressed by a linear heat equation has been replaced by the anisotropic diffusion since the classical model was proposed by Perona and Malik in 1990.
It is well known that the PM model yields very favorable results for eliminating image noise while retaining object edges and reducing oscillations. However, it also possesses some unsatisfactory properties in some cases; that is, it easily loses the contrast and texture details and produces staircase effects and so on [6, 7]. To remedy the loss of contrast and texture details, some scholars proposed an iterative regularization method for image restoration [8], and norm in the fidelity term in (1) was replaced by norm in [9–12]. Staircase effect is one of the unsatisfactory properties of the PM model. On the other hand, in [13–17], some models for image noise removal have been introduced to prevent staircase effect, which are based on higher order derivatives to energy norm or integrate high-order derivatives to original model. In addition, a total generalized variation (TGV) [18, 19] has been presented, which contains order greater than or equal to 2, unlike TV, which only considers first-order derivatives. It has some effects on image denoising, which not only leads to piecewise polynomial intensities but also obtains almost perfect edges while efficiently avoiding staircase effect.
Recently, geometric attributes of images have been incorporated into the image denoising processing. The most commonly used local geometrical features are the curvature of the image level curve [20, 21], the Gauss curvature of the image surface [22–25], and the mean curvature of the image surface [26]. Although the aforementioned models have some effects on retaining image contrast, textures, and removing staircase effects, the results are still unsatisfactory. In short, the previously proposed models are integer-order differential based model in essence, they may produce blur edges, and their effect for textures preservation is not well. In particular, the variational method based high-order derivatives of the denoised image often lead to speckle artifacts, that is, isolated black and white speckles.
Thus, in order to preserve textures and achieve a good tradeoff between eliminating speckle artifacts and restraining staircase effect, as a compromise between second-order PDEs and high-order PDEs, fractional-order PDEs based methods have been investigated and applied in the fields of computer vision. Unlike our familiar derivative operator, that is the integer-order derivative operator, the fractional-order derivative operator expands the integer-order derivative operator; it possesses the “nonlocal” property because the fractional-order derivative at a point depends upon the characteristics of the entire function and not just the values in the vicinity of the point [27–30], which is beneficial to improve the performance of textures preservation. Bai and Feng (BF) [31] present an anisotropic diffusion equation based on fractional-order space derivatives for image denoising. Their proposed model is interpolated between the PM equation and the fourth-order anisotropic diffusion equation, which is modeled bywhere and and are the adjoints of and , respectively. The numerical results implemented fractional developmental equation in frequency domain and showed that both of the staircase effect and the speckle artifacts can be eliminated effectively. Guidotti and Lambers [32] extended gradient operator of energy norm from first-order to fractional-order. It has some effects on image denoising. Zhang et al. [33, 34] generalized the classical integer-order TV regularization model for image denoising, which is based on the Grünwald-Letnikov fractional-order derivative. Numerical examples show that the new models can eliminate both staircase effect and speckle artifacts when the parameters are chosen appropriately. Moreover, they are able to retain textures better than the TV model. For the case of calculation, Pu et al. [35] scheduled a class of fractional-order differential masks and illustrated that fractional-order differentiation can handle well fine structures such as texture. With regard to the tomographic reconstruction of the blurry and noisy binary images, a variational approach based on adapted fractional-order Hilbert spaces has been proposed by Bergounioux and Trélat [36]. A class of fully fractional anisotropic diffusion models was presented by Janev et al. [37] for noise removal, where spatial as well as time fractional-order derivatives are employed. Cuesta et al. [38] deduced a Volterra matrix-valued equation for image noise removal from a generalization of the linear fractional integral equation. However, they still have many great potential drawbacks. They simply used the first-order gradient operator of the energy norm to displace fractional-order one and the algorithms do not consider the fractional power of the energy norm, so the loss of image contrasts is not kept well after denoising.
In this paper, we aim to use two new techniques, fractional calculus and difference curvature to image denoising. In comparison to the Gauss curvature and the mean curvature, the difference curvature can more effectively distinguish between ramps and edges. The fractional-order derivatives are able to better deal with the textures of the image and achieve a good tradeoff between eliminating speckle artifacts and restraining staircase effect on the entire image. We combine the difference curvature and fractional calculus frames for one diffusion. Near object edges, the difference curvature is large and the total variation term approximates to the norm of fractional-order gradient flow for the image in order to preserve the edges; away from the edges (in ramp and flat regions), the difference curvature is small and the total variation term approximates to the norm of fractional-order gradient flow for the image in order to prevent staircase effect. It is confirmed that the difference curvature is able to distinguish effectively between ramps and edges. Therefore the proposed model can effectively distinguish between ramps and edges, while improving the performance of texture preservation and preventing staircase effect and speckle artifacts.
The remainder of this paper is organized as follows. In Section 2, the fractional-order derivatives and the new edge indicator (difference curvature) are briefly reviewed. In Section 3, based on fractional-order derivatives and the difference curvature the proposed model and algorithm are presented. Section 4 develops a numerical algorithm to implement our proposed model. Section 5 is devoted to experimental examples. Comparing the simulation results obtained by our proposed denoising model with the baseline denoising methods reports a noticeable improvement in the quality of the denoised images evaluated qualitatively and quantitatively. Section 6 is a conclusion.
2. Preliminaries
2.1. Fractional-Order Derivatives
As a generalization form of the integer-order derivative, the fractional-order derivative has a long history, but it has no uniform definition until now. Three class definitions, Grümwald-Letnikov definition, Riemann-Liouville definition, and Caputo definition [29], are generally carried out by some investigations in recent literatures. Let us introduce the representation of an -fold integral for as follows:Replacing in (3), we obtain a definition of a fractional-order integralwhere is Euler’s continuous gamma function. In particular, if , then . Let , the operator defined asfor is called the Riemann-Liouville (R-L) differential operator with order . However, in this paper, we apply the frequency domain definition to extend the step from integer to fractional, because it is easier for calculating. Therefore, the fractional-order partial derivatives can be formulated as It is natural that the 2D fractional-order gradient operator is able to be formulated as , with the associated magnitude .
The definition of fractional-order Hilbert space goes by using the Fourier transform as follows [31, 36]. For arbitrary number , defineendowed with the normwhere is the Fourier transform of . We denote by the dual of .
2.2. Difference Curvature
It is generally appreciated that the gradient of an image cannot effectively distinguish between the image edges and the image ramps. It is well known that the second-order derivatives are able to effectively distinguish edges and ramps in the case of one-dimensional signal [14]. Similarly, in the two-dimensional case of images, in order to effectively distinguish edges and ramps, we will introduce a new feature descriptor called difference curvature (DC). We definewhere and represent the second-order derivatives in the direction of the gradient and in the direction perpendicular to , respectively, and the operator denotes the absolute value. Here The value of is large at the object edges, while it is small in flat and ramp regions [39], so the edges can be distinguished from flat and ramp regions when the is computed for an image. The property analysis of the new feature descriptor is as follows.(a)At image edges, is large and is small, so is large.(b)In flat and ramp regions, and are both small, so is small.(c)At isolated noise, and are both large but they are almost equal, so is small.
3. Difference Curvature-Based Fractional Anisotropic Diffusion Model
3.1. Description of Our Proposed Model
Let us recall a classic example of noise removal scheme driven by curvature, which is an isotropic diffusion scheme, where the smoothing process is formulated by a partial differential equation. Let denote an image intensity function, is the time (scale parameter), and denotes the speed function which is dependent on the curvature of the level curve; then the diffusion equation can be formulated as If , (11) is the standard curvature evolution equation [21].
With introducing a new coordinate and assigning the grey level of the image to the coordinate, we obtain a level function . Then the image can be viewed as its zero level set corresponds to the surface . In [23], the mean curvature of image surface evolution for noise removal is proposed, which is modeled byA modified mean curvature motion for image smooth is proposed in [24]. Chan et al. [25] developed a diffusion equation which uses the mean curvature as the conductance term, that is, the mean curvature driven diffusion (MCDD) in which the amount of diffusion is dependent on the mean curvature, which can be written aswhere is the mean curvature.
Lee and Seo [26] proposed a Gauss curvature driven diffusion (G-CDD) model in which the amount of diffusion is determined by the Gauss curvature; that is,where is a decreasing and continuous function, with , and is the Gauss curvature of the image surface .
It is well known that the mean curvature and the Gauss curvature are the average and product of and , respectively, which can be defined aswhere MC is mean curvature and GC is Gauss curvature.
Let us discuss the performance comparisons of gradient and three curvatures, that is, mean curvature, Gauss curvature, and difference curvature. Compared with the flaw of anisotropic diffusion which has a tendency to estimate the false edges during the evolution, the difference curvature driven diffusion for images is more appropriate than the gradient driven diffusion to detect the intensity oscillations because the difference curvature can distinguish the edges and ramp regions; however the gradient difference between edges and ramps is not significant.
Although the mean curvature and the difference curvature can effectively distinguish edges and ramps, the mean curvature cannot distinguish edges and isolated noise in a noisy image. For the Gaussian curvature, especially when the image is with isolated noise, the values at corners are large, so the denoising based on the Gaussian curvature driven diffusion cannot retain corners. Meanwhile, the gradient and three curvatures both have a flaw of forming piecewise constant images and indistinctness of texture details because they use the integer-order derivative.
Inspired by fractional-order anisotropic diffusion and the (mean/Gauss) curvature driven diffusion by introducing fractional calculus instead of the (mean/Gauss) curvature, we propose a novel diffusion model for image denoising which uses the difference curvature in the conductance term to describe the intensity variations in images, which is named difference curvature driven fractional-order anisotropic diffusion (DCFAD) and modeled bywhere and are the adjoints of and , respectively, is conductance (diffusivity) function based on Laplacian kernel [40], and is a threshold constant. As shown in Figure 1, compared with two conductance functions of PM diffusion model, the diffusion speed of conductance function based on Laplacian kernel is between them at object edges (where DC is large), and it is smaller than them in flat and ramp regions (where DC is large).
4. Numerical Implementation
In this section, to solve (16), we introduce the frequency domain definition for the fractional-order derivative because it is easy to implement. We suppose the original discrete image has pixels and is uniformly sampled from its continuous version starting at . It is a finite and discrete 2D signal, which stands for the actual image pixels in resolution. We employ the 2D discrete Fourier transform (2D DFT) to approximate the spatial fractional-order difference. Note that we usually handle an image on a bounded domain, while the Fourier domain fractional-order derivative is defined on the whole space. Thus, the image should be prolonged onto the whole space in practical computations.
Fractional central partial difference formulation (see [31]) can have the following form:where are DFT frequencies that correspond to and , respectively. is the inverse 2D DFT operator.
The adjoints of operators and are calculated by the following formulation, respectively:where is the complex conjugation.
Operator has the form , where stands for a diagonal matrix and represents the matrix multiplication. The adjoint operator of is calculated by the following formulation:where
The same is obtained for and . As is proven in [31], when the image size is an odd integer, then is real valued. When is an even integer, the domain of function is not symmetric with respect to , so usually there is a complex component in . However, the image size is large enough (e.g., , as used in the experiments); this complex component is very small due to the fact that the symmetry is minimally disturbed so that it can be ignored in computations.
To solve (16), the finite difference schemes are applied as follows. Let be the time step; we calculate the iterative scheme of the initial image and work in the Fourier domain, only returning to the spatial domain when calculating and , where .
We suppose that the image is defined on , where is the space step size. We use the notation . We denote for or .
The discrete approximation is defined as follows:Our noise removal method for solving (16) is as follows.
Step 1. Let the input image be and set , , , , and ; calculate the 2D DFT of .
Step 2. Calculate the -order central partial differences and using (18).
Step 3. Calculate and in the spatial domain.
Calculate .
Step 4. Calculate and set ; if , calculate the 2D IDFT of and then stop; else go to Step 2.
5. Experimental Results
To analyze and show the good denoising capability of our proposed difference curvature driven fractional-order anisotropic diffusion (DCFAD) for image denoising, some experiments are carried out to compare the denoising result of the proposed method with that of several baseline denoising methods, including the PDE-based methods, the PM model [4], fourth-order PDEs [13], mean curvature-based variational model proposed by Zhu and Chan (ZC) [22], Gauss curvature driven diffusion (G-CDD) [26], the fractional anisotropic diffusion proposed by Bai and Feng (BF) [31], and the other class methods: nonlocal means filtering (NLMF) [41], the wavelet-based method (Bayes least-squares estimate of Gaussian scale mixture, BLS-GSM) [42], and the total generalized variation (TGV) [19]. The proper selection of the fractional-order of our proposed model depends on the noisy image. We have tested the order from 0.1 to 3.0 with step size of 0.2 on every testing image and we find that can get the best results on the “Barbara” image, can get the best results on the “Bicycle” image, and can get the best results on the “Baboon” image. The time step is set to be [35] and . For our model, the peak signal to noise ratio (PSNR) and the mean structural similarity (MSSIM) [43] are adopted to illustrate denoising performance and structure preservation in terms of objective criterion and visual fidelity. We adopt a uniform stopping criterion for all the baseline denoising methods and the proposed method; that is, or the iteration counter becomes larger than 5000 otherwise, where denotes the denoised image in the th iteration. The results of all methods are obtained by performing under MATLAB R2011a running on a desktop with Intel Core i3-2120 CUP 3.30 GHz and 8 GB of memory under 64-bit Windows 7.
Figure 2 shows the results of comparison of some denoising methods on a digital image “Barbara.” The values of PSNR and MSSIM on image “Barbara” are presented in Table 1. In Figure 2(a), the original image is presented and the noisy image with stand variance is shown in Figure 2(b). Figures 2(c)–2(g) show the denoising results with the PM model, fourth-order PDEs, ZC, G-CDD, BF, NLMF, BLS-GSM, TGV, and our proposed DCFAD model, respectively. In order to perform better than the other models, we present the denoising results of a small patch of the original noisy image in Figure 3.
(a)
(b)
(c)
(d)
(e)
(f)
(g)
(h)
(i)
(j)
(k)
(a)
(b)
(c)
(d)
(e)
(f)
(g)
(h)
(i)
(j)
From Figure 3 and Table 1, we can see the following performances. A large scale cartoon part (face) shows staircase effect in the PM model because the flat regions have been smoothed very strongly. The denoising result of the fourth-order PDEs shows speckle artifacts. Small scale texture part (scarf) cannot be kept for the PM and fourth-order PDEs. As seen from Figures 3(d) and 3(i), the noise of the ZC model at the high-frequency edges and textures regions (the stripes on Barbara’s scarf) is reserved; that is, it cannot distinguish edges and isolated noise in the noisy image, while it cannot keep small scale texture when the image includes a high level of noise. In Figure 3(e), the result of the G-CDD model is able to distinguish object edges and ramps but retain the noise at the textures and corners, which will lead to the loss of small scale textures and corners when the image includes a high level of noise. The results of Figures 3(f) and 3(i) exhibit that the BF model yields the false edges because the gradient indicator is sensitive to a high level of noise at the oscillations. As shown in Figure 3(g), the denoising result obtained by the NLMF model is very good but smoothes the image so strongly that it leads to the blur of small scale textures. From Figure 3(h), the wavelet-based method BLS-GSM yields artificial ringing artifacts when the model is reconstructed by wavelet. Although we can see that the result of TGV model avoids blocky effects (staircase effects) near discontinuities in smooth regions and leads to a nature piecewise smooth, the TGV model cannot distinguish the edges and ramp regions because the gradient difference between edges and ramps is not significant. Meanwhile, the TGV contains order greater than or equal to two, which is not able to avoid speckle effects. The values of MSSIM on image “Barbara” imply that our proposed model possesses good local structure preservation capability for the original image. The denoising capability and the local structure-preserving capability of our proposed method can be improved when the fractional-order derivatives and the difference curvature are computed on the predenoised image. Comparing with the denoising results of the other eight models, our proposed model has a good performance in terms of objective criterion and visual fidelity, which obtains better PSNR and is able to retain important features such as edges, ramps, corners, and textures.
Figure 4 shows the results of comparison of some denoising methods on “bicycle” image with many geometric shapes including disks and squares that possess significant amount of highly oscillatory regions (stripes of various thicknesses in the right upper corner, i.e., textures). The noisy image with stand variance is shown in Figure 4(b). Table 1 presents the PSNR and MSSIM of Figure 4. In order to further study the performance of all models, we present the denoising results for a small patch of the original noisy image in Figure 5. It can be confirmed from the enlargements that our proposed model does not introduce other artifacts and has more distinct edges enhancement effect on texture regions than other methods. Thus the difference curvature calculated on the noisy image can improve the denoising performance. For the image denoising results, as it can be seen, PM method preserves edges and corners but has oversmoothed the stripes of various thicknesses in the right upper corner. The problem for fourth-order PDEs model is that speckle effects are not able to be avoided. The ZC model exhibits better visual results on the stripes of various thicknesses in the right upper corner in comparison to the first two methods, but some blurry regions appear in stripes of various thicknesses in the right upper corner in comparison to the proposed method. In Figure 5, the denoising result shows the G-CDD model can preserve the ramp but the denoising capability in the corners is much worse. Figure 5(f) shows the BF model loses a large number of high-frequency edges. As shown in Figure 5(g), the NLMF model smoothes strongly the stripes of various thicknesses in the right upper corner. The wavelet-based method BLS-GSM tends to show some artifacts in Figure 5(h). The TGV model exhibits some blurry regions on the stripes of various thicknesses in the left upper corner compared with the proposed method.
(a)
(b)
(c)
(d)
(e)
(f)
(g)
(h)
(i)
(j)
(k)
(a)
(b)
(c)
(d)
(e)
(f)
(g)
(h)
(i)
(j)
Figure 6 shows the results on a “Baboon” image of size . The noisy image with stand variance is shown in Figure 6(b). Figure 6(c) indicates that the PM model brings in a false edge at the nasolabial groove. The fourth-order PDEs method is able to effectively overcome the staircase effect, but speckle effects are introduced. The ZC model shows the corners preserving capability is worse. It can be observed that the G-CDD model cannot preserve edges at the nasolabial groove. Some blurry edges at the nasolabial groove appear in result of the BF model. The result obtained by the NLMF model smoothes strongly the flat regions of Baboon, which leads to the blur of the hairs on Baboon. The result from wavelet-based method BLS-GSM exhibits good denoising capability but creates some artifacts. As seen, our proposed model is able to better save the Baboon’s nasal surface in smooth areas and nonlinearly keep the bridge of the nose in the regions where the gray level changes greatly and also enhance Baboon’s hairs in those regions where the gray level does not change obviously. The TGV model exhibits speckle effects on Baboon’s nasal surface in smooth areas.
(a)
(b)
(c)
(d)
(e)
(f)
(g)
(h)
(i)
(j)
(k)
In order to perform better than the other models to preserve edges, ramps, corners, and textures, denoising results of a small patch of the original noisy image are demonstrated in Figure 7. To show the practical applicability of our model, we test our proposed method on a CT Scan medical image corrupted by the random white Gaussian additive noise with stand variance . In Figure 8, we show our model can retain the important features, such as curvy edges, straight edges, ramps, corners, and small scale textures.
(a)
(b)
(c)
(d)
(e)
(f)
(g)
(h)
(i)
(j)
(a)
(b)
(c)
(d)
(e)
(f)
(g)
(h)
(i)
(j)
(k)
To illustrate the does not depend on image size, we run our proposed algorithm on three benchmarked images, “Barbara,” “Lena,” and “Peppers,” with two different resolutions, and . Table 2 reports the corresponding results on the MSSIM with . It can be confirmed from Table 2 that the MSSIM values obtained by the testing images with resolution of are close to that of , respectively. It indicates that in our adoptive algorithm does not depend on image size.
6. Conclusions
In this paper, we present an introduction of two new techniques, fractional calculus and difference curvature, to the field of image processing and the implementation of fractional-order PDE, which can effectively distinguish between ramps and edges while preserving important structures such as edges, ramps, corners, and textures and preventing staircase effect and speckle artifacts. The experimental results show that our proposed method is capable of preserving the low frequency contour feature in the smooth area, nonlinearly maintaining the high-frequency edge and texture details in those areas where gray level changes greatly, and also without smoothing the important fine details in those areas where gray level has little changed while preserving corners and effectively distinguishing between ramps and edges.
Conflict of Interests
The authors declare that there is no conflict of interests regarding the publication of this paper.
Acknowledgments
This work was supported by the Fundamental Research Funds for the Central Universities, no. 106112014CDJZR188801, and by Frontier and Application Foundation Research Program of CQ CSTC, no. cstc2014jcyjA40037.