Many applications of detecting point-like objects in astronomical images require, in a subsequent phase, to restore the background. Here we mean the ``flat" parts which will usually consist in an extended object and a true sky background. Applying the criteria discussed in the above section allows to generate a mask or several masks which are combined into a single, final one. The ideal case of point-like objects with Gaussian shapes would lead to a mask composed of perfect circular areas at the location of each ``star". Defects or saturations in the real images introduce distortions of those areas. They may even present a hole if the image is severely saturated and will appear as rings which may further be partly open. These problems can be alleviated by a simple morphological operation which consists in a slight geometrical extension, a ``dilatation", of the above areas. Usually, an extension by one or two pixels is adequate and allows to generate the final mask. It is therefore a binary image where the valid pixels of the background are set to 1 and the invalid pixels of the point-like objects to be removed are set to 0. The purpose of the next and final phase is to replace patches of invalid pixels in the original intensity image by a local interpolation, so-called pyramidal interpolation, using a non-linear multiresolution method which reminds of the general lines of the Haar wavelet analysis and synthesis.
This method is preferred to a global surface fitting procedure because it preserves the original intensity values of the background where it is valid, allowing both better local and global adjustments due to the multiresolution process. Incidentally, a strictly local interpolation is inapplicable because the remaining valid points may appear as a set of unconnected domains.
The procedure starts by shrinking the original 
pixels image to a
image. Each pixel in the latter corresponds to an area
of 
pixels in the former. Their value is a mean of valid pixels
in the related area. If 3 or 4 pixels in this area are invalid,
the mean is classified as invalid. After this first phase, we obtain an
image of means and a corresponding binary image or mask of valid or invalid
means. Repeating this procedure builds another couple of images
(means and mask) of lower resolution. The image of means is called the
analysis image.
The procedure is repeated l times until either the extinction of patches of invalid pixels is obtained, or Nl+1 or Ml+1 are too small (that is equal to or less than 4 pixels by side). If some final pixels are invalid, they are replaced by values estimated with a global adjustment of a quadric surface over the remaining valid pixels. After this, all pixels are valid at the lowest resolution level and the synthesis process can start.
The synthesis is done step by step as follows:
,
The procedure stops when the initial resolution is recovered. The method restores, in each patch of invalid pixels, the low resolution estimate from neighbouring valid pixels.
This procedure needs a minimal percentage of valid pixels (
to
obtain meaningful results. When the number of valid pixels is very low,
the threshold to obtain a significant mean from a area in an
analysis step can be reduced to one valid pixel.
Turning to practical considerations, the iterative nature of this algorithm
allows an easy implementation, in spite of its complex appearance.
Due to the pyramidal nature, it uses less dynamic memory (i.e. for an
initial 
pixels image, we need a temporary storage equal to
,
and its complexity is 
operations,
far better than spline 2-D interpolation methods. And last but not least,
it is very stable even with ill-conditionned data (for instance, very irregular
spatial distributions of valid data).
Figure 3: The logarithm of the intensity image of comet P/Halley
after specific corrections described in the text
Figure 4: The bidimentional histogram of the principal curvatures
(k1, k2 diagram)