# Blind Image Deblurring based on Kernel Mixture††thanks: This work has been
submitted to the IEEE for possible publication. Copyright may be transferred
without notice, after which this version may no longer be accessible.

Sajjad Amrollahi Biyouki
Department of Industrial and Systems Engineering
The University of Tennessee
Knoxville, TN 37996
<EMAIL_ADDRESS>
&Hoon Hwangbo
Department of Industrial and Systems Engineering
The University of Tennessee
Knoxville, TN 37996
<EMAIL_ADDRESS>

###### Abstract

Blind Image deblurring tries to estimate blurriness and a latent image out of
a blurred image. This estimation, as being an ill-posed problem, requires
imposing restrictions on the latent image or a blur kernel that represents
blurriness. Different from recent studies that impose some priors on the
latent image, this paper regulates the structure of the blur kernel. We
propose a kernel mixture structure while using the Gaussian kernel as a base
kernel. By combining multiple Gaussian kernels structurally enhanced in terms
of scales and centers, the kernel mixture becomes capable of modeling nearly
non-parametric shape of blurriness. A data-driven decision for the number of
base kernels to combine makes the structure even more flexible. We apply this
approach to a remote sensing problem to recover images from blurry images of
satellite. This case study shows the superiority of the proposed method
regulating the blur kernel in comparison with state-of-the-art methods that
regulates the latent image.

_K_ eywords Blind deconvolution $\cdot$ Gaussian kernel $\cdot$ Mixture of
kernels $\cdot$ Remote Sensing $\cdot$ Image Restoration

##  1 Introduction

Image restoration is widely used when images do not provide desired quality in
terms of clarity and contrast of target objects as a consequence of noisy
disturbance so-called blurriness. Such poor-quality images are often observed
in remote sensing [23],[30], underwater objects detection [12],[19], and
healthcare applications [9],[28]. Due to great need in many real-world
applications, a wide spectrum of methods have been developed to restore a
“best-quality” latent image from a blurred image [8],[7],[22].

In general, a blur process is modeled as the convolution of a latent image and
blurriness represented by a kernel. The deconvolution method, the inverse
process of the convolution, is the process of extracting the latent image by
deconvolving a blurred image into the latent image and the blur kernel. If
some prior knowledge is given for the blur kernel so if a specific kernel can
be assumed a priori, the deconvolution process becomes straightforward as it
only requires estimating the latent image. This type of approaches that use a
pre-defined kernel are called non-blind deconvolution [21],[5]. The kernels
are, however, unknown or only partially known in real-world problems. With the
rapid increase in the usage of images for an analysis of various systems,
assuming a pre-defined kernel is often too restrictive. In this regard, the
blind deconvolution [1],[11] (also referred to as blind image deblurring) that
estimates both the latent image and the blur kernel has attracted great
attentions in the last two decades.

The blind deconvolution process requires solving an ill-posed problem. This is
because, for a given, single blurred image, there can be an infinite number of
solutions for the latent image and the blur kernel that satisfy the system of
equations defined for the blur process model:

$\mathbf{B}=\mathbf{I}*\mathbf{K}+\mathbf{n}$ (1)

where $*$ is the convolution operator, $\mathbf{B}$ is a degraded (blurred)
image, $\mathbf{I}$ is a recovered (latent) image, $\mathbf{K}$ is a blur
kernel, and $\mathbf{n}$ is an additive noise. Therefore, to derive a unique
solution for the blind deconvolution problem, additional restrictions need to
be imposed. For this purpose, some priors and regularization terms have been
used to maintain and intensify image edges (restricting latent images) or to
diminish harmful degradation and noise (restricting blur kernels). Low-rank
prior [20], dark channel prior [18], and graph-based prior [4] have been
developed recently to restrict the latent image, and total variation (TV)
regularization [6] and gradient prior [26] have been applied to either the
blur kernel or the latent image.

Recently, the approaches restricting the latent image, e.g., through
statistical priors [8],[17] or image related priors [20],[18],[4], have been a
mainstream in blind image deblurring research. However, image characteristics
can vary by different images and applications, and there is no guarantee that
the restrictions justified for a specific latent image work well for other
blind deconvolution problems. On the other hand, blurriness is typically
originated from common sources, including atmospheric turbulence (mostly for
remote sensing), out-of-focus, camera shake or object motion [27],[10]. This
suggests that restricting the blur kernel rather than the latent image can be
more effective and generalized more easily to broad applications of the blind
deconvolution. In this paper, we propose a novel kernel structure to model the
blur kernel, which can be leveraged for any blurry image no matter what image
prior (on the latent image) is appropriate.

A Gaussian kernel has been used to model the blur kernel, especially to
estimate atmospheric turbulence blurs [27]. In general, a circular shape
modeled by a simple Gaussian blur kernel is not sufficient to represent the
underlying blurriness, and there could be multiple sources of blurriness
making the overall shape complicated. To address this problem, we propose
using a mixture of several Gaussian kernels as a blur kernel while allowing
different scales, centers, and rotations of each individual kernel. This
kernel mixture is capable of modeling a complex shape of blurriness from
symmetric to asymmetric and from circular to linear without significant
limitation. Consequently, the kernel mixture shows flexible behaviors as if it
were estimated nonparametrically even with its parameterized structure. All
the decisions to define the kernel structure, including how many base kernels
to combine, are data-driven, so the proposed kernel mixture is adaptive and
can be applied to any blind deconvolution problem.

The main contributions of this paper can be summarized as:

  * •

This paper develops a novel kernel structure, a kernel mixture, that can be
applied to a broad class of blind image deblurring problems, independent of
the characteristics of latent images.

  * •

The parametric structure of the proposed kernel, induced by Gaussian base
kernels, restricts and characterizes the blur kernel effectively producing a
good solution for the ill-posed blind deconvolution problem.

  * •

The proposed kernel is flexible in modeling blurriness; with different scales,
centers, and rotations, the Gaussian kernels become capable of modeling
various shapes of blurriness.

  * •

The proposed kernel is adaptive to given images since the determination of its
structure is data-driven.

The rest of the paper is organized as follows. Section 2 reviews the related
works in the blind image deblurring domain. Section 3 describes the proposed
method in detail elaborating the development of the kernel mixture, the
optimization of associated parameters, and the overall process of blind
deconvolution based on the kernel mixture. Section 4 presents a case study of
deblurring noisy satellite images, discusses the dataset and experimental
settings, and compares the proposed method with other state-of-the-art
benchmark methods. Section 5 discusses future research directions and
concludes the paper.

##  2 Related Works

