1. Introduction

In this paper, we consider the image deblurring whose task is to nd the true image f from the observed image g which is blurred and noised as

g = Pf + n, (1)

where P is the linear transform that determines the blurring process and n represents a mean zero Gaussian noise. We denote pixel sets of f and g by and, respectively.

To denote image pixels, we will use the single index without the paren-thesis for some cases and the double index with the parenthesis for other cases. For instance, fv = f(v1,v2) is the intensity of the image f at the pixel v = (v1; v2). With the single pixel index, we will treat images g and f as one dimensional vectors. In such case, P in (1) will be regarded as an Λ |×| Ω matrix. We will also use following convention throughout this paper: boldface alphabets for images or matrices and normal alphabets for their intensities or entries, as in f = (fv), g = (gb), and P = (Pb,v).

The image deconvolution has many applications in science and engi-neering elds, and many methods have been proposed for it [1]. The image deconvolution, however, often produces unsatisfactory results due to various obstacles. The paper [2] characterizes obstacles of the image deconvolution into 4 categories; noise, insucient deconvolution, boundary artifacts, and incorrect blurring model. In this paper, we will discuss the suppression of boundary artifacts as the main concern.

The blurring process makes some near-boundary pixels of the observed image g to be influenced by `unseen pixels' (pixels in Ω Λ). It is often regarded that the existence of unseen pixels causes boundary artifacts. To suppress boundary artifacts by removing unseen pixels, various boundary condition methods have been proposed. Among them, periodic, reactive, and anti-reective BCs haven been most popular. The periodic boundary condition extends the true image to satisfy the periodic condition on unseen pixels across the boundary. Similarly, reflective and anti-reflective boundary conditions extend the true image to satisfy corresponding conditions across the boundary. For details, see, e.g., [3,4]. These boundary condition methods, however, do not suppress boundary artifacts effectively in cases when imposed conditions are greatly mismatched to characteristics of images to be recovered [5]. Considering this fact, in this paper, we will not impose any restriction on unseen image pixels for the reduction of boundary artifacts. Such approaches are called free boundary condition methods in [5].

In this paper, we assume that the blurring matrix P satises

P b,v ≥0 for  all  b∈Λ  and  v∈Ω, (2)

∑ v∈Ω P b,v =1 for  all  b∈Λ, (3)

and

0< ∑ b∈Λ P b,v ≤1 for  all  v∈Λ (4)

Described three conditions (2), (3), and (4) are very mild assumptions on image deblurring problems. For example, the condition (2) must hold for any kind of photon diusion related deblurring processes. The second condition simply indicates that measurement sensitivities on are uniform. The third condition implies that in image deblurring we will consider image pixels only that give some contribution to the observation. The third condition denes normalization coecients.

w v = ∑ b∈Λ P b,v , (5)

which measures the contribution of the image pixel v ∈ Ω to the observation on the image pixel set . In this paper we will call the

Q = W^-1P^t (6)

normalized backprojector, where W is the diagonal matrix dened by (Wx)v = w_vx_v for each image x on defined on Ω and P^t is the transpose matrix. Here the term `normalized' comes from the fact that Q1_Λ = 1_Ω ,where 1_A is the all-one image dened on A.

