1 Introduction

The image restoration problem consists in the reconstruction of the best estimate of an object "x'' from the knowledge of a blurred image $\tilde y$. In the general case, the transformation suffered by "x'' is described by a Fredholm integral equation of the first kind:

$$\tilde y ( r ) =\int h(r,\alpha)x(\alpha){\rm d}\alpha.$$ (1)

The function $h(r,\alpha)$ is the kernel of the integral equation. For the deconvolution problem, the kernel h is space invariant and is called the Point Spread Function (PSF). The above relation then becomes:

$$\tilde y ( r )=\int h(r-\alpha)x(\alpha){\rm d}\alpha.$$ (2)

Experimentally, the main difficulty in recovering x arises from the fact that we can only access a noise-corrupted version y of $\tilde y$, for which the convolution relationship does not apply. The general problem is then to obtain a solution x minimizing a distance D between the model and the observation, thus we can write the problem in the form:

$${\rm Minimize}/x:\ D \lfloor y,\int h(r-\alpha)x(\alpha){\rm d}\alpha\rfloor.$$ (3)

The particular definition of the distance D must be established in accordance with the statistics of the noise. This problem belongs to the general class of ill-posed problems in the sense of Hadamar, extensively studied in the literature (Bertero 1989; Bertero et al. 1980, 1995; Demoment 1986). The main difficulties lie in the instabilities appearing in the solution due to noise amplification during the reconstruction process. To avoid such effects and to obtain stable and physically suitable results, a regularization of the problem is necessary. A prior knowledge of the properties of the admissible solutions is required; however, classical constraints for a non-negative and spatially limited solution do not allow avoiding the instabilities. To obtain stable solutions, some regularization by a smoothness constraint is required. An explicit regularization implies the choice of the smoothing operator and the strength of the constraint, that is a compromise between the consistency with the data and the amount of details wanted in the reconstructed image. For iterative algorithms, an elementary regularization can be made by stopping the iterations. It is what we have done in the present paper.

The mathematical developments we present in this paper are written in the algebraic form corresponding to the discrete problem. The above relations then write in the form of a linear system $\tilde y=Hx$. The solution x is obtained by minimizing the distance D(y,Hx). The operator H takes a particular form in the case of the convolution we examine here. This assumption makes the numerical simulations easier, but most of the results obtained in this paper remain valid for the general case of Eq. (1).

The deconvolution problem can be solved using a likelihood (or Log likelihood) maximization (Lane 1996; Zaccheo & Gonsalves 1996). We analyze it for two different aspects of the noise process (Gaussian and Poisson processes), leading to the Image Space Reconstruction Algorithm (ISRA) and the Richardson-Lucy Algorithm (RLA).

The paper is organized as follows. In Sect. 2, after a brief recall on the classical unconstrained Gaussian case leading to the generalized inverse solution, we write ISRA in additive form, to make it appear as the solution of the least squares problem under non-negativity constraint. Then we proceed in the same way for RLA that corresponds to the Log-Likelihood maximization with non-negativity constraint for the Poisson case. The additive form makes the comparison between these two algorithms easier, and evidences interesting features related to variances of the noise. Next, we compare the behavior of the two algorithms on simulated data and analyze the spectral extension due to the non-linearity of the algorithms. In Sect. 3, we propose a new method to determine the iteration number corresponding to the best-restored image, by comparison with the results of the Wiener filter. This is a slight modification of a procedure presented in a recent ESO/NOAO Conference (Lanteri et al. 1998). Two appendices are used for particular computations and conclusions are given in Sect. 4.