Integration time adjusted completeness

Dean Keithly; Dmitry Savransky; Corey Spohn

doi:10.1117/1.JATIS.7.3.037002

1 September 2021 Integration time adjusted completeness

Dean Keithly, Dmitry Savransky, Corey Spohn

Journal of Astronomical Telescopes, Instruments, and Systems, Vol. 7, Issue 3, 037002 (September 2021). https://doi.org/10.1117/1.JATIS.7.3.037002

Abstract

Future, large-scale, exoplanet direct-imaging missions will be capable of discovering and characterizing Earth-like exoplanets. These mission designs can be evaluated using completeness, the fraction of planets from some population that are detectable by a telescope at an arbitrary observation time. However, the original formulation of completeness uses instrument visibility limits and ignores additional integration time and planetary motion constraints. Some of the sampled planets used to calculate completeness may transit in and out of an instrument’s geometric and photometric visibility limits while they are being observed, thereby causing the integration time agnostic calculation to overestimate completeness. We present a method for calculating completeness that accounts for the fraction of planets that leave the visibility limits of the telescope during the integration time period. We define completeness using the aggregate fraction of an orbital period during which planets are detectable, calculated using the specific times that planets enter and leave an instrument’s visibility limits and the integration time. To perform this calculation, we derive analytical methods for finding the planet-star projected separation extrema, times past periastron that these extrema occur, and times past periastron that the planet-star projected separation intersects a specific separation circle. We also provide efficient numerical methods for calculating the planet-star difference in magnitude extrema and times past periastron corresponding to specific values Δmag. Our integration time adjusted completeness shows that, for a planned star observation at 25 pc with 1-day and 5-day integration times, integration time adjusted completeness of Earth-like planets is reduced by 1% and 5% from the integration time agnostic completeness, respectively. Integration time adjusted completeness calculated in this manner also provides a computationally inexpensive method for finding dynamic completeness—the completeness change on subsequent observations.

1. Introduction

Direct imaging blind search mission schedules can be optimized¹ by maximizing completeness²—the fraction of exoplanets from an assumed planet population that are detectable by a particular instrument at an arbitrary observation time. Completeness is typically parameterized by a limiting planet-star brightness difference ( $Δ {mag}_{\lim}$ ), the inner working angle (IWA) of the instrument, and the outer working angle (OWA) of the instrument. The original Monte Carlo approach developed by Brown² involved creating a cloud of synthetic planets by sampling the underlying Keplerian orbital elements (KOE) and physical parameters of a planet population and determining the fraction of those individually simulated planets within the visible limits of the instrument. Multiplying completeness by the exoplanet occurrence rate gives the expected exoplanet yield for observing a given star. Although completeness is a good metric for predicting instrument performance, the calculation described above only captures an instant in time and does not include whether the time that a planet is within the visible region of the telescope is sufficient to actually make a detection.

The original method of calculating completeness was developed by Brown.² Our adaptation of this method in EXOSIMS³ is extensively outlined in Ref. 1 and involves randomly generating planet KOE, planet geometric albedo ( $p$ ), and planetary radius ( $R$ ) sampled from NASA’s Exoplanet Program Analysis Group Study Analysis Group 13 (SAG13)¹^,⁴^,⁵ probability density functions, generating a finely binned 2D histogram of $Δ mag$ versus $s$ (planet-star separation), and fitting a 2D spline to the histogram bins. In all current methods, completeness is calculated as the double integral of $s$ and $Δ mag$ over the joint probability density function as in Eq. (17) of Ref. 1 and Eq. (7) in Ref. 6. The limits of integration define a detectable planet as one in which $Δ mag < Δ {mag}_{\lim}$ and $d_{i} IWA < s < d_{i} OWA$ , where $d_{i}$ is the distance of the host star from the spacecraft. While the $Δ {mag}_{\lim}$ in this paper is used to describe a general upper limit of integration for calculating completeness, this $Δ {mag}_{\lim}$ can be formulated as a function of integration time as in Eq. (12) of Ref. 1 by making assumptions about a multitude of instrument parameters and external noise sources. This approach, as conventionally implemented, is a good estimator for completeness,¹ but it does not take planet motion in time into account.

Searches like the Gemini Planet Imager and Nancy Grace Roman Space Telescope (Roman) Coronagraphic Instrument (CGI)⁷ are sensitive to larger planets with larger planet-star separations and longer orbital periods.¹ They make use of Brown completeness to plan blind searches. However, future telescope concepts like the Habitable Exoplanet Observatory (HabEx)⁸ and the Large UV/Optical/IR Surveyor (LUVOIR)⁹ seek to find smaller Earth-like exoplanets, with shorter orbital periods around stars farther away. Their search criteria will result in Brown completeness based yield overestimations due to planet motion out of instrument visibility limits. Figure 1 shows a schematic of a direct imaging observation; it demonstrates how a planet can be within the photometric and geometric visibility limits of the instrument and still not be detected. The red regions of the planet’s orbit indicate where the planet is within the visible limits but not detectable because the planet will move out of the region in less time than it takes to detect the planet [case (b) in this figure]. Since there will be some Earth-like exoplanets that are counted toward the completeness score but are not actually within the visibility constraints of the instrument long enough to be directly imaged, we need a new method that only counts targets within the instrument visibility long enough to be observed or characterized.

Fig. 1

Schematic of a direct imaging observation. The line represents the projection of a planet’s orbit about its parent star (yellow) into the plane of the sky as seen by a distant observer. The arrow indicates the planet’s direction of motion. Blue portions of the orbit (a) indicate times when the planet is detectable. Red portions of the orbit (b) indicate times when the planet is within the visibility limits of the instrument, but it is not detectable due to integration time constraints. The dashed portions of the orbit (c) indicate times when the planet is too faint to be observed. The shaded gray regions represent the projected inner and OWAs of the instrument. The planet is unobservable on portions of the orbit intersecting these gray regions. Green dots indicate transitions in and out of instrument visibility limits.

Figure 2 demonstrates this phenomenon in the separation versus $Δ mag$ phase space commonly used to define instrument contrast curves and completeness. Figure 2(b) shows a short-period, fast-moving, Mars-like planet that transitions into and out of the assumed instrument’s visible limits before it can be detected. Figure 2(c) shows the effect for a Uranus-like planet that crosses through and exits the visible region in less time than it takes for it to be detected. In these cases, Brown’s completeness calculation would include simulated planets in these regions, thus overestimating the overall detection yield. Additionally, Fig. 2(d) shows that the last possible moment a Neptune-like planet could be imaged before it leaves the visible limits of the telescope. This demonstrates a portion of the orbit where the planet is within the visible limits of the telescope but is not detectable.

Fig. 2

The $Δ mag$ versus $s$ curve for Neptune (blue) and Mars (red). The instrument’s visible region (white) is bounded by 0.066 and 0.517 arc sec, and $Δ mag = 31$ at $d_{i} = 22.87$ pc (we selected instrument parameters to make a demonstration on solar system planets about a Sun-like star, but any instrument will have planets with similar effects). Any planet in the grayed region is not visible by the instrument. Dots indicate calculated intersection points between the planet and the visible region bounds. The black dashed portions of planet orbits indicate portions of the orbit where the planet is not detectable. (b) A portion of Mars’s orbit where the planet enters and exits the visible limits of the instrument in less time than the required integration time.¹ (c) A portion of Neptune’s orbit where the planet enters and exits the visible region of the instrument in less time than the required integration time. (d) Neptune’s detectable region adjusted for integration time. If the observation of Neptune in (d) is started when it is left of the blue × and in the visible region, then it will be detected; if it is started right of the blue ×, the planet will not be detected because the required integration time is larger than the time the planet will be in the visible region of the instrument.

This paper presents a method, implemented in the EXOSIMS modeling software,³ for calculating integration time adjusted completeness that allows us to investigate exactly how much the original completeness definition overestimates planet yields. EXOSIMS is an exoplanet direct imaging mission modeling software used to simulate the Roman CGI⁷ and future mission concepts such as HabEx⁸ and LUVOIR.⁹ EXOSIMS simulates populations of planets, observatories, instruments, and underlying dynamics to create full, end-to-end mission simulations. We derive distributions of potential mission yield from ensembles of these mission simulations. EXOSIMS utilizes completeness as a heuristic for target selection and observation scheduling,¹ implementing both Brown’s Monte Carlo approach to completeness² as well as the analytical completeness formulation from Ref. 6. However, EXOSIMS determines whether actual detections occur by evaluating the effective signal-to-noise ratio (SNR) of individual simulated observations on simulated planets. Keithly et al.¹ discussed the full optimization process for a single-visit blind search and showed that the ideal yield from an EXOSIMS simulation ensemble (when mission constraints are discounted) for a limited blind search mission of Roman is equivalent to the yield calculated via Brown Completeness.

Comparisons between Brown completeness-based exoplanet yield estimations and Monte Carlo of full mission simulations have been done in Ref. 10 and the Standard Definitions and Evaluation Team Final Report.¹¹ The yields calculated by completeness are higher than the average yield from a Monte Carlo of full missions simulated in EXOSIMS as noted in Table 6 of Ref. 11. These works attribute the difference in yields to the additional mission constraints captured by the full mission simulation approach, inefficiencies in the optimization algorithms, and inefficiencies in the scheduling algorithms implemented in EXOSIMS.¹² The unmentioned assumption is that Brown completeness-based yield estimates are accurate, but this paper demonstrates how the omission of planet motion and integration time reduces the resulting completeness yield estimates.

To know the amount of time that an individual planet spends within the observable limits of an instrument, we need to calculate the times in a planet’s orbit when a planet enters or exits these limits. This means that we need methods for calculating when a planet has a given planet-star separation and $Δ mag$ . Once we know all of the times when a planet enters or exits an instrument’s visibility limits, we know the fraction of time that the planet is able to be detected with that instrument. Averaging over these visibility fractions for a large number of samples evaluates completeness in a fundamentally different way from its original formulation. These same methods also allow for the evaluation of integration time adjusted completeness as well as a new method of calculating dynamic completeness.¹³

Section 2 presents the derivation of this chain of calculations and discusses practical aspects of their implementation. Section 3 provides detailed validation of the methods and presents results of the various calculations enabled by them. Finally, in Section 4, we discuss various aspects of the algorithms and results, and we lay out future applications for this methodology.

2. Methods

In this section, we present our detailed process for calculating integration time adjusted completeness. We do this by finding the time windows in which a planet is within the separation and $Δ mag$ visibility limits of an instrument, discounting each time window by the integration time. The general overview of this process is as follows.

1. Calculate locations of apparent intersections between the projection of the 3D orbit into the plane of the sky and the $s_{WA}$ circle about the star.
- (a) Parameterize the 3D orbit as an ellipse in the plane of the sky.
- (b) Formulate the planet-star separation equation and its derivative.
- (c) Solve for the magnitude and locations of planet-star separation extrema.
- (d) Identify and assign the subset of algebraic solutions that are separation extrema.
- (e) Classify the expected number of intersections between the projected ellipse and the $s_{WA}$ circle.
- (f) Solve for intersections between the projection of the orbital ellipse and the $s_{WA}$ circle.
- (g) Identify and assign the subset of algebraic solutions that form intersections.
2. Calculate locations and magnitudes of $Δ mag$ intersections.
- (h) Formulate an expression isolating the $Δ mag$ and $ν$ terms and its derivative.
- (i) Express this formulation as a polynomial and solve for the algebraic solutions.
- (j) Identify and assign the subset of algebraic solutions that are $Δ mag$ extrema.
- (k) Classify the expected number of $Δ mag$ and $Δ {mag}_{\lim}$ intersections.
- (l) Solve for $Δ mag$ and $Δ {mag}_{\lim}$ intersections.
- (m) Identify and assign the subset of algebraic solutions that form the intersections.
3. Calculate $ν$ from X and Y for each intersection and extrema.
4. Calculate $t$ from $ν$ of each intersection and extrema.
5. Combine times of $s_{WA}$ and $Δ {mag}_{\lim}$ intersections to create time windows between intersections.
6. Identify time windows in which the planet is visible or not visible.
7. Calculate integration time adjusted completeness averaging the orbital fraction of time that a planet is visible discounted by the integration time.

The general equations for $Δ mag$ and $s$ used in these derivations are

Eq. (1)

Δ mag = - 2.5 \log_{10} (p {(\frac{R}{| {\underset{̲}{r}}_{k / i} |})}^{2} Φ (β)),

and

Eq. (2)

s = ‖ {\underline{r}}_{k / i} - ({\underline{r}}_{k / i} \cdot {\underset{̲}{\hat{r}}}_{i / SC}) {\underset{̲}{\hat{r}}}_{i / SC} ‖ .

Here

R

is the planet radius,

p

is the geometric albedo of the planet, and

Φ (β)

is the planet phase function. The other variables are defined in Fig. 3, where

β

is the star-planet-observer angle (also called the phase angle),

{\underline{r}}_{k / i}

is the vector from star

i

to planet

k

, and

{\underset{̲}{\hat{r}}}_{i / SC}

is the unit vector from the spacecraft (SC) to the star (

\underset{̲}{\hat{r}}

is the unit vector of

\underline{r}

). The plane of the sky for a given observation lies in the

\underset{̲}{\hat{x}}

and

\underset{̲}{\hat{y}}

plane, where

{\underset{̲}{\hat{r}}}_{i / SC}

defines

\underset{̲}{\hat{z}}

of the target system. For our purposes, the direction of

\underset{̲}{\hat{x}}

is arbitrary, but is typically taken to be a well-defined, inertially fixed direction, such as the ICRS mean equinox or pole direction. Here and throughout this paper,

i

in a subscript refers to the

i

’th target star and

i

(not subscripted) refers to the orbit inclination.

Fig. 3

The orbital path of planet $k$ (black) in a general $X Y Z$ Cartesian coordinate system. ${\underline{r}}_{k / i}$ (straight blue arrow) describes the location of the planet $k$ (blue circle) relative to star $i$ (yellow circle). The plane of the sky is noted by the translucent red parallelogram entirely in the $\underline{\hat{x}}$ and $\underline{\hat{y}}$ plane. The dashed ellipse (red) is the projection of the planet orbit (black) onto the plane of the sky. $\underline{\hat{h}}$ is the orbit angular momentum vector. $\underline{\hat{e}}$ is the orbit eccentricity vector. $\underline{\hat{n}}$ is the line of nodes. The star $*$ is generally referred to in subscripts as $*$ and often referred to in subscripts as the $i$ ’th star. This is not to be confused with the variable $i$ , which is the inclination of the planet orbit. $ν$ is the true anomaly of the planet. $ω$ is the argument of periapsis of the planet. $Ω$ is the longitude of the ascending node of the planet.

2.1.

Projected Orbit and Separation Intersection

In this section, we derive an analytical expression for the true anomaly ( $ν$ ) of $s$ -orbit intersection points between a circle in the plane of the sky and the projection of a 3D Keplerian orbit on the plane of the sky. A practical example of a circle in the plane of the sky is the projected inner or OWA of the instrument, equal to $IWA d_{i}$ (or $OWA d_{i}$ ) for target star distance $d_{i}$ . An orbit defined by KOE takes the shape of an ellipse in the orbital plane. The perpendicular projection of a 3D elliptical orbit onto the plane of the sky is given by the $\underset{̲}{\hat{x}}$ and $\underset{̲}{\hat{y}}$ components of the 3D orbit. However, the resulting expression is not easily solvable for the true anomalies at intersection points. Instead, we can express the orbit projection in the plane of the sky as another ellipse. This simplifies the difficult intersection problem into analytically solvable subproblems and yields a relatively simple solution to the intersection between a circle in the plane of the sky and the projected orbit in the plane of the sky. Here we include the brief outline of this process, with the full procedure detailed in Appendix B.

1. Project the 3D orbital ellipse into a 2D projected ellipse and formulate the separation equation.
2. Find the expected number of $s$ intersection points.
3. Calculate $s$ intersection point coordinates.

2.1.1.

Projection of the elliptical orbit

Given the KOE of a 3D orbit and a plane to project it on, we calculate the semimajor and semiminor axes of the projected orbit ( $a_{p}$ and $b_{p}$ , respectively) as well as the angle between the projected semimajor axis and $\underset{̲}{\hat{x}}$ , which we call ( $θ$ ). We first convert the KOE into orbital radius components in the Cartesian $X Y Z$ coordinate system as in Fig. 3:

Eq. (3)

X = r (\cos (Ω) \cos (ω + ν) - \sin (Ω) \sin (ω + ν) \cos (i)) Y = r (\sin (Ω) \cos (ω + ν) + \cos (Ω) \sin (ω + ν) \cos (i)) Z = r \sin (i) \sin (ω + ν),

where

r

is the orbital radius (magnitude of

{\underline{r}}_{k / i}

) given by

Eq. (4)

r = \frac{a (1 - e^{2})}{1 + e \cos (ν)} .

Ω

is the longitude of the ascending node,

ω

is the argument of periapsis,

ν

is the true anomaly, and

e

is the eccentricity of the planet’s orbit.

Figure 4 shows a schematic view of the projection of the orbit onto the plane of the sky. There are two particularly important points associated with the projected orbit ellipse. First, $F$ is the filled focus (star location) of the orbit and retains the same coordinates in the plane of the sky. When observing a star, $F$ is the center of all working angle circles and the origin of the $X Y Z$ coordinate system. Note that in Fig. 4 the projected ellipse is located well below the original 3D ellipse for clarity, but the points $F$ indicated by the orange circle and $F$ indicated by the orange × are in fact coincident.

Fig. 4

The original 3D elliptical orbit of a planet (black line) containing points $A$ , $B$ , $C$ , and $D$ ; the gray endpoints of the purple semimajor axis $\overline{A B}$ and semiminor axis $\overline{C D}$ . The planet orbits about the star, which is located at focus $F$ (orange circle). The orange dot and orange × denoted as $F$ are the same point in space, but the red projection is shown as offset from the original orbit for clarity. The projection of the original 3D elliptical orbit onto the $X Y$ plane is given by the red ellipse containing points $A^{'}$ , $B^{'}$ , $C^{'}$ , and $D^{'}$ (blue diamonds; perpendicular projections of $A$ , $B$ , $C$ , and $D$ ). Point $O$ is the geometric center of the 3D elliptical orbit and projects to $O^{'}$ in the plane of the sky. Point $P$ (green circle) is any arbitrary point along the original 3D ellipse and maps to the semimajor axis and semiminor axis components $H$ and $K$ , respectively (green ×). $P^{'}$ (pink circle) is the perpendicular projection of $P$ , and $H^{'}$ and $K^{'}$ (pink ×) are projections of $H$ and $K$ , respectively. The components of these perpendicular projections preserve the ratios of their values to the semimajor and semiminor axes, which, given the equation for an ellipse, can be used to prove the projection of an ellipse is itself an ellipse. The blue lines $\overline{A^{'} B^{'}}$ and $\overline{C^{'} D^{'}}$ form conjugate diameters of the red ellipse and are the projection of the semimajor axis and semiminor axis of the original 3D ellipse onto the plane of the sky. These conjugate diameters can then be used to find the semimajor axis and semiminor axis of the red projected ellipse $\overline{I R}$ and $\overline{S T}$ , respectively. ${\underline{x}}_{d r}$ , ${\underline{y}}_{d r}$ , ${\underline{z}}_{d r}$ are the components of the derotated reference frame ( $d r$ ) as defined in the text.

Second, $O$ is the geometric center of the 3D orbit and is found most efficiently by averaging the $X Y Z$ locations of the planet at apoapsis and periapsis:

Eq. (5)

\bar{F O} = \frac{1}{2} ({\underline{r}}_{k / i} (ν = 0) + {\underline{r}}_{k / i} (ν = π)),

where

\overline{F O}

is the line segment (equivalently Euclidean vector) from

F

to

O

and

{\underline{r}}_{k / i} (ν)

is the evaluation of Eq. (3) for the given value of

ν

. As

\underset{̲}{\hat{x}}

and

\underset{̲}{\hat{y}}

define the plane of the sky, the

X Y Z

