OSA's Digital Library

Optics Express

Optics Express

  • Editor: C. Martijn de Sterke
  • Vol. 20, Iss. 14 — Jul. 2, 2012
  • pp: 15196–15208
« Show journal navigation

Boundary diffraction wave integrals for diffraction modeling of external occulters

Eric Cady  »View Author Affiliations


Optics Express, Vol. 20, Issue 14, pp. 15196-15208 (2012)
http://dx.doi.org/10.1364/OE.20.015196


View Full Text Article

Acrobat PDF (884 KB)





Browse Journals / Lookup Meetings

Browse by Journal and Year


   


Lookup Conference Papers

Close Browse Journals / Lookup Meetings

Article Tools

Share
Citations

Abstract

An occulter is a large diffracting screen which may be flown in conjunction with a telescope to image extrasolar planets. The edge is shaped to minimize the diffracted light in a region beyond the occulter, and a telescope may be placed in this dark shadow to view an extrasolar system with the starlight removed. Errors in position, orientation, and shape of the occulter will diffract additional light into this region, and a challenge of modeling an occulter system is to accurately and quickly model these effects. We present a fast method for the calculation of electric fields following an occulter, based on the concept of the boundary diffraction wave: the 2D structure of the occulter is reduced to a 1D edge integral which directly incorporates the occulter shape, and which can be easily adjusted to include changes in occulter position and shape, as well as the effects of sources—such as exoplanets—which arrive off-axis to the occulter. The structure of a typical implementation of the algorithm is included.

© 2012 OSA

1. Introduction

One of the major goals of the coming decades is to directly image a terrestrial exoplanet. Direct imaging allows spectral information about the planet’s atmosphere to be extracted from the spectrum of reflected light, including the potential presence of biomarkers, gases whose presence and abundance may indicate the presence of life. Detecting planets directly is a difficult task, however, for two reasons: intensity ratio (contrast) and angular resolution. An Earth-twin which orbits a Sun-twin at 1AU emits 1010 times less flux than the parent star in the visible band; if the system was located at 10 parsecs from Earth, the angular separation of the two objects would be 100mas, which for most proposed space telescopes puts the separation under 4λ/D. Specialized methods are required to remove this flux at these small working angles.

One proposed method of suppressing the starlight is an occulter, a spacecraft with a shaped edge flown in front of the telescope to block the starlight before it arrives at the telescope. The occulter size (tens of meters) and distance (tens of thousands of kilometers) are chosen so the angular extent of the occulter is smaller than some desired angle, on the order of 100 milliarcseconds, so exoplanets in orbit about the star will still be visible. The edge is shaped to control diffraction, with the form chosen to suppress the light to a factor of 1010 across the telescope aperture and over a wide spectral passband.

A major advantage of this approach, compared to an internal coronagraph, is that it significantly loosens the tolerances on the optical quality and thermal stability of the telescope optics, as well as eliminating the need for a wavefront control system for active correction. Conversely, the occulter itself must be built, deployed, and flown in formation to its own set of tolerances, and demonstrating this capability is an area of active research [1

1. N. J. Kasdin, D. N. Spergel, R. J. Vanderbei, D. Lisman, S. Shaklan, M. Thomson, P. Walkemeyer, V. Bach, E. Oakes, E. Cady, S. Martin, L. Marchen, B. Macintosh, R. E. Rudd, J. Mikula, and D. Lynch, “Advancing technology for starlight suppression via an external occulter,” Proc. SPIE 8151, 81510J (2011). [CrossRef]

]. These tolerances include limits on permanent errors in the occulter shape due to manufacturing and deployment errors; transient shape deformations from thermal fluctuations and vibration; and position and orientation changes of the occulter relative to the telescope-target axis from formation flying [2

2. S. B. Shaklan, M. C. Noecker, T. Glassman, A. S. Lo, P. J. Dumont, N. J. Kasdin, E. J. Cady, R. Vanderbei, and P. R. Lawson, “Error budgeting and tolerancing of starshades for exoplanet detection,” Proc. SPIE 7731, 77312G (2010). [CrossRef]

, 3

3. T. Glassman, A. Johnson, A. Lo, D. Dailey, H. Shelton, and J. Vogrin, “Error analysis on the NWO starshade,” Proc. SPIE 7731, 773150 (2010). [CrossRef]

].

The scale of an occulter—and the distance it must maintain with respect to its target—precludes optical testing on the ground. All of these effects must be modeled in order to build an error budget and verify that an occulter can satisfy these requirements. The capture of the large dynamic range in the resulting field—105 or greater in amplitude—drives the accuracy of the propagator used for the models, while the modeling of time-varying thermal shape deformations and closed-loop formation flying in broadband light drive the speed of the calculation. Some of the most in-depth simulated observations can take from hours to days to run, and the ability to quickly and accurately model propagation past an occulter is thus a major enabling factor in the characterization of occulter system performance. Section 2 gives an overview of the current techniques; the new technique presented in this work is described in Sections 3 and 3.1, and a suggested implementation scheme, along with a pseudocode overview of the algorithm and some computational results, are given in Section 4.

