## Abstract

Conventional centroid estimation fits a template point spread function (PSF) to image data. Because the PSF is typically not known to high accuracy, systematic errors exist. Here, we present an accurate centroid displacement estimation algorithm by reconstructing the PSF from Nyquist-sampled images. In absence of inter-pixel response variations, this method can estimate centroid displacement between two 32×32 images to sub-micropixel accuracy. Inter-pixel response variations can be calibrated in Fourier space by using laser metrology. The inter-pixel variations of Fourier transforms of the pixel response functions can be conveniently expressed in terms of powers of spatial wavenumbers. Calibrating up to the third-order terms in the expansion, the displacement estimation is accurate to a few micro-pixels. This algorithm is applicable to a new mission concept of performing mirco-arcsecond level relative astrometry using a 1 m telescope for detecting terrestrial exoplanets and high-precision photometry missions.

## 1. Introduction

Consider an *N*×*N* detector array. Let *I*_{mn} be the intensity measured by pixel in the *m*-th row and *n*-th column, which we shall use ordered pair (*m*,*n*),*m*,*n*=1,2,…,*N* to represent. A straightforward estimate of centroid is
1.1where (*x*_{mn},*y*_{mn}) is the coordinate of the centre of pixel (*m*,*n*) at the detector plane and the summation is over all the pixels within a window relevant to the centroid estimation. Because the pixels at large distances from the centre of the image detect very few photons, but are heavily weighted by their large coordinates, we restrict the window of pixels to a small size to ensure a good signal-to-noise ratio (SNR). Therefore, estimation (1.1) suffers systematic errors for using a relatively small window of pixels. The point spread function (PSF) fitting algorithms (Stone 1989) supersede this straightforward centroid estimation because they avoid the amplification of noise from multiplying large coordinates. The main challenge for the PSF fitting approach is the knowledge of the PSF. Computing PSF using a diffraction model requires precise knowledge of the optical system and wavefront aberrations that are usually not available. In the past, Gaussian functions have been popularly used (Stone 1989). However, in order to achieve micro-pixel level accuracy centroiding, a more precise PSF is needed. Mighell (2005) has done work on PSF fitting using the digital images based on 21-point damped sinc interpolation. In this paper, we work along the same line to reconstruct the PSF from the pixelated images using the fact that the PSF is a bandwidth-limited signal. Theoretically, a bandwidth-limited signal can be reconstructed to any accuracy as long as the sampling is above the Nyquist frequency, and the number of the sample is sufficiently large. Using 32×32 images, the truncation error is less than one micro-pixel for centroid displacement estimation.

Owing to the complex microstructure of detectors and the charge diffusion effect between pixels, the pixel detection response varies over the physical detection area of a pixel (Jorden *et al.* 1994), which is referred to as *intra*-pixel variation. The pixel counts is a convolution of the pixel response function (PRF) with the photon energy flux function. If each pixel has identical response function, pixel counts still represent values of a bandwidth-limited function, which we call an effective PSF, sampled at the pixel grid. The sampling theorem is applicable for reconstructing the effective PSF. Past measurements found dominant intra-pixel detection variations are indeed common to pixels (Kavaldjiev & Ninkov 1998). However, for micro-pixel centroid estimation, it is necessary to take into account the pixel-to-pixel differences of the response functions, which is referred to as *inter*-pixel variations. The inter-pixel variations make the pixel counts no longer precisely represent the sample values of a strictly bandwidth-limited function. To characterize inter-pixel variations, we use laser metrology to measure the pixel responses in Fourier space. It is convenient to parametrize the Fourier transforms of the PRFs in terms of powers of spatial wavenumbers using the Taylor series expansions. The leading order effect of inter-pixel variations is the average pixel response or flat-field response, which can be measured as response to a uniform E-field. The first-order correction is an effective geometric pixel location shift from a regular grid for each pixel. The second- and third-order corrections are quadratic and cubic polynomials of the spatial wavenumbers. As we include more terms in the expansion, the model and the centroid estimation become more accurate.

