Äîêóìåíò âçÿò èç êýøà ïîèñêîâîé ìàøèíû. Àäðåñ îðèãèíàëüíîãî äîêóìåíòà : http://www.stsci.edu/stsci/meetings/irw/proceedings/schulzt.ps.gz Äàòà èçìåíåíèÿ: Mon May 23 18:44:30 1994 Äàòà èíäåêñèðîâàíèÿ: Sun Dec 23 05:19:46 2007 Êîäèðîâêà: IBM-866 Ïîèñêîâûå ñëîâà: ï ï ï ï ï ï ï ï ï ï ï ï ï ï ï ï

The Restoration of HST Images and Spectra II
Space Telescope Science Institute, 1994
R. J. Hanisch and R. L. White, eds.
Simultaneous Phase Retrieval and Deblurring for the Hubble Space Telescope
Timothy J. Schulz and Stephen C. Cain
Michigan Technological University, Houghton, MI 49931‘1295
Abstract. Most methods proposed for restoring images acquired by the Hubble Space
Telescope rely on prior knowledge of the telescope's pointíspread function; however, for
many images, this function is not known precisely and must be inferred from the noisy
measured data. In this paper, we address this problem and discuss a maximumílikelihood
estimation technique for simultaneously determining the nature of the aberrations and for
recovering the underlying object from a noisy, degraded image.
1. Introduction
Various methods have been proposed to restore images acquired by the Hubble Space Telescope
(HST). Many of these techniques rely on prior knowledge of the telescope's pointíspread function
(PSF). However, as noted by R. White (White 1993),
Most HST deconvolutions are carried out using observed PSFs, because theoretically
computed PSFs (using Tiny TIM software, for example) are usually not a very good
match to the observations. The noise in observed PSFs presents a problem, however:
when using a noisy PSF, the restoration algorithm ought to account for the fact that
both the image and the PSF are noisy.
This problem # the simultaneous estimation of both the object and the PSF # has been referred to
as blind deconvolution. Unfortunately, the naive application of a blind deconvolution procedure to
HST data will surely lead to difîculties since the trivial solution in which the PSF is estimated as
a point or impulse and the object is estimated as the data will produce a perfect match to the data.
Therefore, some type of constraints must be used to avoid both this solution and the equally trivial
one in which the object is estimated as a point source and the PSF is estimated as the data.
For a related problem in groundíbased astronomy, Schulz (1993) has recently proposed a
technique for processing a sequence of images degraded by turbulenceíinduced phase errors. With
this technique, the trivial solutions mentioned in the previous paragraph are avoided by using the
prior knowledge that the pointíspread functions are determined by phase errors distributed across
the telescope's aperture. In this paper, this approach is used to derive an image restoration technique
applicable to data acquired by the HST.
2. Data Model
We describe the CCD image data collected by the HST by using the model discussed in detail by
Snyder, Hammoud, and White (1993) in which the data collected at pixel y of the CCD detector,
d(y), is described as
d(y) = n obj (y) + n 0 (y) + g(y); y 2 Y ; (1)
where n obj (y) is the number of objectídependent photoíelectrons,n 0 (y) is the number of backgroundí
dependent photoíelectrons, g(y) is CCD readíout noise, y is a twoídimensional index, and Y is the
support region for the detector array. We assume that the random variables n obj (y), n 0 (y), and g(y)
206

