Alternating method based on framelet l₀-norm and TV regularization for image restoration

ABSTRACT

This paper presents an efficient alternating method for image deblurring and denoising, it is based on our new model with TV norm to denoise the deblurred image which is penalized by $l_0$ norm of the framelet transform. Our combinational algorithm performs deblurring and denoising alternately. We apply the proximity algorithm to solve the TV regularization for denoising part, and the mean doubly augmented Lagrangian (MDAL) method is used to solve $l_0$ minimization in the analysis-based sparsity for the deblurring part. Experimental results show that the proposed alternating minimization method is robust for a different type of blur and noise. We also do some comparisons with related existing methods to demonstrate that our method has a more significant effect.

KEYWORDS:

• Image restoration
• augmented lagrangian method
• $l_0$ regularization
• framelet transform
• proximity operator

MSC SUBJECT CLASSIFICATIONS:

• 65N22
• 68U10
• 65N21

1. Introduction

Digital image restoration is an important problem in image processing which is widely applied to medical and astronomical imaging, film restoration, image and video coding [1–8] and so on. The image restoration task can be formulated as a linear inverse problem: (1) $f=Au+ϵ,$(1) where $f\in {\mathit{R}}^{{\mathit{N}}^{\mathit{2}}}$ is the observed degraded image, $A\in {\mathit{R}}^{{\mathit{N}}^{\mathit{2}}\times {\mathit{N}}^{\mathit{2}}}$ represents a spatial-invariant blur matrix which is constructed by the blur kernel, $u\in {\mathit{R}}^{{\mathit{N}}^{\mathit{2}}}$ is the original image to be restored, $ϵ\in {\mathit{R}}^{{\mathit{N}}^{\mathit{2}}}$ denotes an additive Gaussian white noise with some kind of deviation σ.

To obtain a reasonable approximated solution u of this problem, various regularization-based methods have been proposed in the past decades. One of the most popular is the Rudin–Osher–Fatemi (ROF) model with total variation (TV) regularization, proposed by Rudin, Osher and Fatemi [9]. A discrete version of the unconstrained TV regularization is given by (2) $\underset{u}{min}\frac{1}{2}\parallel Au-f{\parallel }_2^2+\alpha \parallel u{\parallel }_{\mathrm{TV}},$(2) where $\parallel \cdot {\parallel }_2$ denotes the Euclidean norm, α is a positive regularization parameter which measures the trade-off between a good fit and a regularized solution, and $\parallel \cdot {\parallel }_{TV}$ is the TV regularization term. When A is an identity matrix, the ROF model produces the denoising case. The model is initially solved by the partial differential equation, which has high computational complexity.

Difficulty in minimizing the functional appearing in (2(2) $\underset{u}{min}\frac{1}{2}\parallel Au-f{\parallel }_2^2+\alpha \parallel u{\parallel }_{\mathrm{TV}},$(2) ) lies in the nondifferentiability of the total variation semi-norm and the large dimension of the underlying images. A number of numerical methods have been proposed to address these issues. The most popular is direct optimization [9–11] by solving the associated Euler-Lagrangian equations [12–16]. Other methods including subgradient descents, second-order cone programming, interior point methods, graph-based approaches, Nesterov's first-order explicit schemes [17–21], the split Bregman iteration methods and its variations [22–25], which introduces an intermediate variable to make the ROF model separable so that the minimization of the resulting functional is easier to solve.

Micchelli, Shen and Xu [26] proposed the proximity algorithm to solve the ROF model for image denosing which is simple and efficient. The total-variation can be viewed as the composition of a convex function (the $l_1$ norm for the anisotropic total-variation or the $l_2$ norm for the isotropic total-variation) with a linear transformation, then use the fixed-point algorithm for computing the proximity operator of the composition of the convex function with a linear transformation.

A fast TV minimization method is proposed by Huang, Ng and Wen [27] for image restoration. They modified the above unconstrained TV problem (2(2) $\underset{u}{min}\frac{1}{2}\parallel Au-f{\parallel }_2^2+\alpha \parallel u{\parallel }_{\mathrm{TV}},$(2) ) as follow: (3) $\begin{array}{rl}& \underset{u,v}{min}\frac{1}{2}\parallel Au-f{\parallel }_2^2+{\alpha }_1\parallel u-v{\parallel }_2^2+{\alpha }_2\parallel v{\parallel }_{\mathrm{TV}}\\ & =\underset{v}{min}\underset{u}{min}\{\parallel Au-f{\parallel }_2^2+{\alpha }_1\parallel u-v{\parallel }_2^2\}+{\alpha }_2\parallel v{\parallel }_{\mathrm{TV}},\end{array}$(3) where ${\alpha }_1,{\alpha }_2$ are positive regularization parameters. In (3(3) $\begin{array}{rl}& \underset{u,v}{min}\frac{1}{2}\parallel Au-f{\parallel }_2^2+{\alpha }_1\parallel u-v{\parallel }_2^2+{\alpha }_2\parallel v{\parallel }_{\mathrm{TV}}\\ & =\underset{v}{min}\underset{u}{min}\{\parallel Au-f{\parallel }_2^2+{\alpha }_1\parallel u-v{\parallel }_2^2\}+{\alpha }_2\parallel v{\parallel }_{\mathrm{TV}},\end{array}$(3) ), we can interpret the TV minimization scheme to denoise the deblurred image u. The TV norm for the image restoration process is able to preserve edges very well.

In [28], Dong and Zhang proposed a framelet-based method which the image term u had a sparse representation in some redundant transform domain, and penalized by $l_0$ norm which is widely researched in [29–33]. Synthesis-based sparsity and analysis-based sparsity priors are two famous sparse approximations. Because analysis-based approach has a smoother solution than that of synthesis-based approach, the focus of the paper is the analysis-based approach, the detail information about the two sparse representation can be referred to [34,35].

It has been proved in [27] that it is a good idea to bring in auxiliary variable v to denoise and deblur image alternately and use the TV minimization scheme in denoising step. Since the $l_0$ regularization based on the framelet transform has the outstanding performance on image deblurring, we apply the analysis-based approach for the deblurring step. Motivated by the mentality of [27] and the excellent effect of [28], we propose a new optimization model: (4) $\underset{u,v}{min}\mathit{J}\mathit{(}\mathit{u}\mathit{,}\mathit{v}\mathit{)}\equiv \underset{\mathit{u}\mathit{,}\mathit{v}}{min}\frac{\mathit{1}}{\mathit{2}}\parallel \mathit{A}\mathit{u}\mathit{-}\mathit{f}{\parallel }_{\mathit{2}}^{\mathit{2}}\mathit{+}\parallel \lambda \cdot \mathit{W}\mathit{u}{\parallel }_{\mathit{0}}\mathit{+}\frac{\alpha }{\mathit{2}}\parallel \mathit{u}\mathit{-}\mathit{v}{\parallel }_{\mathit{2}}^{\mathit{2}}\mathit{+}{\lambda }_{\mathit{1}}\parallel \mathit{v}{\parallel }_{\mathrm{TV}}\mathit{,}$(4) where W is a matrix based on framelet transform, this framelet form tight wavelet frame, detail description can be found in [36], $\lambda ,\alpha ,{\lambda }_1$ are positive regularization parameters, α measures the trade-off between a deblurred image u and a denoised image v, λ and ${\lambda }_1$ measure the amount of regularization to a deblurring image u and a denoising image v, respectively.

