EDP Sciences
Free Access
Issue
A&A
Volume 569, September 2014
Article Number A28
Number of page(s) 10
Section Astronomical instrumentation
DOI https://doi.org/10.1051/0004-6361/201423680
Published online 12 September 2014

© ESO, 2014

1. Introduction

Almost twenty years ago Mayor & Queloz (1995) made the first detection of an exoplanet, and more than a thousand are now recognized1. External occulters will be formidable instruments for performing detailed spectralanalysis of exo-Earths that have previously been detected by other methods.

The first use of an external occulter dates back to Evans (1948), who used an occulting disk set in front of a very small coronagraph of Lyot (1939). Very soon solar astronomers have identified diffraction effects of a circular occulter and considered various kinds of apodized, shaped, or multiple occulters for space missions, as well described in the review paper of Koutchmy (1988). Because of implementation difficulties, toothed or multiple disks are preferred over radially graded ones, which have been described in Purcell & Koomen (1962) and Newkirk Jr & Bohlin (1965). Next to the successful Solar and Heliospheric Observatory mission (Brueckner et al. 1995), future solar coronagraphy envisages formation-flying spacecrafts with an occulter of diameter 1.5 m at about 150 m in front of a small telescope (Vives et al. 2009).

The parameters are completely different for exoplanets,. In rounded numbers, the occulting angle is about 0.2 arcsec against 2000 arcsec for the Sun, and a large 4 m diameter telescope is needed with an occulter of 50 m diameter at a distance of 80 000 km. The positive aspect of exoplanet experiments is that diffraction phenomena are less difficult to simulate numerically than in the solar case. This is due to two reasons. There are fewer than 20 Fresnel zones compared with 8000 in the solar case, and the main source of diffracted light is a point source instead of a huge 1/2° extended object. Nevertheless, the problem remains difficult because of the high expected dynamicrange.

Progression of ideas for exoplanets first followed developments for the Sun, with little interaction, but now the theoretical calculations and instrumental developments have resulted in much more developments than in the solar case. In his review of space astronomy Spitzer (1962) mentioned the technique of apodization envisaged by solar astronomers. Elaborated petal-shaped occulters can already be found in Marchal (1985), and Copi & Starkman (2000) proposed a circular apodizing transmission on a square mask. Studies of shaped and apodized occulters have been reported in many publications, such as Cash (2006), Arenberg et al. (2007), Vanderbei et al. (2007), Soummer et al. (2010), Cash (2011) and Wasylkiwskyj & Shiri (2011), to cite just a few. Although partially transparent petal shaped occulters are envisaged by Shiri & Wasylkiwskyj (2013), the most advanced technological developments are for petaled star shades. The two-step procedure leading to an optimal shaped occulter is well described in Kasdin et al. (2013). These authors explained that they first seek for the optimal variable transmission of an apodized occulter and then chose a sufficient number of petals for the shaped occulter to best approach the theoretical result given by the variable transmission. We here focus only on the first part of this approach, that is, on the search for an ideal transmission, leaving the final design of a shaped occulter to a future work.

Our purpose is to show that the occulter can be optimized advantageously from the telescope focal plane, which will minimize the level of light that is diffracted by the star in a zone of interest of the focal plane for the exoplanet. The resulting illumination on the telescope aperture appears to be apodized from center to edge. The integrated flux on the telescope pupil is no longer minimized. The observation is nevertheless improved because most of the residual flux is trapped behind the geometrical image of the mask, in a blind area for the observation of the exoplanet. Differences in transmission with an occulter that minimizes the flux over the aperture are small, but the gain in the residual background light at the level of the exoplanet range from 1.9 to 50 depending on the spectral bandwidth, which might correspond to a gain in integration time of a factor 3.6 to 2500 for certain observing conditions, inphoton counting mode.

The paper is organized as follows: the fundamental relations that describe the complex amplitudes diffracted by the occulter in the aperture and focal planes are given in Sect. 2. The procedure for optimizing the transmission of a circularly symmetric apodized occulter to minimize the flux on the telescope aperture or for selected regions of the focal plane are described in Sect. 3. Section 4 contains results on the optimization problem and a discussion. Additional information is given in two appendices based on an analytic approach of the Fresnel diffraction of a pure circular disk and of a linearly apodized disk.

2. Expressions of complex amplitudes and intensities in the aperture and focal planes

We denote with the overall diameter of the occulter. For an occulting mask, it is convenient to write its transmission as t(r) = 1 − f(r), with the attenuation function f(r) constrained by0 ≤ f(r) ≤ 1 for | r | ≤ Ω, and f(r) = 0 for | r | > Ω. For a wave of unit amplitude, the Fresnel diffraction at the telescope aperture at a distance z from the occulter (see Fig. 1) can be written as (1)where τz(r) = exp(r2/λz) is a quadratic phase-term corresponding to a diverging lens, and J0(r) is the Bessel function of the first kind.