We use simulations to demonstrate the capability of calibrating the inter-pixel variations for performing micro-pixel level centroid estimations by including up to third-order terms in the Taylor series expansion of the Fourier transforms of the PRFs. This algorithm can be used by a new mission concept (Malbet *et al.* 2010) to perform micro-arcsecond level relative astrometry with a 1 m telescope and thus detect terrestrial exo-planets in the habitable zone.

In addition to precise astrometry, our pixel calibration technique is also applicable to high-precision photometry, which requires characterizing the PRFs. With the calibration of pixel responses and the reconstruction of PSF, we can remove photometric errors owing to pointing errors and detector response variations, which limit the performance. Because we measure the pixel responses in Fourier space, it is especially convenient for studying photometry when the images are not too under-sampled because we only need to measure the pixel response to a modest spatial frequency range to cover the bandwidth of the PSF.

## 2. Model and algorithm description

### (a) Pixel intensity model

To simplify the formulation, we assume that the image is stable over a sampling period, so that the temporal integration by the detector is simply an overall factor of the duration of the sampling period, which can be absorbed into an overall normalization factor.

A model describing the photo-electron counts recorded by pixel (*m*,*n*) is expressed as:
2.1where *Q*_{mn}(*x*,*y*) is the PRF of pixel (*m*,*n*) to a point illumination at (*x*,*y*) in the detector plane (Kavaldjiev & Ninkov 1998, 2001), *I*(*x*,*y*) is an input intensity function, and (*x*_{c},*y*_{c}) is the location of the centroid of the image. For a point source, which we will consider exclusively, *I*(*x*,*y*) is the PSF and proportional to the square of the focal plane E-field *E*(*x*,*y*):
2.2By Fourier optics (Goodman 1998), *E*(*x*,*y*) is related to the E-field at the input pupil of the telescope *E*_{i}(*x*,*y*) via a Fourier transform^{1}
2.3where is a normalization factor, *λ* is the wavelength of the light, *f* is the focal length of the telescope and *P*(*x*,*y*) is the aperture function, whose value is 1 inside the aperture and 0 outside the aperture. We will focus on a monochromatic light case and briefly discuss the polychromatic case in §5. The polychromatic case requires one extra integral over the wavenumber, which makes the formulation slightly complicated. Because the integration over the wavenumber is a linear operation, the main steps of derivation remain the same.

Performing a change of variable (*k*_{x},*k*_{y})=2*π*/(*λf*)(*x*^{′},*y*^{′}), we rewrite expression (2.3) as a two-dimensional Fourier integral for *E*(*x*,*y*):
2.4where (*k*_{x},*k*_{y}) represents a two-dimensional spatial frequency vector. Because *P*(*x*,*y*) vanishes outside the aperture of the telescope, *E*(*x*,*y*) is a two-dimensional bandwidth-limited signal. By the theorem of convolution, the Fourier transform of *I*(*x*,*y*) is the Fourier transform of *E*(*x*,*y*) convolved with the Fourier transform of the complex conjugate of *E*(*x*,*y*). Therefore, *I*(*x*,*y*) is also a bandwidth-limited function, with bandwidth being twice that of *E*(*x*,*y*) from the process of convolution. Let *D* be the diameter of the aperture of the telescope, the bandwidth of *I*(*x*,*y*) is then limited by |*k*_{x}|<2*πD*/(*λf*),|*k*_{y}|<2*πD*/(*λf*).^{2} Assuming the pixel size is smaller than the Nyquist sampling spacing *fλ*/(2*D*), by sampling theorem, the intensity function *I*(*x*,*y*) can be precisely reconstructed from a pixelated image given an infinitely large detection array. In practice, this is complicated by two things. First of all, we cannot have an infinite number of samples. Therefore, we have to truncate the reconstruction process at some finite size, which introduces truncation error. We found that for 32×32 array size (including about the 7th Airy ring), the truncation errors are small enough for micro-pixel centroid estimation. Secondly, the pixel intensity model involves PRF *Q*_{mn}(*x*,*y*), which depends on pixels in general. If *Q*_{mn}(*x*,*y*) does not depend on pixels, the pixel intensities correspond to sampled values of an effective PSF, which is the convolution of *I*(*x*,*y*) with the common PRF. Because convolution does not alter the bandwidth, the effective PSF has the same bandwidth as *I*(*x*,*y*). Thus, it is bandwidth limited and can be reconstructed from a pixelate image. In fact, experiments found that the leading order intra-pixel variation is common for all the pixels (Kavaldjiev & Ninkov 1998). For micro-pixel level centroiding, however, we cannot ignore the inter-pixel variations of *Q*_{mn}(*x*,*y*). The pixel counts no longer strictly correspond to sampled values of a bandwidth-limited function. Fortunately, the inter-pixel variation is only a small fraction of the total response, whose Fourier transform can be characterized by laser metrology fringe measurements.

