# Measurability of the Heliocentric Momentum Enhancement from a Kinetic Impact: The Double Asteroid Redirection Test (DART) Mission

Rahil Makadia, Steven R. Chesley, Davide Farnocchia, Shantanu P. Naidu, Damya Souami, Paolo Tanga, Kleomenis Tsiganis, Masatoshi Hirabayashi, Siegfried Eggl

► To cite this version:

Rahil Makadia, Steven R. Chesley, Davide Farnocchia, Shantanu P. Naidu, Damya Souami, et al. Measurability of the Heliocentric Momentum Enhancement from a Kinetic Impact: The Double Asteroid Redirection Test (DART) Mission. The Planetary Science Journal, 2024, 5, 10.3847/PSJ/ad1bce. obspm-04753309

HAL Id: obspm-04753309

https://hal-obspm.ccsd.cnrs.fr/obspm-04753309v1

Submitted on 30 Oct 2024

**HAL** is a multi-disciplinary open access archive for the deposit and dissemination of scientific research documents, whether they are published or not. The documents may come from teaching and research institutions in France or abroad, or from public or private research centers.

L'archive ouverte pluridisciplinaire **HAL**, est destinée au dépôt et à la diffusion de documents scientifiques de niveau recherche, publiés ou non, émanant des établissements d'enseignement et de recherche français ou étrangers, des laboratoires publics ou privés.

<figure><img src="image_1.png" /></figure>

Distributed under a Creative Commons Attribution 4.0 International License
<figure><img src="image_2.png" /></figure>

# Measurability of the Heliocentric Momentum Enhancement from a Kinetic Impact: The Double Asteroid Redirection Test (DART) Mission

Rahil Makadia¹, Steven R. Chesley², Davide Farnocchia², Shantanu P. Naidu², Damya Souami³,⁴,⁵,⁶, Paolo Tanga⁴, Kleomenis Tsiganis⁷, Masatoshi Hirabayashi⁸, and Siegfried Eggl¹,⁹,¹⁰

¹Department of Aerospace Engineering, University of Illinois at Urbana-Champaign, Urbana, IL 61801, USA; makadia2@illinois.edu

²Jet Propulsion Laboratory, California Institute of Technology, Pasadena, CA 91109, USA

³LESIA, Observatoire de Paris—Section Meudon, F-92195 Meudon Cedex, France

⁴Université Côte d'Azur, Observatoire de la Côte d'Azur, CNRS, Laboratoire Lagrange, F-06304 Nice Cedex 4, France

⁵Department of Astronomy and Earth & Planetary Science, 501 Campbell Hall, University of California, Berkeley, Berkeley, CA 94720, USA

⁶naXys, Department of Mathematics, University of Namur, Rue de Bruxelles 61, B-5000 Namur, Belgium

⁷Department of Physics, Aristotle University of Thessaloniki, GR 54124 Thessaloniki, Greece

⁸Daniel Guggenheim School of Aerospace Engineering, Georgia Institute of Technology, Atlanta, GA 30332, USA

⁹Department of Astronomy, University of Illinois at Urbana-Champaign, Urbana, IL 61801, USA

¹⁰National Center for Supercomputing Applications, University of Illinois at Urbana-Champaign, Urbana, IL 61801, USA

Received 2023 September 15; revised 2023 December 12; accepted 2024 January 4; published 2024 February 14

## Abstract

The NASA Double Asteroid Redirection Test (DART) has demonstrated the capability of successfully conducting kinetic impact-based asteroid deflection missions. The changes in the Didymos–Dimorphos mutual orbit as a result of the DART impact have already been measured. To fully assess the heliocentric outcome of deflection missions, the heliocentric momentum enhancement parameter, $\beta_{\odot}$, needs to be determined and disentangled from other nongravitational phenomena such as the Yarkovsky effect. Here we explore the measurability of $\beta_{\odot}$ resulting from DART, which we estimate simultaneously with nongravitational accelerations using a least-squares filter. Results show that successful stellar occultation measurements of the Didymos system in the second half of 2024 in addition to the ones in the 2022–2023 campaigns can achieve a statistically significant estimate of $\beta_{\odot}$, with an uncertainty slightly above 20% for an assumed $\beta_{\odot} = 3$. Adding additional occultation measurements and pseudorange measurements from the Hera spacecraft operations at Didymos starting in 2027 decreases this relative uncertainty to under 6%. We find that pre-impact occultation observations combined with post-impact occultations would have yielded substantially higher signal-to-noise ratios on the heliocentric deflection. Additionally, pre-impact occultations would also have enabled a statistically significant $\beta_{\odot}$ estimate using only one additional occultation in 2023 September. Therefore, we conclude that future asteroid deflection missions would greatly benefit from both pre- and post-deflection occultation measurements to help assess the resulting orbital changes.

*Unified Astronomy Thesaurus concepts: Asteroid dynamics (2210); Near-Earth objects (1092); Impact phenomena (779)*

## 1. Introduction

NASA's Double Asteroid Redirection Test (DART; Rivkin et al. 2021) mission successfully changed the mutual orbit of Dimorphos, the secondary asteroid in the (65803) Didymos binary asteroid system, around Didymos on 2022 September 26. The change in Dimorphos's orbital period using observations from the first few weeks after impact was estimated to be $-33 \pm 0.3$ minutes (Thomas et al. 2023). DART was the world's first full-scale planetary defense technology demonstration mission, and it successfully showed that the trajectory of a sub-kilometer-sized asteroid can be altered through a kinetic impactor (KI). Deflecting near-Earth objects (NEOs) using a KI involves a hypervelocity collision between the spacecraft and the target (Ahrens & Harris 1994). The idea is that, given sufficient lead time, even a small deflection can push a potentially hazardous object out of Earth's path.

A characteristic feature of the KI method is that materials from the surface of the target asteroid, such as fine particles and pebbles, get kicked off in significantly complex ways. For

DART, even boulders were kicked off the surface of the asteroid as a result of the impact (Jewitt et al. 2023). A fraction of the particles in this ejecta cloud possess enough momentum to escape the gravitational influence of the target body. The escaping ejecta impart additional momentum onto the target. Therefore, the overall deflection of the body is a result of the combined momenta delivered by the spacecraft and this escaping ejecta. This entire exchange is characterized by the momentum enhancement parameter, $\beta$ (Holsapple & Housen 2012). More specifically, $\beta$ is the ratio of the momentum imparted to the target asteroid to the momentum of the impactor at the time of impact.

