# Nonlocal Adaptive Image Denoising Model

Authors: Xiaoli Sun; Chen Xu; Min Li
Journal: Mathematical Problems in Engineering
 (2013)

---

## Abstract

When denoising with the method of Weickert’s anisotropic diffusion equation, the textures and details will be compromised. A fidelity term is added to Weickert’s equation, and the coefficient of fidelity term will vary adaptively with the instant image, which makes  the diffusion term and the fidelity term come to a better compromise. Otherwise, when deciding the edge directions, because of the strong smoothness of linear Gaussian function, a few other edge directions hiding in the main direction will be lost. To preserve these detailed edge directions, Gaussian kernel is substituted for nonlinear wavelet threshold. In addition, in order to preserve the textures and details as much as possible, a nonlocal diffusion tensor was introduced, and the two eigenvalues are reset by combining the two methods: edge-enhancing diffusion and coherence-enhancing diffusion. Experiments show that the new model has an obvious effect on preserving textures and details.


---

## Body

## 1. Introduction

Partial differential equation has been widely used in image processing. PDE-based image processing can retrospect to Cabor and, Jain, but it is really founded by  Koenderink and Wikin  [1]; they proposed the concept of scale space and their contributions made up of the foundation of PDE-based image processing. The scale space of image can be obtained through evolving image by classical heat diffusion equation, but it is not the only way to get the scale space of image. In the late 1980s, Hummel proposed the criterion of constructing scale space. Based on the criterion, Perona and Malik proposed their anisotropic diffusion equation.However, Perona-Malik model only takes the position information into account, and the diffusion coefficient of the position is decided by the magnitude of gradient modulus. Weickert introduced the concept of structure tensor and proposed his anisotropic diffusion equation. Weickert’s diffusion equation takes local variations of the gradient orientation into account, in some sense, it is really anisotropic. Two mean models of this type are edge-enhancing model and coherence-enhancing model. But the two models have limitations.In order to overcome the defects, a fidelity term is added to Weickert’s diffusion equation, and the coefficient of fidelity term can adaptively vary with the instant image. In addition, to preserve more detailed edge directions hiding in the main direction, we substitute the nonlinear wavelet threshold operator for the Gaussian kernel. In order to preserve the textures and details as much as possible, a nonlocal diffusion tensor was introduced, and the two-eigenvalues are reset by combining the two methods: edge-enhancing diffusion and coherence-enhancing diffusion. Experiments show that the new method has an obvious improvement in vision effect.


## 2. Weickert’s Anisotropic Diffusion Equation