Image deconvolution methods can be grouped into two general categories: the
non-blind deconvolution where the kernel information is known and the blind
deconvolution where the kernel is also unknown and subject to estimation. In
the non-blind deconvolution domain, Wiener filter [25] and Richardson-Lucy
algorithm [21] are the most well-known methods among earlier works, and they
are still in use for the image restoration problems. The major shortcoming of
these methods is in their noise sensitivity that leads to ringing artifacts in
the recovered image. In addition, these methods require assuming a specific
kernel, but it is hard to find a proper kernel that works well for different
images/applications. Albeit more difficult to implement, the blind
deconvolution has been used more broadly with better capability of image
recovery in general. Earlier studies mostly focused on removing motion blur
[7],[8],[22] caused by dynamic movement of an object while an image is taken.
Recent papers also considered other types of blurriness stemmed from various
sources such as atmosphere turbulence and camera shake [4],[27].To solve a
blind deconvolution problem, some Bayesian approaches, specifically Maximum a
Posteriori (MAP) estimation, and other optimization techniques have been used.

A few decades ago, Likas and Galatsanos [13] proposed using a hierarchical
Bayesian modeling for blind deblurring. They used Gaussian probability
distribution for modeling image prior, blur kernel, and hyperparameters of the
priors. They employed the variational Expectation Maximization (EM) to obtain
the MAP estimates from their Bayesian model. Fergus et al. [8] used Miskin’s
ensemble learning [15] for a variational Bayesian approach while assuming a
Gaussian mixture prior for the latent image. Inspired by this study, other
researchers proposed more efficient approaches by considering various priors
[2],[3],[16]. Babacan et al. [2] employed the variational Bayesian approach
while assuming sparse priors for the latent image (super-Gaussian priors).
Babacan et al. [3] imposed total variation prior on the latent image and
assumed a Gaussian blur kernel. Molina et al. [16] proposed simultaneous
autoregressions as priors for both the latent image and blur kernel and used
gamma distributions to model the hyperparameters of the priors. The major
weakness of this type of approaches based on the MAP estimation is their
strong dependency on the choice of priors and the lack of generality as a
consequence. It has been shown that the blind deconvolution tends to estimate
a trivial unblur image when an MAP approach is applied [2],[11]. In addition,
when sparse priors are used, the computational performance of an MAP
estimation exacerbates as the objective function for the estimation becomes
non-convex [26].

Others have used optimization techniques to solve a blind deconvolution
problem. The main idea is to solve an individual optimization problem for each
of the latent image and blur kernel while keeping one constant in the
estimation of the other and iteratively updating the estimates [29], [6]. You
and Kaveh [29] introduced such an alternating optimization problem in which
they regularized both the latent image and blur kernel by using the Laplace
operator. Chan and Wong [6] took advantage of the alternating optimization
structure and proposed the usage of total variation regularization, which can
improve recovering the edges of an image. Since then, the alternating
optimization-based approaches have been evolved in two different ways: i)
introducing novel image priors and ii) developing new blur kernel structures.

Most recent works have developed more sophisticated image priors, including
$l_{0}$-norm prior [17], low-rank prior [20], dark channel prior [18], and
most recently, reweighed graph-based prior [4]. Pan et al. [17] proposed using
the $l_{0}$-norm that regulates the number of nonzero pixels as it can
distinguish a clear image from a blurred image based on their opposite
behaviors in terms of nonzero intensities. Regulating the $l_{0}$-norm prior,
however, makes the objective function for estimating the latent image and blur
kernel non-convex, so the estimation problem becomes hard to solve [4]. Ren et
al. [20] enforced image prior to be low-rank approximation of a degraded image
by using a weighted nuclear norm minimization and combined the image prior
with the gradient map of the image. However, this requires solving pixel-based
singular value decomposition that has $O(N^{3})$ complexity [4]. Pan et al.
[18] observed that the dark channel of blurred images should be less sparse
than that of clear images due to the nature of a typical blur process. Based
on this observation, they proposed regulating the dark channel of the latent
image and making it sparse. This process though involves computing nonlinear
dark channel and its $l_{0}$-norm, so implementing the approach is
computational intensive [4]. Most recently, Bai et al. [4] incorporated a
graph structure to model a blurred image in nodal domain and used the graph
structure as the image prior. This was the first effort mapping pixels to a
domain other than frequency and real domains. All the discussed methods with
novel image priors have improved the quality of image recovery. However, they
may not be effective in extracting a latent image if the characteristics
assumed for the priors are not so obvious in a given image. In other words,
their quality can vary substantially depending on the images given.

In the domain of blur kernel regularization, Xu and Jia [26] utilized
influential edges of an image (image gradients) to create an initial kernel
and refined the kernel by using an iterative support detection. Although this
approach regularizes the blur kernel, the kernel estimation strongly relies on
the construction of edges of which characteristics vary by images. More
recently, other studies have imposed specific structures on the blur kernel by
combining multiple kernels of known structures. Mai and Liu [14] fused blur
kernels estimated from other methods by using Gaussian Conditional Random
Fields. Their approach outperformed the other methods from which the
individual kernels were extracted. However, this requires implementing all the
other methods, and the entire estimation process becomes computationally
expensive. Furthermore, its performance depends on the quality of individual
kernels. If none of the individual kernels is capable of capturing specific
properties of the true blurriness, the kernel fusion will also fail to model
such a property. Later, Lee and Hwang [10] modeled a blur kernel as a linear
combination of basic two-dimensional patterns. To construct the patterns, they
used “one-dimensional” Gaussian density with a different scale parameter for
each pattern. Since this blur kernel relies on simple Gaussian densities, it
cannot represent various shapes of blurriness. In this paper, we develop a
kernel mixture that combines structure-enhanced Gaussian kernels, which is
capable of modeling almost any shape of blurriness.

##  3 Methodology

In this section, we develop the kernel mixture of structure-enhanced Gaussian
kernels and present its capability of modeling general shapes of blurriness.
We also propose a blind image deblurring formulation that properly regulates
the latent image and blur kernel in the context of the kernel mixture. Then,
we describe the alternating optimization problem used to estimate the blur
kernel and latent image.

###  3.1 Mixture of Structure-Enhanced Gaussian Kernels

The major objective of this section is to develop a blur kernel that is
flexible in modeling blurriness while maintaining a certain parametric
structure. We use Gaussian kernels as base kernels to impose parametric
structures and improve the structures in terms of scales, centers, and
rotations to grant various shape characteristics when they are combined.

The general model of our proposed blur kernel $\mathbf{K}$ is defined as a
mixture of $N$ Gaussian kernels as

$\mathbf{K}=\sum_{t=1}^{N}K_{G,t}$ (2)

where $K_{G,t}$ for $t=1,\ldots,N$ is a base Gaussian kernel. The base kernels
will have different shape complexity which is determined by the underlying
blurriness of a given image. We present possible base kernel structures and
the shape of their mixture from the simplest to the most complicated ones with
their distinguished specifications.