Any iterative deblurring approaches to (1) have two types of matrixvector multiplications; one is of the form Px for some image x dened on Ω and the other is of the form P^ty form some image y defined on Λ. We call the former projection of the image x and the latter backprojection of the image y. We also call the matrix P the projector and the transpose matrix P^t the standard backprojector. Here the term `standard' is used to distinguish the normalized backprojector Q from P^t.

In this paper, we will show that a simple replacement of the standard backprojector P^t with the normalized backprojector Q can reduce boundary artifacts eciently. We will also show why such a simple leads to the efficient reduction of boundary artifacts. For these purpose, we will consider two iterative deblurring methods; the Landweber iteration and the conjugate gradient method.

In this work, we distinguish boundary artifacts from ringing artifacts which are often noticeable near strong edges in image deblurring. Ringing artifacts are caused by the fact that pixel values in smooth image region have faster convergence trend than pixel values near edges in iterative de-blurring methods. Ringing artifacts make ripples near edges in deblurred images [1]. Those ripples by ringing artifacts sometimes look similar to rip-ples by boundary artifacts. But, they are different in at least following two aspects: (a) Ripples by ringing artifacts show lines that resemble strong edge lines nearby, but ripples by boundary artifacts show lines that resem- ble boundary lines, thus, for rectangular boundaries, ripples by boundary artifacts are always exactly vertical or horizontal straight lines. (b) Ringing artifacts are caused by non-uniformity in convergence trend and, in some degree, inevitable in all iterative deblurring methods, while boundary artifacts are caused by mistreatment on image pixel values and hence preventable if pixel values are treated correctly, for instance, by using the normalized back projector proposed in this paper. Obviously, the last statement is the main claim of this paper.

This paper is outlined as follows. In Section 2 we explain why the use of the normalized backprojector is a natural choice in image deblurring applications. In Section 3 we present simulation results of the proposed method. Finally, we have discussion and conclusion in Section 4.

2. Methods

To solve (1), we consider the Tikhonov regularization formulated by

min f ( ‖g−pf ‖ 2 +λ‖f ‖ w 2 ) (7)

Where || f ||_w is the norm dened by the weighted inner product

<x,z > w = ∑ v∈Ω x v z v w v (8)

The use of the w-weighted inner product is motivated by the intention to treat image pixels proportionally to their contribution to the observation on . Simulation studies in Section 3 will show that the use of the w-weighted norm leads to better boundary artifact removal than the standard norm in(7).

In cases when = 0, the Tikhonov regularization (7) becomes the least square problem. The least square problem can be approximated by the Landweber iteration that, starting from arbitrary f⁰, takes

fⁿ⁺¹ = fⁿ −βP^t (Pfⁿ − g) (9)

with 0 < β < 2б_max , where б_max is the largest singular value of P. As a special case of gradient descent, the Landweber iteration is derived from the dierentiation of the variational form, i.e.,

∇|| g-Pf||² = P^t (pfⁿ -g) (10)

Here we emphasize that the dierential is computed with respect to the standard inner product. To nd out what happens if the differentiation is computed with respect to the inner product < .;. >U dened by a positive definite matrix U, i.e., U=Σ_v∈ΩΣ_v'∈Ω,X_vU_v,v'Z_v' we generalize the denition of the gradient of a real-valued differentiable function ϕ at the image x with respect to the inner product <.;.>U, denoted by ∇uϕ (X) the unique image dened on satisfying

lim h→0 |φ(x+h)−φ(x)−<α,h > U | <h,h > U 1/2 =0 (11)

A simple consideration shows that

∇uϕ(X)=U^-1∇ϕ(X) (12)

This result indicates that if the differentiation of the least square functional is computed with respect to the w-weighted inner product in (8), then the corresponding Landweber iteration becomes

f ⁿ⁺¹ = f ⁿ −β Q(Pfⁿ -g) (13)

where Q is the normalized backprojector in (6). Simulation studies in Section 3 will show that this normalized Landweber iteration (13) removes boundary artifacts more eectively than the standard Landweber iteration (9).

With the same argument, the Tikhonov regularization (7) depends on the inner product used in the dierentiation. For example, if the standard inner product is used, then the Tikhonov regularization (7) leads to the normal equation of the form

P^tg = (P^tP + λW)f (14)

In Section 3, we will approximate the solution of (14) by the conjugate gradient iteration. Similarly, if the w-weighted inner product is used, then the Tikhonov regularization (7) leads to the normal equation of the form

Qg = (QP + λI) (15)

where I is the identity matrix. In Section 3, we will approximate the solution of (15) by the normalized conjugate gradient iteration. Here the normalized conjugate gradient iteration means that the iteration is performed with respect to the w-weighted inner product. Simulation studies in Section 3 will show that the normalized conjugate gradient iteration applied to (15) removes boundary artifacts more eectively than the standard conjugate gradient iteration applied to (14).

Before we close this section, we note that the normalized backprojector Q is the transpose of the projector P with respect to the w-weighted inner product in a sense that for all images x on

[y;Px] = y^tPx = (Qy)^tWx =< Qy; x >w (16)

Ω and all images y on , where [; ] is the standard inner product of images on . We also note that the governing matrix QP+I in (15) is symmetric with respect to the w-weighted inner product, in a sense that for all images x; z on.

< z; (QP + I)x >_w=< (QP + I)z; x >_w (17)

This is the reason why the normalized conjugate gradient iteration can can performed for (15).

It is also worth to note that algebraically, (14) and (15) are identical. Thus, if conjugate gradient iterations are completely executed, then final results must be identical. In practical use of conjugate gradient iterations, however, incompletely iterated conjugate gradient approximates are preferred. In such case, results are very dierent.

3. Simulation Studies

We conducted simulations to compare normalized versions of Landweber and conjugate gradient iterations with standard version to evaluate the effect of the normalized backprojector. For this purpose, we used test images (a) `couple' and `boat' in Figure 1.