To illustrate this, we rewrite the pixel intensity model as:
2.5where *a* is the spacing between adjacent pixels and we have used the following relations for Fourier transforms and :
2.6and
2.7The PRF-related effects are all captured by the Fourier transforms .

For the case that all pixels have identical response function,
2.8 no longer depends on indices *m*,*n*,
2.9Putting this in the expression (2.5), we obtain
2.10Defining an effective PSF
2.11*I*_{mn}(*x*_{c},*y*_{c}) corresponds to values of function at grid points (*x*,*y*)=((*m*+1/2)*a*−*x*_{c},(*n*+1/2)*a*−*y*_{c}), *m*,*n*=0,±1,±2,…. Because vanishes for |*k*_{x}|≥2*πD*/(*λf*) or |*k*_{y}|≥2*πD*/(*λf*), is bandwidth limited and can be reconstructed from *I*_{mn}(*x*_{c},*y*_{c}) by the sampling theorem.

However, response function *Q*_{mn}(*x*,*y*) depends on (*m*,*n*) in general. It is useful to define an average response function
2.12and its Fourier transform
2.13Because the common intra-PRF does not affect the property of the pixel intensities to be sampled values of a bandwidth-limited signal, it is convenient to factor out the common pixel response to parametrize the pixel-dependent terms as
2.14where *q*_{mn} represents the flat-field response of pixel (*m*,*n*), Δ*x*_{mn} and Δ*y*_{mn} represent effective geometric pixel location shifts along *x*- and *y*-directions for pixel (*m*,*n*) to deviate from a regular pixel grid location ((*m*+1/2)*a*,(*n*+1/2)*a*). *α*_{mn},*β*_{mn},*γ*_{mn} are parameters specifying quadratic behaviour in the amplitude of *Q*_{mn}(*k*_{x},*k*_{y}) that are pixel dependent. *a*_{mn},*b*_{mn},*c*_{mn},*d*_{mn} are third-order term coefficients. This parametrization is based on a Taylor series expansion of in terms of powers of *k*_{x} and *k*_{y}.

### (b) Centroid displacement estimation algorithm

We now turn to present the algorithm for estimating the centroid displacement between two images of the same PSF. The main idea is to resample one of the two images at a grid shifted from the default grid by an offset. When the resampled image best matches the second image, the offset at which we resampled the first image is the estimated centroid displacement. We call the first image (the image to be resampled) the *reference* image.

Based on pixel grids, we use the following discrete frequencies
2.15In terms of discrete frequencies, model (2.5) is expressed as:
2.16where we have defined
2.17Pixel response calibration measures and thus the ratio at a set of spatial frequencies covering the bandwidth of PSF. We use expression (2.14) to parametrize in low-order terms. If keeping the terms up to the third power, we have then parameters *q*_{mn}, Δ*x*_{mn},Δ*y*_{mn}, *α*_{mn},*β*_{mn},*γ*_{mn}, *a*_{mn},*b*_{mn},*c*_{mn},*d*_{mn}, which can be estimated by fitting expansion (2.14) to measure at a set of spatial frequencies pre-determined by the metrology system for pixel response calibration. With *q*_{mn},Δ*x*_{mn},Δ*y*_{mn},… estimated, it is straightforward to use expansion (2.14) to evaluate the ratios .