This relation is equivalent to Eq. (4) of Vanderbei et al. (2007), neglecting here the term of propagation of a plane wave. In general, this integral does not admit an analytical solution, except for f(r) = 1, as shown by Aime (2013) and outlined here in Appendix A. Therefore, Eq. (1), which appears as the Hankel transform of f(r)τz(r), must be evaluated numerically. Computations of these transforms are delicate, as discussed by Lemoine (1994) for example, who recommended a sampling of the function based on zeroes of the Bessel function J0. In practice, after several tests described in Appendix A and B and below, we directly used the function NIntegrate of Mathematica, which performs very well for this kind of calculation, thanks to recent improvements of numerical integration.

This wavefront ψ(r) arrives onto the aperture of the telescope of transmission P(r), and we denote with ϕ(r) = P(r)ψ(r) the complex amplitude of the wave going through the telescope aperture. For the sake of simplicity, we assume in the following a perfectly circular telescope of diameter 2R, centered on the optical axis of the system, so that P(r) = Π(r/R), defining here for convenience that the box distribution Π(r) equals 1 for | r | < 1 and 0 elsewhere. Note that a more realistic telescope with a central obscuration can be easily inserted into the calculations. Using the circular symmetry of the problem, the complex amplitude of the wave in the focal plane can be written as a Hankel transform of ϕ(r) of the form: (2)where F is the focal length of the telescope. Note that here the phase term τF(r) is compensated for by the lens phase function of the telescope, meaning that a simple 2D Fourier transform or Hankel transform denoted by substitutes the more complex Fresnel propagation of Eq. (1).

It is moreover convenient to calibrate the focal plane in terms of angular units θ = r/F on the sky, and the intensity in the focal plane of the telescope can be written (3)instead of just the Airy pattern (4)for the direct observation without external occulter with a perfect telescope of circular aperture of diameter 2R.

Limiting the analysis to only geometrical aspects, the occulter prohibits the planet observation for an angular radius smaller than Ω /z to the star, and the light coming from the planet passes over the occulter for an angle , a variant of the geometric inner working angle Ω /z of Kasdin et al. (2013), taking into account the telescope size. Note that all the wave and intensity functions defined in this section are illustrated in Fig. 1.

thumbnail Fig. 1

Schematic diagram of an external occulter coronagraph with the principal notations used in the paper. The observation, occulter, aperture, and focal telescope planes are shown.

Open with DEXTER

3. Procedure used for optimizing f(r) for integrated intensities

The basic idea of an external occulter coronagraph is to avoid letting the stellar light enter the telescope while leaving the light coming from the exoplanet unchanged. In a first approximation, if the angular position of the planet is larger than θ0, we can expect that the flux coming from the planet is almost unaffected by the occulter, and a fair measurement of the efficiency of the external occulter can be the normalized residual integrated intensity over the telescope aperture, for example the quantity Γ defined by (5)where Ψ(r) = | ψ(r) | 2 is the wavelength-dependent intensity and Δλ is the spectral bandwidth of the experiment, simply taken equal to λMλm here. In this analysis, the same weight is assigned to each wavelength, independently of the brightness of the source, the quantum efficiency of the detector, or the presence of acquisition filters (Soummer et al. 2010). It would be easy to perform a differently weighted average for a more realistic analysis, if necessary. Our numerical computation has been made from λm = 380 nm to λM = 750 nm and with R = 2 m.

To find an optimal shape for the external occulter, a natural approach would be to minimize Γ with respect to f. Note that this approach has been used in Wasylkiwskyj & Shiri (2011), but this is not the only way to optimize the shape of an occulter. For instance, Vanderbei et al. (2007) and Kasdin et al. (2013) proposed to minimize the upper bound c on the real and imaginary part of ψ(r),∀r ∈ [ 0,R ]. This approach will ensure that the intensity in the aperture plan will be below 2c2 and also minimize the flux. We emphasizethat while we decided here to minimize the flux, our approach can be readily adapted to the minimization of the upper bound instead of the whole flux Γ.

Because of the conservation of energy, the flux is conserved from aperture to focal planes. Nevertheless, we are in fact interested in the level of background light produced by the star in the area where we observe the planet. A quantity representative of the residual flux may be given by the residual light γ over a zone of the focal plane, such as (6)where θ1 and θ2 define the region of interest of the observation (Fig. 1), corresponding for example to the habitable zone around the star, or the possible excursion of a known planet over a star, with probably θ1 close to θ0. We therefore have γ ≤ Γ, the equality holding in the limit θ1 = 0 and θ2 = ∞.

The two measures of residual light Γ and γ are based on different points of view. If one manages to completelyshut offthe light in the aperture, then the residual light in the observation area will also shut off, but the reverse is false. In other words, minimizing Γ or γ with respect to f will lead to different optimal solutions. We believe that minimizing γ is better when the objective is to observe an exoplanet in a known area of the focal plane.

3.1. Decomposition of f(r) on a basis of functions fk(r)