The simplest structure of a two-dimensional base kernel is an isotropic zero-
mean Gaussian kernel with a scale parameter $\sigma$:

$K_{G}^{\text{simple}}(x,y;\sigma^{2})\propto\exp{(-\frac{x^{2}+y^{2}}{2\sigma^{2}})}$
(3)

where $x$ and $y$ are the locations of a pixel being evaluated, relative to a
target pixel location in the horizontal and vertical axes of an image,
respectively. For example, suppose pixel A in Fig. 1 is a target pixel. When
evaluating blurriness propagating from the target pixel to others, the
location of B relative to A is calculated by $(x,y)=(-1,1)$ as it moves to the
left and to the top each by one pixel. Similarly, the location of C and D are
$(1,0)$ and $(0,-1)$, respectively.

Figure 1: Relative locations of pixels

In Eq. (3), the kernel $K_{G}^{\text{simple}}$ is determined by the single
scale parameter $\sigma$. When using a single parameter to describe the scale
toward all directions, the propagation of blurriness is modeled as a circular
shape. Even the mixture of $K_{G}^{\text{simple}}$ will remain as a concentric
circular shape. As such, this kernel is too simple to represent general
blurriness.

(a) Blurriness shape with individual scales

(b) Scale-enhanced base kernel

(c) Mixture of $K_{G}^{\text{scale}}$

Figure 2: Kernel mixture with scale enhancement

To improve the simple structure, we allow a different scale parameter for each
axis and construct a covariance matrix to account for the varying scales of
different axes. A scale-enhanced Gaussian kernel is expressed in matrix form
as

$K_{G}^{\text{scale}}(\mathbf{p};\bm{\Sigma})\propto\exp{(-\frac{1}{2}\mathbf{p}^{T}\bm{\Sigma}^{-1}\mathbf{p})}$
(4)

where $\mathbf{p}=\begin{bmatrix}x&y\end{bmatrix}^{T}$,
$\bm{\Sigma}=\text{diag}(\sigma_{x}^{2},\sigma_{y}^{2})$, and $\sigma_{x}^{2}$
and $\sigma_{y}^{2}$ are the variances associated with each axis. By allowing
different scales on each axis, this kernel can be represented as an elliptic
shape as shown in Fig. 2(a). In an extreme condition, as one of
$\sigma_{x}^{2}$ or $\sigma_{y}^{2}$ gets close to zero, the elliptic shape
can converge to a linear shape since the kernel is defined on the discretized
grid formed by the pixels. A combination of $K_{G}^{\text{scale}}$ kernels
produces a centralized multi-layer cross-line shape as shown in Fig. 2(c).

(a) Blurriness shape with a non-zero center

(b) Center-enhanced base kernel

(c) Mixture of $K_{G}^{\text{center}}$

Figure 3: Kernel mixture with center enhancement

(a) Blurriness shape with rotation

(b) Rotation-enhanced base kernel

(c) Mixture of $K_{G}^{\text{rotation}}$

Figure 4: Kernel mixture with rotation enhancement

Still, the blur kernel constructed by the mixture of $K_{G}^{\text{scale}}$
cannot model asymmetric blurriness or multi-source blurriness because every
base kernel is symmetric and centered at zero. To overcome this weakness, we
propose using a non-zero center for each base kernel which is formulated as

$\displaystyle
K_{G}^{\text{center}}(\mathbf{p};\bm{\Sigma},\bm{\mu})\propto\exp{\left(-\frac{1}{2}(\mathbf{p}-\bm{\mu})^{T}\bm{\Sigma}^{-1}(\mathbf{p}-\bm{\mu})\right)}$
(5)

where $\bm{\mu}=\begin{bmatrix}\mu_{x}&\mu_{y}\end{bmatrix}^{T}$ and $\mu_{x}$
and $\mu_{y}$ are the relative locations of the center on each axis. Fig. 3
displays the structure of a Gaussian kernel with a non-zero center as well as
a mixture of the same kind. With a non-zero center, the
$K_{G}^{\text{center}}$ can be located anywhere within the range where the
kernel is defined. By combining multiple of them, the blur kernel as a whole
can present an asymmetric shape of blurriness as shown in Fig. 3(c). The
mixture is also capable of modeling a sparse blur kernel by combining multiple
base kernels with small scales. With this great flexibility, the ultimate
shape of the blur kernel will be determined by the underlying shape of
blurriness inherent in a blurred image. As such, this blur kernel exhibits
almost nonparametric behaviors.

One limitation of the mixture of $K_{G}^{\text{center}}$ is that the shapes of
the base kernels should be parallel to horizontal and vertical axes. We relax
this condition and further consider rotations of the base kernels by using a
rotation matrix $\mathbf{R}$ and a rotation angle $\theta$:

$K_{G}^{\text{rotation}}(\mathbf{p};\bm{\Sigma},\bm{\mu},\theta)\propto\exp{\left(-\frac{1}{2}(\mathbf{Rp}-\bm{\mu})^{T}\bm{\Sigma}^{-1}(\mathbf{Rp}-\bm{\mu})\right)},\mathbf{R}=\begin{bmatrix}\cos\theta&-\sin\theta\\\
\sin\theta&\cos\theta\end{bmatrix}.$ (6)

This formulation provides the most enhanced kernel structure, which leverages
all the aforementioned structural features. When $\theta=0$,
$K_{G}^{\text{rotation}}$ is equivalent to $K_{G}^{\text{center}}$. If
$\bm{\mu}$ is also a zero vector, the kernel becomes $K_{G}^{\text{shape}}$.
In fact, the shape complexity of each base kernel is determined by a given
image if we use the most advanced structure as a base kernel. The mixture of
$K_{G}^{\text{rotation}}$ can model nearly any shape of blurriness thanks to
its high flexibility as shown in Fig. 4(c).

###  3.2 Blind Deblurring Formulation

For blind image deblurring, both the latent image and blur kernel need to be
estimated. In a Bayesian framework, the MAP estimates of $\mathbf{K}$ and
$\mathbf{I}$ are determined by maximizing the joint posterior distribution,
$p(\mathbf{I,K|B})$, as

$\mathbf{I}^{*},\mathbf{K}^{*}=\underset{\mathbf{I,K}}{\mbox{arg\\!
max}}\hskip 2.84526ptp(\mathbf{I,K|B})=\underset{\mathbf{I,K}}{\mbox{arg\\!
max}}\hskip 2.84526ptp(\mathbf{B|I,K})p(\mathbf{K})p(\mathbf{I})$ (7)

while assuming independence between $\mathbf{K}$ and $\mathbf{I}$. This MAP
estimation requires specifying $p(\mathbf{K})$ and $p(\mathbf{I})$, the prior
distributions of the blur kernel and latent image, respectively.

