A New Extrapolation Economy Cascadic Multigrid Method for Image Restoration Problems ()
1. Introduction
In the process of image formation, transmission, recording and processing, it is easy to be affected by equipment, environment, human factors and so on. It is very important to recover the noisy and blurred image. Therefore, it is very meaningful to construct a fast and effective image restoration algorithm, which has become one of the research hotspots in the field of image processing.
Consider the following two-dimensional image degradation model,
(1)
where
is an impulse response function, which is called the point spread function in the spatial domain,
represents the two-dimensional convolution operation,
represents the input image to be processed,
represents the additive noise, and
represents the degraded image.
Ignoring the noise
, Equation (1) can be expressed as
. (2)
In practical problems, it is usually necessary to restore the blur- and noise-contaminated image. Therefore, we usually consider the following problem,
(3)
where
represents the blur- and noise-contaminated image.
Equation (3) is an ill-posed problem. In the case of direct solution, a small disturbance of the right-hand side may lead to great changes, so the recovered image is usually of low quality. The Tikhonov regularization method can eliminate this instability by solving the following model
(4)
where the first term is the fidelity term,
is the regularization term, also known as constraints, and α is the regularization parameter, which is used to balance the fidelity term and regularization term. The discrete regularization term is denoted by
, and the selection of the regularization operator
is crucial for the quality of image restoration. The most representative one is the nonlinear anisotropic diffusion equation proposed by Perona and Malik in [1] , and the regularization operator is expressed as follows
(5)
the diffusion coefficients are chosen as follows,
(6)
where k is the threshold.
For nonlinear Equation (4), [2] and [3] propose the explicit and semi-implicit discretization schemes, which have good preserving effect on the edges of the image and can provide high quality restoration of images. However, whether it is the explicit or semi-implicit discretization schemes its computation cost is very large. In addition, the selection of the regularization parameters and the time intervals are not easy.
Equation (2) is discretized to obtain the following linear system
(7)
where H represents the fuzzy operator, which is usually a real symmetric block Teoplitiz matrix, u is the original image, and g represents the blur and noise-free image.
The blur- and noise-contaminated image g is obtained by adding random normally distributed noise e with mean 0 and standard deviation 1, such as
(8)
Thus, Equation (3) is discretized to obtain the following linear system
(9)
we set the solution is
.
For the linear Equation (9), the matrix H is usually a large-scale singular matrix. The system of linear equations is usually an ill-posed problem. Although we can solve Equation (9) by some iterative method, due to the loss of high-frequency information, these iterative methods may produce ringing effect and cannot accurately recover edges, which seriously affects the quality of image restoration.
In order to solve the problems lots of effective methods have been proposed.
In [4] , a standard iterative method is proposed to compute several approximate solutions and then extract a new candidate solution from the linear subspace expanded by these approximate solutions. This method can be used for large-scale problems. In [5] , LSQR algorithm and RRGMRES algorithm are proposed to solve the large-scale linear ill-posed problems, and corresponding stopping criteria are proposed, which can obtain meaningful approximate solutions. In [6] , the regularization properties of the global GMRES method are analyzed, and a regularized global GMRES method is proposed for solving linear discrete ill-posed problems with multiple right-hand sides contaminated by errors. In [7] , NL-means algorithm is proposed to correct the ladder effect.
Multigrid method is one of the most effective methods for solving elliptic boundary value problems. The cascadic multigrid method is proposed in [8] . The advantage of this method lies in its simplicity and efficiency. In [9] presents an economical cascadic multigrid method. Compared with the conventional cascadic multigrid method, the proposed method can get the same efficient solutions with fewer iterations. [10] [11] [12] have proposed some extrapolation cascadic multigrid methods. The methods use the new extrapolation with quadratic interpolation, which can provide a better initial value for the fine grid level.
Based on the advantages of these methods, many scholars have applied them to image restoration. [13] [14] have proposed multigrid method for image restoration. The step effect is reduced by introducing wavelet soft threshold denoising and smoothing. In [15] , a cascadic multilevel method is proposed to solve a linear discrete ill-posed problem with error contamination on the right side. Numerical experiments show that the algorithm is effective in image restoration. In [16] , a class of nonlinear edge-preserving continuation operators are proposed to be used as prolongation operators for cascadic multigrid methods. This method can reduce the ringing effect and restore good image edges. In [17] , the selection of fuzzy operator and the smooth method are further studied, and it is shown that the cascadic multilevel method has good efficiency for image restoration.
Combining the advantages of the new extrapolation method and the cascadic multigrid method, a new extrapolation economy cascadic multigrid method is proposed for image restoration in this paper. The experimental results show that this method not only improves computational efficiency but also ensures good recovery quality. The new method has the following advantages:
1) The new extrapolation formula and the quadratic interpolation are combined to provide more accurate values for the fine mesh level.
2) Fewer iterations are costed on each grid level without affecting the accuracy.
3) The edge preserving denoising operator can preserve the edge of image while denoising.
4) The local smoothing operator reduces the ladder effect and ensures good recovery quality.
The structure of this paper is as follows. In Section 2, introduces the basic framework of the new extrapolation economy cascadic multigrid. In Section 3, the linear prolongation operator, the new extrapolation formula prolongation operator combined with quadratic interpolation and the smoothing operator generated by the weighted least squares approximation is introduced. In Section 4, the edge preserving denoising operator is introduced. In Section 5, some numerical results are given to illustrate the effectiveness of the new method. In Section 6, a summary is presented.
2. New Extrapolation Economy Cascadic Multigrid Method
Selecting a Template
Let
, we define the following norm
(10)
For the Equation (9), we can solve it by using the conjugate gradient method (CG) or the minimal residual method (MR). Let the initial iterative value
, here l represents the grid level. We use the following stopping rule to terminate the iteration.
Let
, for a given
, and
is a constant and independent of
, when
, stop the iteration.
By using CG or MR to solve the linear model directly, although it cannot effectively restore the edge which seriously affects the quality of image restoration, it is simple and the time is short. In view of these two advantages, we can use CG or MR as the smoothing operator in the new extrapolation economic cascadic multigrid method.
For the sake of discussion, we set u as an
column vector and
,
is a nested subspace. The subscript l indicates the number of the grid level,
denotes the coarsest grid level, and
indicates the finest grid level. Accordingly, we have
. When N is odd,
, when N is even,
.
In the cascadic multigrid algorithm, the following linear models need to be solved iteratively at each grid level
(11)
where
is the representation of the fuzzy operator H in the
.
Set
be the kth element of
. Let
(12)
We define the restriction operator
, which satisfies
(13)
When
,
.
is generated in the following way,
(14)
where
is the adjoint of
, when
,
.
In [16] , a cascadic multiresolution method is proposed. Numerical examples show that the method has good performance in image restoration. It is explained that for Equation (9), lots of iterations will affect the recovery effect. Therefore, the number of iterations is particularly critical, and the cascadic multigrid method proposed by [9] is of reference significance. We combine the advantages of the algorithms of [9] and [16] propose the following improved economic cascadic multigrid method.
Algorithm 1: Improved Economical Cascadic Multigrid Method for Image Restoration (IECMG)
Step 1: Set
,
.
Step 2: When
, do:
1) Prolongation:
or
.
2) Local smoothing:
.
3) Preserving edge and removing noise:
.
4) Smoothing:
,
.
Step 3:
.
In [9] [10] [11] [12] , the economic cascadic multigrid method and the new extrapolation cascadic multigrid method are proposed respectively. The new extrapolation interpolation method can provide better initial values for the next level and has better convergence. Therefore, we propose the following new extrapolation economy cascadic multigrid method.
Algorithm 2: New Extrapolation Economy Cascadic Multigrid method for Image Restoration (EECMG)
Step 1: Set
,
,
.
Step 2: When
, do:
1) New extrapolation:
.
2) Quadratic interpolation prolongation:
.
3) Local smoothing:
.
4) Preserving edge and removing noise:
.
5) Smoothing: Let
,
Step 3:
The number of iterations
on each level is selected as follows [9] .
1) If
, then
.
2) If
, then
Where
,
,
,
,
,
,
,
represents the smallest positive integer not less than t.
Note: In the above algorithms,
denotes smoothing
times by using CG or MR. In Sections 3 and Sections 4, we will introduce the linear prolongation operator
, nonlinear prolongation operator
, new extrapolation operator
, local smoothing operator
and edge preserving denoising operator
in Algorithm 1 and Algorithm 2 in detail.
3. Prolongation Operator and Local Smoothing Operator
3.1. Linear Prolongation Operator
Let
denote the pixel value in
on the level l. The linear operator
can be defined as follows,
(15)
Because the accuracy of the image on the fine grid level by piecewise linear interpolation is not high, a nonlinear prolongation operator, combined a new extrapolation formula with quadratic interpolation, is constructed to solve this problem.
3.2. New Extrapolation Quadratic Interpolation Prolongation Operator
Let us take the one-dimensional triplet grid
as an example, the pixel is denoted by
, and the corresponding pixel value is denoted by
. Let
be the coarsest grid level and contains pixels
. The pixels contained in
and
are represented as follows,
(16)
the corresponding pixel value are expressed as,
(17)
We combine the new extrapolation formula with quadratic interpolation, which is called the new extrapolation quadratic interpolation method. This can provide better initial values on the fine grid level
. The details are as follows (see Figure 1).
1) New extrapolation
,
2) Quadratic interpolation
,
Next, we apply this idea to the two-dimension problems. As shown in Figure 2, the hollow circles represent the pixels on the
, the solid black dots represent the pixel values on the
grid level, and the rectangular boxes represent the pixel values on the
grid level. We use the new extrapolation formula to get the image value at the black node and the quadratic interpolation to obtian the image value at the rectangular box.
Figure 1. New extrapolation quadratic interpolation method for the 1-dimensional.
Figure 2. New extrapolation quadratic interpolation method for the 2-dimensional.
1) New extrapolation
,
2) Quadratic interpolation
,
3.3. Local Smoothing Operator
In order to get a smoother image, the local weighted least squares approximation is used to smooth the image in [2] . As shown in Figure 3, the pixel value
is related to the image value of the surrounding eight points. The process is represented by the operator
as follows.
First, we introduce the weighting function,
(18)
then we solve the locally weighted least squares problem,
(19)
Set
represent the solutions of the above problem, then the updated pixel value is
.
4. Edge Preserving Denoising Operator
The isotropic diffusion equation diffuses equally along the tangential direction and the normal direction of the image edge, so it cannot protect the edge, nor can it preserve the original fine structure of the image well, so that the image becomes blurred. In order to deal with these shortcomings, an anisotropic diffusion equation is proposed. One is the nonlinear anisotropic diffusion equation proposed by Perona and Malik, abbreviated as P-M equation [1] , which is as follows,
(20)
where the
is the image of the t moment,
is the diffusion coefficient. The more classical diffusion coefficient is as follows,
(21)
Set
, discretization of (20) is as follows,
(22)
Taking
, we have,
(23)
where
denotes the set of four neighboring points centered on
. So according to the above equation we have,
(24)
Equation (24) can be expressed in matrix form as,
(25)
where the elements of
are as follows,
(26)
5. Numerical Examples
We use numerical examples to illustrate the effectiveness of the new extrapolation economy cascadic multigrid method for restoring blur- and noise-contaminated image. The elements in the fuzzy operator H are defined by the following function.
where
denotes the Kronecker product. σ represents the variance of the Gaussian point diffusion function. The “band” represents the half bandwidth of h. The degree of blurring is denoted by
, and the following cases are considered in this paper,
Defining the noise levels as follows,
the noise level v is set to the following values,
In this experiment, the blur and noise level of the image is denoted as
(
). For example,
denotes that
,
,
.
The image restoration effect is measured by the following peak signal ratio,
where u represents blur- and noise-free image,
represents the recovered image.
Figure 4(a) is the original image, Figures 4(b)-(e) are the images with four different levels of blur and noise.
Firstly, CG algorithm and MR algorithm are used for image restoration. When the stopping rule is met, the iteration times of both algorithms are 4. The restoration effect of four different levels of blur- and noise-contaminated images are shown in Figure 5.
Table 1 and Figure 5 show that the MR algorithm has a faster computational speed than the CG algorithm in the same case. For b1v1 and b2v2, the recovery effect of CG and MR is almost the same. For b3v3 and b4v4, the MR algorithm has better recovery effect.
Table 1. Numerical results of CG and MR.
In the following experiments, with different smoothers, Algorithm 1 and Algorithm 2 are abbreviated as shown in the Table 2. Algorithm 1 and Algorithm 2 with different smoother. Next, we compare the numerical results of restoration Algorithm 1 and Algorithm 2.
From Table 3 and Table 4, Figure 6 and Figure 7, it can be found that,
1) The quadratic interpolation can obtain better accuracy initial value. So IECMG-P is better than IECMG-L.
2) New extrapolation with quadratic interpolation can improve the accuracy of the initial values. So EECMG is better than IECMG.
3) With MR as smoother, EECMG (MR) and IECMG (MR) have the best restoration quality and spend the least recovery time.
Moreover, Figure 5, Figure 8 and Figure 9 show that EECMG (MR) is best, whether PSNR or CPU time.
Finally, we show the edge preserving ability of the new extrapolation cascadic multigrid method.
Figure 10(a) is the original image, Figure 10(b) is the third-level polluted image. Figure 10(c) is the restored image by using EECMG (MR), which PSNR = 27.9617. The good effect of edge preservation is shown in Figure 10(d).
Table 2. Algorithm 1 and Algorithm 2 with different smoother.
(a)(b)
Figure 6. PSNR and Time of IECMG (CG) and EECMG (CG). (a) PSNR; (b) Time.
Table 3. Numerical results of IECMG (CG) and EECMG (CG).
(a)(b)
Figure 7. PSNR and TIME of IECMG (MR) and EECMG (MR). (a) PSNR; (b) Time.
Table 4. Numerical results of IECMG (MR) and EECMG (MR).
(a)(b)
Figure 9. Numerical results of CG, MR, IECMG and EECMG. (a) PSNR; (b) Time.
6. Conclusions
Combining the advantages of the new extrapolation method and the economical cascadic multigrid method, this paper proposes a new extrapolation economical cascadic multigrid method for image restoration, which not only improves the computational efficiency but also ensures the good quality of restoration. Numerical experiments show that,
1) The edge preserving denoising operator can achieve the effect of preserving image edges and removing noise;
2) The weighted least squares approximation smooths the image locally and reduces the staircase effect;
3) The new extrapolation formula and quadratic interpolation construct a nonlinear continuation operator, which can provide more accurate pixel values for the fine grid level and ensure good recovery quality;
4) The iteration number selection method controls the iteration number of each grid level, which greatly reduces the computational work.
Although this paper puts forward several kinds of algorithms for linear model of image restoration problem and has made some research achievements, there are still many works worth further research, such as:
1) Whether the algorithm proposed in this paper can also have good restoration effect on other image restoration models?
2) Whether the algorithm in this paper can restore large-size images?
3) Can the algorithm in this paper be applied to the restoration of three-dimensional images?
Acknowledgements
This work was supported by Natural Science Foundation of China (11661027), the Guangxi Natural Science Foundation (2020GXNSFAA159143), and National Natural Science Foundation of China (12161027).