As already presented in Jacquinot & Roizen-Dossier (1964) for the systematic search for apodizing properties, a common way to optimize the shape of a function is to force this function to be a weighted sum of basis functions. We here aim to optimize the weights αk such that (7)Jacquinot & Roizen-Dossier (1964) have shown that, if the expansion contains an infinite number K of terms, “the absolute optimal function according to the criterion chosen” is obtained regardless of the bases used. This result was obtained in the different context of Fraunhofer diffraction, but given the linearity of the equations, the result still holds for our purpose. In practice, however, the number K of elements of the bases must be limited and the choice of basis functions becomes important.

thumbnail Fig. 2

Left: illustration of the basis function used to represent the attenuation function f(r) with an example of the combination of three of these functions (right). For clarity of representation, the values of Ωm and Δ are not to scale of those used for the fk(r) basis functions.

Open with DEXTER

After various tests, we used trapezoidal functions for the fk(r), equal to 1 for r ≤ Ωm + (k − 1)Δ and equal to 0 for r ≥ Ωm + kΔ. The functions linearly decrease from 1 to 0 between Ωm + (k − 1)Δ and Ωm + kΔ. These function, illustrated in Fig. 2, are very similar to the binary disk for small Δ. It is the same for their Fresnel diffraction pattern, as shown in Appendix B, but this small difference is sufficient to lead to much better numerical results. Note that these basis functions lead to a transmission that is piece-wise linear as opposed to the use of binary disk that leads to piece-wise constant transmission function. For the numerical analysis, we set Δ equal to 5 cm and used K = 300 functions spreading between Ωm = 10 m to ΩM = Ω = 25 m. An example of such functions is given in Fig. 2. The reason for the existence of an innermost opaque section (i.e. t(r) = 0 or f(r) = 1 for r ≤ Ωm) is of technical origin linked with the spacecraft, as indicated for example by Vanderbei et al. (2007).

Note that expressing function f(r) as a sum of basis function leads to a more general model potentially at the cost of more complex numerical computation. For instance, one might choose to use a more elaborated basis of functions such as prolate functions. This would potentially lead to very good results because those functions are known to provide efficient apodizers, and might also be excellent occulters. Nevertheless, the physical constraints on f(r) as discussed in the following would be more difficult to express. In this particular case, it would be advantageous to compute these constraint on a finite sampling of the function, as in Vanderbei et al. (2007).

3.2. Constraints on f(r) and fk(r)

The finalattenuation function f(r) has to follow a number of physical constraints. On the one hand, some of those constraints can be easily enforced through a wise choice of the basis functions. For instance the smallest radius of apodization is set by Ωm, while the largest radius is selected by Ω = Ωm + KΔ. On the other hand, the attenuation function f(r) has to be equal to 1 for r< Ωm, and since all fk(r) functions are equal to 1 for r ≤ Ωm, we have to enforce the linear constraint . This constraint can also be expressed as a linear equality of the form 1α = 1, where 1 is a vector of ones.

Another more complicated constraint is the fact that 0 ≤ f(r) ≤ 1, r. These block constraints, thanks to the simple structure of the fk functions, can be expressed as a linear inequality constraints Aαb with Another constraint that has been used in previous works is the smoothness of the transmission function. This is typically enforced using derivatives of the function, but in our case smoothness can be readily enforced using a classical quadratic regularization term αα. Indeed, minimizing the euclidean norm of the α vector under the constraint 1α = 1 tends to promote similar values on all αk, leading to a smooth transmission function (linear apodization for an infinite regularization). This kind of regularization is equivalent to the constraint on the curvature of f described in Kasdin et al. (2013).

An additional constraint that has been enforced in previous literature is to force f to be monotonic, which in this case is a decreasing function w.r.t. the radius r. This constraint, also enforced in Kasdin et al. (2013), has an extremely simple form in our formulation, in effect it is the positivity constraint α ≥ 0. To have a grasp of the performance loss due to this additional constraint, we performed an optimization with and without this constraint in the numerical experiments.

Note that in this work, we made the design choice of using linear apodization functions as basis functions for the transmission. These choices lead to the above-mentioned physical constraints on the weights α. An equivalent formulation has been proposed in Vanderbei et al. (2007) and Kasdin et al. (2013), where the basis function are trapezoidal functions overlapping such that the final transmission is also piece-wise linear. Copi & Starkman (2000), Cash (2011) and Wasylkiwskyj & Shiri (2011), on the other hand, used offset high-order polynomial or hyper-Gaussian functions.

Below, we present the optimization problems derived from minimizing the intensity Γ in the telescope aperture and minimizing the intensity γ in a region of the focal plane. The optimal functions are denoted as fΓ(r) and fγ(r) for the optimal attenuation w.r.t. the aperture plane and w.r.t. a selected region of the focal plane, respectively.

3.3. Minimizing the integrated starlight in the telescope aperture plane: fΓ(r)