2. Current occulter modeling techniques

Occulter designs take the general form of a solid central disk surrounded by Np identical tapering structures (petals) around the edge, giving the whole structure the general appearance of a flower. (An example is shown in Fig. 1.) This shape is a result of their provenance as (0,1)-valued approximations to apodizers; the general method to design an occulter is to determine a smooth apodization profile A(r) which can provide the necessary starlight suppression, and then convert it into a binary shape with a sufficient number of petals so that performance is not undermined. The actual form of A(r) can be determined by optimization [4

4. R. J. Vanderbei, E. J. Cady, and N. J. Kasdin, “Optimal occulter design for finding extrasolar planets,” Astrophys. J. 665, 794–798 (2007). [CrossRef]

, 5

5. E. Cady, L. Pueyo, R. Soummer, and N. J. Kasdin, “Performance of hybrid occulters using apodized pupil Lyot coronagraphy,” Proc. SPIE 7010, 70101X (2008). [CrossRef]

] or selecting an existing functional form [6

6. C. J. Copi and G. D. Starkman, “The Big Occulting Steerable Satellite [BOSS],” Astrophys. J. 532, 581–592 (2000). [CrossRef]

, 7

7. W. Cash, “Detection of earth-like planets around nearby stars using a petal-shaped occulter,” Nature 442, 51–53 (2006). [CrossRef] [PubMed]

].

Fig. 1 An example of a typical occulter with 20 petals, with the associated Cartesian and polar coordinate systems.

For an occulter designed in this manner, the region Ω in which the occulter is opaque can be written directly:
Ω={(r,θ):0rR,θΘ(r)}
(1)
whereΘ(r)=n=0Np1[2πnNpπNpA(r),2πnNp+πNpA(r)],
(2)
with (r, θ) being polar coordinates in the plane of the occulter. Each of the [...,...] is a set of points which defines the outline of a single petal; ∪ denotes the union of these sets, which sets the boundary for the entire occulter. (We denote this boundary as Ω, as it will be used later.)

For a typical occulter, the occulter-telescope distance z = 5 × 107m; the wavelength λ = 5 × 10−7m; the occulter radius R = 25m; and the maximum excursions in the plane of the telescope |ξ| and |η| ≤ 3m. These values place the occulter well within the Fresnel regime; that is, the region in which it is valid to approximate the exponent in the propagation integral with the first two terms of a power series (the Fresnel approximation). This approximation is suitable when [8

8. J. W. Goodman, Introduction to Fourier Optics (McGraw-Hill, 1996).

].
z3π4λ[(xξ)2+(yη)2]max2.
(3)

Further, a standard propagation integral computes the electric field past a finite aperture; the occulter, on the other hand, is a finite obstruction with an effectively infinite aperture. Here, we can turn to Babinet’s theorem [9

9. M. Born and E. Wolf, Principles of Optics (Cambridge University Press, 1999).

], which notes that field propagated through an aperture and the field propagated past its complement sum to the field that would have resulted were no obstruction present at all.

If a plane wave with amplitude A and wavelength λ is normally-incident on the occulter, as would be the case for an occulter properly aligned with a target star, the electric field at a distance z past the occulter would be written with the Fresnel approximation and Babinet’s principle as
U(ξ,η)=Aexp(2πizλ)(11iλzΩexp{πiλz[(rcosθξ)2+(rsinθη)2]}rdrdθ)=Aexp(2πizλ)(11iλz)Ωexp{πiλz[(xξ)2+(yη)2]}dxdy),
(4)
where ξ and η are Cartesian coordinates at the downstream location. The “1–” term before the integral results from the use of Babinet’s principle. Errors in (for example) occulter shape will modify Ω, but Eq. (4) will hold regardless. Computing the integral for points in the plane of a telescope aperture is the key step in modeling the performance of an occulter.

One could attempt this integral by brute force: discretizing the shape of the occulter over a sufficiently fine grid to capture all of the necessary structure, and propagating this to the telescope aperture with direct integration or a matrix Fourier transform[10

10. R. Soummer, L. Pueyo, A. Sivaramakrishnan, and R. J. Vanderbei, “Fast computation of Lyot-style coronagraph propagation,” Opt. Express 15(24), 15935–15951, (2007). [CrossRef] [PubMed]

]. This approach is not an easy task, as the occulter will be tens of meters across, with tolerances at the level of tens of microns[2

2. S. B. Shaklan, M. C. Noecker, T. Glassman, A. S. Lo, P. J. Dumont, N. J. Kasdin, E. J. Cady, R. Vanderbei, and P. R. Lawson, “Error budgeting and tolerancing of starshades for exoplanet detection,” Proc. SPIE 7731, 77312G (2010). [CrossRef]

], requiring a grid up to 106 × 106 or more to capture the full shape of the profile, not including array padding.