coordinates of the geometric center of the projected ellipse (

O^{'}

) are given by

Eq. (6)

O^{'} = ⟨ \bar{F O} \cdot \underset{̲}{\hat{x}}, \bar{F O} \cdot \underset{̲}{\hat{y}}, 0 ⟩ .

Appendix J provides a proof that any generic 3D ellipse projects to a 2D ellipse, as shown graphically in Fig. 4. The projection of the semimajor axis and semiminor axis from the 3D orbit onto the 2D plane of the sky form conjugate diameters of the projected ellipse. Any two diameters of an ellipse are conjugate diameters if and only if the tangent line to the ellipse at an endpoint of one diameter is parallel to the other diameter.¹⁴ Each pair of conjugate diameters of an ellipse has a corresponding tangent parallelogram, sometimes called a bounding parallelogram. In Fig. 4, the pair of blue lines ( $\overline{A^{'} B^{'}}$ and $\overline{C^{'} D^{'}}$ ) are specific examples of conjugate diameters. The pair of purple lines ( $\bar{I R}$ and $\overline{S T}$ ) are also examples of conjugate diameters, and they form the semimajor and semiminor axes of the projected ellipse.

The semimajor and semiminor axes of the projected ellipse may be found from any conjugate diameters. We take the conjugate diameters $\overline{A^{'} B^{'}}$ and $\overline{D^{'} C^{'}}$ and draw line $\overline{Q Q^{'}}$ through $B^{'}$ , perpendicular to $\overline{D^{'} C^{'}}$ as shown by the gray line in Fig. 5. Points $Q$ and $Q^{'}$ are chosen such that $\bar{B^{'} Q} = \bar{B^{'} Q^{'}} = \bar{O^{'} D^{'}}$ . The principal axes $\overline{I R}$ and $\overline{S T}$ lie on the bisectors of the angles formed by lines $\overline{O^{'} Q}$ and $\overline{O^{'} Q^{'}}$ . From this construction, we calculate the projected ellipse semimajor axis, semiminor axis, and angular offset of the semimajor axis from $\underset{̲}{\hat{x}}$ .¹⁶ $\overline{I R}$ and $\overline{S T}$ are given by

Eq. (7)

\bar{I R} = \bar{O^{'} Q^{'}} + \bar{O^{'} Q},

Eq. (8)

\bar{T S} = \bar{O^{'} Q^{'}} - \bar{O^{'} Q} .

Defining the angle between

\bar{O^{'} B^{'}}

and

\bar{O^{'} D^{'}}

as

ϕ

, the cosine rule applied to triangles

O^{'} B^{'} Q

and

O^{'} B^{'} Q^{'}

yielding

Eq. (9)

{| \bar{O^{'} Q} |}^{2} = {| \bar{O^{'} B^{'}} |}^{2} + {| \bar{O^{'} D^{'}} |}^{2} - 2 | \bar{O^{'} B^{'}} | | \bar{O^{'} D^{'}} | \sin ϕ,

Eq. (10)

{| \bar{O^{'} Q^{'}} |}^{2} = {| \bar{O^{'} B^{'}} |}^{2} + {| \bar{O^{'} D^{'}} |}^{2} + 2 | \bar{O^{'} B^{'}} | | \bar{O^{'} D^{'}} | \sin ϕ,

where

| \overline{O^{'} B^{'}} |

denotes the length of

\bar{O^{'} B^{'}}

. Inserting these into Eq. (7), we obtain

Eq. (11)

\bar{O^{'} R} \cdot \bar{O^{'} S} = | \bar{O^{'} B^{'}} | | \bar{O^{'} D^{'}} | \sin ϕ,

Eq. (12)

{| \bar{O^{'} R} |}^{2} + {| \bar{O^{'} S} |}^{2} = {| \bar{O^{'} B^{'}} |}^{2} + {| \bar{O^{'} D^{'}} |}^{2} .

| \overline{I R} |

must be twice the semimajor axis of the projected ellipse (

a_{p}

), so

Eq. (13)

a_{p} = \frac{| \bar{I R} |}{2} = | \bar{O^{'} R} |,

and

| \overline{T S} |

must be twice the semiminor axis of the projected ellipse, so

Eq. (14)

b_{p} = \frac{| \bar{T S} |}{2} = | \bar{O^{'} S} | .

The angle of the semimajor axis of the projected ellipse from

\underset{̲}{\hat{x}}

is calculated using the average of the angles between

\overline{O^{'} Q}

and

\underset{̲}{\hat{x}}

and

\overline{O^{'} Q^{'}}

and

\underset{̲}{\hat{x}}

:

Eq. (15)

θ = \frac{1}{2} (\tan^{- 1} (\frac{\bar{O^{'} Q} \cdot \underset{̲}{\hat{y}}}{\bar{O^{'} Q} \cdot \underset{̲}{\hat{x}}}) + \tan^{- 1} (\frac{\bar{O^{'} Q^{'}} \cdot \underset{̲}{\hat{y}}}{\bar{O^{'} Q^{'}} \cdot \underset{̲}{\hat{x}}})) .

Appendix E provides the full expressions for

a_{p}

,

b_{p}

, and

θ

via full expansions of Eqs. (13)–(15), respectively.

Fig. 5

The projected ellipse (red) has semimajor axis $\overline{I R}$ and semiminor axis $\overline{T S}$ (purple). Both are calculable from the conjugate diameters that are the projections of the semimajor and semiminor axes of the original 3D ellipse ( $\overline{A^{'} B^{'}}$ , $\overline{C^{'} D^{'}}$ , blue). The line $\overline{Q Q^{'}}$ is drawn such that it is perpendicular to the smaller conjugate diameter ( $\overline{C^{'} D^{'}}$ ) and is bisected by $B^{'}$ : $| \overline{B^{'} Q} | = | \overline{B^{'} Q^{'}} | = | \overline{O^{'} C^{'}} |$ . The semimajor axis of the projected ellipse is the angular bisector of $\overline{O^{'} Q}$ and $\overline{O^{'} Q^{'}}$ .¹⁵ Finally, the semiminor axis $\overline{T S}$ of the projected ellipse is perpendicular to the semimajor axis.

Now that we know all of the parameters necessary to describe the projected ellipse, we can standardize this ellipse into a simpler form to simplify subsequent calculations. We define a new frame ( $d r$ as in Fig. 4) as the derotation and geometric centering of the projected orbit such that the semimajor axis of the projected ellipse ( $\bar{O^{'} I}$ ) is aligned with ${\underset{̲}{\hat{x}}}_{d r}$ , the semiminor axis of the projected ellipse ( $\overline{O^{'} S}$ ) is aligned with ${\underset{̲}{\hat{y}}}_{d r}$ , and $O^{'}$ is the origin of the $d r$ coordinate system. The $d r$ coordinates of the star location $F$ ( $x_{*}$ , $y_{*}$ ) are given by a simple rotation of the projection of $\overline{O^{'} F}$ onto the $\underset{̲}{x}, \underset{̲}{y}$ plane by angle $θ$ :

Eq. (16)

{[\begin{matrix} x_{*} \\ y_{*} \end{matrix}]}_{dr} = - [\begin{matrix} \cos (θ) & \sin (θ) \\ - \sin (θ) & \cos (θ) \end{matrix}] {[\begin{matrix} \bar{O^{'} F} \cdot \underset{̲}{\hat{x}} \\ \bar{O^{'} F} \cdot \underset{̲}{\hat{y}} \end{matrix}]}_{X Y Z} .

2.1.2.

Global and local extrema of planet-star separation

Before solving for the true anomalies where the orbit’s projected separation $s$ is equal to $s_{WA}$ (general working angle separation $WA d_{i}$ ), we first need to know how many of these $s_{WA}$ -orbit intersections we are looking for. We find the expected number of solutions by finding the $s$ extrema throughout the orbit. We do so by solving for the roots of the derivative of the projected planet-star separation.

We start with the general equation for an ellipse:

Eq. (17)

{(\frac{x_{e}}{a_{p}})}^{2} + {(\frac{y_{e}}{b_{p}})}^{2} = 1,

where

x_{e}

and

y_{e}

are the coordinates of any point on the ellipse and

a_{p}

,

b_{p}

are the semimajor and semiminor axes, respectively. We rewrite this in terms of

x_{e}

, giving

Eq. (18)

y_{e} = b_{p} \sqrt{1 - \frac{x_{e}^{2}}{a_{p}^{2}}} .

The projected separation is given by

Eq. (19)

s^{2} = {(- x_{*} + x_{e})}^{2} + {(- y_{*} + y_{e})}^{2} .

Taking the derivative of Eq. (19) with respect to

x_{e}

, substituting Eq. (18) (and its derivative), and setting it equal to zero, we have

Eq. (20)

0 = \frac{δ s^{2}}{δ x_{e}} = - 2 x_{*} + 2 x_{e} + \frac{2 b_{p} x_{e} y_{*}}{a_{p} \sqrt{a_{p}^{2} - x_{e}^{2}}} - \frac{2 b_{p}^{2} x_{e}}{a_{p}^{2}} .

Isolating the square root term to one side and squaring both sides of the equation gives

Eq. (21)

{(- 2 x_{*} + 2 x_{e} - \frac{2 b_{p}^{2} x_{e}}{a_{p}^{2}})}^{2} = {(\frac{2 b_{p} x_{e} y_{*}}{a_{p} \sqrt{a_{p}^{2} - x_{e}^{2}}})}^{2},

which is expanded with coefficients of

x_{e}

, to get the polynomial expression

Eq. (22)

0 = x_{e}^{4} + \frac{- 8 a_{p}^{2} x_{*} + 8 b_{p}^{2} x_{*}}{(4 a_{p}^{4} - 8 a_{p}^{2} b_{p}^{2} + 4 b_{p}^{4}) / a_{p}^{2}} x_{e}^{3} + \frac{- 4 a_{p}^{4} + 8 a_{p}^{2} b_{p}^{2} + 4 a_{p}^{2} x_{*}^{2} - 4 b_{p}^{4} + 4 b^{2} y_{*}^{2}}{(4 a_{p}^{4} - 8 a_{p}^{2} b_{p}^{2} + 4 b_{p}^{4}) / a_{p}^{2}} x_{e}^{2} + \frac{8 a_{p}^{4} x_{*} - 8 a_{p}^{2} b_{p}^{2} x_{*}}{(4 a_{p}^{4} - 8 a_{p}^{2} b_{p}^{2} + 4 b_{p}^{4}) / a_{p}^{2}} x_{e} + \frac{- 4 a_{p}^{4} x_{*}^{2}}{(4 a_{p}^{4} - 8 a_{p}^{2} b_{p}^{2} + 4 b_{p}^{4}) / a_{p}^{2}} .

This expression is a fourth-order polynomial, which we write in standard quartic form as

Eq. (23)

0 = x_{e}^{4} + A_{0} x_{e}^{3} + B_{0} x_{e}^{2} + C_{0} x_{e} + D_{0},

where

A_{0}

,

B_{0}

,

C_{0}

, and

D_{0}

are the constants and functions of

a_{p}

,

b_{p}

,

x_{*}

, and

y_{*}

, defined in Appendix F. We now apply the analytical solutions of the general quartic expression in standard form given in Appendix K.

Solving this quartic gives us a set of four $x_{e}$ solutions ( $x_{e}$ ) corresponding to two global extrema and two local extrema (if the latter exist for a particular orbit). Figure 6 shows a schematic representation of an orbit with four extrema. We use the imaginary components of the solutions, geometry of the ellipse, and magnitude of the higher order terms in the quartic solutions to identify which solutions belong to which extrema. We first use the magnitude of the imaginary components of solutions to determine how many extrema there are and filter out solutions that are not extrema. The algebraic solutions of the quartic polynomial all give values in the first quadrant (quadrants 1 to 4 are numbered counter-clockwise such that quadrant 1 has strictly positive coordinates, see the four quadrants of Fig. 6), so we must use geometry of the problem to determine the proper sign of each extremum’s true coordinates $x_{e, g}$ and $y_{e, g}$ [ $g$ references solutions 0 to 3 Eq. (91)–Eq. (94)].

Fig. 6

Planet orbit (black) in the $d r$ frame and planet-star separation extrema. The minimum separation (cyan) always occurs in the same quadrant as the star in the $d r$ frame (orange ×). The maximum separation (red) always occurs in the quadrants opposite the star. The local minimum separation (magenta) and local maximum separation (gold) occur in the same half-plane about the $y$ axis as the star, but in opposite half-planes about the $x$ axis as the star (fourth quadrant). This configuration applies to the vast majority of orbits. A small number of edge cases exist, including circular orbits and some edge on orbits. The blue dots represent the foci of the projected ellipse in the $d r$ frame, the orange dot is the ellipse center, and the purple dashed lines are the semimajor and semiminor axes.

All numerical solutions to the quartic have some degree of imaginary component due to accumulation of numerical errors. Quartic solution sets with only two extrema will have two real solutions (solutions with small imaginary components only due to numerical error) and two solutions with large imaginary components (algebraic solutions that are artifacts to be thrown away). The majority of KOE have only two $s$ extrema. The only case in which no extrema exist is in a circular face-on orbit. In this specific case, all solutions to the quartic will be nearly identical and have large imaginary components, and the resulting $s$ extrema will be identical ( $s_{\min} = s_{\max}$ ). We assume that all solutions are real if $| I (x_{e, g}) | < 10^{- 5} \forall g \in {0,1, 2,3}$ . We assume only two solutions are real if $| I (x_{e}) | < 10^{- 5}$ for only two solutions. We define a new ordered set containing either two or four elements depending on magnitude of the imaginary components as

Eq. (24)

x_{R} = {| R (x_{e, g}) | : | I (x_{e, g}) | < 10^{- 5} \forall g \in 0..3},

where the absolute value is due to the algebraic solutions of the quartic being only defined in the first quadrant.

We now know the number of expected solutions from the dimension of $x_{R}$ but need to leverage the quartic algebraic solution and geometry to determine which components belong to which quadrant. For orbits with $x_{R}$ containing only two solutions, the first two solutions to the quartic ( $x_{R, 0}$ and $x_{R, 1}$ ) produce the largest and smallest $x_{R, g} - x_{*}$ , respectively. Due to the shape of an ellipse, these must necessarily produce $s_{\min}$ and $s_{\max}$ . We define the set $y_{R}$ as the application of Eq. (18) to each element of $x_{R}$ . We define four separation quantities from the possible sign combinations of the coordinate magnitudes in $x_{R}$ and $y_{R}$ as

Eq. (25)

s_{- \pm 0} = \sqrt{{(x_{R, 0} - x_{*})}^{2} + {(y_{R, 0} \pm y_{*})}^{2}},

Eq. (26)

s_{+ \pm 0} = \sqrt{{(x_{R, 0} + x_{*})}^{2} + {(y_{R, 0} \pm y_{*})}^{2}},

Eq. (27)

s_{- \pm 1} = \sqrt{{(x_{R, 1} - x_{*})}^{2} + {(y_{R, 1} \pm y_{*})}^{2}},

Eq. (28)

s_{+ \pm 1} = \sqrt{{(x_{R, 1} + x_{*})}^{2} + {(y_{R, 1} \pm y_{*})}^{2}} .

Typically, $s_{-- 1}$ is smaller than $s_{-- 0}$ , but in $≪ 0.01 %$ of cases, $s_{-- 0}$ can be the smallest separation due to numerical error. (In this case, we ensure $| s_{-- 1} - s_{-- 0} | < 10^{- 8}$ and swap values.) In addition to $s_{-- 1} < s_{-- 0}$ , we also know that $s_{- + 0} < s_{- + 1}$ , $s_{+ - 1} < s_{+ - 0}$ , and $s_{++ 1} < s_{++ 0}$ for all cases. Using this knowledge, we reduce the number of comparisons that we need to find the minimum planet-star separation to

Eq. (29)

s_{\min} = {\begin{cases} s_{-- 1} where (s_{-- 1} < s_{- + 0}) and (s_{-- 1} < s_{+ - 1}) and (s_{-- 1} < s_{++ 1}), \\ s_{- + 0} where (s_{- + 0} < s_{-- 1}) and (s_{- + 0} < s_{+ - 1}) and (s_{- + 0} < s_{++ 1}), \\ s_{+ - 1} where (s_{+ - 1} < s_{- + 0}) and (s_{+ - 1} < s_{-- 1}) and (s_{+ - 1} < s_{++ 1}), \\ s_{++ 1} where (s_{++ 1} < s_{- + 0}) and (s_{++ 1} < s_{+ - 1}) and (s_{++ 1} < s_{-- 1}) . \end{cases}

The maximum separation similarly is found as

Eq. (30)

s_{\max} = {\begin{matrix} s_{-- 0}, where (s_{-- 0} > s_{- + 0}) and (s_{-- 0} > s_{+ - 1}) and (s_{-- 0} > s_{++ 0}), \\ s_{- + 0}, where (s_{- + 0} > s_{-- 0}) and (s_{- + 0} > s_{+ - 1}) and (s_{- + 0} > s_{++ 0}), \\ s_{+ - 1}, where (s_{+ - 1} > s_{- + 0}) and (s_{+ - 1} > s_{-- 0}) and (s_{+ - 1} > s_{++ 0}), \\ s_{++ 0}, where (s_{++ 0} > s_{- + 0}) and (s_{++ 0} > s_{+ - 1}) and (s_{++ 0} > s_{-- 0}) . \end{matrix}

We are able to find

x_{e}

and

y_{e}

of the minimum and maximum separation diamonds drawn in Fig. 6 using the same logic as for finding

s_{\min}

and

s_{\max}

.

For the orbits with four solutions in $x_{R}$ , the first two elements will always be the global extrema, and the last two elements will be the local extrema, which are given by

Eq. (31)

s_{\min} = s_{-- 1},

Eq. (32)

s_{\max} = s_{++ 0},

Eq. (33)

s_{lmin} = {\begin{cases} s_{- + 3}, & s_{- + 2} > s_{- + 3} \\ s_{- + 2}, & else \end{cases},

Eq. (34)

s_{lmax} = {\begin{cases} s_{- + 2}, & s_{- + 2} > s_{- + 3} \\ s_{- + 3}, & else \end{cases},

where

s_{- + 2}

and

s_{- + 3}

are calculated in the same manner as

s_{- + 0}

and

s_{- + 0}

. We are able to find the coordinates of all four using the same logic as for finding

s_{\min}

,

s_{\max}

,

s_{lmin}

, and

s_{lmax}

. This procedure yields the coordinates of all existing extrema in the

d r

frame. To find their locations on the projection of orbit in the plane of the sky, we apply the inverse of Eq. (16).

2.1.3.

Intersections between a circle and an ellipse

We find $s_{WA}$ -orbit intersections by formulating the circle-ellipse intersections as another quartic, solving this, and assigning the algebraic solutions to intersections. We assign solutions using the number of $s$ extrema, the size of the $s_{WA}$ intersecting circle relative to these $s$ extrema, and the ellipse geometry. To formulate the $s_{WA}$ -orbit intersections as a quartic, we start with Eq. (19) and substitute in Eq. (18). This gives us a separation equation solely as a function of $x_{e}$ and star location. We expand and transform this into a general polynomial of $x_{e}$ with a general $s$

Eq. (35)

0 = (\frac{a_{p}^{4} - 2 a_{p}^{2} b_{p}^{2} + b_{p}^{4}}{a_{p}^{4}}) x_{e}^{4} + (\frac{- 4 a_{p}^{2} x + 4 b_{p}^{2} x_{*}}{a_{p}^{2}}) x_{e}^{3} + (\frac{2 a_{p}^{2} b_{p}^{2} - 2 a_{p}^{2} s^{2} + 6 a_{p}^{2} x_{*}^{2} + 2 a_{p}^{2} y_{*}^{2} - 2 b_{p}^{4} + 2 b_{p}^{2} s^{2} - 2 b_{p}^{2} x_{*}^{2} + 2 b_{p}^{2} y_{*}^{2}}{a_{p}^{2}}) x_{e}^{2} + (- 4 b_{p}^{2} x + 4 s^{2} x - 4 x_{*}^{3} - 4 x_{*} y_{*}^{2}) x_{e} + (b_{p}^{4} - 2 b_{p}^{2} s^{2} + 2 b_{p}^{2} x_{*}^{2} - 2 b_{p}^{2} y_{*}^{2} + s^{4} - 2 s^{2} x_{*}^{2} - 2 s^{2} y_{*}^{2} + x_{*}^{4} + 2 x_{*}^{2} y_{*}^{2} + y_{*}^{4}) .

As in Sec. 2.1.2, we divide by the leading coefficient to convert this to the general quartic form:

Eq. (36)

0 = x_{e}^{4} + A_{1} x_{e}^{3} + B_{1} x_{e}^{2} + C_{1} x_{e} + D_{1},

with

A_{1}

,

B_{1}

,

C_{1}

, and

D_{1}

given in Appendix G. We solve this using the general quartic solution as given in Appendix K. This results in a solution for the

x_{e}

of intersections in the first quadrant of the

d r

frame.

We always have four algebraic $x_{e}$ solutions that may or may not correspond to actual $s_{WA}$ -orbit intersections. The number of $s$ extrema (either two or four) and projected separation relative to these extrema determine how many intersections will occur. (If we expect two intersections, then two of the four algebraic solutions must be real solutions, and the other two are some combination of repeated roots or non-physical imaginary solutions.) In each of these cases, we must handle the assignment of $x_{e}$ and $y_{e}$ solutions to the correct quadrants. Figure 7 shows a schematic representation of two cases corresponding to four total intersections.

Fig. 7

Diagrams of orbits where the projected ellipse in the $d r$ frame produces four intersections (green dots) with the $s_{WA}$ separation circle (green circle). The general orbit’s projected ellipse (black) is centered at the origin, and the projected ellipse axes (purple) define the ${\bar{\hat{x}}}_{d r}$ and ${\overline{\hat{y}}}_{d r}$ axes. (Note that the projected ellipse semimajor axis has been derotated such that the star is always located in the first quadrant.) The separation circle center (orange dot) is the star’s location relative to the orbit’s projected ellipse. The $x_{h}$ points are the quartic solutions, where subscripts are reordered in ascending distance from the star’s $x$ position. This ordering means that $x_{2}$ must always occur in the fourth quadrant, $x_{0}$ may occur in either the first or second quadrants, and $x_{1}$ may occur in either the third or fourth quadrants. For 99.992% of KOE sampled from the SAG13 population that produce four intersections, the $x_{3}$ intersection occurs in the first quadrant (a). The other 0.008% of KOE result in the $x_{3}$ intersection occurring in the fourth quadrant (b).

For KOE with four extrema and $s_{\min} < s_{WA} < s_{lmin}$ , we know that there will be two intersections on the same $y$ side of the ellipse as the star (quadrants 1 or 2 in the $d r$ frame). Of the four quartic solutions to Eq. (36) that we have to choose from, we know that $x_{0}$ is one of them. The other solution could either be $x_{1}$ if $I (x_{1}) < 10^{- 9}$ or $x_{3}$ if $I (x_{1}) > 10^{- 9}$ .

For KOE with four extrema and $s_{lmax} < s_{WA} < s_{\max}$ , we know that there will be two intersections on the opposite $x$ side of the ellipse as the star (quadrants 2 and 3 in the $d r$ frame). In all cases, $x_{0}$ and $x_{1}$ are the intersection solutions. This is because the first two solutions have the largest term in the quartic solution [Eqs. (91) and (92)]. $x_{1}$ is slightly smaller than $x_{0}$ because it subtracts the second largest term. The relative magnitudes of $x_{0}$ and $x_{1}$ determine that $y_{0}$ occurs in the same side of the ellipse as the star and $y_{1}$ must occur on the opposite side of the ellipse as the star. Therefore, ( $x_{0}$ , $y_{0}$ ) occurs in quadrant 2, and ( $x_{1}$ , $y_{1}$ ) occurs in quadrant 3.

For KOE with four extrema and $s_{lmin} < s_{WA} < s_{lmax}$ , we know that there will be four intersections. Unlike in Sec. 2.1.2 and the rest of this paper where $x_{0}$ through $x_{3}$ are ordered as in Appendix K, we order $x_{h}$ based off $Δ x_{h} = | x_{h} - x_{*} |$ . We order the quartic solutions from $x_{0}$ to $x_{3}$ such that $x_{0}$ is where $\min ({Δ x_{h} \forall h \in {0,1, 2,3}})$ and $x_{3}$ is where $\max ({Δ x_{h} \forall h \in {0,1, 2,3}})$ . These newly ordered $Δ x_{h}$ correspond to those shown in Fig. 7. The ( $x_{2}, y_{2}$ ) intersection occurs in either quadrant 1 or 2, but it is always above and to the left of the star in the $d r$ frame and has the second largest $Δ x_{h}$ component. The ( $x_{3}$ , $y_{3}$ ) intersection occurs in either quadrant 1 or 4, but it always has the largest $Δ x_{h}$ component. We resolve the sign of $y_{3}$ by testing it in both quadrants. In $> 99.992 %$ of cases, we assign $y_{3}$ to quadrant 1 as in Fig. 7(a), but for a minority of cases, its correct assignment is quadrant 4 as in Fig. 7(b). The ( $x_{1}, y_{1}$ ) intersection always occurs in quadrant 4. The ( $x_{0}, y_{0}$ ) intersection occurs in either quadrant 3 or quadrant 4, but it is always below and to the left of the star in the $d r$ frame and has the smallest $Δ x_{h}$ component.

For KOE with two extrema and $s_{\min} < s_{WA} < s_{\max}$ , we know that there will be two circle-ellipse intersections. The KOE determine where the star is located in the $d r$ frame. The location of the star in the $d r$ frame relative to the vertices of the projected ellipse [ $(0, b_{p})$ , $(0, - b_{p})$ , $(- a_{p}, 0)$ , and $(a_{p}, 0)$ ] determines the star-vertex separation ordering and subsequently to which quadrants the two intersection solutions belong. Instead of calculating the star-vertex separation for each orbit, we divide the first quadrant into four regions that specify the star location type [types 0 to 3 as indicated in Fig. 8(a)]. This means that any KOE with the star in location type 2 has the associated star-vertex separation ordering.

Fig. 8

The regions identifying the star type and which vertices are closest to the host star. (a) The first quadrant of a projected ellipse (black curve) with the semimajor and semiminor axes (purple) and three dashed lines dividing the quadrant into four regions defining the separation ordering. The pink dashed line represents the line of points equidistant from $(0, b_{p})$ and $(a_{p}, 0)$ . The gray dashed line represents the line of points equidistant from $(- a_{p}, 0)$ and $(0, - b_{p})$ . The turquoise dashed line represents the line of points equidistant from $(0, - b_{p})$ and $(a_{p}, 0)$ . The yellow dashed line represents the line of points equidistant from $(- a_{p}, 0)$ and $(0, b_{p})$ . We say that the star in (a) is a type 2 star and has the associated separation ordering. (b)–(e) Color each region of the ellipse, identifying which vertex is closest, second closest, third closest, and fourth closest. The star in (a) is type 2 and has the top vertex as the closest as seen in (b). These plots are based off the projected ellipse of a planet with $a = 0.40 AU$ , $e = 0.23$ , $i = 0.69 rad$ , $Ω = 3.49 rad$ , and $ω = 5.64 rad$ .

Using the equidistant lines between ellipse vertices, we divide the first quadrant into four regions (types 0 to 3) as shown in Fig. 8(a). Regions 0 and 2 are divided by the line defined by

Eq. (37)

y_{a_{p} - x, b_{p} - y} = \frac{a_{p}}{b_{p}} x + \frac{a_{p}^{2}}{2 b_{p}} - \frac{b_{p}}{2},

where a projected ellipse and star is of type 0 if

y_{*} > y_{a_{p} - x, b_{p} - y} (x_{*})

and one of types 1, 2, or 3 otherwise depending on the other equidistant lines. Similarly, regions 1 and 2 are divided by the line defined by

Eq. (38)

y_{a_{p} + x, b_{p} - y} = - \frac{a_{p}}{b_{p}} x + \frac{a_{p}^{2}}{2 b_{p}} - \frac{b_{p}}{2},

regions 2 and 3 are divided by the line defined by

Eq. (39)

y_{a_{p} + x, b_{p} + y} = \frac{a_{p}}{b_{p}} x - \frac{a_{p}^{2}}{2 b_{p}} + \frac{b_{p}}{2},

and the star-ellipse classification type is calculated similar to type 0.

An example star is shown in Fig. 8(a). This star location is type 2 because it is in the region bounded by the three dashed lines (gray, pink, and teal). This type 2 star has the separation ordering indicated in Fig. 8(a) and Table 1. Let us consider a $s_{WA}$ separation circle such that $s_{a_{p} - x_{*}, y_{*}} < s_{WA} < s_{x_{*}, b_{p} + y_{*}}$ . We, therefore, know that the two intersections must occur in quadrant 2 and quadrant 4 of the $d r$ frame.

Table 1

Separation order from smallest to largest by star location type.

Type	Condition	Separation order from smallest to largest
Type	Condition	First	Second	Third	Fourth
0	$s_{a_{p} + x_{}, y_{}} < s_{x_{}, b_{p} + y_{}}$	$s_{x_{}, b_{p} - y_{}}$	$s_{a_{p} - x_{}, y_{}}$	$s_{a_{p} + x_{*}, b_{p}}$	$s_{x_{}, b_{p} + y_{}}$
1	$s_{x_{}, b_{p} + y_{}} < s_{a_{p} - x_{}, y_{}}$	$s_{x_{}, b_{p} - y_{}}$	$s_{x_{}, b_{p} + y_{}}$	$s_{a_{p} - x_{}, y_{}}$	$s_{a_{p} + x_{}, y_{}}$
2	$s_{a_{p} - x_{}, y_{}} < s_{x_{}, b_{p} + y_{}}$ $s_{x_{}, b_{p} + y_{}} < s_{a_{p} + x_{}, y_{}}$ $s_{x_{}, b_{p} - y_{}} < s_{a_{p} - x_{}, y_{}}$	$s_{x_{}, b_{p} - y_{}}$	$s_{a_{p} - x_{}, y_{}}$	$s_{x_{}, b_{p} + y_{}}$	$s_{a_{p} + x_{}, y_{}}$
3	$s_{a_{p} - x_{}, y_{}} < s_{x_{}, b_{p} - y_{}}$	$s_{a_{p} - x_{}, y_{}}$	$s_{x_{}, b_{p} - y_{}}$	$s_{x_{}, b_{p} + y_{}}$	$s_{a_{p} + x_{}, y_{}}$

The distances $s_{x_{*}, b_{p} + y_{*}}$ , $s_{x_{*}, b_{p} - y_{*}}$ , $s_{a_{p} + x_{*}, y_{*}}$ , and $s_{a_{p} - x_{*}, y_{*}}$ in Fig. 8(a) are the distances of the star to each of the ellipse vertices and are calculated by

Eq. (40)

s_{x_{*}, b_{p} \pm y_{*}} = \sqrt{x_{*}^{2} + {(b_{p} \pm y_{*})}^{2}},

Eq. (41)

s_{a_{p} \pm x_{*}, y_{*}} = \sqrt{{(a_{p} \pm x_{*})}^{2} + y_{*}^{2}} .

Figure 8(a) and Table 1 define orbit types sorted by vertex distances from smallest to largest based on star location type in the first quadrant of the

d r

frame.

Finally, after determining the correct number of intersection solutions and the proper quadrants to which these solutions belong, we rerotate and translate the intersection solution locations back into the 2D projection of the 3D orbit.

2.2.

$Δ mag$ Intersections

In this section, we present our method for calculating the values of $ν$ on a planet’s orbit where the planet has a specific value of $Δ {mag}_{\lim}$ , called $Δ mag$ intersections. As in Sec. 2.1, to compute these solutions, we first need to calculate all $Δ mag$ extrema ( $Δ {mag}_{\min}$ , $Δ {mag}_{\max}$ , $Δ {mag}_{lmin}$ , and $Δ {mag}_{lmax}$ ). The process for calculating $Δ mag$ extrema is detailed in Appendix C, but briefly included as follows:

1. express $Δ mag$ as a polynomial in $\cos (ν)$ ;
2. find the values of $ν$ and $Δ mag$ for all polynomial roots; and
3. remove invalid and duplicate solutions.

The general process for calculating the true anomalies where a planet has $Δ {mag}_{\lim}$ is similar to the extrema-finding process but with some minor modifications. After calculating the $Δ mag$ extrema, we determine how many $Δ mag$ intersections a given orbit should have with a particular $Δ {mag}_{\lim}$ value. If $Δ {mag}_{\lim} < Δ {mag}_{\min}$ or $Δ {mag}_{\lim} > Δ {mag}_{\max}$ , there are no intersections. If $Δ {mag}_{\min} < Δ {mag}_{\lim} < Δ {mag}_{lmin}$ or $Δ {mag}_{lmax} < Δ {mag}_{\lim} < Δ {mag}_{\max}$ , then there are exactly two intersections. If $Δ {mag}_{lmin} < Δ {mag}_{\lim} < Δ {mag}_{lmax}$ , then there are exactly four intersections. When the orbit does not contain local $Δ mag$ extrema, there are just two $Δ mag$ intersections if $Δ {mag}_{\min} < Δ {mag}_{\lim} < Δ {mag}_{\max}$ .

Knowing the number of solutions to expect, we follow the same steps as above: finding a governing polynomial equation, solving for its roots (represented by the ordered set $x$ ), and filtering out the relevant solutions. We additionally throw out any solutions with large errors from the input $Δ {mag}_{\lim}$ . The full process outline is included in Appendix D and discussed in depth below.

We start with the definition of $Δ mag$ given in Eq. (1). Although this expression contains multiple terms that are fully defined by an orbit’s KOE, it is also a function of other planet properties, including the planet’s radius, geometric albedo, and phase function. To make the mathematical development presented below tractable, we assume that the quasi-Lambert phase function¹⁷ is a sufficient approximation of any planet’s phase function. In general, the quasi-Lambert phase function is a better representation of the Earth’s phase function than the Lambert phase function.

This quasi-Lambert phase function is given by

Eq. (42)

Φ_{L} (β) = \cos^{4} (\frac{β}{2}) .

We take advantage of this function’s form by substituting the half-angle formula:

Eq. (43)

\cos (\frac{β}{2}) = \sqrt{\frac{1 + \cos (β)}{2}} .

At the same time, from the orbital geometry defined in Fig. 3 and Eq. (3), we write

Eq. (44)

β = \cos^{- 1} (\sin (i) \sin (ν + ω)),

where the final term is expanded by the angle addition formula as

Eq. (45)

\sin (ν + ω) = \sin (ν) \cos (ω) + \cos (ν) \sin (ω) .

Substituting in

β

from Eq. (44) expanded with Eq. (45) into Eq. (43) allows us to reduce the order of the fully substituted Eq. (42). Note that this expression for

β

depends on making the approximation that the observer-star and observer-planet vectors are parallel, which introduces minor error given the large distances to even the nearest stars.¹⁸

After making these substitutions, as well as taking Eq. (4) for the planet-star distance term $| {\underset{̲}{r}}_{k / i} |$ , simplifying, and collecting all of the non-orbital planet parameters on the left side, we find

Eq. (46)

\frac{a^{2} {(1 - e^{2})}^{2}}{p R^{2}} 10^{- 2.5 Δ mag} = \frac{1}{4} {(e \cos (ν) + 1)}^{2} {(1 + \sin (i) \cos (ω) \sin (ν) + \sin (i) \sin (ω) \cos (ν))}^{2} .

We now have an expression that isolates terms containing

ν

and can be decomposed into a numerically solvable polynomial. To improve solving efficiency and determine which subset of planets should have

ν

solutions for a given

Δ mag

, we need to first calculate the

Δ mag

extrema over the full orbit. The outline of the process for calculating these

Δ mag

extrema is included in Appendix C.

To calculate these $Δ mag$ extrema, we first find the derivative of Eq. (46) and multiply by two (for simplification purposes) to get

Eq. (47)

0 = - e^{2} \sin^{2} (i) \sin^{3} (ν) \cos (ν) \cos^{2} (ω) - 3 e^{2} \sin^{2} (i) \sin^{2} (ν) \sin (ω) \cos^{2} (ν) \cos (ω) - 2 e^{2} \sin^{2} (i) \sin (ν) \sin^{2} (ω) \cos^{3} (ν) + e^{2} \sin^{2} (i) \sin (ν) \cos^{3} (ν) \cos^{2} (ω) + e^{2} \sin^{2} (i) \sin (ω) \cos^{4} (ν) \cos (ω) - 2 e^{2} \sin (i) \sin^{2} (ν) \cos (ν) \cos (ω) - 3 e^{2} \sin (i) \sin (ν) \sin (ω) \cos^{2} (ν) + e^{2} \sin (i) \cos^{3} (ν) \cos (ω) - e^{2} \sin (ν) \cos (ν) - e \sin^{2} (i) \sin^{3} (ν) \cos^{2} (ω) - 4 e \sin^{2} (i) \sin^{2} (ν) \sin (ω) \cos (ν) \cos (ω) - 3 e \sin^{2} (i) \sin (ν) \sin^{2} (ω) \cos^{2} (ν) + 2 e \sin^{2} (i) \sin (ν) \cos^{2} (ν) \cos^{2} (ω) + 2 e \sin^{2} (i) \sin (ω) \cos^{3} (ν) \cos (ω) - 2 e \sin (i) \sin^{2} (ν) \cos (ω) - 4 e \sin (i) \sin (ν) \sin (ω) \cos (ν) + 2 e \sin (i) \cos^{2} (ν) \cos (ω) - e \sin (ν) - \sin^{2} (i) \sin^{2} (ν) \sin (ω) \cos (ω) - \sin^{2} (i) \sin (ν) \sin^{2} (ω) \cos (ν) + \sin^{2} (i) \sin (ν) \cos (ν) \cos^{2} (ω) + \sin^{2} (i) \sin (ω) \cos^{2} (ν) \cos (ω) - \sin (i) \sin (ν) \sin (ω) + \sin (i) \cos (ν) \cos (ω),

which is a function of the

\sin (ν)

and

\cos (ν)

terms. We now fully expand this expression and substitute in

Eq. (48)

\sin (ν) = \sqrt{1 - \cos^{2} (ν)}

to get an expression in

\cos (ν)

only. We define

x = \cos (ν)

, isolate the

\sqrt{1 - x^{2}}

term, square both sides, and expand to get an eighth degree polynomial in

x

of the form

Eq. (49)

0 = A_{2} x^{8} + B_{2} x^{7} + C_{2} x^{6} + D_{2} x^{5} + E_{2} x^{4} + F_{2} x^{3} + G_{2} x^{2} + H_{2} x + I_{2} .

The coefficients of this expression are included in Eq. (75) in Appendix H.

Although eighth degree polynomials do not have analytical solutions, a numerical root solver can determine the eight roots of this function denoted as the ordered set $x$ . We discard any $x_{p} > 1 \forall p \in {0 \dots 8}$ , $x_{p} < - 1 \forall p \in {0 \dots 8}$ , or solutions with large imaginary components. We then calculate the remaining true anomaly solutions by $ν_{0} = \cos^{- 1} (x)$ and $ν_{1} = 2 π - ν_{0}$ . Of the remaining valid solutions, we identify whether each solution is an extremum by evaluating whether $Δ mag (ν_{0} \pm δ ν)$ are both larger or both smaller than $Δ mag (ν_{0})$ . If identified as a potential extremum, the smallest and largest extrema are assigned to $Δ {mag}_{\min}$ and $Δ {mag}_{\max}$ . The remaining extrema are checked for duplicates, which are identified by solutions with ( $ν$ , $Δ mag$ ) values close to existing extrema. If any solutions are remaining, an additional check and assignment is then made for the local extrema. Through this process, a solution identified in $ν_{0}$ has the associated solution in $ν_{1}$ removed.

To calculate the $Δ mag$ intersections, we apply the same process for turning Eq. (46) into a polynomial as in the $Δ mag$ extrema calculation. The full outline for calculating the $Δ mag$ intersections is given in Appendix D. By following this process, we arrive at

Eq. (50)

ξ = A_{3} x^{8} + B_{3} x^{7} + C_{3} x^{6} + D_{3} x^{5} + E_{3} x^{4} + F_{3} x^{3} + G_{3} x^{2} + H_{3} x + I_{3},

where

ξ

is the collection of constants on the left side of Eq. (46). The coefficients of this polynomial are given in Eq. (76) of Appendix I.

We again use a numerical root solver to find the roots of this function $x$ and again discard all $x_{p} > 1$ $\forall p \in {0 \dots 8}$ , $x_{p} < - 1 \forall p \in {0 \dots 8}$ , and solutions with large imaginary components. We then calculate the true anomalies of all remaining solutions by $ν_{0} = \cos^{- 1} (x)$ and $ν_{1} = 2 π - ν_{0}$ . We then evaluate ${Δ mag}_{0} = Δ mag (ν_{0})$ and ${Δ mag}_{1} = Δ mag (ν_{1})$ and remove solutions where ${| Δ mag}_{0} - Δ {mag}_{\lim} | > 0.01$ and $| {Δ mag}_{1} - Δ {mag}_{\lim} | > 0.01$ . (Note that 0.01 is $< \pm 0.08 %$ error on $Δ mag$ .) We do an iterative process of selecting ( $ν$ , $Δ mag$ ), which are unique (not duplicate solutions) and are closest to the expected $Δ {mag}_{\lim}$ , until we have the expected number of solutions. In the majority of cases in which solutions exist, there are generally only two viable solutions, the assignment of which is simple. In some cases, the solution selection is more ambiguous as double roots are possible.

2.3.

$ν$ from $X$ and $Y$

Sections 2.1.2, 2.1.3, and 2.2 give locations of extrema at intersections in the plane of the sky of the form $(x, y)$ , but we need to know the true anomalies of the orbit where these intersections occur.

We start with $X$ and $Y$ in Eq. (3) and solve for $(1 + e \cos (ν)) / (a (1 - e^{2}))$ , resulting in

Eq. (51)

\frac{1 + e \cos (ν)}{a (1 - e^{2})} = \frac{1}{Y} [\sin Ω \cos (ω + ν) + \cos Ω \sin (ω + ν) \cos i],

and

Eq. (52)

\frac{1 + e \cos (ν)}{a (1 - e^{2})} = \frac{1}{X} [\cos (Ω) \cos (ω + ν) - \sin (Ω) \sin (ω + ν) \cos (i)] .

In both of these expressions, we substitute angle addition formulas of

\sin (ω + ν)

and

\cos (ω + ν)

and subsequently set the two equations equal to each other.

We set Eqs. (51) and (52) equal to each other, expand, isolate a $\cos (ν)$ term and $\sin (ν)$ term, solve for $ν$ , and rearrange to get

Eq. (53)

ν = \tan^{- 1} [\frac{- \frac{X}{Y} \sin Ω \cos ω - \frac{X}{Y} \cos Ω \cos i \sin ω + \cos Ω \cos ω - \sin Ω \cos i \cos ω}{- \frac{X}{Y} \sin Ω \sin ω + \frac{X}{Y} \cos Ω \cos i \cos ω + \cos Ω \sin ω + \sin Ω \cos i \cos ω}] .

We now have an analytical expression for

ν

solely as a function of the KOE,

X

, and

Y

of a particular point on the orbital ellipse.

By combining the original $X$ and $Y$ equations to solve for $ν$ analytically, we have created two potential solutions at $ν$ and $ν \pm π$ . One of these is the correct $ν$ value and the other is not. To calculate the correct intersection point, we calculate the separations at both $ν$ and $ν \pm π$ to find absolute error in $s$ and use the smaller error of the two.

2.4.

Calculate $t$ from $ν$

We have $ν$ for the locations where $s$ extrema, $Δ mag$ extrema, $s_{WA}$ -orbit intersections, and $Δ mag$ intersections occur, but we need them in terms of time. We calculate eccentric anomaly $E$ of these events directly from $ν$ as

Eq. (54)

E = \tan^{- 1} (\frac{\sqrt{1 - e^{2}} \sin ν}{e + \cos ν}),

which gives the corresponding time

Eq. (55)

t = \frac{E - e \sin (E) T}{2 π},

where

T

is the orbital period of a general planet.

2.5.

Converting from Star to Star

The planet-star intersection points and visibility ranges can be saved as either true anomalies or as specific times for a reference star. It makes the most sense to store the values as times and use $1 M_{⊙}$ as the reference time. The times are scaled to any star mass $M$ by

Eq. (56)

T = T_{⊙} \sqrt{\frac{M_{⊙}}{M}} .

2.6.

Calculating Completeness

We calculate integration time adjusted completeness using the aggregated fraction of time that planets are detectable. Using the methods described in Sec. 2 on the instrument’s photometric visibility limit ( $Δ {mag}_{\lim}$ ) and astrometric visibility limits (IWA and OWA), we calculate the specific times that the planet enters or exists the instruments visibility limits. By collecting these times and inspecting intermediate test points, we identify the time windows in which any planet is detectable by the instrument. Given an integration time ( $t_{\max}$ ) required to reach some SNR consistent with a clear detection (typically a value > 5), we throw out all visibility time windows less than $t_{\max}$ and discount all other visibility windows by $t_{\max}$ . Dividing the sum of all integration time discounted visibility windows for a planet by its orbital period gives us the fraction of time that a planet is detectable by the instrument. Aggregating the fractions of time that planets are detectable by the instrument gives us integration time adjusted completeness. In mathematical form, this is given as

Eq. (57)

C_{t_{\max}} = \frac{1}{N} (\sum_{\forall k} \frac{\sum_{\forall j} (δ t_{j, k} - t_{\max} U_{k})}{T_{k}}) .

Here

δ t_{j, k}

is the

j

’th time window larger than

t_{\max}

in which the

k

’th planet is visible,

U_{k}

is a boolean that indicates that the planet is always visible 0 or at least sometimes visible 1,

N

is the total number of planets, and

T_{k}

is the orbital period for the

k

’th planet.

We also calculate the completeness for an Earth-like exoplanet of the Earth-like exoplanets subpopulation:

Eq. (58)

C_{\oplus, t_{\max}} = \frac{1}{N_{\oplus}} (\sum_{\forall k \in \oplus} \frac{\sum_{\forall j} (δ t_{j, k} - t_{\max} U_{k})}{T_{k}}),

where

N_{\oplus}

is the total number of Earth-like planets simulated.

Similarly, this capability can be extended to the entirety of the Kopparapu classification scheme.⁴ Our method uses orders of magnitude fewer exoplanets and an equivalent memory but fills in the $Δ mag$ space using an order of magnitude better accounting with strategically calculated true anomalies.

2.7.

Calculating Dynamic Completeness

The methods presented in this paper can also be used to calculate dynamic completeness.¹³ Dynamic completeness extends the original concept of completeness by considering the fraction of planets observable on subsequent observations of the same target. For the original formulation of completeness, which relies on simulating a cloud of planets, this requires propagating every simulated planet along its orbits. Dynamic completeness is frequently used to compute the fraction of planets in a population that are initially undetectable but become detectable upon a second observation some time later.

We define $P$ as the unordered set of all planets, $P_{detected, 1}$ as the set of planets detected at observation time one, and $P_{undetected, 1}$ as the set of planets not detected at observation time one. To find the sets of detected planets at the first observation, we start by generating a random observation time for the $k$ ’th planet ( $t_{start, k}$ ) between 0 and $T_{k}$ . We then determine which planets are within the instrument’s visibility limits at $t_{start, k}$ to create $P_{detected, 1}$ and its complement $P_{undetected, 1}$ .

For a subsequent observation some time later ( $t_{wait}$ ), we find the time past periastron at which the observation occurs by [ $\mod (t_{start} + t_{wait}, T_{k})$ ] and determine which planets are within the instrument’s visibility limits at that time. This defines the unordered set of planets detected at observation time two as $P_{detected, 2}$ . Its complement is $P_{undetected, 2}$ , the set of planets not detected at observation time two.

Dynamic completeness for the second observation of a target star is given by the fraction of planets undetected in the initial observation but detected on the second one:

Eq. (59)

C_{2} = \frac{1}{N} n (P_{undetected, 1} \cap P_{detected, 2}),

where

n (P)

is a function giving the number of elements in set

P

and

N

is the total number of planets.

Dynamic completeness of the $m$ ’th visit is similarly given by the fraction of previously undetected planets that are detected on that visit:

Eq. (60)

C_{m} = \frac{1}{N} n ((⋂_{j = 0}^{m - 1} P_{undetected, j}) \cap P_{detected, m}) .

3. Results

3.1.

$ν$ from $s$

We apply the methods developed Sec. 2 to calculate the planet-star separation extrema, the true anomalies at which they occur, the times past periastron at which the extrema occur, the true anomalies where the separation circle intersects with the projected orbital ellipse, and the times at which these intersections occur for a planet orbit. Figure 9 shows the projected separation of a sample orbit with global and local extrema calculated via Sec. 2 methods. This figure also demonstrates the ability of our methods to find all true anomalies (and times) when the projected separation takes a specific value (in this case $s_{WA} = 1 AU$ ).

We now determine the accuracy of the methods implemented. To do this, we calculate the true anomalies of $s_{WA}$ -orbit intersections for $10^{5}$ orbits using the method in Sec. 2.1.3 and $s_{WA} = 1 AU .$ We then calculate the planet-star separation by substituting these true anomalies back into Eq. (2) and plot the error between these separations and the input $s_{WA}$ in Fig. 10. Of the $10^{5}$ planets orbits simulated, $\sim 25, 000$ orbits produced two or more intersections. The largest intersection error observed is $< 10^{- 6}$ . Machine precision error is $\sim 10^{- 16}$ , and the square root of this is $10^{- 8}$ , which indicates that we have approached machine precision in our method.

3.2.

$ν$ from Δmag

We apply the methods developed Sec. 2 to calculate the $Δ mag$ extrema, the true anomalies at which they occur, the times past periastron at which the extrema occur, the true anomalies where the $Δ mag$ intersections occur, and the times at which these $Δ mag$ intersections occur for a planet orbit. Figure 11 shows the $Δ mag$ of a sample orbit with global and local extrema calculated via Sec. 2 methods. This figure also demonstrates the ability to find all true anomalies (and times) when the $Δ mag$ takes a specific value (in this case $Δ {mag}_{\lim} = 25.0$ ).

To check the error in true anomaly produced by the $Δ mag$ intersection method in this paper, we compare the $ν$ of $Δ mag$ intersections calculated with alternative numerical solving methods. The method presented in this paper finds $Δ mag$ intersections of $10^{6}$ planets within $\sim 6.3 s$ . The first error checking method uses a cubic spline fit to 300 ( $ν$ , $Δ mag$ ) points along each planet’s orbit. We then subtract $Δ {mag}_{\lim} = 29$ (our test point) from the spline and find the roots. This univariate spline root solving method is capable of executing in 419 s on the 510,120 planets in the population that produced two intersections, a rate of $8.1 \times 10^{- 4} s$ per planet. The univariate spline root solving method is $100 \times$ slower than the $Δ mag$ intersection method in Sec. 2.2. The separation error between the univariate spline method and the $Δ mag$ intersection method is $< 6 \times 10^{- 5}$ for $> 99.98 %$ of planets and $< 10^{- 4}$ for the other 0.02%. These high-error targets are planets with low variation in $Δ mag (ν)$ orbits. The second-error checking method solves for the ( $ν$ , $Δ mag$ ) intersections using a numerical minimization method on a random subset of $10^{4}$ planets. This method is far more inefficient, taking 1646 s on $10^{4}$ planets with two intersections at a rate of 0.1646 s per planet. Using a minimization function allows us to determine the error in $Δ mag$ of an intersection point to within $10^{- 8}$ . This numerical minimization method independently confirms the accuracy of solutions to the univariate spline roots method on the limited number of targets tested. We plotted the normalized frequency of true anomaly error of both the numerical minimization method and univariate spline root solving method compared against the $Δ mag$ intersection method in Fig. 12. Note that the total number of incidences of a given error is normalized by the bin width and total number of targets, so the non-normalized frequency of the $3 \times 10^{- 11}$ bin is incredibly large compared with the $5 \times 10^{- 4}$ bin. Both methods indicate that the resulting true anomalies of the $Δ mag$ intersections are within $10^{- 4} rad$ of each other.

Although the $Δ mag$ intersection method featured in this paper is orders of magnitude faster than the univariate-spline-roots method and just as accurate, the quasi-Lambert phase function required for its use is not the best phase function for all planets; however, it does fit Earth-like planet phase function quite well. At a substantial time cost, the univariate-spline-roots method can be used for planets with any phase function.

3.3.

Convergence and Validation

We test convergence of our method for calculating completeness by repeatedly calculating completeness on $10^{5}$ planets. We test convergence of Brown’s method by repeatedly calculating completeness over a set of a logarithmically increasing number of planets. We use the same SAG13 planet population parameters from Ref. 1 but make the substitution of a quasi-Lambert phase function to define the underlying planet population. Since this technique will be most relevant to future highly capable telescopes, we are using the HabEx IWA of 0.045 arc sec, OWA of 6 arc sec, and upper $Δ mag$ limit of 25 on a star that is 10 pc away.⁸ This resulting instrument visibility limits and planet population cover a well populated region of the $Δ mag$ versus $s$ joint probability density function [see Fig. 1 in Ref. 1 for a comparable example of the joint probability density function]. Convergence of the two completeness methods is demonstrated in Fig. 13.

Both methods presented in Fig. 13 demonstrate similar percent convergence to their own converged means. Completeness repeatedly calculated using the method in Sec. 2.6 with $10^{5}$ planets has a mean of 0.25785 with a standard deviation of 0.0010 for $4.28 \times 10^{8}$ planets. In comparison, the Brown completeness method converges to 0.25783. This error in the converged mean of the two methods is consistent with the completeness error identified in Fig. 14. We expect that completeness calculated using the method in this paper will have comparatively better convergence in sparsely populated regions of the $Δ mag$ versus $s$ joint probability density function (see Fig. 1 in Ref. 1).

Since completeness for an individual planet is the fraction of time that a planet is visible to the instrument, we validate individual planet visibility windows using test points. To calculate a ground truth fraction of time that a planet is visible, we create $10^{5}$ test points, evenly spaced in time, and calculate both the planet-star separation and planet-star $Δ mag$ at each of theose test points. We then determine whether each point is within the visibility limits of the instrument. The fraction of visible points is the fraction of time that a planet is visible. We repeated this for 25,000 planets due to computational cost of the test point method. The histogram of error in completeness for individual planets calculated using the test point method and method in this paper is shown in Fig. 14. The maximum error in completeness is observed to be near $5 \times 10^{- 5}$ , which is consistent with the converged error between Brown completeness and integration time adjusted completeness.

The test point method that we used to validate integration time adjusted completeness calculations took over 52 h on only 25,000 planets (less than a quarter of the planets needed to compute completeness using integration time adjusted completeness). The computation time for the integration time adjusted completeness calculation on $10^{5}$ planets is 21.34 s with a standard deviation of 0.4 s over 1000 calculations. Brown completeness, including the joint probability density function generation, has a mean execution time of 3.19 s and standard deviation of 0.13 s over 1000 calculations. This execution time of all of these methods scales linearly with the number of planets. Some room for optimization exists in integration time adjusted completeness.

3.4.

Completeness Versus Integration Time

To evaluate the effects of integration time on completeness, we calculate integration time adjusted completeness for various integration times, star distances, and planet populations (SAG13 and Earth-like planet population). Figure 15 shows the integration time adjusted completeness for the SAG13 planet population¹ and the Earth-like planet subpopulation defined in Appendix L. The decrease in completeness for longer integration times is most prominent for nearby stars for both populations. As expected, longer integration times decrease completeness.

For the assumed observatory parameters, the Brown completeness of an Earth-like population for a target 5 pc away is 0.583. As we increase integration time to 1, 2, and 5 days, the associated integration time adjusted completeness decreases by 0.63%, 1.27%, and 3.15%, respectively. Recalculating Brown completeness for a target at 25 pc in Fig. 15(b) gives 0.583. As we increase integration time to 1, 2, and 5 days on this 25 pc target, the associated integration time adjusted completeness decreases by 0.989%, 1.97%, and 4.92%, respectively. Our integration time adjustment of completeness shows Brown completeness overestimates exoplanet yields for any observation. Brown completeness applied to stars farther away results in a substantial overestimation of exoplanet yield.

Technically, the $Δ {mag}_{\lim}$ used as the upper limit for calculating completeness is a function of integration time and will approach a theoretical upper limit as the noise floor is reached. By calculating $Δ {mag}_{\lim}$ using Eq. (12) of Ref. 1, we can calculate completeness at each integration time to get the completeness versus integration time curve. We created a 4-m telescope as in Appendix M and evaluated completeness at each integration time to create an example demonstrating how Brown completeness and integration time adjusted completeness vary with integration time in Fig. 16. We can see that the completeness of both methods tracks very closely until they diverge. Brown completeness approaches its asymptotic limit as expected, whereas integration time adjusted completeness approaches a maximum and begins decreasing past $\sim 0.2 day$ .

3.5.

Dynamic Completeness

We compute dynamic completeness for the example test case in Fig. 1 of Ref. 13 using the method described in this paper and replicating the approach of the original work. This computation is only for a second epoch. We replicated the original work by sampling a large number of planets, finding the planets initially visible, propagating these planets to some time past initial observation, and determining what fraction of these planets were not detected the first time but were detected the second time. We replicated Brown’s work using both the Lambert phase function (orange line in Fig. 17) to match the original work as well as the quasi-Lambert phase function (red line in Fig. 17).

We compared the dynamic completeness computation time for 1000 dynamic completeness calculations testing 1000 individual points in time past the initial observation. Brown’s dynamic completeness had an average execution time of 34.62 s with a standard deviation of 0.99 s. Calculating dynamic completeness using the method presented in this work resulted in an average execution time of 4.757 s with a standard deviation of 0.036 s. The dynamic completeness calculation presented in this paper is $\approx 7 \times$ faster than the traditional Brown completeness method, enabling its use in optimization.

4. Discussion

4.1.

Convergence

Both methods of calculating completeness suffer from sparse sampling of the exoplanet population parameter space. Brown completeness is most affected as it needs planets sampled over $Ω$ , $ω$ , $a$ , $e$ , $i$ , $ν$ , $p$ , and $R$ . The current implementation of our method needs the same parameters except for $ν$ . We are therefore surprised that Brown’s method and our method share such similar completeness convergence. We hypothesize that the similar convergence results are because the integration bounds used for calculating completeness are over a dense region of the joint probability density function. If we examined a more sparse region of the joint probability density function such as high $Δ mag$ or larger separation (where the larger, lower occurrence rate planets occur), we expect our completeness methods to have marginally better convergence for low numbers of planets. Regardless, integration time adjusted completeness easily reaches completeness errors below $10^{- 3}$ (translates to $\approx 0.8 %$ where $C = 0.12$ , the largest completeness of an observation in Ref. 1), which is a sufficient error for use in optimization.

4.2.

Reducing Parameter Spaces

The integration time adjusted completeness method implemented in EXOSIMS is not cached for computational efficiency like the Brown completeness is. Although we are able to store the joint probability density function of $Δ mag$ versus $s$ with Brown’s method, integration time completeness requires the KOE of each simulated planet to be stored. The curse of dimensionality prohibits us from finely sampling the entire subspace of exoplanets, making the caching of integration time adjusted completeness prohibitive, but not all our methods require all 6 KOE, $p$ , and $R$ . So, in the future, it may be possible to create a representative subsampling of planets weighted by occurrence rate based on the parameters needed for different calculations. For example, $s$ extrema, $s_{WA}$ -orbit intersections, and $Δ mag$ extrema calculations only require $ω$ , $a$ , $e$ , and $i$ (Note that the $ν$ locations of the $Δ mag$ extrema are independent of $p$ and $R$ , but the magnitude of $Δ mag$ extrema are dependent upon them). The photometric property $p \times R$ could then be sampled and independently combined with individual ( $ω$ , $a$ , $e$ , $i$ ) combinations. The primary benefit is that, given some instrument parameters and star properties, we could determine the subset of parameters that can physically be observed prior to calculating completeness, thus making the per-star completeness calculations more efficient.

4.3.

Limitations

The approach we implemented in this paper did not vary the $Δ {mag}_{\lim}$ as we varied the integration time. If properly executed, $Δ {mag}_{\lim}$ should increase with increased integration time. However, calculating the planet true anomaly intersections with a given $Δ {mag}_{\lim}$ is the most expensive calculation, so repeatedly calculating this to convergence is not desirable.

The $Δ {mag}_{\lim}$ of coronagraph designs are working angle dependent and therefore separation dependent, but we assume that this limit is a fixed quantity. If there is substantial variation in the limiting $Δ mag$ , it is better to use the conservative value with this method. If a higher level of predictive accuracy is desired, these methods could be used to ascertain if any intersection with the $Δ {mag}_{\lim}$ exists and subsequently another method could be used on the subset of planets with known intersections.

4.4.

Impact of the Integration Time Adjusted Completeness

As we showed in this paper, integration time adjusted completeness and Brown completeness converge to the same value to within 0.00002 when evaluated at $t_{\max} = 0$ day. The integration time adjusted completeness for an observation of an Earth-like planet population on a 1- $M_{⊙}$ star 5 pc away and integration time of $\approx 1$ day is 0.64% lower than the comparable Brown completeness. For reference, the Roman target list in Table 9 of Ref. 1 has a maximal integration time of 1.71 days but most are $< 0.6$ day. This means that we could expect an average reduction in yield below 0.64% when observing a population of Earth-like planets around Sun-like stars. This integration time completeness adjustment is within the 3.19%, $3 σ$ , margin of error from 1000 Monte Carlo simulations of the cycle 6 Roman in Ref. 1 and cannot be considered statistically significant. However, a $\sim 4$ -day integration time on a Sun-like star 25 pc away observing an Earth-like population has a decrease in completeness above the Ref. 1 $3 σ$ margin of error. However, the threshold of statistical significance should not deter the widespread use of integration time adjusted completeness as the adoption of this method can shore up the differences between completeness-based yields and simulation-based yields. We can also say that integration time has a muted effect when observing stars that are farther away. This can most likely be attributed to the increase in $s_{\min}$ and resulting decrease in the total time-fraction that smaller period planets spend within the visible limits of the telescope. Big planets with the orbital radius of Jupiter will move more slowly and be less affected by integration times.

In Appendix M, we optimized an exoplanet direct imaging mission maximizing single-visit Brown completeness yield and included the resulting Design Reference Mission (DRM) in Table 2. The resulting DRM has a Brown completeness yield of 387.16 exoplanet detections in a single-visit detection survey on average, but the integration time adjusted completeness yield expects only 354.67 exoplanets on average. This means that Brown completeness for the particular instrument parameters optimized over the SAG13 planet population with the particular mission parameters overpredicts the actual exoplanet yield by 32.49 exoplanets on average, 9.16% more than the integration time adjusted completeness yield. When calculating completeness, we are careful to scale the completeness limits of integration by the star’s luminosity. For the integration time adjusted completeness method, we additionally scale the planet’s orbital periods based off the mass of the host star. The percentage difference above and beyond that expected from Fig. 15 can be found by looking at these two adjustments applied to each star and the target list in Table 2. The average star in the target list has a larger mass and a brighter luminosity than that of the Sun. The brighter luminosity results in a smaller $s_{\min}$ and $s_{\max}$ , which serves to incorporate more smaller semimajor axis planets into the completeness calculations, and the larger star masses result in shorter periods, meaning that the planet visibility windows all decrease in duration.

Integration time adjusted completeness is crucial for determining the ability of an Earth-like planet to be spectrally characterized. A spectral characterization with a coronagraph could take between a few days and 60 days. Because of the long integration time, the planners of a HabEx use the “characterization completeness” of 10% and maximum integration time of 60 days as a filter on stars to consider observing.⁸ The calculation of this critical filter could be substantially improved by considering integration time adjusted completeness if the spectrum of the planet must be taken all at once (i.e., observation spanning multiple weeks) and not spread across multiple epochs. Just because a planet can be detected does not mean that it can be spectrally characterized.

4.5.

Dynamic Completeness and Computation Cost

The greatest benefit of using the methods in this paper to calculate completeness is the marginal additional cost of calculating dynamic completeness. Generally, dynamic completeness requires the computation of true anomalies from time, which is prohibitive for $> 10^{5}$ orbits at $> 100$ different times in the future. With our methods, calculating visibility windows allows us to use boolean operations to compute dynamic completeness about $7 \times$ faster than Brown’s dynamic completeness method. The computational cost of calculating completeness using the method described in this paper is an order of magnitude larger than Ref. 2, but it can use orders of magnitude fewer planets to do so. Unlike Brown’s method, within the computation time of our method, we also get additional desirable information about the detectable planets such as their $s$ extrema and $Δ mag$ extrema.

4.6.

Revisiting the Same Exoplanet

Another limitation present in the planning of exoplanet direct-imaging missions is the telescope keep-out angles. Figure 3 of Ref. 1 contains a keep-out map for a subset of the DRM created in that paper. The smallest percent of time that a target star is visible for the Roman is nominally 28%. Due to symmetry of the keep-out region, this translates into two separate time windows of visibility of $\sim 51$ days. Instead of considering the integration time an integration time input, it can also be considered a revisit time input for determining the probability of being able to observe a planet twice in the same target-star visibility window.

4.7.

Exoplanet Classification

The underlying methods in this paper are used to find locations along a planet’s orbit where $s$ and $Δ mag$ intersections occur. The methods in this paper, with some modification, can also be used to probabilistically classify an exoplanet subtype.⁴ If an exoplanet is detected with a particular ( $s$ , $Δ mag$ ) and an uncertainty region of $s \pm σ_{s}$ and $Δ mag \pm σ_{Δ mag}$ , then the methods in this paper can be applied to each of these four bounding lines. By finding the average orbital time-fraction that exoplanets of each type spend in the bounding uncertainty box, we can find the probability that the exoplanet belongs to a specific exoplanet subtype. This requires additional work beyond the scope of this paper.

5. Conclusion

We have demonstrated an accurate method for calculating integration time adjusted completeness and its adaptation to calculating dynamic completeness. In the process, we also created fast and accurate methods for calculating the true anomalies where a planet’s orbit has specific values of projected separation $Δ mag$ and their extrema. We demonstrated how to use these methods to calculate a more accurate integration time adjusted completeness using the fractions of time that a planet is detectable by an instrument. We demonstrate that traditional methods of calculating completeness overestimate the number of observable planets because they do not subtract the integration time used in observing the target. For a Sun-like star at 25 pc with 1-day and 5-day integration times, integration time adjusted completeness of Earth-like planets is reduced by 1% and 5%, respectively. We applied integration time adjusted completeness to a target list optimized using the Brown completeness method and found that Brown completeness overestimated yields by 9.61%. We also demonstrated that our methods can be used to quickly calculate dynamic completeness for determining the optimal time to revisit a target star.

6. Appendix A: Common Notation

There are many common variable forms used in this paper that are simply summarizable. Any variable $\underline{x}$ refers to a 3D vector in $X$ , $Y$ , and $Z$ coordinates of something. We use $x$ to indicate an array of variables, specifically used when referring to multiple roots of a polynomial. $| \overline{A B} |$ refers to the length of line segment $\overline{A B}$ . Line segment $\overline{A B}$ is treated as a vector. $C^{'}$ is a projection of point $C$ .

There are multiple different subscripts with different meanings used in this paper. A subscript with $x_{i}$ refers to the index of a host star (out of some whole, non-descript, star catalog). A subscript $x_{k}$ refers to an individual planet out of the whole large set of planets. Subscripts of $x_{\min}$ , $x_{\max}$ , $x_{lmin}$ , and $x_{lmax}$ are descriptors on the individual variable $x$ that indicate that the variable is associated with the minimum, maximum, local minimum, or local maximum, respectively. When solving for the distance between points inside an ellipse and the vertices or co-vertices of that ellipse with semimajor axis $a$ , we use the shorthand notation of the form $s_{a + x, y}$ . [This is specifically the distance between point ( $x$ , $y$ ) and ( $- a$ , 0).]

We denote coefficients of the polynomials of the four methods that solve the quartic in this paper as $A_{#}$ .

In Fig. 4, we reference many points on the 3D elliptical orbit and 2D projection of this orbit into the plane of the sky. Points on the original 3D ellipse are labeled $B$ , $H$ , $P$ , $C$ , $K$ , $O$ , $D$ , $A$ , $B^{'}$ , $H^{'}$ , $P^{'}$ , $C^{'}$ , $K^{'}$ , $O^{'}$ , $D^{'}$ , and $A^{'}$ .

In Appendix K, we use $p_{0}$ to $p_{11}$ as intermediate constants for simplifying the full quartic expression. $P$ , $D$ , and $Δ$ (by itself and only in this section) are intermediate constants derived from quartic coefficients for determining the sign of the quartic solutions.

$K_{0}$ to $K_{9}$ are intermediate simplifying coefficients used in representing Eqs. (13)–(15) in Sec. 2.1.1. The full expansion of these equations is included in Appendix E.

Fig. 9

Planet-star separation of the planet is plotted in black. The separation extrema are indicated by diamonds, where the maximum is red, local maximum is yellow, local minimum is magenta, and minimum is teal. Lines are drawn at the minimum and maximum separations. For the input separation of $s = 1 AU$ , the green dots are the analytically calculated orbit intersections. (a) The separation versus true anomaly and (b) the separation versus time past periastron for a planet orbit with $a = 1.82 AU$ , $e = 0.09$ , $Ω = 3.37 rad$ , $ω = 4.86 rad$ , and $i = 1.25 rad$ .

Fig. 10

$s_{WA}$ -orbit intersection error histogram calculated for a separation of $s_{WA} = 1 AU$ for $10^{5}$ planet orbits generated from the SAG13 planet population. Of these $10^{5}$ planet orbits, 6201 orbits have two intersections with the $s_{WA}$ circle, and 17,952 have four intersections with the $s_{WA}$ circle. This results in a total of 84,210 planet-star intersections. After calculating and identifying true anomalies of intersections using the methods described in this paper, we evaluated the planet-star separation of each orbit at the true anomalies and calculated error from the input $s_{WA}$ .

Fig. 11

(a) The $Δ mag$ of a planet (black line) plotted against $ν$ and (b) time past periastron. $Δ mag$ extrema are indicated by diamonds, where the maximum is red, local maximum is yellow, local minimum is magenta, and minimum is teal. Separation minimum and maximum are indicated by horizontal lines. For an input $Δ {mag}_{\lim} = 25$ (green line), we calculated the specific true anomalies (green dots) where the planet’s $Δ mag$ intersects this line. The specific planet KOE are $a = 5.36 AU$ , $e = 0.56$ , $ω = 5.06 rad$ , $Ω = 0.69 rad$ , $i = 0.81 rad$ ), $p = 0.3$ , and $R = 4 R_{\oplus}$ .

Fig. 12

A histogram of $ν (Δ mag)$ error between the root solving method described in this paper and a root solver of a cubic spline fit method (purple) for 510, 120 planet orbits that produce two intersections (totaling 1,020,240 datapoints). We also compare the error between the root solving method described in this paper and an optimization method (blue) for $10^{4}$ planet subset of the planet orbits that produce intersections. The optimization method investigates fewer planets because it is orders of magnitude more computationally expensive, but it is a fundamentally different approach to validating our method than a root solver.

Fig. 13

Convergence of % error in completeness for varying numbers of planets used in its computation for Brown’s method (green) and the method presented in this work with $t_{\max} = 0$ day (purple). Completeness at the maximum number of planets is assumed to be the converged value of completeness of the respective methods. The converged mean of Brown’s method is 0.25783 compared with the converged mean of the method in this paper of 0.25785. The standard deviation of the method in this paper is 0.0010.

Fig. 14

Per planet completeness error histogram between the integration time adjusted completeness method with $t_{\max} = 0$ day and the test point method for 25,000 planets.

Fig. 15

Integration time adjusted completeness values of (a) SAG13 planet population and (b) Earth-like planet population with IWA = 0.045 arc sec, OWA = 6 arc sec, and $Δ mag = 25$ for varying integration times and star distances.

Fig. 16

Brown completeness and integration time adjusted completeness versus integration time for HIP 32,279 assuming a mass of 1.564 $M_{⊙}$ , luminosity of 7.12 $L_{⊙}$ , and distance of $d_{i} = 3.51$ pc. The assumed telescope is a 4-m monolith with IWA = 0.045 arc sec and OWA = 6 arc sec and $Δ {mag}_{\lim}$ computed using the instrument noise model as in Ref. 1.

Fig. 17

Dynamic completeness assuming the same parameters as in Ref. 13 ( $Δ mag = 26$ , $M_{HIP 29271} = 1.103 M_{⊙}$ , $d_{i} = 10.215$ pc, OWA = 600 arc sec, IWA = 0.075 arc sec, $0.7 \sqrt{L} \leq a \leq 1.5 \sqrt{L}$ , $0 \leq e \leq 0.35$ , $p = 0.26$ , $R = 1 R_{\oplus}$ , and $L = 0.83 L_{⊙}$ ). In this work (blue), we assume a quasi-Lambert phase function, but in Ref. 13, a Lambert phase function is assumed. Work was done to replicate Brown’s original work (Brown Lambert, orange line) and replicate Brown’s work using the quasi-Lambert phase function (Brown quasi-Lambert, red line).

7. Appendix B: True Anomalies of $s$ Intersection Process

The process outlined here is used in Sec. 2.1.2.

1. Define KOE.
2. ${\underset{̲}{r}}_{k / i}$ from KOE.
3. Orbiting foci $F$ from KOE.
4. Calculate 3D orbit ellipse center.
5. Project 3D ellipse to 2D ellipse.
- (a) Prove that 3D ellipses project to 2D ellipses.
- (b) Project 3D ellipse center to 2D ellipse center.
- (c) Project ellipse orbiting foci into plane of the sky.
- (d) Project semiminor axis line of 3D ellipse to 2D ellipse.
- (e) Project semimajor axis line of 3D ellipse to 2D ellipse.
- (f) Calculate $\overline{Q Q^{'}}$ construction line from semiminor axis and semiminor axis.
- (g) Calculate projected ellipse semimajor axis and semiminor axis from $\overline{O^{'} Q^{'}}$ and $\overline{O^{'} Q}$ .
- (h) Calculate projected semimajor axis angular offset from $X$ axis of plane of the sky.
6. Derotate the projected ellipse.
7. Center the projected ellipse.
8. Find minimum, maximum, local minimum, and local maximum of projected planet-star separation.
- (i) Formulate projected planet-star separation equation.
- (j) Reformat and set $δ s^{2} / δ x_{e}$ equal to 0.
- (k) Combine coefficients into standard quartic form.
- (l) Use standard general quartic solutions.
- (m) Take the absolute value and only real component of $x_{e}$ solutions.
- (n) Calculate $y_{e}$ associated with each solution.
- (o) Assign solutions to minimum, maximum, local minimum, and local maximum planet-star separations.
  - i. For all solutions (all real and only 2 real):
    A. Calculate $s_{\min}$ from $x_{1}$ and assign $(x_{\min}, y_{\min})$ of quadrant 1.
    B. Calculate $s_{\max}$ from $x_{0}$ and assign $(x_{\max}, y_{\max})$ of quadrant 1.
  - ii. All real solutions (where all solutions have $I (y) < 10^{- 5}$ , additionally assign).
    C. The larger of $s_{(x_{2})}$ and $s_{(x_{3})}$ becomes $s_{lmax}$ and assign $(x_{lmax}, y_{lmax})$ of quadrant 1.
    D. The smaller of $s_{(x_{2})}$ and $s_{(x_{3})}$ becomes $s_{lmin}$ and assign $(x_{lmin}, y_{lmin})$ of quadrant 1.
  - iii. Assign solution signs to proper quadrants.
9. Find circle and projected derotated centered ellipse intersection points.
- (p) Formulate circle-ellipse intersection equation.
- (q) Reformat into standard quartic form.
- (r) Use standard general quartic solutions.
- (s) Classify intersection solutions.
  - iv. Inside outer separation bounds, $s_{\min} < s < s_{\max}$ .
  - v. Inside local min/max separation bounds, $s_{lmin} < s < s_{lmax}$ .
  - vi. Outside outer separation bounds (no intersections).
- (t) Calculate intersection solutions for planets with $s_{\min}$ , $s_{\max}$ , $s_{lmin}$ , and $s_{lmax}$ .
  - vii. Two intersections on same $y$ side of ellipse as the star, $s_{\min} < s < s_{lmin}$ $x = {\begin{matrix} x_{3}, I (x_{1}) > 10^{- 10} \\ x_{1}, I (x_{1}) < 10^{- 10} \\ x_{0} \end{matrix}$ .
  - viii. Four intersections where, $s_{lmin} < s < s_{lmax}$
    • Calculate $Δ x = x_{k} - | x |$ , ( $k \in {0,1, 2,3}$ ).
    • Sort $Δ x_{k}$ from min to max and rearrange the associated $x_{k}$ to match that order.
    • Same $X Y$ points is $x_{3}$ .
    • Same $X$ opposite $Y$ is $x_{2}$ .
    • Opposite $X$ Same $Y$ is $x_{0}$ .
    • Opposite $X$ opposite $Y$ is $x_{1}$ .
  - ix. Two intersections on opposite $x$ side of ellipse as the star, $s_{lmax} < s < s_{\max}$ .
    • $x_{0}$ and $x_{1}$ .
    • $y_{1}$ is opposite sign of $| y |$ .
- (u) Calculate intersection solutions for planets with only $s_{\min}$ and $s_{\max}$ .
  - x. Calculate projected ellipse quadrant separation bounds $s_{x, b + y}$ , $s_{x, b - y}$ , $s_{a + x, y}$ , and $s_{a - x, y}$ .
  - xi. Identify star location type and quadrant order.
    • Type 0 occurs where $s_{x + a, y} < s_{x, y + b}$ . Smallest to largest order: $s_{x, b - y}$ , $s_{a - x, y}$ , $s_{x + a, y}$ , $s_{x, y + b}$ .
    • Type 1 occurs where $s_{x, y + b} < s_{a - x, y}$ . Smallest to largest order: $s_{x, b - y}$ , $s_{x, y + b}$ , $s_{a - x, y}$ , $s_{x + a, y}$ .
    • Type 2 occurs where $s_{a - x, y} < s_{x, y + b}$ and $s_{x, y + b} < s_{x + a, y}$ and $s_{x, b - y} < s_{a - x, y}$ . Smallest to largest order: $s_{x, b - y}$ , $s_{a - x, y}$ , $s_{x, y + b}$ , $s_{x + a, y}$ .
    • Type 3 occurs where $s_{a - x, y} < s_{x, b - y}$ . Smallest to largest order: $s_{a - x, y}$ , $s_{x, b - y}$ , $s_{x, y + b}$ , $s_{x + a, y}$ .

8. Appendix C: $ν$ from $Δ mag$ Extrema process

Here we present the full outline of our process for calculating true anomalies where the planet’s orbit has a $Δ mag$ extrema. These steps are discussed in Sec. 2.2.

1. Substitute into Eq. (1) components and expand until all terms are functions of KOE.
2. Manipulate the $Δ mag$ equation until all $ν$ terms are isolated on the right side.
3. Take the derivative.
4. Replace $\sin (ν)$ with $\sqrt{1 - \cos {(ν)}^{2}}$ .
5. Replace all $\cos (ν)$ with $x$ .
6. Isolate all square root terms to general form $A^{2} \sqrt{1 - x^{2}} = B^{2}$ .
7. Square both sides.
8. Find polynomial coefficients of $x$ .
9. Solve with root solver.
10. Filter out $x$ solutions with imaginary components ( $> 10^{- 7}$ ) ( $10^{- 8}$ is approximately $\sqrt{10^{- 16}}$ , machine precision).
11. Filter out invalid solutions, $1 < x$ or $- 1 > x$ .
12. Compute $ν_{0} = \cos^{- 1} (x)$ and $ν_{1} = 2 π - ν_{0}$ .
13. Compute associated ${Δ mag}_{0}$ and ${Δ mag}_{1}$ .
14. Remove solutions that are not extrema.
15. Assign solutions $Δ {mag}_{\max} = \min ({Δ mag}_{0}, {Δ mag}_{1})$ .
16. Assign solutions $Δ {mag}_{\min} = \max ({Δ mag}_{0} {, Δ mag}_{1})$
17. Remove duplicate solutions.
18. Remove assigned solutions.
19. Check if solution is extrema.
20. Assign $Δ {mag}_{lmin}$ and $Δ {mag}_{lmax}$ .

9. Appendix D: $ν$ from ${Δ mag}_{\lim}$ Intersection Process

Here we present the full outline of our process for calculating true anomalies where the planet’s orbit has a specified $Δ mag$ , which we call an intersection. These steps are discussed in Sec. 2.2.

1. Find $Δ mag$ extrema.
2. Find planets where $Δ mag < Δ {mag}_{\max}$ , $Δ mag > Δ {mag}_{\min}$ , and $Δ mag > Δ {mag}_{lmax}$ or $Δ mag < Δ {mag}_{lmin}$ (these planets have only two intersections).
3. Substitute into $Δ mag$ equation components and expand until all terms are functions of KOE.
4. Manipulate the equation until all $ν$ terms are isolated on the right side.
5. Expand.
6. Replace $\sin (ν)$ with $\sqrt{1 - \cos {(ν)}^{2}}$ .
7. Replace all $\cos (ν)$ with $x$ .
8. Isolate all square root terms to general form $A^{2} \sqrt{1 - x^{2}} = B^{2}$ .
9. Square both sides.
10. Expand.
11. Find polynomial coefficients of $x$ .
12. Solve with root solver.
13. Filter out solutions with imaginary components ( $> 10^{- 7}$ ).
14. Filter out solutions $1 < x_{k}$ or $- 1 > x_{k}$ .
15. Take arrays $ν_{0} = \cos^{- 1} (x)$ and $ν_{1} = 2 π - ν_{0}$ .
16. Solve for associated arrays $Δ {mag}_{0}$ and $Δ {mag}_{1}$ .
17. Remove solutions where $| Δ {mag}_{0} - Δ mag | > 0.01$ and $| Δ {mag}_{1} - Δ mag | > 0.01$ (Note that 0.01 is $< \pm 0.08 %$ error on dmag).
18. Remove duplicate solutions.
19. Verify that there are only two viable solutions.
20. Assign solutions to $ν$ .

10. Appendix E: Projected Ellipse Semimajor Axis and Semiminor Axis, and $θ$

We write the analytical expressions for $a_{p}$ , $b_{p}$ , and $θ$ using full expansions of Eqs. (13)–(15). These expressions are far too long to be practically conveyed, so we use the following intermediary parameters here and only here to simplify these expressions:

Eq. (61)

K_{0} = e (1 - e^{2}) K_{1} = \sin (Ω) \cos (ω) + \sin (ω) \cos (Ω) \cos (i) K_{2} = \sqrt{\frac{e + 1}{1 - e}} K_{3} = a (1 - e^{2}) K_{4} = - \sin (Ω) \sin (ω) \cos (i) + \cos (Ω) \cos (ω) K_{5} = \sqrt{K_{1}^{2} K_{3}^{2} + K_{4}^{2} K_{3}^{2} + K_{3}^{2} \sin^{2} (i) \sin^{2} (ω)} K_{6} = K_{3} (- \sin (Ω) \cos (ω) - \sin (ω) \cos (Ω) \cos (i)) / (1 - e) K_{7} = K_{3} (\sin (Ω) \sin (ω) \cos (i) - \cos (Ω) \cos (ω)) / (1 - e) K_{8} = K_{3} (\sin (Ω) \cos (ω + 2 \tan^{- 1} (K_{2})) + \sin (ω + 2 \tan^{- 1} (K_{2})) \cos (Ω) \cos (i)) / (e \cos (2 \tan^{- 1} (K_{2})) + 1) K_{9} = K_{3} (- \sin (Ω) \sin (ω + 2 \tan^{- 1} (K_{2})) \cos (i) + \cos (Ω) \cos (ω + 2 \tan^{- 1} (K_{2}))) / (e \cos (2 \tan^{- 1} (K_{2})) + 1) .

The simplified expression for the semimajor axis of the projected ellipse is

Eq. (62)

a_{p} = \frac{1}{2} \sqrt{{(\frac{K_{0} a^{2} K_{1}}{K_{5}} - \frac{K_{0} a^{2} K_{4}}{K_{5}} + K_{6} - K_{9})}^{2} + {(\frac{K_{0} a^{2} K_{1}}{K_{5}} + \frac{K_{0} a^{2} K_{4}}{K_{5}} + K_{7} + K_{8})}^{2}} + \frac{1}{2} \sqrt{{(\frac{K_{0} a^{2} K_{1}}{K_{5}} - \frac{K_{0} a^{2} K_{4}}{K_{5}} - K_{7} + K_{8})}^{2} + {(\frac{K_{0} a^{2} K_{1}}{K_{5}} + \frac{K_{0} a^{2} K_{4}}{K_{5}} + K_{6} + K_{9})}^{2}} .

The simplified expression for the projected ellipse semiminor axis is

Eq. (63)

b_{p} = - \frac{1}{2} \sqrt{{| \frac{K_{0} a^{2} K_{1}}{K_{5}} - \frac{K_{0} a^{2} K_{4}}{K_{5}} + K_{6} - K_{9} |}^{2} + {| \frac{K_{0} a^{2} K_{1}}{K_{5}} + \frac{K_{0} a^{2} K_{4}}{K_{5}} + K_{7} + K_{8} |}^{2}} + \frac{1}{2} \sqrt{{| \frac{K_{0} a^{2} K_{1}}{K_{5}} - \frac{K_{0} a^{2} K_{4}}{K_{5}} - K_{7} + K_{8} |}^{2} + {| \frac{K_{0} a^{2} K_{1}}{K_{5}} + \frac{K_{0} a^{2} K_{4}}{K_{5}} + K_{6} + K_{9} |}^{2}} .

The simplified expression for the angle between the projected ellipse semimajor axis and

X

axis is

Eq. (64)

θ = \frac{1}{2} [\tan^{- 1} (\frac{\frac{K_{0} a^{2} K_{1}}{K_{5}} - \frac{K_{0} a^{2} K_{4}}{K_{5}} + K_{6} - K_{9}}{\frac{K_{0} a^{2} K_{1}}{K_{5}} + \frac{K_{0} a^{2} K_{4}}{K_{5}} + K_{7} + K_{8}}) + \tan^{- 1} (\frac{\frac{K_{0} a^{2} K_{1}}{K_{5}} + \frac{K_{0} a^{2} K_{4}}{K_{5}} + K_{6} + K_{9}}{- \frac{K_{0} a^{2} K_{1}}{K_{5}} + \frac{K_{0} a^{2} K_{4}}{K_{5}} + K_{7} - K_{8}})] .

11. Appendix F: $s$ Extrema Polynomial

We reduced Eq. (22) from Sec. 2.1.2 into the standard reduced form of a quartic expression:

Eq. (65)

0 = x_{e}^{4} + \frac{- 2 a_{p}^{2} x_{*} + 2 b_{p}^{2} x_{*}}{(a_{p}^{4} - 2 a_{p}^{2} b_{p}^{2} + b_{p}^{4}) / a_{p}^{2}} x_{e}^{3} + \frac{- a_{p}^{4} + 2 a_{p}^{2} b_{p}^{2} + a_{p}^{2} x_{*}^{2} - b_{p}^{4} + b_{p}^{2} y_{*}^{2}}{(a_{p}^{4} - 2 a_{p}^{2} b_{p}^{2} + b_{p}^{4}) / a_{p}^{2}} x_{e}^{2} + \frac{2 a_{p}^{4} x_{*} - 2 a_{p}^{2} b_{p}^{2} x_{*}}{(a_{p}^{4} - 2 a_{p}^{2} b_{p}^{2} + b_{p}^{4}) / a_{p}^{2}} x_{e} + \frac{- a_{p}^{4} x_{*}^{2}}{(a_{p}^{4} - 2 a_{p}^{2} b_{p}^{2} + b_{p}^{4}) / a_{p}^{2}} .

The coefficients of the quartic expression are

Eq. (66)

A_{0} = \frac{- 2 a_{p}^{2} x_{*} + 2 b_{p}^{2} x_{*}}{(a_{p}^{4} - 2 a_{p}^{2} b_{p}^{2} + b_{p}^{4}) / a_{p}^{2}},

Eq. (67)

B_{0} = \frac{- a_{p}^{4} + 2 a_{p}^{2} b_{p}^{2} + a_{p}^{2} x_{*}^{2} - b_{p}^{4} + b_{p}^{2} y_{*}^{2}}{(a_{p}^{4} - 2 a_{p}^{2} b_{p}^{2} + b_{p}^{4}) / a_{p}^{2}},

Eq. (68)

C_{0} = \frac{2 a_{p}^{4} x_{*} - 2 a_{p}^{2} b_{p}^{2} x_{*}}{(a_{p}^{4} - 2 a_{p}^{2} b_{p}^{2} + b_{p}^{4}) / a_{p}^{2}},

Eq. (69)

D_{0} = \frac{- a_{p}^{4} x_{*}^{2}}{(a_{p}^{4} - 2 a_{p}^{2} b_{p}^{2} + b_{p}^{4}) / a_{p}^{2}} .

12. Appendix G: $s$ Intersection Polynomial

After dividing by the coefficient of $x_{e}^{4}$ in Eq. (35) from Sec. 2.1.3, we get

Eq. (70)

0 = x_{e}^{4} - \frac{4 a_{p}^{2} x_{*}}{a_{p}^{2} - b_{p}^{2}} x_{e}^{3} + \frac{2 a_{p}^{2} (a_{p}^{2} b_{p}^{2} - a_{p}^{2} s^{2} + 3 a_{p}^{2} x_{*}^{2} + a_{p}^{2} y_{*}^{2} - b_{p}^{4} + b_{p}^{2} s^{2} - b_{p}^{2} x_{*}^{2} + b_{p}^{2} y_{*}^{2})}{a_{p}^{4} - 2 a_{p}^{2} b_{p}^{2} + b_{p}^{4}} x_{e}^{2} + \frac{4 a_{p}^{4} x_{*} (- b_{p}^{2} + s^{2} - x_{*}^{2} - y_{*}^{2})}{a_{p}^{4} - 2 a_{p}^{2} b_{p}^{2} + b_{p}^{4}} x_{e} + \frac{a_{p}^{4} (b_{p}^{4} - 2 b_{p}^{2} s^{2} + 2 b_{p}^{2} x_{*}^{2} - 2 b_{p}^{2} y_{*}^{2} + s^{4} - 2 s^{2} x_{*}^{2} - 2 s^{2} y_{*}^{2} + x_{*}^{4} + 2 x_{*}^{2} y_{*}^{2} + y_{*}^{4})}{a_{p}^{4} - 2 a_{p}^{2} b_{p}^{2} + b_{p}^{4}} .

The coefficients of the general quartic expression in Eq. (36) are

Eq. (71)

A_{1} = - \frac{4 a_{p}^{2} x_{*}}{a_{p}^{2} - b_{p}^{2}},

Eq. (72)

B_{1} = \frac{2 a_{p}^{2} (a_{p}^{2} b_{p}^{2} - a_{p}^{2} s^{2} + 3 a_{p}^{2} x_{*}^{2} + a_{p}^{2} y_{*}^{2} - b_{p}^{4} + b_{p}^{2} s^{2} - b_{p}^{2} x_{*}^{2} + b_{p}^{2} y_{*}^{2})}{a_{p}^{4} - 2 a_{p}^{2} b_{p}^{2} + b_{p}^{4}},

Eq. (73)

C_{1} = \frac{4 a_{p}^{4} x_{*} (- b_{p}^{2} + s^{2} - x_{*}^{2} - y_{*}^{2})}{a_{p}^{4} - 2 a_{p}^{2} b_{p}^{2} + b_{p}^{4}},

and

Eq. (74)

D_{1} = \frac{a_{p}^{4} (b_{p}^{4} - 2 b_{p}^{2} s^{2} + 2 b_{p}^{2} x_{*}^{2} - 2 b_{p}^{2} y_{*}^{2} + s^{4} - 2 s^{2} x_{*}^{2} - 2 s^{2} y_{*}^{2} + x_{*}^{4} + 2 x_{*}^{2} y_{*}^{2} + y_{*}^{4})}{a_{p}^{4} - 2 a_{p}^{2} b_{p}^{2} + b_{p}^{4}} .

13. Appendix H: $Δ mag$ Extrema Polynomial

Here we include the coefficients of Eq. (49) from Sec. 2.2 as

Eq. (75)

A_{2} = e^{4} \sin^{4} (i), B_{2} = 3 e^{3} (e \sin (ω) + \sin (i)) \sin^{3} (i), C_{2} = \frac{1}{4} e^{2} (4 e^{2} \sin^{2} (i) \sin^{2} (ω) - 8 e^{2} \sin^{2} (i) + 8 e^{2} \sin^{2} (ω) + 5 e^{2} + 34 e \sin (i) \sin (ω) + 13 \sin^{2} (i)) \sin^{2} (i), D_{2} = \frac{1}{2} e [2 e^{3} \sin^{2} (i) \sin^{3} (ω) - 8 e^{3} \sin^{2} (i) \sin (ω) + 3 e^{3} \sin (ω) + 5 e^{2} \sin^{3} (i) \sin^{2} (ω) - 11 e^{2} \sin^{3} (i) + 10 e^{2} \sin (i) \sin^{2} (ω) + 7 e^{2} \sin (i) + 17 e \sin^{2} (i) \sin (ω) + 3 \sin^{3} (i)) \sin (i), E_{2} = - \frac{5 e^{4} \sin^{4} (i) \sin^{2} (ω)}{4} + \frac{5 e^{4} \sin^{4} (i)}{4} - \frac{7 e^{4} \sin^{2} (i) \sin^{2} (ω)}{4} - \frac{3 e^{4} \sin^{2} (i)}{2} + \frac{e^{4}}{4} + \frac{3 e^{3} \sin^{3} (i) \sin^{3} (ω)}{2} - 10 e^{3} \sin^{3} (i) \sin (ω) + \frac{7 e^{3} \sin (i) \sin (ω)}{2} + 2 e^{2} \sin^{4} (i) \sin^{2} (ω) - \frac{21 e^{2} \sin^{4} (i)}{4} + 4 e^{2} \sin^{2} (i) \sin^{2} (ω) + \frac{7 e^{2} \sin^{2} (i)}{2} + \frac{7 e \sin^{3} (i) \sin (ω)}{2} + \frac{\sin^{4} (i)}{4}, F_{2} = - \frac{3 e^{4} \sin^{3} (i) \sin^{3} (ω)}{2} + \frac{3 e^{4} \sin^{3} (i) \sin (ω)}{2} - \frac{3 e^{4} \sin (i) \sin (ω)}{2} - \frac{e^{3} \sin^{4} (i) \sin^{4} (ω)}{2} - \frac{5 e^{3} \sin^{4} (i) \sin^{2} (ω)}{2} + 3 e^{3} \sin^{4} (i) - 5 e^{3} \sin^{2} (i) \sin^{2} (ω) - \frac{7 e^{3} \sin^{2} (i)}{2} + \frac{e^{3}}{2} - \frac{e^{2} \sin^{3} (i) \sin^{3} (ω)}{2} - 8 e^{2} \sin^{3} (i) \sin (ω) + \frac{5 e^{2} \sin (i) \sin (ω)}{2} + \frac{e \sin^{4} (i) \sin^{2} (ω)}{2} - 2 e \sin^{4} (i) + e \sin^{2} (i) \sin^{2} (ω) + \frac{3 e \sin^{2} (i)}{2} + \frac{\sin^{3} (i) \sin (ω)}{2}, G_{2} = - \frac{e^{4} \sin^{4} (i) \sin^{4} (ω)}{4} + \frac{e^{4} \sin^{4} (i) \sin^{2} (ω)}{2} - \frac{e^{4} \sin^{4} (i)}{4} - \frac{e^{4} \sin^{2} (i) \sin^{2} (ω)}{2} + \frac{e^{4} \sin^{2} (i)}{2} - \frac{e^{4}}{4} - \frac{7 e^{3} \sin^{3} (i) \sin^{3} (ω)}{2} + \frac{7 e^{3} \sin^{3} (i) \sin (ω)}{2} - \frac{7 e^{3} \sin (i) \sin (ω)}{2} - \frac{5 e^{2} \sin^{4} (i) \sin^{4} (ω)}{4} - e^{2} \sin^{4} (i) \sin^{2} (ω) + \frac{9 e^{2} \sin^{4} (i)}{4} - 5 e^{2} \sin^{2} (i) \sin^{2} (ω) - \frac{5 e^{2} \sin^{2} (i)}{2} + \frac{e^{2}}{4} - \frac{3 e \sin^{3} (i) \sin^{3} (ω)}{2} - 2 e \sin^{3} (i) \sin (ω) + \frac{e \sin (i) \sin (ω)}{2} - \frac{\sin^{4} (i)}{4} + \frac{\sin^{2} (i)}{4}, H_{2} = - \frac{e^{3} \sin^{4} (i) \cos^{4} (ω)}{2} + e^{3} \sin^{2} (i) \cos^{2} (ω) - \frac{e^{3}}{2} + \frac{5 e^{2} \sin^{3} (i) \sin (ω) \cos^{2} (ω)}{2} - \frac{5 e^{2} \sin (i) \sin (ω)}{2} - \frac{3 e {(\cos 2 (i) - 1)}^{2} (\cos (4 w) - 1)}{64} + \frac{e \sin^{4} (i) \cos^{4} (ω)}{2} - \frac{5 e \sin^{2} (i) \sin^{2} (ω)}{2} - \frac{e \sin^{2} (i) \cos^{2} (ω)}{2} - \frac{\sin^{3} (i) \sin^{3} (ω)}{2}, and I_{2} = - \frac{e^{2} \sin^{4} (i) \cos^{4} (ω)}{4} + \frac{e^{2} \sin^{2} (i) \cos^{2} (ω)}{2} - \frac{e^{2}}{4} + \frac{e \sin^{3} (i) \sin (ω) \cos^{2} (ω)}{2} - \frac{e \sin (i) \sin (ω)}{2} - \frac{{(\cos (2 i) - 1)}^{2} (\cos (4 w) - 1)}{128} - \frac{\sin^{2} (i) \sin^{2} (ω)}{4} .

14. Appendix I: $Δ mag$ Intersection Polynomial

Here we include the coefficients of Eq. (50) from Sec. 2.2 as

Eq. (76)

A_{3} = \frac{e^{4} \sin^{4} (i)}{16}, B_{3} = \frac{1}{4} e^{3} (e \sin (ω) + \sin (i)) \sin^{3} (i), C_{3} = \frac{1}{8} e^{2} (e^{2} \cos^{2} (i) \cos^{2} (ω) - 3 e^{2} \cos^{2} (ω) + 3 e^{2} + 8 e \sin (i) \sin (ω) - 3 \cos^{2} (i) + 3) \sin^{2} (i), D_{3} = \frac{1}{4} e [- e^{3} \sin^{3} (ω) \cos^{2} (i) + e^{3} \sin^{3} (ω) + e^{3} \sin (ω) \cos^{2} (i) - 2 e^{2} \sin^{3} (i) \cos^{2} (ω) - 4 e^{2} \sin (i) \cos^{2} (ω) + 6 e^{2} \sin (i) - 6 e \sin (ω) \cos^{2} (i) + 6 e \sin (ω) + \sin^{3} (i)) \sin (i), E_{3} = \frac{e^{4} \sin^{4} (i) \sin^{4} (ω)}{16} - \frac{e^{4} \sin^{4} (i) \sin^{2} (ω)}{8} + \frac{e^{4} \sin^{4} (i)}{16} + \frac{e^{4} \sin^{2} (i) \sin^{2} (ω)}{8} - \frac{e^{4} \sin^{2} (i)}{8} + \frac{e^{4}}{16} + e^{3} \sin^{3} (i) \sin^{3} (ω) - e^{3} \sin^{3} (i) \sin (ω) + e^{3} \sin (i) \sin (ω) - e^{2} ξ \sin^{2} (i) \sin^{2} (ω) + \frac{e^{2} ξ \sin^{2} (i)}{2} + \frac{3 e^{2} \sin^{4} (i) \sin^{2} (ω)}{4} - \frac{3 e^{2} \sin^{4} (i)}{4} + \frac{3 e^{2} \sin^{2} (i) \sin^{2} (ω)}{2} + \frac{3 e^{2} \sin^{2} (i)}{4} + e \sin^{3} (i) \sin (ω) + \frac{\sin^{4} (i)}{16}, F_{3} = \frac{e^{3} \sin^{4} (i) \sin^{4} (ω)}{4} - \frac{e^{3} \sin^{4} (i) \sin^{2} (ω)}{2} + \frac{e^{3} \sin^{4} (i)}{4} + \frac{e^{3} \sin^{2} (i) \sin^{2} (ω)}{2} - \frac{e^{3} \sin^{2} (i)}{2} + \frac{e^{3}}{4} - e^{2} ξ \sin (i) \sin (ω) + \frac{3 e^{2} \sin^{3} (i) \sin^{3} (ω)}{2} - \frac{3 e^{2} \sin^{3} (i) \sin (ω)}{2} + \frac{3 e^{2} \sin (i) \sin (ω)}{2} - 2 e ξ \sin^{2} (i) \sin^{2} (ω) + e ξ \sin^{2} (i) + \frac{e \sin^{4} (i) \sin^{2} (ω)}{2} - \frac{e \sin^{4} (i)}{2} + e \sin^{2} (i) \sin^{2} (ω) + \frac{e \sin^{2} (i)}{2} + \frac{\sin^{3} (i) \sin (ω)}{4}, G_{3} = - \frac{e^{2} ξ \sin^{2} (i) \cos^{2} (ω)}{2} - \frac{e^{2} ξ}{2} + \frac{3 e^{2} \sin^{4} (i) \cos^{4} (ω)}{8} - \frac{3 e^{2} \sin^{2} (i) \cos^{2} (ω)}{4} + \frac{3 e^{2}}{8} - 2 e ξ \sin (i) \sin (ω) - e \sin^{3} (i) \sin (ω) \cos^{2} (ω) + e \sin (i) \sin (ω) - \frac{ξ \sin^{2} (i) \sin^{2} (ω)}{2} + \frac{ξ \sin^{2} (i) \cos^{2} (ω)}{2} + \frac{{(\cos (2 i) - 1)}^{2} (\cos (4 w) - 1)}{256} - \frac{\sin^{4} (i) \cos^{4} (ω)}{8} + \frac{3 \sin^{2} (i) \sin^{2} (ω)}{8} + \frac{\sin^{2} (i) \cos^{2} (ω)}{8}, H_{3} = - e ξ \sin^{2} (i) \cos^{2} (ω) - e ξ + \frac{e \sin^{4} (i) \cos^{4} (ω)}{4} - \frac{e \sin^{2} (i) \cos^{2} (ω)}{2} + \frac{e}{4} - ξ \sin (i) \sin (ω) - \frac{\sin^{3} (i) \sin (ω) \cos^{2} (ω)}{4} + \frac{\sin (i) \sin (ω)}{4}, and I_{3} = ξ^{2} - \frac{ξ \sin^{2} (i) \cos^{2} (ω)}{2} - \frac{ξ}{2} + \frac{\sin^{4} (i) \cos^{4} (ω)}{16} - \frac{\sin^{2} (i) \cos^{2} (ω)}{8} + \frac{1}{16} .

15. Appendix J: Proof that Ellipses Project to Ellipses

Let us suppose that we have some ellipse on an arbitrary plane indicated by the black ellipse in Fig. 4. We will say $\overline{A B}$ and $\overline{C D}$ are the principal axes of some this 3D ellipse. $\overline{A B}$ and $\overline{C D}$ intersect at point $O$ , the geometric center of this 3D ellipse. We will say point $P$ is any point on the ellipse with $H$ being the projection of point $P$ onto axis $\overline{A B}$ and $K$ being the projection of point $P$ onto axis $\overline{C D}$ . By the definition of an ellipse, we have

Eq. (77)

\frac{{\bar{O H}}^{2}}{{\bar{O B}}^{2}} + \frac{{\bar{O K}}^{2}}{{\bar{O D}}^{2}} = 1 .

Now we let $A^{'}$ , $B^{'}$ , $C^{'}$ , $D^{'}$ , $O^{'}$ , $K^{'}$ , $H^{'}$ , and $P^{'}$ be the perpendicular projections of points $A$ , $B$ , $C$ , $D$ , $O$ , $K$ , $H$ , and $P$ onto any given plane. Since perpendicular projections preserve the ratios of segments on a line, we can now say

Eq. (78)

\frac{{\bar{O^{'} H^{'}}}^{2}}{{\bar{O^{'} B^{'}}}^{2}} + \frac{{\bar{O^{'} K^{'}}}^{2}}{{\bar{O^{'} D^{'}}}^{2}} = 1 .

This equation means that point

P^{'}

belongs to the ellipse having

A^{'} B^{'}

and

C^{'} D^{'}

as conjugate diameters. For an ellipse, two diameters are conjugate if and only if the tangent line to the ellipse at an endpoint of one diameter is parallel to the other diameter.

16. Appendix K: Quartic Solution

The quartic function has a known analytical solution.¹⁹ We start with a set of useful simplifying terms

Eq. (79)

p_{0} = {(- \frac{3 A^{2}}{8} + B)}^{3},

Eq. (80)

p_{1} = {(A (\frac{A^{2}}{8} - \frac{B}{2}) + C)}^{2},

Eq. (81)

p_{2} = - A (A (\frac{3 A^{2}}{256} - \frac{B}{16}) + \frac{C}{4}) + D,

Eq. (82)

p_{3} = - \frac{3 A^{2}}{8} + B,

Eq. (83)

p_{4} = 2 A (\frac{A^{2}}{8} - \frac{B}{2}),

Eq. (84)

p_{5} = - \frac{p_{0}}{108} - \frac{p_{1}}{8} + \frac{p_{2} p_{3}}{3},

Eq. (85)

p_{6} = \sqrt[3]{\frac{p_{0}}{216} + \frac{p_{1}}{16} - \frac{p_{2} p_{3}}{6} + \sqrt{\frac{p_{5}^{2}}{4} + \frac{{(- p_{2} - \frac{p_{3}^{2}}{12})}^{3}}{27}}},

Eq. (86)

p_{7} = \frac{A^{2}}{4} - \frac{2 B}{3},

Eq. (87)

p_{8} = \frac{2 p_{2} + \frac{p_{3}^{2}}{6}}{3 p_{6}},

Eq. (88)

p_{9} = \sqrt{- 2 \sqrt[3]{p_{5}} + p_{7}},

Eq. (89)

p_{10} = \sqrt{2 p_{6} + p_{7} + p_{8}},

and

Eq. (90)

p_{11} = \frac{A^{2}}{2} - \frac{4 B}{3} .

The solutions are four piecewise functions:

Eq. (91)

x_{0} = {\begin{cases} - \frac{A}{4} + \frac{p_{9}}{2} - \frac{\sqrt{p_{11} + 2 \sqrt[3]{p_{5}} + \frac{- 2 C - p_{4}}{p_{9}}}}{2} & for p_{2} + \frac{p_{3}^{2}}{12} = 0 \\ - \frac{A}{4} - \frac{p_{10}}{2} - \frac{\sqrt{p_{11} - 2 p_{6} + p_{9} + \frac{2 C + p_{4}}{p_{10}}}}{2} & otherwise \end{cases},

Eq. (92)

x_{1} = {\begin{cases} - \frac{A}{4} + \frac{p_{9}}{2} + \frac{\sqrt{p_{11} + 2 \sqrt[3]{p_{5}} - \frac{2 C + p_{4}}{p_{9}}}}{2} & for p_{2} + \frac{p_{3}^{2}}{12} = 0 \\ - \frac{A}{4} - \frac{p_{10}}{2} + \frac{\sqrt{p_{11} - 2 p_{6} + p_{9} + \frac{2 C + p_{4}}{p_{10}}}}{2} & otherwise \end{cases},

Eq. (93)

x_{2} = {\begin{cases} - \frac{A}{4} - \frac{p_{9}}{2} - \frac{\sqrt{p_{11} + 2 \sqrt[3]{p_{5}} - \frac{- 2 C - p_{4}}{p_{9}}}}{2} & for p_{2} + \frac{p_{3}^{2}}{12} = 0 \\ - \frac{A}{4} + \frac{p_{10}}{2} - \frac{\sqrt{p_{11} - 2 p_{6} + p_{9} + \frac{- 2 C - p_{4}}{p_{10}}}}{2} & otherwise \end{cases},

and

Eq. (94)

x_{3} = {\begin{cases} - \frac{A}{4} - \frac{p_{9}}{2} + \frac{\sqrt{p_{11} + 2 \sqrt[3]{p_{5}} - \frac{- 2 C - p_{4}}{p_{9}}}}{2} & for p_{2} + \frac{p_{3}^{2}}{12} = 0 \\ - \frac{A}{4} + \frac{p_{10}}{2} + \frac{\sqrt{p_{11} - 2 p_{6} + p_{9} + \frac{- 2 C - p_{4}}{p_{10}}}}{2} & otherwise \end{cases} .

Theoretically, we always have four solutions. Theoretically, we can use equations of coefficients to determine how many roots are real, complex, and double roots. The general form of the quartic defines the expressions $Δ$ , $P$ , $D_{2}$ , $R$ , and $Δ_{0}$ , which are

Eq. (95)

Δ = 256 D^{3} - 192 A C D^{2} - 128 B^{2} D^{2} + 144 B C^{2} D - 27 C^{4} + 144 A^{2} B D^{2} - 6 A^{2} C^{2} D - 80 A B^{2} C D + 18 A B C^{3} + 16 B^{4} D - 4 B^{3} C^{2} - 27 A^{4} D^{2} + 18 A^{3} B C D - 4 A^{3} C^{3} - 4 A^{2} B^{3} D + A^{2} B^{2} C^{2},

Eq. (96)

P = 8 B - 3 A^{2},

Eq. (97)

D_{2} = 64 D - 16 B^{2} + 16 A^{2} B - 16 A C - 3 A^{4},

Eq. (98)

R = A^{3} + 8 C - 4 A B,

Eq. (99)

Δ_{0} = B^{2} - 3 A C + 12 D .

After evaluation of these constants, we can determine

Eq. (100)

solution types = {\begin{cases} 2 distinct R roots & 2 complex conjugate roots & Δ < 0 \\ 4 R distinct roots & Δ > 0 & P < 0 & D < 0 \\ 2 pair complex conjugate roots & Δ > 0 & (P > 0 | D > 0) \end{cases} .

However, the use of this theoretical classification does not work when numerical errors are introduced. The accumulation of numerical errors causes the solutions to be improperly classified for a sufficiently large number of cases to prohibit the use of these classifications. The

Δ

expression rarely evaluates to

Δ = 0

when it should do so quite frequently. Numerical rounding errors can also frequently result in the evaluation

Δ < 0

when in reality

Δ > 0

and vice versa.

17. Appendix L: Earth-Like Subpopulation

We use a definition of Earth-like exoplanet similar to that used in Ref. 8. We use a more conservative planetary radius range of $0.9 R_{\oplus} \leq R \leq 1.4 R_{\oplus}$ . We use the flux-at-planet range of $0.3586 \leq L_{plan} \leq 1.1080$ . To classify a planet, we use the time-averaged incident flux on the planet calculated by

Eq. (101)

L_{plan} = \frac{L_{*}}{{(a + \frac{a e^{2}}{2})}^{2}} .

18. Appendix M: 4-m Monolith Design Reference Mission

We calculated the Brown completeness and integration time adjusted completeness for a 4-m monolith telescope with a contrast of $10^{- 10}$ on the SAG13 planet population discussed in Ref. 1. The $Δ mag$ and separations of integration are scaled by the luminosity of each candidate target star. The orbital periods are similarly scaled by the mass of each target star in the case of integration time adjusted completeness. We optimized the mission schedule using the sequential least squares quadratic programming method discussed in Ref. 1 with the Brown completeness method. Our approach included the adoption of the use of the local zodiacal light minimum and filtering of target stars.

We assume a total telescope time of 213.67 days ( $\sim 7$ months) and that each observation requires 0.1 day of overhead time and 0.042 day settling time. The coronagraph parameters used are from a vortex charge six coronagraph at a wavelength of 500 nm with an IWA of 0.045 arc sec and OWA of 2.127 arc sec and a working angle of integration consistent with that used in Ref. 1.

The resulting Brown completeness optimized integration times for a DRM are included in Table 2. The summed Brown completeness for this DRM is 68.89, whereas the same summed integration time adjusted completeness is 63.11, a difference of 5.78. We can multiply either of these numbers by the SAG13 exoplanet occurrence rate of 5.62 to get the expected number of exoplanets detected in a single-visit blind search. These are 387.16 and 354.67 on average, respectively. This means that 32.49 planets are lost in a mission by not taking into account integration times when optimizing the mission schedule.

Table 2

4-m monolith DRM.

Name	Dist (pc)	Int. time (day)	Δmag	Cbrown	CIAC
HIP 439	4.34	0.399	19.557	0.049	0.041
HIP 522	25.71	0.181	23.476	0.038	0.033
HIP 544	13.67	0.6	23.403	0.16	0.141
HIP 746	16.78	0.066	27.142	0.247	0.24
HIP 910	18.75	0.383	24.916	0.16	0.146
HIP 950	21.28	0.341	24.504	0.106	0.095
HIP 1292	17.5	0.338	22.181	0.054	0.047
HIP 1392 B	15.19	0.293	21.362	0.039	0.033
HIP 1475	3.57	0.625	20.237	0.105	0.089
HIP 1499	23.25	0.139	21.736	0.018	0.016
HIP 1599	8.59	0.374	25.793	0.47	0.438
HIP 1803	20.86	0.241	22.209	0.036	0.031
HIP 2021	7.46	0.179	26.942	0.55	0.525
HIP 2072	23.81	0.129	25.773	0.103	0.097
HIP 2081	25.97	0.053	26.432	0.082	0.078
HIP 2711	25.48	0.193	23.752	0.045	0.039
HIP 3093	11.06	0.745	23.731	0.242	0.213
HIP 3170	24.96	0.192	23.257	0.038	0.033
HIP 3419	29.53	0.035	26.545	0.048	0.046
HIP 3497	22.06	0.184	21.928	0.025	0.022
HIP 3505	26.75	0.169	24.171	0.044	0.039
HIP 3583	15.16	0.567	23.919	0.17	0.151
HIP 3765	7.45	0.861	23.974	0.38	0.339
HIP 3810	23.45	0.27	24.446	0.078	0.07
HIP 3909	15.75	0.515	24.685	0.207	0.187
HIP 4148	14.17	0.293	21.192	0.039	0.033
HIP 4151	18.74	0.353	25.07	0.166	0.152
HIP 5336	7.55	0.66	24.859	0.459	0.418
HIP 5862	15.11	0.46	25.001	0.238	0.217
HIP 6379	16.81	0.217	21.188	0.028	0.024
HIP 6706	25.21	0.207	23.99	0.051	0.046
HIP 6813	28.62	0.13	24.535	0.037	0.033
HIP 7235	19.05	0.164	21.17	0.021	0.018
HIP 7339	20.66	0.236	22.106	0.035	0.03
HIP 7513	13.49	0.296	25.765	0.309	0.288
HIP 7601	27.38	0.129	23.04	0.023	0.02
HIP 7734	21.37	0.17	21.643	0.022	0.019
HIP 7918	12.74	0.521	24.978	0.295	0.269
HIP 7978	17.43	0.483	24.302	0.153	0.137
HIP 7981	7.53	0.733	24.611	0.44	0.398
HIP 8102	3.65	0.323	26.364	0.745	0.704
HIP 8362	10.07	0.753	24.254	0.318	0.285
HIP 8433	27.86	0.119	23.024	0.021	0.019
HIP 8497	23.19	0.241	25.078	0.098	0.09
HIP 8796	19.42	0.141	26.116	0.18	0.17
HIP 8903	17.99	0.08	26.814	0.217	0.209
HIP 9007	17.85	0.199	25.89	0.207	0.194
HIP 9236	22.01	0.069	26.565	0.138	0.132
HIP 9884	20.18	0.069	26.753	0.172	0.165
HIP 10138	10.78	0.746	23.544	0.231	0.204
HIP 10306	28.87	0.115	23.779	0.026	0.023
HIP 10644	10.78	0.528	25.092	0.359	0.329
HIP 10723	24.35	0.218	23.513	0.047	0.041
HIP 10798	12.67	0.605	22.975	0.146	0.128
HIP 11029	29.66	0.097	23.551	0.02	0.018
HIP 11783	26.7	0.179	24.768	0.055	0.049
HIP 12114	7.18	0.88	23.866	0.379	0.337
HIP 12186	25.78	0.173	23.307	0.034	0.03
HIP 12444	21.76	0.299	23.526	0.067	0.059
HIP 12623	24.19	0.249	24.661	0.076	0.068
HIP 12653	17.17	0.488	24.473	0.167	0.15
HIP 12777	11.13	0.324	25.829	0.381	0.356
HIP 12843	14.22	0.352	25.471	0.28	0.259
HIP 13402	10.35	0.771	23.572	0.245	0.216
HIP 13665	26.71	0.159	23.616	0.035	0.031
HIP 14146	27.17	0.112	25.539	0.062	0.058
HIP 14150	20.64	0.162	21.435	0.021	0.018
HIP 14632	10.54	0.328	25.855	0.402	0.375
HIP 14954	22.58	0.299	24.511	0.09	0.08
HIP 15330	12.01	0.669	24.437	0.278	0.25
HIP 15371	12.03	0.592	24.779	0.303	0.275
HIP 15457	9.14	0.567	25.073	0.413	0.378
HIP 15510	6.04	0.426	25.751	0.581	0.541
HIP 16245	21.68	0.264	25.16	0.121	0.111
HIP 16537	3.21	0.377	26.121	0.768	0.722
HIP 16852	13.96	0.339	25.542	0.289	0.268
HIP 17378	9.04	0.286	26.162	0.466	0.438
HIP 17420	13.95	0.326	21.324	0.044	0.038
HIP 17439	16.03	0.286	21.475	0.039	0.033
HIP 17651	17.63	0.25	25.64	0.206	0.192
HIP 18859	18.83	0.429	24.311	0.132	0.118
HIP 19076	16.94	0.466	23.433	0.114	0.101
HIP 19205	27.6	0.147	24.021	0.036	0.032
HIP 19335	21.0	0.341	23.973	0.09	0.08
HIP 19849	4.98	0.504	25.516	0.625	0.579
HIP 19859	21.33	0.244	22.41	0.039	0.034
HIP 19884	13.04	0.19	20.397	0.024	0.02
HIP 19893	20.46	0.206	25.608	0.152	0.141
HIP 19921	18.24	0.322	25.261	0.183	0.168
HIP 19990	28.94	0.122	24.383	0.033	0.03
HIP 21421	20.43	0.039	27.351	0.17	0.166
HIP 21547	29.43	0.107	23.929	0.025	0.022
HIP 21770	20.17	0.24	25.449	0.153	0.141
HIP 21861	28.67	0.127	24.318	0.034	0.03
HIP 22263	13.28	0.635	24.402	0.243	0.218
HIP 22336	26.42	0.154	23.183	0.029	0.026
HIP 22449	8.07	0.218	26.617	0.516	0.491
HIP 23311	8.71	0.832	23.26	0.258	0.226
HIP 23482	26.08	0.185	24.074	0.047	0.042
HIP 23693	11.65	0.434	25.351	0.346	0.319
HIP 23783	26.29	0.189	24.587	0.055	0.049
HIP 23875	27.4	0.044	26.48	0.067	0.064
HIP 23941	25.46	0.208	24.364	0.057	0.051
HIP 24186	3.91	0.377	19.298	0.044	0.037
HIP 24786	25.03	0.183	23.045	0.034	0.029
HIP 24813	12.63	0.458	25.207	0.311	0.285
HIP 25110	20.89	0.343	24.749	0.12	0.108
HIP 25278	14.39	0.517	24.847	0.246	0.224
HIP 25544	19.2	0.183	21.347	0.024	0.02
HIP 25623	13.02	0.163	20.148	0.02	0.017
HIP 25878	5.68	0.483	20.34	0.068	0.058
HIP 26394	18.32	0.446	24.074	0.128	0.114
HIP 26779	12.28	0.645	23.215	0.171	0.15
HIP 27072	8.93	0.256	26.326	0.474	0.448
HIP 27288	21.61	0.108	26.165	0.141	0.134
HIP 27321	19.44	0.155	26.012	0.178	0.167
HIP 27435	15.18	0.544	23.539	0.145	0.128
HIP 27628	26.73	0.081	25.936	0.071	0.066
HIP 27887	13.0	0.376	21.413	0.053	0.045
HIP 27890	26.25	0.19	24.754	0.058	0.053
HIP 28103	14.88	0.2	26.143	0.282	0.266
HIP 28360	24.87	0.031	27.217	0.098	0.096
HIP 28908	25.62	0.143	22.489	0.022	0.019
HIP 28954	15.27	0.327	21.575	0.046	0.039
HIP 29271	10.2	0.582	24.95	0.369	0.337
HIP 29295	5.75	0.401	20.008	0.051	0.043
HIP 29525	17.95	0.337	22.3	0.055	0.048
HIP 29568	16.72	0.423	22.752	0.083	0.073
HIP 29650	20.82	0.358	24.425	0.109	0.097
HIP 29800	19.25	0.394	24.751	0.145	0.131
HIP 29860	19.25	0.396	23.768	0.102	0.09
HIP 30503	21.88	0.228	22.401	0.036	0.031
HIP 31592	19.75	0.225	25.568	0.163	0.151
HIP 32349	2.63	0.025	30.177	0.872	0.861
HIP 32362	18.0	0.145	26.221	0.21	0.199
HIP 32366	24.38	0.207	23.238	0.041	0.036
HIP 32439	17.87	0.465	24.414	0.152	0.136
HIP 32480	16.72	0.498	24.565	0.18	0.163
HIP 32607	29.6	0.051	26.098	0.046	0.043
HIP 32765	25.26	0.215	24.454	0.061	0.055
HIP 32984	8.71	0.777	22.753	0.205	0.179
HIP 33094	25.89	0.154	22.824	0.026	0.023
HIP 33277	17.24	0.471	23.721	0.125	0.111
HIP 33690	18.33	0.245	21.638	0.033	0.028
HIP 33817	14.65	0.425	22.086	0.071	0.061
HIP 34017	19.13	0.372	23.259	0.082	0.072
HIP 34834	21.43	0.231	25.368	0.13	0.12
HIP 35136	16.89	0.503	24.178	0.156	0.139
HIP 35296 A	14.59	0.444	22.214	0.077	0.067
HIP 36366	18.05	0.245	25.626	0.197	0.183
HIP 36399	27.24	0.132	23.032	0.024	0.021
HIP 36439	20.24	0.38	24.354	0.114	0.102
HIP 36795	25.3	0.184	25.172	0.076	0.069
HIP 37279	3.51	0.063	28.816	0.796	0.778
HIP 37349	14.21	0.281	21.13	0.037	0.031
HIP 37606	24.67	0.233	24.476	0.067	0.06
HIP 37826	10.36	0.073	27.756	0.442	0.43
HIP 38784	17.19	0.358	22.262	0.059	0.051
HIP 38908	16.2	0.532	24.279	0.174	0.155
HIP 39157	16.77	0.216	21.173	0.028	0.024
HIP 39342	17.31	0.195	21.121	0.025	0.022
HIP 39757	19.48	0.083	26.621	0.184	0.177
HIP 39780	23.29	0.266	24.017	0.068	0.06
HIP 39903	19.98	0.305	25.155	0.148	0.135
HIP 40035	22.38	0.292	23.955	0.075	0.066
HIP 40693	12.49	0.674	23.689	0.205	0.181
HIP 40702	19.56	0.196	25.75	0.171	0.159
HIP 40706	28.63	0.124	25.177	0.045	0.041
HIP 40843	18.27	0.438	24.634	0.156	0.141
HIP 41211	26.62	0.159	23.493	0.033	0.029
HIP 41484	22.25	0.216	22.368	0.033	0.029
HIP 41926	12.21	0.621	22.933	0.151	0.132
HIP 42438	14.36	0.607	24.195	0.204	0.182
HIP 42808	11.14	0.631	22.574	0.141	0.123
HIP 43587	12.34	0.673	23.552	0.197	0.173
HIP 43726	17.39	0.442	23.35	0.105	0.092
HIP 43797	24.15	0.229	23.68	0.052	0.046
HIP 44075	21.03	0.323	23.609	0.077	0.068
HIP 44127	14.51	0.129	26.677	0.3	0.288
HIP 44143	26.42	0.182	24.346	0.049	0.044
HIP 44897	19.19	0.377	23.372	0.086	0.076
HIP 45170	20.36	0.232	21.979	0.033	0.028
HIP 45333	19.57	0.397	24.602	0.134	0.12
HIP 45343	6.11	0.566	20.808	0.09	0.076
HIP 45617	17.27	0.131	20.537	0.016	0.014
HIP 46509	17.33	0.342	25.261	0.201	0.184
HIP 46580	12.91	0.281	20.919	0.036	0.031
HIP 46733	23.82	0.112	25.913	0.105	0.099
HIP 46843	17.79	0.161	20.93	0.021	0.018
HIP 46853	13.48	0.158	26.543	0.327	0.312
HIP 47080	11.37	0.689	24.424	0.295	0.265
HIP 47592	15.01	0.457	25.025	0.241	0.22
HIP 48113	18.37	0.425	24.727	0.159	0.143
HIP 48331	11.16	0.237	20.399	0.03	0.025
HIP 48833	28.15	0.136	24.138	0.034	0.031
HIP 49081	15.05	0.571	24.431	0.206	0.185
HIP 49593	28.24	0.137	25.074	0.047	0.043
HIP 49669	24.31	0.029	27.69	0.106	0.105
HIP 49809	27.74	0.141	23.924	0.034	0.03
HIP 49908	4.87	0.989	22.892	0.356	0.312
HIP 50075	22.85	0.198	22.331	0.03	0.026
HIP 50384	22.81	0.257	23.348	0.054	0.047
HIP 50505	20.24	0.198	21.651	0.026	0.022
HIP 50564	21.37	0.309	24.921	0.119	0.108
HIP 50954	16.22	0.219	25.922	0.244	0.229
HIP 51459	12.78	0.464	25.177	0.305	0.28
HIP 51502	21.5	0.333	24.57	0.105	0.094
HIP 51523	21.81	0.298	24.914	0.112	0.102
HIP 51814	26.53	0.179	24.362	0.049	0.044
HIP 51933	25.08	0.195	23.413	0.04	0.035
HIP 53229	29.09	0.108	25.297	0.043	0.04
HIP 53253	29.13	0.101	25.391	0.044	0.041
HIP 53721	14.06	0.509	24.905	0.257	0.234
HIP 53910	24.45	0.039	26.951	0.103	0.1
HIP 54035	2.54	1.0	21.604	0.308	0.263
HIP 54182	28.99	0.125	24.696	0.037	0.033
HIP 54211	4.86	0.263	19.093	0.027	0.023
HIP 54745	21.93	0.213	22.22	0.032	0.027
HIP 54872	17.91	0.074	26.903	0.219	0.212
HIP 55691 A	13.71	0.174	20.377	0.022	0.018
HIP 55779	27.22	0.158	24.189	0.041	0.037
HIP 56445	27.23	0.133	23.075	0.024	0.021
HIP 56452	9.56	0.799	24.001	0.31	0.276
HIP 56802	26.73	0.158	23.596	0.034	0.03
HIP 56997	9.61	0.688	24.618	0.365	0.329
HIP 56998	12.4	0.178	20.179	0.022	0.018
HIP 57443	9.22	0.55	25.115	0.413	0.378
HIP 57507	17.47	0.358	22.357	0.061	0.053
HIP 57632	11.0	0.088	27.443	0.417	0.404
HIP 57757	10.93	0.259	26.148	0.398	0.374
HIP 57939	9.09	0.782	22.992	0.22	0.192
HIP 58001	25.5	0.037	26.884	0.089	0.087
HIP 58345	10.16	0.549	21.707	0.093	0.08
HIP 58576	12.76	0.676	24.144	0.236	0.21
HIP 58803	25.32	0.212	24.371	0.058	0.052
HIP 58910	5.52	0.78	21.563	0.167	0.143
HIP 59072	19.76	0.204	25.689	0.166	0.154
HIP 59199	14.94	0.253	25.855	0.273	0.256
HIP 59750	22.35	0.244	22.816	0.043	0.038
HIP 59774	24.69	0.072	26.26	0.097	0.092
HIP 60965	26.63	0.051	26.407	0.075	0.071
HIP 61053	21.78	0.261	22.834	0.047	0.041
HIP 61084	27.15	0.037	26.765	0.07	0.068
HIP 61174	18.28	0.271	25.475	0.188	0.174
HIP 61291	16.18	0.248	21.274	0.033	0.028
HIP 61317	8.44	0.386	25.755	0.475	0.442
HIP 62145	14.88	0.318	21.442	0.043	0.037
HIP 62207	17.38	0.456	23.567	0.116	0.102
HIP 62523	16.93	0.418	22.789	0.083	0.073
HIP 62951	7.53	0.378	20.357	0.049	0.041
HIP 62956	25.31	0.03	27.46	0.093	0.091
HIP 63076	29.29	0.112	24.085	0.028	0.025
HIP 63613	27.87	0.107	25.502	0.055	0.051
HIP 63721	4.62	0.538	20.255	0.081	0.068
HIP 64394	9.13	0.381	25.732	0.447	0.416
HIP 64408	20.67	0.337	24.857	0.127	0.115
HIP 64583	18.2	0.38	25.015	0.174	0.159
HIP 64792	17.56	0.467	24.586	0.166	0.15
HIP 64924	8.56	0.553	25.147	0.439	0.403
HIP 65109	18.02	0.076	26.859	0.216	0.209
HIP 65352	15.45	0.25	21.164	0.033	0.028
HIP 65355	16.14	0.162	20.661	0.021	0.017
HIP 65530	21.17	0.228	22.163	0.034	0.029
HIP 65721	17.99	0.426	24.802	0.169	0.153
HIP 66249	22.71	0.095	26.18	0.124	0.117
HIP 66765	15.65	0.293	21.439	0.039	0.034
HIP 67153	19.4	0.236	25.54	0.169	0.157
HIP 67155	5.41	0.319	19.535	0.036	0.03
HIP 67275	15.62	0.345	25.386	0.242	0.224
HIP 67422 A	13.41	0.362	21.413	0.05	0.043
HIP 67927	11.4	0.139	26.884	0.397	0.38
HIP 68030	24.86	0.165	22.587	0.026	0.023
HIP 68184	10.06	0.729	22.959	0.194	0.169
HIP 68682	16.98	0.42	22.833	0.085	0.074
HIP 68933	18.03	0.072	26.912	0.216	0.209
HIP 69671	21.22	0.258	22.541	0.043	0.037
HIP 69701	22.24	0.194	25.491	0.121	0.112
HIP 69713	29.07	0.124	24.773	0.037	0.034
HIP 69965	18.03	0.43	23.529	0.106	0.093
HIP 69972	11.8	0.592	22.482	0.124	0.108
HIP 69989	26.1	0.181	23.903	0.043	0.038
HIP 70319	17.19	0.411	22.815	0.082	0.072
HIP 70497	14.53	0.261	25.85	0.284	0.266
HIP 70857	19.89	0.178	21.426	0.023	0.02
HIP 70873	23.74	0.172	22.267	0.025	0.022
HIP 71075	26.61	0.056	26.311	0.074	0.071
HIP 71181	13.22	0.337	21.254	0.046	0.039
HIP 71284	15.83	0.314	25.495	0.241	0.224
HIP 71855	20.0	0.203	21.649	0.027	0.023
HIP 71957	18.27	0.204	25.823	0.197	0.185
HIP 72567	18.17	0.429	23.644	0.109	0.097
HIP 72603	22.98	0.285	24.398	0.082	0.073
HIP 72622	23.24	0.062	26.537	0.119	0.114
HIP 72848	11.51	0.715	23.602	0.219	0.193
HIP 73100	25.11	0.202	23.679	0.046	0.04
HIP 73165	26.9	0.165	25.043	0.058	0.052
HIP 73996	19.55	0.347	25.001	0.15	0.136
HIP 74273	24.2	0.167	22.366	0.025	0.022
HIP 74537	17.67	0.311	21.982	0.046	0.04
HIP 74605	25.34	0.211	24.386	0.059	0.053
HIP 74702	15.85	0.274	21.367	0.036	0.031
HIP 74975	25.38	0.211	24.352	0.058	0.052
HIP 75181	14.81	0.589	24.092	0.188	0.167
HIP 75718	20.58	0.139	21.214	0.018	0.015
HIP 75809	21.85	0.176	21.803	0.024	0.02
HIP 76219	28.93	0.123	24.381	0.033	0.03
HIP 76267	23.01	0.038	27.087	0.125	0.121
HIP 76829	17.44	0.34	25.26	0.198	0.182
HIP 77052	14.66	0.578	23.722	0.165	0.146
HIP 77070	22.68	0.084	26.306	0.125	0.119
HIP 77257	12.12	0.385	25.5	0.338	0.313
HIP 77358	15.35	0.525	23.372	0.132	0.117
HIP 77760	15.89	0.352	25.336	0.234	0.216
HIP 77801	17.35	0.437	23.278	0.102	0.089
HIP 77952	12.38	0.126	26.908	0.365	0.351
HIP 78072	11.25	0.272	26.058	0.385	0.361
HIP 78459	17.24	0.488	24.429	0.164	0.147
HIP 78527	21.03	0.173	25.748	0.145	0.135
HIP 78775	14.52	0.451	22.244	0.079	0.069
HIP 79190	14.67	0.304	21.331	0.041	0.035
HIP 79248	17.57	0.328	22.113	0.051	0.044
HIP 79537	13.89	0.226	20.787	0.029	0.025
HIP 79578	21.55	0.168	21.662	0.022	0.019
HIP 79672	13.9	0.625	24.242	0.217	0.194
HIP 79822	29.74	0.107	24.346	0.028	0.026
HIP 80179	27.27	0.165	24.698	0.049	0.044
HIP 80331	28.23	0.053	26.202	0.058	0.055
HIP 80337	12.78	0.641	24.472	0.261	0.235
HIP 80686	12.12	0.492	25.116	0.32	0.293
HIP 81300	9.75	0.794	23.977	0.302	0.268
HIP 81935	14.26	0.166	20.386	0.021	0.017
HIP 82003	9.81	0.179	19.732	0.021	0.017
HIP 82020	26.11	0.192	24.816	0.061	0.055
HIP 82396	19.54	0.085	26.601	0.183	0.175
HIP 82587	29.19	0.111	23.885	0.026	0.023
HIP 82588	17.25	0.303	21.827	0.043	0.037
HIP 82621	26.94	0.13	22.785	0.021	0.019
HIP 82860	15.26	0.419	25.138	0.241	0.22
HIP 83000	28.04	0.084	25.743	0.056	0.052
HIP 83389	18.59	0.25	21.73	0.034	0.029
HIP 83541	17.84	0.284	21.815	0.04	0.034
HIP 83591	10.71	0.218	20.197	0.027	0.023
HIP 83601	20.67	0.317	23.214	0.067	0.058
HIP 83609	5.62	0.311	19.561	0.035	0.029
HIP 83990	13.62	0.26	20.929	0.034	0.029
HIP 84143	22.53	0.104	26.117	0.126	0.119
HIP 84478	5.95	0.95	23.071	0.332	0.29
HIP 84862	14.33	0.578	24.57	0.231	0.208
HIP 85042	19.52	0.316	22.628	0.056	0.049
HIP 85235	12.8	0.588	22.853	0.136	0.119
HIP 85295	7.7	0.45	20.704	0.063	0.053
HIP 86201	23.16	0.251	25.007	0.097	0.088
HIP 86400	11.0	0.644	22.624	0.147	0.128
HIP 86486	21.45	0.29	25.037	0.121	0.11
HIP 86614	22.84	0.226	25.235	0.107	0.098
HIP 86620	23.06	0.257	23.538	0.057	0.05
HIP 86742	25.09	0.079	26.132	0.091	0.086
HIP 86796	15.51	0.512	24.744	0.215	0.195
HIP 86974	8.31	0.253	26.393	0.501	0.474
HIP 88175	23.55	0.233	25.073	0.094	0.085
HIP 88635	29.7	0.065	25.833	0.044	0.041
HIP 88694	17.55	0.439	23.382	0.104	0.092
HIP 88771	26.63	0.088	25.85	0.071	0.066
HIP 88972	11.02	0.686	23.01	0.178	0.156
HIP 89042	17.61	0.476	24.3	0.15	0.135
HIP 89348	22.92	0.28	24.79	0.094	0.085
HIP 89474	22.82	0.213	22.54	0.034	0.03
HIP 89962	18.54	0.154	26.107	0.197	0.186
HIP 90496	23.97	0.088	26.129	0.105	0.099
HIP 90790	13.25	0.399	21.587	0.059	0.05
HIP 91262	7.68	0.03	29.135	0.559	0.549
HIP 91768	3.52	0.41	19.289	0.05	0.041
HIP 92024	28.55	0.136	24.804	0.041	0.037
HIP 92043	19.21	0.229	25.603	0.174	0.162
HIP 92161	28.89	0.115	25.259	0.044	0.041
HIP 92549	25.67	0.187	23.697	0.042	0.037
HIP 93747	25.46	0.055	26.449	0.088	0.084
HIP 93858	16.95	0.437	23.018	0.094	0.082
HIP 93966	21.43	0.285	23.049	0.056	0.049
HIP 94083	27.3	0.159	24.403	0.044	0.039
HIP 94376	29.87	0.06	25.889	0.043	0.04
HIP 95149	18.83	0.293	22.147	0.045	0.039
HIP 95319	15.76	0.465	22.822	0.097	0.085
HIP 95447	15.18	0.537	24.662	0.217	0.196
HIP 95501	15.53	0.158	26.357	0.269	0.256
HIP 95995 A	16.96	0.367	22.278	0.061	0.053
HIP 96100	5.75	0.532	25.374	0.576	0.531
HIP 96258	25.11	0.198	23.568	0.043	0.038
HIP 96441	18.34	0.281	25.424	0.186	0.172
HIP 96895	21.08	0.315	23.431	0.07	0.062
HIP 96901	21.21	0.284	22.947	0.054	0.047
HIP 97649	5.12	0.061	28.669	0.689	0.673
HIP 97650	27.87	0.135	23.712	0.03	0.027
HIP 97675	19.19	0.403	24.679	0.143	0.129
HIP 97944	14.05	0.552	23.011	0.128	0.113
HIP 98036	13.7	0.257	25.942	0.308	0.289
HIP 98066	25.99	0.196	24.561	0.057	0.051
HIP 98470	21.25	0.322	23.721	0.078	0.069
HIP 98698	12.86	0.225	20.606	0.029	0.024
HIP 98767	15.86	0.539	23.947	0.159	0.142
HIP 98792	15.77	0.222	21.051	0.029	0.025
HIP 98819	17.77	0.453	23.791	0.122	0.108
HIP 98921	18.79	0.371	23.036	0.076	0.067
HIP 98959	17.73	0.425	23.298	0.098	0.086
HIP 99031	23.95	0.242	23.966	0.061	0.054
HIP 99137	23.43	0.197	22.533	0.031	0.027
HIP 99240	6.11	0.299	26.312	0.598	0.564
HIP 99461	6.01	0.782	24.537	0.499	0.452
HIP 99572	27.7	0.118	22.892	0.02	0.018
HIP 99701	6.2	0.402	20.137	0.052	0.044
HIP 99825	8.91	0.825	23.923	0.323	0.287
HIP 100017	17.57	0.453	23.653	0.117	0.104
HIP 100511	26.25	0.131	22.509	0.02	0.018
HIP 100925	19.52	0.216	21.657	0.029	0.025
HIP 101345	24.4	0.212	23.393	0.044	0.039
HIP 101612	27.79	0.153	24.797	0.046	0.042
HIP 101983	24.66	0.232	24.425	0.066	0.059
HIP 101997	14.38	0.499	22.622	0.101	0.088
HIP 102040	20.95	0.242	22.254	0.037	0.032
HIP 102333	24.17	0.195	25.258	0.09	0.083
HIP 102422	14.27	0.201	26.191	0.299	0.282
HIP 102485	14.68	0.297	25.672	0.275	0.256
HIP 102488	22.29	0.068	26.555	0.133	0.128
HIP 103389	21.97	0.291	23.579	0.067	0.059
HIP 103458	22.14	0.176	21.863	0.024	0.02
HIP 103673	28.01	0.12	23.163	0.022	0.02
HIP 104214	3.5	0.778	24.704	0.657	0.599
HIP 104217	3.5	0.966	23.852	0.562	0.503
HIP 104239 A	17.57	0.13	20.584	0.016	0.014
HIP 105090	3.95	1.019	22.792	0.39	0.341
HIP 105199	15.04	0.077	27.164	0.291	0.282
HIP 105858	9.26	0.361	25.804	0.445	0.415
HIP 106440	4.95	0.288	19.245	0.031	0.026
HIP 106559	27.1	0.14	23.208	0.027	0.023
HIP 106696	14.62	0.283	21.212	0.037	0.032
HIP 107089	21.2	0.174	25.728	0.142	0.132
HIP 107350	17.89	0.43	23.466	0.104	0.092
HIP 107556	11.87	0.142	26.818	0.38	0.364
HIP 107649	15.99	0.539	24.202	0.173	0.154
HIP 107975	27.45	0.143	23.663	0.031	0.028
HIP 108036	26.62	0.173	24.239	0.046	0.041
HIP 108870	3.62	0.625	25.183	0.688	0.633
HIP 109176	11.73	0.252	26.12	0.371	0.349
HIP 109378	21.56	0.164	21.628	0.022	0.018
HIP 109422	18.28	0.399	24.905	0.168	0.153
HIP 109427	28.3	0.071	25.894	0.055	0.052
HIP 109821	22.05	0.24	22.624	0.04	0.035
HIP 110649	20.56	0.369	24.363	0.11	0.098
HIP 111449	22.68	0.294	24.314	0.083	0.074
HIP 112117	23.64	0.22	23.07	0.042	0.037
HIP 112447	16.3	0.286	25.586	0.233	0.217
HIP 112935	27.28	0.155	24.059	0.038	0.034
HIP 113137	26.75	0.131	22.721	0.021	0.019
HIP 113283	7.61	0.849	22.821	0.242	0.211
HIP 113357	15.61	0.554	24.306	0.186	0.167
HIP 113368	7.7	0.065	28.225	0.554	0.54
HIP 113421	19.86	0.31	22.688	0.055	0.048
HIP 113576	8.22	0.273	20.018	0.033	0.028
HIP 113860	29.42	0.109	24.09	0.027	0.024
HIP 114046	3.28	1.0	21.937	0.308	0.265
HIP 114430	27.65	0.123	23.011	0.022	0.019
HIP 114570	24.59	0.192	25.221	0.085	0.077
HIP 114622	6.54	0.811	24.377	0.46	0.415
HIP 114924	20.5	0.362	24.042	0.098	0.087
HIP 114948	20.54	0.356	23.94	0.094	0.083
HIP 114996	23.06	0.154	25.676	0.112	0.104
HIP 116085	16.93	0.296	21.711	0.041	0.036
HIP 116584	26.41	0.133	25.433	0.068	0.063
HIP 116727	14.1	0.186	26.3	0.305	0.29
HIP 116745	11.42	0.44	21.427	0.065	0.056
HIP 116771	13.71	0.317	25.66	0.3	0.279
HIP 120005	6.11	0.543	20.704	0.083	0.07

Acknowledgments

This work was funded by the Science Investigation Team of the Nancy Grace Roman Space Telescope under NASA (Grant No. NNX15AB40G). This work was an expansion of the work done in Ref. 20.

References

1.

D. R. Keithly et al., “Optimal scheduling of exoplanet direct imaging single-visit observations of a blind search survey,” J. Astron. Telesc. Instrum. Syst., 6 027001 (2020). https://doi.org/10.1117/1.JATIS.6.2.027001 Google Scholar

2.

R. A. Brown, “Single-visit photometric and obscurational completeness,” Astrophys. J., 624 1010 –1024 (2005). https://doi.org/10.1086/429124 ASJOAB 0004-637X Google Scholar

3.

D. Savransky, C. Delacroix and D. Garrett, “EXOSIMS: exoplanet open-source imaging mission simulator,” (2017). Google Scholar

4.

R. K. Kopparapu et al., “Exoplanet classification and yield estimates for direct imaging missions,” Astrophys. J., 856 122 (2018). https://doi.org/10.3847/1538-4357/aab205 ASJOAB 0004-637X Google Scholar

5.

D. Garrett, “Exoplanet direct imaging detection metrics and exoplanet popu- lations,” Cornell University, (2018). Google Scholar

6.

D. Garrett and D. Savransky, “Analytical formulation of the single-visit completeness joint probability density function,” Astrophys. J., 828 20 (2016). https://doi.org/10.3847/0004-637X/828/1/20 ASJOAB 0004-637X Google Scholar

7.

N. J. Kasdin et al., “The Nancy Grace Roman Space Telescope Coronagraph Instrument (CGI) technology demonstration,” Proc. SPIE, 11443 114431U (2020). https://doi.org/10.1117/12.2562997 Google Scholar

8.

Science and Technology Definition Team, “HabEx final report,” (2019). Google Scholar

9.

Science and Technology Definition Team, “LUVOIR final report,” (2019). Google Scholar

10.

R. Morgan et al., “Exoplanet yield estimation for decadal study concepts using EX,” in 227th Meeting Am. Astron. Soc., (2020). Google Scholar

11.

R. Morgan et al., “The standard definitions and evaluation team final report a common comparison of exoplanet yield,” (2019). Google Scholar

12.

Jet Propulsion Laboratory, “Habitable exoplanet observatory final report,” (2019). Google Scholar

13.

R. A. Brown and R. Soummer, “New completeness methods for estimating exoplanet discoveries by direct detection,” Astrophys. J., 715 122 –131 (2010). https://doi.org/10.1088/0004-637X/715/1/122 ASJOAB 0004-637X Google Scholar

14.

B. Spain, Analytical Conics, Dover, London (2007). Google Scholar

15.

W. H. Besant, Conic Sections, Treated Geometrically, 9th ed.London(1895). Google Scholar

16.

T. A. A. Broadbent, L. Kuipers and R. Timman, Handbook of Mathematics, 5th ed.Springer, Berlin (1970). Google Scholar

17.

E. Agol, “Rounding up the wanderers: optimizing coronagraphic searches for extrasolar planets,” Mon. Not. R. Astron. Soc., 374 (4), 1271 –1289 (2007). https://doi.org/10.1111/j.1365-2966.2006.11232.x MNRAA4 0035-8711 Google Scholar

18.

D. Savransky, E. Cady and N. J. Kasdin, “Parameter distributions of Keplerian orbiTS,” Astrophys. J., 728 (1), 66 (2011). https://doi.org/10.1088/0004-637X/728/1/66 ASJOAB 0004-637X Google Scholar

19.

M. Abramowitz and I. A. Stegun, Handbook of Mathemaical Functions with Formulas, Graphsm and Mathematical Tables, 9th ed.Dover(1972). Google Scholar

20.

D. Keithly and D. Savransky, “Integration time adjusted completeness,” Proc. SPIE, 11443 1144324 (2020). https://doi.org/10.1117/12.2562636 Google Scholar

Biography

Dean Keithly is a PhD candidate in the Sibley School of Mechanical and Aerospace Engineering at Cornell University. He received his bachelor’s degree in mechanical engineering from Michigan Technological University, where he built the Oculus-ASR satellite thermal control subsystem. He received his master’s degree in systems engineering from Cornell University, working on projects funded by the NASA Institute of Advanced Concepts. His primary research focus is on spacecraft modeling, mission modeling, optimization, and scheduling.

Dmitry Savransky is a professor in the Sibley School of Mechanical and Aerospace Engineering at Cornell University and PI of the Space Imaging and Optical Systems Laboratory. He received his PhD from Princeton University in 2011 and was a postdoctoral researcher at Lawrence Livermore National Laboratory, where he worked on integration, testing and commissioning of the Gemini Planet Imager. His research focuses on the optimization of ground and space-based exoplanet imaging surveys, control of autonomous optical systems, and advanced image processing.

Corey Spohn received a Bachelor of Science in both physics and engineering science and mechanics from Virginia Tech in May 2018. During his time at Virginia Tech he focused heavily on astrophysics while conducting research on anti-frosting surfaces, inter-droplet ice bridging, and flying snakes. For his senior design project, he investigated human-robot dynamics in a crowded panic situation by having humans play the children’s game sharks and minnows with robots in a motion capture studio. He began working towards a PhD in the fall of 2018 and joined the Space Imaging and Optical Systems Lab in December 2018 where he is currently researching exoplanet mission scheduling.

CC BY: © The Authors. Published by SPIE under a Creative Commons Attribution 4.0 Unported License. Distribution or reproduction of this work in whole or in part requires full attribution of the original publication, including its DOI.

Citation Download Citation

Dean Keithly, Dmitry Savransky, and Corey Spohn "Integration time adjusted completeness," Journal of Astronomical Telescopes, Instruments, and Systems 7(3), 037002 (1 September 2021). https://doi.org/10.1117/1.JATIS.7.3.037002

Received: 30 March 2021; Accepted: 16 August 2021; Published: 1 September 2021

Access the abstract

JOURNAL ARTICLE
54 PAGES

DOWNLOAD PAPER SAVE TO MY LIBRARY

GET CITATION

CITATIONS

Cited by 1 scholarly publication.

Explore citations on Lens.org

KEYWORDS

Planets

Stars

Exoplanets

Visibility

Telescopes

Monte Carlo methods

Picosecond phenomena

1.

Introduction

Fig. 1

Fig. 2

2.

Methods

Eq. (1)

Eq. (2)

Fig. 3

2.1.

Projected Orbit and Separation Intersection

2.1.1.

Projection of the elliptical orbit

Eq. (3)

Eq. (4)

Fig. 4

Eq. (5)

Eq. (6)

Eq. (7)

Eq. (8)

Eq. (9)

Eq. (10)

Eq. (11)

Eq. (12)

Eq. (13)

Eq. (14)

Eq. (15)

Fig. 5

Eq. (16)

2.1.2.

Global and local extrema of planet-star separation

Eq. (17)

Eq. (18)

Eq. (19)

Eq. (20)

Eq. (21)

Eq. (22)

Eq. (23)

Fig. 6

Eq. (24)

Eq. (25)

Eq. (26)

Eq. (27)

Eq. (28)

Eq. (29)

Eq. (30)

Eq. (31)

Eq. (32)

Eq. (33)

Eq. (34)

2.1.3.

Intersections between a circle and an ellipse

Eq. (35)

Eq. (36)

Fig. 7

Fig. 8

Eq. (37)

Eq. (38)

Eq. (39)

Table 1

Eq. (40)

Eq. (41)

2.2.

Δmag Intersections

Eq. (42)

Eq. (43)

Eq. (44)

Eq. (45)

Eq. (46)

Eq. (47)

Eq. (48)

Eq. (49)

Eq. (50)

2.3.

ν from X and Y

Eq. (51)

Eq. (52)

Eq. (53)

2.4.

Calculate t from ν

$Δ mag$ Intersections

$ν$ from $X$ and $Y$

Calculate $t$ from $ν$

$ν$ from $s$

$ν$ from Δmag

Appendix B: True Anomalies of $s$ Intersection Process

Appendix C: $ν$ from $Δ mag$ Extrema process

Appendix D: $ν$ from ${Δ mag}_{\lim}$ Intersection Process

Appendix E: Projected Ellipse Semimajor Axis and Semiminor Axis, and $θ$

Appendix F: $s$ Extrema Polynomial