Substituting the decomposition of f(r) of Eqs. (7) into (1), we obtain (8)and the corresponding intensity becomes (9)with with ψ = [ ψ1(r),...,ψK(r) ] a rank one matrix of general coefficient , where ψk(r) is the Fresnel diffraction of a wave of unit amplitude for the basis function occulter fk(r). The flux Γ can then be expressed as (10)where is a K × K matrix integrating all the surface of the aperture and all the wavelength in (λm,λM). The optimization problem is finally defined as (11)where I is the identity matrix and μ is a positive regularization parameter. As discussed in Wasylkiwskyj & Shiri (2011), regularization is important in this case because the matrix KΓ is poorly conditioned and would lead to difficulties when solving the problem. Moreover, as discussed in the previous section, the μ parameter weight the quadratic regularization on the α vector, promoting smoothness in the transmission function. This smoothness has also been shown to be important in practice in Kasdin et al. (2013), where an equivalent smoothness constraint was inserted in the optimization problem. Also note that in practice one can easily express the smoothness constraint of Kasdin et al. (2013) in our optimization problem instead of using a regularization term. Since with the trapezoidal basis, the second-order derivatives are exactly computed using the finite difference between adjacent α components, the corresponding constraints can be obtained by adding lines to the A matrix that contain finite difference operators, and lines to b that contain the constraint parameter σ of Kasdin et al. (2013).

The optimization problem defined Eq. (11) is a constrained quadratic optimization Boyd & Vandenberghe (2004). The linear constraints forbid the use of the closed-form solution as proposed in Wasylkiwskyj & Shiri (2011), but there exist several efficient optimization strategies in the literature, among which the active set approach presented in Vanderbei & Shanno (1999).

3.4. Minimizing the starlight in a region of the telescope focal plane: fγ(r)

To express, the intensity in the focal plane, we used the same decomposition of the function f(r) as a sum of functions fk(r) as above, the complex amplitude in the focal plane can be written as a sum of elementary function φk(r), Fourier transforms of P(r)ψk(r) according to Eqs. (7) and (2). The intensity in the focal plane of Eq. (3) becomes (12)and the intergrated density (13)The resulting optimization problem is exactly the same that of Eq. (11), but with a different metric Kγ.

The optimal functions fΓ(r) and fγ(r) are different because the solutions α are different. Because of conservation of energy between aperture and focal planes, these functions become identical in the special case θ1 = 0,θ2 = ∞, where KΓ = Kγ. Our optimization in the focal plane was performed for θ1 = 0.1 arcsec and θ2 = 0.5 arcsec.

4. Results and discussion

4.1. Numerical computations

In the following, we denote ΨΓ(r) and ΦΓ(θ) the intensities obtained in the aperture and focal plane of the telescope using an external occulter of transmission 1 − fΓ(r) minimizing the integrated intensity Γ and described in Sect. 3.3. Similarly, we denote Ψγ(r) and Φγ(θ) the results obtained using an occulter transmission 1 − fγ(r), as described in Sect. 3.4.

Note that while obtaining the matrices KΓ and Kγ is a computationally intensive problem, we can use heavily parallel computing to estimate as a first step the ψk(r) and φk(θ) functions through numerical integration. This was performed using Mathematica, with K computational jobs computed by 100 parallel processes on eight core Intel processors. In the numerical experiments, both ψk(r) and φk(θ) were regularly sampled on 3001 samples from respectively r ∈ (0,R + 1) and θ ∈ (0,0.644) arcsec. Note that ψk(r) was computed outside of the telescope so that we can used it for the tolerance analysis in Sect. 4.4. Moreover, the behavior of the different optimization strategies outside of the telescope aperture is also interesting. The wavelength λ has also been sampled regularly with 21 samples in the visible light interval (380,750) nm. After computing ψk(r) and φk(θ), a numeric 2D integration can be performed to obtain the matrices KΓ and Kγ. The regularization parameter μ in the optimization problem was chosen adaptively with μ = μ0maxk,l | Kk,l | to be less dependent on the metric matrix K. μ0 was chosen by hand to allow good attenuation performance with a limited dynamic in the solution function f. For the monochromatic study, we set μ0 = 10-8 and for the more complex study on a wide bandwidth, less regularization was necessary with μ0 = 10-10 because of a more complex problem and a better conditioning of the K matrix.

Finally, the code for all the numerical experiments will be freely available on the authors website to promote reproducible research2.

thumbnail Fig. 3

Radial cut of the attenuation of the functions f(r) for aperture and focal optimization at a wavelength of 562 nm. The dashed curves denoted as 0 correspond to an optimization with α0 constraint.

Open with DEXTER

4.2. Optimization for monochromatic light