To reconstruct the PSF, we need to estimate . Let us choose the coordinate so that the centroid of the reference image is at the origin, i.e. *x*_{c}=0, *y*_{c}=0 and its pixel intensities are represented by *I*^{ref}_{mn}. can be solved by inverting reference intensity relation
2.18in the same fashion as inverse Fourier transform. may be viewed as the frequency representation of the image. With , expression (2.16) can be used to resample the reference image. The centroid displacement of the second image *I*_{mn} relative to the reference image can be estimated by solving a weighted nonlinear least-squares fitting as
2.19where *W*_{mn} is a weight factor for pixel (*m*,*n*). The choice of *W*_{mn} affects sensitivity to noise. Here, we simply use equal weights and save optimizing the weights to achieve best noise sensitivity as a future topic.

When we use discrete frequencies (2.15), the signal is parametrized as a periodical function in space with period *Na*. For *N*=32, the image encircles the 7th Airy ring; the PSF is small enough at the boundary, so that the truncation error is negligible.

## 3. Pixel response calibration using laser metrology

We measure the PRFs in Fourier space by observing the response of the detector to a sinusoidal illumination pattern. The sinusoidal pattern is generated by fringes from interference of two distant laser metrology beams so that the wavefront is close to a plane wave over the spatial extension of the detector. To make the identification of the laser fringes easier, we use an acousto-optic modulator (AOM) to offset one laser’s frequency from the other by a few hertz, so that the fringes move across the detector array. Therefore, the illuminating intensity may be modelled as:
3.1where and are the intensities of the two laser beams at the detector plane, is the spatial wavenumber of the laser fringes and Δ*ω* is the angular frequency difference between the two lasers introduced by the AOM. The output counts of pixel (*m*,*n*) is then
3.2where we have used the definition (2.7). The temporal variation of the intensity at each pixel is a sinusoidal function with angular frequency Δ*ω* plus a constant. We estimate the amplitude and phase of the sinusoidal temporal variation to get the complex Fourier transform because the overall intensities of the two metrology beams and can be easily measured. By having different separations between the laser beams or altering the distances between the metrology laser sources and the detector, we can generate fringes of different spatial frequencies and thus measure at a set of values of . Note that it is adequate to measure the Fourier transform of the the PRFs to cover the bandwidth of the PSF. For Nyquist-sampled images, we only need to measure the spatial frequency to *π*/*a*.

## 4. Results using simulated data

We use simulated data to validate our concept of the algorithm. In §4*a*, we first show that for an ideal detector whose pixels have the same response function, we can achieve sub-micro-pixel level accuracy in centroid displacement estimation.

### (a) Results for an ideal detector

In our simulation, the diameter of the telescope is *D*=1 m with focal length *f*=40 m to be consistent with mission concept Nearby Earth Astrometric Telescope (NEAT) (Malbet *et al.* 2010). We focus on a monochromatic source with wavelength *λ*=600 nm. In §5, we will discuss the case of a polychromatic source. We only consider an array of 32×32 pixels, i.e. the dimension *N*=32 because this is sufficient for achieving micro-pixel centroid estimation. The corresponding *λ*/*D* at the focal plane is 24 μm. The pixel size is 10 μm, which samples above the Nyquist frequency (1/(12 μm)). All the pixels have the same PRF, which is displayed in the left plot in figure 1. It is modelled as a Gaussian function multiplied by low-order polynomials (up to fourth order),
4.1where *x*,*y* are coordinates in the detector plane in unit of pixel with origin at the centre of the pixel, *r*_{g}=0.5, and the coefficients *c*_{i} are drawn from Gaussian random number generators with standard deviation being 0.05. The mean values are all 0 except for *c*_{0}, whose mean is 0.8. This PRF is roughly similar to the intra-pixel variation for a backside-illuminated charge-coupled device (CCD) (Piterman & Ninkov 2002). We include *λ*/20 r.m.s. low-order wavefront aberrations parametrized by the first 15 Zernike polynomials, whose amplitudes are randomly generated. The right plot in figure 1 displays the phase aberration over the telescope pupil used in our simulation. The two plots (*a* for *X* centroid and *b* for *Y* centroid) in figure 2 display the centroid displacement estimation errors for a grid of *X* and *Y* offsets within range [−0.5,0.5] pixel between the two images. These errors are owing to truncation and are no more than 0.1 μ-pixel.