Instead of assuming particular distributions for priors, we apply the
alternating optimization technique. To construct an optimization problem, we
modify the objective function of the MAP estimation in Eq. (7). By taking the
negative logarithm, Eq. (7) becomes

$\mathbf{I}^{*},\mathbf{K}^{*}=\underset{\mathbf{I,K}}{\mbox{arg\\!
min}}\hskip 2.84526pt-
log(p(\mathbf{B|I,K}))-log(p(\mathbf{K}))-log(p(\mathbf{I})).$ (8)

A more general optimization problem can be formulated by replacing the log
probability in Eq. (8) by some loss functions:

$\mathbf{I}^{*},\mathbf{K}^{*}=\underset{\mathbf{I,K}}{\mbox{arg\\!
min}}\hskip
2.84526pt\ell(\mathbf{I}*\mathbf{K},\mathbf{B})+\lambda_{1}\ell_{\mathbf{K}}(\mathbf{K})+\lambda_{2}\ell_{\mathbf{I}}(\mathbf{I})$
(9)

where the first term measures the loss occurring from estimating $\mathbf{B}$
as the convolution of the estimates of $\mathbf{I}$ and $\mathbf{K}$,
$\ell_{\mathbf{K}}$ and $\ell_{\mathbf{I}}$ are prior terms calculating some
loss associated with the estimates of $\mathbf{K}$ and $\mathbf{I}$,
respectively. $\lambda_{1}$ and $\lambda_{2}$ regularizes the prior terms
while presenting their relative importance to the first term, the estimation
loss.

Typically, the kernel prior, $\ell_{\mathbf{K}}(\mathbf{K})$ is used to
stabilize the blur kernel, and the image prior,
$\ell_{\mathbf{I}}(\mathbf{I})$, is used to recover the latent image with
sharp edges. These priors play a significant role in blind deblurring. While
this paper focuses on developing a blur kernel with a flexible structure, a
proper prior is still needed to control the level of flexibility in the blur
kernel. To achieve that, we solve the following minimization problem:

$\mathbf{I}^{*},\mathbf{K}^{*}=\underset{\mathbf{I},\mathbf{K}}{\mbox{arg\\!
min}}\hskip
2.84526pt||\mathbf{I}*\mathbf{K}-\mathbf{B}||_{2}^{2}+\lambda_{1}||\mathbf{K}||_{2}^{2}+\lambda_{2}||\bm{\sigma}^{2}||_{2}^{2}+\lambda_{3}(||\nabla\mathbf{I}||_{2}^{2}+||\mathbf{I}||_{2}^{2})$
(10)

where
$\bm{\sigma}^{2}=(\sigma_{x,1}^{2},\ldots,\sigma_{x,N}^{2},\sigma_{y,1}^{2},\ldots,\sigma_{y,N}^{2})^{T}$
is a vector including the diagonal elements
$(\sigma_{x,t}^{2},\sigma_{y,t}^{2})$ of the covariance matrix for each base
kernel $K_{G,t}$ for $t=1,\ldots,N$, $\nabla I$ is the image gradient, and
$||\cdot||_{2}$ is the $l_{2}$-norm. We use $||\mathbf{K}||_{2}^{2}$ and
$||\bm{\sigma}^{2}||_{2}^{2}$ as the kernel prior and
$||\nabla\mathbf{I}||_{2}^{2}+||\mathbf{I}||_{2}^{2}$ as the image prior.
Regulating $||\mathbf{K}||_{2}^{2}$ and $||\mathbf{I}||_{2}^{2}$ induces the
sparsity of $\mathbf{K}$ and $\mathbf{I}$, respectively. The inclusion of
$||\nabla\mathbf{I}||_{2}^{2}$ restricts the estimate of $\mathbf{I}$ by
eliminating tiny gradient segments but keeping large ones only. In addition to
these common priors, we propose to include $||\bm{\sigma}^{2}||_{2}^{2}$, say
a covariance prior, to further regulate the blur kernel.

Input: Degraded image $\mathbf{B}$, number of base kernels $N$, kernel size
$h$

Let $i\leftarrow 0$;

Initialize the latent image, $\mathbf{I}^{0}\leftarrow\mathbf{B}$;

Generate random numbers to initialize the model parameters,
$\bm{\Sigma}_{t}^{0}$ and $\bm{\mu}_{t}^{0}$ for $t=1,\ldots,N$;

Use the model parameters to initialize base kernels $K_{G,t}^{0}$ for
$t=1,\ldots,N$;

Initialize the blur kernel $\mathbf{K}^{0}$ by combining $K_{G,t}^{0}$ for
$t=1,\ldots,N$ according to Eq. (2);

repeat

Update the number of iteration, $i\leftarrow i+1$;

Blur kernel estimation steps:

  * Given $\mathbf{I}^{i-1}$, estimate $\mathbf{K}^{i}$ by optimizing Eq. (11) with respect to $\bm{\Sigma}_{t}^{i}$ and $\bm{\mu}_{t}^{i}$ for $t=1,\ldots,N$;

* Normalize $\mathbf{K}^{i}$;

Latent image restoration steps:

  * Given $\mathbf{K}^{i}$, recover an intermediate latent image $\mathbf{I}^{i}$ by using Eq. (13);

* Update a tuning parameter, $\lambda_{3}\leftarrow\lambda_{3}/1.1$;

until  _
$\frac{||\mathbf{K}^{i}-\mathbf{K}^{i-1}||_{2}}{||\mathbf{K}^{i-1}||_{2}}<\epsilon\quad\text{and}\quad\frac{||\mathbf{I}^{i}-\mathbf{I}^{i-1}||_{2}}{||\mathbf{I}^{i-1}||_{2}}<\epsilon$
_;

Output: Latent image $\mathbf{I}^{*}\leftarrow\mathbf{I}^{i}$

Algorithm 1 Alternating optimization for blind deblurring with the kernel
mixture

As done in typical coefficient shrinkage, e.g., ridge regression, the
covariance prior enforces insignificant $\sigma_{x,t}^{2}$ and
$\sigma_{y,t}^{2}$ to be close to zero and makes the corresponding base kernel
$K_{G,t}$ almost negligible in modeling blurriness. This property is used to
determine the number of base kernels for the kernel mixture. Instead of
predetermining the exact number of base kernels, we include a large enough
number of base kernels in the model. Then, some of them with little impact
will become negligible due to the covariance prior.

###  3.3 Estimation Procedure