For a single monochromatic wave, here λ = 562 nm, we have represented the apodizing curves fΓ(r) and fγ(r) in Fig. 3 and the resulting intensities ΨΓ(r) and Ψγ(r) (top curves), and ΦΓ(θ) and Φγ(θ) (bottom curves) in Fig. 4. We emphasize that these curves were obtained assuming a perfect circular telescope of radius R. The f-curves are only slightly different and are really similar to a linear apodization. Moreover, in the monochromatic case, enforcing a monotonic f leads to the same results for the thickness of the trace, which will not be the case anymore for a wide spectral bandwidth, as we show below.

thumbnail Fig. 4

Top curves: intensity in the telescope aperture plane for an optimization that minimizes the integrated intensity Γ in the aperture plane (blue curve) or γ for optimizing the observation between 0.1 and 0.5 arcsec in the focal plane (red curve), for λ = 562 nm. Note the apodization profile for the γ curve and the inversion of relative position of the curves outside the region of optimization. Bottom curves: focal plane intensity obtained for the two optimizations. An airy pattern corresponding to a 10-12 fainter exoplanet in the observation area is reported in black. The apodization area defined by the smallest and largest radius of the occulter is also reported in green. Note the inversion of behavior of curves with regard to the overall occulter profile.

Open with DEXTER

The results for the Ψ and Φ-curves are quite different. The residual starlight ΨΓ(r) is very well attenuated across the telescope pupil, as found by Vanderbei et al. (2007). In contrast, Ψγ(r) presents a center-to-limb decrease which suggests an apodized structure.It is interesting to note the behavior of the curves beyond the region of optimization. There, the relative position of the curves is reversed very rapidly, the focal plane optimization becoming lower than the pupil plane optimization. This is an aspect that deserves further studies.

The intensity in the focal plane is also very interesting. While fΓ tends to minimize the intensity at all angles θ, fγ only focuses on the region of interest, and a strong residual light is observed in the blind area behind the occulter.

The resulting flux Γ for all the different optimizations is given in the upper part of Table 1. Interestingly, fΓ leads to a flux Γ in the aperture plane 105 times fainter than fγ. Nevertheless, in the focal plane the flux γ in the observation area is 50 times fainter when using fγ. This clearly shows the importance of using focal plane optimization for optimal apodization.

4.3. Optimization for a wide spectral bandwidth

Next we investigated the two optimization strategies on a wide bandwidth. For this we compute a regular sampling in the interval λ (380, 750) nm and computing the average flux across λ. Note that both a sampling of 20 and 100 values were investigated, and while we used in the experiments a sampling of 100 for a better precision, the final results were very similar. In practice the variation of the wavefront is extremely smooth w.r.t.λ and the selected sampling in sufficient.

For multiple wavelengths, λ ∈ (380, 750) nm, the results are very different. First the optimal functions are much more similar, as illustrated in Fig. 5. The functions also tend to oscillate and the monotonic variants are this time different from the free variant of the optimal functions.

The resulting intensities in the aperture and focal plan are plotted in Fig. 6. In the aperture plane, the intensity is still apodized for fΓ, but the flux is much more attenuated than in the monochromatic case. In the focal plane, the intensities retain the same tendencies, but the gain due to the focal plane optimization is lower. Interestingly, enforcing a monotonic function through α ≥ 0 leads to a clear loss in terms of attenuation for both fΓ and fγ.

thumbnail Fig. 5

Radial cut of the attenuation of the functions f(r) for aperture and focal optimization for the whole spectral bandwidth between 380 nm and 750 nm.

Open with DEXTER

Intermsof quantitative performances, Table 1 (lower part) shows that the flux in the observation area of the focal plane is still more attenuated with fγ with a gain of 1.9 in the chromatic case (1.6 with positivity constraint). While the gain in performance is smaller than in the monochromatic case, it is far from negligible when observing faint objects around a star. Interestingly, in the chromatic case, the α ≥ 0 constraint leads to a loss of 2 in attenuation (for both fΓ and fγ), proving that relaxing the constraint can lead to additional gain as long as the apodization is physically realizable.

An additional visualization comparing the performances of fΓ versus fγ in the focal plane is given in Fig. 7. The intensity Φγ(θ) is plotted as a function of ΦΓ(θ) for θ ∈ (0,0.7) arcsec. The red line allows a comparison between the methods; if a curve is below this line, it means that Φγ(θ) < ΦΓ(θ) and vice versa. We clearly see that Φγ(θ) > ΦΓ(θ) for θ<θ1 in the blind area behind the occulter and that Φγ(θ) < ΦΓ(θ)θ ∈ (θ1,θ2) in the observation area of the focal plane (illustrated in green in the figure).

thumbnail Fig. 6

Same representation as in Fig. 4, but for an optimization made for the whole spectral bandwidth between 380 nm and 750 nm. The behavior of the inversion of the curves is similar to the monochromatic case, although it is somewhat reduced by the effect of the chromatism.

Open with DEXTER

Table 1

Resulting light flux in the aperture (Γ) of focal plane (γ) for all the optimization schemes.

4.4. Tolerance analysis to positioning errors of the telescope