### (b) Pixel response calibration and results

For realistic detectors, PRFs vary from pixel to pixel, which we call *inter-pixel* variations. Pixel response calibration measures the PRFs of all the pixels. To concentrate on pixel response calibration, we assume that the optical system and thus the PSF are stable and defer the study of wavefront stability to §4*d*. The Fourier transforms of the PRFs *Q*_{mn}(*x*,*y*) of all pixels are measured at a set of spatial frequencies covering the bandwidth of the PSF. To do this, we illuminate the detector using a sinusoidal intensity pattern, which may be achieved by interfering two metrology lasers at various spatial frequencies. The fringe pattern gives a sinusoidal illumination on the pixel array, and the separation of the two lasers determines the wavelengths of the fringes. AOM is used to offset the frequency of one laser from the other by a few hertz. The temporal variation of the pixel intensities can be used to estimate Fourier transforms as discussed in §3. We use expansion (2.14) to parametrize . For micro-arcsecond accuracy, we nominally keep the terms up to second order or third order in *k*_{x} and *k*_{y}.

To simulate the inter-pixel response variations, we make the low-order coefficients *c*_{i} in equation (4.1) pixel dependent by replacing *c*_{i} with , where are random numbers drawn from zero mean Gaussian random number generators with standard deviation 0.01. We also add 2 per cent white noise to the pixel response as a multiplicative factor as well as random geometric pixel location shifts for all the pixels along both *x*- and *y*-directions with 0.01 pixel r.m.s. The left plot in figure 3 displays typical intra-pixel responses of a 10×10 pixel array at the upper left corner of the 32×32 array. Here, the inter-pixel cross talk from diffusion is not displayed to avoid overlap. The zeroth-order calibration measures the pixel response to a flat field, i.e. Fourier transform of the PRF at spatial wavenumber (*k*_{x},*k*_{y})=(0,0). The right plot in figure 3 displays a typical flat-field calibration results. Figure 4 shows systematic centroid displacement estimation errors with only flat-field response calibration for different centroid offsets between the two images, whose range is [−0.5,0.5] pixel along both *x*- and *y*-directions. The next level of calibration measures the effective pixel locations deviating from a regular grid. The effective location of pixel (*m*,*n*) is estimated by fitting to the phase of the estimated Fourier transform at different values of spatial frequencies available from the metrology system. Figure 5 displays the estimated effective pixel location deviating from a regular grid. Now the errors shown in figure 6 are significantly reduced to tens of micro-pixels after including the effective pixel locations in the estimation. However, the errors are still large for micro-pixel level centroiding. We further include the second-order amplitude terms in expansion (2.14) and estimate coefficients *α*_{mn},*β*_{mn}, and *γ*_{mn} by fitting expansion (2.14) to the Fourier transforms measured by metrology calibration. The second-order coefficients are displayed in figure 7. The two plots in figure 8 display the centroid estimation errors for performing pixel response calibration that estimates the flat-field response, the effective pixel locations and the second-order amplitude corrections. Including the third-order terms in expansion (2.14) in our pixel response calibration enables us to achieve centroid displacement estimation accuracy of a few micro-pixels over [0.5,0.5] pixel displacement range between the two images along both *x*- and *y*-directions. The corresponding coefficients are displayed in figure 9. Figure 10 shows systematic centroid displacement estimation errors for including through the third-order phase terms. The r.m.s. is only a few micro-pixels.