The first advantage of this paper is the new model which alternates deblurring and denoising, our new model bases on TV norm and $l_0$ norm in framelet transform is superior than relative models. The second advantage is designing a combinational algorithm which consists of the mean doubly augmented Lagrangian (MDAL) method and the proximity algorithm. The MDAL method is to solve the $l_0$-norm in framelet transform model (4(4) $\underset{u,v}{min}\mathit{J}\mathit{(}\mathit{u}\mathit{,}\mathit{v}\mathit{)}\equiv \underset{\mathit{u}\mathit{,}\mathit{v}}{min}\frac{\mathit{1}}{\mathit{2}}\parallel \mathit{A}\mathit{u}\mathit{-}\mathit{f}{\parallel }_{\mathit{2}}^{\mathit{2}}\mathit{+}\parallel \lambda \cdot \mathit{W}\mathit{u}{\parallel }_{\mathit{0}}\mathit{+}\frac{\alpha }{\mathit{2}}\parallel \mathit{u}\mathit{-}\mathit{v}{\parallel }_{\mathit{2}}^{\mathit{2}}\mathit{+}{\lambda }_{\mathit{1}}\parallel \mathit{v}{\parallel }_{\mathrm{TV}}\mathit{,}$(4) ). Because the iterative sequence $(u^k,d^k,b^k)$ obtained by augmented Lagrangian (AL) method can not guarantee the convergence of $(u^k,d^k)$, and the doubly augmented Lagrangian (DAL) method may converge slowly and unsteadily, this motivated to consider the mean values of ${\overline{u}}^k=1/(k+1)\sum _{j=0}^k{u}^k$ and ${\overline{d}}^k=1/(k+1)\sum _{j=0}^k{\mathrm{d}}^k$ as the outputs of the deblurring step. The proposed algorithm is more effective for image restoration, since it can improve the deblurring and denoising quality, respectively.

The paper is organized as following. In Section 2, we show the detailed implementation of our combinational algorithm based on the MDAL method and proximity algorithm for solving model (4(4) $\underset{u,v}{min}\mathit{J}\mathit{(}\mathit{u}\mathit{,}\mathit{v}\mathit{)}\equiv \underset{\mathit{u}\mathit{,}\mathit{v}}{min}\frac{\mathit{1}}{\mathit{2}}\parallel \mathit{A}\mathit{u}\mathit{-}\mathit{f}{\parallel }_{\mathit{2}}^{\mathit{2}}\mathit{+}\parallel \lambda \cdot \mathit{W}\mathit{u}{\parallel }_{\mathit{0}}\mathit{+}\frac{\alpha }{\mathit{2}}\parallel \mathit{u}\mathit{-}\mathit{v}{\parallel }_{\mathit{2}}^{\mathit{2}}\mathit{+}{\lambda }_{\mathit{1}}\parallel \mathit{v}{\parallel }_{\mathrm{TV}}\mathit{,}$(4) ), Convergence analysis is given in Section 3, Experimental results are exhibited in Section 4, A brief summary is given in Section 5.

2. Alternating minimization method for deblurring and denoising

In this section, we propose an efficient alternating method for our new model (4(4) $\underset{u,v}{min}\mathit{J}\mathit{(}\mathit{u}\mathit{,}\mathit{v}\mathit{)}\equiv \underset{\mathit{u}\mathit{,}\mathit{v}}{min}\frac{\mathit{1}}{\mathit{2}}\parallel \mathit{A}\mathit{u}\mathit{-}\mathit{f}{\parallel }_{\mathit{2}}^{\mathit{2}}\mathit{+}\parallel \lambda \cdot \mathit{W}\mathit{u}{\parallel }_{\mathit{0}}\mathit{+}\frac{\alpha }{\mathit{2}}\parallel \mathit{u}\mathit{-}\mathit{v}{\parallel }_{\mathit{2}}^{\mathit{2}}\mathit{+}{\lambda }_{\mathit{1}}\parallel \mathit{v}{\parallel }_{\mathrm{TV}}\mathit{,}$(4) ). Deblurring alternates with Denoising. We use MDAL to solve $l_0$ regularization based on the framelet transform for deblurring and employ proximity algorithm based on total variation minimization to denoise the deblurred images, that is deblurring alternates with denoising.

In (4(4) $\underset{u,v}{min}\mathit{J}\mathit{(}\mathit{u}\mathit{,}\mathit{v}\mathit{)}\equiv \underset{\mathit{u}\mathit{,}\mathit{v}}{min}\frac{\mathit{1}}{\mathit{2}}\parallel \mathit{A}\mathit{u}\mathit{-}\mathit{f}{\parallel }_{\mathit{2}}^{\mathit{2}}\mathit{+}\parallel \lambda \cdot \mathit{W}\mathit{u}{\parallel }_{\mathit{0}}\mathit{+}\frac{\alpha }{\mathit{2}}\parallel \mathit{u}\mathit{-}\mathit{v}{\parallel }_{\mathit{2}}^{\mathit{2}}\mathit{+}{\lambda }_{\mathit{1}}\parallel \mathit{v}{\parallel }_{\mathrm{TV}}\mathit{,}$(4) ), there are two unknown variables, one is the deblurred image u, and the other is the denoised image v. Starting from an initial guess $v^{\left(0\right)}$, we can get a series of $u^{\left(n\right)},v^{\left(n\right)}$: $u^{\left(1\right)},v^{\left(1\right)},u^{\left(2\right)},v^{\left(2\right)},\dots ,u^{\left(i\right)},v^{\left(i\right)},\dots$ such that (5) $u^{\left(i\right)}=\arg \underset{u}{min}\frac{1}{2}\parallel Au-f{\parallel }_2^2+\frac{\alpha }{2}\parallel u-v^{(i-1)}{\parallel }_2^2+\parallel \lambda \cdot Wu{\parallel }_0,$(5) (6) $v^{\left(i\right)}=\arg \underset{v}{min}\frac{\alpha }{2}\parallel v-u^{\left(i\right)}{\parallel }_2^2+{\lambda }_1\parallel v{\parallel }_{\mathrm{TV}},$(6) with $i=1,2,\dots n$.

2.1. Deblurring step

The first step of the alternating method is deblurring process which solves the $l_0$-norm-based minimization of (5(5) $u^{\left(i\right)}=\arg \underset{u}{min}\frac{1}{2}\parallel Au-f{\parallel }_2^2+\frac{\alpha }{2}\parallel u-v^{(i-1)}{\parallel }_2^2+\parallel \lambda \cdot Wu{\parallel }_0,$(5) ): (7) $\underset{u}{min}\frac{1}{2}\parallel Au-f{\parallel }_2^2+\frac{\alpha }{2}\parallel u-v^{(i-1)}{\parallel }_2^2+\parallel \lambda \cdot Wu{\parallel }_0.$(7) Here we take the effective replacement Wu=d, rewrite minimization (7(7) $\underset{u}{min}\frac{1}{2}\parallel Au-f{\parallel }_2^2+\frac{\alpha }{2}\parallel u-v^{(i-1)}{\parallel }_2^2+\parallel \lambda \cdot Wu{\parallel }_0.$(7) ) as a constrained optimization problem: (8) $\underset{u,d}{min}\frac{1}{2}\parallel Au-f{\parallel }_2^2+\frac{\alpha }{2}\parallel u-v^{(i-1)}{\parallel }_2^2+\parallel \lambda \cdot d{\parallel }_{0}s.t.d=Wu.$(8) The alternating direction method of multipliers (ADMM) [37] is applied to the augmented Lagrangian (AL) [38] of (7(7) $\underset{u}{min}\frac{1}{2}\parallel Au-f{\parallel }_2^2+\frac{\alpha }{2}\parallel u-v^{(i-1)}{\parallel }_2^2+\parallel \lambda \cdot Wu{\parallel }_0.$(7) ): (9) $\begin{array}{rl}\underset{u,d}{min}\mathcal{S}(u,d,b)& =\underset{u,d}{min}\frac{1}{2}\parallel Au-f{\parallel }_2^2+\frac{\alpha }{2}\parallel u-v^{(i-1)}{\parallel }_2^2+\parallel \lambda \cdot d{\parallel }_0+⟨b,Wu-d⟩\\ & +\frac{\mu }{2}\parallel Wu-d{\parallel }_2^2,\end{array}$(9) where $\mu >0$ controls the weight of the penalty term.

