Single Image Blind Deblurring Using Multi-Scale Latent Structure Prior

TL;DR

This paper introduces a multi-scale latent structure prior for blind image deblurring, enabling efficient coarse-to-fine restoration of sharp images from blurry inputs, with theoretical support and competitive results.

Contribution

The paper proposes a novel multi-scale latent structure prior for blind deblurring, incorporating a coarse-to-fine approach and efficient local self-example matching, kernel estimation, and deblurring.

Findings

Achieves competitive deblurring results with faster speed.

Validates the latent structure prior theoretically.

Extends to non-uniform blind deblurring.

Abstract

Blind image deblurring is a challenging problem in computer vision, which aims to restore both the blur kernel and the latent sharp image from only a blurry observation. Inspired by the prevalent self-example prior in image super-resolution, in this paper, we observe that a coarse enough image down-sampled from a blurry observation is approximately a low-resolution version of the latent sharp image. We prove this phenomenon theoretically and define the coarse enough image as a latent structure prior of the unknown sharp image. Starting from this prior, we propose to restore sharp images from the coarsest scale to the finest scale on a blurry image pyramid, and progressively update the prior image using the newly restored sharp image. These coarse-to-fine priors are referred to as \textit{Multi-Scale Latent Structures} (MSLS). Leveraging the MSLS prior, our algorithm comprises two phases: 1) we first preliminarily restore sharp images in the coarse scales; 2) we then apply a refinement process in the finest scale to obtain the final deblurred image. In each scale, to achieve lower computational complexity, we alternately perform a sharp image reconstruction with fast local self-example matching, an accelerated kernel estimation with error compensation, and a fast non-blind image deblurring, instead of computing any computationally expensive non-convex priors. We further extend the proposed algorithm to solve more challenging non-uniform blind image deblurring problem. Extensive experiments demonstrate that our algorithm achieves competitive results against the state-of-the-art methods with much faster running speed.

Tables

Table 1. TABLE I: Quantitative comparison of all methods over the entire dataset (640 blurry images).

Blind Deblurring Method	Mean error ratio	Worst error ratio	Success rate (r 5)
Cho et al. [56]	28.1	165.0	11.7%
Krishnan et al. [5]	11.7	133.2	24.8%
Levin et al. [4]	6.6	40.9	46.7%
Cho & Lee [46]	8.7	111.1	65.5%
Xu & Jia [50]	3.6	65.3	85.8%
Sun et al. [55]	2.5	30.5	93.4%
Michaeli & Irani [16]	2.6	9.3	95.9%
Lai et al. [9]	2.1	17.9	97.3%
Ours	2.1	14.0	98.4%
Pan et al. [10]	1.6	8.8	99.1%
Ours*	1.6	10.4	99.7%

Table 2. TABLE II: Running time of different algorithms on real images.

Image	[5]	[4]	[16]	[10]	Ours
Flower	170 s	142 min	276 min	14 min	74 s
Picasso	189 s	148 min	292 min	17 min	102 s
Pietro	248 s	153 min	306 min	19 min	103 s
Roma	463 s	667 min	2167 min	41 min	282 s

Equations

$$\mathbf{b}=\mathbf{k}\otimes \mathbf{x}+\mathbf{n},$$

$$\mathbf{b}[i,j]=\mathbf{k}\otimes \mathbf{x}=\iint \mathbf{k}(u,v)\mathbf{x}(i-u,j-v)dudv,$$

$${\mathbf{b}}_{\alpha }[i,j]$$

$$=\left(\mathbf{k}\otimes (\mathbf{x}\otimes \mathbf{h})\right)\downarrow =\left(\mathbf{k}\otimes {\mathbf{x}}_h\right)\downarrow$$

$$=\iint \mathbf{k}(u,v){\mathbf{x}}_h(\alpha i-u,\alpha j-v)dudv,$$

$${\mathbf{b}}_{\alpha }[i,j]$$

$$=\iint {\mathbf{k}}_{\alpha }(t,s)\cdot {\mathbf{x}}_{\alpha }(i-t,j-s)dtds.$$

$${\mathbf{k}}_{\alpha }(u,v)={\alpha }^2\mathbf{k}(\alpha u,\alpha v).$$

$${\mathbf{K}}_{\alpha }(\tilde{u},\tilde{v})$$

$$={\alpha }^2\left(\frac{1}{{\alpha }^2}\mathbf{K}(\tilde{u}/\alpha ,\tilde{v}/\alpha )\right)$$

$$=\mathbf{K}(\tilde{u}/\alpha ,\tilde{v}/\alpha )$$