### (c) Noise sensitivity

So far, we have not included any noise. We now study the sensitivity to photon shot noise, which is the dominant source of random errors. We set the total number of photons for each image to be 10^{8} and use the Poisson random number generator to simulate the photon shot noise. We simulated 1000 images with their centroid positions randomly distributed within ±0.5 pixel relative to a reference image along both *x*- and *y*-directions. The reference image has the same level of photon and shot noise. Figure 11 displays the Allan deviation of the estimated centroid displacements relative to a reference image as function of the total number of photons integrated. The green dashed line shows an empirical sensitivity formula
4.2for uncertainty of centroid estimation using an equal weight in the least-squares fitting, where *N*_{ph} is the total number of photons. To reach 10 mirco-pixel precision, we shall need about 2.2×10^{10} photons.

The Allan deviation shows the variation between 1000 images with respect to the same reference image. How does the noise in the reference image affect the results? It turns out that the noises in both images contribute to the displacement estimate at the same level. For the 1000 centroid displacements, the noises in the reference image lead to mostly an overall offset to all the displacements, i.e. it is an excellent approximation that
4.3where *d*(*A*,*R*_{1}),*d*(*B*,*R*_{1}),*d*(*A*,*R*_{2}) and *d*(*B*,*R*_{2}) represent the displacements of centroids of images A and B relative to reference images *R*_{1} and *R*_{2}, respectively. Figure 12 shows the centroid displacements of first 100 simulated images relative to two different reference images, which are the same image with two different realizations of photon shot noises. It is easy to see that the effect of different reference images is an offset common to all the 100 centroid displacements with less than 1 micro-pixel standard deviation. Letting *B* be *R*_{1} and *A* be *R*_{2}, relation (4.3) reduces to
4.4where we have used *d*(*A*,*A*)=*d*(*R*_{1},*R*_{1})=0. The estimated displacement between images *A* and *R*_{1} is insensitive to the choice of using *R*_{1} or *A* as the reference image. Noises in *A* and *R*_{1} contribute to the random errors in *d*(*A*,*R*_{1}) at the same level.

We are performing an equal weight least-squares fitting, it is possible to optimize weights *W*_{mn} to achieve the best sensitivity to noise, which will be a subject for future study.

### (d) Wavefront stability

In this section, we study the sensitivity of our algorithm to wavefront changes, i.e. the two images are taken with slightly different wavefronts. The left plot in figure 13 shows a randomly generated small low-order wavefront change with *λ*/1350 r.m.s. (by random generation) used in our study. We note that high-order wavefront aberrations generate speckles in the image plane while the low-order wavefront causes a change in the shape of the PSF. For PSF fitting, the centroid estimation is more sensitive to low-order wavefront aberrations. The main effect of a wavefront change between two observations is an overall centroid displacement error common to all the stars in the field. Even though this overall centroid displacement error is not calibrated, it causes no problem for differential astrometry because it is common to the two stars under observation. Our study in this section is to show that this error is indeed common for the two stars.

Because the two stars are imaged at different portions of the detector with inter-pixel response variations and we only simulated a detector of size 32×32 (nominal detector portion), we instantiate a second 32×32 detector (second detector portion) whose PRFs differ from our nominal detector as a result of a different realization of the random parameters in pixel response simulation. We study the effect of the wavefront change on the centroid displacement estimation as well as its dependency on the detect portion as a result from inter-pixel response variation.

Figure 13*b* shows the difference of the pixel responses between the second portion and the nominal portion, which comes from two different instantiations of random low-order polynomials. We first consider the case where the wavefront changes by an amount of *λ*/1350 r.m.s., as shown in the figure 13*a*, between taking the two images whose centroid displacement we are interested to estimate using the nominal detector. Figure 14 displays the centroid displacement errors owing to the *λ*/1350 r.m.s. low-order wavefront change for a grid of centroid offsets between the two images for the nominal portion. The main effect of a wavefront change between the two images is an overall constant independent of the actual displacement between the two images; the variation is only a few micro-pixels.