The AL method solving (7(7) $\underset{u}{min}\frac{1}{2}\parallel Au-f{\parallel }_2^2+\frac{\alpha }{2}\parallel u-v^{(i-1)}{\parallel }_2^2+\parallel \lambda \cdot Wu{\parallel }_0.$(7) ) can be formulated as iterative steps: (10) $\begin{array}{rl}& (u^{k+1},d^{k+1})=\arg \underset{u,d}{min}\mathcal{S}(u,d,b^k),\\ & b^{k+1}=b^k+(W{u}^{k+1}-d^{k+1}).\end{array}$(10) The AL method can only guarantee the convergence of the dual variable $b^k$ but not the primal variables $(u^k,d^k)$, so the doubly augmented Lagrangian (DAL) is introduced for (7(7) $\underset{u}{min}\frac{1}{2}\parallel Au-f{\parallel }_2^2+\frac{\alpha }{2}\parallel u-v^{(i-1)}{\parallel }_2^2+\parallel \lambda \cdot Wu{\parallel }_0.$(7) ), that is we rewrite the minimizing problem as: (11) $\tilde{\mathcal{S}}(u,d,b,\tilde{u},\tilde{d})=\mathcal{S}(u,d,b)+\frac{\gamma }{2}\parallel u-\tilde{u}{\parallel }_2^2+\frac{\gamma }{2}\parallel d-\tilde{d}{\parallel }_2^2.$(11) where $\mathcal{S}(u,d,b)$ is the functional in (9(9) $\begin{array}{rl}\underset{u,d}{min}\mathcal{S}(u,d,b)& =\underset{u,d}{min}\frac{1}{2}\parallel Au-f{\parallel }_2^2+\frac{\alpha }{2}\parallel u-v^{(i-1)}{\parallel }_2^2+\parallel \lambda \cdot d{\parallel }_0+⟨b,Wu-d⟩\\ & +\frac{\mu }{2}\parallel Wu-d{\parallel }_2^2,\end{array}$(9) ), the iterative scheme of the DAL method for (11(11) $\tilde{\mathcal{S}}(u,d,b,\tilde{u},\tilde{d})=\mathcal{S}(u,d,b)+\frac{\gamma }{2}\parallel u-\tilde{u}{\parallel }_2^2+\frac{\gamma }{2}\parallel d-\tilde{d}{\parallel }_2^2.$(11) ) consists of the following three subproblems: (12) $\begin{array}{rl}{u}^{(k+1,i)}& =\arg \underset{u}{min}\frac{1}{2}\parallel Au-f{\parallel }_2^2+\frac{\alpha }{2}\parallel u-v^{(i-1)}{\parallel }_2^2\\ & +\frac{\mu }{2}\parallel Wu-d^{(k,i)}+b^{(k,i)}{\parallel }_2^2+\frac{\gamma }{2}\parallel u-u^{(k,i)}{\parallel }_2^2,\end{array}$(12) (13) $\begin{array}{rl}{d}^{(k+1,i)}& =\arg \underset{d}{min}\parallel \lambda \cdot d{\parallel }_0+\frac{\mu }{2}\parallel d-(W{u}^{(k+1,i)}+b^{(k,i)}){\parallel }_2^2\\ & +\frac{\gamma }{2}\parallel d-d^{(k,i)}{\parallel }_2^2,\end{array}$(13) (14) $b^{(k+1,i)}=b^{(k,i)}+(W{u}^{(k+1,i)}-d^{(k+1,i)}).$(14) As the subproblem of solving u is differentiable, we solve the following positive definite linear system: (15) $(A^{T}A+\alpha +\mu +\gamma )u=A^{T}f+\alpha v^{(i-1)}+{\mu }_1{W}^T(d^{(k,i)}-b^{(k,i)})+\gamma u^{(k,i)},$(15) where $A^T$ is the conjugate transpose of A.

Then we use fast fourier transform (FFT) to perform the main computation: $u^{(k+1,i)}={\mathrm{FFT}}^{-1}\left(\frac{\mathrm{FFT}(A^{T}f+\alpha v^{(i-1)}+\mu W^T(d^{(k,i)}-b^{(k,i)})+\gamma u^{(k,i)})}{\mathrm{FFT}(A^{T}A+\alpha +\mu +\gamma )}\right),$ here A is usually a matrix of block circulant with circulant blocks when periodic boundary conditions.

