Speckle noise and dynamic range in coronagraphic images
Abstract
This paper is concerned with the theoretical properties of high contrast coronagraphic images in the context of exoplanet searches. We derive and analyze the statistical properties of the residual starlight in coronagraphic images, and describe the effect of a coronagraph on the speckle and photon noise. Current observations with coronagraphic instruments have shown that the main limitations to high contrast imaging are due to residual quasistatic speckles. We tackle this problem in this paper, and propose a generalization of our statistical model to include the description of static, quasistatic and fast residual atmospheric speckles. The results provide insight into the effects on the dynamic range of wavefront control, coronagraphy, active speckle reduction, and differential speckle calibration. The study is focused on groundbased imaging with extreme adaptive optics, but the approach is general enough to be applicable to space, with different parameters.
1 Introduction
Direct imaging of faint companions or planets around a bright star is a very difficult task, where the contrast ratio and the angular separation are the observable parameters. The problem consists of detecting a faint source over a bright and noisy background, mainly due to the diffracted stellar light. High contrast ratios and small angular separations correspond to the most difficult case. Typically, for extrasolar giant planets, contrast ratios of about are expected in the near infrared (J,H,K bands), based on models for relatively young objects of about 100 Myr (Chabrier & Baraffe, 2000; Baraffe et al., 2003; Burrows et al., 2004). According to these models, older objects would be an order of magnitude fainter. Terrestrial planets are much fainter than giant planets, about 3 to 4 orders magnitudes fainter depending on the wavelength range.
In the case of groundbased observations with Adaptive Optics (AO), the residual uncorrected aberrations produce random intensity fluctuations of the background, which appear as speckles in the field. In direct noncoronagraphic high quality images, these speckles mainly appear at the position of the diffraction rings of the star. In Fig.1, we show two highquality point spread functions (PSF) obtained with the AO system at Palomar (Troy et al., 2000; Hayward et al., 2001), where the speckles are clearly visible. This phenomenon, also known as “speckle pinning” (Bloemhof et al., 2001), can be explained using an expansion of the point spread function (PSF) (Sivaramakrishnan et al., 2002, 2003; Bloemhof, 2003a; Perrin et al., 2003; Bloemhof, 2004b). An alternative approach using a statistical model (Aime & Soummer, 2004b) enables a deeper understanding of the phenomenon, especially from the statistical point of view, by providing information on the speckle variance. The effect of a coronagraph on speckle noise is well explained with this approach. In particular, we show that static aberrations produce residual speckle pinning after a coronagraph, which has important implications in high contrast imaging. Statistical properties of long exposure AO images were also studied by Fusco & Conan (2004).
The dynamic range of a coronagraphic image corresponds to the faintest companion that can be detected at a given position in the field, at the detection limit. It is usually expressed as a magnitude difference or an intensity ratio, relative to the unocculted central star at some signal to noise ratio (SNR) level. Although the dynamic range is a twodimensional map of the detection sensitivity in the field, it is often represented as a radial profile showing the magnitude difference as a function of the angular separation . A radial profile should be acceptable for most applications, especially when the coronagraph does not have particular asymmetries (Lyot, 1939; Aime et al., 2002; Soummer, 2005). Otherwise, a twodimensional map of sensitivity is required when the coronagraph shows an asymmetric response (Kuchner & Spergel, 2003; Kasdin et al., 2003; Rouan et al., 2000). Understanding, measuring and predicting the dynamic range is still one of the important issues in this field, with implications for instrumentation (design, observing strategies) and data reduction and analysis.
During the design phase of future planet finder instruments using Extreme Adaptive Optics (ExAO) and coronagraphy, for example the Gemini Planet Imager (Macintosh et al., 2004, 2006), or the ESO/VLT SPHERE (Beuzit et al., 2005; Fusco et al., 2005), it is necessary to anticipate what part of the observable parameterspace (contrast vs. separation) can be probed, and link it to the actual physical parameter space (mass vs. semi major axis). Such studies have been carried out by Brown (2004b, a, 2005) in the context of terrestrial planet searches.
When operating an existing highcontrast instrument, like the Lyot project coronagraph (Oppenheimer et al., 2004), the dynamic range has to be measured to evaluate the performance (Soummer et al., 2006; Hinkley et al., 2006). In the case of a detected object, photometry and astrometry (Digby et al., 2006; Sivaramakrishnan & Oppenheimer, 2006; Marois et al., 2006b) are necessary to help determine the objects characteristics. Dynamic range computations are also important in the case of nondetection, to determine which part of the parameter space has been probed by the experiment, and which physical objects can be ruled out.
The dynamic range is limited by the intensity fluctuations close to the star. These are due to several sources: speckle noise, photon noise, detector noise, background noise etc. Speckle noise is known to be the main limitation for high contrast imaging, either in direct images (Marois et al., 2003, 2005; Masciadri et al., 2005; Marois et al., 2006a), or coronagraphic images (Beuzit et al., 1997; Oppenheimer et al., 2001; Boccaletti et al., 2003, 2004; Hinkley et al., 2006). Speckles find their origin in wavefront imperfections (amplitude and phase errors), whether they correspond to uncorrected atmospheric residual errors (residual atmospheric speckles) or slowvarying wavefront caused by mechanical or thermal deformations. The main problem comes from these quasistatic speckles, which can be calibrated either using active prefocal methods Malbet et al. (1995); Give’on et al. (2006); Bordé & Traub (2006) or using postprocessing (Marois et al., 2000; Sparks & Ford, 2002; Marois et al., 2006a). An alternative approach based on nonredundant masking has been achieved by Lloyd et al. (2006). Assuming a good enough calibration of these quasistatic speckles, the physical limit of the system is set by the residual atmospheric aberrations. The noise limitation in high dynamic range images has been studied by several authors, using numerical simulations (Boccaletti, 2004; Cavarroc et al., 2006), or other theoretical or empirical approaches (Angel, 1994; Racine et al., 1999; Perrin et al., 2003; Guyon, 2005; Sivaramakrishnan et al., 2005).
In Aime & Soummer (2004b), we modeled the statistics of AOcorrected, direct images, and discussed the effect in terms of signal to noise ratio on these images qualitatively. In this paper, we examine the effects a coronagraph has on the properties of residual speckles. A semianalytical method to compute the dynamic range can be derived from the statistical properties of the speckle and photon noise. We compare our results with purely numerical simulations. Results in this paper are presented using an Apodized Pupil Lyot Coronagraph (APLC) as an example (Aime et al., 2002; Soummer et al., 2003a; Soummer, 2005), but are valid for any other type of mask coronagraph (Rouan et al., 2000; Kuchner & Traub, 2002; Soummer et al., 2003c; Mawet et al., 2005), pure apodizers or shaped pupils (Jacquinot & RoizenDossier, 1964; Nisenson & Papaliolios, 2001; Soummer et al., 2003b; Kasdin et al., 2003; Aime, 2005b). Furthermore, the statistical model can be modified to include both static and quasistatic aberrations, and we discuss the coherent interaction between residual atmospheric and quasistatic aberrations.
Our theoretical model and results apply for both space and groundbased imaging. However, we illustrate the results with simulations in the case of groundbased Extreme Adaptive Optics (ExAO) and coronagraphy.
2 Propagation through a coronagraph in the presence of aberrations
In this section we assume that all the static aberrations, mainly optical quality (polishing for example) and misalignment errors, can be represented at the entrance pupil of the telescope. The case of out of pupil aberrations is discussed by Marois et al. (2006c) and does not affect our monochromatic statistical model. We consider an instrument with an ExAO system and a generic coronagraph which can describe any type of mask designs (Rouan et al., 2000; Aime et al., 2002; Kuchner & Traub, 2002; Soummer et al., 2003c; Soummer, 2005; Mawet et al., 2005). The formalism also applies to the shaped pupil approach which corresponds to the case of direct imaging with apodization (Kasdin et al., 2003). The ExAO and telescope characteristics are chosen to be consistent with current or future projects on eightmeter telescopes and illustrations are given for an APLC.
These coronagraphs consist of optical filtering in four successive planes denoted 1,2,3,4 hereafter. The first plane corresponds to the entrance aperture (possibly apodized), the second plane is the focal plane where a mask is applied (opaque, graded or phaseshifting), the third plane corresponds to a relay pupil plane where a diaphragm is applied (the Lyot Stop), and finally the fourth plane corresponds to the final focal plane.
In the pupil plane, we model the wavefront using a static phase aberration , a residual random atmospheric phase , and amplitude aberrations , including the eventual apodization. The complex amplitude is:
(1) 
where the function describes the normalized aperture transmission: , and is the coordinate vector, used in both pupil and field. Following the notations of Aime & Soummer (2004b), we can write the wavefront complex amplitude at the entrance pupil as the coherent sum of three terms:
(2) 
where is a deterministic term corresponding to a perfect plane wave, is a deterministic complex term corresponding to the static aberrations, and is a random term with zero mean corresponding to the uncorrected part of the wavefront. The probability density function (PDF) of this complex amplitude is illustrated in Fig.2 without and with AO. Static aberrations are not included in this figure.
is defined as the mean of the complex amplitude, averaged over the pupil:
(3) 
is therefore the Strehl Ratio of the system (Hardy, 1998).
With , Eq.2 implies . Integrating this equation, we obtain:
(4) 
and therefore:
(5) 
Assuming that the phase errors are stationary over the aperture, we obtain:
(6) 
Note that the difference between Eq.5 and Eq.6 comes from the fact that is deterministic and is random: Although can be defined with zeromean over the aperture, each independent realizations of do not necessarily have an exactly zero average over the aperture.
The specific case of quasistatic aberrations is not considered in this section and will be treated in Sec.3.2. The three terms of Eq.2 are illustrated in Fig.3. The length of the vector is arbitrary in the figure, to illustrate that the modulus of is not unity and that the vectors , , and are defined according to the definitions above.
In the first focal plane, a coronagraphic mask is applied at the center of the image of the star. Writing the mask transmission as , allows us to accommodate any type of mask coronagraph, including Lyot, APLC, BandLimited, Phase Masks. For example, a classical hardedged Lyot coronagraph (or APLC), is described using a tophat function for . The complex amplitude of the wave in the focal plane is given by a scaled Fourier Transform (FT) of this pupil amplitude (Goodman, 1996):
(7) 
where the symbol denotes the scaled FT. For clarity, we will omit the wavelengthdependent scaling factors. For the complete chromatic formalism see for example Aime et al. (2002); Soummer et al. (2003a); Aime (2005a). In the next pupil plane, the complex amplitude before the Lyot stop is also the sum of three terms:
(8) 
where is the complex amplitude in the lyot stop plane for a perfect wavefront (Soummer et al., 2003a). The two other terms and correspond to the propagation of the terms and respectively. For example:
(9) 
The perfect coronagraph term and the term are shown in intensity in Fig.4. The coronagraphs rejects most of the starlight outside the image of the aperture in the Lyot plane for the perfect part of the wave, but most of the energy remains inside the aperture for the speckle part.
The Lyot stop is applied in this plane. In the case of an APLC, the Lyot stop is identical to the entrance pupil; in all other cases, the Lyot stop is undersized. With , and with the notations and , we obtain the complex amplitude in the final focal plane:
(10)  
where denotes the focal wave amplitude of the coronagraph in the perfect case, following the notations of Aime et al. (2002). The convolution product in Eq.10 has a negligible effect outside the mask area. Indeed, the spatial extension of is limited to the occulting mask area, and is a rapidly decreasing function (the Airy amplitude in a perfect case) whose characteristic size is . The result of the convolution does not extend much beyond the mask area, which is illustrated in Fig.5, where we show an example of the effect of the convolution term using a numerical simulation. This can also be explained by considering the propagation of phase ripples through a coronagraph and constructing a Bode diagram (Sivaramakrishnan et al., 2007).
Grouping all the deterministic terms together, we introduce
(11) 
and we obtain a similar expression to the case without coronagraph where (Aime & Soummer, 2004b):
(12) 
The wave amplitude after a coronagraph appears as a sum of a deterministic term , and a random term , at each position in the focal plane, outside the mask area. The deterministic term corresponds to the focal plane complex amplitude of the coronagraph in the presence of static aberrations. The random term , associated with the speckles, is identical to the case without coronagraph (Aime & Soummer, 2004b): the coronagraph has a negligible effect on the random part of the wavefront, as illustrated in Fig.4 and Fig.5. Formally, the effect of the coronagraph is to replace the wave amplitude without coronagraph by the coronagraphic amplitude . In the case of pure apodizers (shaped pupils), the direct apodized term is used.
The random term is non stationary in the field. The profile for can be computed from a simulation of the AO system, as we detail in Sec.4.2. Loworder aberrations can also be included in this description, but usually require a specific study, as for example in Sivaramakrishnan et al. (2005, 2007). The profile or the twodimensional map for can be computed independently, considering a perfect coronagraph in the presence of deterministic static aberrations (and normalized to the SR). Even in the case of an ideal coronagraph that cancels all the star light for a perfect wave, the deterministic term still contains the terms due to the static aberrations and will contribute to speckle pinning, as discussed below.
3 Statistical properties of direct or coronagraphic images
3.1 Statistical model and properties of speckles
3.1.1 Complex amplitude distribution
In this section, we discuss the distribution of the complex amplitude in the focal plane. We consider the case of monochromatic direct images for simplicity. The case of coronagraphic images is formally identical to the coronagraphic case, according to the approximations described in the previous paragraph (Eq.12). The focal plane complex amplitude is the Fourier transform of the pupil plane complex amplitude:
(13) 
The complex amplitude in the focal plane is therefore a sum of the random complex term weighted by the Fourier complex phasors. At the center of the image, the Fourier phasors vanish, so a special treatment for this particular case (and the transition region around it) is necessary (Soummer & Ferrari, 2007). Outside the central point of the image, the distribution of the complex amplitude can be derived using known results in Signal Processing, based on reasonable assumptions. We assume that the complex amplitude in the pupil plane can be represented by discrete values (an implicit assumption in any numerical simulation), and that the correlation of the complex amplitude between two points in the pupil plane decreases with distance between them. Under these hypotheses, it can be shown that the distribution of the complex amplitude in the focal plane is asymptotically circular Gaussian (Brillinger, 1981). We remind the reader here that if the real and imaginary parts of a complex number are Gaussian, its distribution is said to be Gaussian. If the real and imaginary parts are independent and have same variance, the distribution is said to be circular Gaussian and denoted . See Fig.2 of Aime & Soummer (2004b) for an illustration of the focal plane PDFs. The circularity of the Gaussian distribution is due to the Fourier phasors mixing the complex amplitudes in the complex plane in the Fourier integral (13), where varies between and . For positions in the focal plane such as , the Fourier phase term therefore varies between 0 and and this circularization occurs. The complex amplitude of the wave in the focal plane follows a circular Gaussian law, decentered by the mean of the amplitude and denoted: .
In Fig.6, we give an illustration of the distribution of the complex amplitude in the four successive planes of the coronagraph. This illustration is based on numerical simulations of a perfect APLC coronagraph and of an ExAO system. In the first pupil plane, we have a decentered crescent (see Fig.2). In the first focal plane, in this example at the top of an Airy ring, has a high absolute value, and the distribution is Gaussian, decentered by this amount. Detailed illustrations of the decentered Gaussian statistics as a function of the position in the field can be found in Aime & Soummer (2004b). In the following Lyot plane the coronagraph almost completely removes the perfect part of the wave (see Fig.4), and the resulting distribution is similar to the initial distribution of the complex amplitude in the pupil, but centered at the origin. Finally in the last focal plane, without static aberrations, , as and the result is a centered circular Gaussian distribution.
3.1.2 Intensity distribution
In this section we derive the Probability Density Function (PDF) of the intensity from that of the wave complex amplitude. Our problem is formally equivalent to the case of laser speckles added to a coherent background, which has been studied extensively (Goodman, 1975, 2006), in particular in the context of holography. We introduce the two intensity terms:
(14) 
Note that and are both functions of r, and that can describe both the direct or coronagraphic case, with and without static aberrations. Following Goodman, the joint PDF for the intensity and phase can be obtained from PDF of the complex amplitude, using the simple cartesianpolar change of variables , where the modulus of the Jacobian of this transformation is and integrating the phase to find the PDF for the intensity.
An alternative derivation of the PDF for the intensity is to consider the properties of Gaussian distributions. As discussed in the previous section, the speckle term is a circular Gaussian distribution . The instantaneous intensity corresponding to the complex amplitude of Eq.12 is simply:
(15)  
where and denote the real and imaginary parts. Using the properties of circular Gaussian distributions, and are independent Gaussian random variables of same variance . We can rewrite the intensity with real and imaginary terms of variance unity:
(16) 
where .
The random variable follows a decentered with two degrees of freedom: , with a decentering parameter , (Johnson et al., 1995, chap. 29).
The probability density function for is therefore:
(17) 
where is the regularized confluent hypergeometric function and the confluent hypergeometric function defined as:
(18) 
Finally, the probability density function of the intensity is:
(19) 
This expression is equivalent ^{1}^{1}1The Mathematica software (Wolfram, 1999) can be used to derive these expressions, and the equivalence between Eq.19 and Eq.20 can be verified easily using the functions Simplify and FunctionExpand. to the modified Rician distribution derived by Goodman (1975) and used by Cagigal & Canales (1998); Canales & Cagigal (1999); Cagigal & Canales (2000); Canales & Cagigal (2001):
(20) 
where denotes the zeroorder modified Bessel function of the first kind. The Rician distribution is illustrated in Fig.7. A comparison between the Rician model and simulation data will be presented in Sec.3.1.4.
An interesting particular case is when the background is zero and only the speckle term is present. Making in Eq.20 (this happens at the zeros of the perfect PSF or using a perfect coronagraph), the PDF reduces to:
(21) 
This PDF corresponds to the wellknown negative exponential density for a fully developed speckle pattern (e.g. laser speckle pattern) (Goodman, 2000). Finally, the distribution at photon counting levels can be obtained performing a PoissonMandel transformation of the high flux PDF in Eq.20. An analytical expression of this PDF has been given in (Aime & Soummer, 2004a).
The mean and variance of the intensity can be obtained by several ways. A first method (Goodman, 1975, 2000), is to express the mean intensity and the second order moment of the intensity as a function of and . The second order moment for the intensity is the fourth order moment for the complex amplitude: (omitting the variables r for clarity), which can be simplified using the properties of Gaussian distributions, with we obtain: . A second method is to derive a general analytical expression for the moments of the Rician distribution. This can be be obtained either from the definition of the moments of Eq.20 (Goodman, 1975), or computing the derivatives of the moment generating function (Aime & Soummer, 2004a). The instantaneous intensity in the focal plane (Eq.15) can be written as:
(22) 
Since (circular Gaussian distribution), the mean intensity is simply the sum of the deterministic diffraction pattern with a halo produced by the average of the speckles: or , respectively for direct and coronagraphic images. The variance also finds a simple analytical expression, and we have:
(23) 
The variance associated with photodetection can be added to this expression to obtain the total variance ,
where is the variance associated to the poisson statistics: .
The total variance is therefore:
(24) 
In the case of direct images, the term corresponds to the perfect Point Spread Function (PSF) scaled to the SR. In the case of coronagraphic images, the focal plane intensity is not invariant by translation, and therefore it is technically not a true PSF. However, we will use the term “coronagraphic PSF” for simplicity and to follow the general usage in the community. The term is a function of the radial distance , which can describe an actual AO halo. These PSFs and halo structures have been studied analytically (Moffat, 1969; Racine, 1996; Racine et al., 1999). It is also possible to determine the halo profile directly from numerical simulations, and an illustration of and is shown in Fig.8. The long exposure PSF profile is the sum of these two contributions, the halo clearing effect for higher SR (Sivaramakrishnan et al., 2001) is clearly visible between the two figures. The shape of the halo is due to the spatially filtered wavefront sensor (Poyneer & Macintosh, 2004) used in this simulation.
In Fig.9 we show the effect of a coronagraph on the term, while the term is left unmodified as explained in Sec.2. In this figure we only consider one of the previous two AO cases. In this example, the coronagraph is good enough to render the constant background term negligible when compared to the speckle term .
3.1.3 Effect of a coronagraph on speckle pinning
In a direct, non coronagraphic image, the term coupling the deterministic and random parts in Eq.22 corresponds to the socalled “speckle pinning”, discussed by several authors (Bloemhof et al., 2001; Sivaramakrishnan et al., 2002; Bloemhof, 2003b; Sivaramakrishnan et al., 2003; Perrin et al., 2003; Bloemhof, 2004b), using an expansion of . In this expansion approach, the PSF consist of a sum of terms, some of which contain a multiplicative factor , corresponding to the diffracted field (the Airy amplitude for a circular aperture). These terms inherit the zeros of and contribute all together to the pinned speckles. Pinned speckles therefore have zero amplitude at the Airy nulls, and nonpinned speckles remain at these locations. The first order PSF expansion term, denoted by in Perrin et al. (2003), corresponds to the pinned speckles and the second order () to nonpinned speckles. Higher order terms contribute to pinned and nonpinned speckles. This expansion approach provides particularly interesting insight into the spatial properties and symmetries of the speckles for each order of the expansion (Sivaramakrishnan et al., 2003; Perrin et al., 2003; Bloemhof, 2004b). The expansion and our decomposition are therefore very similar. In Eq.22, our pinning term includes the contribution of all pinning terms from the infinite expansion. However, our constant term in the decomposition is not 1 and this term can include the effect of static aberrations, as discussed below.
The analysis of the statistical properties of the speckles enables a deeper understanding of the pinning phenomenon. As shown in Eq.22 and Eq.3.1.2, the pinned speckle term of Eq.22 does not contribute to the mean intensity, but only contributes to the variance. A numerical simulation is used in Fig.7 to illustrate the Rician distribution for a direct image, at three different positions in the field: one at the top of an Airy ring (strong pinning effect), one at a PSF zero (no speckle pinning) and one at an intermediate position. Speckle pinning can be well illustrated by the analysis of these PDFs, as speckle intensity and fluctuations are amplified by the term . This can be seen in the PDFs in Fig.7, where the widths increase with . Depending on the amplitude of the Airy pattern at successive rings, the intensity is alternatively large and small and the variance of the speckles is amplified accordingly by the coherent part of the wave amplitude , with corresponding intensity . At the zeroes of the PSF, no amplification occurs and the statistics is equivalent to that of a fully developed speckle pattern (exponential statistics). It is important to note that speckles fluctuations are not fully cancelled there (Fig.7), but simply not amplified and their statistics is that of laser speckles.
The effect of a coronagraph on speckles can be well explained from this statistical modeling. Formally, as shown in Eq.12, the effect of a coronagraph is to replace the telescope PSF by the onaxis coronagraphic PSF after a coronagraph . This term includes the effect of static aberrations if they are included in the model.