In the case of a binary system like Didymos, there are two $\beta$ values to be considered: a local Dimorphos $\beta$ that includes all ejecta that escape Dimorphos's gravitational pull (determined to be between 2.2 and 4.9 by Cheng et al. 2023), and a heliocentric system momentum enhancement parameter, $\beta_{\odot}$, caused by the ejecta that escape the Didymos binary system's gravitational influence. ESA's Hera mission will visit the Didymos system to study the consequences of the DART experiment in detail, which will help constrain the mass of Dimorphos, the shape of the crater caused by DART, and consequently the momentum enhancement parameters (Michel et al. 2022). For a more comprehensive discussion of how this

<figure><img src="image_3.png" /></figure>

Original content from this work may be used under the terms of the [Creative Commons Attribution 4.0 licence](https://creativecommons.org/licenses/by/4.0/). Any further distribution of this work must maintain attribution to the author(s) and the title of the work, journal citation and DOI.
information will aid in further understanding of the mutual orbit dynamics of the binary system, the reader is referred to Richardson et al. (2022, 2024).

Makadia et al. (2022) used statistical estimates of $\beta_{\odot}$ to determine that the DART impact would not put the Didymos system on a collision course with Earth. Given the importance of $\beta_{\odot}$ in both planning and deflection assessment of future KI missions, we set out to investigate whether the current observation data are sufficient to effectively estimate the heliocentric momentum enhancement factor from classical astrometry, derived astrometry (from stellar occultations), and pseudorange measurements of the target asteroid. If not, we explore observation scenarios that allow a statistically significant estimate of $\beta_{\odot}$.

More specifically, we explore the measurability of the heliocentric momentum enhancement caused by the DART impact. The $\beta_{\odot}$ value for kinetic impact-based planetary defense missions can be estimated as part of a least-squares orbit determination filter. This allows $\beta_{\odot}$ to be computed alongside the target body's state and other relevant orbital parameters such as the $A_2$ transverse nongravitational acceleration due to the Yarkovsky effect as defined in Marsden et al. (1973). To reliably estimate $\beta_{\odot}$, one needs to separate its contribution from that of the Yarkovsky effect, which can introduce a semimajor-axis drift on the order of hundreds of meters per year (Farnocchia et al. 2013; Del Vigna et al. 2018; Greenberg et al. 2020; Dziadura et al. 2022). At the same time, the $\beta_{\odot}$ from the DART impact causes a near-instantaneous change in semimajor axis. Therefore, $A_2$ and $\beta_{\odot}$ are highly correlated, and isolating the signal for both parameters is a challenge when it comes to estimating them.

Section 2 presents the methods used to understand the measurability of the heliocentric momentum enhancement, including the setup of the least-squares estimator and the treatment of astrometric data. This section also describes how we validated the orbit propagator and orbit filter used to generate the results in this work. In Section 3 we present our results on measuring $\beta_{\odot}$ for DART using the current observation data, as well as simulated future observations. Those are based on various observation scenarios that allow the $\beta_{\odot}$ value of the DART impact to be determined. This section also discusses the correlation between the $A_2$ Yarkovsky effect parameter and $\beta_{\odot}$. Following this, Section 4 discusses observation scenarios that would enable a solid $\beta_{\odot}$ estimate in future kinetic impact-based planetary defense missions. Finally, Section 5 summarizes our results and provides the key points we would like the reader to take away from this work.

## 2. Methods

If we assume that a kinetic impact causes a (near) instantaneous change in the momentum vector of an asteroid, then $\beta_{\odot}$ is defined as

$$ \Delta \mathbf{p} = \beta_{\odot} \cdot \mathbf{p}_{\text{DART}}, \qquad (1) $$

where $\Delta \mathbf{p}$ is the vector change in momentum of the asteroid, $\beta_{\odot}$ is once more the heliocentric momentum enhancement parameter, and $\mathbf{p}_{\text{DART}}$ is the relative momentum vector of the impactor spacecraft with respect to the target asteroid. This equation assumes that the escaping ejecta momentum travels in the opposite direction as the impactor momentum relative to the target asteroid. The change in the momentum of the asteroid

<figure><img src="image_4.png" /><figcaption>Figure 1. Illustration of pre- and post-DART orbits and their associated uncertainties for the Didymos system.</figcaption></figure>

can be interpreted as a discontinuity separating the pre-impact orbit from a new post-impact orbit. Both pre- and post-impact orbits of the target body have uncertainties that can overwhelm the deflection signal for some time. This leads to an inability of getting a statistically significant estimate for $\beta_{\odot}$, as shown in the results presented in Section 3.1.

Figure 1 shows a schematic diagram of the nominal pre- and post-impact orbits of the Didymos system and the corresponding uncertainties. Once the impact happens, the heliocentric orbit of the target body is altered. New observations following the impact then track this new orbit and reduce its uncertainties. Given enough time post-impact and a sufficient number of high-quality observations, the two orbits become distinctly separate if the uncertainties are small enough. The change between the pre- and post-impact orbits is the result of the heliocentric momentum change of the system, the magnitude of which is parameterized by $\beta_{\odot}$.

A full three-dimensional formulation of the momentum enhancement parameter is discussed in Statler et al. (2022). As mentioned earlier, Equation (1) assumes that the escaping ejecta momentum travels in the opposite direction to the impactor momentum relative to the target asteroid. High-fidelity ejecta dynamics models can instead be used to determine the direction of this ejecta momentum. We will focus on finding the magnitude of $\beta_{\odot}$ first and discuss determining the direction of the net escaping ejecta momentum vector in Section 2.3.

Since the heliocentric momentum change can also be thought of as an impulsive $\Delta V$ to the system, $\beta_{\odot}$ can be directly estimated via differential correction, an iterative orbit determination process that makes use of least-squares fitting of trajectories to observations to determine free parameters. Setting $\beta_{\odot}$ as a solve-for parameter as opposed to a calculated quantity allows us to directly understand its dependencies on the other estimated parameters and allows us to easily add or remove directional constraints on the total impulse delivered by the KI spacecraft.

### 2.1. GRSS Least-squares Filter

This article presents the Gauss-Radau Small-body Simulator (GRSS),<fnref n="11" /> an open-source small-body orbit propagation and orbit determination tool. The orbit determination algorithm in the GRSS library is based on the least-squares filter, the algorithm for which comes from Vallado (2022). The goal of a

<footnote n="11">¹¹ GRSS is publicly available at https://github.com/rahil-makadia/grss</footnote>
least-squares filter is to minimize the cost function,

$$J = \mathbf{b}^T \mathbf{W} \mathbf{b}, \qquad (2)$$

where $J$ is the cost and $\mathbf{b}$ is the $N$-dimensional vector of residuals of each observation. The $i$th element of this vector is

$$\mathbf{b}_i = \mathbf{O}_i - \mathbf{C}_i. \qquad (3)$$

$\mathbf{O}$ is the vector of observables, $\mathbf{C}$ is the vector of computed observations, and $\mathbf{W}$ is the $N \times N$ weight matrix. The observation data used in GRSS come from two sources: (1) optical astrometry from the International Astronomical Union's Minor Planet Center (MPC),<fnref n="12" /> and (2) radar astrometry from the JPL Small-Body Radar page.<fnref n="13" /> Additional measurements, such as stellar occultations and spacecraft measurements, can be ingested using the MPC's ADES format.<fnref n="14" />

The computed measurements are acquired by propagating an initial $N_{par}$-dimensional vector of solve-for parameters (where $N_{par}$ is the number of fitted parameters), $\mathbf{X}_0$ consisting of the cometary orbital elements (Milani & Gronchi 2009), any nongravitational parameters according to the Farnocchia et al. (2013) transverse acceleration model, and, in our case, $\beta_\odot$ to the epoch of each observation. This initial guess is updated through a process called differential corrections. The correction at every iteration is calculated using the equation

$$\Delta \mathbf{X} = \mathbf{P} \mathbf{A}^T \mathbf{W} \mathbf{b}, \qquad (4)$$

where $\mathbf{P} = (\mathbf{A}^T \mathbf{W} \mathbf{A})^{-1}$ is the covariance matrix of the solution and $\mathbf{A}$ is the gradient of the observations with respect to the solve-for parameters. GRSS calculates this gradient using first-order central differences. This corrected state (calculated as $\mathbf{X}_0 + \Delta \mathbf{X}$) is computed until the change in nominal solution is within a prescribed tolerance. The solution is said to have converged at that point. This classical process of nonlinear weighted least-squares estimation can be traced back to the work of C. F. Gauss (Gauss 1809).

The GRSS library uses the IAS15 algorithm, which is a 15th-order integrator based on Gauss–Radau spacings, which also gives the library its name (Everhart 1985; Rein & Spiegel 2014). The propagator submodule within GRSS uses a force model that incorporates the following:

1. Newtonian point-mass interactions (Newton 1687);

2. Einstein–Infeld–Hoffman point-mass relativistic effects (Einstein et al. 1938; Moyer 2003);

3. Zonal harmonics accelerations from the Sun and Earth (Laplace 1798; Kaula 1966); and

4. $A_1, A_2, A_3$ nongravitational acceleration terms in the heliocentric radial, transverse, and normal directions, respectively (Marsden et al. 1973).

The Newtonian accelerations are calculated for the Sun, planets, the Moon, Pluto, and the Big16 main-belt asteroid perturbers from the DE441 ephemerides. On the other hand, the relativistic accelerations are only calculated for the Sun, planets, the Moon, and Pluto. Marsden et al. (1973) define the nongravitational acceleration terms. However, in this work, only the acceleration in the transverse direction is calculated. This is done using an inverse square expression of the form

$A_T = A_2(1 \text{ au/r})^d$, where $r$ is the heliocentric distance and $d=2$ is chosen to match the level of absorbed solar radiation (Farnocchia et al. 2013). It is possible to vary this exponent for a better orbit fit as was done by Chesley et al. (2014). However, in the case of Didymos, retaining $d=2$ already provides an excellent orbit fit to all the observations.

Before the optical astrometry is passed through the orbit determination code, it needs to be treated to account for star catalog biases and weighted appropriately to ensure statistical reliability of the orbit solution given by the least-squares estimator. We use the astrometric catalog debiasing scheme from Eggl et al. (2020) to account for biases in the star catalogs used to reduce the optical astrometry reported to the MPC.

In addition to the debiasing, the optical astrometry acquired from the MPC also needs to be assigned weights before it is passed through the filter. We use the Vereš et al. (2017) weighting scheme to assign weights to observations depending on the observatory and the star catalog used for each observation. We also adopt the Vereš et al. (2017) rule of inflating the weights by a scaling factor of $\sqrt{N_{obs}/4$ for each observation if the number of observations from the same observatory on the same night is $N_{obs} > 4$. Once these optical observations are debiased and properly weighted, they are ready for the orbit filter. The radar astrometry uncertainty acquired from the JPL radar API is accepted at face value.

During the filtering process, one more algorithm is applied to the residuals. The GRSS orbit filter also has the capability to perform automatic outlier rejection according to the work of Carpino et al. (2003). This algorithm rejects and readmits optical measurements in the orbit fit based on a rejection threshold $\chi_{rej}$ and a recovery threshold $\chi_{rec}$. We set $\chi_{rej} = 3$ and $\chi_{rec} = 2.8$ for GRSS. The $\chi$ values for each optical measurement pair are calculated as

$$\chi_i^2 = \mathbf{b}_i \mathbf{\Gamma}_{\mathbf{b}_i}^{-1} \mathbf{b}_i^T, \qquad (5)$$

where $\mathbf{b}_i$ is the vector of residuals for the $i$th optical observation pair (R.A. and decl.) and $\mathbf{\Gamma}_{\mathbf{b}_i}$ is the residual covariance matrix of that pair of residuals, calculated as

$$\mathbf{\Gamma}_{\mathbf{b}_i} = \mathbf{\Gamma}_i \pm A_i \mathbf{P} \mathbf{A}_i^T. \qquad (6)$$

Here $A_i$ is the $2 \times N_{par}$ array of the gradient for the $i$th observation (the first row corresponds to the partial derivative of the R.A. with respect to the fitted parameters, and the second row is for the decl. measurement), and $\mathbf{\Gamma}_i$ is the $2 \times 2$ covariance matrix of the same optical observation. The second term on the right-hand side of the equation has a negative sign when the observation is included in the fit and a positive sign when it is rejected. In our implementation, the initial guess given to the filter is allowed to converge without any outlier rejection on the first pass. Once a converged solution is found, outlier rejection is turned on, and the filter will successively remove and/or recover any optical measurements with too high of a $\chi$ value. The weighted residuals are once again allowed to converge with this outlier rejection scheme turned on, and the converged solution after this stage is the output of the GRSS orbit filter.

In addition to the data from MPC and the JPL radar astrometry database, GRSS also has the capability to generate synthetic observations for the purposes of covariance analysis. These synthetic data are generated from a reference trajectory

<footnote n="12">¹² https://www.minorplanetcenter.net</footnote>

<footnote n="13">¹³ https://ssd.jpl.nasa.gov/sb/radar.html</footnote>

<footnote n="14">¹⁴ https://minorplanetcenter.net/iau/info/ADES.html</footnote>
<table><caption>**Table 1**
JPL Solution 204 and Its Associated $1\sigma$ Uncertainties for Didymos at 2023 January 12.5 TDB in Heliocentric J2000 Ecliptic Coordinates</caption><thead><tr><th>Element</th><th>Value</th><th>Uncertainty</th></tr></thead><tbody><tr><td>e (·)</td><td>0.3832947329867846</td><td>1.33837E-10</td></tr><tr><td>q (au)</td><td>1.012947259273837</td><td>2.75166E-10</td></tr><tr><td>t<sub>p</sub> (MJD TDB)</td><td>59873.790 058 610 495</td><td>9.56587E-08</td></tr><tr><td>Ω (deg)</td><td>72.99221760945375</td><td>2.18542E-07</td></tr><tr><td>ω (deg)</td><td>319.5559920731379</td><td>2.51502E-07</td></tr><tr><td>i (deg)</td><td>3.414168800300961</td><td>1.61880E-08</td></tr><tr><td>A<sub>2</sub> (au day<sup>-2</sup>)</td><td>-1.042309666988643e-14</td><td>6.49716E-16</td></tr></tbody></table>

**Note.** *e*: eccentricity; *q*: perihelion distance; *t*<sub>p</sub>: time of perihelion passage; Ω: longitude of ascending node; ω: argument of perihelion; *i*: Inclination; A<sub>2</sub>: transverse nongravitational acceleration parameter.

and can be generated with or without any constant biases or Gaussian noise around the reference trajectory. The synthetic observation feature is heavily used in this work to investigate the number and type of observations needed to measure β<sub>☉</sub>.

## 2.2. GRSS Validation

Before the GRSS orbit propagator and orbit filter can be used to estimate β<sub>☉</sub>, they must be validated against tools that have been demonstrably successful in the past. We individually tested the propagator and orbit determination code with the heliocentric trajectory of the (65803) Didymos system barycenter. This was a natural choice, since the Didymos system was also the target of the DART mission. Additional testing with different asteroids was also done and is shown in the examples available on the GRSS documentation.<sup>15</sup>

To validate the propagator, we took the latest Didymos orbital solution from the JPL Small-Body Database<sup>16</sup> (Solution 204 at the time of writing) and propagated it for 3000 days into the future. The initial solution is shown in Table 1. The difference in the final propagated position when compared to JPL's internal small-body propagation code was 2.4 m. Another longer-term propagation of 250 yr exhibited differences on the order of 1 km, which can be attributed to differences in the numerical scheme used to integrate the equations of motion.

Testing the orbit determination code is inherently more comprehensive, as it accumulates errors from both the integration and the observation model. The observation data set for Solution 204 includes a diverse set of observation types. There are more than 3000 optical measurements from ground-based observatories that include optical astrometry from telescopes around the world, radar delay measurements from the Goldstone Solar System Radar facility and Arecibo Observatory, optical measurements from the DART spacecraft before impact, and ground-based stellar occultation measurements. These stellar occultation measurements provide high-quality R.A.-decl. astrometry, and therefore very accurate observation models must be used when fitting these data.

Table 2 shows the orbit solution produced using GRSS. The GRSS solution and the JPL Solution 204 show excellent agreement with each other, with all orbital elements within $1\sigma$ of each other. Additionally, the Bhattacharya distance, which is a measure of the distance between two statistical distributions,

<table><caption>**Table 2**
GRSS Orbit Solution and Its Associated $1\sigma$ Uncertainties for Didymos at 2023 January 12.5 TDB in Heliocentric J2000 Ecliptic Coordinates</caption><thead><tr><th>Element</th><th>Value</th><th>Uncertainty</th></tr></thead><tbody><tr><td>e (·)</td><td>0.38329473307756273</td><td>1.32739E-10</td></tr><tr><td>q (au)</td><td>1.01294725920442</td><td>2.86507E-10</td></tr><tr><td>t<sub>p</sub> (MJD TDB)</td><td>59873.790058557846</td><td>7.75412E-08</td></tr><tr><td>Ω (deg)</td><td>72.99221758229581</td><td>2.00783E-07</td></tr><tr><td>ω (deg)</td><td>319.555992043250 06</td><td>2.39655E-07</td></tr><tr><td>i (deg)</td><td>3.4141688210436767</td><td>1.59297E-08</td></tr><tr><td>A<sub>2</sub> (au day<sup>-2</sup>)</td><td>-1.009633133758832e-14</td><td>6.09611E-16</td></tr></tbody></table>

**Note.** *e*: eccentricity; *q*: perihelion distance; *t*<sub>p</sub>: time of perihelion passage; Ω: longitude of ascending node; ω: argument of perihelion; *i*: Inclination; A<sub>2</sub>: transverse nongravitational acceleration parameter.

<figure><img src="image_5.png" /><figcaption>**Figure 2.** Postfit residuals of the optical observations used in the Didymos orbit solution produced by GRSS. The rejected observations are marked with a red circle.</figcaption></figure>

is only 0.060. The complement to this distance is the Bhattacharya coefficient, which measures the similarity between two statistical distributions, which is 0.942 in this case. These two metrics indicate that the two solutions are exceedingly close. Given a similar input, orbit solutions computed by the GRSS library are providing results of comparable quality to those of the JPL solar system dynamics group. Figure 2 contains further support for this statement. Here, the postfit residuals of the optical measurements used by GRSS are shown. The red circles around the scatter points indicate that those residuals were automatically removed by the outlier rejection scheme. We see that most observations with high residuals are automatically rejected. The few that remain are weighted such that, even with the high residuals, they do not have a meaningful effect on the final solution.

The remaining differences in the solutions can be explained by the fact that the JPL solutions sometimes incorporate manual weighting of the observations on top of the Vereš et al. (2017) rule. Additionally, the JPL solution assumes a 1s uncertainty on the timing of all observations. This is not modeled by GRSS, thus resulting in a lower uncertainty than if time uncertainty were taken into account. There are also some differences introduced by the fact that JPL propagates analytical partials along with the state and GRSS uses numerical

<footnote n="15"><sup>15</sup> https://rahil-makadia.github.io/grss/</footnote>

<footnote n="16"><sup>16</sup> https://ssd.jpl.nasa.gov/tools/sbdb_lookup.html#/?sstr=didymos</footnote>
<figure><img src="image_6.png" /><figcaption>**Figure 3.** Postfit correlation values from the covariance matrix of the Didymos orbit solutions produced by JPL (left) and GRSS (right).</figcaption></figure>

derivatives. This leads to small differences in the covariance matrix between the two solutions. However, Tables 1 and 2 show that the uncertainties are quite close to each other. Additionally, Figure 3 shows the correlation factors for each of the two orbit solutions. The two sets of correlation factors are in good agreement with each other. This indicates that it is not just the $1\sigma$ uncertainties on the diagonals of the covariance matrices that agree well, but it is both covariance matrices that are very similar to each other.

GRSS does not just perform well at fitting optical data. Figure 4 shows the postfit $\chi$ values of each observable in the fit. We see that the radar delay values (shown in green) are closely packed around 0, indicating that the GRSS observation model does not have trouble fitting radar data either. The radar data from 2003 have lower rms values than the ones in 2022 because they were regenerated using a shape model and therefore have lower uncertainties on the measurement (Naidu et al. 2020). Due to the lower uncertainties, the filter assigns a tighter weight to the 2003 delay points, which leads to better unnormalized residuals.

## 2.3. DART Kinetic Impact Geometry

A consequence of the DART impact on the Didymos system barycenter is a momentum transfer of magnitude $\beta_\odot$ times the momentum of the DART spacecraft at impact. Figure 5 depicts the direction of this impulse and its associated formal $1\sigma$ uncertainty using the observation scenario from Section 3.5, as well as the directions of the velocity vector of the DART spacecraft and the Didymos system barycenter at the time of impact. The blue vector represents the direction of the velocity of the DART spacecraft (relative to the Didymos system barycenter) at the time of impact, and the red arrow represents the direction of the velocity of the system barycenter relative to the Sun at the same epoch. The uncertainty of the delta-v imparted by DART is shown by the green ellipsoid. More details on how this ellipsoid was calculated are presented at the end of Section 3.5.

The minor axis of the ellipsoid is significantly smaller than the major and intermediate axes. Additionally, the minor axis is oriented along nearly the same direction as the DART spacecraft velocity at impact. In this work, we therefore make the assumption that the resulting impulse is instantaneously applied in the same direction as the DART spacecraft's velocity relative to the Didymos system barycenter at impact. Fixing the direction of the impulse in a direction that cuts through the smallest axis of the ellipsoid helps the filter when estimating $\beta_\odot$.

Since this work focuses on measurability of $\beta_\odot$, rather than on an actual estimate of $\beta_\odot$, this assumption of alignment of the

<figure><img src="image_7.png" /><figcaption>**Figure 4.** Postfit weighted residual values of the observations used in the Didymos orbit solution produced by GRSS. The rejected observations are marked with a red circle.</figcaption></figure>

<figure><img src="image_8.png" /><figcaption>**Figure 5.** Visualization of the DART impact circumstances and the resulting impulse on the Didymos system barycenter. The uncertainty ellipsoid is shown in an edge-on (left) and face-on (right) orientation. Note that the direction vectors have been scaled to a magnitude of $1E-7$ km s$^{-1}$. The ellipsoid reflects the observational circumstances of Section 3.5.</figcaption></figure>

DART velocity vector and the impulse to the barycenter is justified. Measurement of $\beta_\odot$ from the DART impact would need a better estimate of this impulse direction. This direction would then be informed by high-fidelity ejecta dynamics simulations and observations of the ejecta (Fahnestock et al. 2022; Li et al. 2023). Conversely, accurately estimating the three-dimensional $\beta_\odot$ properties can help better constrain the ejecta plume as well.

## 3. $\beta_\odot$ Estimation Scenarios for DART

Can we detect $\beta_\odot$ for DART? Could additional observations help us reach that point? In this section, we investigate a set of scenarios that aim to answer these questions and enable the detection of $\beta_\odot$ for the DART impact. Constants and parameters used in this section come from Daly et al. (2023), the DART Design Reference Asteroid (DRA) in Chabot et al. (2024), and the SPICE kernels for the DART mission.<fnref n="17" /> These are listed in Table 3. Moreover, it is important to note that all the $\beta_\odot$ standard deviations presented in this section are the formal $1\sigma$ uncertainties from the least-squares estimation.

<footnote n="17">¹⁷ https://naif.jpl.nasa.gov/pub/naif/pds/pds4/dart/dart_spice/</footnote>
<table><caption>**Table 3**
DART Impact Circumstances and Didymos System Constants</caption><thead><tr><th>Constant</th><th>Value</th></tr></thead><tbody><tr><td>UTC time</td><td>2022 Sep 26 23:14:24.183</td></tr><tr><td>Spacecraft mass</td><td>579.4 kg</td></tr><tr><td>Relative velocity [ICRF, x]</td><td>3.573 km s<sup>-1</sup></td></tr><tr><td>Relative velocity [ICRF, y]</td><td>-4.642 km s<sup>-1</sup></td></tr><tr><td>Relative velocity [ICRF, z]</td><td>-1.856 km s<sup>-1</sup></td></tr><tr><td>Didymos system mass</td><td>6.05e11 kg</td></tr></tbody></table>

### 3.1. Current Observation Data

The first question that needs to be answered when looking at $\beta_{\odot}$ measurability is whether it can be estimated with the already-existing observation data. Therefore, we took all of the Didymos astrometry up until the time of writing (2023 August) and attempted to estimate $\beta_{\odot}$. This data set included all the observations in the JPL Solution 204 data set and more than 2000 additional ground-based optical astrometry points from 2022 August to 2023 June from the MPC.

Results in this case showed a $\beta_{\odot}$ estimate of $-4.8 \pm 16.5$, indicating a signal-to-noise ratio (S/N) of about 0.3. This shows that the existing data are not enough to give a statistically significant estimate. Since $\beta_{\odot}$ is the link between the pre- and post-impact orbits of the Didymos system, both of them need to be well constrained to obtain a significant $\beta_{\odot}$ estimate. Since the post-impact observation arc is less than 9 months long at the time of writing this article, a $\beta_{\odot}$ estimate is not possible at this time. Therefore, we look at the future observation scenarios needed to get a reliable $\beta_{\odot}$ estimate.

### 3.2. Case 1: New Occultations in 2024

In the first observation scenario, we consider the possibility of a handful of new stellar occultation measurements at a monthly cadence between 2024 June and October. The 2024 June start time was chosen since there is a good chance of getting successful occultation measurements in the second half of 2024. This prediction is based on the density of star background that would enable an occultation measurement, regardless of the duration and magnitude of the occulted star. Additionally, the uncertainty of the current orbit of Didymos propagated to the beginning of its observability period in 2024 is only ≈1 Didymos body radius. The 2024 October cutoff was chosen because that is the planned launch date of ESA's Hera mission (Michel et al. 2022), and we wanted to explore the possibility of estimating $\beta_{\odot}$ before the launch.

These synthetic future stellar occultation measurement uncertainties were assumed to be at 2.5 mas. This is slightly higher than the average uncertainty of the existing Didymos occultation measurements (Chesley et al. 2023), to ensure a conservative estimate. These observations were generated from the JPL Solution 204 for Didymos as a reference trajectory. A reference value of $\beta_{\odot} = 3$ was set as the “true” value to assess the accuracy of the estimates returned by the orbit filter. This assumption is justified because the nominal $\beta$ value (when not accounting for Dimorphos mass uncertainties) for the DART impact is about 3.6. $\beta_{\odot}$ must be lower than the $\beta$ value since only a subset of ejecta that escaped Dimorphos could have escaped the system's Hill sphere.

Additionally, even though it is not statistically significant, the $\beta_{\odot}$ estimate from Section 3.1 is less than 0.5σ away from

this assumption. This implies that $\beta_{\odot} = 3$ is not an unreasonable assumption to make. It is, however, important to stress that this is not the $\beta_{\odot}$ value from the DART impact, but a reference value to assess how well we could measure $\beta_{\odot}$ under different scenarios. The $\beta_{\odot}$ estimate in this case was $3.2 \pm 0.7$, at an S/N ~ 4.9. Therefore, occultations in the second half of 2024 already give a long-enough observation arc for a reliable $\beta_{\odot}$ estimate before the Hera launch.

### 3.3. Case 2: New Occultations in 2024 and 2027

After analyzing the effects of the 2024 occultations on the $\beta_{\odot}$ S/N, we decided to add more occultations in 2027, when the Hera spacecraft will be at the Didymos system. Therefore, five additional monthly occultations from 2027 February to June were added to the observation data set for Case 1. The $\beta_{\odot}$ estimate in this case was $3.0 \pm 0.4$, at an S/N ~ 8. Therefore, the 2027 occultations increased the S/N by a factor >1.5. However, all five of these additional occultations might not be possible in early 2027, so we decided to explore another case.

### 3.4. Case 3: New Hera Measurements in 2027

After the Hera mission's planned launch in 2024 and a 2 yr cruise, the spacecraft would arrive at the Didymos system in late 2026 and would begin operations in early 2027. As part of these operations, the spacecraft would regularly communicate with radar ground stations on Earth. Some of these tracking passes could then be converted into pseudorange measurements of the asteroid barycenter, as was done with the OSIRIS-REx mission at Bennu (Farnocchia et al. 2021). Case 3 considers converting some of these Hera tracking measurements into pseudorange measurements of the Didymos system barycenter at a monthly cadence from 2027 February to June.

These pseudorange measurements are added to existing observations without any new occultations. The expected uncertainties are ~5 m, based on detailed numerical simulations of the Hera Radio Science experiment (Zannoni et al. 2018). However, all synthetic pseudorange measurements of the barycenter of the Didymos system in this study are weighted at an uncertainty of 15 m to ensure a conservative estimate. This is also conservative compared to the 2 m uncertainties of the Farnocchia et al. (2021) OSIRIS-REx pseudorange points and should not produce overly optimistic S/N values. In this case, the $\beta_{\odot}$ estimate is $3.0 \pm 0.2$, with S/N ~ 17. This increase in S/N is due to the fact that the Hera pseudorange measurements provide extremely accurate astrometry, and these measurements are taken at a point where the post-impact arc is almost 5 yr long.

### 3.5. Case 4: New Occultations in 2024 and 2027 with New Hera Measurements in 2027

The last observation scenario considered for estimating the $\beta_{\odot}$ for DART is the superset of the three future observation scenarios considered up to this point. These data include all existing observations and assume that new occultation measurements in 2024 and 2027 will be successful and that the Hera mission will provide pseudorange measurements for the Didymos system barycenter. In this scenario, the estimate of $\beta_{\odot}$ is $3.0 \pm 0.2$, with an S/N ~ 17. This is the best estimate of $\beta_{\odot}$ (negligibly better S/N than Case 3), as it is based on the largest set of high-precision observations.
<table><caption>**Table 4**
DART $\beta_{\odot}$ and $A_2$ Estimation Summary for Various Simulated Future Observation Scenarios</caption><thead><tr><th>Scenario</th><th>$\beta_{\odot}$</th><th>$\sigma_{\beta_{\odot}}$</th><th>S/N<sub>$\beta_{\odot}$</sub></th><th>$A_2$ (au day<sup>-2</sup>)</th><th>$\sigma_{A_2}$ (au day<sup>-2</sup>)</th><th>S/N<sub>$A_2$</sub></th><th>Corr ($\beta_{\odot}$, $A_2$)</th></tr></thead><tbody><tr><td>Case 1</td><td>3.21</td><td>0.66</td><td>4.85</td><td>-1.03E-14</td><td>7.21E-16</td><td>14.21</td><td>0.57</td></tr><tr><td>Case 2</td><td>3.03</td><td>0.39</td><td>7.84</td><td>-1.03E-14</td><td>6.99E-16</td><td>14.75</td><td>0.67</td></tr><tr><td>Case 3</td><td>3.01</td><td>0.18</td><td>16.75</td><td>-1.04E-14</td><td>4.19E-16</td><td>24.78</td><td>1.00</td></tr><tr><td>Case 4</td><td>3.03</td><td>0.18</td><td>17.01</td><td>-1.03E-14</td><td>4.23E-16</td><td>24.43</td><td>1.00</td></tr></tbody></table>

Note. The target value to be retrieved in these simulations is $\beta_{\odot} = 3$, $A_2 = -1.04E-14$ au day$^{-2}$.

The S/N value in this case has another point of validation. In Figure 5, the ellipsoid of the delta-$\nu$ uncertainty is calculated using this observation scenario. A separate capability of GRSS that can be used to estimate unconstrained $\Delta V$ vectors was used to obtain the 3 × 3 subcovariance matrix corresponding to the DART $\Delta V$. The eigenvalues of this covariance matrix yield the uncertainty of this $\Delta V$ (and therefore the dimensions of the ellipsoid), and the eigenvectors give us the orientation of the uncertainty ellipsoid.

Using the estimate of this unconstrained $\Delta V$ vector and its uncertainty, we find that the theoretical maximum S/N value of the $\Delta V$ magnitude for this observation scenario is 20.3. We consider the S/N of the magnitude of the $\Delta V$ vector instead of the components since this is similar to estimating the scalar $\beta_{\odot}$ value. (The S/N for the components of the $\Delta V$ vector in such a scenario would be much smaller owing to a lack of directional constraint.) Additionally, there is a small 16° offset between the eigenvector corresponding to the maximum S/N and the assumed direction of the escaping ejecta momentum in this work. Due to this offset, the maximum achievable S/N for this observation scenario reduces to 19.5. This maximum possible S/N is only 15% higher than the one we achieve with the assumption about the direction of the escaping ejecta momentum, and therefore it serves as a sanity check for the $\beta_{\odot}$ S/N values we get from the least-squares estimator. Therefore, the assumption about the direction of the escaping ejecta is not significantly affecting the measurability of $\beta_{\odot}$.

### 3.6. Summary and Discussion

Table 4 shows the summary of the simulated estimates for $\beta_{\odot}$ and the $A_2$ transverse nongravitational acceleration parameter for all four cases described in this section. We see the expected progression of the S/N as more observations of the Didymos system are added from Case 1 to Case 4. Additionally, we present the fitted values and S/Ns for the fit of the $A_2$ parameter along with $\beta_{\odot}$ for each of the cases. This is presented because the two parameters are highly correlated (as shown by the last column of Table 4). The Yarkovsky effect causes a secular semimajor-axis drift that has been detected in many other asteroids, with Bennu having the best detection owing to the OSIRIS-REx mission (Farnocchia et al. 2021). On the other hand, the DART impact impulsively changed the orbit of the Didymos system barycenter. So when estimating both parameters simultaneously, the filter has to disentangle the semimajor-axis drift from the Yarkovsky effect, $\dot{a}$, from the semimajor-axis change $\Delta a$ due to the DART impact.

To further illustrate this correlation, Figure 6 shows the 3$\sigma$ covariance ellipses generated using the 2 × 2 subcovariance matrix for the four cases. We see the uncertainty ellipses get smaller as the S/N on both $\beta_{\odot}$ and $A_2$ increases when more observations are added in each scenario. The correlation

<figure><img src="image_9.png" /><figcaption>**Figure 6.** $\beta_{\odot}$-$A_2$ covariance mapping for all four future observation scenarios considered in this work.</figcaption></figure>

between the two parameters in Cases 1 and 2 is 0.568 and 0.670, respectively. Once the Hera pseudorange measurements are added to the estimate in Cases 3 and 4 (highlighted in Figure 7), the correlation between $\beta_{\odot}$ and $A_2$ jumps to 0.996 in both scenarios.

Additionally, Figure 8 shows the complete correlation matrix of the orbit fit for Case 4. The strong correlation between $\beta_{\odot}$ and $A_2$ is even more apparent in this depiction. This highlights the importance of having an accurate picture of the acceleration due to the Yarkovsky effect when it comes to estimating $\beta_{\odot}$ for kinetic impact-based missions. In the case of Didymos, as well as many other asteroids, such as Apophis, an accurate estimate of the Yarkovsky acceleration has been enabled by either radar delay measurements or stellar occultations, or both. It is important to note that all near-Earth asteroids that have had successful occultation campaigns at the time of writing have been observed first by radar facilities. For this reason, when considering $\beta_{\odot}$ estimation scenarios for a future KI mission, we consider these two types of measurements both before and after impact.

### 4. $\beta_{\odot}$ Estimation Scenarios for Future KI Missions

A key takeaway from the previous sections is that the orbit uncertainties for the Didymos system barycenter were reduced significantly with the aid of the successful occultation
<figure><img src="image_10.png" /><figcaption>**Figure 7.** $\beta_{\odot}-A_2$ covariance mapping for Cases 3 and 4 (zoomed-in view of Figure 6).</figcaption></figure>

<figure><img src="image_11.png" /><figcaption>**Figure 8.** Postfit correlation values from the covariance matrix of the orbit solution used to estimate $\beta_{\odot}$ in Case 4.</figcaption></figure>

campaigns. Combined with the excellent radar astrometry provided by radar delay measurements, occultations could aid the success of future planetary defense missions, e.g., through a successful measurement of $\beta_{\odot}$.

Therefore, here we consider a future DART-like KI mission at Didymos. An interesting exercise in this respect is to investigate how pre-impact stellar occultation measurements would have affected our ability to determine the heliocentric outcome of the DART experiment. The weights of the occultation measurements in this scenario correspond to 25% of the angular diameter of Didymos at the time of the synthetic occultation. This is similar to the average weight at 23% of the angular diameter of Didymos for existing occultation

measurements. This exercise also assumes that the pre-impact orbit is good enough to allow for successful occultation observations.

In this scenario, we retain Case 1 from Section 3 but assume that the first 8 out of 16 occultations used in the fit for Solution 204 occurred on the first day of the 8 months leading up to impact, i.e., eight occultations on the first day of the month from 2022 February to September. The time interval over which these eight occultations occur can be compressed or extended as appropriate for an occultation opportunity, but what matters here is the number of high-accuracy occultation measurements. In this scenario, with the five additional occultation measurements in 2024 corresponding to Case 1, the S/N of $\beta_{\odot}$ increases by more than a factor of 2 to S/N ~ 12.

In another observation scenario, the pre-impact occultations can also help to determine how early we can get a meaningful (S/N ≥ 3) estimate of $\beta_{\odot}$. In this scenario, we still retain the assumption of splitting the first 8 out of 16 occultations to happen before the impact from 2022 February to September. Only one additional occultation in 2023 September would already provide an estimate of $\beta_{\odot} = 3.0 \pm 0.9$, with an S/N ~ 3. This means that pre-impact occultations would have reduced the amount of time it would have taken to obtain a statistically significant estimate of $\beta_{\odot}$ by more than a year.

However, all asteroids with successful occultation measurements have first been observed by ground-based radar stations. (The dependence on radar data to predict occultation by NEOs is expected to change as new high-accuracy asteroid astrometry data come in from Gaia measurements,<fnref n="18" /> but it remains true for now.) This is because radar measurements can reduce the uncertainty of the orbit solution to a point where successful occultation campaigns can take place. Therefore, to minimize the time between impact and a significant $\beta_{\odot}$ estimate (which is enabled by stellar occultation measurements), the target asteroid should be observed using radar (which can enable successful stellar occultations in the first place). This means that it is a combination of the complementary pre- and post-impact radar and stellar occultation measurements that are necessary to get a $\beta_{\odot}$ estimate for a KI mission as quickly as possible. Therefore, future KI missions should plan a combination of radar and occultation measurements to increase the likelihood of a rapid and reliable assessment of mission success.

## 5. Conclusions

We explored the measurability of the heliocentric momentum enhancement parameter ($\beta_{\odot}$) for kinetic impact-based asteroid deflection missions for the specific case of the DART mission. Making use of the least-squares estimator in the publicly available GRSS library, we found that the heliocentric momentum enhancement factor $\beta_{\odot}$ does not show a signal with the observation set in hand as of 2023 August. However, if future high-precision observations are added, such as occultation measurements in the second half of 2024, estimates for $\beta_{\odot}$ from the DART impact on Dimorphos would become statistically significant. Adding pseudorange measurements of the Didymos system using the Hera spacecraft in 2027 would result in detections of $\beta_{\odot}$ with an S/N of 17.

<footnote n="18">¹⁸ https://www.cosmos.esa.int/web/gaia/news; Gaia Data Releases 4 and 5 are expected in 2025 and 2030, respectively.</footnote>
We find that pre-impact occultation measurements would have had a significant impact on how quickly a significant $β_⊙$ estimate could have been achieved with DART. Successful occultation observation campaigns for the Didymos system in the year before the impact could have resulted in a solid detection of $β_⊙$ as early as 2023 September. Therefore, we conclude that for future KI missions, stellar occultation measurements and radar observations that enable the occultation measurements should be planned for and considered an integral part of the observation campaign before and after impact. This would enable a rapid and statistically significant assessment of the momentum transferred during a kinetic impact. In turn, this would help us to conclude whether Earth is safe by judging whether a future impact event has been safely averted.

## Acknowledgments

R.M. and S.E. acknowledge funding from a NASA Space Technology Graduate Research Opportunities (NSTGRO) award, NASA contract No. 80NSSC22K1173. This work was supported in part by the DART mission, NASA contract No. 80MSFC20D0004 to JHU/APL. The work of S.R.C., D.F., and S.P.N. was carried out at the Jet Propulsion Laboratory, California Institute of Technology, under a contract with the National Aeronautics and Space Administration (No. 80NM0018D0004). The authors also acknowledge the ACROSS Collaboration for predicting occultation events and conducting successful campaigns that led to high-accuracy astrometric measurements. The ACROSS project is supported under the OSIP ESA contract No. 4000135299/21/NL/GLC/ov. This study makes use of data obtained by the Observing Working Group of the DART Investigation Team. D.S. thanks Action Fédératrice Gaia of the Paris Observatory for financial support and acknowledges support as a Fulbright Visiting Scholar (2022–2023) at UC Berkeley. Part of this work was supported by the Programme National de Planétologie (PNP) of CNRS-INSU cofunded by CNES and the BQR program of the Observatoire de la Côte d’Azur. The authors also thank two anonymous reviewers for their comments and suggestions.

## ORCID iDs

Rahil Makadia https://orcid.org/0000-0001-9265-2230

Steven R. Chesley https://orcid.org/0000-0003-3240-6497

Davide Farnocchia https://orcid.org/0000-0003-0774-884X

Shantanu P. Naidu https://orcid.org/0000-0003-4439-7014

Damya Souami https://orcid.org/0000-0003-4058-0815

Paolo Tanga https://orcid.org/0000-0002-2718-997X

Kleomenis Tsiganis https://orcid.org/0000-0003-3334-6190

Masatoshi Hirabayashi https://orcid.org/0000-0002-1821-5689

Siegfried Eggl https://orcid.org/0000-0002-1398-6302

## References

Ahrens, T. J., & Harris, A. W. 1994, in Hazards Due to Comets and Asteroids, ed. T. Gehrels, M. Matthews, & A. Schumann (Tucson, AZ: Univ. Arizona Press), 897

Carpino, M., Milani, A., & Chesley, S. R. 2003, Icar, 166, 248

Chabot, N. L., Rivkin, A. S., Cheng, A. F., et al. 2024, PSJ, in press

Cheng, A. F., Agrusa, H. F., Barbee, B. W., et al. 2023, Natur, 616, 457

Chesley, S. R., Farnocchia, D., Nolan, M. C., et al. 2014, Icar, 235, 5

Chesley, S. R., Naidu, S. P., Bellerose, J., et al. 2023, in Asteroids, Comets, Meteors (Flagstaff, AZ), https://www.hou.usra.edu/meetings/acm2023/pdf/2314.pdf

Daly, R. T., Ernst, C. M., Barnouin, O. S., et al. 2023, Natur, 616, 443

Del Vigna, A., Faggioli, L., Milani, A., et al. 2018, A&A, 617, A61

Dziadura, K., Oszkiewicz, D., & Bartczak, P. 2022, Icar, 383, 115040

Eggl, S., Farnocchia, D., Chamberlin, A. B., & Chesley, S. R. 2020, Icar, 339, 113596

Einstein, A., Infeld, L., & Hoffmann, B. 1938, AnMat, 39, 65

Everhart, E. 1985, IAUC, 83, 185

Fahnestock, E. G., Cheng, A. F., Ivanovski, S., et al. 2022, PSJ, 3, 206

Farnocchia, D., Chesley, S., Vokrouhlický, D., et al. 2013, Icar, 224, 1

Farnocchia, D., Chesley, S. R., Takahashi, Y., et al. 2021, Icar, 369, 114594

Gauss, C. F. 1809, in Theoria Motus Corporum Coelestium in Sectionibus Conicis Solem Ambientium, ed. F. Perthes & I. H. Besser

Greenberg, A. H., Margot, J. L., Verma, A. K., Taylor, P. A., & Hodge, S. E. 2020, AJ, 159, 92

Holsapple, K. A., & Housen, K. R. 2012, Icar, 221, 875

Jewitt, D., Kim, Y., Li, J., & Mutchler, M. 2023, ApJL, 952, L12

Kaula, W. M. 1966, Theory of Satellite Geodesy (Waltham, MA: Blaisdell Publishing Company)

Laplace, P. S. 1798, in Traité de Mécanique Céleste, ed. J. B. M Duprat (Paris: L'Imprimerie de Crapelet)

Li, J. Y., Hirabayashi, M., Farnham, T. L., et al. 2023, Natur, 616, 452

Makadia, R., Raducan, S. D., Fahnestock, E. G., & Eggl, S. 2022, PSJ, 3, 184

Marsden, B. G., Sekanina, Z., & Yeomans, D. K. 1973, AJ, 78, 211

Michel, P., Küppers, M., Bagatin, A. C., et al. 2022, PSJ, 3, 160

Milani, A., & Gronchi, G. 2009, Theory of Orbit Determination (Cambridge: Cambridge Univ. Press)

Moyer, T. D. 2003, Formulation for Observed and Computed Values of Deep Space Network Data Types for Navigation, JPL Deep-Space Communications and Navigation Series (New York: Wiley)

Naidu, S., Benner, L., Brozovic, M., et al. 2020, Icar, 348, 113777

Newton, I. 1687, Philosophiae Naturalis Principia Mathematica (London: Jussu Societatis Regiae ac Typis Josephi Streater)

Rein, H., & Spiegel, D. S. 2014, MNRAS, 446, 1424

Richardson, D. C., Agrusa, H. F., Barbee, B., et al. 2022, PSJ, 3, 157

Richardson, D. C., Agrusa, H. F., Barbee, B., et al. 2024, PSJ, submitted

Rivkin, A. S., Chabot, N. L., Stickle, A. M., et al. 2021, PSJ, 2, 173

Statler, T. S., Raducan, S. D., Barnouin, O. S., et al. 2022, PSJ, 3, 244

Thomas, C. A., Naidu, S. P., Scheirich, P., et al. 2023, Natur, 616, 448

Vallado, D. A. 2022, Fundamentals of Astrodynamics and Applications (5th ed.; El Segundo, CA: Microcosm Press)

Vereš, P., Farnocchia, D., Chesley, S. R., & Chamberlin, A. B. 2017, Icar, 296, 139

Zannoni, M., Tommei, G., Modenini, D., et al. 2018, AdSpR, 62, 2273