We now examine the effect of the same wavefront change on the second detector portion. Figure 15 displays the centroid errors for using the second detector. To avoid common error cancellation, we on purpose centre the reference image at [0.25, 0.25] pixel instead of [0,0] pixel as for the nominal detector portion. Again, we can see that the dominant effect of wavefront change between the two images is an overall offset in the centroid estimation and this offset differs from the offset for the nominal detector portion by only a few micro-pixels. Because the overall offset is not sensitive to the difference between two detector portions, the overall differential displacement caused by the wavefront difference cancels between two detector portions leaving a residual of a few micro-pixel r.m.s. as shown in figure 16. This shows even though different stars are centred at portions of the detector with different pixel responses, the effect owing to a wavefront change of *λ*/1350 r.m.s. leads to a common error for all the stars in the displacement estimation up to micro-pixel errors. We note again that the common error between stars cancels for differential astrometry.

## 5. Polychromatic effect

In this section, we discuss the polychromatic effect. For white light source-like stars, image intensity model (2.5) needs an extra integration over the photon wavelengths weighted with the source spectrum,
5.1where *S*(*λ*) is the source spectral energy density function and we have included the wavelength dependencies in Fourier transforms of the monochromatic PSF function and the pixel detection response function . As a leading order approximation, we ignore the spectral dependency in the pixel detection function. The broadband image model has the same expression as the monochromatic model with the following replacement
5.2Because the integral over the photon wavelength is a linear operation, the polychromatic signal is still a bandwidth-limited signal, whose bandwidth is determined by the shortest wavelength. As far as the pixelated images are Nyquist sampled, our algorithms works as for the case of monochromatic source. Figure 17 displays the centroid offset estimation errors for simulated chromatic pixelated images pairs displaced by various offsets along *x*- and *y*-directions assuming an ideal tophat pixel response detection.

The errors are at sub-micro pixel. A slight complication comes from the dependency of pixel response on the wavelength, which requires pixel response calibration using metrology at multiple wavelengths. The results in Piterman & Ninkov (2002) for backside illuminated CCD show weak dependency on the wavelengths especially at the longer wavelength. Because only the inter-pixel variations of the spectral dependency in pixel response lead to error and only the variations of these errors between different observations matter for differential astrometry, we expect this effect to be small in general and can be calibrated using laser metrology at a few wavelengths. We defer the study of this to a future work.

## 6. Conclusions and future work

We have presented a systematic frame work for accurate centroid displacement estimation and detector characterization. For an ideal detector with all the pixels having the same PRF, the PSF can be accurately reconstructed without any detector calibration. For realistic detectors whose PRFs share a dominant common portion and small pixel-to-pixel variations, it is possible to measure the pixel responses functions in Fourier space using laser metrology fringes. Measuring the pixel response in Fourier space is especially convenient for well-sampled images. Keeping a few low-order terms (e.g. third order) in the Taylor series expansion of the Fourier transform in wavenumbers, we can achieve centroid estimations accurate to a few micro-pixels assuming the pixel response does not depend on wavelength with wavefront stability better than one-thousandth of a wave in r.m.s.

In the future, we will include the wavelength and polarization dependency in the pixel responses and study the centroid estimation for different source spectral types using metrology to calibrate the detector at multiple wavelengths and different polarizations. Effects of response nonlinearity and observation strategy for avoiding detector saturation are also important for completing this work.

This technology will enable micro-arcsecond level differential astrometry using 1 m telescopes and thus detect earth-like exo-planets in the habitable zone. This method is also applicable to precise photometry.

## Acknowledgements

This work was prepared at the Jet Propulsion Laboratory, California Institute of Technology, under a contract with the National Aeronautics and Space Administration. ©2011 California Institute of Technology. Government sponsorship acknowledged.

## Footnotes

- Received April 22, 2011.
- Accepted July 5, 2011.

- This journal is © 2011 The Royal Society