In the case of static aberrations in the pupil (due to polishing and alignment errors for example), these aberrations propagate through the coronagraph, according to the description given in Sec.2, or in Sivaramakrishnan et al. (2007). The deterministic response of the coronagraph (Eq.11) therefore includes the effect of these static aberrations. As a result, static aberrations leaking through a coronagraph contribute to speckle pinning, even in the case of a perfect coronagraph. Such effect can be produced for example by dead actuators on the deformable mirror (Sivaramakrishnan et al., 2006).
Following Aime & Soummer (2004b), we can illustrate the effect of a coronagraph on speckle noise, by breaking down the total variance (Eq.24) into two contributions, and . The term is the part of the variance that can be removed by a coronagraph by changing , thus affecting speckle pinning. In Fig.10 and Fig.11 we illustrate these variance contributions for direct and coronagraphic images. In the case of the direct image, the ringing aspect of the variance corresponds to speckle pinning (where the variance is amplified). The coronagraphic variance profile is smooth (no amplification of the speckles is seen at the location of the diffraction rings). Details of the numerical simulation are given in Sec.4.2.
3.1.4 Test of the Rician distribution with numerical simulations
Tests of the Rician statistics on real data using the Lick observatory AO system have been carried out by Fitzgerald & Graham (2006), showing that the Rician model is consistent with the data. Several complications exist in real data, so we test the Rician distribution on numerical simulations to determine whether the model is acceptable in a simple case consistent with the hypothesis used in the physical model, without any additional complicating circumstances or noise. We used PAOLA to generate 10000 AOcorrected instantaneous phase screens corresponding to an ExAO system, on a 8m telescope (we used a telescope geometry compatible with Gemini or VLT). The parameters chosen for this simulation include 44 actuators across the pupil, an integration time and time lag of 0.5ms for a magnitude star, and observations in the Hband. The atmosphere include a typical profile for CerroPachon and the seeing is assumed to be 1.4 arcec. The Strehl Ratio of these simulated images is .
In each image, we extract the intensities values along a radius to construct 50 intensity series in the focal plane at these 50 pixel locations. An example of the first 200 intensity values at an arbitrary location is given in Fig.12.
For each of these 50 points, we performed a Maximum Likelihood estimation of the parameters and , assuming a Rician distribution for the data. The Likelihood has been computed as function of the two parameters for the unbinned data and maximized using optimization routines of Mathematica. We then perform the and KolomogorovSmirnoff test statistics on these results. We use ten identical identical bins for the test. In Fig.13, we show two examples of binned data with error bars (here due to the Poisson statistics), superimposed with the Rician distribution for the estimated and parameters.
Our implementation of the KolmogorovSmirnoff test is based on a MonteCarlo estimation of the distribution. Indeed, when parameters are estimated from the data (here and ), the distribution of the KS test is not known analytically and must be estimated from MonteCarlo simulations. For that, we calculated a distribution of KS test values for 100,000 random drawings of random Rician data where the parameters have to be estimated from the data. This empirical distribution of the KS test values is used to generate the empirical right tail values for the test statistic.
The results of the two tests (right tail values) are given in Fig.14, for each of the 50 points along the radial axis (between 0 and ). It is interesting to note that for the first two points (onaxis intensity distribution), both tests conclude that the Rician model is incompatible with the data at the 5% level. Although the statistical properties at the central point are not directly relevant for evaluating the detection limits in high contrast images, this question is important in itself as it corresponds to the distribution of Strehl Ratio. This problem has been studied by Soummer & Ferrari (2007) where it is shown that the distribution at the center of the image is given by a “reversed” noncentral Gamma distribution. This problem was studied independently by Gladysz et al. (2006); Gladysz (2006); Christou et al. (2006) with similar results.
Outside of the center, most pixels pass the two tests except for about 2 pixel locations, which are rejected, by one test or the other. This result is consistent with the 5% level we chose, given the total number of locations (50 pixels) and we conclude that the Rician model is compatible with the simulated data. We have verified on a few sets of simulations that the location of these pixel failing the test is not relevant and simply due to the statistics. On the contrary, the central pixel fails the tests systematically.
It is interesting to compare the estimated parameters with the actual parameters known from the simulation. In Fig.15 we compare the actual profile known from the simulation data and compare it to the reconstructed profile from the estimation of the Rician statistics. Both curves show a very good agreement, even somewhat surprisingly at the center where the Rician model is wrong. This result is compatible with the satisfactory results of the test statistics. Fitzgerald & Graham (2006) discussed this type of approach as a possible way to reconstruct a telescope perfect PSF from the statistics of the PSFs. They note that this procedure might be difficult to implement with real data, but is at least verified in simulations in this simple case.
3.2 Effect of quasistatic aberrations on the statistical properties
Observations at high contrast have shown that the main dynamic range limitations are due to longlived quasistatic speckles (Beuzit et al., 1997; Oppenheimer et al., 2001; Marois et al., 2003; Boccaletti et al., 2003, 2004; Soummer et al., 2006). Practical solutions have been proposed to overcome the quasistatic speckles and will be implemented on the next generation of high contrast imagers (GPI, Sphere) (Wallace et al., 2006; Marois et al., 2006a) or considered for space projects (Bordé & Traub, 2006; Give’on et al., 2006). In ExAO coronagraphic images, Hinkley et al. (2006) have studied residual speckle lifetimes, and evidenced two types of speckles, one with lifetimes of a few seconds and the other with lifetimes of a few hundred seconds. Although there is no reason to believe the universality of these particular values, this problem of quasistatic speckles is general enough to need further understanding and theoretical insight. Quasistatic speckles are produced by the slow deformations of the optics (thermal and mechanical) that happen for example as the telescope slews to follow the star. The typical time scale for these slow variations of aberrations are tens of minutes.
We can generalize the approach of Sec.2 to include two types of aberrations with different time scales. We will assume that the error wavefront at the entrance pupil consists of the coherent addition of two terms: a fastvarying wavefront related to atmospheric residuals, and a quasistatic aberration term. We assume that quasistatic aberrations can be described as the sum of a deterministic static aberration and a slowly evolving, zeromean complex aberration: the static aberrations correspond to the polishing errors and aberrations produced by the optical design and actual misalignments in the system. Quasistatic aberrations corresponds to slow zeromean fluctuations around this static (and therefore deterministic) contribution. Under these hypothesis, we can add a fourth term to the wavefront decomposition of Eq.2:
(25) 
where, as previously, A corresponds to the perfect wave and to the static aberrations (deterministic). The error terms and correspond respectively to the AOcorected and quasistatic aberrations. We assume that the lifetime of and are respectively and . We consider the longest lifetime as the time unit, and we denote by the number of fastspeckle realizations during a slowspeckle lifetime:
(26) 
where and are the number of realizations for each processes during a long exposure of duration . Following the results of the previous sections, we can write the field amplitude in the focal plane during successive short time intervals as:
(27) 
As previously, where corresponds to the focal plane field (direct of coronagraphic) in the perfect case without aberrations, and to the static aberration (e.g. polishing errors or actual optical design aberrations). corresponds to the AOcorrected atmospheric aberrations and to the quasistatic aberrations. Outside of the central point in the focal plane, and are independent circular Gaussian distributions with zero mean: and with respective lifetimes and . Also, the spatial properties of and are different as they result from different aberrations sources, and their corresponding variances and can be generated considering the appropriate spatial power spectra.
During a time unit corresponding to the lifetime of a slow speckle, the intensity is:
(28)  
The expected value of is , where follows the Gaussian distribution . Therefore, the mean intensity is:
(29) 
Reminding that . This expression is consistent with the simple case with one type of speckle (Eq.3.1.2), as it is expressed here for an exposure .
The variance of the intensity is defined as:
(30)  
The covariance can be easily computed using the expansion of the high order moments of a complex Gaussian vector as a function of its high order cumulants (McCullagh, 1987; Ferrari, 2006). The covariance finds a simple expression:
(31) 
With and , we obtain immediately the covariance and the variance (for k=l):
(32) 
Finally, the variance of the intensity (Eq.30), for an exposure , becomes:
(33)  
where the variance for a short exposure :
(34) 
is consistent with the case with one type of speckle (Eq.3.1.2), if .