A better approach is to take advantage of the structure of the occulter: in the absence of errors, it has an Np-fold symmetry in the azimuthal direction, and most errors introduce small perturbations to that. Taking advantage of the symmetry leads to three main approaches to modeling propagation past an external occulter:
  1. Integrals in r: The Bessel function expansion [4

    4. R. J. Vanderbei, E. J. Cady, and N. J. Kasdin, “Optimal occulter design for finding extrasolar planets,” Astrophys. J. 665, 794–798 (2007). [CrossRef]

    ] does the integral over θ explicitly, reducing the 2D integral to a series of single integrals in r whose contributions fall off exponentially fast near the optical axis, which makes computation extremely quick in this region. The series can be modified to incorporate some shape errors [11

    11. E. Cady, “Design, tolerancing, and experimental verification of occulters for finding extrasolar planets,” PhD thesis, Princeton University, 2010.

    ] at the expense of slower convergence. Unfortunately, it is slower at locations far from the optical axis, and requires a good deal of effort to incorporate many types of errors, such as occulter tilt, into the model. Aside from the propagation modeling, the first, radially-independent term in this series is used in optimization approaches [4

    4. R. J. Vanderbei, E. J. Cady, and N. J. Kasdin, “Optimal occulter design for finding extrasolar planets,” Astrophys. J. 665, 794–798 (2007). [CrossRef]

    , 5

    5. E. Cady, L. Pueyo, R. Soummer, and N. J. Kasdin, “Performance of hybrid occulters using apodized pupil Lyot coronagraphy,” Proc. SPIE 7010, 70101X (2008). [CrossRef]

    ] to design occulter shapes.
  2. Integrals in θ: The Dubra-Ferrari integral [12

    12. A. Dubra and J. A. Ferrari, “Diffracted field by an arbitrary aperture,” Am. J. of Phys. 67, 87–92 (1999). [CrossRef]

    ] starts one step back from the Fresnel form of Eq. (4), at the Raleigh-Sommerfeld diffraction integral, and evaluates the integral over r to produce a pair of single integrals in angular coordinates, one of which holds for points inside the geometric extent of the starshade, and one which holds outside. This approach can incorporate most errors straightforwardly and spends the same amount of time to determine the field any point in the telescope aperture plane. Because of the manner in which the integrals are segmented, however, the integrands are more difficult to determine correctly for points physically outside the extent of the starshade.
  3. Perturbations to a nominal shape: The slit approach [13

    13. P. Dumont, S. Shaklan, E. Cady, J. Kasdin, and R. Vanderbei, “Analysis of external occulters in the presence of defects,” Proc. SPIE 7440, 744008 (2009). [CrossRef]

    ] assumes the electric field for an occulter with no errors has been calculated by one of the above methods, and approximates the difference between the shape with and without errors as a series of long, thin boxes whose Fresnel transforms can be calculated quickly and exactly. This method is very fast but necessarily approximate.

In this paper, we propose a fourth method—the boundary-diffraction-wave approach—which seems to out-perform the three extant methods for general purpose use. It maintains the constant-time property of the Dubra-Ferrari integral with a 2–3× improvement in computational speed; see Sec. 4.1 for further discussion. It incorporates all errors on the shape of the occulter natively, and can be written straightforwardly to include position and orientation errors, as well as model waves from off-axis sources, such as an exoplanetary system or the finite diameter of the star. It can also take advantage of the structure of the integral to evaluate the field at multiple wavelengths with little additional computation. The approach of defining the occulter solely in terms of its boundary is a good fit for tolerancing the sorts of errors that are expected for an occulter. Manufacturing errors can be applied by specifically deforming desired edge regions; errors in petal placement can be included by translating and rotating only the boundary points associated with that petal.

3. The boundary diffraction wave formulation

Most diffraction approaches to propagation through a finite aperture find their basis in the Huygens-Fresnel principle, which expresses the electric field at a point following the aperture as the sum over spherical waves emitting from every point within the aperture. However, a different representation was introduced by Maggi and Rubinowicz which splits the field into two parts, a part which propagates according to geometrical optics and a part which incorporates the effect of diffraction from the boundary. This boundary-diffraction-wave (BDW) representation was codified by the work of Miyamoto and Wolf [14

14. K. Miyamoto and E. Wolf, “Generalization of the Maggi-Rubinowicz theory of the boundary diffraction wave part I,” J. Opt. Soc. Am 52, 615–622 (1962). [CrossRef]

, 15

15. K. Miyamoto and E. Wolf, “Generalization of the Maggi-Rubinowicz theory of the boundary diffraction wave part II,” J. Opt. Soc. Am 52, 626–636 (1962). [CrossRef]

], who showed that the boundary integral could be represented for any incident field entirely by a line integral around a vector potential. We will review their analysis briefly here.

