WIPE| is a Fourier synthesis method recently developed in radio imaging and optical interferometry (Lannes et al. 1994, 1996). The name of WIPE is associated with that of ÇLEAN|, a deconvolution method intensively used in astronomy (Högbom 1974; Schwarz 1978). (WIPE| can also be used as a deconvolution technique.) As a matter of fact, WIPE| was devised, quite independently, on the grounds of well known properties in harmonic analysis and band-limited interpolation (Lannes et al. 1987a). There also exists a probabilistic interpretation of this method (Maréchal & Lannes 1997). Nevertheless, to some extent, WIPE can equally well be regarded as an updated version of ÇLEAN|. This paper is aimed precisely at clarifying the relationship between these two methods.
In Sect. 2,
we first show
that CLEAN lies in the family
of "matching pursuit" techniques
(Mallat & Zhang 1993).
In other words, the principle of the "iterative
harmonic analysis'' of ÇLEAN|, as it is
exhibited in Schwarz (1978),
can be regarded as a special case of a more general approach.
We then present the reconstruction methods that can be implemented
for solving the related optimization problems,
first, without any control of the propagation of
errors (Sect. 3),
and afterwards, with a complete control of the stability
parameters (Sect. 4).
The first weak point of CLEAN
can thus be exhibited:
in situations of astrophysical interest (for example,
when observing extended complex objects with dilute apertures),
the matching process of CLEAN is .
As shown in Sect. 5,
the regularization principle of WIPE|,
which can be applied to CLEAN as it is,
remedies this lack of robustness.
Unfortunately,
CLEAN has another weak point:
the basis functions used in the matching pursuit process,
the `clean beams,'
are not well suited for reconstructing
the boundaries of the structuring entities of the images
at the selected resolution level.
As a result, CLEAN must be interrupted before the best possible
fit is reached.
As indicated in Sect. 5,
the basic version of WIPE overcomes this difficulty
in a simple and efficient manner.
The multiresolution extension of WIPE|,
underlying the main aspects
of this paper,
reinforces the arguments justifying
this methodological choice.
The guiding idea of our analysis is based on the theoretical framework presented in Sect. 1.1. (The reader who is not familiar with the related basic notions is invited to consult Appendices 1 to 4 of Lannes et al. 1987a). As CLEAN was first used as a Fourier synthesis technique, it was natural to compare the principles of CLEAN and WIPE| in this context. Such a comparative study can only be started from a common statement of the problem. The corresponding formulation, which makes the theoretical framework more attractive, is presented in Sect. 1.2. We then also specify the general conditions of the numerical simulations in support of our analysis.
In many inverse problems,
the "reconstructed object" is described in
an "object representation space" E
generated by m vectors
selected among a family
of M vectors
.
The latter may be regarded as the "atoms" of the object representation
under consideration.
The linear space generated by all the vectors of
,
,
is a subspace of some
space:
the "object space"
.
The "data vector"
lies in another
space: the "data space"
.
(When the data are complex quantities,
it is always possible to work in the real linear space
underlying the complex linear space directly involved in the analysis.)
To a first approximation,
is related to the object (or image)
to be reconstructed
via a linear operator A from
into
.
The basic problem is to choose the vectors
in
,
and thereby to construct the object representation space E .
The solutions are then obtained
by minimizing on E the quadratic functional:
Clearly
is the norm on
;
the scalar products
on
and
, as well as the norms,
will be distinguished
by the subscripts o and d, respectively.
Let F be the "image" of E by A
(the space of the 's,
spanning E ),
be the operator from E onto F induced by A ,
and
be the projection of
onto F
(see Fig. 1 (click here)).
The vectors
minimizing q on E ,
the solutions of the problem,
are such that
.
They are identical up to a vector lying in the kernel of
.
Note that
.
(By definition,
the kernel of A ,
,
also referred to as the null space of A ,
is the set of vectors
such that
).
The condition
is a necessary condition for
to be reduced to
.
As
is orthogonal to F ,
the solutions
of the problem are characterized by the property:
On denoting by the adjoint of A ,
this property can also be written in the form
where r is regarded as a residue.
The solutions
are therefore characterized by the fact that the corresponding residue
is orthogonal to E , i.e.,
This condition is of course equivalent to
,
where
is the projection (operator)
of
onto E .
The solutions of the problem are therefore
the solutions of the equation on E :
For any and any
,
we have
hence
.
This explicitly shows that Eq. (3) is the "normal equation":
When the problem is well posed,
is a one-to-one map (
);
the solution is then unique:
there exists only one vector
such that
.
This vector,
, is said to be the least-squares solution
of the equation
.
In this case,
let be a variation of
in F , and
be the corresponding variation
of
(see Fig. 1 (click here)).
As shown in Appendix 1,
the robustness of the reconstruction
process is then governed by the inequality
where is the condition number of
:
and
respectively denote the greatest lower bound and
the least upper bound of
for the
's with norm unity in E:
The closer to 1 is the condition number, the easier and the more robust
is the reconstruction process.
As
and
are the smallest and largest eigenvalues
of
, respectively.
Figure 1: Uniqueness of the solution and robustness
of the reconstruction process.
Operator A is an operator from the object space
into the data space
.
The object representation space E
is a particular subspace of
(see text).
The image of E by A , the range of
,
is denoted by F .
In this Euclidean representation,
is the projection
of the data vector
onto F .
The inverse problem must be stated so that
is
a one-to-one map from E onto F ,
the condition number
having a reasonable value (see Eq. (5))
In the analysis of the different methods that can be devised and implemented for solving the problem, three main aspects must be examined:
1) the precise definition of , A and E ;
2) the selected minimization technique;
3) the robustness of the reconstruction process.
The first point is essential for the interpretation of the results.
It is related to the statement of the problem
and, if need be, to its "regularization;"
is not necessarily the vector of the experimental data.
In such a case, the definition of A must of course take into account
the corresponding manipulations.
As regards E ,
which may also be involved in the regularization of the problem,
it is important to note that this space
may be constructed, in a global manner or step by step, interactively or
automatically.
As suggested by the statement of point 2,
many different techniques can be used for minimizing q on E ;
some of these are certainly more efficient than others, but
in this case, the choice is not crucial.
The last point concerns the study of the propagation of errors.
The part played by inequality (5) in the development of the
corresponding error analysis
shows that a good reconstruction procedure must also provide, in particular,
the condition number
.
In the problems of Fourier synthesis
encountered in astronomy,
the "object function" of interest, ,
is a real-valued function of
an angular position variable
.
The "object model variable"
lies in some
object space
whose vectors, the functions
,
are defined at a very high level of resolution.
In this sense, one may say that
emulates
the Hilbert space of square integrable
real-valued functions
.
More precisely,
is characterized by two key parameters:
the extension
of its field,
and its resolution scale
;
is therefore much smaller than
the resolution limit that can be reasonably expected
at the end of the reconstruction process:
,
where N is some power of 2.
(The larger is N , the more oversampled is the object field.)
To define the object space more explicitly,
we introduce the finite grid (see Fig. 2 (click here)):
Let be an integer vector lying
in
:
.
Clearly, the functions
where
form an orthogonal set
with
In our comparative analysis of CLEAN
and WIPE|,
is the
space generated by
the basis vectors
,
spanning G
(see Fig. 2 (click here)).
The functions
lying in this space
can therefore be expanded in the form:
In the general context of this paper,
the functions ,
which play the role of interpolation or scaling functions,
can be regarded as the "elementary particles"
of the object representation.
Evidently,
other orthogonal sets of
scaling functions can be used
(cf. Sect. 4.3 of Lannes et al. 1994).
Let us finally note that
is
a subspace of
;
its inner product can be explicitly expressed in the form
(cf. Eqs. (8)-(11)):
Thus,
Figure 2: Object grid
for N = 8
The Fourier transform of
is defined by the relationship
where is a two-dimensional
angular spatial frequency:
.
According to Eqs. (11), (8) and (9),
we therefore have:
where
with
and
The "experimental data" are blurred
values of
on a finite list of frequencies in the Fourier domain:
.
As is a
function,
it is natural to define
as the complex conjugate of
.
The `experimental frequency list'
is defined consequently:
if
, then
(except for
:
in the convention adopted in the sequel,
either the null frequency does not lie in
,
or there exists only one occurrence of this point).
The experimental frequency coverage
generated by
is therefore centrosymmetric
(see Fig. 3 (click here)a).
The experimental data vector
lies in the "experimental data space"
.
By definition,
is the real Euclidean space underlying the space
of
functions on
such that
.
This space, whose dimension is equal to
,
is endowed with the scalar product
where , and
W is a given weighting function
(see for instance Eqs. (17) and (18) further on).
Note that this scalar product is real:
.
In ÇLEAN|,
the data space
reduces to
,
whereas in WIPE|,
is a natural extension of
(cf. Sect. 5.1).
Let be the Fourier domain:
.
In Fourier synthesis, the frequency coverage to be synthesized
is a centrosymmetric region
.
In Fig. 3 (click here)a, this region is the disc centred on the origin.
CLEAN and WIPE share a common objective, that of
the image to be reconstructed.
This image,
,
is defined so that its Fourier transform
is quadratically negligible outside
.
More explicitly,
is defined
by a convolution relation of the form:
The "synthetic beam" s is a function
resulting from the choice of :
the "clean beam" in ÇLEAN|, the "neat beam" in WIPE
(see Fig. 3 (click here)b).
The Fourier data corresponding to
are then defined by the relationship:
Clearly,
lies in
.
The neat beam can be regarded as a sort of
optimal clean beam:
the optimal apodized point spread function
that can be designed
within the limits of the principle.
(Apodization theory is concerned with the
control of the
distribution of light over the exit pupil of an optical system
in order to achieve a suppression of side lobes of the
diffraction pattern
-- Slepian 1965).
More precisely,
the neat beam s is a
centrosymmetric function lying in the object space
,
and satisfying the following properties:
where
.
Figure 3: General conditions of the simulation;
a) experimental frequency coverage and frequency coverage
to be synthesized;
b) representation of the corresponding neat beam (for
).
The experimental frequency list corresponding to the
experimental frequency coverage shown in a)
includes 211 frequency points
(
).
The Fourier grid is of the form
where
with
.
The neat beam s represented in b)
corresponds to the frequency coverage to be synthesized
.
It is centred in the object grid
where
.
The resolution limit of the reconstruction process
(the full width of s at half maximum)
is of the order of
where
is the diameter of
In the terminology adopted in this paper,
an atom such as s
is a finite linear combination of elementary
particles .
The integers
involved in this linear combination
lie in a subset
of G
(see Figs. 2 (click here) and 3 (click here)b).
As explicitly shown in Sect. 2
of Lannes et al. (1996),
the computation of the corresponding coefficients does not raise
any particular difficulty.
In a first approach to Fourier synthesis, Eq. (13) suggests that the basis functions of the object representation space E should be translated versions of s : a finite number of atoms centred on the nodes of grid G . This very natural idea, which is exploited as it is, in the matching pursuit principle of CLEAN (cf. Sect. 2.2), may be completely relaxed in WIPE|. More precisely, the matching pursuit process may be performed at the level of the elementary particles. The resolution limit of the reconstruction process is then kept under control thanks to an appropriate regularization principle (cf. Sect. 5.1).
The simulations presented in this paper
correspond to the general conditions defined in Fig. 3 (click here).
The "object function"
was defined as a set of "
functions"
centred on the nodes of grid
.
The image to be reconstructed (
)
then lies in
.
Figure 4 (click here) gives an idea of the complexity of this image.
Figure 4: Image to be reconstructed, ,
at the resolution level defined in Fig. 3 (click here).
Note that the portion of the field represented here
is twice as large as that defined in Fig. 3 (click here)b
The Fourier data were blurred by adding a Gaussian noise.
More precisely, for all the
(see Fig. 3 (click here)a),
the standard deviation of
was set equal to 5% of
the total flux of the object:
.
The weighting function W introduced in Eq. (12) was
explicitly defined by the formula
where