Patrizia Boccacci

Introduction to inverse problems in imaging

Part 1 Image deconvolution: some mathematical tools examples of image blurring the ill-posedness of image deconvolution regularization methods iterative regularization methods statistical methods. Part 2 Linear inverse problems: examples of linear inverse problems singular value decomposition (SVD) inversion methods revisited Fourier based methods for specific problems comments and concluding remarks. Part 3 Mathematical appendices. References.

Image deblurring with Poisson data: From cells to galaxies.

A semiconductor device including an N-type semiconductor substrate which includes arsenic as an impurity, a first electrode formed on a main surface of the N-type semiconductor substrate, a ground surface formed on another surface of the N-type semiconductor substrate, a second electrode formed on the ground surface and ohmically-contacted with the N-type semiconductor substrate, a semiconductor element formed in the N-type semiconductor substrate and flowing current between the first electrode and the second electrode during ON-state thereof. The device has a reduced ON-resistance thereof.

Super-resolution in computational imaging.

Super-resolution is a word used in different contexts but mainly in the case of methods aimed at improving the resolution of an optical instrument beyond the diffraction limit. Such a result may be achieved by means of specific instrumental techniques (such as, for instance, interferometry) or by means of a suitable processing of a digital image; in the latter case we will use the expression computational super-resolution (CSR). In this paper we describe the basic concepts underlying CSR without using the mathematics required for establishing its theoretical validity. The aim is to introduce a wide audience to this topic, to specify the situations where CSR is feasible and to emphasize the point that unlimited CSR is not possible.

Efficient gradient projection methods for edge-preserving removal of Poisson noise

Several methods based on different image models have been proposed and developed for image denoising. Some of them, such as total variation (TV) and wavelet thresholding, are based on the assumption of additive Gaussian noise. Recently the TV approach has been extended to the case of Poisson noise, a model describing the effect of photon counting in applications such as emission tomography, microscopy and astronomy. For the removal of this kind of noise we consider an approach based on a constrained optimization problem, with an objective function describing TV and other edge-preserving regularizations of the Kullback–Leibler divergence. We introduce a new discrepancy principle for the choice of the regularization parameter, which is justified by the statistical properties of the Poisson noise. For solving the optimization problem we propose a particular form of a general scaled gradient projection (SGP) method, recently introduced for image deblurring. We derive the form of the scaling from a decomposition of the gradient of the regularization functional into a positive and a negative part. The beneficial effect of the scaling is proved by means of numerical simulations, showing that the performance of the proposed form of SGP is superior to that of the most efficient gradient projection methods. An extended numerical analysis of the dependence of the solution on the regularization parameter is also performed to test the effectiveness of the proposed discrepancy principle.

On the Recovery and Resolution of Exponential Relaxation Rates from Experimental Data: A Singular-Value Analysis of the Laplace Transform Inversion in the Presence of Noise

The problem of numerical inversion of the Laplace transform is considered when the inverse function is of bounded, strictly positive support. The recent eigenvalue analysis of McWhirter and Pike for infinite support has been generalized by numerical calculations of singular values. A priori knowledge of the support is shown to lead to increased resolution in the inversion, and the number of exponentials that can be recovered in given levels of noise is calculated.

Image restoration in chirp-pulse microwave CT (CP-MCT)

Chirp-pulse microwave computed tomography (CP-MCT) is a technique for imaging the distribution of temperature variations inside biological tissues. Even if resolution and contrast are adequate to this purpose, a further improvement of image quality is desirable. In this paper, we discuss the blur of CP-MCT images and we propose a method for estimating the corresponding point spread function (PSF). To this purpose we use both a measured and a computed projection of a cylindrical phantom. We find a good agreement between the two cases. Finally the estimated PSF is used for deconvolving data corresponding to various kinds of cylindrical phantoms. We use an iterative nonlinear deconvolution method which assures nonnegative solutions and we demonstrate the improvement of image quality which can be obtained in such a way.

A discrepancy principle for Poisson data

In applications of imaging science, such as emission tomography, fluorescence microscopy and optical/infrared astronomy, image intensity is measured via the counting of incident particles (photons, γ-rays, etc). Fluctuations in the emission-counting process can be described by modeling the data as realizations of Poisson random variables (Poisson data). A maximum-likelihood approach for image reconstruction from Poisson data was proposed in the mid-1980s. Since the consequent maximization problem is, in general, ill-conditioned, various kinds of regularizations were introduced in the framework of the so-called Bayesian paradigm. A modification of the well-known Tikhonov regularization strategy results in the data-fidelity function being a generalized Kullback–Leibler divergence. Then a relevant issue is to find rules for selecting a proper value of the regularization parameter. In this paper we propose a criterion, nicknamed discrepancy principle for Poisson data, that applies to both denoising and deblurring problems and fits quite naturally the statistical properties of the data. The main purpose of the paper is to establish conditions, on the data and the imaging matrix, ensuring that the proposed criterion does actually provide a unique value of the regularization parameter for various classes of regularization functions. A few numerical experiments are performed to demonstrate its effectiveness. More extensive numerical analysis and comparison with other proposed criteria will be the object of future work.