Miyamoto and Wolf [14

14. K. Miyamoto and E. Wolf, “Generalization of the Maggi-Rubinowicz theory of the boundary diffraction wave part I,” J. Opt. Soc. Am 52, 615–622 (1962). [CrossRef]

] begin by writing the Kirchhoff formulation (see e.g. [9

9. M. Born and E. Wolf, Principles of Optics (Cambridge University Press, 1999).

]) for the diffraction integral from the aperture. Consider a volume of space, bounded by a surface S; the field at any point P within the volume may be expressed as an integral over all surface points Q on S:
U(P)=14πS(U(Q)Qexp(iks)sexp(iks)sQU(Q))ndS14πSV(Q,P)ndS.
(5)
Here we let s be the distance between P and Q, n be a vector normal to S at point Q, and ∇Q is the gradient evaluated at point Q. In the case of the occulter, the volume is the space z ≥ 0, and S stretches over the z = 0 plane, as well as having a component at infinity which is assumed to vanish. (Showing that this term vanishes requires some minor additional assumptions; see Sec. 8.3.2 of [9

9. M. Born and E. Wolf, Principles of Optics (Cambridge University Press, 1999).

]).

It can then be shown [14

14. K. Miyamoto and E. Wolf, “Generalization of the Maggi-Rubinowicz theory of the boundary diffraction wave part I,” J. Opt. Soc. Am 52, 615–622 (1962). [CrossRef]

] that for V of this form, there exists a vector potential W, such that V = ∇ × W, and for plane waves incident on an aperture, that associated vector potential is
W(Q,p)=14πAexp(ikpr)exp(iks)ss^×p1+s^p,
(6)
U(P)=S×W(Q,P)ndS.
(7)
Here r is the vector from the origin to P; Q is a point within the aperture on S, and r′ is the vector from the origin to Q. k = 2π/λ is the wavenumber and λ is the wavelength under consideration. The various s-variables become
srr,
(8)
s=s,
(9)
s^=ss.
(10)
and we assume the form of the plane wave in the vector direction p to be
U0(P)=Aexp(ikpr),
(11)
with amplitude A. A coordinate diagram is shown in Fig. 2.

Fig. 2 A diagram of the points, coordinate systems, and vectors used in this paper. The left grid is in the plane of the occulter, and the right grid in the plane of the telescope aperture.

Next, we consider for the moment the case of an opaque screen with an opening over the region Ω, a subregion of S with boundary Ω. We recall here that S is the entire bounding surface of the half-plane z ≥ 0, including the z = 0 plane; Ω is the section of it on which the occulter lies. (We are not considering the occulter case at the moment, but its complement; we will return to the occulter shortly.) In this case, the electric field from Eq. (7) across the aperture is:
Uap(P)=Ω×W(Q,P)ndΩ,
(12)
Our approach is to reduce this integral to a line integral, using Stokes’ Theorem. To do this, we consider singularities of the vector potential.

In general, W will have singularities somewhere on S for a given P, as otherwise the field in the half-space past the aperture will be zero. In the specific case of the plane wave, the lone singularity occurs at 1 + ŝ· p = 0, which occurs only when ŝ is parallel to the propagation direction of the plane wave; thus, the singularity will fall in the aperture only for points P which fall in the beam as dictated by geometric optics. (Fig. 3 shows an example of this.) Miyamoto and Wolf [15

15. K. Miyamoto and E. Wolf, “Generalization of the Maggi-Rubinowicz theory of the boundary diffraction wave part II,” J. Opt. Soc. Am 52, 626–636 (1962). [CrossRef]

] show that when Stokes’ Theorem is applied, we get two distinct cases for the field at point P, depending on the singularity location:
Uap(P)=Aexp(ikpr)+U(B)(P)forpointsinsidetheextentoftheaperture=U(B)(P)otherwise,
(13)
where U(B)(P) is the counterclockwise line integral about the edge of Ω:
U(B)(P)=ΩW(Q,P)d.
(14)
Here, is a unit vector in the direction tangent to the edge at any point, and d is a differential element of the boundary.

Fig. 3 A diagram showing the projected extent of the aperture given by geometric optics. A plane wave incident on the aperture in plane A with vector direction p produces a shifted aperture in plane B. Two angles which parameterize the direction of the plane wave, ψ1 and ψ2, are shown as well.

