James G. Nagy

Restoring images degraded by spatially-variant blur

Restoration of images that have been blurred by the effects of a Gaussian blurring function is an ill-posed but well-studied problem. Any blur that is spatially invariant can be expressed as a convolution kernel in an integral equation. Fast and effective algorithms then exist for determining the original image by preconditioned iterative methods. If the blurring function is spatially variant, however, then the problem is more difficult. In this work we develop fast algorithms for forming the convolution and for recovering the original image when the convolution functions are spatially variant but have a small domain of support. This assumption leads to a discrete problem involving a banded matrix. We devise an effective preconditioner and prove that the preconditioned matrix differs from the identity by a matrix of small rank plus a matrix of small norm. Numerical examples are given, related to the Hubble Space Telescope (HST) Wide-Field/Planetary Camera. The algorithms that we develop are applicable to other ill-posed integral equations as well.

Iterative methods for image deblurring: a Matlab object-oriented approach ∗

In iterative image restoration methods, implementation of efficient matrix vector multiplication, and linear system solves for preconditioners, can be a tedious and time consuming process. Different blurring functions and boundary conditions often require implementing different data structures and algorithms. A complex set of computational methods is needed, each likely having different input parameters and calling sequences. This paper describes a set of Matlab tools that hide these complicated implementation details. Combining the powerful scientific computing and graphics capabilities in Matlab, with the ability to do object-oriented programming and operator overloading, results in a set of classes that is easy to use, and easily extensible.

Restoration of atmospherically blurred images by symmetric indefinite conjugate gradient techniques

We consider an ill-posed deconvolution problem from astronomical imaging with a given noise-contaminated observation, and an approximately known convolution kernel. The limitations of the mathematical model and the shape of the kernel function motivate and legitimate a further approximation of the convolution operator by one that is self-adjoint. This simplifies the reconstruction problem substantially because the efficient conjugate gradient method can now be used for an iterative computation of a (regularized) approximation of the true unblurred image. Since the constructed self-adjoint operator fails to be positive definite, a symmetric indefinite conjugate gradient technique, called MR-II is used to avoid a breakdown of the iteration. We illustrate how the L-curve method can be used to stop the iterations, and suggest a preconditioner for further reducing the computations.

FFT-based preconditioners for Toeplitz-block least squares problems

Discretized two-dimensional deconvolution problems arising, e.g., in image restoration and seismic tomography, can be formulated as least squares computations,

Iterative image restoration using approximate inverse preconditioning

\min \| {b - Tx} \|_2

Kronecker product and SVD approximations in image restoration

, where T is often a large-scale rectangular Toeplitz-block matrix. The authors consider solving such block least squares problems by the preconditioned conjugate gradient algorithm using square nonsingular circulant-block and related preconditioners, constructed from the blocks of the rectangular matrix T. Preconditioning with such matrices allows efficient implementation using the one-dimensional or two-dimensional fast Fourier transform (FFT). Two-block preconditioners, related to those proposed by T. Chan and J. Olkin for square nonsingular Toeplitz-block systems, are derived and analyzed. It is shown that, for important classes of T, the singular values of the preconditioned matrix are clustered around one. This extends the authors’ earlier work on preconditioners for Toeplitz least squares iterations for one-dimensiona...

Quasi-Newton approach to nonnegative image restorations

Removing a linear shift-invariant blur from a signal or image can be accomplished by inverse or Wiener filtering, or by an iterative least-squares deblurring procedure. Because of the ill-posed characteristics of the deconvolution problem, in the presence of noise, filtering methods often yield poor results. On the other hand, iterative methods often suffer from slow convergence at high spatial frequencies. This paper concerns solving deconvolution problems for atmospherically blurred images by the preconditioned conjugate gradient algorithm, where a new approximate inverse preconditioner is used to increase the rate of convergence. Theoretical results are established to show that fast convergence can be expected, and test results are reported for a ground-based astronomical imaging problem.

Numerical methods for coupled super-resolution

Abstract Image restoration applications often result in ill-posed least squares problems involving large, structured matrices. One approach used extensively is to restore the image in the frequency domain, thus providing fast algorithms using fast Fourier transforms (FFTs). This is equivalent to using a circulant approximation to a given matrix. Iterative methods may also be used effectively by exploiting the structure of the matrix. While iterative schemes are more expensive than FFT-based methods, it has been demonstrated that they are capable of providing better restorations. As an alternative, we propose an approximate singular value decomposition (SVD), which can be used in a variety of applications. Used as a direct method, the computed restorations are comparable to iterative methods but are computationally less expensive. In addition, the approximate SVD may be used with the generalized cross validation method to choose regularization parameters. It is also demonstrated that the approximate SVD can be an effective preconditioner for iterative methods.

Optimal Kronecker Product Approximation of Block Toeplitz Matrices

Abstract Image restoration, or deblurring, is the process of attempting to correct for degradation in a recorded image. Typically the blurring system is assumed to be linear and spatially invariant, and fast Fourier transform (FFT) based schemes result in efficient computational image restoration methods. However, real images have properties that cannot always be handled by linear methods. In particular, an image consists of positive light intensities, and thus a nonnegativity constraint should be enforced. This constraint and other ways of incorporating a priori information have been suggested in various applications, and can lead to substantial improvements in the reconstructions. Nevertheless, such constraints are rarely implemented because they lead to nonlinear problems which require demanding computations. We suggest efficient implementations for three nonnegatively constrained restorations schemes: constrained least squares, maximum likelihood and maximum entropy. We show that with a certain parameterization, and using a Quasi-Newton scheme, these methods are very similar. In addition, our formulation reveals a connection between our approach for maximum likelihood and the expectation–maximization (EM) method used extensively by astronomers. Numerical experiments illustrate that our approach is superior to EM both in terms of accuracy and efficiency.

Covariance-Preconditioned Iterative Methods for Nonnegatively Constrained Astronomical Imaging

The process of combining, via mathematical software tools, a set of low-resolution images into a single high-resolution image is often referred to as super-resolution. Algorithms for super-resolution involve two key steps: registration and reconstruction. Most approaches proposed in the literature decouple these steps, solving each independently. This can be effective if there are very simple, linear displacements between the low-resolution images. However, for more complex, nonlinear, nonuniform transformations, estimating the displacements can be very difficult, leading to severe inaccuracies in the reconstructed high-resolution image. This paper presents a mathematical framework and optimization algorithms that can be used to jointly estimate these quantities. Efficient implementation details are considered, and numerical experiments are provided to illustrate the effectiveness of our approach.