A study of the tolerance sensitivity to positioning errors for a 4 m circular telescope out of its on-axis position was carried out up to an offset of 1 m.

Because of this offset, the resulting pupil plane and focal plane intensities are no longer circularly symmetric, and the computation of the final focal plane image can no longer use the Hankel transformation of Eq. (2). The focal plane image must be computed taking the modulus squared of the two-dimensional Fourier transform of the complex amplitude of the wave on the telescope aperture.

In practice, the study was made using a discrete Fast Fourier Transform of an array of 1024 × 1024 points corresponding to an overall zone of 20 × 20 m, inside which the telescope of diameter 2 m was defined by 1 or 0 transmitting pixels over about 205 points in diameter. As a result, altogether 33 081 points were set equal to 1. This percentage of zero-padding is enough to obtain a satisfactorily sampled focal plane image.

Positioning errors were simulated by moving the complex amplitude obtained in the shadow of the occulter over the telescope aperture in one direction, by steps of 1 pixel or about 2 cm. The two-dimensional array was filled using an interpolation of the one-dimensional complex amplitude of the wavefront computed with Eq. (1).

The focal plane image was obtained taking the modulus squared of the Fourier transform of the wave on the telescope aperture. The effect of wavelength was taken into account afterward. The diffraction pattern increases in size with the wavelength. For that, a two-dimensional interpolation of the result was used to resample the image. Moreover, a scaling factor in intensity inversely proportional to the square of the wavelength was applied, as in Eq. (4). As a result, the integrated flux in the focal plane becomes wavelength independent and equal to that crossing through the telescope aperture. As before, for the sake of simplicity, we assumed that the received flux is constant over the spectral bandwidth.

Figures 8 and 9 present results obtained for transverse displacements of the telescope up to 1 m from its nominal position. We recall that the occulter shapes optimize the measurement when the telescope is at its nominal on-axis position, either for the integrated intensity on the telescope aperture, or for a region in the focal plane between 0.1 and 0.5 arcsec. In both cases the computation was made for the whole bandwidth of λ ∈ (380, 750) nm.

thumbnail Fig. 7

Parametric plot of the observed intensities in the focal plane: x-axis, aperture plane optimization, y-axis, focal plane optimization. The curve is drawn for the wide spectral bandwidth without positivity constraints. From top-right to bottom-left the curve follows the increase of θ, clearly showing the concentration of light behind the occulter for the focal plane optimization. The line y = x clearly shows that the γ-curves always give better results than the Γ-curves, especially in the observation area (reported in green).

Open with DEXTER

Three-dimensional representations of the normalized intensities in the telescope aperture were used to enhance differences of intensities in the results of the two optimizations, and the apodization-like pattern obtained for the focal plane optimization. Both curves were clipped to make the central parts of the images well visible. As expected, there is much less light for the aperture plane minimization when the telescope is at its nominal position. Since the optimization was calculated for a 4 m diameter telescope, there is a rise in the intensity at the edges for this large offset. This increases becomes somewhat surprisingly in favors of the focal plane optimization for a large off-set, as already discussed in Fig. 6.

In Fig. 8 focal plane intensities are also represented in color levels, with a common color scale for the telescope on-axis and off-axis of 50 cm and 1 m. The two circles of radius θ1 = 0.1 arcsec and θ2 = 0.5 arcsec define the angular area inside which the residual intensity is integrated. The focal plane optimization is found to always be better than the aperture plane optimization. From comparing these curves, it is clear that the residual diffracted light remains concentrated longer inside the first circle for the focal plane optimization, which is especially clear for the 50 cm offset images. This behavior is confirmed by Fig. 9 which represents the integrated intensity in this region with a logarithmic scale.

thumbnail Fig. 8

Top curves: pupil plane intensities for an offset of 1 m for λ ∈ (380, 750) nm. Density curves below represent the focal plane intensities (log scale) for telescope positions: on-axis (upper curves) and offsets of 50 cm (intermediate curves) and 1 m (lower curves). Left: aperture optimization, right: focal plane optimization.

Open with DEXTER

thumbnail Fig. 9

Average residual light level γ for the focal plane a) and pupil plane b) optimization, as a function of the telescope offset in centimeters for λ ∈ (380, 750) nm. The average is computed over the region from 0.1 to 0.5 arcsec delimited in Fig. 8 by two circles.

Open with DEXTER

5. Conclusion

We addressed the problem of optimizing of the shape of a circularly symmetric apodized occulter. To this end, we proposed a general expression for the attenuation function, as a weighted sum of basis functions. We expressed the optimization problem as the minimization of the flux, first in the pupil plane, as is commonly done in the literature, and then in the focal plane of a 4 m diameter telescope, for a region between 0.1 to 0.5 arcsec in which the exoplanet is supposed to be observed.

