Method

The method adopted might be described as image restoration incorporating designated sources having specified shape and known positions. For astronomical images, the designated sources are naturally taken to be point sources (stars, QSOs, AGNs), so the specified shape is the delta function. Clearly, an observer can commonly, with high confidence, identify many objects in an image as point sources. Accordingly, it makes sense to develop codes that use this extra information and, at least for these designated point sources, thereby eliminate the ringing associated with failed attempts to recover delta functions.

In 1-D vector notation, the adopted model for , the intensity distribution on the celestial sphere, is

Here represents distributed emission and the designated point sources. This latter representation is achieved by defining to be zero everywhere except at the locations () of point sources.

This model of the object can be mapped onto the image plane using the equation of image formation, and the resulting prediction for the intensity distribution in the image plane can then be compared to the observed distribution . The posed image restoration problem is thus to determine and in order to achieve an optimum fit of to .

An obvious first thought is to adopt the R-L procedure and determine and with an iterative scheme that asymptotically yields Maximum Likelihood estimates subject to the appropriate non-negativity and normalization constraints. But this results in indeterminacy since an individual star can be represented either in the vector or by the corresponding single-pixel peak added to the vector, or indeed by any linear combination of these extremes.

To eliminate this indeterminacy, we must impose a constraint on that excludes single-pixel peaks. In other words, we must impose a resolution limit on that is somewhat larger than the limit implied by the numerical discretization. This is achieved by introducing regularization into the optimization problem for determining and . Of the many possible regularization terms (see, e.g., Titterington 1985), we choose the entropic form

where and where the floating default (Horne 1985, Lucy 1994)

Here the resolution kernel is a bell-shaped function whose width is the required resolution limit. Structure in on a smaller scale results in a decrease in the entropy and so is disfavored.