Figure 1: Test images. (a)(left) `peppers' (b)(right) `boat'.

For simulation studies, we consider translation-invariant blurring models. In such cases, the blurring process Pf can be described by the truncated convolution of a point spread function (PSF) k and the image f . To be specific,

(pf) b = ∑ v∈Ω P b,v f v = ∑ v∈(b− supp k ) k b − f v v . (18)

Here supp_k, the support of k, is {u | K_u > 0}g. The conditions described in (2), (3), and (4) imply that the PSF k is nonnegative, its components have sum 1, and the point (0; 0) Є supp_k ; Pf is dened on Λ, where b Є Λ if and only if b - supp_k С Ω. We used the horizontally and vertically symmetrical Gaussian PSF k, and regard the index of the center of the Gaussian PSF as (0; 0). For the diagonal PSF, we assumed that the upper-left corner pixel has the index (0; 0) and the PSF decays in the diagonal direction from the upper-left corner pixel to the lower-right corner pixel. Figure 2 shows 15 × 15 (a)(left) Gaussian and (b)(right) diagonal PSFs.

Figure 2: PSFs. (a)(left) Gaussian (b)(right) diagonal.

Figure 3 shows observed images, which are blurred by (a) Gaussian and (b) diagonal PSFs and then noised mean zero Gaussian noises with the standard deviation σ = 0:5% of means of blurred images.

Figure 3: Observed images, which are blurred by 15 × 15 (a)(left) Gaussian and (b)(right) diagonal PSFs and then noised mean zero Gaussian noises with the stan-dard deviation σ = 0:5% of means of blurred images.

In simulations, we chose the 1000-th iterates for Landweber iterations and 100-th iterates for conjugate gradient iterations.

Figure 4 shows deblurred images from Figure 3(a) by (a) standard and (b) normalized Landweber iterations, and Figure 5 shows deblurred images from Figure 3(b) by (a) standard and (b) normalized Landweber iterations. Visual comparison in Figure 4 and 5 clearly shows that the normalized backprojector removes boundary artifacts more eectively in the Landwe- ber iteration. Figure 4 and 5 show artifacts in recovering near boundary pixels. This is well expected phenomenon since near boundary pixels give less contribution to the observation, recovering them is subject to produce more artifacts than recovering near center image pixels. The main problem of boundary artifacts is the propagation of errors as in Figure 4(a) and 5(a).

Figure 4: Deblurred images from Figure 3(a) by (a)(left) standard and (b) (right) normalized Landweber iterations.

Figure 5: Deblurred images from Figure 3(b) by (a)(left) standard and (b) (right) normalized Landweber iterations.

Figure 6 shows deblurred images from Figure 3(a) by (a) standard and (b) normalized conjugate gradient iterations, and Figure 7 shows deblurred images from Figure 3(b) by (a) standard and (b) normalized conjugate gra-dient iterations. Visual comparison in Figure 6 and 7 clearly shows that the normalized backprojector removes boundary artifacts more effectively in the conjugate gradient iteration.