This formulation gives the field following an aperture, but in our case we have the opposite: the occulter is opaque, and the surroundings are open. Here we return to Babinet’s theorem. As we started with the plane wave in Eq. (11), this implies our occulter field U(P) satisfies
U0(P)=U(P)+Uap(P),
(15)
or
U(P)=ΩW(Q,P)dforpointsinthegeometricshadowoftheocculter,=AeikprΩW(Q,P)dotherwise.
(16)
The part of Eq. (16) which depends on geometric optics is termed U(G)(P):
U(G)(P)=0forpointsinsidetheextentoftheaperture=U0(P)otherwise.
(17)
This equation is the general form of the occulter field, in the boundary-diffraction-wave formulation.

3.1. Vector potentials for occulter propagation

From this suitably general representation, we can make appropriate substitutions for the specific case of an occulter; these will simplify computation of electric fields later. We define the vector from the origin to P as r = (ξ, η, z); r′, the vector from the origin to Q, is defined as r′ = (x, y, 0). (See again Fig. 2.) We can write the s-variables explicitly as
srr=(xξ,yη,z),
(18)
s=s=[z2+(xξ)2+(yη)2]1/2,
(19)
s^=ss.
(20)

We consider first the case where the plane wave is normally-incident on the occulter, as starlight would be when the system is correctly aligned. Given this, our vector terms can be written as
p=(0,0,1),
(21)
pr=0,
(22)
s^×p=1s(yη,(xξ),0),
(23)
1+s^p=1zs,
(24)
and the potential in Eq. (6) becomes
exp(iks)=exp(ikz)exp[ik(sz)],
(25)
W(Q,P)=14πAexp(ikz)exp[ik(sz)]szs^×p.
(26)

More generally, an off-axis source at an angle of ψ1 off-axis and ψ2 from the x-axis will have
p=(sinψ1cosψ2,sinψ1sinψ2,cosψ1),
(27)
pr=xsinψ1cosψ2+ysinψ1sinψ2.
(28)
Before the expression for the vector potential becomes too complicated, we will define three intermediate variables (f, g, and h) which we can substitute into Eq. (19); these choices will greatly simplify subsequent computations.
f=zsinψ1+(xξ)cosψ2cosψ1+(yη)sinψ2cosψ1,
(29)
g=(xξ)sinψ2+(yη)cosψ2,
(30)
h=zcosψ1(xξ)cosψ2sinψ1(yη)sinψ2sinψ1sp.
(31)
Substituting these variables gives:
s=[f2+g2+h2]1/2,
(32)
1+s^p=1hs,
(33)
s^×p=1s(fsinψ2+gcosψ2cosψ1,fcosψ2+gsinψ2cosψ1,gsinψ1)
(34)
exp(ikh)=exp[ikzcosψ1ik(xξ)cosψ2sinψ1ik(yη)sinψ2sinψ1]=exp(ikpr)exp(ikzcosψ1)exp[ik(ξcosψ2+ηsinψ2)sinψ1],
(35)
W(Q,P)=14πAexp(ikzcosψ1)exp[ik(ξcosψ2+ηsinψ2)sinψ1]exp(ik(sh))shs^×p,
(36)
a similar relation to Eq. (26). This equation is still exact; for ψ1 ≪ 1, the usual case for objects in an exoplanetary system, we can note that
1s(sh)=1f2+g2(1+1[1+f2+g2h2]1/2)2f2+g2,
(37)
shf2+g22h,
(38)
and so the potential becomes
W(Q,P)12πAexp(ikzcosψ1)exp[ik(ξcosψ2+ηsinψ2)sinψ1]exp(ikf2+g22h)f2+g2v(f,g),
(39)
v(f,g)=(fsinψ2+gcosψ2cosψ1,fcosψ2+gsinψ2cosψ1,gsinψ1)=s×p.
(40)
This form proves to be particularly useful, as (f2 + g2) and v · can be calculated without loss of precision, as neither will be of order z for the small ψ1 case. A typical occulter might have z = 5×107m, R = 25m, and |ξ| and |η| ≤ 3m, and k ≈ 107m−1; a representative exoplanet target for this occulter might have ψ1 = 5 × 10−7 rad. We note that exponential terms of order kz have been separated from smaller terms; these should be calculated independently, as kz ∼ 1014, and combining terms will lose precision in evaluating the exponent.

4. Efficient implementation

Evaluating the field at each downstream point requires the determination of two parts: the geometric field and the boundary field. We assume, to begin, that we are given a set of N points representing the edge of the occulter which form a simple closed curve, as well as a list of M downstream points and a list of L wavelengths at which the field is to be determined. We also assume that the edge is sufficiently well-sampled that the section of the edge between explicitly-specified points is well-modeled by a linear segment, as this allows us to use the midpoint rule for numerically approximating the integral. Given the slowly-changing shape over the majority of the petal, this is a good assumption; small regions such as petal tips may be specified more finely.