Numerical experiments show that the focal plane optimization leads to a better attenuation of the starlight than the pupil plane optimization. One interesting result is that while the focal plane optimization leads to a strong flux in the pupil plane, most of this light is concentrated in the area behind the occulter in the focal plane, which is not a problem for a perfect telescope. Finally, the robustness of our approach to lateral misalignment of the occulter was investigated and showed that the focal plane optimization also yields a better robustness.

As already indicated, the optimization study was made for a 4 m telescope assumed to be at its on-axis position. A simple procedure to make the experience less sensitive to defect position would be to perform the optimization for a larger aperture than is actually used. A better procedure would be to estimate the statistics of pointing errors and include them in the optimization procedure. Then a better tolerance to pointing errors would be obtained, probably at the expenses of a tolerable loss in the starlight rejection. Nevertheless, our analysis showed that telescope offsets of a few decimeters will not strongly reduce the efficiency of the occulter.

This conclusion was obtained for the condition that the telescope has a perfectly circular aperture and another study is required to model a telescope with central obscuration. Moreover, it would be interesting to study the effect of a a more realistic telescope with small defaults of phase and amplitude. Deviations from this perfection could lead to somewhat different results, but the exigency of quality is limited to a maximum of 0.5 arcsec, or about twenty Airy rings for a 4 m diameter telescope. Finally, with the knowledge of the exact PSF of the telescope, one can readily adapt our optimization approach. A comparison of pupil vs focal optimization for the maximum-magnitude minimization as in Kasdin et al. (2013) would also be interesting.


1

See for example http://exoplanet.eu/

Acknowledgments

We want to thank the referee Jeremy Kasdin for his very constructive comments, and in particular for his suggestion to study the effect of telescope misalignment. This work was granted access to the HPC and visualization resources of Centre de Calcul Interactif hosted by Université Nice Sophia Antipolis (http://calculs.unice.fr/en).

References

Appendix A: Quality test of the numerical computation using an analytical expression for a sharp-edged circular occulter

thumbnail Fig. A.1

From top to bottom: real, imaginary, squared modulus and unwrapped phase of ψ(r) for a circular occulter of radius Ω = 25 m (the vertical line) at z = 80 000 km for λ = 0.55 μm. The continuous line represents the Lommel series, dots are the result of a direct numerical computation of Eq. (1) using NIntegrate of Mathematica.

Open with DEXTER

Aime (2013) has shown, using an approach similar to Born & Wolf (1999), that an analytic expression using Lommel series could be derived for the Fresnel diffraction of a sharp circular occulter of the form f(r) = f(r,Ω) = Π(r/ Ω), defining here Π(r) as the function equal to 1 for | r | ≤ 1 and 0 otherwise. At the wavelength λ, the complex amplitude ψ(r,Ω) of the wave at the distance z for

  • (i)

    inside the geometrical umbra (r< Ω) and,

  • (ii)

    outside (r> Ω), can be written using two series: (A.1)At r = Ω, both series converge to the same value given by (A.2)

These series are alternative series for the real and imaginary terms, and according to the Leibniz estimate, an upper bound of the error for the sum limited to n = N terms is given by the absolute value of the N + 1 term. The convergence is easy when r/ Ω is smaller than 1, for the first sum for r< Ω inside the geometrical umbra, or the inverse outside. The main error occurs near the transition zone r = Ω.

In Fig. A.1 we compare the analytical and numerical results for an example of a circular occulter of 50 m diameter set at 80 000 km, and for a wavelength of 550 nm. The Lommel series were computed for 100 terms, and there is no difference between the two methods. Results were drawn for the real, imaginary, modulus, and unwrapped phase of the Fresnel diffraction. Note that this is just an example, in particular, the real and imaginary part of the diffraction pattern are highly sensitive to a small variation of any of the parameters. The modulus of the intensity pattern clearly shows the visible Poisson spot at the center of the pattern. The unwrapped phase presents a strong phase variation in the center of the umbra and a quieter structure outside the geometric umbra, that is for a region utilized by the telescope for the planet observation. The consistency between the two analytic and numerical results is recognized here as a proof of the quality of the numerical calculation.

Appendix B: Numerical calculation of linearly apodized basis functions fk(r)

thumbnail Fig. B.1

Differences of moduli of Fresnel diffractions for raw and linearly apodized functions for Ω values of 10 m (thick curve) and 10 m (thin curve).

Open with DEXTER

thumbnail Fig. B.2

Representation for r between 0 and 2 m of the moduli of ψk(r) of linearly apodized occulters of radii 10 m and 10.05 m a). The two curves are almost over imposed, and curve b) represents 100 times their differences.

Open with DEXTER

As indicated in the body of the paper, we tried several forms for the set of apodized basis functions fk(r), and decided to use the simple trapezoidal functions described in the body of the paper because they affect the final solution the least. We did notsucceed in deriving an analytic form for the Fresnel diffraction of these functions, but comforted by the results obtained on the Fresnel diffraction of Π(r/ Ω), we used the direct numerical calculation of Mathematica.