To solve Eq. (10), we use alternating optimization as described in Alg. 1. The
ultimate purpose of Alg. 1 is to recover the latent image $\mathbf{I}$ by
modeling a proper blur kernel $\mathbf{K}$, given the degraded image
$\mathbf{B}$. To construct $\mathbf{K}$ as a kernel mixture, the information
about the number of base kernels, $N$, and the kernel size, $h$, is also
needed. The kernel size does not need to be the same with the size of
$\mathbf{B}$, but it will be sufficient if the blur kernel $\mathbf{K}$ is
large enough to model the blurriness in the degraded image $\mathbf{B}$. Once
the inputs are ready, we use $\mathbf{B}$ to initialize the latent image,
$\mathbf{I}^{0}$. Then, we generate random numbers to set the initial values
of model parameters, initialize the base kernels, and combine them to form the
initial blur kernel, $\mathbf{K}^{0}$. After the initialization, we keep
updating the blur kernel and latent image forming an iterative loop. Within a
loop, we optimize $\mathbf{K}$ given the latent image at hand,
$\mathbf{I}^{i-1}$, and use the resulting optimal solution to update
$\mathbf{K}^{i}$. The updated blur kernel, $\mathbf{K}^{i}$, is then used to
obtain a new estimate of the latent image, $\mathbf{I}^{i}$. This procedure
will be repeated until the estimates of the latent image and blur kernel
converge. In fact, the procedure in the loop consists of two sub-problems, one
for blur kernel estimation and another for latent image restoration, which is
further elaborated in the following sections.

####  3.3.1 Blur kernel

Under the alternating optimization framework, a latent image estimate is given
for the blur kernel estimation. Then, the latent image is assumed constant, so
the objective function in Eq. (10) reduces to

$E(\mathbf{K})=||\mathbf{I}*\mathbf{K}-\mathbf{B}||_{2}^{2}+\lambda_{1}||\mathbf{K}||_{2}^{2}+\lambda_{2}||\bm{\sigma}^{2}||_{2}^{2}.$
(11)

For the blur kernel estimation, we minimize the energy function,
$E(\mathbf{K})$, with respect to the model parameters that compose
$\mathbf{K}$, by using the Conjugate Gradient (CG) method. Once a new estimate
of $\mathbf{K}$ is obtained, it is normalized to ensure
$\sum_{i}\sum_{j}\mathbf{K}_{ij}=1$.

####  3.3.2 Latent image

The goal of this sub-problem is to recover the latent image given a blur
kernel estimate. Based on Eq. (10), the energy function is formulated as

$E(\mathbf{I})=||\mathbf{I}*\mathbf{K}-\mathbf{B}||_{2}^{2}+\lambda_{3}(||\nabla\mathbf{I}||_{2}^{2}+||\mathbf{I}||_{2}^{2}).$
(12)

For the minimization of Eq. (12), the closed-form solution exists. By using
the Fast Fourier Transform (FFT), the solution is formulated [26] as

$\mathbf{I}^{*}=\mathcal{F}^{-1}(\frac{\overline{\mathcal{F}(\mathbf{K})}\mathcal{F}(\mathbf{B})}{\overline{\mathcal{F}(\mathbf{K})}\mathcal{F}(\mathbf{K})+\lambda_{3}[\overline{\mathcal{F}(\bm{\partial_{x}})}\mathcal{F}(\bm{\partial_{x}})+\overline{\mathcal{F}(\bm{\partial_{y}})}\mathcal{F}(\bm{\partial_{y}})]+\lambda_{3}})$
(13)

where $\mathcal{F}(\cdot)$, $\mathcal{F}^{-1}(\cdot)$ are the FFT and inverse
FFT operators, respectively, and $\overline{\mathcal{F}(\cdot)}$ denote the
complex conjugate operator of FFT. $\bm{\partial_{x}}$ and $\bm{\partial_{y}}$
are the horizontal and vertical partial differential operators, respectively.
To achieve a better estimate of $\mathbf{I}$, we keep updating $\lambda_{3}$
at each iteration as has been done in [4],[17].

##  4 Comparative Experiments

In this section, we use satellite image data and compare the performance of
our proposed method with other state-of-the-art methods. We describe the
dataset and experimental settings used for the implementation of the methods
and discuss the comparison results.

###  4.1 Dataset

The dataset used in this study includes a simulated image of satellite that is
convolved with unknown kernels of various noises to generate a wide spectrum
of blurred noisy images. These synthetic images are distorted by a two-layer
wind model whose strength turbulence parameters create different types of
noisy images, characterized as $Dr=10$ (Dr10) and $Dr=20$ (Dr20) [24]. The
dataset consists of 200 images in total, including 100 images for each
distortion type of Dr10 or Dr20. The images are in RGB (three channels) format
and their patch size is 365$\times$365\. The authors of Swindle et al. [24]
generously provided the 200 blurred images and the original (simulated)
satellite image for the analysis of this study. As such, we are not aware of
what kernels were used for the convolution and how the final blurry images
were created.

Fig. 5(a) and 5(b) presents a sample of Dr10 and Dr20 image, respectively.
Dr20 images are noisier and more degraded, and it is even hard to distinguish
the object of satellite. Fig. 5(c) shows the simulated satellite image without
any degradation, which will be referred to as the real image. The different
images demonstrate the distortion severity of the dataset, especially for
Dr20.

(a) Blurred image at Dr10

(b) Blurred image at Dr20

(c) Simulated (real) Image

Figure 5: Sample images and the original image without degradation

###  4.2 Experimental Settings

Our proposed method is evaluated in comparison with other state-of-the-art
methods in both quantitative and qualitative fashion. The benchmark methods
include Xu and Jia [26], Ren et al. [20], Pan et al. [18], and Bai et al. [4].
These papers are not only recent but also have shown their superior
performance to others. In addition, their methods are publicly available
online as software or Matlab source codes, so the implementation of these
methods can be accurate and simple. Xu and Jia [26] can be implemented via
software, but for its implementation, a specific kernel size needs to be
supplied. There are three options of small, medium, and large kernels. The
medium kernel provided the best quality of deblurred images when applied Xu
and Jia [26] to our data. For a fair comparison, we consider the same size of
kernel for all the methods including ours, i.e., $h=31\times 31$ (the value
for the medium size kernel). Our method is implemented in Python on an i7-8700
CPU system.

To fully specify our model and estimation method, we need to determine a set
of parameters other than the model parameters. This includes the number of
base kernels, $N$, and the regularization parameters of
$\lambda_{1},\lambda_{2}$, and $\lambda_{3}$.