Application of the projected Landweber method to the estimation of the source time function in seismology

The empirical Green function (EGF) model, which is used in this paper for the analysis of the waveforms of low-energy earthquakes, consists in assuming that the propagating medium and the recording instrument can be treated as a linear system and that the impulse response function of the system can be approximated by the waveform of a very small earthquake. The deconvolution of the Green function event from the waveform of a larger one, located at approximately the same position, provides information about the source time function (STF) of the latter. Linear inversion methods do not yield satisfactory estimations of the STF which must be positive and causal. Moreover, an estimate of the duration (support) of the STF should be desirable. In this paper we apply to this problem the so-called projected Landweber method, which is an iterative nonlinear method allowing for the introduction of constraints on the solution. The implementation of the method is easy and efficient. We first validate the method by means of synthetic data, generated by the use of waveforms of a seismic swarm that occurred in the Ligurian Alps (north-western Italy) during July 1993. Then, taking into account the indications provided by the simulations, the method has been applied to the inversion of real data, yielding satisfactory results also in the case of quite complex events.

Resolution in Diffraction-limited Imaging, a Singular Value Analysis

In a previous paper, methods of singular function expansions have been applied to the analysis of coherent imaging when the object and image domains are allowed to differ. In this paper the method is extended to incoherent illumination, restricting the analysis to the aberration-free case. While singular functions and singular values for coherent imaging are related in a simple way to the prolate spheroidal functions and their eigenvalues, such relations do not exist for the incoherent imaging case. In spite of this difficulty many properties of singular functions and singular values are derived in this paper and asymptotic estimates are obtained in the limit of large space-bandwidth product. For small values of the space-bandwidth product, the singular values are computed numerically and by means of these results it is shown that super-resolution, in the sense of improving on previous criteria in the presence of noise, can be achieved.

Efficient deconvolution methods for astronomical imaging: algorithms and IDL-GPU codes

Context. The Richardson-Lucy method is the most popular deconvolution method in astronomy because it preserves the number of counts and the non-negativity of the original object. Regularization is, in general, obtained by an early stopping of Richardson-Lucy iterations. In the case of point-wise objects such as binaries or open star clusters, iterations can be pushed to convergence. However, it is well-known that Richardson-Lucy is an inefficient method. In most cases and, in particular, for low noise levels, acceptable solutions are obtained at the cost of hundreds or thousands of iterations, thus several approaches to accelerating Richardson-Lucy have been proposed. They are mainly based on Richardson-Lucy being a scaled gradient method for the minimization of the Kullback-Leibler divergence, or Csiszar I-divergence, which represents the data-fidelity function in the case of Poisson noise. In this framework, a line search along the descent direction is considered for reducing the number of iterations. Aims. A general optimization method, referred to as the scaled gradient projection method, has been proposed for the constrained minimization of continuously differentiable convex functions. It is applicable to the non-negative minimization of the Kullback-Leibler divergence. If the scaling suggested by Richardson-Lucy is used in this method, then it provides a considerable increase in the efficiency of Richardson-Lucy. Therefore the aim of this paper is to apply the scaled gradient projection method to a number of imaging problems in astronomy such as single image deconvolution, multiple image deconvolution, and boundary effect correction. Methods. Deconvolution methods are proposed by applying the scaled gradient projection method to the minimization of the Kullback-Leibler divergence for the imaging problems mentioned above and the corresponding algorithms are derived and implemented in interactive data language. For all the algorithms, several stopping rules are introduced, including one based on a recently proposed discrepancy principle for Poisson data. To attempt to achieve a further increase in efficiency, we also consider an implementation on graphic processing units. Results. The proposed algorithms are tested on simulated images. The acceleration of scaled gradient projection methods achieved with respect to the corresponding Richardson-Lucy methods strongly depends on both the problem and the specific object to be reconstructed, and in our simulations the improvement achieved ranges from about a factor of 4 to more than 30. Moreover, significant accelerations of up to two orders of magnitude have been observed between the serial and parallel implementations of the algorithms. The codes are available upon request.