The solution of d subproblem can be described as follows: (16) $d^{(k+1,i)}={\mathcal{H}}_{\lambda ,\mu ,\gamma }(W{u}^{(k+1,i)}+b^{(k,i)},d^{(k,i)}),$(16) where the operator $\mathcal{H}$ is a generalized hard-thresholding operator defined componentwisely as (17) $\left({\mathcal{H}}_{\lambda ,\mu ,\gamma }\right(x,y){)}_i=\left\{\begin{array}{ll}0,& if|\frac{\mu x_i+\gamma y_i}{\mu +\gamma }|<\sqrt{\frac{2{\lambda }_i}{\mu +\gamma }},\\ \frac{\mu x_i+\gamma y_i}{\mu +\gamma },& \text{otherwise.}\end{array}\right.$(17) However, the DAL method may converge slowly and unsteadily, so we use the mean values as the outputs of deblurring, instead the sequence $(u^{(k,i)},d^{(k,i)})$ as $({\overline{u}}^{(k,i)},{\overline{d}}^{(k,i)})$, where ${\overline{u}}^{k,i}=1/(k+1)\sum _{j=0}^k{u}^{(k,i)}$ and ${\overline{d}}^{k,i}=1/(k+1)\sum _{j=0}^k{d}^{(k,i)}$. The new sequence $({\overline{u}}^{(k,i)},{\overline{d}}^{(k,i)})$ does well converge and satisfies the constraint Wu=d asymptotically, this is the mean doubly augmented Lagrangian method.

2.2. Denoising step

The second step is to denoise the image u generated by the previous deblurring step, this is the minimizer of the optimization problem (6(6) $v^{\left(i\right)}=\arg \underset{v}{min}\frac{\alpha }{2}\parallel v-u^{\left(i\right)}{\parallel }_2^2+{\lambda }_1\parallel v{\parallel }_{\mathrm{TV}},$(6) ): (18) $\underset{v}{min}\frac{\alpha }{2}\parallel v-u^{\left(i\right)}{\parallel }_2^2+{\lambda }_1\parallel v{\parallel }_{TV}.$(18) Here we employ the proximity algorithm in this denoising step. Total variation can be viewed as the composition of a convex functional,the $l_1$ norm for the anisotropic total-variation or the $l_2$ norm for the isotropic total-variation. We use fixed-point algorithm for computing the proximity operator, this method is mathematically sound simple and efficient for solving (19(18) $\underset{v}{min}\frac{\alpha }{2}\parallel v-u^{\left(i\right)}{\parallel }_2^2+{\lambda }_1\parallel v{\parallel }_{TV}.$(18) ).

So we consider the following isotropic total-variation problem: (19) $\underset{v}{min}\frac{\alpha }{2}\parallel v-u^{\left(i\right)}{\parallel }_2^2+\varphi \circ B\left(v\right),$(19) where we denote $\varphi ={\lambda }_1\parallel \cdot {\parallel }_2$. The vector that achieves the minimum in (20(19) $\underset{v}{min}\frac{\alpha }{2}\parallel v-u^{\left(i\right)}{\parallel }_2^2+\varphi \circ B\left(v\right),$(19) ), denoted by $pro{x}_{\varphi \circ B}{u}^{\left(i\right)}$, exists and is unique. for convenience we introduce the affine transformation $T:R^m\to R^m$ defined, for all $y\in R^m$ by $Ty:=Bx+(I-{\mu }_{1}B{B}^t)y,$ and the operator $H:R^m\to R^m$ given, for a fixed $x\in R^d$ by $H:=(I-pro{x}_{1/\left({\mu }_1\right)\varphi })\circ T,$ where ${\mu }_1>0$ satisfies specified condition.

Micchelli, Shen and Xu [26] have proved that the solution of the minimization problem (20(19) $\underset{v}{min}\frac{\alpha }{2}\parallel v-u^{\left(i\right)}{\parallel }_2^2+\varphi \circ B\left(v\right),$(19) ) is just a fixed point of the operator H which involves the proximity operator of the convex function $1/{\mu }_1\varphi$ and the linear transformation B. Based on the above analysis, we focus on finding the fix-points of the operator H.

Now we consider the Picard iterates of H. For a given initial point $v^0\in R^m$, the Picard sequence of H is defined by $v^{n+1}:=H\left(v^n\right),n\in N.$ The k-averaged operator $H_k$ of H with $k\in (0,1)$ is defined as following: (20) $H_k:=kI+(1-k)H.$(20) To speed up convergence of the iteration, we can also use Gaussian-Seidel variation of the fixed-point algorithm based on the proximity operator for model (20(19) $\underset{v}{min}\frac{\alpha }{2}\parallel v-u^{\left(i\right)}{\parallel }_2^2+\varphi \circ B\left(v\right),$(19) ).

Let D denote the $N\times N$ matrix defined by the following: $D:=\left(\begin{array}{cccc}0& & & \\ -1& 1& & \\ & \ddots & \ddots & \\ & & -1& 1\end{array}\right),$ and choose B to be a $2N^2\times N^2$ matrix given by (21) $B:=\left(\begin{array}{c}{I}_N\otimes D\\ D\otimes I_N\end{array}\right),$(21) where $I_N$ is the $N\times N$ identity matrix and the notation ⊗ denotes the Kronecker product.

To present a Gauss-Seidel variation of the fixed-point algorithm, we specialize the above linear transformation T to our current setting. $T:R^{2N^2\to R^{2N^2}}$ for a given $u^{\left(i\right)}\in R^{N^2}$ is defined by $Tv=B{u}^{\left(i\right)}+(I-{\mu }_{1}B{B}^t)v,\mathrm{\forall }v\in R^{2N^2}$.

By (22(21) $B:=\left(\begin{array}{c}{I}_N\otimes D\\ D\otimes I_N\end{array}\right),$(21) ), we have (22) $Tv=\left(\begin{array}{c}(I\otimes D)u^{\left(i\right)}+v_u-{\mu }_1\left(\right(I\otimes D{D}^t)v_u+(D^t\otimes D\left)v_l\right)\\ (D\otimes I)u^{\left(i\right)}+v_l-{\mu }_1\left(\right(D\otimes D^t)v_u+(D{D}^t\otimes I\left)v_l\right)\end{array}\right).$(22) With this notations, we have the following iteration scheme: (23) $\begin{array}{rl}{V}_U^{n+1}& =k{V}_U^n+(1-k)(I-pro{x}_{{\lambda }_1/(\alpha \ast {\mu }_1)}\parallel \cdot {\parallel }_2)\left(p\right),\\ V_L^{n+1}& =k{V}_L^n+(1-k)(I-pro{x}_{{\lambda }_1/(\alpha \ast {\mu }_1)}\parallel \cdot {\parallel }_2)\left(q\right),\end{array}$(23) where (24) $\begin{array}{rl}p& =D{U}^{\left(i\right)}+V_U^n-{\mu }_1(D{D}^t{V}_U^n+D{V}_L^{n}D),\\ q& =U^{\left(i\right)}{D}^t+V_L^n-{\mu }_1(D^t{V}_U^n{D}^t+V_L^{n}D{D}^t).\end{array}$(24) are two components of (23(22) $Tv=\left(\begin{array}{c}(I\otimes D)u^{\left(i\right)}+v_u-{\mu }_1\left(\right(I\otimes D{D}^t)v_u+(D^t\otimes D\left)v_l\right)\\ (D\otimes I)u^{\left(i\right)}+v_l-{\mu }_1\left(\right(D\otimes D^t)v_u+(D{D}^t\otimes I\left)v_l\right)\end{array}\right).$(22) ).

In the end, we obtain the minimization of (19(18) $\underset{v}{min}\frac{\alpha }{2}\parallel v-u^{\left(i\right)}{\parallel }_2^2+{\lambda }_1\parallel v{\parallel }_{TV}.$(18) ): (25) $pro{x}_{\varphi \circ B}{U}^{\left(i\right)}=U^{\left(i\right)}-{\mu }_1(D^t{V}_U+V_{L}D).$(25) Similar Gaussian–Seidel variation of the fixed-point algorithm for can be found in [26].

2.3. Summary of the alternating algorithm Fl0pro

The alternating algorithm based on $l_0$ norm in framelet transform for deblurring and proximity algorithm for denoising is named as Fl0pro.

The algorithm is summarized as follows.

3. Convergence analysis

Convergence of iteration for problem (4(4) $\underset{u,v}{min}\mathit{J}\mathit{(}\mathit{u}\mathit{,}\mathit{v}\mathit{)}\equiv \underset{\mathit{u}\mathit{,}\mathit{v}}{min}\frac{\mathit{1}}{\mathit{2}}\parallel \mathit{A}\mathit{u}\mathit{-}\mathit{f}{\parallel }_{\mathit{2}}^{\mathit{2}}\mathit{+}\parallel \lambda \cdot \mathit{W}\mathit{u}{\parallel }_{\mathit{0}}\mathit{+}\frac{\alpha }{\mathit{2}}\parallel \mathit{u}\mathit{-}\mathit{v}{\parallel }_{\mathit{2}}^{\mathit{2}}\mathit{+}{\lambda }_{\mathit{1}}\parallel \mathit{v}{\parallel }_{\mathrm{TV}}\mathit{,}$(4) ) can be decomposed into two parts, denoising and deblurring, we focus on the existence of optimal solutions to problem (5(5) $u^{\left(i\right)}=\arg \underset{u}{min}\frac{1}{2}\parallel Au-f{\parallel }_2^2+\frac{\alpha }{2}\parallel u-v^{(i-1)}{\parallel }_2^2+\parallel \lambda \cdot Wu{\parallel }_0,$(5) ), the convergence of the fixed-point algorithm for (6(6) $v^{\left(i\right)}=\arg \underset{v}{min}\frac{\alpha }{2}\parallel v-u^{\left(i\right)}{\parallel }_2^2+{\lambda }_1\parallel v{\parallel }_{\mathrm{TV}},$(6) ) is obvious.

The penalty problem (5(5) $u^{\left(i\right)}=\arg \underset{u}{min}\frac{1}{2}\parallel Au-f{\parallel }_2^2+\frac{\alpha }{2}\parallel u-v^{(i-1)}{\parallel }_2^2+\parallel \lambda \cdot Wu{\parallel }_0,$(5) ) can be in general formulated as: (26) $\underset{u}{min}\{\parallel \lambda \cdot Wu{\parallel }_0+{\phi }_1(\parallel Au-f{\parallel }_2)+\alpha {\phi }_2(\parallel u-v{\parallel }_2\left)\right\},$(26) where the functions ${\phi }_1,{\phi }_2:R_{+}\to R_{+}$ are nondecreasing continuous functions.

For any positive integer k, denote ${\mathbb{N}}_k:=1,2,\dots ,k\mathrm{and}{\mathbb{N}}_k^0:=0,1,2,\dots ,k.$ Given ${\phi }_1,{\phi }_2$ satisfying the above conditions, the integer L, the sets $\{s_i^1,s_i^2,s_i\in \mathbb{N}:i\in {\mathbb{N}}_L^0\}$, $\{{\rho }_i^1,{\rho }_i^2\ge 0:i\in {\mathbb{N}}_L^0\}$, $\{{\mathrm{\Omega }}_i^1,{\mathrm{\Omega }}_i^2,{\mathrm{\Omega }}_i\subseteq {\mathbb{R}}^n:i\in {\mathbb{N}}_L^0\}$ are defined by the following iteration $\begin{array}{rl}{\rho }_{i+1}^1& =min\left\{{\phi }_1\right(\parallel Au-f{\parallel }_2),\parallel Wu{\parallel }_0\le s_i-1,u\in {\mathbb{R}}^n\},\\ {\rho }_{i+1}^2& =min\left\{{\phi }_2\right(\parallel u-v{\parallel }_2),\parallel Wu{\parallel }_0\le s_i-1,u\in {\mathbb{R}}^n\},\\ s_{i+1}^1& =min\{\parallel Wu{\parallel }_0:{\phi }_1(\parallel Au-f{\parallel }_2)={\rho }_{i+1}^1,u\in {\mathbb{R}}^n\},\\ s_{i+1}^2& =min\{\parallel Wu{\parallel }_0:{\phi }_2(\parallel u-v{\parallel }_2)={\rho }_{i+1}^2,u\in {\mathbb{R}}^n\},\\ & s_{i+1}=min\{s_{i+1}^1,s_{i+1}^2\},\\ {\mathrm{\Omega }}_{i+1}^1& =\arg min\{\parallel Wu{\parallel }_0:{\phi }_1(\parallel Au-f{\parallel }_2)={\rho }_{i+1}^1,u\in {\mathbb{R}}^n\},\\ {\mathrm{\Omega }}_{i+1}^2& =\arg min\{\parallel Wu{\parallel }_0:{\phi }_2(\parallel u-v{\parallel }_2)={\rho }_{i+1}^2,u\in {\mathbb{R}}^n\},\\ & {\mathrm{\Omega }}_{i+1}={\mathrm{\Omega }}_{i+1}^1\cap {\mathrm{\Omega }}_{i+1}^2\end{array}$ with initial values $\begin{array}{rl}{\rho }_0^1& =min\left\{{\phi }_1\right(\parallel Au-f{\parallel }_2):u\in {\mathbb{R}}^n\},\\ {\rho }_0^2& =min\left\{{\phi }_2\right(\parallel u-v{\parallel }_2):u\in {\mathbb{R}}^n\},\\ s_0^1& =min\{\parallel Wu{\parallel }_0:{\phi }_1(\parallel Au-f{\parallel }_2)={\rho }_0^1,u\in {\mathbb{R}}^n\},\\ s_0^{2}7& =min\{\parallel Wu{\parallel }_0:{\phi }_2(\parallel u-v{\parallel }_2)={\rho }_0^2,u\in {\mathbb{R}}^n\},\\ & s_0=min\{s_0^1,s_0^2\},\\ {\mathrm{\Omega }}_0^1& =\arg min\{\parallel Wu{\parallel }_0:{\phi }_1(\parallel Au-f{\parallel }_2)={\rho }_0^1,u\in {\mathbb{R}}^n\},\\ {\mathrm{\Omega }}_0^2& =\arg min\{\parallel Wu{\parallel }_0:{\phi }_2(\parallel u-v{\parallel }_2)={\rho }_0^2,u\in {\mathbb{R}}^n\}\\ & {\mathrm{\Omega }}_0^{=}{\mathrm{\Omega }}_0^1\cap {\mathrm{\Omega }}_0^2.\end{array}$ It is well defined and there the following properties are straight forward.

1. $s_0>s_1>\cdots >s_L=0.$

2. $0\le {\rho }_0^1<{\rho }_1^1\cdots <{\rho }_L^1={\phi }_1(\parallel f{\parallel }_2).$

3. $0\le {\rho }_0^2<{\rho }_1^2\cdots <{\rho }_L^2={\phi }_2(\parallel v{\parallel }_2).$

4. ${\mathrm{\Omega }}_L=\left\{00_n\right\}$, ${\mathrm{\Omega }}_i\ne \mathrm{\varnothing }$ and ${\mathrm{\Omega }}_i\bigcap {\mathrm{\Omega }}_j=\mathrm{\varnothing }$, for $i\ne j,i,j\in {\mathbb{N}}_L^0.$

5. ${\mathrm{\Omega }}_i=\{u:\parallel Wu{\parallel }_0=s_i,{\phi }_1(\parallel Au-f{\parallel }_2)={\rho }_i^1,{\phi }_2(\parallel u-v{\parallel }_2)={\rho }_i^2,u\in {\mathbb{R}}^n\}$ for $i\in {\mathbb{N}}_L^0.$

where $00_n$ denotes the n-dimensional vector with 0 as its components.

In this way, we can see that in order to establish the existence of optimal solutions to problem (27(26) $\underset{u}{min}\{\parallel \lambda \cdot Wu{\parallel }_0+{\phi }_1(\parallel Au-f{\parallel }_2)+\alpha {\phi }_2(\parallel u-v{\parallel }_2\left)\right\},$(26) ), it suffices to prove the existence of optimal solutions to problem: $\begin{array}{rl}& \underset{u}{min}\{\parallel \lambda \cdot Wu{\parallel }_0+{\phi }_1(\parallel Au-f{\parallel }_2)+\alpha {\phi }_2(\parallel u-v{\parallel }_2):\parallel Wu{\parallel }_0\ge s_0,u\in {\mathbb{R}}^n\},\\ & \underset{u}{min}\{\parallel \lambda \cdot Wu{\parallel }_0+{\phi }_1(\parallel Au-f{\parallel }_2)+\alpha {\phi }_2(\parallel u-v{\parallel }_2):s_{i+1}\le \parallel Wu{\parallel }_0\le s_i-1,u\in {\mathbb{R}}^n\}\end{array}$ for all $i\in {\mathbb{N}}_{L-1}^0$.

Similarly as [32], we can prove the following results.

Theorem 1

Given ${\phi }_1,{\phi }_2$ satisfying the above conditions, let $L,s_i^1,s_i^2,s_i$ and ${\rho }_i^1,{\rho }_i^2,{\mathrm{\Omega }}_i^1,{\mathrm{\Omega }}_i^2,{\mathrm{\Omega }}_i$ for $i\in {\mathbb{N}}_L^0$ be defined by the above. Then, the optimal solution set of problem (27(26) $\underset{u}{min}\{\parallel \lambda \cdot Wu{\parallel }_0+{\phi }_1(\parallel Au-f{\parallel }_2)+\alpha {\phi }_2(\parallel u-v{\parallel }_2\left)\right\},$(26) ) is nonempty for any $\lambda ,\alpha >0$. Moreover, we have the following statements:

1. The optimal value of problem (27(26) $\underset{u}{min}\{\parallel \lambda \cdot Wu{\parallel }_0+{\phi }_1(\parallel Au-f{\parallel }_2)+\alpha {\phi }_2(\parallel u-v{\parallel }_2\left)\right\},$(26) ) is $min\left\{b_i\right(\lambda ,\alpha ):i\in {\mathbb{N}}_L^0\}$,

2. ${\cup }_{i\in {\wedge }^{\ast }}{\mathrm{\Omega }}_i$ is the optimal solution set of problem (27(26) $\underset{u}{min}\{\parallel \lambda \cdot Wu{\parallel }_0+{\phi }_1(\parallel Au-f{\parallel }_2)+\alpha {\phi }_2(\parallel u-v{\parallel }_2\left)\right\},$(26) ), where ${\wedge }^{\ast }:=\arg min\left\{b_i\right(\lambda ,\alpha ):i\in {\mathbb{N}}_L^0\}$ and $b_i(\lambda ,\alpha )=\lambda s_i+{\rho }_i^1+\alpha {\rho }_i^2$.

Details are omitted here.

We present the convergence error of DAL method and MDAL method for the image ‘Barbara’ and ‘Kitten’ in Figure 1, it shows that our MDAL has better performance than DAL.

Figure 1. Relatively errors for the DAL method (dash line) and the MDAL method (solid line).

4. Numerical experiments

In this section, we test four different images which are ‘Barbara’ of size $200\times 200$, ‘Kitten’ and ‘Bridge’ of size $256\times 256$, ‘Mandrill’ of size $512\times 512$ shown in Figure 2 to show the effectiveness of our proposed method and algorithm for image restoration problems. In addition, we show some blur kernels of our tests added in the figures in Figure 3. In the following examples, Signal-to-noise ratio (SNR), peak signal-to-noise ratio (PSNR), and structural similarity (SSIM) are used to measure the quality of the restored images, and blurred signal to noise ratio (BSNR) for indicating the noise size of the observed images. They are defined as follows: $\begin{array}{rl}SNR(\tilde{u},u)& =10{\log }_{10}\frac{\parallel u{\parallel }^2}{\parallel \tilde{u}-u{\parallel }^2},\\ PSNR(\tilde{u},u)& =10{\log }_{10}\frac{{255}^{2}p}{\parallel \tilde{u}-u{\parallel }^2},\\ BSNR(f,u)& =10{\log }_{10}\frac{\parallel f{\parallel }^2}{\parallel ϵ{\parallel }^2},\end{array}$ where $u,f,ϵ,\tilde{u}$ denote the original image, the observed image, the noise vector added in the test, and the restored image, p is the size of the image, respectively. Mainly, we have the larger SNR, PSNR, SSIM values mean better restored image, and the smaller BSNR value means bigger noise added to the image.

Figure 2. Original images for ‘Barbara’,‘Kitten’,‘Bridge’ and ‘Mandrill’, respectively.

Figure 3. Blur kernels: GB(9,2), AB(9,9), MB(20,0), respectively.

We compare the proposed method which based on $l_0$ norm in framelet transform for deblurring and proximity algorithm for denoising (Fl0pro) with the fast TV method (Fast-TV) in (3(3) $\begin{array}{rl}& \underset{u,v}{min}\frac{1}{2}\parallel Au-f{\parallel }_2^2+{\alpha }_1\parallel u-v{\parallel }_2^2+{\alpha }_2\parallel v{\parallel }_{\mathrm{TV}}\\ & =\underset{v}{min}\underset{u}{min}\{\parallel Au-f{\parallel }_2^2+{\alpha }_1\parallel u-v{\parallel }_2^2\}+{\alpha }_2\parallel v{\parallel }_{\mathrm{TV}},\end{array}$(3) ), and the MDAL method for $l_0$ minimization (Fl0) [28]. In all experiments, the stop criterion $\parallel u^{k+1}-u^k{\parallel }_2/\parallel u^k{\parallel }_2<5e-04$ is used to terminate iteration. Moreover, we set the maximum iterative number as 100. All the experiments are performed under Windows 7 and MATLAB R2014a running on a desktop (Intel(R) Core(TM) i7-4790 CPU @ 3.60  GHz).

4.1. Parameters selection

In this part, we will mainly discuss some basic guidelines on how to choose the parameters appeared in our algorithms. We determine the best values of parameters from their tested values such that the SNR, PSNR, SSIM values are the relative maximum.

For our proposed model (4(4) $\underset{u,v}{min}\mathit{J}\mathit{(}\mathit{u}\mathit{,}\mathit{v}\mathit{)}\equiv \underset{\mathit{u}\mathit{,}\mathit{v}}{min}\frac{\mathit{1}}{\mathit{2}}\parallel \mathit{A}\mathit{u}\mathit{-}\mathit{f}{\parallel }_{\mathit{2}}^{\mathit{2}}\mathit{+}\parallel \lambda \cdot \mathit{W}\mathit{u}{\parallel }_{\mathit{0}}\mathit{+}\frac{\alpha }{\mathit{2}}\parallel \mathit{u}\mathit{-}\mathit{v}{\parallel }_{\mathit{2}}^{\mathit{2}}\mathit{+}{\lambda }_{\mathit{1}}\parallel \mathit{v}{\parallel }_{\mathrm{TV}}\mathit{,}$(4) ), we need to consider the parameters $(\alpha ,\lambda ,\mu ,\gamma ,{\mu }_1,{\lambda }_1)$. Here α measures the trade-off between the deblurred image u and the denoised image v, we find for all tested cases that α can be set to be candidate values, i.e. $\{1e-4,2e-4,3e-4,1e-3,5e-3\}$. $\lambda ,\mu ,\gamma$ are related to the MDAL method. We choose smaller or larger λ when BSNR is equal to 40 or 30 dB. If the parameter λ is too big, it may result in the restored image over-smooth. If λ is too small, it may induce incompleteness of denoising. μ is only related to the convergence rate of the augmented Lagrangian, setting different values of μ does not change the outcome. The parameter γ controls the regularity of the sequence $\left(\right(u^{\left(k\right)},d^{\left(k\right)})$, here it has two candidate values $\{1e-5,2e-5\}$. According to the computational results, the choose of ${\mu }_1$ should satisfy the condition of ${\mu }_1<1/(4{\sin }^2/(N-1)\pi \left(2N\right))$, here we choose ${\mu }_1=\frac{1}{4},\frac{1}{8},\frac{1}{16}$. For ${\lambda }_1$, we select ${\lambda }_1=1e-2,1e-1$. From our trials, we know that too large ${\lambda }_1$ leads to too smooth of denoised image, and too small ${\lambda }_1$ may cause noise residue.

4.2. Numerical results

In this part, we present numerical results for three kinds of cases, the noise is the same Gaussian,and the blur is three kinds of different types.

Case 1. Gaussian blur and Gaussian noise

In this section, we demonstrate the efficiency of our proposed method for different images with Gaussian blur and Gaussian noise. We do the comparison of the three algorithms by the SNR, PSNR, SSIM values, number of iterations and CPU running time when the maximum iterative number is the same. Here the Gaussian blur is $3\times 3,5\times 5$, $7\times 7$, $9\times 9$ pixels with standard deviation ${\sigma }^2=2$ and the additive Gaussian noise is 40 or 30 dB.

In Table 1, the ‘Blur/ BSNR’ means that the degraded images have the type and level of blur and the size of noise, the SNR values of the degraded image are showed in the column of ‘Blur/ BSNR’. According to the values we see that the SNR, PSNR, SSIM values of the restored images by the proposed method (Fl0pro) is almost bigger than those by the fast TV method (Fast-TV) in (3(3) $\begin{array}{rl}& \underset{u,v}{min}\frac{1}{2}\parallel Au-f{\parallel }_2^2+{\alpha }_1\parallel u-v{\parallel }_2^2+{\alpha }_2\parallel v{\parallel }_{\mathrm{TV}}\\ & =\underset{v}{min}\underset{u}{min}\{\parallel Au-f{\parallel }_2^2+{\alpha }_1\parallel u-v{\parallel }_2^2\}+{\alpha }_2\parallel v{\parallel }_{\mathrm{TV}},\end{array}$(3) ), and the MDAL method based on $l_0$ minimization (Fl0). The proposed Fl0pro method is generally much better than Fl0 in computation efficiency, but consumes more running time in seconds and number of iterations than Fast-TV in general.

Table 1. SNR, PSNR(dB), SSIM values, iteration numbers, and CPU time of different methods for the different images corrupted by Gaussian blur and Gaussian noise.

Image	Blur/ BSNR	Method	SNR/PSNR	SSIM	Ites	Time
Barbara	$\mathrm{GB}(3,2)/4010.28$	Fl0pro	19.93/33.04	0.9497	12	0.78
		Fl0	18.23/31.33	0.9423	15	1.14
		Fast-TV	18.68/31.79	0.9434	12	0.09
	$\mathrm{GB}(5,2)/408.42$	Fl0pro	15.91/29.02	0.8889	11	0.73
		Fl0	15.26/28.36	0.8755	15	1.08
		Fast-TV	15.11/28.22	0.8852	14	0.09
	$\mathrm{GB}(7,2)/408.29$	Fl0pro	11.82/24.92	0.7846	10	0.61
		FL0	11.52/24.62	0.7647	21	1.35
		Fast-TV	10.90/24.01	0.7728	17	0.12
	$\mathrm{GB}(9,2)/308.10$	Fl0pro	9.34/22.44	0.6844	13	0.80
		FL0	9.10/22.20	0.6822	20	1.28
		Fast-TV	9.23/22.34	0.6512	19	0.13
Kitten	$\mathrm{GB}(3,2)/4013.42$	Fl0pro	21.62/33.55	0.9523	9	0.87
		FL0	20.38/32.30	0.9503	8	0.74
		Fast-TV	21.09/33.02	0.9444	11	0.12
	$\mathrm{GB}(5,2)/4010.17$	Fl0pro	17.10/29.03	0.9046	8	0.77
		FL0	16.46/28.39	0.8977	11	1.01
		Fast-TV	16.91/28.84	0.9023	13	0.14
	$\mathrm{GB}(7,2)/409.12$	Fl0pro	14.14/26.07	0.8348	13	1.22
		FL0	13.70/25.63	0.8237	14	1.23
		Fast-TV	13.81/25.74	0.8358	17	0.18
	$\mathrm{GB}(9,2)/308.75$	Fl0pro	11.78/23.70	0.7378	12	1.11
		FL0	10.77/22.70	0.7395	24	2.13
		Fast-TV	11.66/23.59	0.7135	20	0.21
Bridge	$\mathrm{GB}(3,2)/4010.53$	Fl0pro	18.71/32.54	0.9532	9	0.84
		FL0	18.38/32.21	0.9491	8	0.72
		Fast-TV	18.54/32.37	0.9559	11	0.12
	$\mathrm{GB}(5,2)/408.45$	Fl0pro	14.42/28.25	0.8803	11	1.04
		FL0	14.25/28.08	0.8650	11	0.97
		Fast-TV	13.92/27.75	0.8658	18	0.19
	$\mathrm{GB}(7,2)/407.67$	Fl0pro	11.84/25.67	0.7827	17	1.59
		FL0	11.65/25.49	0.7738	15	1.36
		Fast-TV	11.25/25.08	0.7476	26	0.27
	$\mathrm{GB}(9,2)/307.37$	Fl0pro	9.73/23.56	0.6710	13	1.16
		FL0	9.25/23.08	0.6671	24	2.08
		Fast-TV	9.71/23.54	0.6567	28	0.29
Mandrill	$\mathrm{GB}(3,2)/407.90$	Fl0pro	16.88/32.48	0.9454	8	3.15
		FL0	16.59/32.19	0.9434	6	2.25
		Fast-TV	16.44/32.04	0.9483	12	0.58
	$\mathrm{GB}(5,2)/406.15$	Fl0pro	12.49/28.09	0.8731	9	3.53
		FL0	12.25/27.85	0.8593	9	3.37
		Fast-TV	11.64/27.24	0.8632	13	0.62
	$\mathrm{GB}(7,2)/405.69$	Fl0pro	9.53/25.13	0.7630	10	3.87
		FL0	9.48/25.08	0.7622	11	4.13
		Fast-TV	8.61/24.22	0.7299	15	0.72
	$\mathrm{GB}(9,2)/305.54$	Fl0pro	7.15/22.75	0.6275	6	2.32
		FL0	6.83/22.43	0.6261	16	5.92
		Fast-TV	7.04/22.64	0.6260	24	1.14

Here we show the image of ‘Barbara’ for GB(5,2)/40, the ‘Kitten’ and ‘Bridge’ for GB(7,2)/40, and the ‘Mandrill’ for GB(9,2)/30, the others have the similar effect. From Figure 4, the visual quality of restored images by using our proposed method (Fl0pro) are quite competitive with those restored by using the other two existent methods. We can find the images restored by the fast TV method appearing stair-casing, especially obvious for higher blur and noise level. Also, the MDAL method based on $l_0$ minimization cannot wipe off the noise completely, the images restored by the Fl0 still exist noise.

Figure 4. (a1), (b1), (c1), (d1) blurred images degraded by GB(5)/40, GB(7)/40, GB(7)/40, GB(9)/30, respectively; (a2), (b2), (c2), (d2) corresponding restored images by Fast-TV; (a3), (b3), (c3), (d3) corresponding restored images by Fl0; (a4), (b4), (c4), (d4) corresponding restored images by Fl0pro.

Case 2. Average blur and Gaussian noise Here we act on different images with Average blur and Gaussian noise to demonstrate the efficiency of our proposed methods. We also do the comparison of the three algorithms by the SNR, PSNR, SSIM values, number of iterations and CPU running time. Here the Average blur is $5\times 5$, $7\times 7$, $9\times 9$ pixels and the additive Gaussian noise is 40 or 30 dB.

From Table 2, we see that the SNR, PSNR, SSIM values of the restored images by our proposed method (Fl0pro) are bigger than those by the fast TV method (Fast-TV) in (3(3) $\begin{array}{rl}& \underset{u,v}{min}\frac{1}{2}\parallel Au-f{\parallel }_2^2+{\alpha }_1\parallel u-v{\parallel }_2^2+{\alpha }_2\parallel v{\parallel }_{\mathrm{TV}}\\ & =\underset{v}{min}\underset{u}{min}\{\parallel Au-f{\parallel }_2^2+{\alpha }_1\parallel u-v{\parallel }_2^2\}+{\alpha }_2\parallel v{\parallel }_{\mathrm{TV}},\end{array}$(3) ) and the MDAL method based on $l_0$ minimization (Flo). However, the proposed Fl0pro method needs more running time in seconds and number of iterations than Fast-TV in general.

Table 2. SNR, PSNR(dB), SSIM values, iteration numbers, and CPU time(second) of different methods for the different images corrupted by average blur and Gaussian noise.

Image	Blur/ BSNR	Method	SNR/PSNR	SSIM	Ites	Time
Barbara	$\mathrm{AB}(5,5)/407.84$	Fl0pro	16.12/29.22	0.8895	11	0.67
		Fl0	15.15/28.26	0.8756	13	0.83
		Fast-TV	15.34/28.44	0.8811	14	0.10
	$\mathrm{AB}(7,7)/407.67$	Fl0pro	13.65/26.75	0.8385	9	0.55
		Fl0	13.01/26.11	0.8105	18	1.13
		Fast-TV	12.81/25.92	0.8301	15	0.11
	$\mathrm{AB}(9,9)/307.11$	Fl0pro	10.43/23.54	0.7158	12	0.97
		FL0	10.25/23.36	0.6937	21	1.38
		Fast-TV	10.18/23.29	0.7020	20	0.14
Kitten	$\mathrm{AB}(5,5)/409.31$	Fl0pro	17.37/29.30	0.9024	10	0.94
		FL0	16.16/28.09	0.8842	10	0.94
		Fast-TV	16.98/28.91	0.8975	12	0.13
	$\mathrm{AB}(7,7)/407.73$	Fl0pro	15.12/27.05	0.8559	13	1.21
		FL0	14.62/26.54	0.8456	11	1.03
		Fast-TV	14.88/26.81	0.8597	14	0.15
	$\mathrm{AB}(9,9)/306.88$	Fl0pro	11.58/23.50	0.7186	11	1.06
		FL0	11.34/23.26	0.7015	17	1.49
		Fast-TV	11.33/23.26	0.7123	20	0.22
Bridge	$\mathrm{AB}(5,5)/407.74$	Fl0pro	15.17/29.01	0.8945	11	1.01
		FL0	14.83/28.66	0.8783	11	0.98
		Fast-TV	15.10/28.93	0.8944	14	0.15
	$\mathrm{AB}(7,7)/406.44$	Fl0pro	13.48/27.31	0.8447	13	1.25
		FL0	13.25/27.08	0.7998	15	1.31
		Fast-TV	13.14/26.97	0.8351	15	0.16
	$\mathrm{AB}(9,9)/305.60$	Fl0pro	10.39/24.22	0.6990	16	1.50
		FL0	10.28/24.11	0.6522	20	1.77
		Fast-TV	10.20/24.03	0.6851	24	0.25
Mandrill	$\mathrm{AB}(5,5)/405.57$	Fl0pro	13.36/28.96	0.8856	9	3.39
		FL0	13.09/28.69	0.8821	7	2.57
		Fast-TV	12.50/28.10	0.8842	13	0.63
	$\mathrm{AB}(7,7)/404.88$	Fl0pro	11.44/27.04	0.8345	9	3.33
		FL0	11.19/26.79	0.8164	8	2.96
		Fast-TV	10.35/25.95	0.8147	14	0.68
	$\mathrm{AB}(9,9)/304.51$	Fl0pro	8.08/23.68	0.6631	11	4.05
		FL0	7.96/23.56	0.6370	13	4.73
		Fast-TV	8.02/23.63	0.6590	16	0.78

Here we show the experimental results about the image of ‘Barbara’ for AB(5)/40, the ‘Kitten’ and ‘Bridge’ for AB(7)/40, and the ‘Mandrill’ for AB(9)/30, the other different images with different blur and different noise have the similar effect which do not show exhaustively. From Figure 5, We can find the images restored by the fast TV method appearing stair-casing. Also, the MDAL method based on $l_0$ minimization cannot wipe off the noise completely. From the experiments we conducted, we obtain that the validity and effectiveness of our method and algorithm.

Figure 5. (a1), (b1), (c1), (d1) blurred images degraded by AB(5)/40, AB(7)/40, AB(7)/40, AB(9)/30, respectively; (a2), (b2), (c2), (d2) corresponding restored images by Fast-TV; (a3), (b3), (c3), (d3) corresponding restored images by Fl0; (a4), (b4), (c4), (d4) corresponding restored images by Fl0pro.

Case 3. Motion blur and Gaussian noise Here we do experiments to show the efficiency of our proposed methods for different images with Motion blur and Gaussian noise. We do the comparison of the three algorithms by the SNR, PSNR, SSIM values, number of iterations and CPU running time. Here the Motion blur is the length of len=9,20 with an angle $\theta =0$ in a counterclockwise direction and the additive Gaussian noise is 40 or 30 dB.

In Table 3, we see that the SNR, PSNR, SSIM values of the restored images by the proposed method Fl0pro are better than those by the fast TV method (Fast-TV) in (3(3) $\begin{array}{rl}& \underset{u,v}{min}\frac{1}{2}\parallel Au-f{\parallel }_2^2+{\alpha }_1\parallel u-v{\parallel }_2^2+{\alpha }_2\parallel v{\parallel }_{\mathrm{TV}}\\ & =\underset{v}{min}\underset{u}{min}\{\parallel Au-f{\parallel }_2^2+{\alpha }_1\parallel u-v{\parallel }_2^2\}+{\alpha }_2\parallel v{\parallel }_{\mathrm{TV}},\end{array}$(3) ), and MDAL method based on $l_0$ minimization (Flo). The computational efficiency of our proposed method is lower than Fast-TV method.

Table 3. SNR, PSNR(dB), SSIM values, iteration numbers, and CPU time(second) of different methods for the different images corrupted by motion blur and Gaussian noise.

Image	Blur/ BSNR	Method	SNR/PSNR	SSIM	Ites	Time
Barbara	$\mathrm{MB}(9,0)/407.54$	Fl0pro	17.92/31.02	0.9050	12	0.78
		Fl0	17.18/30.29	0.8920	16	1.08
		Fast-TV	16.71/29.82	0.9043	21	0.15
	$\mathrm{MB}(20,0)/305.77$	Fl0pro	9.59/22.69	0.6427	16	0.98
		Fl0	9.19/22.30	0.6361	23	1.53
		Fast-TV	9.09/22.20	0.6383	38	0.27
Kitten	$\mathrm{MB}(9,0)/409.13$	Fl0pro	20.13/32.06	0.9372	13	1.21
		FL0	19.24/31.17	0.9215	14	1.24
		Fast-TV	19.83/31.76	0.9298	12	0.13
	$\mathrm{MB}(20,0)/306.52$	Fl0pro	13.03/24.96	0.7742	13	1.21
		FL0	12.58/24.51	0.7781	25	2.18
		Fast-TV	12.42/24.34	0.7661	28	0.29
Bridge	$\mathrm{MB}(9,0)/408.19$	Fl0pro	18.25/32.08	0.9415	13	1.25
		FL0	17.94/31.77	0.9322	15	1.33
		Fast-TV	17.76/31.59	0.9387	21	0.23
	$\mathrm{MB}(20,0)/305.93$	Fl0pro	11.08/24.91	0.7422	21	1.90
		FL0	10.92/24.75	0.7357	26	2.27
		Fast-TV	11.03/24.86	0.7408	23	0.24
Mandrill	$\mathrm{MB}(9,0)/406.49$	Fl0pro	17.35/32.95	0.9407	7	2.81
		FL0	17.14/32.74	0.9405	8	2.97
		Fast-TV	16.28/31.88	0.9406	19	0.92
	$\mathrm{MB}(20,0)/304.64$	Fl0pro	10.04/25.64	0.7655	11	4.15
		FL0	9.83/25.43	0.7629	13	4.79
		Fast-TV	9.24/24.84	0.7419	19	0.91

Here we show the image of ‘Barbara’ and ‘Kitten’ for MB(9,0)/40, and ‘Bridge’, ‘Mandrill’ for MB(20,0)/30, the others have the similar effect. From Figure 6, the visual quality of restored images by using our proposed method (Fl0pro) are quite competitive with those restored by using the other two existent methods. We can find the images restored by the fast TV method and the MDAL method based on $l_0$ minimization do have the disadvantages as the above two blur cases, which our proposed method (Fl0pro) can avoid in an excellent way. Then we can say that our proposed method and algorithm are valid and effective.

Figure 6. (a1), (b1), (c1), (d1) blurred images degraded by MB(9)/40, MB(9)/40, MB(20)/30, MB(20)/30, respectively; (a2), (b2), (c2), (d2) corresponding restored images by Fast-TV; (a3), (b3), (c3), (d3) corresponding restored images by Fl0; (a4), (b4), (c4), (d4)corresponding restored images by Fl0pro.

5. Conclusion

In this paper, we introduce fast alternating minimization methods for image deblurring and denoising. Motivated by the mentality of [27], the excellent effect of [28], and the simple and valid proximity algorithm proposed in [26], we obtain an alternating method which bases TV norm to denoise the deblurred image using analysis-based approach. We use MDAL method to solve $l_0$ minimization based on framelet transform for the deblurring part and use proximity algorithm to solve isotropic total variation for the denoising part. The experimental outcomes show that the proposed alternating minimization method is robust for additive noise and different blur kernels comparing with relative method existence.

Disclosure statement

No potential conflict of interest was reported by the authors.

Additional information

Funding

The authors are partially supported by Science Challenge project TZ2016002, NSFC [No. 11171154,11671050,11771055] and 3D numerical simulation platform TB14-1 of China Academy of Engineering Physics.