In this generalized expression of the variance, the pinning term , still exist. However, it now includes a term corresponding to the pinning of quasistatic speckle. The coefficient (the ratio of the speckle lifetimes) can be very large, which makes this pinning particularly efficient: a quasistatic halo which is times lower than produces the same pinning effect in direct images. This pinning contribution is still directly affected by a coronagraph, which reduces (or cancels) the term. In the case of a real, imperfect coronagraph, although no atmospheric pinning might be present, residual pinning of quasistatic speckles may occur due to the amplification by the factor if the quasistatic aberrations are not sufficiently small.

The cross term corresponds to speckle pinning of the atmospheric by quasistatic aberrations (coherent amplification of the atmospheric speckles by the quasistatic speckles). If the system is dominated by the quasistatic speckles ( large), this term is negligible compared to the halo .
It is possible to specify the performance of a coronagraph to avoid speckle pinning, using this expression for the variance. Pinning terms are negligible compared to the halo terms if:
(35) 
assuming the cross term negligible. In the case of only one type of speckles, or when the system is dominated by quasistatic speckles, this condition simply becomes (Aime & Soummer, 2004b). The necessary rejection of the coronagraph can be defined as a ratio of the coronagraphic PSF to the direct PSF . A similar approach was used by Bloemhof (2004a) to quantify the coronagraph effect on pinned speckles using the expansion approach.
4 Dynamic range in coronagraphic images
In this paper we adopt a definition of the dynamic range based on an expression of the Signal to Noise Ratio (S/R), assuming a typical level detection limit. Based on the knowledge of the noise statistics, more advanced methods are possible, where the dynamic range can be defined in terms of probability of detection and false alarm (Michel & Ferrari, 2003; Ferrari et al., 2006). In particular, because of the nonGaussian statistics (exponential distribution for a perfect coronagraph), confidence levels have to be carefully defined, as they do not correspond to the usual values for a Gaussian distribution, C. Marois et al. (2007, in preparation). In this section we define analytically the dynamic range in a coronagraphic experiment, using the variance expressions obtained above. This approach enables a semianalytical method. We present some simulations results for the dynamic range.
4.1 Signal to Noise ratio and dynamic range
For simplicity in the calculation of the variance (Sec.3.2), we expressed the intensity for a time unit as a sum of short exposures (Eq.28). In the expression of the variance we derived (Eq.33), the terms , , correspond implicitly to a number of photons during an exposure . However, it is convenient in practice to define normalized intensity terms for a single photon at the entrance aperture and scale them using the star flux . A normalization term is therefore necessary in order to use such normalized intensities. Recalling that a long exposure consists of exposures of durations , we obtain the final expression of the variance for a long exposure :
(36)  
The photon noise is obtained by applying the same normalization, multiplying Eq.29 by :
(37) 
and introducing for consistency with the case with one type of speckle. The signal from a companion is , where is the planet flux and is a function that describes the offaxis response of the coronagraph to a companion. can be computed independently. For simplicity, we use the maximum intensity of the offaxis companion for the signal, but other approaches such as matched filtering can also be used (Soummer et al., 2006). The signal to noise ratio (S/N) is:
(38) 
This expression for the signal to noise ratio can be immediately converted into an expression for the dynamic range, assuming a given S/N level. The dynamic range, denoted is simply the intensity ratio between an offaxis companion and the star that produces a given S/N level: . Since , we have:
(39) 
Other noise contributions can be easily added in this expression if necessary, but a detailed study of any specific instrument is beyond the scope of this paper.
4.2 Semianalytical method based on the statistical model
The dynamic range (Eq.39), can be computed directly from numerical simulations of the terms , , and . We follow the equations used in the model (Eq.3 to Eq.12) to construct numerical estimates for and (or , , if applicable). The term is obtained from its definition (Eq.3). We consider the average over the aperture and over a few independent realizations, generated with PAOLA (Jolissaint et al., 2006) and including both amplitude and phase screens . As described in Sec.2, the term corresponds to the Strehl Ratio, and is used to generate the profiles: direct and coronagraphic normalized PSFs are calculated in the perfect case without aberrations, and multiplied by .
The term is generated following the notation of Eq.10: , where is the zeromean complex amplitude at the aperture. The term is obtained averaging independent realizations of . It is easy to verify that the total intensity remains correctly normalized with this procedure.
Finally, radial profiles are generated by averaging these results azimuthally. Satisfactory profiles are obtained with only a few realizations (typically a few tens), where purely numerical simulations require large numbers of independent realizations to estimate the same variance. Once the various and terms have been generated, they can be readily combined according to Eq.36 or Eq.38 to produce variance or detection plots. Multiple instrumental cases can be easily studied this way, combining direct or various coronagraphic images with different cases of AO is done rapidly using our analytical expressions (for example to study the effect of stellar magnitude on the dynamic range). This method also enables the study of quasistatic speckles, which is otherwise computationally expensive.
We compared the variance obtained with the analytical model with the variance obtained from a purely numerical simulation, computing a large number of independent realizations. This process is timeconsuming because of the wavefront propagation through the coronagraph for each phase screen. Excellent agreement with a full numerical simulation of the variance is obtained with a few realizations , using the semianalytical method (Fig.16). We note that the model gives an incorrect estimate of the variance behind the focal plane occulting mask. This is expected given the approximations described in Eq.11 and Eq.12. This is not a problem when studying dynamic range, as these expressions are weighted by the offaxis transmission of the coronagraph, and results are only relevant outside the inner working angle of the instrument.
4.3 Results and comparison with numerical simulations
We present some dynamic range simulation results, using the semianalytical method based on Eq.39. We consider an eightmeter telescope, an AO system with 44 actuators across the aperture, a spatially filtered wavefront sensor (Poyneer & Macintosh, 2004), and an open loop frequency of 2.5 kHz. The atmospheric simulation uses typical seeing values and profiles at Mauna Kea. A few cases have been simulated, with stellar magnitudes ranging from V=4 to V=8. Strehl Ratios are obtained between and . Atmospheric scintillation is also simulated as amplitude screens with PAOLA, which implements the method of Roddier (1981). We assumed an instrumental throughput of and an arbitrary exposure time of . The atmospheric lifetime is assumed to be for observations in the Hband and the stellar magnitude is V=4. The coronagraph in the simulation is an APLC similar to the one under study for GPI (Soummer, 2005). With a detection level for this simulation, we illustrate the dynamic range (including photon and speckle noise only) in Fig.17. In this example, photon and speckle noise are approximately at the same level inside the AO control region. Outside the control region, the dynamic range is set by the speckle noise level.
It is interesting to note that the relative position of the speckle and photon contributions do not depend on exposure time (Eq.39). Indeed, only the speckle lifetime and stellar flux have an effect on the relative contributions of speckle and photon noise to the dynamic range. This is also true when both fast and quasistatic speckles are present.
In terms of the contributions of and to the dynamic range as defined in Aime & Soummer (2004b), the APLC coronagraph is almost perfect in this simulation and the part of the variance corresponding to speckle pinning, , has been completely removed. The case illustrated in Fig.10 and Fig.11 corresponds to the same simulation and shows the effect of the APLC coronagraph canceling the pinning contribution completely.
Finally, we performed a simulation including both fast and quasistatic aberrations. We considered for that a simple approach where the power spectrum of the quasistatic aberrations follows a power law , with a lifetime of 600s (10min). This choice is consistent with the assumption that quasistatic aberrations result from polishing errors, usually well described by such power laws. We considered both photon and speckle noise, according to Eq.39. In this example, although the quasistatic wavefront error is very small (10nm RMS), it dominates the error budget and limits the dynamic range inside the control radius of the AO system. This result depends directly on the power spectrum chosen to generate these aberrations. This could explain the results obtained with the Lyot project coronagraph where the dynamic range plots do not show the expected haloclearing region within the control radius of the AO system Hinkley et al. (2006). A slightly shallower power spectrum, or a higher level of quasistatic aberrations would easily reproduce the observed dynamic range by filling the cleared region completely. This effect could be due mainly to the presence of broken actuators (Sivaramakrishnan et al., 2006) in the AO system.
4.4 Discussion on the effects of speckle reduction techniques
The noise budgets expressed in Eq.36 and Eq.37 enable a general understanding of high contrast imaging for exoplanet detection. In the case of direct imaging without coronagraph, the main limitation is set by speckle pinning which appears as an amplification of the variance levels at the position of the maxima of the diffracted pattern. In the absence of static aberrations, pinning happens at the top of the Airy rings. In the presence of static aberrations and even with a perfect coronagraph, residual pinning exists at the position of the maxima of the coronagraphic PSF corresponding to the propagated static aberrations through the coronagraph. For example, in the case of loworder static aberrations described by Zernike polynomials, the location of the maxima of the PSF (direct or coronagraphic) dictates where pinning happens. These locations are no longer at the top of the perfect Airy pattern in this case. Quasistatic aberrations also add to the problem by creating additional pinning terms (Eq.34). This phenomenon was observed at the AEOS telescope with the Lyot project coronagraph, where pinned speckles have been associated with dead actuators on the deformable mirror (Sivaramakrishnan et al., 2006), or AO waffle mode (Makidon et al., 2005).
In reality, coronagraphs are not perfect and suffer from various sources of errors, including chromaticity, manufacturing imperfections, alignment errors, etc. The cancellation of speckle pinning in the final error budget provides a constraint on the coronagraph performance. This can be achieved based on the comparison between the terms and described above.
With current state of the art AO systems and coronagraphs alone, it is clear that the residual noise terms (even if the pinning term is well controlled) would not enable the detection of giant planets with typical contrasts of at 0.5 arcsec, or Earthlike planets with contrast of at 0.1 arcsec. Additional speckle reduction techniques have been developed (Bordé & Traub, 2006; Give’on et al., 2006) to improve the contrast performance further in combination with wavefront control and coronagraphs. These techniques affect the term by modifying the pupil aberrations to create a dark zone in the search region of the field. In the case of groundbased high contrast imagers, for example with the GPI speckle calibration system (Wallace et al., 2006), the speckle reduction will only affect the term, as the sensing involved in the technique is slow. Some speckle reduction techniques may also affect the residual speckle pinning in some region of the field, as they create a static wavefront aberration which creates a particularly asymmetric PSF presenting a dark zone (Codona et al., 2006; Serabyn et al., 2006).
The combination of coronagraphy and speckle nulling will therefore remove most of the terms contributing to the noise variance . The limit for the association of coronagraphy and speckle nulling is set by the level of residual fast atmospheric aberrations. Also, by decreasing the term , the speckle reduction system puts a stronger constraint on the coronagraph, whose performance must suppress pinned speckles to a level below other speckles in the dark region produced by speckle nulling.
Finally, the remaining parts of the speckle noise which have not been removed by a coronagraph or by the speckle reduction system, can be removed further using postprocessing differential methods, combining multiwavelength information from dualband imaging, or an integral field spectrograph Marois et al. (2000); Sparks & Ford (2002); Marois et al. (2005), or using differential rotation Marois et al. (2006a). These differential methods affect all the residual speckle variance contributions (pinning and halo) in the same way. The combination of these three stages (coronagraph, speckle reduction and speckle processing) is necessary to reach the necessary dynamic range regime for exoplanet imaging. This approach has been chosen for the design of GPI (Macintosh et al., 2006), and have been implemented in the laboratory by Trauger & Traub (2007).
5 Conclusion
This paper develops theoretical insight in high contrast imaging, using a statistical model to understand the performance limitations of coronagraphic experiments. We discuss the application of the statistical model (Aime & Soummer, 2004b) to the case of coronagraphic images, and show that the statistical model is valid outside of the focal plane mask occulting area. We show some limitations of the statistical model at the central point of the image, and a modification of the model is proposed by Soummer & Ferrari (2007) to solve this problem.
It has been confirmed by several observations that the main limitation for high contrast imaging is due to the presence of quasistatic speckles. These do not average out over time, and are difficult to calibrate. We presented the first theoretical attempt to understand the effect of quasistatic speckles and their interaction with other aberrations (static and residual atmospheric). We obtained an analytical expression for the variance of the intensity, which includes speckle and photon noise in the presence of static, quasistatic and rapidlyvarying aberrations. This result enables the use of a semianalytical method, which has been compared successfully with purely numerical simulations. This semianalytical method enables fast simulations and the possibility to compare and combine easily various types of parameters for instrument design studies. It also enables a breakdown of the noise into specific contributions (pinning, halo, speckle, photon). This method can be used to model real observations and extract information on relative noise contributions from various sources (residual atmospheric noise, quasistatic noise).
The model also provides an understanding of the speckle pinning effect in high Strehl images, where bright speckles appear at the location of the diffraction rings. This effect, explained by a coherent amplification of the speckle noise, and can be cancelled by a coronagraph. The theory provides insight on the required level of performance for the coronagraph to effectively cancel speckle pinning. In the presence of quasistatic aberrations, we showed that pinning also exists between the perfect part of the wave and quasistatic aberrations, but is weighted by the ratio of the lifetimes (which can be a large number of the order of a hundred or more). This puts strong requirements on the level of static and quasistatic aberrations and may explain the results currently obtained with the Lyot project coronagraph, the first instrument of its kind to use an extreme AO system. Detail modeling of the Lyot project data, and specific simulations for the future generation of high contrast imager (GPI), will be done in the future.
This model was developed for the monochromatic case. Generalization to wide band is not straightforward, as the speckle noise is highly correlated with wavelength, which makes this generalization difficult, but can be used for multiwavelength speckle subtractions. Also, a generalization of the model to the case of subtraction residuals in multiwavelength imaging would be particularly interesting and is a topic for future study.
Although this presentation is focused on the groundbased case, the same formalism can be directly used for spacebased observations, with similar results, but different orders of magnitudes (in particular for the speckle lifetimes). An identical expression of the variance of the noise was established independently by Shaklan et al. (2005) for the TPF error budget.
Once a coronagraph has removed the speckle pinning contribution, the dynamic range can be improved further by use of speckle reduction techniques (speckle nulling) or speckle subtraction (ADI, multiwavelength or polarization). The effect of speckle nulling on the noise budget is now well understood from the expression of the variance. The model can be used to draw the requirements on the relative performance of the coronagraph and speckle reduction system. This theory brings an understanding of the effects on the dynamic range of the various elements in a high contrast instrument, including wavefront control, coronagraphy, speckle reduction and differential calibration techniques. All these stages are crucial and necessary to reach high contrast imaging, and they are going to be implemented in the next generation of high contrast imagers, currently in development.
References
 Aime (2005a) Aime, C. 2005a, PASP, 117, 1012
 Aime (2005b) —. 2005b, A&A, 434, 785
 Aime & Soummer (2004a) Aime, C. & Soummer, R. 2004a, in EAS Publications Series, ed. C. Aime & R. Soummer, 89–101
 Aime & Soummer (2004b) Aime, C. & Soummer, R. 2004b, ApJ, 612, L85
 Aime et al. (2002) Aime, C., Soummer, R., & Ferrari, A. 2002, A&A, 389, 334
 Angel (1994) Angel, J. R. P. 1994, Nature, 368, 203
 Baraffe et al. (2003) Baraffe, I., Chabrier, G., Barman, T. S., Allard, F., & Hauschildt, P. H. 2003, A&A, 402, 701
 Beuzit et al. (2005) Beuzit, J.L., Mouillet, D., Dohlen, K., & Puget, P. 2005, in SF2A2005: Semaine de l’Astrophysique Francaise, ed. F. Casoli, T. Contini, J. M. Hameury, & L. Pagani, 215–+
 Beuzit et al. (1997) Beuzit, J.L., Mouillet, D., Lagrange, A.M., & Paufique, J. 1997, A&AS, 125, 175
 Bloemhof (2003a) Bloemhof, E. E. 2003a, ApJL, 582, L59
 Bloemhof (2003b) —. 2003b, ApJL, 582, L59
 Bloemhof (2004a) —. 2004a, ApJ, 610, L69
 Bloemhof (2004b) —. 2004b, Optics Letters, 29, 2333
 Bloemhof et al. (2001) Bloemhof, E. E., Dekany, R. G., Troy, M., & Oppenheimer, B. R. 2001, ApJ, 558, L71
 Boccaletti (2004) Boccaletti, A. 2004, in EAS Publications Series, ed. C. Aime & R. Soummer, 165–176
 Boccaletti et al. (2003) Boccaletti, A., Chauvin, G., Lagrange, A.M., & Marchis, F. 2003, A&A, 410, 283
 Boccaletti et al. (2004) Boccaletti, A., Riaud, P., Baudoz, P., Baudrand, J., Rouan, D., Gratadour, D., Lacombe, F., & Lagrange, A.M. 2004, PASP, 116, 1061
 Bordé & Traub (2006) Bordé, P. J. & Traub, W. A. 2006, ApJ, 638, 488
 Brillinger (1981) Brillinger, D. 1981, Time Series : Data Analysis and Theory (McGrawHill)
 Brown (2004a) Brown, R. A. 2004a, ApJ, 610, 1079
 Brown (2004b) —. 2004b, ApJ, 607, 1003
 Brown (2005) —. 2005, ApJ, 624, 1010
 Burrows et al. (2004) Burrows, A., Sudarsky, D., & Hubeny, I. 2004, ApJ, 609, 407
 Cagigal & Canales (1998) Cagigal, M. P. & Canales, V. F. 1998, Optics Letters, 23, 1072
 Cagigal & Canales (2000) —. 2000, JOSA, 17, 1312
 Canales & Cagigal (1999) Canales, V. F. & Cagigal, M. P. 1999, AO, 38, 766
 Canales & Cagigal (2001) —. 2001, Optics Letters, 26, 737
 Cavarroc et al. (2006) Cavarroc, C., Boccaletti, A., Baudoz, P., Fusco, T., & Rouan, D. 2006, A&A, 447, 397
 Chabrier & Baraffe (2000) Chabrier, G. & Baraffe, I. 2000, ARA&A, 38, 337
 Christou et al. (2006) Christou, J. C., Gladysz, S., Redfern, M., Bradford, L. W., & Roberts, L. C. J. 2006, in Proc AMOS technical conference, Maui, Hawaii, 1014 September 2006, 528–537
 Codona et al. (2006) Codona, J. L., Kenworthy, M. A., Hinz, P. M., Angel, J. R. P., & Woolf, N. J. 2006, in Groundbased and Airborne Instrumentation for Astronomy. Edited by McLean, Ian S.; Iye, Masanori. Proceedings of the SPIE, Volume 6269, pp. (2006).
 Digby et al. (2006) Digby, A. P., Hinkley, S., Oppenheimer, B. R., Sivaramakrishnan, A., Lloyd, J. P., Perrin, M. D., Roberts, L. C. J., Soummer, R., Brenner, D., Makidon, R. B., Shara, M., Kuhn, J., Graham, J. R., Kalas, P. G., & Newburgh, L. 2006, Submited to ApJ
 Ferrari (2006) Ferrari, A. 2006, in EAS Publications Series, 85–101
 Ferrari et al. (2006) Ferrari, A., Carbillet, M., Serradel, E., Aime, C., & Soummer, R. 2006, in Proceedings IAU Colloquium No. 200, ed. C. Aime & J. Vakili, F. Darcourt
 Fitzgerald & Graham (2006) Fitzgerald, M. P. & Graham, J. R. 2006, ApJ, 637, 541
 Foo et al. (2005) Foo, G., Palacios, D. M., & Swartzlander, Jr., G. A. 2005, Optics Letters, 30, 3308
 Fusco & Conan (2004) Fusco, T. & Conan, J.M. 2004, Journal of the Optical Society of America A, 21, 1277
 Fusco et al. (2005) Fusco, T., Rousset, G., Mouillet, D., & Beuzit, J.L. 2005, in SF2A2005: Semaine de l’Astrophysique Francaise, ed. F. Casoli, T. Contini, J. M. Hameury, & L. Pagani, 219–+
 Give’on et al. (2006) Give’on, A., Kasdin, N. J., Vanderbei, R. J., & Avitzour, Y. 2006, JOSA A, 23
 Gladysz (2006) Gladysz, S. 2006, PhD thesis, National University of Ireland, Galway
 Gladysz et al. (2006) Gladysz, S., Christou, J. C., & Redfern, M. 2006, in Advances in Adaptive Optics II. Edited by Ellerbroek, Brent L.; Bonaccini Calia, Domenico. Proceedings of the SPIE, Volume 6272, pp. (2006).
 Goodman (1975) Goodman, J. 1975, in topics in applied physics : laser speckle and related phenomena, Dainty Ed. (springer verlag berlin)
 Goodman (1996) Goodman, J. 1996, Introduction to Fourier Optics (Mac Graw Hill)
 Goodman (2000) Goodman, J. W. 2000, Statistical Optics (Wiley Classics Library), wiley classics library ed. edn., Wiley classics library (New York: WileyInterscience), 576
 Goodman (2006) —. 2006, Speckle Phenomena in Optics (Englewood, Colo.: Roberts and Company Publishers), 384
 Guyon (2005) Guyon, O. 2005, ApJ, 629, 592
 Hardy (1998) Hardy, J. W. 1998, Adaptive Optics for Astronomical Telescopes (Oxford University Press, USA), 448
 Hayward et al. (2001) Hayward, T. L., Brandl, B., Pirger, B., Blacken, C., Gull, G. E., Schoenwald, J., & Houck, J. R. 2001, PASP, 113, 105
 Hinkley et al. (2006) Hinkley, S., Oppenheimer, B. R., Soummer, R., Sivaramakrishnan, A., Roberts, L. C. J., Kuhn, J., Makidon, R. B., Perrin, M. D., Lloyd, J. P., Kratter, K., & Brenner, D. 2006, Submited to ApJ
 Jacquinot & RoizenDossier (1964) Jacquinot, P. & RoizenDossier, B. 1964, Progress in Optics, Vol. 3 (Wolf, E.)
 Johnson et al. (1995) Johnson, N., Kotz, S., & Balakrishnan, N. 1995, Continuous Univariate Distributions, Vol. 2 (John Wiley & Sons Inc.)
 Jolissaint et al. (2006) Jolissaint, L., Véran, J.P., & Conan, R. 2006, Journal of the Optical Society of America A, 23, 382
 Kasdin et al. (2003) Kasdin, N. J., Vanderbei, R. J., Spergel, D. N., & Littman, M. G. 2003, ApJ, 582, 1147
 Kuchner & Spergel (2003) Kuchner, M. J. & Spergel, D. N. 2003, ApJ, 594, 617
 Kuchner & Traub (2002) Kuchner, M. J. & Traub, W. A. 2002, ApJ, 570, 900
 Lloyd et al. (2006) Lloyd, J. P., Martinache, F., Ireland, M. J., Monnier, J. D., Pravdo, S. H., Shaklan, S. B., & Tuthill, P. G. 2006, ApJ, 650, L131
 Lyot (1939) Lyot, B. 1939, MNRAS, 99, 580
 Macintosh et al. (2006) Macintosh, B., Graham, J., Palmer, D., Doyon, R., Gavel, D., Larkin, J., Oppenheimer, B., Saddlemyer, L., Wallace, J. K., Bauman, B., Evans, J., Erikson, D., Morzinski, K., Phillion, D., Poyneer, L., Sivaramakrishnan, A., Soummer, R., Thibault, S., & Veran, J.P. 2006, in Advances in Adaptive Optics II. Edited by Ellerbroek, Brent L.; Bonaccini Calia, Domenico. Proceedings of the SPIE, Volume 6272, pp. (2006).
 Macintosh et al. (2004) Macintosh, B. A., Bauman, B., Wilhelmsen Evans, J., Graham, J. R., Lockwood, C., Poyneer, L., Dillon, D., Gavel, D. T., Green, J. J., Lloyd, J. P., Makidon, R. B., Olivier, S., Palmer, D., Perrin, M. D., Severson, S., Sheinis, A. I., Sivaramakrishnan, A., Sommargren, G., Soummer, R., Troy, M., Wallace, J. K., & Wishnow, E. 2004, in Advancements in Adaptive Optics. Edited by Domenico B. Calia, Brent L. Ellerbroek, and Roberto Ragazzoni. Proceedings of the SPIE, Volume 5490, pp. 359369 (2004)., ed. D. Bonaccini Calia, B. L. Ellerbroek, & R. Ragazzoni, 359–369
 Makidon et al. (2005) Makidon, R. B., Sivaramakrishnan, A., Perrin, M. D., Roberts, L. C., Oppenheimer, B. R., Soummer, R., & Graham, J. R. 2005, PASP, 117, 831
 Malbet et al. (1995) Malbet, F., Yu, J. W., & Shao, M. 1995, PASP, 107, 386
 Marois et al. (2000) Marois, C., Doyon, R. ., Racine, R. ., & Nadeau, D. 2000, PASP, 112, 91
 Marois et al. (2005) Marois, C., Doyon, R., Nadeau, D., Racine, R., Riopel, M., Vallée, P., & Lafrenière, D. 2005, PASP, 117, 745
 Marois et al. (2003) Marois, C., Doyon, R., Nadeau, D., Racine, R., & Walker, G. A. H. 2003, in EAS Publications Series, ed. C. Aime & R. Soummer, 233–243
 Marois et al. (2006a) Marois, C., Lafrenière, D., Doyon, R., Macintosh, B., & Nadeau, D. 2006a, ApJ, 641, 556
 Marois et al. (2006b) Marois, C., Lafrenière, D., Macintosh, B., & Doyon, R. 2006b, ApJ, 647, 612
 Marois et al. (2006c) Marois, C., Phillion, D. W., & Macintosh, B. 2006c, in Groundbased and Airborne Instrumentation for Astronomy. Edited by McLean, Ian S.; Iye, Masanori. Proceedings of the SPIE, Volume 6269, pp. (2006).
 Masciadri et al. (2005) Masciadri, E., Mundt, R., Henning, T., Alvarez, C., & Barrado y Navascués, D. 2005, ApJ, 625, 1004
 Mawet et al. (2005) Mawet, D., Riaud, P., Absil, O., & Surdej, J. 2005, ApJ, 633, 1191
 McCullagh (1987) McCullagh, P. 1987, Tensor methods in statistics, Monographs on statistics and applied probability (Wiley, New York)
 Michel & Ferrari (2003) Michel, O. & Ferrari, A. 2003, in EAS Publications Series, ed. C. Aime & R. Soummer, 129–146
 Moffat (1969) Moffat, A. F. J. 1969, A&A, 3, 455
 Nisenson & Papaliolios (2001) Nisenson, P. & Papaliolios, C. 2001, ApJ Letters, 549
 Oppenheimer et al. (2004) Oppenheimer, B. R., Digby, A. P., Newburgh, L., Brenner, D., Shara, M., Mey, J., Mandeville, C., Makidon, R. B., Sivaramakrishnan, A., Soummer, R., Graham, J. R., Kalas, P., Perrin, M. D., Roberts, L. C., Kuhn, J. R., Whitman, K., & Lloyd, J. P. 2004, in Advancements in Adaptive Optics. Edited by Domenico B. Calia, Brent L. Ellerbroek, and Roberto Ragazzoni. Proceedings of the SPIE, Volume 5490, pp. 433442 (2004)., ed. D. Bonaccini Calia, B. L. Ellerbroek, & R. Ragazzoni, 433–442
 Oppenheimer et al. (2001) Oppenheimer, B. R., Golimowski, D. A., Kulkarni, S. R., Matthews, K., Nakajima, T., CreechEakman, M., & Durrance, S. T. 2001, AJ, 121, 2189
 Perrin et al. (2003) Perrin, M. D., Sivaramakrishnan, A., Makidon, R., Oppenheimer, B. R., & Graham, J. R. 2003, ApJ, 596, 702
 Poyneer & Macintosh (2004) Poyneer, L. A. & Macintosh, B. 2004, Journal of the Optical Society of America A, vol. 21, Issue 5, pp.810819, 21, 810
 Racine (1996) Racine, R. 1996, PASP, 108, 699
 Racine et al. (1999) Racine, R., Walker, G. A. H., Nadeau, D., Doyon, R., & Marois, C. 1999, PASP, 111, 587
 Roddier (1981) Roddier, F. 1981, Prog. Optics, Volume 19, p. 281376, 19, 281
 Rouan et al. (2000) Rouan, D., Riaud, P., Boccaletti, A., Clénet, Y., & Labeyrie, A. 2000, PASP, 112, 1479
 Serabyn et al. (2006) Serabyn, E., Wallace, J. K., Troy, M., Mennesson, B., Haguenauer, P., Gappinger, R. O., & Bloemhof, E. E. 2006, in Advances in Adaptive Optics II. Edited by Ellerbroek, Brent L.; Bonaccini Calia, Domenico. Proceedings of the SPIE, Volume 6272, pp. (2006).
 Shaklan et al. (2005) Shaklan, S. B., Marchen, L., Green, J. J., & Lay, O. P. 2005, in Techniques and Instrumentation for Detection of Exoplanets II. Edited by Coulter, Daniel R. Proceedings of the SPIE, Volume 5905, pp. 110121 (2005)., ed. D. R. Coulter, 110–121
 Sivaramakrishnan et al. (2003) Sivaramakrishnan, A., Hodge, P. E., Makidon, R. B., Perrin, M. D., Lloyd, J. P., Bloemhof, E. E., & Oppenheimer, B. R. 2003, in HighContrast Imaging for ExoPlanet Detection. Edited by Alfred B. Schultz. Proceedings of the SPIE, Volume 4860, pp. 161170 (2003)., 161–170
 Sivaramakrishnan et al. (2001) Sivaramakrishnan, A., Koresko, C. D., Makidon, R. B., Berkefeld, T., & Kuchner, M. J. 2001, ApJ, 552, 397
 Sivaramakrishnan et al. (2002) Sivaramakrishnan, A., Lloyd, J. P., Hodge, P. E., & Macintosh, B. A. 2002, ApJL, 581, L59
 Sivaramakrishnan & Oppenheimer (2006) Sivaramakrishnan, A. & Oppenheimer, B. R. 2006, ApJ, 647, 620
 Sivaramakrishnan et al. (2006) Sivaramakrishnan, A., Oppenheimer, B. R., Perrin, M. D., Roberts, L. C., Makidon, R. B., Soummer, R., Digby, A. P., Bradford, L. W., Skinner, M. A., Turner, N. H., & Ten Brummelaar, T. A. 2006, in IAU Colloq. 200: Direct Imaging of Exoplanets: Science & Techniques, ed. C. Aime & F. Vakili, 613–616
 Sivaramakrishnan et al. (2007) Sivaramakrishnan, A., Soummer, R., Oppenheimer, B., & Pueyo, L. 2007, Submited to ApJ
 Sivaramakrishnan et al. (2005) Sivaramakrishnan, A., Soummer, R., Sivaramakrishnan, A. V., Lloyd, J. P., Oppenheimer, B. R., & Makidon, R. B. 2005, ApJ, 634, 1416
 Soummer (2005) Soummer, R. 2005, ApJ, 618, L161
 Soummer et al. (2003a) Soummer, R., Aime, C., & Falloon, P. E. 2003a, A&A, 397, 1161
 Soummer et al. (2003b) Soummer, R., Aime, C., Ferrari, A., & Falloon, P. E. 2003b, in HighContrast Imaging for ExoPlanet Detection. Edited by Alfred B. Schultz. Proceedings of the SPIE, Volume 4860, pp. 211220 (2003)., ed. A. B. Schultz, 211–220
 Soummer et al. (2003c) Soummer, R., Dohlen, K., & Aime, C. 2003c, A&A, 403, 369
 Soummer & Ferrari (2007) Soummer, R. & Ferrari, A. 2007, ApJL, in press
 Soummer et al. (2006) Soummer, R., Oppenheimer, B. R., Hinkley, S., Sivaramakrishnan, A., Makidon, R. B., Oppenheimer, B. R., Digby, A. P., Brenner, D., Kuhn, J., Perrin, M. D., Roberts, L. C. J., & Kratter, K. 2006, in EAS Publications Series, ed. C. Aime & M. Carbillet
 Sparks & Ford (2002) Sparks, W. B. & Ford, H. C. 2002, ApJ, 578, 543
 Trauger & Traub (2007) Trauger, J. T. & Traub, W. A. 2007, Nature, 446, 771
 Troy et al. (2000) Troy, M., Dekany, R. G., Brack, G., Oppenheimer, B. R., Bloemhof, E. E., Trinh, T., Dekens, F. G., Shi, F., Hayward, T. L., & Brandl, B. 2000, in Proc. SPIE Vol. 4007, p. 3140, Adaptive Optical Systems Technology, Peter L. Wizinowich; Ed., 31–40
 Wallace et al. (2006) Wallace, J. K., Bartos, R., Rao, S., Samuele, R., & Schmidtlin, E. 2006, in Advances in Adaptive Optics II. Edited by Ellerbroek, Brent L.; Bonaccini Calia, Domenico. Proceedings of the SPIE, Volume 6272, pp. (2006).
 Wolfram (1999) Wolfram, S. 1999, The Mathematica book, Fourth Edition (Cambridge University Press)