The geometric field is still a plane wave of form Eq. (11) outside the aperture, and is 0 everywhere inside; finding the geometric field thus becomes a problem of checking whether a downstream point falls behind it or not. We do this by treating the edge as a polygon and finding the winding number of the downstream point (ξ0,η0); this can be done efficiently with an O(N) routine such as polywind in Numerical Recipes [16

16. W. H. Press, S. A. Teukolsky, W. T. Vetterling, and B. P. Flannery, Numerical Recipes. The Art of Scientific Computing (Cambridge University Press, 2007).

], and holds regardless of the complexity of the occulter shape. This approach also speeds up multiband calculations, as the geometric extent of the shadow determined this way is wavelength-independent. Lateral errors in the occulter-telescope position may be included by adding a Δξ and Δη to all points at the telescope plane, and that should be done first. Note that in the case of an off-axis source, the geometric shadow is shifted laterally as well; in this case, the winding number should be calculated around (ξ0zsinψ1 cosψ2,η0zsinψ1 sinψ2). We then have:
U(G)(ξ0,η0,z)=Aexp(ikzcosψ1)exp[ik(ξ0cosψ2+η0sinψ2)sinψ1],forpointsoutsidetheshadow=0otherwise.
(41)
In some cases—for example, when the telescope aperture remains close to the center of the occulter and the perturbations to the ideal shape of the occulter are small—the winding number may be able to be determined independently. This determination should be done if possible, as even an efficient algorithm can increase runtime significantly for boundaries containing a large number of points.

The boundary field can be evaluated efficiently with midpoint-rule quadrature running around the edge. First, for each of the N line segments running counterclockwise around the occulter edge, we derive the midpoint of line segment j going from (xj, yj, 0) to (xj+1, yj+1, 0):
(xj(m),yj(m),0)=([xj+xj+1]/2,[yj+yj+1]/2,0),
(42)
and the vector pointing along the segment
(xj(),yj(),0)j=(xj+1xj,yj+1yj,0).
(43)
with (xj+1, yj+1) equal to (x1, y1) when j = N. This derivation need only be done once at the beginning, as it holds for all M points and L wavelengths. (We note that, while the midpoint rule has been used for the quadrature, this is certainly not a requirement, and higher order methods could be used; Eq. (39) will still hold. Using a higher-order method may come at the cost of runtime.)

For each point in the downstream field, we calculate f, g, and h for all N segments following Eq. (29):
fj=zsinψ1+(xj(m)ξ0)cosψ2cosψ1+(yj(m)η0)sinψ2cosψ1,
(44)
gj=(xj(m)ξ0)sinψ2+(yj(m)η0)cosψ2,
(45)
hj=zcosψ1(xj(m)ξ0)cosψ2sinψ1(yj(m)η0)sinψ2sinψ1,
(46)
and we then need only two intermediate terms:
T1j=fj2+gj2,
(47)
T2j=v(fj,gj)d=xj()(zsinψ1sinψ2+(yj(m)η0)cosψ1)yj()(zsinψ1cosψ2+(xj(m)ξ0)cosψ1),
(48)
to give us:
U(B)(ξ0,η0,z)=12πAexp(ikzcosψ1)exp[ik(ξ0cosψ2+η0sinψ2)sinψ1]×j=1Nexp(ikT1j2hj)T1jT2j.
(49)
Note that all of the steps prior to Eq. (49) are wavelength-independent; different values of k may be iterated over in the final step. It is advisable to keep terms of order kz (e.g. kzcosψ1 or kh) separate when evaluating to maintain numerical precision. A typical implementation might take the form shown in Algorithm 1. The for-loops over j in particular are well-suited to vectorization where this is supported.

Algorithm 1:. Pseudocode representation of BDW calculation

table-icon
View This Table

4.1. Computational results

In the standard formulation of the Dubra-Ferrari integral (DF), points outside the geometric shadow are calculated with a different set of integrals which require searches for nearest- and further-neighbor points along the boundary for each point at the telescope aperture. To simplify direct comparison between the Dubra-Ferrari approach and the boundary-diffraction-wave approach for this paper, we choose to select a region which lies entirely in the geometric shadow of the occulter. This selection also fixes the winding number to 1 at all points under consideration, letting us bypass an explicit calculation. This choice is not an unreasonable assumption; the telescope will generally be entirely in the geometric shadow for a properly aligned occulter-telescope system.

Both approaches have similar formulations, with the primary difference coming in the fact that the Dubra-Ferrari integral uses angular coordinates θ centered about (ξ0, η0), and so the angular distance between occulter edge points must be rederived for each point at the telescope aperture. This calculation turns out to be one of the major sources of overhead during the DF calculation.