We recall that these functions are linearly apodized disks of radii varying from 10 m to 25 m in steps of 5 cm, the apodization applies for 5 cm at the edge. The resulting ψk(r) are very similar to those obtained for the functions Π(r/ Ω). Differences between moduli are lower than 1%, as shown in Fig. B.1 for two extreme cases, but essential to obtain proper results on diffraction patterns. In Fig. B.2 we focus on the central part of the diffraction pattern, for r ≤ 2 m. Curves are drawn for a raw and linearly apodized occulter of 10 m and 10.05 m.

All Tables

Table 1

Resulting light flux in the aperture (Γ) of focal plane (γ) for all the optimization schemes.

All Figures

thumbnail Fig. 1

Schematic diagram of an external occulter coronagraph with the principal notations used in the paper. The observation, occulter, aperture, and focal telescope planes are shown.

Open with DEXTER
In the text
thumbnail Fig. 2

Left: illustration of the basis function used to represent the attenuation function f(r) with an example of the combination of three of these functions (right). For clarity of representation, the values of Ωm and Δ are not to scale of those used for the fk(r) basis functions.

Open with DEXTER
In the text
thumbnail Fig. 3

Radial cut of the attenuation of the functions f(r) for aperture and focal optimization at a wavelength of 562 nm. The dashed curves denoted as 0 correspond to an optimization with α0 constraint.

Open with DEXTER
In the text
thumbnail Fig. 4

Top curves: intensity in the telescope aperture plane for an optimization that minimizes the integrated intensity Γ in the aperture plane (blue curve) or γ for optimizing the observation between 0.1 and 0.5 arcsec in the focal plane (red curve), for λ = 562 nm. Note the apodization profile for the γ curve and the inversion of relative position of the curves outside the region of optimization. Bottom curves: focal plane intensity obtained for the two optimizations. An airy pattern corresponding to a 10-12 fainter exoplanet in the observation area is reported in black. The apodization area defined by the smallest and largest radius of the occulter is also reported in green. Note the inversion of behavior of curves with regard to the overall occulter profile.

Open with DEXTER
In the text
thumbnail Fig. 5

Radial cut of the attenuation of the functions f(r) for aperture and focal optimization for the whole spectral bandwidth between 380 nm and 750 nm.

Open with DEXTER
In the text
thumbnail Fig. 6

Same representation as in Fig. 4, but for an optimization made for the whole spectral bandwidth between 380 nm and 750 nm. The behavior of the inversion of the curves is similar to the monochromatic case, although it is somewhat reduced by the effect of the chromatism.

Open with DEXTER
In the text
thumbnail Fig. 7

Parametric plot of the observed intensities in the focal plane: x-axis, aperture plane optimization, y-axis, focal plane optimization. The curve is drawn for the wide spectral bandwidth without positivity constraints. From top-right to bottom-left the curve follows the increase of θ, clearly showing the concentration of light behind the occulter for the focal plane optimization. The line y = x clearly shows that the γ-curves always give better results than the Γ-curves, especially in the observation area (reported in green).

Open with DEXTER
In the text
thumbnail Fig. 8

Top curves: pupil plane intensities for an offset of 1 m for λ ∈ (380, 750) nm. Density curves below represent the focal plane intensities (log scale) for telescope positions: on-axis (upper curves) and offsets of 50 cm (intermediate curves) and 1 m (lower curves). Left: aperture optimization, right: focal plane optimization.

Open with DEXTER
In the text
thumbnail Fig. 9

Average residual light level γ for the focal plane a) and pupil plane b) optimization, as a function of the telescope offset in centimeters for λ ∈ (380, 750) nm. The average is computed over the region from 0.1 to 0.5 arcsec delimited in Fig. 8 by two circles.

Open with DEXTER
In the text
thumbnail Fig. A.1

From top to bottom: real, imaginary, squared modulus and unwrapped phase of ψ(r) for a circular occulter of radius Ω = 25 m (the vertical line) at z = 80 000 km for λ = 0.55 μm. The continuous line represents the Lommel series, dots are the result of a direct numerical computation of Eq. (1) using NIntegrate of Mathematica.

Open with DEXTER
In the text
thumbnail Fig. B.1

Differences of moduli of Fresnel diffractions for raw and linearly apodized functions for Ω values of 10 m (thick curve) and 10 m (thin curve).

Open with DEXTER
In the text
thumbnail Fig. B.2

Representation for r between 0 and 2 m of the moduli of ψk(r) of linearly apodized occulters of radii 10 m and 10.05 m a). The two curves are almost over imposed, and curve b) represents 100 times their differences.

Open with DEXTER
In the text

Current usage metrics show cumulative count of Article Views (full-text article views including HTML views, PDF and ePub downloads, according to the available data) and Abstracts Views on Vision4Press platform.

Data correspond to usage on the plateform after 2015. The current usage metrics is available 48-96 hours after online publication and is updated daily on week days.

Initial download of the metrics may take a while.