Fig. 3(c) and 4(c) implies that the number of kernels is one of the most
critical factors defining the overall shape of blur kernel. As its impact
being so crucial, we do not choose the exact number of kernels with any prior
knowledge, but we make the selection adaptive to the given image and
underlying blurriness structure therein. This is achieved by modeling the
covariance prior in Eq. (10) while leaving the value of $N$ as a sufficiently
large number. From our preliminary study, we found that $N=9$ could provide
enough level of flexibility for the kernel mixture. For the regularization
parameters, we sampled a few images and selected $\lambda_{1}=10^{-4}$ and
$\lambda_{2}=10^{-2}$ based on the quantitative metrics described in Section
4.3 and visual inspection of recovered images. On the other hand, as described
in Alg. 1, $\lambda_{3}$ is kept updated in each iteration of the estimation
procedure, which is initialized to $10^{-2}$.

One thing to note is, for the base kernels, we use $K_{G}^{\text{center}}$
instead of $K_{G}^{\text{rotation}}$ to form the kernel mixture, $\mathbf{K}$.
Although $K_{G}^{\text{rotation}}$ is the most flexible and capable of
modeling the most complicated shape of blurriness, it adds additional $N$
rotation parameters to estimate, one for each base kernel. Then, the total
number of parameters to estimate becomes $5N$ which is 45 when $N=9$. From the
perspective of the nonlinear optimization in Eq. (11), these additional
parameters require significantly more computations. At least for the images
used in this study, the additional parameters do not add much benefit in terms
of the image recovery (also see the similarity between Fig. 3(c) and 4(c)) but
increase the variance of the final estimate producing a poorer quality of a
recovered image.

###  4.3 Performance Measures

We use two performance measures to quantify and compare the quality of various
blind deblurring methods. The first measure is root mean square error (RMSE)
that is widely used to evaluate the estimation (or prediction) accuracy of a
model not only in the blind deblurring domain but also in general machine
learning applications. The second measure is peak signal-to-noise ratio (PSNR)
that is often used in computer vision and image processing applications to
quantify the quality of the reconstructed images.

The RMSE and PSNR can be calculated as follows:

$\displaystyle\text{RMSE}(\mathbf{I},\hat{\mathbf{I}})=\sqrt{\frac{1}{n_{h}n_{v}}\sum_{i=1}^{n_{h}}\sum_{j=1}^{n_{v}}(\mathbf{I}_{ij}-\hat{\mathbf{I}}_{ij})^{2}}$
(14)
$\displaystyle\text{PSNR}(\mathbf{I},\hat{\mathbf{I}})=\frac{\max_{i,j}\hat{\mathbf{I}}_{i,j}-\min_{i,j}\hat{\mathbf{I}}_{i,j}}{\sum_{i=1}^{n_{h}}\sum_{j=1}^{n_{v}}(\mathbf{I}_{ij}-\hat{\mathbf{I}}_{ij})^{2}/n_{h}n_{v}}$
(15)

where $\mathbf{I}_{ij}$ and $\hat{\mathbf{I}}_{ij}$ for $i=1,\ldots,n_{h}$ and
$j=1,\ldots,n_{v}$ are the values of $i$th horizontal and $j$th vertical pixel
from the real and recovered image, respectively. Since the size of the images
used in this study is $365\times 365$, $n_{h}=n_{v}=365$. The RMSE evaluates
how different the real and recovered images are whereas the PSNR quantifies
how much (peak) variation is captured by the recovered image relative to the
average of remaining variation. As such, the smaller the RMSE is and the
larger the PSNR is, the better the quality of a recovered image is. Because
they measure different aspects of the quality, we use both of them to compare
the proposed method with others.

###  4.4 Comparison Results

We randomly sample 10 images from each of the Dr10 and Dr20 datasets to
compare the performance.

####  4.4.1 Dr10 Results

Table 1 shows the RMSE and PSNR values for the ten recovered Dr10 images. In
all cases, the proposed method outperforms Xu and Jia [26] and Ren et al.
[20]. In addition, our method performs better in most cases than the most
sophisticated method, the graph-based blind deblurring [4]. Still, the graph-
based method consistently shows descent performance while producing the best
recovered image in one case. In couple other cases, the dark channel prior
method [18] performs the best, but its performance is worse than the graph-
based method on average. Even though the dark channel method can provide
really great quality for certain images, it can also suffer from considerably
poor performance; see the results of Image 3, 9, and 10 where its RMSE and
PSNR significantly deviate from the best values. On the contrary, our proposed
method not only has the most best cases but also remains close to the best
performance whenever it loses the first place. All these observations are
applied to both RMSE and PSNR measures. Overall, the proposed method has the
lowest RMSE and the highest PSNR on average. Fig. 6 visualizes the relative
performance of all methods, demonstrating the superiority of our method.

Table 1: Performance measures for randomly sampled Dr10 images; the boldface
highlights the best value for each image.

RMSE PSNR Image Xu et al. [26] Ren et al. [20] Pan et al. [18] Bai et al. [4]
Ours Xu et al. [26] Ren et al. [20] Pan et al. [18] Bai et al. [4] Ours 1
16.89 19.47 14.26 14.55 14.90 24.98 23.75 26.46 26.29 26.08 2 23.51 23.25
20.31 20.18 19.20 22.12 22.21 23.39 23.44 23.87 3 22.38 18.12 16.81 14.11
13.27 22.54 24.38 25.03 26.55 27.08 4 19.43 18.44 14.62 14.95 14.61 23.77
24.22 26.24 26.05 26.25 5 17.33 20.54 15.56 17.19 16.20 24.76 23.29 25.70
24.84 25.35 6 16.00 14.72 13.77 14.96 14.84 25.45 26.18 26.76 26.04 26.11 7
18.29 22.33 19.40 19.27 17.73 24.30 22.56 23.78 23.84 24.57 8 18.78 20.29
14.68 15.34 15.71 24.10 23.40 26.21 25.83 25.62 9 18.90 20.99 15.27 14.50
11.96 24.01 23.10 25.86 26.31 27.99 10 19.41 13.50 15.88 12.49 13.19 23.78
24.72 25.52 27.61 27.14 Mean 19.10 19.16 16.05 15.75 15.16 23.98 23.78 25.50
25.68 26.00

(a) RMSE results

(b) PSNR results

Figure 6: Visualization of relative performance - Dr10

(a)

(b)

(c)

(d)

(e)

(f)

Figure 7: Dr10 - Results output (Image 3): (a) Blurred images, (b) Xu and Jia
[26], (c) Ren et al. [20], (d) Pan et al. [18], (e) Bai et al. [4], (f)
Proposed method

(a)

(b)

(c)

(d)

(e)

(f)

Figure 8: Dr10 - Results output (Image 9): (a) Blurred images, (b) Xu and Jia
[26], (c) Ren et al. [20], (d) Pan et al. [18], (e) Bai et al. [4], (f)
Proposed method