Figure 6: Deblurred images from Figure 3(a) by (a)(left) standard and (b)(right) normalized conjugate gradient iterations.

Figure 7: Deblurred images from Figure 3(b) by (a)(left) standard and (b) (right) normalized Landweber iterations.

As mentioned earlier, the use of the w-weighted inner product in the Tikhonov regularization (7) is motivated by the intention to treat image pixels proportionally to their contribution to the observation on . To find out whether such intension leads to improved boundary artifact removal or not, we considered the Tikhonov regularization formulated by

min f ( ‖g-pf ‖ 2 +λ‖f ‖ 2 ) (19)

As shown in Section 2, this variational problem depends on the inner product used in the dierentiation. The use of the standard inner product leads to the normal equation of the form

P^tg = (P^tP + λI)f, (20)

where I is the identity matrix, and the use of the w-weighted inner product leads to

Qg = (QP + λW^-1)f ; (21)

As in (7), we approximated solution of (20) and (21) by standard and normalized conjugate gradient iterations, respectively.

Figure 8 shows deblurred images from Figure 3(a) by (a) standard and (b) normalized conjugate gradient iterations applied to (20) and (21), respectively, and Figure 9 shows deblurred images from Figure 3(b) by (a) standard and (b) normalized conjugate gradient iterations applied to (20) and (21), respectively. Again, visual comparison in Figure 8 and 9 clearly shows that the normalized backprojector removes boundary artifacts more effectively in the conjugate gradient iteration.

Figure 8: Deblurred images from Figure 3(a) by (a)(left) standard and (b)(right) normalized conjugate gradient iterations.

Figure 9: Deblurred images from Figure 3(b) by (a)(left) standard and (b)(right) normalized conjugate gradient iterations.

The comparison between Figure 8(b) and 6(b), however, shows that the boundary artifacts removal in Gaussian deblurring by the normalized conjugate gradient iteration applied to (21) is not as good as that applied to (15). For detailed comparison, see Figure 8(d) and 6(d) and RSE results of Figure 8(b) and 6(b). Similarly, the comparison between Figure 9(b) and 7(b) shows that the boundary artifacts removal in diagonal deblurring by the normalized conjugate gradient iteration applied to (21) is not as good as that applied to (15), either. For detailed comparison, see Figure 9(d) and 7(d) and RSE results of Figure 9(d) and 7(d). Based on these two comparisons, we can conclude that the use of the w-weighted inner product improves boundary artifact removal, as compared with the use of the standard inner product, in the normalized conjugate gradient iterations.

On the other hand, Figure 8(a) and 6(a) produced almost identical results. Similarly, Figure 9(a) and 7(a) produced almost identical results. Based on these comparisons, we can conclude that the use of the w-weighted inner product does not improve boundary artifact removal in the standard conjugate gradient iterations. This result also shows that the use of the normalized backprojector instead of the standard backprojector really removes boundary artifacts, not the use of the w-weighted norm itself in the Tikhonov regularization.

4. Conclusion and Discussion

In this paper, we suggest the normalized backprojector Q = W^-1P^t as a replacement of the standard backprojector P^t for the suppression of boundary artifacts related to the deblurring problem g = Pf + n in (1).

Simulation results showed that the proposed normalized backprojector can remove boundary artifacts eectively in Landweber and conjugate gradient iterations.

The success of the proposed normalized backprojector can be viewed as a benet of using unmatched projector/backprojector pairs, i.e., (P;Q) instead of (P;P^t). The use of unmatched projector/backprojector pairs in medical imaging area has been gaining the popularity since it often provides faster approximation, faster computation, or some specic effect such as ring artifact removal [6,7]. In this work we used the unmatched projector/backprojector pair in image deblurring problem to remove the boundary artifacts. We expect that it is also possible to accelerate the convergence or reduce the computational burden, by using unmatched projector/backprojector pair in image deblurring problems.

Competing Interests

The authors declare that they have no competing interests.