Simultaneous Phase Retrieval and Deblurring 207
are statistically independent of each other and of n obj (y 0 ), n 0 (y 0 ), and g(y 0 ), for y 6= y 0 , and that the
objectídependent photoelectrons fn obj (y); y 2 Yg are Poisson distributed with the mean function
i(y) = fi(y)
X
x2X
h(y; x; ff)o(x); y 2 Y ; (2)
where fo(x); x 2 X g is the object's intensity function, ffi(y); y 2 Yg is the #atíîeld response
function, fh(y; x; ff); y; x 2 Y \Theta X g is the PSF, and ff is a collection of parameters that are
determined by the optical pathílength error (Goodman 1968) due to the spherical aberration, focusing
error, or other aberrations. Additionally, we assume that the backgroundídependent photoelectrons
are Poisson distributed with a known mean function fi 0 (y); y 2 Yg, and that the CCD readout
noise fg(y); y 2 Yg are identicallyídistributed Gaussian randomívariables with known mean m
and variance oe 2 .
Because of the parameterization by ff, the PSF cannot take an arbitrary functional form, and the
trivial solutions previously discussed are avoided. In the monochromatic, space invariant situation,
the parameterized PSF is obtained by spatially sampling the continuous function
h c (v; z; ff) =
fi fi fi fi
Z
A(u)e j 2‹
Ö W (u;ff) e j 2‹
Ö u\Delta(v\Gammaz) du
fi fi fi fi
2
; (3)
where A(u) is the telescope's pupil function, Ö is the wavelength of the observed light, W (u; ff)
is the optical pathílength error due to the optical aberrations, and u is a twoídimensional spatial
variable that indexes over the effective telescope pupil. The pathílength error is commonly speciîed
through either a pointíbyípoint description
W (u; ff) = ff(u) ; (4)
or a polynomial expansion
W (u; ff) =
N
X
n=1
ff n z n (u) ; (5)
where the polynomials fz n (u)g are typically chosen as the Zernike polynomials, orthogonal over
an annular region (Burrows 1990). For a system such as the HST, the discrete PSF can be created
by averaging the continuous PSF over detector regions. Mathematically, this procedure is described
as
h(x; y; ff) =
Z Z
w(x; y; v; z)h c (v; z; ff)dvdz ; (6)
where w(x; y; v; z) is a jointly discreteícontinuous function that describes the averaging and samí
pling performed by the CCD detector array.
3. Problem Statement
Following Snyder, Hammoud, and White (1993), we simplify the noise model by deîning the
modiîed data
~
d(y) = d(y) + oe 2 \Gamma m ; (7)
and by considering the situation in which the CCD readíout noise variance is large. In this case, the
modiîed data are approximately equal in distribution to a Poisson process whose mean function is
~ i(y) = i(y) + i 0 (y) + oe 2 : (8)
Therefore, the pseudoílogílikelihood takes the form
L(ff; o) = \Gamma
X
y2Y
~ i(y) +
X
y2Y
~
d(y) ln ~ i(y) ; (9)

208 Schulz & Cain
where terms that do not depend on the unknown object o or parameters ff have been omitted.
The maximum likelihood estimator of these parameters will maximize L for a particular data set
f ~
d(y); y 2 Yg, subject to the physical constraint that the object estimate be a nonnegative function.
An explicit, closed form solution to this problem is not known. Therefore, in the next section we
propose a technique for producing estimates numerically by using the expectationímaximization
(EM) algorithm (Dempster et al. 1977).
4. Solution
When ff is known, the PSF is known and the numerical procedure is the modiîcation of the
RichardsoníLucy iterations (Richardson 1972, Lucy 1974) described by Snyder (1990) and Snyder,
Hammoud, and White (1993):
o new (x) =
1
ï
fi(x)
o old (x)
X
y2Y
h(y; x; ff) fi(y) ~
d(y)
~ i old (y)
; (10)
where
ï
fi(x) =
X
y2Y
fi(y)h(y; x; ff) ;
and
~ i old (y) = fi(y)
X
x2X
h(y; x; ff)o old (x) + i 0 (y) + oe 2 :
When ff is unknown, the problem is more difîcult. For our initial investigation, we consider
the simpliîed situation for which the #atíîeld response function is a constant fi(y) = fi, the PSF
is spatially invariant, and the light is monochromatic. For this case, the EM algorithm for updating
both the object and the aberration parameters becomes,
o new (x) =
1
H 0
o old (x)
X
y2Y
h(y \Gamma x; ff old )
~
d(y)
~ i old (y)
; (11)
and
ff new
= arg max
ff
'' X
x2X
s(x) ln h(x; ff)
#
; (12)
where
H 0 =
X
y
h(y; ff) ; (13)
and where
s(x) = h(x; ff old
)
X
y2Y
o old
(y \Gamma x)
~
d(y)
~ i old (y)
: (14)
Notice that H 0 is independent of the aberration parameters ff. The optimization problem deîned in
Eq. 12 can be performed by using standard optimization techniques; however, the maximum need
not be found. Instead, if ff new is chosen such that
X
x2X
s(x) ln h(x; ff new
) Ö
X
x2X
s(x) ln h(x; ff old
) ; (15)
then the iterations will be a generalized EM algorithm (Dempster et al. 1977) and still have the
desirable property
L(ff new ; o new ) Ö L(ff old ; o old ) : (16)
Eq. 15 can be satisîed, for instance, by performing a few iterations of any numerical optimization
procedure.