Fig. 7 and 8 illustrate some deblurring outcomes from all methods as well as
the blurred image that has been used as an input (Image 3 and 9,
respectively). As expected, the visual outcomes are in accordance with the
results derived from the quantitative measures. By taking a closer look into
the deblurred images, we observe that the dark channel method (Fig. 7(d)) and
the graph-based method (Fig. 7(e)) create non-smooth edges of the object,
especially around the wings of the satellite. On the contrary, the recovered
image from our method (Fig. 7(f)) shows fine edges around the object. The
difference in the continuity (or smoothness) of nearby pixel values is
reflected in the RMSE and PSNR calculation where our method attains better
values than the other two. Fig. 7(b) and 7(c) show crude reconstruction of the
images without much details on the surface of the satellite, which causes
higher RMSE and lower PSNR for the corresponding methods.

Fig. 7(d) shows a higher contrast than Fig. 7(e) and 7(f) though its contrast
level is still lower than that of Fig. 7(b) and 7(c). This is likely because
the dark channel method regularizes the dark channel of the image producing
more of darker pixels. By having a relatively high contrast, the image in Fig.
7(d) looks clearer at distance. However, the high contrast, when it is not
appropriate, can result in poorer performance in the image recovery (see Table
1 for Image 3 and 9).

The similar argument can be made for Fig. 8. The recovered images from Xu and
Jia [26] and Ren et al. [20] lack detail. The dark channel method and the
graph-based method generate even more corrupted edges of the object. The image
recovered from the dark channel method has a higher contrast than the images
from the graph-based method and our method. All the results in terms of RMSE,
PSNR, and visual inspection demonstrate that the proposed method is
competitive with and superior to other state-of-the-art methods in recovering
the latent images.

Table 2: Performance measures for randomly sampled Dr20 images; the boldface
highlights the best value for each image.

RMSE PSNR Image Xu et al. [26] Ren et al. [20] Pan et al. [18] Bai et al. [4]
Ours Xu et al. [26] Ren et al. [20] Pan et al. [18] Bai et al. [4] Ours 1
23.42 23.80 21.67 22.64 21.08 22.15 22.01 22.82 22.44 23.06 2 24.80 26.35
26.17 24.23 23.34 21.65 21.13 21.18 21.85 22.18 3 24.02 25.49 23.45 20.57
20.37 21.93 21.42 22.14 23.28 23.36 4 23.28 26.19 24.18 22.78 21.18 22.20
21.18 21.87 22.39 23.02 5 21.06 24.32 22.56 22.14 18.84 23.07 21.82 22.47
22.64 24.04 6 23.85 21.93 23.69 23.29 21.26 21.99 22.72 22.05 22.20 22.99 7
23.42 25.05 24.56 24.09 22.61 22.15 21.57 21.74 21.90 22.45 8 22.60 21.80
19.97 21.24 19.90 22.46 22.77 23.53 23.00 23.56 9 19.75 22.59 20.28 19.77
19.77 23.63 22.46 23.40 23.62 23.62 10 20.47 21.92 19.39 17.96 18.17 23.32
22.73 23.79 24.45 24.35 Mean 22.67 23.94 22.59 21.87 20.65 22.45 21.98 22.49
22.77 23.26

(a) RMSE results

(b) PSNR results

Figure 9: Visualization of relative performance - Dr20

####  4.4.2 Dr20 Results

As discussed earlier, Dr20 images involve more blurriness, so it is natural to
observe worse performance of the methods compared to Dr10 results. Table 2
shows the RMSE and PSNR calculations for randomly chosen ten Dr20 images. For
these images with a higher level of noises, the results show the overwhelming
superiority of the proposed method. taking advantage of explicitly modeling
the structure of the blur kernel. Our method attains the lowest RMSE and the
highest PSNR in all but one cases. Even for the case where the graph-based
method achieves the best quality, our method is the only method producing a
competitive result to the best value. This superior performance is likely
because our method explicitly models the structure of the blur kernel that
represents blurriness.

##  5 Concluding Remarks

In this paper, we propose a novel blind image deblurring method that imposes a
specific structure on the blur kernel and achieves great flexibility in
modeling blurriness. To this end, we develop structure-enhanced Gaussian
kernels and form a mixture of the kernels to model the blur kernel. While the
behavior of the resulting blur kernel is regulated within a parametric
structure, it can represent various shapes of blurriness. Still, the modeling
capability and flexibility of the blur kernel depend on how many kernels to
incorporate in the kernel mixture. To address this issue, we reformulate the
optimization framework for kernel estimation and let the optimization process
itself decide the number of kernels through a covariance prior of the blur
kernel.

Our experimental results based on satellite image data show that the proposed
method outperforms descent state-of-the-art methods that employed some complex
image priors, in both quantitative and qualitative manners. Our method attains
the lowest RMSE and the highest PSNR, and this superiority becomes more
apparent when the noise level gets higher. From visual inspection, while other
methods suffer from discontinuity in nearby pixel values, the proposed method
recovers images with smooth edges without any corruption. Yet, the purpose of
this method regularizing the blur kernel is not in replacing image prior-based
methods but to advance general blind deconvolution practices through a proper
combination with the image prior-based methods.

In this paper, we combined several kernels in an additive way which was the
simplest form of integrating kernels. Applying a more advanced method for the
kernel fusion such as Gaussian Conditional Random Fields would be an
interesting research topic. In addition, by improving the optimization process
of kernel estimation, the proposed method can be applied to much broader
deblurring problems of significant importance. Therefore, finding better
heuristics and deriving a closed form solution based on proper assumptions
would be a valuable research task.

## Acknowledgment

The authors would like to thank Ryan Swindle, Douglas Hope, Michael Hart, and
Stuart Jefferies for providing the dataset used in this study.