In PM diffusion equation, the diffusion degree is estimated by the magnitude of gradient modulus. In order to take the local variations of the gradient orientation into account, Weickert defined a more general descriptor, and structure tensor is proposed. The anisotropic diffusion becomes as follows:(1)∂u∂t=div(D(Jρ(∇uσ))∇u),u(0,x)=u0(x),
where uσ=u(t,•)*Kσ,  Kσ is a Gaussian kernel, which avoids the influence of noise in estimating the   gradient;  Jρ(∇uσ)=Kρ*(∇uσ∇uσ⊥), it is called a structure tensor [2, 3]; the local information is averaged by convolving component wise ∇uσ∇uσ⊥ with a Gaussian kernel Kρ, which avoids the fact that neighbouring gradients with same orientation but opposite sign offset one another. D(Jρ) is called diffusion tensor; it has the same orthonormal basis of eigenvectors as Jρ, and the choice of corresponding eigenvalues of D(Jρ) depends on the demand. Suppose that the two orthonormal bases and corresponding eigenvalues are v1, v2 and μ1, μ2, v1 indicates the orientation maximizing the gray-value fluctuations, while v2 is perpendicular to v1. The eigenvalues μ1 and μ2 convey the shape information. Isotropic structures are characterized by μ1≅μ2≅0, that is, the smooth areas, linelike edges by μ1≫μ2≈0, and corners by μ1≥μ2≫0.In edge-enhancing diffusion model [4], J(∇uσ)=(∇uσ∇uσ⊥). Supposing that the two orthonormal bases and corresponding eigenvalues are v1,v2 and μ1,μ2, the eigenvalues of diffusion tensor D(Jρ) are
(2)λ1={1,ifμ1=0,1-exp(-1μ14)otherwise,λ2=1,
where λ1 is the eigenvalue that corresponds to v1 and λ2 is the eigenvalue that corresponds to v2. We can notice that  λ2 is large and λ1 decreases with the increase of μ1; this process can preserve edges, while, J(∇uσ) cannot get right edge direction. If σ is too small, the noises cannot be fully avoided, and the decided direction of edge  will be  fluctuated greatly; if σ is too large, the neighbouring gradients with the same orientation but opposite sign will offset one another, and we cannot get right edge directions. So, edge-enhancing diffusion is only applicable for images only with large structures.Different from edge-enhancing diffusion model, in coherence-enhancing model [5], Jρ(∇uσ)=Kρ*(∇uσ∇uσ⊥). Supposing that the two orthonormal bases and corresponding eigenvalues are v1,v2 and μ1,μ2, the eigenvalues of diffusion tensor D(Jρ) are
(3)λ1=α,λ2={α,ifμ1=μ2,α+(1-α)exp(-1(μ1-μ2)4),otherwise,
where α∈(0,1), which keeps the diffusion tensor uniformly positive definite. Because of the act of large-scale Gaussian kernel, a few edge directions hiding in the main direction will be lost. In addition, because of the small α, even if in the smooth area, λ2 is much larger than λ1, and many artificial edges are formed in the smooth area. So, coherence-enhancing model is only applicable for images with only linelike edges such as fingerprint image.


## 3. The Decision of Edge Directions Based on Nonlinear Wavelet Threshold

Based on the previous analysis, we know that some detailed edge directions will be lost because of the strong smoothness of Gaussian kernel [6]. So, we substitute the Gaussian kernel for the nonlinear wavelet threshold. The new structure tensor based on nonlinear wavelet threshold is as follows:
(4)JSn(∇uσ)=Sn(∇uσ∇uσ⊥)=[Sn((uσ,x)2)Sn(uσ,xuσ,y)Sn(uσ,xuσ,y)Sn((uσ,y)2)],
where Sn is the nonlinear wavelet threshold operator. The comparison of edge directions of a remote image, based on Sn and Gaussian operator Kρ, is illustrated in Figure 1. The edge direction is depicted by grey values: horizontal direction is depicted in black and vertical direction in white. We choose Sn as the soft wavelet threshold operator.The comparison of edge directions for remote sensing image.
(a)
(b)
(c)Figure1(a)    is   the result of a  noise remote sensing image, and   Figure 1(b)   is the result of edge directions found by Kρ. We can see that some detailed directions have been lost.   Figure 1(c)   is the result of edge directions found by Sn. We can even see that more detailed edge directions have been detected, and it get more accurate edge directions.


## 4. The New Model Based on Nonlocal Diffusion Tensor

Suppose that the image functionu(X)∈L2(Ω),  Ω∈R2 is a bounded open region, and the notion of derivatives was extended to a nonlocal framework by the following definition [7]:
(5)∂yu(X)=u(Y)-u(X)d~(X,Y),
where X=(x1,x2)∈Ω, Y=(y1,y2)∈Ω, 0<d~(X,Y)≤∞ is a positive measure defined between the points X and Y. To keep up with standard notations related to graphs [8], we define the weights as follows:
(6)ω(X,Y)=d~-2(X,Y).
Thus, 0<ω(X,Y)≤∞. We assume, that the weights are symmetric; that is, ω(X,Y)=ω(Y,X).The nonlocal derivative can be written as(7)∂yu(X)=[u(Y)-u(X)]ω(X,Y).The nonlocal gradient∇ωu(X):Ω→Ω×Ω is defined as the vector of all partial derivatives:
(8)∇ωu(X)=[u(Y)-u(X)]ω(X,Y).We denote vectors asv=v(X,Y)∈Ω×Ω. The standard L2 inner product is used for functions
(9)〈u1,u2〉=∫Ωu1(X)u2(X)dX.
For vectors, we define a dot product
(10)(v1·v2)(X)=∫Ωv1(X,Y)v2(X,Y)dY
and an inner product
(11)(v1·v2)=〈v1·v2,1〉=∫Ω×Ωv1(X,Y)v2(X,Y)dXdY.
The magnitude of a vector is
(12)|v|(X)=v1·v1=∫Ωv2(X,Y)dY.With the above inner products, for vectorv=v(X,Y), the nonlocal divergence is defined as
(13)divωv(X)=∫Ω[v(X,Y)-v(Y,X)]ω(X,Y)dY.The Laplacian can now be defined by(14)Δωu(X)=12divω(∇ωu(X))=∫Ω[u(Y)-u(X)]ω(X,Y)dY.Note that in order to get the standard Laplacian definition which relates to the graph Laplacian, we need a factor of1/2.We can also formulate nonlocal curvatures as follows:(15)Kω=divω(∇ωu|∇ωu|)=∫Ω(u(Y)-u(X))ω(X,Y)×(1|∇ωu(X)|+1|∇ωu(Y)|)dY,
where
(16)|∇ωu|(Q)=∫Ω(u(Z)-u(Q))2ω(Q,Z)dZ.Based on the analysis in Section3 and the definition of nonlocal operators mentioned previously, we extend the local tensor diffusion to nonlocal one; that is,
(17)∂u∂t=divω(D(Jsn(∇ωuσ))∇ωu),u(0,x)=u0(x),
where ∇ω is the nonlocal gradient and divω is the nonlocal divergence. D is a symmetric positive definite nonlocal diffusion tensor. Starting from image u0(x), series of images u(t,x) are obtained, the noises in which will be gradually reduced. The final denoised image should preserve textures and details of the original image as much as possible. Next, how to determine the weights ω(X,Y) and the nonlocal diffusion tensor D will be discussed.


### 4.1. The Definition of Weightsω(X,Y)

The calculation of weights was divided into two steps. First, a distance measure of two points should be given, and then the weights can be defined according to the distance measure.Taking into account the similarity and Euclidean distance, the distance measure ofX=(x1,x2) and Y=(y1,y2) is given [9] as follows:
(18)dα(X,Y)=∫ΩGα(Z)|f(X+Z)-f(Y+Z)|2dZ+∥X-Y∥,
where Gα(Z) is a Gauss function with standard deviation α, the first integral term measures the similarity of X and Y,  and ∥X-Y∥ is the Euclidean distance between X and Y. Only both the gray value and the Euclidean distance of the two pixels are similar; the two pixels are the similar pixels.Now, the weightsω(X,Y) can be defined according to the distance measure dα(X,Y). For model (17), in the diffusion process, we hope the similar pixels interact with each other, so the weights of similar pixels should be larger; that is, the weights ω(X,Y) should be a decreasing function of distance dα(X,Y).   For each point we define the following set A of area |A|=γ (a parameter) within a search neighborhood  s(X) around X (where A⊂S(X)⊆Ω, |A|≪|s(X)|).
(19)A(x)=argminA∈s(X)[∫Adα(X,Y)dY,|A|=γ≪|s(X)|],
where |A| is the area of region A. In the discrete case, |A| expresses the number of pixels in region A and  γ is an area parameter. The significance of (19) is that in a larger neighborhood s(X)∈Ω of X, γ pixels whose gray values are the most similar to X are searched; the γ pixels constitute a set recorded as A(x). The weights of pixel Y∈Ω to X are defined as
(20)ω(X,Y)={1,Y∈A(X)orX∈A(Y)0,other.


### 4.2. The Design of Diffusion TensorD

In the area of isotropic structures, strong diffusion is needed in two directions to denoise efficiently, and the two eigenvalues should be large; in linelike edges, strong diffusion is needed in edge direction to enhance the edges, little diffusion is needed in the direction perpendicular to edges to preserve edges, the eigenvalue in edge direction should be large, and the one perpendicular to edges should be small; in the corner, little diffusion is needed to preserve the corner and the two eigenvalues should be small.Based on the above analysis, the two eigenvalues are reset by combining edge-enhancing diffusion model and coherence-enhancing model. Lettings=μ12+μ22, the new eigenvalues are as follows:
(21)λ1*=c(x,y)λ11+(1-c(x,y))λ12,λ2*=c(x,y)λ21+(1-c(x,y))λ22,
where λ1 and  λ2 are as follows:
(22)λ1:{λ11=1-exp(-1μ14),λ21=1,λ2:{λ12=α,λ22=1-exp(-(μ1-μ2)42),
where c(x,y) is a weight function, and it decreases with s=μ12+μ22 and satisfies 0≤c(x,y)≤1. In this paper, we choose c(x,y)=k/(k+s).  c(x,y) can be other operators such as exponent form.Now, we analyze the eigenvalues in different areas. In isotropic-structure areas,μ1≅μ2≅0, so λ1*→1; λ2*→1, the two eigenvalues are large, which assure the strong diffusion in two directions and the effective de-noising in smooth areas.. In linelike edges,  μ1≫μ2≈0. Now λ1* is smaller, and λ2* is larger, which assures the preservation and enhancement of edges. In the corner, μ1≥μ2≫0, λ1* and λ2* are both smaller, which assures the little diffusion in two directions and the preservation of corners.Similar to local structure tensor proposed above, the nonlocal structure tensor of each pixel is defined as follows:(23)JSn(∇ωuσ)=Sn(∇ωuσ∇ωuσ⊥),
where uσ=Gσ*u, Sn is the nonlinear wavelet threshold operator [10].Supposing that the two orthonormal bases and the corresponding eigenvalues ofJSn(∇ωuσ) are v1, v2 and  μ1, μ2, the two orthonormal bases of diffusion tensor D are same  as   JSn(∇ωuσ), and the two corresponding eigenvalues of diffusion tensor D are reset as above.After determining the weightsω(X,Y) and the nonlocal diffusion tensor D, the nonlocal tensor diffusion model (17) is obtained.


## 4.1. The Definition of Weightsω(X,Y)

The calculation of weights was divided into two steps. First, a distance measure of two points should be given, and then the weights can be defined according to the distance measure.Taking into account the similarity and Euclidean distance, the distance measure ofX=(x1,x2) and Y=(y1,y2) is given [9] as follows:
(18)dα(X,Y)=∫ΩGα(Z)|f(X+Z)-f(Y+Z)|2dZ+∥X-Y∥,
where Gα(Z) is a Gauss function with standard deviation α, the first integral term measures the similarity of X and Y,  and ∥X-Y∥ is the Euclidean distance between X and Y. Only both the gray value and the Euclidean distance of the two pixels are similar; the two pixels are the similar pixels.Now, the weightsω(X,Y) can be defined according to the distance measure dα(X,Y). For model (17), in the diffusion process, we hope the similar pixels interact with each other, so the weights of similar pixels should be larger; that is, the weights ω(X,Y) should be a decreasing function of distance dα(X,Y).   For each point we define the following set A of area |A|=γ (a parameter) within a search neighborhood  s(X) around X (where A⊂S(X)⊆Ω, |A|≪|s(X)|).
(19)A(x)=argminA∈s(X)[∫Adα(X,Y)dY,|A|=γ≪|s(X)|],
where |A| is the area of region A. In the discrete case, |A| expresses the number of pixels in region A and  γ is an area parameter. The significance of (19) is that in a larger neighborhood s(X)∈Ω of X, γ pixels whose gray values are the most similar to X are searched; the γ pixels constitute a set recorded as A(x). The weights of pixel Y∈Ω to X are defined as
(20)ω(X,Y)={1,Y∈A(X)orX∈A(Y)0,other.


## 4.2. The Design of Diffusion TensorD

In the area of isotropic structures, strong diffusion is needed in two directions to denoise efficiently, and the two eigenvalues should be large; in linelike edges, strong diffusion is needed in edge direction to enhance the edges, little diffusion is needed in the direction perpendicular to edges to preserve edges, the eigenvalue in edge direction should be large, and the one perpendicular to edges should be small; in the corner, little diffusion is needed to preserve the corner and the two eigenvalues should be small.Based on the above analysis, the two eigenvalues are reset by combining edge-enhancing diffusion model and coherence-enhancing model. Lettings=μ12+μ22, the new eigenvalues are as follows:
(21)λ1*=c(x,y)λ11+(1-c(x,y))λ12,λ2*=c(x,y)λ21+(1-c(x,y))λ22,
where λ1 and  λ2 are as follows:
(22)λ1:{λ11=1-exp(-1μ14),λ21=1,λ2:{λ12=α,λ22=1-exp(-(μ1-μ2)42),
where c(x,y) is a weight function, and it decreases with s=μ12+μ22 and satisfies 0≤c(x,y)≤1. In this paper, we choose c(x,y)=k/(k+s).  c(x,y) can be other operators such as exponent form.Now, we analyze the eigenvalues in different areas. In isotropic-structure areas,μ1≅μ2≅0, so λ1*→1; λ2*→1, the two eigenvalues are large, which assure the strong diffusion in two directions and the effective de-noising in smooth areas.. In linelike edges,  μ1≫μ2≈0. Now λ1* is smaller, and λ2* is larger, which assures the preservation and enhancement of edges. In the corner, μ1≥μ2≫0, λ1* and λ2* are both smaller, which assures the little diffusion in two directions and the preservation of corners.Similar to local structure tensor proposed above, the nonlocal structure tensor of each pixel is defined as follows:(23)JSn(∇ωuσ)=Sn(∇ωuσ∇ωuσ⊥),
where uσ=Gσ*u, Sn is the nonlinear wavelet threshold operator [10].Supposing that the two orthonormal bases and the corresponding eigenvalues ofJSn(∇ωuσ) are v1, v2 and  μ1, μ2, the two orthonormal bases of diffusion tensor D are same  as   JSn(∇ωuσ), and the two corresponding eigenvalues of diffusion tensor D are reset as above.After determining the weightsω(X,Y) and the nonlocal diffusion tensor D, the nonlocal tensor diffusion model (17) is obtained.


## 5. Anisotropic Diffusion Equation with Adaptive Fidelity Term

### 5.1. Anisotropic Diffusion Equation with Scalar-Parameter Fidelity Term

New eigenvalues are obtained by the method of Section4, and the new diffusion tensor, expressed as D(Jsn), is determined. In order to preserve the detailed structures as much as possible and make the final denoised image faithful to initial noise image, a fidelity term with parameter λ is added to the equation, and the new diffusion equation is as follows:
(24)∂u∂t=divω(D(Jsn(∇ωuσ))∇ωu)+λ(u0-u).Suppose that the noise can be approximated by Gaussian white noise, and the variance isζ [11], when the equation reaches its stable condition, and  u does not vary with the time t. So the equation is
(25)divω(D(Jsn(∇ωuσ))∇ωu)+λ(u0-u)=0.After denoising, we hope that the residual part contains all noises; that is,(26)1Ω∫Ω(u-u0)2dxdy=ζ2.Multiplying (25) by (u-u0) and integrating  over Ω,  we get
(27)∫Ωdivω(D(Jsn(∇ωuσ))∇ωu)•(u-u0)dxdy=λ∫Ω(u-u0)2dxdy.Substituting (26) into (27), we get
(28)λ=1ζ2|Ω|∫Ωdivω(D(Jsn(∇ωuσ))∇ωu)•(u-u0)dxdy.After each step of iteration, substitute the instant imageu into (28), and we get the new parameter λ of fidelity term, which is proper for the instant time scale. Next iteration is processed by substituting the new λ into (24). Experiments show that diffusion equation with fidelity term is better in preserving edges and details than  the one  without fidelity term.


### 5.2. Anisotropic Diffusion Equation with Adaptive Fidelity Term

Let the initial noise image beu0; our model consists of three components: u0=uC+uNC+un,  uC is the Cartoon approximation, and uNC is the remainder noncartoon part; it typically consists of textures, small-scale details, thin lines, and so forth. The only assumption we make is that it has zero mean. They should be preserved during the course of denoising. un is an additive noise. In the paper, the residue of an image is uR=u0-u.Now, we define the power of image  as  its variance [10]:
(29)pz=var(uz).The markz, as mentioned above, denotes the different parts of an image. The power of noise is defined as its variance, which is known or can be estimated and denoted by pn.Let us define first a measure which we refer to as the local power [11]:
(30)Pz(x,y)=1|Ω|∫Ω(uz(x~,y~)-E(uz))2ωx,y(x~,y~)dx~dy~,
where ωx,y(x~,y~)=ω(|x~-x|,|y~-y|) is a normalized (∫Ωωx,y(x~,y~)dx~dy~=1) and radially symmetric smoothing window and  E(uz) is the expected value. From the definition of the local power it follows that
(31)∫ΩPz(x,y)dxdy=pz≡var(uz).Different contents need different extent of denoising. In the position of large structures, denoising efficiently is needed. After denoising, we hope that the power of the residual part is equal to the power of noise. While in the position of detailed structures, such as textures, in order to preserve the details, as a compromise, parts of the noise should be left in the denoised image; we hope that the power of the residual part is less than the power of noise. Now, suppose thats(x,y)=PR(x,y)≥0, and it is prior information; the definition of this prior information will be given in next chapter.The scalarλ was replaced by an adaptive function λ-(x,y), and the diffusion equation is as follows:
(32)∂u∂t=divω(D(Jsn(∇ωuσ))∇ωu)+λ-(x,y)(u0-u),
whose stable condition is
(33)λ-(x,y)(u-u0)-divω(D(Jsn(∇ωuσ))∇ωu)=0.
We multiply (33) by (u-u0) and integrate over as follows:
(34)∫Ω((D(Jsn(∇ωuσ))∇ωu)λ-(x,y)(u-u0)2-(u-u0)divω(D(Jsn(∇ωuσ))∇ωu))dxdy=0.
In order to compute λ-(x,y), we denote λ-(x,y) as the convolution of function λ(x,y) with window ωx,y(x,y) as follows:
(35)λ-(x,y)=∫Ωλ(x~,y~)ωx,y(x~,y~)dx~dy~.
Substituting (35) into (34), after a change in the order of integrals in the λ(x,y) we get
(36)∫Ω(λ(x,y)s(x,y)-Q(x,y))dxdy=0,
where Q(x,y)=(u-u0)divω(D(Jsn(∇ωuσ))∇ωu). A sufficient condition is
(37)λ(x,y)=Q(x,y)s(x,y).
Substituting (37) into (35), we can get the value of λ-(x,y).Now, the definition ofs(x,y) will be discussed. In (37), in order to get the value of  λ-(x,y), the problem reduces to finding s(x,y).  s(x,y) should adaptively vary with the content of position (x,y). When (x,y) lies in the large-structure area, s(x,y) should be large and fidelity becomes small, and then the diffusion equation can denoise efficiently in this position. When (x,y) lies in the area which contains details, s(x,y) should be small and fidelity term becomes large, and then the diffusion equation can reject the denoising and preserving details.For images with only large structures (so-called cartoon-type image), it requires full denoising, and the power of the residue should approximate the power of noise; that ispR≈ζ2. In fact, many images are noncartoon types; they typically contain textures, small-scale details, and so forth. If it is still fully denoised, these details will be lost or blurred. To overcome the defect that the denoised image will distort because of the loose of details, denoising and detail preserving should come to a better compromise. So the non-cartoon images require less denoising, and some noises still need to be in the denoised image. Now, the power of the residue should be less than the power of noise; that is pR<ζ2.After the previous analysis, we use the scalar-parameter method mentioned in   Section5.1 to choose the first metastable scale. In our paper, we choose pR=1.5ζ2 (this condition is actually implicit as there is no stable scale with residue power below the noise level). Thus, after the process of scalar parameter, cartoon part of the image uC and the noise and texture parts  uNC+un are differentiated. A majority of details and textures are left in the residual part.Substituting1.5ζ2 for ζ2 into (28), we can get the residue uR by iterating (24). Substituting uR into (27), we can get the local power PR(x,y), we assign
(38)s(x,y)=ζ4PR(x,y).Now, let us analyse (38). When (x,y) lies in the large-structure area, the local power of the residue is almost constant (PR(x,y)≈ζ2), and hence s(x,y)≈ζ2, and this scheme degenerates to the scalar-parameter one mentioned in   Section 5.1. We get a high quality denoising process. When (x,y) lies in the area which contains details, the local power of the residue contains two parts: PNC(x,y)+Pn(x,y), so PR(x,y)≈PNC(x,y)+ζ2 and PR(x,y)≫ζ2, and now, s(x,y) becomes small and λ-(x,y) becomes large. We get a low quality denoising process to preserve detailed structures.Using the scalar-parameter method mentioned in   Section5.1 to choose the first metastable scale, after each step of iteration, substitute the instant image u into (35). We get the new parameter λ- of fidelity term, which is proper for the instant time scale. Next iteration is processed by substituting the new λ- into (32). Experiments show that the diffusion equation with adaptive fidelity term is better in preserving edges and details than that with scalar-parameter fidelity term.


## 5.1. Anisotropic Diffusion Equation with Scalar-Parameter Fidelity Term

New eigenvalues are obtained by the method of Section4, and the new diffusion tensor, expressed as D(Jsn), is determined. In order to preserve the detailed structures as much as possible and make the final denoised image faithful to initial noise image, a fidelity term with parameter λ is added to the equation, and the new diffusion equation is as follows:
(24)∂u∂t=divω(D(Jsn(∇ωuσ))∇ωu)+λ(u0-u).Suppose that the noise can be approximated by Gaussian white noise, and the variance isζ [11], when the equation reaches its stable condition, and  u does not vary with the time t. So the equation is
(25)divω(D(Jsn(∇ωuσ))∇ωu)+λ(u0-u)=0.After denoising, we hope that the residual part contains all noises; that is,(26)1Ω∫Ω(u-u0)2dxdy=ζ2.Multiplying (25) by (u-u0) and integrating  over Ω,  we get
(27)∫Ωdivω(D(Jsn(∇ωuσ))∇ωu)•(u-u0)dxdy=λ∫Ω(u-u0)2dxdy.Substituting (26) into (27), we get
(28)λ=1ζ2|Ω|∫Ωdivω(D(Jsn(∇ωuσ))∇ωu)•(u-u0)dxdy.After each step of iteration, substitute the instant imageu into (28), and we get the new parameter λ of fidelity term, which is proper for the instant time scale. Next iteration is processed by substituting the new λ into (24). Experiments show that diffusion equation with fidelity term is better in preserving edges and details than  the one  without fidelity term.


## 5.2. Anisotropic Diffusion Equation with Adaptive Fidelity Term

Let the initial noise image beu0; our model consists of three components: u0=uC+uNC+un,  uC is the Cartoon approximation, and uNC is the remainder noncartoon part; it typically consists of textures, small-scale details, thin lines, and so forth. The only assumption we make is that it has zero mean. They should be preserved during the course of denoising. un is an additive noise. In the paper, the residue of an image is uR=u0-u.Now, we define the power of image  as  its variance [10]:
(29)pz=var(uz).The markz, as mentioned above, denotes the different parts of an image. The power of noise is defined as its variance, which is known or can be estimated and denoted by pn.Let us define first a measure which we refer to as the local power [11]:
(30)Pz(x,y)=1|Ω|∫Ω(uz(x~,y~)-E(uz))2ωx,y(x~,y~)dx~dy~,
where ωx,y(x~,y~)=ω(|x~-x|,|y~-y|) is a normalized (∫Ωωx,y(x~,y~)dx~dy~=1) and radially symmetric smoothing window and  E(uz) is the expected value. From the definition of the local power it follows that
(31)∫ΩPz(x,y)dxdy=pz≡var(uz).Different contents need different extent of denoising. In the position of large structures, denoising efficiently is needed. After denoising, we hope that the power of the residual part is equal to the power of noise. While in the position of detailed structures, such as textures, in order to preserve the details, as a compromise, parts of the noise should be left in the denoised image; we hope that the power of the residual part is less than the power of noise. Now, suppose thats(x,y)=PR(x,y)≥0, and it is prior information; the definition of this prior information will be given in next chapter.The scalarλ was replaced by an adaptive function λ-(x,y), and the diffusion equation is as follows:
(32)∂u∂t=divω(D(Jsn(∇ωuσ))∇ωu)+λ-(x,y)(u0-u),
whose stable condition is
(33)λ-(x,y)(u-u0)-divω(D(Jsn(∇ωuσ))∇ωu)=0.
We multiply (33) by (u-u0) and integrate over as follows:
(34)∫Ω((D(Jsn(∇ωuσ))∇ωu)λ-(x,y)(u-u0)2-(u-u0)divω(D(Jsn(∇ωuσ))∇ωu))dxdy=0.
In order to compute λ-(x,y), we denote λ-(x,y) as the convolution of function λ(x,y) with window ωx,y(x,y) as follows:
(35)λ-(x,y)=∫Ωλ(x~,y~)ωx,y(x~,y~)dx~dy~.
Substituting (35) into (34), after a change in the order of integrals in the λ(x,y) we get
(36)∫Ω(λ(x,y)s(x,y)-Q(x,y))dxdy=0,
where Q(x,y)=(u-u0)divω(D(Jsn(∇ωuσ))∇ωu). A sufficient condition is
(37)λ(x,y)=Q(x,y)s(x,y).
Substituting (37) into (35), we can get the value of λ-(x,y).Now, the definition ofs(x,y) will be discussed. In (37), in order to get the value of  λ-(x,y), the problem reduces to finding s(x,y).  s(x,y) should adaptively vary with the content of position (x,y). When (x,y) lies in the large-structure area, s(x,y) should be large and fidelity becomes small, and then the diffusion equation can denoise efficiently in this position. When (x,y) lies in the area which contains details, s(x,y) should be small and fidelity term becomes large, and then the diffusion equation can reject the denoising and preserving details.For images with only large structures (so-called cartoon-type image), it requires full denoising, and the power of the residue should approximate the power of noise; that ispR≈ζ2. In fact, many images are noncartoon types; they typically contain textures, small-scale details, and so forth. If it is still fully denoised, these details will be lost or blurred. To overcome the defect that the denoised image will distort because of the loose of details, denoising and detail preserving should come to a better compromise. So the non-cartoon images require less denoising, and some noises still need to be in the denoised image. Now, the power of the residue should be less than the power of noise; that is pR<ζ2.After the previous analysis, we use the scalar-parameter method mentioned in   Section5.1 to choose the first metastable scale. In our paper, we choose pR=1.5ζ2 (this condition is actually implicit as there is no stable scale with residue power below the noise level). Thus, after the process of scalar parameter, cartoon part of the image uC and the noise and texture parts  uNC+un are differentiated. A majority of details and textures are left in the residual part.Substituting1.5ζ2 for ζ2 into (28), we can get the residue uR by iterating (24). Substituting uR into (27), we can get the local power PR(x,y), we assign
(38)s(x,y)=ζ4PR(x,y).Now, let us analyse (38). When (x,y) lies in the large-structure area, the local power of the residue is almost constant (PR(x,y)≈ζ2), and hence s(x,y)≈ζ2, and this scheme degenerates to the scalar-parameter one mentioned in   Section 5.1. We get a high quality denoising process. When (x,y) lies in the area which contains details, the local power of the residue contains two parts: PNC(x,y)+Pn(x,y), so PR(x,y)≈PNC(x,y)+ζ2 and PR(x,y)≫ζ2, and now, s(x,y) becomes small and λ-(x,y) becomes large. We get a low quality denoising process to preserve detailed structures.Using the scalar-parameter method mentioned in   Section5.1 to choose the first metastable scale, after each step of iteration, substitute the instant image u into (35). We get the new parameter λ- of fidelity term, which is proper for the instant time scale. Next iteration is processed by substituting the new λ- into (32). Experiments show that the diffusion equation with adaptive fidelity term is better in preserving edges and details than that with scalar-parameter fidelity term.


## 6. Experiment Results and Analysis

Because “barb,”  “tiger”,  and “remote sensing image” have more detailed structures, we choose these three images as the experiment images. To highlight the treatment of the details, we intercept texture parts: the legs of BARB.We compare the processing results of six models: edge-enhancing model, coherence-enhancing model, local-tensor model without fidelity term, local tensor model with adaptive fidelity term, nonlocal tensor model without fidelity term, and nonlocal-tensor model with adaptive fidelity term. In our experiments,σ=0.5, three-layer wavelet decomposition is made, and time interval τ=0.01.Figures2, 3, and 4 are the denoising results of these six methods. (a) is the original image, (b) is the noise image, and (c) is the results of edge-enhancing model. Because it cannot get the right edge direction, the textures of legs in “woman” have disappeared and the details and fine edges in “Tiger” and remote sensing image have been blurred. (d) is the results of coherence-enhancing model. Although the textures in the leg can be vaguely seen, the details in “Tiger” and remote sensing image are even enhanced; many artificial edges appear in smooth areas. (e) is the results of local-tensor model without fidelity term, (f) is the results of local-tensor model with fidelity term, (g) is the result of nonlocal-tensor model without fidelity term, and (h) is the result of nonlocal-tensor model with fidelity term. We can notice that the models based on nonlocal tensor are better than  the ones  based on local tensor, and the models with fidelity term are better than   the ones  without fidelity term.The comparison of these methods for “barb.”
(a)
(b)
(c)
(d)
(e)
(f)
(g)
(h)The comparison of these methods for “tiger.”
(a)
(b)
(c)
(d)
(e)
(f)
(g)
(h)The comparison of these methods for “remote sensing image.”
(a)
(b)
(c)
(d)
(e)
(f)
(g)
(h)Tables1, 2, and 3 are the performance parameters of these three images, respectively. In order to emphasize the preserving of detailed structures, we choose both parts with more textures and details to process; thus, the two contradictions, preserving details and denoising, should come to a better compromise, so PSNR increases a little, but the vision effect has improved much.Table 1
The comparison of performance parameter of Figure2.
Models
PSNR
Edge-enhancing model
23.3522
Coherence-enhancing model
24.6876
Local tensor model without fidelity term
29.0260
Local tensor model with fidelity term
29.7065
Nonlocal tensor model without fidelity term
29.0712
Nonlocal tensor model with fidelity term
29.9415Table 2
The comparison of performance parameter of Figure3.
Models
PSNR
Edge-enhancing model
24.1498
Coherence-enhancing model
25.3245
Local-tensor model without fidelity term
25.5654
Local-tensor model with fidelity term
26.0946
Nonlocal-tensor model without fidelity term
25.8712
Nonlocal-tensor model with fidelity term
26.8410Table 3
The comparison of performance parameter of Figure3.
Models
PSNR
Edge-enhancing model
26.1514
Coherence-enhancing model
26.8410
Local-tensor model without fidelity term
27.1208
Local-tensor model with fidelity term
27.1335
Nonlocal-tensor model without fidelity term
27.1374
Nonlocal-tensor model with fidelity term
27.2302


## 7. Conclusions

A fidelity term, adaptively varying with the time scale, is added to Weickert’s anisotropic diffusion equation. The conflict, preserving details and denoising, will come to a better compromise. More detailed structures will be preserved during the course of denoising. Otherwise, to preserve a few other edge directions hiding in the main edge direction, Gaussian kernel is substituted for nonlinear wavelet threshold operator. In addition, in order to preserve the textures and details as much as possible, a nonlocal diffusion tensor was introduced and the two eigenvalues are reset by combining the two methods: edge-enhancing diffusion and coherence enhancing diffusion. Experimental results show that the nonlocal adaptive image denoising model has an obvious effect on preserving textures and details.



---

*Source: 605409-2013-05-22.xml*