To compare the runtimes of the Dubra-Ferrari and BDW implementations, we consider a closed occulter boundary specified by 164000 individual points, corresponding to the occulter shown in Fig. 1. We select subsets of these points, and run Matlab implementations of both algorithms on this subset. Figure 4 shows the time to implement the propagation as a function of number of boundary points; the BDW algorithm remains 2–3× faster over a wide range of polygon sizes. Figure 5 shows the intensity profile in a region around the telescope aperture, plotted on a base-10 logarithmic scale. The plots of this intensity profile created with Dubra-Ferrari and BDW are not visually distinguishable; within the telescope aperture, the maximum difference in intensity is 4.8 × 10−15.

Fig. 4 A comparison of the time to run the two propagation algorithms, as a function of number of points along the edge of the entrance aperture.
Fig. 5 Intensity at the plane of the telescope aperture following the occulter, plotted on a base-10 logarithmic scale. The telescope aperture is shown by a dashed circle.

5. Conclusion

We present a new approach to modeling fields following planet-finding occulters, based on the boundary-diffraction-wave formulation laid out by Miyamoto and Wolf. It has proved to be simple to modify the shape and recalculate with this approach to include various types of occulter errors, as well as fast to evaluate, improving on the current best method by 2–3× with no loss in precision. An efficient implementation has been laid out, and we hope it proves useful to others in the field.

Acknowledgments

This work was carried out as a part of the Imaging Performance Study funded by NASA’s Exoplanet Exploration Program; the research was performed at the Jet Propulsion Laboratory, California Institute of Technology, under a contract with the National Aeronautics and Space Administration.

References and links

1.

N. J. Kasdin, D. N. Spergel, R. J. Vanderbei, D. Lisman, S. Shaklan, M. Thomson, P. Walkemeyer, V. Bach, E. Oakes, E. Cady, S. Martin, L. Marchen, B. Macintosh, R. E. Rudd, J. Mikula, and D. Lynch, “Advancing technology for starlight suppression via an external occulter,” Proc. SPIE 8151, 81510J (2011). [CrossRef]

2.

S. B. Shaklan, M. C. Noecker, T. Glassman, A. S. Lo, P. J. Dumont, N. J. Kasdin, E. J. Cady, R. Vanderbei, and P. R. Lawson, “Error budgeting and tolerancing of starshades for exoplanet detection,” Proc. SPIE 7731, 77312G (2010). [CrossRef]

3.

T. Glassman, A. Johnson, A. Lo, D. Dailey, H. Shelton, and J. Vogrin, “Error analysis on the NWO starshade,” Proc. SPIE 7731, 773150 (2010). [CrossRef]

4.

R. J. Vanderbei, E. J. Cady, and N. J. Kasdin, “Optimal occulter design for finding extrasolar planets,” Astrophys. J. 665, 794–798 (2007). [CrossRef]

5.

E. Cady, L. Pueyo, R. Soummer, and N. J. Kasdin, “Performance of hybrid occulters using apodized pupil Lyot coronagraphy,” Proc. SPIE 7010, 70101X (2008). [CrossRef]

6.

C. J. Copi and G. D. Starkman, “The Big Occulting Steerable Satellite [BOSS],” Astrophys. J. 532, 581–592 (2000). [CrossRef]

7.

W. Cash, “Detection of earth-like planets around nearby stars using a petal-shaped occulter,” Nature 442, 51–53 (2006). [CrossRef] [PubMed]

8.

J. W. Goodman, Introduction to Fourier Optics (McGraw-Hill, 1996).

9.

M. Born and E. Wolf, Principles of Optics (Cambridge University Press, 1999).

10.

R. Soummer, L. Pueyo, A. Sivaramakrishnan, and R. J. Vanderbei, “Fast computation of Lyot-style coronagraph propagation,” Opt. Express 15(24), 15935–15951, (2007). [CrossRef] [PubMed]

11.

E. Cady, “Design, tolerancing, and experimental verification of occulters for finding extrasolar planets,” PhD thesis, Princeton University, 2010.

12.

A. Dubra and J. A. Ferrari, “Diffracted field by an arbitrary aperture,” Am. J. of Phys. 67, 87–92 (1999). [CrossRef]

13.

P. Dumont, S. Shaklan, E. Cady, J. Kasdin, and R. Vanderbei, “Analysis of external occulters in the presence of defects,” Proc. SPIE 7440, 744008 (2009). [CrossRef]

14.

K. Miyamoto and E. Wolf, “Generalization of the Maggi-Rubinowicz theory of the boundary diffraction wave part I,” J. Opt. Soc. Am 52, 615–622 (1962). [CrossRef]

15.

K. Miyamoto and E. Wolf, “Generalization of the Maggi-Rubinowicz theory of the boundary diffraction wave part II,” J. Opt. Soc. Am 52, 626–636 (1962). [CrossRef]

16.

W. H. Press, S. A. Teukolsky, W. T. Vetterling, and B. P. Flannery, Numerical Recipes. The Art of Scientific Computing (Cambridge University Press, 2007).