## References

  * Almeida and Almeida [Aug. 2009] Mariana SC Almeida and Luís B Almeida.  Blind and semi-blind deblurring of natural images.  _IEEE Trans. Image Process._ , 19(1):36–52, Aug. 2009.
  * Babacan et al. [2012] S Derin Babacan, Rafael Molina, Minh N Do, and Aggelos K Katsaggelos.  Bayesian blind deconvolution with general sparse image priors.  In _Proc. Eur. Conf. Comput. Vis._ , pages 341–355. Springer, 2012\.
  * Babacan et al. [Nov. 2008] S Derin Babacan, Rafael Molina, and Aggelos K Katsaggelos.  Variational bayesian blind deconvolution using a total variation prior.  _IEEE Trans. Image Process._ , 18(1):12–26, Nov. 2008.
  * Bai et al. [Oct. 2018] Yuanchao Bai, Gene Cheung, Xianming Liu, and Wen Gao.  Graph-based blind image deblurring from a single photograph.  _IEEE Trans. Image Process._ , 28(3):1404–1418, Oct. 2018.
  * Bar et al. [Dec. 2006] Leah Bar, Nahum Kiryati, and Nir Sochen.  Image deblurring in the presence of impulsive noise.  _Int. J. Comput. Vis_ , 70(3):279–298, Dec. 2006\.
  * Chan and Wong [Mar. 1998] Tony F Chan and Chiu-Kwong Wong.  Total variation blind deconvolution.  _IEEE Trans. Image Process._ , 7(3):370–375, Mar. 1998.
  * Cho and Lee [2009] Sunghyun Cho and Seungyong Lee.  Fast motion deblurring.  _ACM Trans. Graph._ , 28(5):1–8, 2009.
  * Fergus et al. [2006] Rob Fergus, Barun Singh, Aaron Hertzmann, Sam T Roweis, and William T Freeman.  Removing camera shake from a single photograph.  _ACM Trans. Graph._ , 25(3):787–794, 2006.
  * Jiang et al. [Jul. 2003] Ming Jiang, Ge Wang, Margaret W Skinner, Jay T Rubinstein, and Michael W Vannier.  Blind deblurring of spiral ct images.  _IEEE Trans. Med. Imag._ , 22(7):837–845, Jul. 2003.
  * Lee and Hwang [May. 2017] Chia-Chen Lee and Wen-Liang Hwang.  Mixture of gaussian blur kernel representation for blind image restoration.  _IEEE Trans. Comput. Imag._ , 3(4):783–797, May. 2017.
  * Levin et al. [Jul. 2011] Anat Levin, Yair Weiss, Fredo Durand, and William T Freeman.  Understanding blind deconvolution algorithms.  _IEEE Trans. Pattern Anal. Mach. Intell._ , 33(12):2354–2367, Jul. 2011.
  * Li et al. [Sep. 2016] Chong-Yi Li, Ji-Chang Guo, Run-Min Cong, Yan-Wei Pang, and Bo Wang.  Underwater image enhancement by dehazing with minimum information loss and histogram distribution prior.  _IEEE Trans. Image Process._ , 25(12):5664–5677, Sep. 2016.
  * Likas and Galatsanos [Jul. 2004] Aristidis C Likas and Nikolas P Galatsanos.  A variational approach for bayesian blind image deconvolution.  _IEEE Trans. Signal Process._ , 52(8):2222–2233, Jul. 2004.
  * Mai and Liu [2015] Long Mai and Feng Liu.  Kernel fusion for better image deblurring.  In _Proc. IEEE Conf. Comput. Vis. Pattern Recognit._ , pages 371–380, June 2015.
  * Miskin and MacKay [2000] James Miskin and David JC MacKay.  Ensemble learning for blind image separation and deconvolution.  In _Advances in Independent Component Analysis_ , pages 123–141. Springer, 2000.
  * Molina et al. [Nov. 2006] Rafael Molina, Javier Mateos, and Aggelos K Katsaggelos.  Blind deconvolution using a variational approach to parameter, image, and blur estimation.  _IEEE Trans. Image Process._ , 15(12):3715–3727, Nov. 2006.
  * Pan et al. [Apr. 2016] Jinshan Pan, Zhe Hu, Zhixun Su, and Ming-Hsuan Yang.  $l_{0}$-regularized intensity and gradient prior for deblurring text images and beyond.  _IEEE Trans. Pattern Anal. Mach. Intell._ , 39(2):342–355, Apr. 2016.
  * Pan et al. [Jun. 2016] Jinshan Pan, Deqing Sun, Hanspeter Pfister, and Ming-Hsuan Yang.  Blind image deblurring using dark channel prior.  In _Proc. IEEE Conf. Comput. Vis. Pattern Recognit._ , pages 1628–1636, Jun. 2016.
  * Peng and Cosman [Feb. 2017] Yan-Tsung Peng and Pamela C Cosman.  Underwater image restoration based on image blurriness and light absorption.  _IEEE Trans. Image Process._ , 26(4):1579–1594, Feb. 2017.
  * Ren et al. [May. 2016] Wenqi Ren, Xiaochun Cao, Jinshan Pan, Xiaojie Guo, Wangmeng Zuo, and Ming-Hsuan Yang.  Image deblurring via enhanced low-rank prior.  _IEEE Trans. Image Process._ , 25(7):3426–3437, May. 2016.
  * Richardson [1972] William Hadley Richardson.  Bayesian-based iterative method of image restoration.  _JoSA_ , 62(1):55–59, 1972.
  * Shan et al. [2008] Qi Shan, Jiaya Jia, and Aseem Agarwala.  High-quality motion deblurring from a single image.  _ACM Trans. Graph._ , 27(3):1–10, 2008.
  * Shen et al. [May. 2012] Huanfeng Shen, Lijun Du, Liangpei Zhang, and Wei Gong.  A blind restoration method for remote sensing images.  _IEEE Geosci. Remote Sens. Lett._ , 9(6):1137–1141, May. 2012.
  * Swindle et al. [Sep. 2018] Ryan Swindle, Douglas Hope, Michael Hart, and Stuart Jefferies.  High-resolution space situational awareness imaging using carbon fiber telescopes.  _J. Appl. Remote Sens_ , 12(4):042406, Sep. 2018\.
  * Wiener [1964] Norbert Wiener.  _Extrapolation, Interpolation, and Smoothing of Stationary Time Series_.  The MIT Press, 1964.  ISBN 0262730057.
  * Xu and Jia [2010] Li Xu and Jiaya Jia.  Two-phase kernel estimation for robust motion deblurring.  In _Proc. Eur. Conf. Comput. Vis._ , pages 157–170. Springer, 2010\.
  * Yan and Shao [Feb. 2016] Ruomei Yan and Ling Shao.  Blind image blur estimation via deep learning.  _IEEE Trans. Image Process._ , 25(4):1910–1921, Feb. 2016.
  * You et al. [Jun. 2019] Chenyu You, Guang Li, Yi Zhang, Xiaoliu Zhang, Hongming Shan, Mengzhou Li, Shenghong Ju, Zhen Zhao, Zhuiyang Zhang, Wenxiang Cong, et al.  Ct super-resolution gan constrained by the identical, residual, and cycle learning ensemble (gan-circle).  _IEEE Trans. Med. Imag._ , 39(1):188–203, Jun. 2019.
  * You and Kaveh [Mar. 1996] Yu-Li You and Mostafa Kaveh.  A regularization approach to joint blur identification and image restoration.  _IEEE Trans. Image Process._ , 5(3):416–428, Mar. 1996.
  * Zhang et al. [Jan. 2019] Shuo Zhang, Guanghui He, Hai-Bao Chen, Naifeng Jing, and Qin Wang.  Scale adaptive proposal network for object detection in remote sensing images.  _IEEE Geosci. Remote Sens. Lett._ , 16(6):864–868, Jan. 2019.