Next: Results
Up: Adaptive Regularized Restoration Algorithms
Previous: Introduction
For the algorithms discussed here, the restoration approaches have
been formulated to handle the problem of an image degraded by a
spatially invariant or spatially variant blur operator with additive
Gaussian noise. This description represents the degradation present
in the HST WFPC images. The image
degradation process can be modeled according to
where a lexicographical ordering of the observed digital image, 
,
the original image, 
, and the additive noise, 
, is used. 
represents the degradation operator of the imaging system. The image
restoration problem calls for obtaining an estimate of given ,
, and some characterization of the noise process, .
Iterative approaches may be used to solve this inverse problem very
effectively. The primary advantages of iterative techniques are
(Schafer et al. 1981, Katsaggelos 1989): (i)
there is no need to explicitly implement the inverse of an operator;
(ii) knowledge about the solution may be directly incorporated into
the restoration process; (iii) the process may be monitored as it
progresses; (iv) the effect of noise may be controlled with certain
constraints; and (v) parameters determining the solution can be
updated as the iteration progresses.
The iterative techniques applied here to HST data are developed
through a set theoretic approach. In this approach, prior constraints
on the solution are imposed by (Katsaggelos et al. 1985, Katsaggelos 1989,
Katsaggelos et al. 1991)
and
where the solution belongs to both ellipsoids described by
Eqs. (2) and (3). In Eq. (2), the operator 
represents a highpass filter which bounds the high frequency energy of
the restored image. If the bounds 
and 
are
known, and the intersection of 
and 
is not empty, the
solution may be found by solving
where 
, the regularization parameter, is equal to
, which controls the trade-off between
fidelity to the data, and smoothness of the solution.
The first algorithm we tested for HST data was the generalized
iterative adaptive algorithm which simultaneously restores the image
and determines a single regularization parameter based on the restored
image at each iteration. This algorithm does not depend on the
initial conditions. In other words, the smoothing functional to be
minimized can be shown to have a unique minimizer. An accurate
estimate of the noise variance is also obtained at convergence. The
approach we have proposed (Kang and Katsaggelos 1993) extracts the
properties of the original image from the partially restored image at
each iteration step.
In order to solve for an optimal restored image and
regularization parameter at each
iteration, we use
as the functional to be minimized. The weighting matrices 
and
are included to make the restoration algorithm spatially
adaptive. The necessary condition for a minimum is that the gradient
of 
with respect to be equal to zero. This results
in
When the noise and the high pass filtered image are stationary,
the weighting matrices and become
symmetric, and thus the equation can be rewritten as
and it becomes
since 
with a proper choice of 
.
When the noise and the high pass filtered image are uncorrelated with
themselves, even though they are nonstationary, and
become diagonal, and therefore we obtain the solution for Eq.
(7). When the noise and the high pass filtered image are
white, then the weighting matrices become constant identity matrices
multiplied by the constant variances. The highly nonlinear term
can be removed by the proper
choice of the regularization functional, based on the global convexity
of the smoothing functional. Since
Eq. (8) is nonlinear, we can not solve
for in a direct way, but we can use an iterative
technique.
The regularization parameter is defined here as a function of the
original image (but in practice becomes a function of an estimate of
the original image). The form of the smoothing functional to be
minimized is of great importance since it preserves convexity and
exhibits only a global minimizer. After investigating the desirable
properties for the regularization functional to satisfy, the following
two forms have been shown to provide optimal solutions (Kang and
Katsaggelos 1993)
where 
,
and
where 
controls convergence and convexity.
Given an optimal choice for 
, we can solve Eq. (8)
by the method of successive approximations with (Schafer et al.
1981, Katsaggelos et al. 1985, Katsaggelos 1989, Katsaggelos et al.
1991):
Using either choice of regularization functional (Eq. (9) or
Eq. (10)), the iterative algorithm does not depend on the
initial condition despite the nonlinearity of the iteration. This is
due to the convexity of the functional and the convergence criteria
satisfied by the globally optimal iteration. The algorithm
is applicable to any type of degradation and stabilizing matrix
(both of which may be spatially varying). So, there is no
requirement that and be block circulant matrices.
It is also important to note that no knowledge of the noise variance, or
of the bound which determines the ellipsoid that expresses the
smoothness of the image, is assumed.
The second iterative algorithm examined here is a nonlinear frequency
domain algorithm in which the regularization parameter is frequency
dependent, and updated at each iteration step. In this case, the algorithm
considers that is a block circulant matrix representing a
spatially invariant blur. Also, the set-theoretic formulation is
constructed in a weighted space, such that
and
where 
and 
are both block circulant weighting matrices. These
matrices are chosen to maximize the speed of convergence at every
frequency component, and to compensate for the near-singular frequency
components of the iteration. The solution which belongs to the
intersection of the ellipsoids given by Eqs. (12) and (13)
is given by
We define 
, 
and 
. The block circulant
matrix 
is the spatial domain representation of the relaxation
parameter which will be called 
in the frequency domain,
and 
is the block circulant spatial domain representation of the
regularization parameter which is shown next in the frequency domain as
. Again, successive approximation may be used to
solve Eq. (14). Since all of the matrices in this equation
are block-circulant, the iteration may be written in the discrete
frequency domain as
where 
represents a single 2-D frequency component. In this
case,
Here, the 
term we use is defined by
where 
, and
and
The optimized parameter 
is given at each frequency component
by (Strand 1974)
According to this iteration, since and 
are frequency
dependent the convergence of the iteration can be accelerated, making
this an attractive algorithm where speed is a concern.
In the course of testing our algorithms, we have made appropriate
consideration of the optimization of these algorithms for HST data.
In particular, we have included a positivity constraint at each step
of the iteration. This imposes the condition that there should be
no negative flux in our image source. The algorithms were developed
with additive Gaussian noise as the observation noise in the model.
So, these algorithms are well equipped to deal with the read-out noise
problem of the WFPC.
Next: Results
Up: Adaptive Regularized Restoration Algorithms
Previous: Introduction