OCIS Codes
(050.1940) Diffraction and gratings : Diffraction
(050.1970) Diffraction and gratings : Diffractive optics
(350.6090) Other areas of optics : Space optics
(120.6085) Instrumentation, measurement, and metrology : Space instrumentation
(070.7345) Fourier optics and signal processing : Wave propagation

ToC Category:
Diffraction and Gratings

History
Original Manuscript: May 9, 2012
Revised Manuscript: June 1, 2012
Manuscript Accepted: June 5, 2012
Published: June 21, 2012

Citation
Eric Cady, "Boundary diffraction wave integrals for diffraction modeling of external occulters," Opt. Express 20, 15196-15208 (2012)
http://www.opticsinfobase.org/oe/abstract.cfm?URI=oe-20-14-15196


Sort:  Author  |  Year  |  Journal  |  Reset  

References

  1. N. J. Kasdin, D. N. Spergel, R. J. Vanderbei, D. Lisman, S. Shaklan, M. Thomson, P. Walkemeyer, V. Bach, E. Oakes, E. Cady, S. Martin, L. Marchen, B. Macintosh, R. E. Rudd, J. Mikula, and D. Lynch, “Advancing technology for starlight suppression via an external occulter,” Proc. SPIE8151, 81510J (2011). [CrossRef]
  2. S. B. Shaklan, M. C. Noecker, T. Glassman, A. S. Lo, P. J. Dumont, N. J. Kasdin, E. J. Cady, R. Vanderbei, and P. R. Lawson, “Error budgeting and tolerancing of starshades for exoplanet detection,” Proc. SPIE7731, 77312G (2010). [CrossRef]
  3. T. Glassman, A. Johnson, A. Lo, D. Dailey, H. Shelton, and J. Vogrin, “Error analysis on the NWO starshade,” Proc. SPIE7731, 773150 (2010). [CrossRef]
  4. R. J. Vanderbei, E. J. Cady, and N. J. Kasdin, “Optimal occulter design for finding extrasolar planets,” Astrophys. J.665, 794–798 (2007). [CrossRef]
  5. E. Cady, L. Pueyo, R. Soummer, and N. J. Kasdin, “Performance of hybrid occulters using apodized pupil Lyot coronagraphy,” Proc. SPIE7010, 70101X (2008). [CrossRef]
  6. C. J. Copi and G. D. Starkman, “The Big Occulting Steerable Satellite [BOSS],” Astrophys. J.532, 581–592 (2000). [CrossRef]
  7. W. Cash, “Detection of earth-like planets around nearby stars using a petal-shaped occulter,” Nature442, 51–53 (2006). [CrossRef] [PubMed]
  8. J. W. Goodman, Introduction to Fourier Optics (McGraw-Hill, 1996).
  9. M. Born and E. Wolf, Principles of Optics (Cambridge University Press, 1999).
  10. R. Soummer, L. Pueyo, A. Sivaramakrishnan, and R. J. Vanderbei, “Fast computation of Lyot-style coronagraph propagation,” Opt. Express15(24), 15935–15951, (2007). [CrossRef] [PubMed]
  11. E. Cady, “Design, tolerancing, and experimental verification of occulters for finding extrasolar planets,” PhD thesis, Princeton University, 2010.
  12. A. Dubra and J. A. Ferrari, “Diffracted field by an arbitrary aperture,” Am. J. of Phys.67, 87–92 (1999). [CrossRef]
  13. P. Dumont, S. Shaklan, E. Cady, J. Kasdin, and R. Vanderbei, “Analysis of external occulters in the presence of defects,” Proc. SPIE7440, 744008 (2009). [CrossRef]
  14. K. Miyamoto and E. Wolf, “Generalization of the Maggi-Rubinowicz theory of the boundary diffraction wave part I,” J. Opt. Soc. Am52, 615–622 (1962). [CrossRef]
  15. K. Miyamoto and E. Wolf, “Generalization of the Maggi-Rubinowicz theory of the boundary diffraction wave part II,” J. Opt. Soc. Am52, 626–636 (1962). [CrossRef]
  16. W. H. Press, S. A. Teukolsky, W. T. Vetterling, and B. P. Flannery, Numerical Recipes. The Art of Scientific Computing (Cambridge University Press, 2007).

Cited By

Alert me when this paper is cited

OSA is able to provide readers links to articles that cite this paper by participating in CrossRef's Cited-By Linking service. CrossRef includes content from more than 3000 publishers and societies. In addition to listing OSA journal articles that cite this paper, citing articles from other participating publishers will also be listed.

Figures

Fig. 1 Fig. 2 Fig. 3
 
Fig. 4 Fig. 5
 

« Previous Article  |  Next Article »

OSA is a member of CrossRef.

CrossCheck Deposited