$$\delta (i,j)=\{\begin{array}{ll}1,& i=0andj=0,\\ 0,& i\ne 0orj\ne 0.\end{array}$$

$${\mathbf{b}}_{\alpha }[i,j]\to {\mathbf{x}}_{\alpha }\otimes \delta ={\mathbf{x}}_{\alpha },$$

$$\underset{\mathbf{k},\mathbf{x}}{\mathrm{argmin}}\Phi (\mathbf{k}\otimes \mathbf{x}-\mathbf{b})+{\lambda }_1{\Psi }_1(\mathbf{k})+{\lambda }_2{\Psi }_2(\mathbf{x})$$

$$d(i,j)=\Vert {\mathbf{P}}_i{\mathbf{x}}_l-{\mathbf{P}}_j{\mathbf{x}}_{pr}{\Vert }_2,j\in \mathcal{N}(i{}^{\prime })$$

$$w(i)=\theta$$

$$W(i,j)$$

$$\mathbf{b}=\mathbf{k}\otimes (\hat{\mathbf{x}}+{\mathbf{x}}_c)+\mathbf{n}.$$

$$\underset{\mathbf{k},{\mathbf{x}}_c}{\mathrm{argmin}}\frac{1}{2}\Vert \mathbf{k}\otimes \nabla (\hat{\mathbf{x}}+{\mathbf{x}}_c)-\nabla \mathbf{b}{\Vert }_2^2+{\lambda }_1\Vert \mathbf{k}{\Vert }_2^2+{\lambda }_2\Vert {\mathbf{x}}_c{\Vert }_1$$

$$\underset{\mathbf{k},\mathbf{v}}{\mathrm{argmin}}\frac{1}{2}\Vert \mathbf{k}\otimes \nabla \hat{\mathbf{x}}+\mathbf{v}-\nabla \mathbf{b}{\Vert }_2^2+{\lambda }_1\Vert \mathbf{k}{\Vert }_2^2+{\lambda }_2\Vert \mathbf{v}{\Vert }_1$$

$$\mathbf{k}={\mathcal{F}}^{-1}\left(\frac{\overline{\mathcal{F}(\nabla \hat{\mathbf{x}})}\mathcal{F}(\nabla \mathbf{b})-\overline{\mathcal{F}(\nabla \hat{\mathbf{x}})}\mathcal{F}(\mathbf{v})}{\overline{\mathcal{F}(\nabla \hat{\mathbf{x}})}\mathcal{F}(\nabla \hat{\mathbf{x}})+2{\lambda }_1}\right)$$

$$v_i=sgn(z_i)\cdot \max \left\{0,|z_i|-{\lambda }_2\right\}$$

$$\underset{\mathbf{x}}{\mathrm{argmin}}\frac{1}{2}\Vert \hat{\mathbf{k}}\otimes \mathbf{x}-\mathbf{b}{\Vert }_2^2+\mu \Vert \nabla \mathbf{x}{\Vert }_1$$

$$\underset{\mathbf{k},\mathbf{v}}{\mathrm{argmin}}\frac{1}{2}\Vert \mathbf{k}\otimes {\nabla }_{\ast }\hat{\mathbf{x}}+{\mathbf{v}}_{\ast }-$$

$$+$$

$$\mathbf{b}=\sum _i{k}_i{\mathbf{H}}_i\mathbf{x}+\mathbf{n}$$

$$\mathbf{b}={\mathbf{A}}_{\mathbf{k}}\mathbf{x}+\mathbf{n}={\mathbf{B}}_{\mathbf{x}}\mathbf{k}+\mathbf{n}$$

$$\underset{\mathbf{k},{\mathbf{x}}_c}{\mathrm{argmin}}\frac{1}{2}\Vert {\mathbf{B}}_{\nabla (\hat{\mathbf{x}}+{\mathbf{x}}_c)}\mathbf{k}-$$

$$+$$

$$({\mathbf{B}}_{\nabla \hat{\mathbf{x}}}^T{\mathbf{B}}_{\nabla \hat{\mathbf{x}}}+2{\lambda }_1\mathbf{I})\mathbf{k}={\mathbf{B}}_{\nabla \hat{\mathbf{x}}}^T\nabla \mathbf{b}-{\mathbf{B}}_{\nabla \hat{\mathbf{x}}}^T\mathbf{v}.$$

$$\underset{\mathbf{x}}{\mathrm{argmin}}\frac{1}{2}\Vert {\mathbf{A}}_{\hat{\mathbf{k}}}\mathbf{x}-\mathbf{b}{\Vert }_2^2+\mu \Vert \nabla \mathbf{x}{\Vert }_1,$$

$$r=\frac{\Vert \mathbf{x}-{\mathbf{x}}_{\hat{\mathbf{k}}}{\Vert }^2}{\Vert \mathbf{x}-{\mathbf{x}}_{\mathbf{k}}{\Vert }^2}$$

Full text

I Introduction

Blur is one of the most common artifacts of digital images. It happens mainly because the camera does not locate in focus, or the camera is held unsteadily, or the objects in the scene are moving quickly over the period of exposure time, resulting in that image sensor accumulates light not only from a single point but also from neighbouring regions. Consequently, a blurry image is captured with unclear edges and details, which greatly degrade the visual quality of the image. As a fundamental degradation model, the blur process is usually assumed to be shift-invariant, and can be represented as the convolution of a sharp image and a blur kernel,

[TABLE]

where ${\mathbf{b}}$ is the observed blurry image, ${\mathbf{x}}$ is the latent sharp image, ${\mathbf{k}}$ is the blur kernel, $\otimes$ is the convolution operation, and ${\mathbf{n}}$ is the additive Gaussian noise.

The blind image deblurring problem is to recover both the latent sharp image ${\mathbf{x}}$ and the blur kernel ${\mathbf{k}}$, given only the observed blurry image ${\mathbf{b}}$ [1, 2]. It is a highly ill-posed problem because the solution is not only unstable but also non-unique. As a convolution in the spatial domain is equivalent to the point-wise product in the frequency domain, the observation ${\mathbf{b}}$ can be factorized into infinite feasible solution pairs of ${\mathbf{x}}$ and ${\mathbf{k}}$. It is difficult to achieve the unique correct pair without sufficient prior knowledge about ${\mathbf{x}}$ and ${\mathbf{k}}$.

Previous methods [3, 4, 5, 6, 7, 8, 9, 10, 11, 12] treated the blind image deblurring as a joint optimization problem and introduced sophisticated image priors, which can penalize image blurriness and promote image sharpness, such as the mixture of Gaussian functions that fits the heavy-tailed prior of natural images [3, 4], normalized sparse prior [5], framelet based prior [6], $l_{0}$-norm based priors [7, 8], color line prior [9], dark channel prior [10], low-rank prior [11] and graph based prior [12], etc. However, these image priors are usually computationally expensive, resulting in complicated optimization algorithms.

In this paper, we propose a novel multi-scale image prior to tackle the challenging blind image deblurring problem, inspired by the well-known self-example prior in image super-resolution (SR) [13, 14, 15]. In image SR, self-example prior is obtained through image down-sampling. An observed image is down-sampled to generate a coarse-scale image, and the self-example prior is modeled as mappings from patches in the coarse-scale image to the corresponding patches in the fine-scale observed image. Although the self-example prior has been successfully applied in image SR, when an observed image is blurred by an arbitrary blur kernel, such as a motion blur, there is no straightforward conclusion that the mappings can still model the image deblurring process. In [16], Michaeli and Irani first analyzed the effect of image down-sampling in the blind image deblurring problem. They found that image down-sampling increased the internal patch recurrence. With this observation, they proposed to promote the internal patch recurrence for blind image deblurring. However, iterative promotion of patch recurrence requires enormous complicated patch matchings [17], which are extremely slow in practice. Different from [16], in this paper, we focus on the extreme case of image down-sampling. We keep down-sampling an observed blurry image to a very coarse level and observe that: the coarse enough image down-sampled from a blurry observation is approximately a low-resolution version of the latent sharp image. We theoretically verify this phenomenon on 2D signals under a general down-sampling operator, and define the coarse enough image as a latent structure prior of the latent sharp image. We further introduce a multi-scale strategy to work in tandem with latent structure prior on an image pyramid. Starting from the coarsest scale, we employ the prior image to blindly restore a finer scale sharp image. The newly restored sharp image is then employed as the latent structure prior of the next finer scale. The strategy works gradually from the coarsest scale to the finest scale. Thus, we name our prior Multi-Scale Latent Structure (MSLS) prior.

Tailored to the MSLS prior, we propose an efficient algorithm to solve the blind image deblurring problem, which is much faster than the aforementioned methods [3, 4, 5, 6, 7, 8, 9, 10, 11, 16, 12]. The proposed algorithm comprises two phases, i.e., we first preliminarily restore sharp images gradually from the coarsest scale to the finest scale, and then perform a refinement in the finest scale. In each scale, we introduce a joint optimization to restore both the latent sharp image and the blur kernel. Specifically, to achieve lower complexity, we alternately apply a sharp image reconstruction with fast local self-example matching inspired by image SR [15], an accelerated kernel estimation with error compensation, and a fast non-blind image deblurring. Furthermore, beyond uniform image blur, we extend the proposed algorithm to solve a more challenging non-uniform blind image deblurring problem. Experimental results demonstrate that in both cases our algorithm achieves the state-of-the-art deblurring performance.

The contributions of this paper can be summarized as follows:

1. We introduce a Multi-Scale Latent Structure (MSLS) prior of latent sharp image for deblurring. We theoretically verify this prior under a general down-sampling operation.

2. We propose a powerful and efficient blind image deblurring algorithm based on the MSLS prior, which achieves the state-of-the-art deblurring performance and is much faster than the previous methods with computationally expensive image priors.

3. We further demonstrate that the proposed algorithm can be extended to tackle more complex non-uniform motion blur, apart from simply uniform blur.

The outline of the paper is as follows. In Section II, we review related work in image deblurring. In Section III, we provide a theoretical analysis of our observation and introduce the MSLS prior. In Section IV, we present the proposed blind deblurring algorithm in detail. Experiments and final conclusion are provided in Section V and Section VI, respectively.

II Related Work

Image deblurring is an ill-posed inverse problem, which is to recover a latent sharp image ${\mathbf{x}}$ from a blurry observation ${\mathbf{b}}$. Based on whether the blur kernel ${\mathbf{k}}$ is known or not, the problem has been divided into two categories, i.e., non-blind image deblurring and blind image deblurring.

II-A Non-blind Image Deblurring

For non-blind image deblurring, the blur kernel is given and the problem is to recover the latent sharp image from the blurry observation with the blur kernel. Non-blind image deblurring is a very unstable process and is easily disturbed by noise, even a small amount of noise will lead to severe distortions in the estimation. This phenomenon is called ill-posedness. Many researches have been done to deal with the ill-posed problem including classical filtering methods [18, 19, 20] and regularization based methods [21, 22, 23, 24, 25, 26]. In the regularization based methods, regularization terms are usually proposed based on the prior of the latent sharp image, such as Total Variation (TV) prior [21] or sparse priors [23, 24]. By combining an $l_{2}$-norm or an $l_{1}$-norm [27] constrained data fidelity term with regularization terms, image deblurring can be modeled as an optimization problem. A latent sharp image can be recovered by solving the optimization problem. Besides Gaussian noise, non-blind deblurring is also disturbed by outliers, e.g., pixel saturation. Interested readers can refer to [28, 29, 30] for more details.

II-B Blind Image Deblurring

For blind image deblurring, the blur kernel is unknown and the problem is to recover both the latent sharp image and the blur kernel from only a blurry observation. The blind image deblurring is much more challenging than the non-blind problem because the solution for the problem is not only unstable but also non-unique.

In order to deal with this highly ill-posed blind image deblurring problem, there have been many pioneering works which attempted to approximate the problem in different perspectives. In early works, blind image deblurring has been processed like image enhancing. Image diffusion and shock filtering methods [31, 32, 33] were proposed to enhance the edges of blurry images blurred by out-of-focus or Gaussian kernel. Another kind of approximation is to introduce parameterized forms of blur kernels or to impose specific constraints, such as centrosymmetry, on the blur kernels [34, 35, 36, 37, 38]. These methods can handle specific blur kernels, but cannot handle general motion blur which is usually irregular. There are also numerous works focusing on the multi-image deblurring [39, 40, 41, 42, 43, 44]. Multiple images can provide more data constraints for kernel estimation and help to solve the blind image deblurring problem. Nevertheless, multiple relevant images in their assumptions are sometimes unavailable, which limits their practicability.

Recently, with the fast progress of regularization and optimization theory, many sophisticated image priors [3, 4, 5, 6, 7, 8, 9, 11, 10, 12] were proposed to handle the single image blind image deblurring problem with general blur kernels, such as the mixture of Gaussians prior that fits the heavy-tailed prior of natural images [3, 4], normalized sparse prior [5], framelet based prior [6], $l_{0}$-norm based priors [8, 7], color line prior [9], dark channel prior [10], low rank prior [11] and graph based prior [12], etc. The priors promote image sharpness and penalize image blurriness, which work as a regularizer in the optimization model guiding the solver to converge to the latent sharp image. Nevertheless, these image priors are usually non-convex, which result in very computationally expensive optimization algorithms.

III Observation and Multi-scale Latent Structure Prior

In this section, we introduce in detail the proposed multi-scale latent structure prior.

III-A Observation and Theoretical Analysis

Michaeli and Irani [16] first analyzed the effect of image down-sampling in the blind image deblurring problem and found that image down-sampling increased the internal patch recurrence. Differently, in this work, we focus on the extreme case of image down-sampling and down-sample a blurry image to a very coarse level. We observe that a coarse enough image down-sampled from a blurry observation is approximately a low-resolution version of the latent sharp image. We provide a theoretical proof on 2D signals to verify this observation under a general down-sampling operation. In the proof, we assume a noise-free blurred image. We empirically demonstrate that the observation is also robust with noise, by deblurring noisy and blurry images in the experiments in Section V.

Claim 1**.**

A coarse enough image down-sampled from a blurry observation is approximately a low-resolution version of the latent sharp image.

Proof: According to the degradation model of (1), a noise-free blurry image ${\mathbf{b}}$ can be represented as

[TABLE]

where ${\mathbf{x}}$ is the latent sharp image and ${\mathbf{k}}$ is the blur kernel.

When ${\mathbf{b}}$ is $\alpha$-times down-sampled, its down-sampled version is given by

[TABLE]

where ${\mathbf{b}}_{\alpha}[i,j]$ is a down-sampled blurry image, ${\mathbf{h}}$ is a general low-pass filter to avoid aliasing, ${\mathbf{x}}_{h}$ is a filtered band-limited image, $\downarrow$ represents the pixel extraction process and $\alpha$ is the down-sampling factor. By replacing $u$, $v$ with $\alpha t$, $\alpha s$, (3) can be further written as

[TABLE]

where ${\mathbf{x}}_{\alpha}$ and ${\mathbf{k}}_{\alpha}$ are the down-sampled sharp image and the down-sampled blur kernel, respectively. From (4), the down-sampled blur kernel is:

[TABLE]

The normalization condition $\iint{\mathbf{k}}_{\alpha}=1$ is satisfied if the assumption $\iint{\mathbf{k}}=1$ is imposed. It follows that the support of down-sampled blur kernel ${\mathbf{k}}_{\alpha}$ is $\alpha$ times smaller than that of blur kernel ${\mathbf{k}}$ in both horizontal and vertical directions.

We further transform ${\mathbf{k}}_{\alpha}$ into the frequency domain:

[TABLE]

where ${\mathcal{F}}(\cdot)$ denotes the Fourier transform and ${\mathbf{K}}(\tilde{u},\tilde{v})={\mathcal{F}}({\mathbf{k}}(u,v))$. ${\mathbf{K}}_{\alpha}$ is expanded by $\alpha$ in the frequency domain. Though ${\mathbf{K}}_{\alpha}$ is expanded, ${\mathbf{x}}_{h}$ and ${\mathbf{x}}_{\alpha}$ are low-pass filtered and band-limited. Thus, aliasing is avoided while down-sampling but only the information within the band-limit is preserved. As the down-sampling becomes deeper, i.e., $\alpha$ is larger, ${\mathbf{K}}_{\alpha}$ in the band-limit tends to constant $1$ and its corresponding spatial response will converges to Dirac delta function $\delta$ [45].

The proof so far works in the continuous domain, i.e., ${\mathbf{k}}$ and ${\mathbf{x}}$ are two continuous signals. For discrete signals in practise, we sample the continuous signals and do normalization on ${\mathbf{k}}$ afterwards to ensure $\sum_{uv}{\mathbf{k}}_{\alpha}(u,v)=1$. Approximately, ${\mathbf{k}}_{\alpha}(u,v)$ will converge to a delta kernel in the discrete domain:

[TABLE]

which satisfies $\sum_{i}\sum_{j}\delta(i,j)=1$. An illustrative experiment in the discrete domain is shown in Fig. 1. Finally, we have

[TABLE]

i.e., the down-sampled blurry image converges to the low-resolution version of the latent sharp image.

III-B Multi-Scale Latent Structure Prior

Claim 1 is very useful for blind image deblurring, since it gets rid of the effect of the unknown blur kernel and provides a low-resolution prior of the latent sharp image, referred as the latent structure prior in our algorithm.

Based on Claim 1, we further introduce a multi-scale strategy to work in tandem with latent structure prior on an image pyramid. Starting from the coarsest scale, the prior image is employed to blindly restore a finer scale sharp image. The newly restored sharp image is then employed as the latent structure prior of the next finer scale. The strategy works gradually from the coarsest scale to the finest scale, thus we name the coarse-to-fine prior—Multi-Scale Latent Structure (MSLS) prior. In the next section, we will elaborate the blind image deblurring algorithm tailored to the proposed MSLS prior.

IV Blind Image Deblurring Algorithm

In this section, we introduce a joint optimization scheme to estimate both the latent sharp image and the blur kernel based on the proposed MSLS prior. The algorithm comprises two phases, i.e., preliminary restoration in coarse scales in Sec. IV-A and refined restoration in the finest scale in Sec. IV-B. We further extend the proposed algorithm to solve non-uniform blind image deblurring problem in Sec. IV-C.

IV-A Preliminary Restoration in Coarse Scales

Blind image deblurring is modeled to solve the following optimization problem,

[TABLE]

where $\Phi({\mathbf{k}}\otimes{\mathbf{x}}-{\mathbf{b}})$ is the data fidelity term, $\Psi_{1}({\mathbf{k}})$ and $\Psi_{2}({\mathbf{x}})$ are the regularizers of blur kernel ${\mathbf{k}}$ and latent sharp image ${\mathbf{x}}$, respectively. Previous methods, such as [3, 4, 5, 6, 7, 8, 9, 10, 11, 16, 12], introduced sophisticated image priors, which are usually either non-convex or computationally expensive.

To tackle the challenging problem mentioned above, we take advantage of the proposed MSLS prior and introduce an efficient local self-example matching strategy to substitute for complex non-convex regularization of ${\mathbf{x}}$ in sharp image reconstruction. In each coarse scale, our algorithm comprises three steps, as sketched in Algorithm 1.

Specifically, in step (a), we blindly reconstruct a sharp image with fast local self-example matching. Step (b) and (c) correspond to a novel kernel estimation with error compensation and a non-blind image deblurring, which are both convex optimization and can be efficiently solved with Fast Fourier Transform (FFT) acceleration. The prior image is first initialized as the coarsest scale image ${\mathbf{x}}_{\alpha}$. After solving one scale, we update the prior image with the newly deblurred result and continue to restore the next finer scale. The diagram of preliminary restoration is illustrated in Fig. 2. In the following, we introduce in detail these three steps:

Step (a): Sharp Image Reconstruction with Local Self-Example Matching: The finer scale latent sharp image is reconstructed with a fast local self-example matching strategy, inspired by [15]. The latent sharp image ${\mathbf{x}}_{l}$ is first initialized as the blurry image ${\mathbf{b}}$. We divide the latent image ${\mathbf{x}}_{l}$ into overlapped small patches (usually of size $5\times 5$ pixels with $50\%$ overlap). Then, each patch $i$ finds its nearest neighbor (NN) in the corresponding local area of the coarser-scale sharp prior image ${\mathbf{x}}_{pr}$, measured by (10),

[TABLE]

where $\mathbf{P}_{i}$ is the patch extraction operation at position $i$. $i{{}^{\prime}}$ is a position on the sharp prior image ${\mathbf{x}}_{pr}$. If ${\mathbf{x}}_{l}$ is down-sampled to the size of ${\mathbf{x}}_{pr}$, the position $i$ on ${\mathbf{x}}_{l}$ will be projected on the position $i^{\prime}$ on ${\mathbf{x}}_{pr}$. $\mathcal{N}(i{{}^{\prime}})$ represents the local neighborhood of the position $i{{}^{\prime}}$. $d(i,j)$ denotes the distance between the patch $i$ and the patch $j$. The prior image ${\mathbf{x}}_{pr}$ is $\beta$-times smaller than ${\mathbf{x}}_{l}$, as shown in Fig. 2.

Finally, we fuse the searched NN patches together by computing weighted average in the overlapped areas, in order to avoid blocking artifacts. The weights are defined by Hamming window111https://en.wikipedia.org/wiki/Window_function. The 2D Hamming window function is defined as follows:

[TABLE]

where $\theta=0.54$ and $\gamma=0.46$ are fixed parameters, $N$ is the size of window, $w$ and $W$ are 1D and 2D weights of Hamming window, $i$ and $j$ are the indices. An illustrative example of a single-scale sharp image reconstruction is shown in Fig. 3.

In our algorithm, local patch matching works not as an approximation of the global patch matching but a better choice for sharp image reconstruction, which avoids the patch over-fitting and considerably reduces the computational complexity of the NN search:

IV-A1 Avoid Patch Over-fitting

In motion blurred images with different depths of field, the blurry patch in the foreground tends to find an over-fitting blurry NN in the background area of the prior image, if the searching distance is not constrained. This misleads the sharp image reconstruction and causes the incorrect kernel estimation in the followed step. We report a single-scale illustrative experiment in Fig. 4. When the search range is constrained by the locality property, more wanted sharp patches are found and the kernel estimation result is also improved.

IV-A2 Accelerate Searching

The computational complexity of local patch matching is $O(C)$, where $C$ is a constant that accounts for the search range. The search range is set to $10\times 10$ for $5\times 5$ sized patches in our algorithm. There are $(10-2\times\lfloor\frac{5}{2}\rfloor)\times(10-2\times\lfloor\frac{5}{2}\rfloor)=36$ times mean square error (MSE) computations at each location.

Step (b): Kernel Estimation with Error Compensation: Given the sharp image reconstructed with MSLS prior, we estimate the blur kernel afterwards. Since the finer-scale sharp image is reconstructed from the coarser-scale prior image, there is an inevitable error between the reconstructed image and the ground-truth.

In this step, we propose a new image blur model (13) for our kernel estimation, in order to compensate for the error:

[TABLE]

where ${\hat{\mathbf{x}}}$ is the reconstructed image and ${\mathbf{x}}_{c}$ is the compensatory layer. With the adjustable ${\mathbf{x}}_{c}$ layer, ${\hat{\mathbf{x}}}+{\mathbf{x}}_{c}$ can be equivalent to the latent sharp image ${\mathbf{x}}$ and the estimated kernel ${\hat{\mathbf{k}}}$ can be more accurate.

Given ${\hat{\mathbf{x}}}$ and model (13), computing ${\hat{\mathbf{k}}}$ and ${\mathbf{x}}_{c}$ results in solving an optimization problem (14). The model (13) is used as a fidelity term to constrain that the ${\hat{\mathbf{k}}}$ and ${\hat{\mathbf{x}}}_{c}$ follow the image convolution process. Instead of directly using (13), we transform it into gradient domain, in order to reduce ringing artifacts [46]. The regularizer of ${\mathbf{k}}$ is penalized by an $l_{2}$-norm and the regularizer of ${\mathbf{x}}_{c}$ is penalized by an $l_{1}$-norm to ensure the spatial sparsity of ${\mathbf{x}}_{c}$:

[TABLE]

where $\lambda_{1}$ and $\lambda_{2}$ are two trade-off parameters to control the regularization strength.

At first glance, (14) is a non-convex optimization problem, since there exists a ${\mathbf{k}}\otimes\nabla{\mathbf{x}}_{c}$. Here we apply a variable substitution to make it convex, i.e., replacing ${\mathbf{k}}\otimes\nabla{\mathbf{x}}_{c}$ with a new variable ${\mathbf{v}}$. Note that ${\mathbf{x}}_{c}$ is assumed sparse and the size of ${\mathbf{k}}$ is much smaller than the size of ${\mathbf{x}}_{c}$, so ${\mathbf{k}}\otimes\nabla{\mathbf{x}}_{c}$ can also be considered sparse. We thus substitute $\|{\mathbf{v}}\|_{1}$ for $\|{\mathbf{x}}_{c}\|_{1}$. Now the objective function (14) can be reformulated as a convex optimization problem (15):

[TABLE]

We alternately compute ${\mathbf{k}}$ and ${\mathbf{v}}$ by fixing one and addressing another. The $\nabla=\{\nabla_{x},\nabla_{y}\}$ are derivative operators in horizontal and vertical directions respectively. The ${\mathbf{v}}=\{{\mathbf{v}}_{x},{\mathbf{v}}_{y}\}$ are the corresponding layers of horizontal and vertical directions and are initialized to zeros.

In (15), terms related to ${\mathbf{k}}$ are all quadratic, so a closed-form solution can be derived for ${\mathbf{k}}$. Moreover, to fast solve ${\mathbf{k}}$, we compute the solution in the frequency domain,

[TABLE]

where $\mathcal{F}(\cdot)$ and $\mathcal{F}^{-1}(\cdot)$ are Fourier and inverse Fourier transform implemented by FFT. The $\overline{\mathcal{F}(\cdot)}$ means the complex conjugate of fourier transformed values.

As each element in ${\mathbf{v}}$ is independent, we use proximal mapping [47] to solve the $l_{1}$-norm optimization problem. The solver works as a fast soft-thresholding:

[TABLE]

where $i$ means the $i$-th pixel in ${\mathbf{v}}_{x}$ or ${\mathbf{v}}_{y}$, $z_{i}=\nabla b_{i}-({\mathbf{k}}\otimes\nabla{\hat{\mathbf{x}}})_{i}$.

An illustrative deblurring example is illustrated in Fig. 5. Fig. 5a is a blurry image. Fig. 5b and Fig. 5c show two blind deblurred results without and with the compensatory layer ${\mathbf{x}}_{c}$. Fig. 5d, Fig. 5e and Fig. 5f are their close-ups of the license plate, respectively. Solving with ${\mathbf{x}}_{c}$ achieves clearer deblurred result. The corresponding estimated kernel is also much better, in which less outliers appear beside the true trajectory.

We further analyze the role of the auxiliary variable ${\mathbf{v}}$ in optimization of (15). Fig. 5g is the close-up of the license plate in ${\hat{\mathbf{x}}}$, which is reconstructed from the coarser-scale prior image. It can be found that ${\hat{\mathbf{x}}}$ contains sharp edges but loses some details, e.g., the license number. To overcome this issue, we introduce the ${\mathbf{x}}_{c}$ to compensate for the error. In this case, ideally, ${\mathbf{x}}_{c}$ should be the license number. To solve both ${\mathbf{k}}$ and ${\mathbf{x}}_{c}$ in (14), we substitute a new variable ${\mathbf{v}}$ for ${\mathbf{k}}\otimes\nabla{\mathbf{x}}_{c}$ and solve (15). Fig. 5h and Fig. 5i are the ${\mathbf{v}}_{x}$ and ${\mathbf{v}}_{y}$ obtained by solving (15). It can be found that the restored ${\mathbf{v}}_{x}$ and ${\mathbf{v}}_{y}$ are the derivatives of the license number $\nabla{\mathbf{x}}_{c}$ convolving with the kernel ${\mathbf{k}}$, which demonstrates that the reconstructed result matches the theoretical analysis.

Step (c): Non-blind Image Deblurring: After the kernel ${\hat{\mathbf{k}}}$ is computed, we conduct the non-blind deblurring to update the latent image ${\mathbf{x}}_{l}$. We apply the Total Variation (TV) regularized non-blind deblurring algorithm, which is formulated as follows:

[TABLE]

where $\mu$ is a regularization parameter. We apply the alternating direction method of multipliers (ADMM) [48, 49] to lead to a fast solver in the frequency domain.

IV-B Refined Restoration in the Finest Scale

Preliminary restoration can fast estimate sharp images in the coarse scales but possibly results in artifacts. In order to achieve more accurate restoration for complex cases, we propose to perform a refined restoration in the finest scale, inspired by [50]. Given the preliminary restored image as initialization, the diagram of refined restoration in the finest scale is shown in Fig. 6.

Since the refined restoration is employed in the finest scale, there is no prior image updated from coarse scales. Thus, we update the prior image by filtering and down-sampling the restored image in each iteration. We use the edge-preserving guided filter [51] to reduce the possible artifacts without introducing any extra blur. Besides, we also employ higher-order derivatives in the data fidelity term to add more constraints for accurate kernel estimation [52]:

[TABLE]

where \nabla_{*}$$={$\nabla_{x}$, $\nabla_{y}$, $\nabla_{xx}$, $\nabla_{yy}$, $\nabla_{xy}$} and {\mathbf{v}}_{*}$$={${\mathbf{v}}_{x}$, ${\mathbf{v}}_{y}$, ${\mathbf{v}}_{xx}$, ${\mathbf{v}}_{yy}$, ${\mathbf{v}}_{xy}$}. The algorithm is outlined in Algorithm 2.

Fig. 7 illustrates an example of our refined restoration for a challenging blur kernel with fine details. Through preliminary estimation and refinement, the kernel becomes more accurate and the artifacts are reduced in the deblurred image. As a comparison, Fig. 7d shows the result of [16], which converges to a sub-optimal solution.

IV-C Extension to Non-Uniform Deblurring

The non-uniform image deblurring is a more challenging problem than the uniform one, which includes complicated camera motions, such as rotation and translation. Specifically, by representing the camera motion as a more general projective transformation [53, 54], the blur process can be modeled as:

[TABLE]

where ${\mathbf{x}}$ and ${\mathbf{b}}$ are vectors of the latent sharp image and the blurry observation, respectively. ${\mathbf{H}}_{i}$ is the $i$-th projective transformation matrix of the camera, $k_{i}$ is the weight of the $i$-th transformation in the blur kernel ${\mathbf{k}}$, ${\mathbf{n}}$ is the additive noise. (20) is a linear transformation, which can be further represented as:

[TABLE]

where matrix ${\mathbf{A}}_{\mathbf{k}}=\sum_{i}k_{i}{\mathbf{H}}_{i}$ and $col({\mathbf{B}}_{\mathbf{x}})_{i}={\mathbf{H}}_{i}{\mathbf{x}}$. $col(\cdot)_{i}$ denotes the $i$-th column.

Our algorithm can be extended to non-uniform deblurring by assuming that non-uniform image blur is locally uniform. Based on this assumption, the MSLS prior still works and we apply the same sharp image reconstruction with local self-example matching. For the kernel estimation in each scale, we replace the convolution model (1) with (21) and modify (14) to the following formulation:

[TABLE]

where $col({\mathbf{B}}_{\nabla({\hat{\mathbf{x}}}+{\mathbf{x}}_{c})})_{i}={\mathbf{H}}_{i}(\nabla({\hat{\mathbf{x}}}+{\mathbf{x}}_{c}))$.

We replace ${\mathbf{B}}_{\nabla{\mathbf{x}}_{c}}{\mathbf{k}}$ with ${\mathbf{v}}$ and replace $\|{\mathbf{x}}_{c}\|_{1}$ with $\|{\mathbf{v}}\|_{1}$ like (15). Then, to solve ${\mathbf{k}}$ results in a linear equation:

[TABLE]

where $col({\mathbf{B}}_{\nabla{\mathbf{x}}_{c}})_{i}={\mathbf{H}}_{i}(\nabla{\mathbf{x}}_{c})$, $col({\mathbf{B}}_{\nabla{\hat{\mathbf{x}}}})_{i}={\mathbf{H}}_{i}(\nabla{\hat{\mathbf{x}}})$, and $\mathbf{I}$ is an identity matrix. Unlike the fast solver in the frequency domain (16), we use conjugate gradient method to solve (23), because non-uniform blur no longer follows the convolution theorem in the frequency domain. We apply the same soft-thresholding to solve ${\mathbf{v}}$ by replacing $z_{i}=\nabla b_{i}-({\mathbf{k}}\otimes\nabla{\hat{\mathbf{x}}})_{i}$ with $z_{i}=\nabla b_{i}-({\mathbf{B}}_{\nabla{\hat{\mathbf{x}}}}{\mathbf{k}})_{i}$ in (17).

For TV regularized non-blind deblurring, (18) becomes:

[TABLE]

which can also be solved with ADMM. Similar modification can be done in the refined restoration of the finest scale.

V Experiments and Discussions

In this section, we have done extensive experiments to demonstrate the performance of our algorithm on the artificial dataset of Sun et al. [55], real uniformly blurred images and non-uniformly blurred images. Our experimental platform is a Windows 7 desktop computer with Intel i5 CPU. For uniform deblurring, 8G memory is more than enough. For non-uniform deblurring, we increase the memory to 24G. We use Matlab R2014b to run and test all the matlab codes. We tune the parameters on the dataset of Sun et al. and find that they generate satisfactory results for all other images. The parameters are set as follows. That is, $\lambda_{1}=\lambda_{3}=5$, $\lambda_{2}=\lambda_{4}=0.05$, $\mu=0.01$. Besides, the down-sampling factor $\beta$ is set to $\log_{2}3$ for the image pyramid construction. We downsample a blurry image scale by scale until its corresponding kernel is a delta function of size $1\times 1$. The $max\_iteration=3$ in the preliminary restoration and the refined restoration as a trade-off between accuracy and speed.

V-A Artificial Blurry Dataset

In this experiment, we compare our algorithm with previous blind image deblurring algorithms on a large dataset introduced by Sun et al. [55]. This dataset includes 640 blurry images (typically $1024\times 768$), which are made by performing convolutions between 80 sharp nature images and 8 blur kernels offered by [4]. Each image is then added $1\%$ white Gaussian noise. The blur kernels are assumed unknown and the kernel sizes are set to $51\times 51$ for all test cases.

In Fig. 8 and 9, we show the visual results deblurred by our algorithm and eight competing algorithms, i.e., Krishnan et al. [5], Lai et al. [9], Pan et al. [10], Michaeli & Irani [16], Cho & Lee [46], Xu & Jia [50], Sun et al. [55], Cho et al. [56]. For Pan et al. [10], the results are generated with their code. For other algorithms, the results are provided by their authors or by Sun et al. in the dataset. As shown in Fig. 8 and 9, our algorithm performs visually better than other algorithms. The edges and details are well restored and the artifacts in our results are much less than those in others.

Besides the visual quality assessment, we also measure the quality of the deblurred results using error ratio, which was introduced by [2]. All the algorithms estimate blur kernels and use the same non-blind deblurring algorithm to restore the latent sharp images. The error ratios are then computed by

[TABLE]

where ${\mathbf{x}}_{{\hat{\mathbf{k}}}}$ is the deblurred result computed using the estimated kernel ${\hat{\mathbf{k}}}$, ${\mathbf{x}}_{{\mathbf{k}}}$ is the deblurred result computed using the ground-truth kernel ${\mathbf{k}}$, ${\mathbf{x}}$ is the ground-truth sharp image. If $r=1$, it means that the blindly deblurred result is as good as the result restored using the ground-truth kernel. The smaller $r$ is, the better the blindly deblurred result is. We also agree with [16] that the visual quality is satisfactory when $r\leq 5$ and consider it as the threshold to decide the success. For all the competing methods except [10], of which the results are provided by their authors or by Sun et al., they estimated the blur kernels and used non-blind deblurring algorithm [25] to restore sharp images. For Pan et al. [10], we generate deblurred images with their code. Their code estimates the blur kernels and uses the non-blind method proposed in [8] to restore sharp images. In order to fairly compare with all the algorithms, we run our algorithms on both settings (Using [25] and [8] for the last-step non-blind deblurring, respectively). Fig. 10 reports the cumulative performance of error ratio. Each curve represents the fractions of 640 images that can be deblurred under different error ratios. For error ratios of results restored using [25] (solid lines), the curves of Sun et al. [55] and Lai et al. [9] are slightly higher than ours when $error\ ratio<2.2$. Since all results in this range are very well recovered, the visual differences between their results and ours are hardly noticed. When $error\ ratio\geq 2.2$, our curve achieves the highest and the improvement is obvious. For error ratios of results restored using [8] (dash lines), the two curves are very close at the beginning and our curve soon becomes higher than that of Pan et al. [10].

In Table I, we further report three statistical measures, i.e., mean error ratio, worst error ratio and success rate, to further compare all the blind deblurring algorithms. According to the table, the mean errors of our algorithm are the smallest and the success rates of our algorithm are the highest among all the algorithms on both settings. Although our worst error ratios are larger than Michaeli & Irani [16] and Pan et al. [10], our second worst results are better on both settings.

In Fig. 11, we show an example of an unsatisfactory blind image deblurring result. In our approach, we require that there are objects with profiles, of which the sizes are much larger than that of the blur kernel. However, in this case, smooth sky accounts for most of the image. Only a windmill and foothills in the foreground provide profiles. The small windmills in the background will be removed after down-sampling. Thus, there aren’t enough constraints in the data fidelity term for the complex blur kernel restoration, leading to the unsatisfactory blind deblurring result.

V-B Real Blurry Images

In the second part of this section, we apply our algorithms on real blurry images and compare with algorithms, of which the implementations are available, i.e., Krishnan et al. [5], Levin et al. [4], Michaeli & Irani [16] and Pan et al. [10]. In Fig. 12, 13, 14 and 15, we blindly deblur an image with different depths of field, a human portrait with a complex background, an image full of humans, and an image with very severe motion blur, respectively. Our algorithm estimates the blur kernels robustly and results in less artifacts in the recovered images.

Table II reports the running time of the competing algorithms and ours on real blurry images in Fig. 12, 13, 14 and 15. All the algorithms run on the Matlab platform. According to Table II, our algorithm runs much faster than other competing algorithms. The efficiency of our algorithm can be further improved with C++ and GPU implementation.

V-C Non-Uniform Blurry Images

In Fig. 16 and 17, we demonstrate the performance of our algorithm on the non-uniform blind deblurring problem. Directly using uniform deblurring algorithm fails to restore the correct sharp image, since the existing blur does not follow the convolution model. as shown in Fig. 16b. Extended to non-uniform blur model (21), our algorithm can tackle the spatial-variant blurry cases. We compare our algorithm with recent non-uniform blind deblurring algorithms [57, 53, 7, 8], of which the results are available. Our results are visually comparable or better than those of the competing methods.

VI Conclusion

In this paper, we propose a Multi-Scale Latent Structure prior, which is derived from down-sampling a blurry image. With this prior, we design an efficient algorithm to jointly restore both the latent sharp image and the blur kernel from a single blurry observation. The qualitative and quantitative experiments demonstrate that our algorithm is competitive against the state-of-the-art methods in blindly recovering uniform and non-uniform blurry images.

Limitation: In our approach, we require that there are objects with profiles, of which the sizes are much larger than that of the blur kernel. If the precondition is disobeyed, objects will be removed after down-sampling. In this case, the data fidelity term of the optimization function fails to work and will converge to an uncertain kernel. Nevertheless, to the best of our knowledge, it is also difficult for other blind deblurring algorithms to deal with the blurry images, in which all objects are smaller than blur kernels.

Recently, there have been methods focusing on specific object deblurring problems, for example, text or face deblurring. These methods can work better by considering the properties of specific objects. In the future, we would like to incorporate the advantage of these techniques with our algorithm to make a more robust unified framework. Besides, although the parameter kernel size in our algorithm is not sensitive and can be set relatively large, it is still more preferred to have it automatically set according to the individuality of the blurry images. In the future, we would like to investigate how to infer the optimal kernel size from a blurry image to make our algorithm more practical.