Simultaneous Phase Retrieval and Deblurring 209
(a) (b)
Figure 1. Simulated data used in the example. (a) 256 \Theta 256 array of data. (b) 64 \Theta 64
subarray of data taken from crowded region in the upper left of (a). The data in (b) are
spatially magniîed by a factor of 4.
5. Example
To demonstrate the use of this technique, we performed a restoration for the simulated star cluster
provided by the Space Telescope Science Institute. The simulated image is shown in Fig. 1. The
crowded region of the true object and of the restored object are shown in Fig. 2. The original
point spread function (PSF), a noisy estimate of the PSF taken from the lower right region of the
data, and the PSF estimate are shown in Fig. 3. The PSF estimate was obtained by estimating
the pointíbyípoint phase errors due to the telescope's aberrations and using these errors to create a
predicted PSF.
6. Conclusions
In this paper we have discussed a technique for the simultaneous estimation of the unknown object
and system aberrationíparameters from data acquired by the Hubble Space Telescope. We have
presented an imageírestoration algorithm very similar to the modiîed RichardsoníLucy technique,
but generalized to allow for the estimation of the unknown system parameters, and, as such,
applicable to the blindídeconvolution problem. In contrast with other approaches in which the
spatial values of the pointíspread function are estimated, with our approach the aberration parameters
are estimated. This makes the technique directly applicable to situations in which the aberrationí
induced pointíspread function is spatially varying. Although the computational burden will surely
increase, the extension of the model to handle this situation is straightforward. The use of this
technique on simulated data has produced promising results; however, the technique should be
tested on real data, and its use with spatially varying, nonímonochromatic pointíspread functions
should be investigated.

210 Schulz & Cain
(a) (b)
Figure 2. (a) Original object. (b) Object estimate after 2000 iterations.
(a) (b) (c)
Figure 3. (a) True PSF. (b) Noisy PSF taken from the measured data. (c) PSF obtained
from the estimated optical pathílength error.

Simultaneous Phase Retrieval and Deblurring 211
References
Burrows, C. 1990, Hubble Space Telescope Optical Telescope Assembly Handbook,Space Telescope
Science Institute, Baltimore
Dempster, A. D., Laird, N. M., & Rubin, D. B. 1977, J. R. Stat. Soc. B, 39, 1
Lucy, L. B. 1974, AJ, 79, 745
Richardson, W. H. 1972, J. Opt. Soc. Am., 62, 55
Schulz, T. J. 1993, J. Opt. Soc. Am. A, 10, 1064
Snyder, D. L. 1990, in The Restoration of HST Images and Spectra, R. L. White & R. J. Allen, eds.,
Space Telescope Science Institute, Baltimore, 56
Snyder, D. L., Hammoud, A. M., & White, R. L. 1993, J. Opt. Soc. Am. A, 10, 1014
White, R. L. 1993, in Image Restoration Newsletter, R. Hanisch, ed., Space Telescope Science
Institute, Baltimore, 1, 11