Autonomous optical navigation using nanosatellite-class instruments: a Mars approach case study

This paper examines the effectiveness of small star trackers for orbital estimation. Autonomous optical navigation has been used for some time to provide local estimates of orbital parameters during close approach to celestial bodies. These techniques have been used extensively on spacecraft dating back to the Voyager missions, but often rely on long exposures and large instrument apertures. Using a hyperbolic Mars approach as a reference mission, we present an EKF-based navigation filter suitable for nanosatellite missions. Observations of Mars and its moons allow the estimator to correct initial errors in both position and velocity. Our results show that nanosatellite-class star trackers can produce good quality navigation solutions with low position (<300m\documentclass[12pt]{minimal} \usepackage{amsmath} \usepackage{wasysym} \usepackage{amsfonts} \usepackage{amssymb} \usepackage{amsbsy} \usepackage{mathrsfs} \usepackage{upgreek} \setlength{\oddsidemargin}{-69pt} \begin{document}$$<300\,\text {m}$$\end{document}) and velocity (<0.15m/s\documentclass[12pt]{minimal} \usepackage{amsmath} \usepackage{wasysym} \usepackage{amsfonts} \usepackage{amssymb} \usepackage{amsbsy} \usepackage{mathrsfs} \usepackage{upgreek} \setlength{\oddsidemargin}{-69pt} \begin{document}$$<0.15\,\text {m/s}$$\end{document}) errors as the spacecraft approaches periapse.


Introduction
Spurred by scores of successful Earth-orbiting applications, designers are beginning to launch nanosatellites on interplanetary missions. As an example, the upcoming 2018 Mars Insight Although other researchers have derived similar navigation filters and applied them to comparable scenarios, there are a number of distinct characteristics that set this study apart. We base our optical design (sensitivity, accuracy, etc.) on the performance that can be expected from nanosatellite-class star trackers. The MDIS-WAC examined by Christian and Lightsey (2010) has a collecting area similar to our baseline design but most other studies differ significantly. Deep Space 1 optics had a narrow field of view and a threshold apparent magnitude of 9 Riedel et al. (1990) and Stastny and Geller (2008) use a value of 9.5 or brighter. Both of these factors suggest fairly long focal lengths and often large collecting areas as well. In turn, the design necessitates frequent re-orientations of the sensor to sequentially image different targets. Such maneuvers may not be possible near periapse when the spacecraft attitude may be constrained by science or orbit injection requirements.
Over the following sections, we develop our model for autonomous navigation during a Mars approach scenario. We start by defining the approach geometry and deriving our primary performance metrics. Our analysis assumes that the accurate state knowledge at periapse is of primary importance. We then introduce the Hikes et al. measurement model (Hikes et al. 2017) and examine how known star tracker noise characteristics can be used to estimate the effective noise encountered in LOS vectors and other measured quantities. Using Monte Carlo simulation, we evaluate the filter performance during the reference scenario. We assess the filter convergence properties, the overall navigational accuracy, and the relative importance of the available measurements.

Reference mission
We begin the development of our estimator with a detailed discussion of the reference mission. This includes definitions of the basic reference frames used in our derivations, an analysis of the approach geometry, as well as the configuration of the spacecraft and the navigation sensors. In particular, we explore how the optical design parameters typical of nanosatellite star trackers-such as focal length, f -number, field of view (FOV) size, detector noiseinteract with OpNav measurement models.

Approach geometry and frames of reference
Our reference scenario concerns the terminal phase of a direct interplanetary transfer from Earth to Mars. Somewhat arbitrarily, we begin the scenario as the spacecraft crosses the sphere of influence (SOI) boundary and approaches Mars on a planeto-centric hyperbolic trajectory (see Fig. 1). The basic frame of reference for analysis is an inertial frame aligned with the theoretical periapse direction at SOI entry. This is designated as the P-frame. The origin is located at the Mars center of mass. The P 1 axis points in the direction of periapse and the P 3 axis lies in the direction of the orbit angular momentum. The P 2 axis completes the right-handed coordinate frame. The definition of the P-frame as a periapse frame is only strictly true at the start of the scenario. Any subsequent perturbations from two-body motion will move the actual periapse. Nevertheless, this frame will remain a convenient reference frame for the remainder of the simulation time.
As the inbound spacecraft approaches on a nominally hyperbolic trajectory, the state of the spacecraft can be defined by the inertial, Cartesian position, r P , and velocity, v P . In our notation, the subscripts indicate the reference frame for a vector-the P-frame in this case. Although our estimator (and the truth simulation) uses Cartesian state representations internally, calculating the corresponding orbital elements of the approach trajectory offers insight into spacecraft behavior. The shape of the hyperbola is defined by the semimajor axis, a, and the eccentricity e. See Appendix C for details on calculating these quantities from the position and velocity values.
The hyperbolic asymptotes will be found at an angle of η ∞ with respect to the P 1 axis, where: We select the negative solution to ensure that the spacecraft lies on the approach leg of the hyperbola. Assuming that the Mars approach occurs at the aphelion of the heliocentric transfer ellipse, then the angle to the Sun, η Sun will be given by: Before further describing the mission geometry, we define several frames of reference that are used throughout this paper. These frames of reference are shown in Table 1. The B-frame definition will depend on the specific pointing rule used on the spacecraft. The sensor frame is denoted S, but if the spacecraft is equipped with more than one star tracker, then this notation can be subscripted, i.e., writing S i to refer to the ith sensor.
We can disregard Frame-I and assume that all inertial measurements and ephemeris predictions are made with respect to Frame-P. We assume that nominal approach trajectory lies in the ecliptic. This allows us to write: where α M is the right ascension of Mars in Frame-E, and C j represent principal axis rotations. (See "Appendix A" for the sense of these rotations.) The in-plane angle η Sun depends on the initial conditions and orients the periapse frame within the ecliptic. The rotation from Frame-P to Frame-O (see Fig. 2) is strictly a function of the spacecraft state. We can write this transform down directly in terms of basis vectors (rather than principal axis rotations): where these basis vectors are defined by where a × is the skew-symmetric matrix of vector a (see "Appendix A"). The last vector is then found from O 2 P = O × 3 P O 1 P . The transform, C B O , captures the effects of pointing rules on spacecraft attitude. In our nominal scenario, we consider two pointing rules: a Mars-pointing mode that keeps the limb of Mars a specified angular distance, ξ , off-center from the sensor boresight, and a velocityaligned mode that keeps the B 3 frame aligned with the velocity vector.
This definition of a Mars-pointing attitude gives the following sequence of rotations: where ρ is the angular radius of Mars, i.e., For velocity-aligned pointing, the transform is most easily expressed by concatenating the basis vectors of B and expressing them in Frame-O: where we define B 1 using the angular momentum and B 3 from the linear velocity: The last set of transforms involve the rotations between S and T . If the center of Mars is at an angle of ρ c from the boresight and in an azimuthal direction θ c , then we can define:

System dynamics
The implementation of our navigation filter relies on two pieces of software: the simulator predicts the true value of the spacecraft state and the synthetic sensor outputs; the estimator implements the navigation filter itself. Both the simulator and the estimator must propagate position, velocity, and attitude throughout the scenario. For the translational motion, we assume that the spacecraft undergoes unpowered flight, perturbed by acceleration disturbances. We also assume attitude motion is actively controlled and follows a positiondependent reference pointing trajectory. We do not explicitly model this control system but instead present a simple model of the attitude error dynamics. The orbital position of the spacecraft is propagated in the periapse frame. (The E-frame would also be suitable.) We use a Cartesian formulation for the dynamics equations.
v P = a Mars,P + δa i,P .
Herev P is the net spacecraft acceleration in the P-frame, a Mars,P is the gravitational acceleration due to Mars, and the δa i,P terms represent other perturbations. The primary gravitational term is given by: a Mars,P = − μ Mars In this study, we introduce third-body gravitational effects from the Sun and Jupiter. These perturbations take the form: where Δr DS is the vector from the disturbing body to the spacecraft and Δr DM is the vector from the disturbing body to Mars (for simplicity we have omitted the P-frame subscript from this expression). The third-body perturbations are present in the truth model for the simulation, but not included in the online estimator. We recognize that including Sun and Jupiter perturbations in the onboard estimator would be prudent for a fielded system. However, we adopt this particular division between truth and estimator dynamics to capture the worst-case modeling error. If the estimator performs well in this case, it will likely improve further with third-body perturbations included into the onboard estimator leaving only higher-order effects unmodeled.
The attitude error dynamics employed in the truth simulation is an approximation to the performance of a real system. We assume that the spacecraft attitude C B P can closely track a reference attitude trajectoryC B P . Representing the attitude error as an error rotation vector, φ (see "Appendix B"), we can write The reference attitude can be decomposed into two separate rotations,C B P = C B O C O P , the former an expression of the pointing rule and the latter a function of the orbit position and velocity. Although this form offers insight into the construction of the overall transformation, it obscures the fact thatC B P ≡C B P ( x).
We model the attitude error dynamics as a first-order system driven by random noise: where w φ as a zero-mean Gaussian random variable with known covariance, Q φ .

Spacecraft state
The state of the spacecraft during the simulations is derived from a number of parameters that govern dynamic motion and spacecraft configuration. Translational motion is captured by estimating position, r P , and velocity, v P -both measured in the periapse frame. The reference attitude is expressed using theC B O rotation, and attitude error is captured with an error rotation vector, φ. The framework incorporates measurements from one or more star trackers, each of which has both a known, nominal, body-frame orientation,C S B , and an unknown mounting error ψ i (one for each sensor). Together the state vector is:

Sensor models
Our sensor model focuses on measurements that can be extracted from the partially illuminated disk of the target planet. Limb (or horizon) measurements have served as the basis of some types of OpNav systems for many years ;Owen Jr et al. (2008) provides a historical perspective and Christian (2015) and Bu et al. (2017) highlight more recent approaches. In addition to standard measurements of the planet's apparent size and position in the field of view, we have identified a number of additional measurements that can be extracted from star tracker images. These include: -Sensor-referenced directions to Phobos and Deimos -Attitude with respect to an inertial frame (Frame-P) -Phase angle between the spacecraft, the Sun, and Mars (by measuring the illuminated fraction of Mars) -The orientation angle of the terminator on the surface of Mars.
Our measurement are based on a pinhole camera model of image formation. We recognize that this is an oversimplification and that real cameras must contend with higher-order aberrations. However, if cameras are well calibrated, the resulting images can be transformed to remove the aberrations and restore the geometric properties of an ideal pinhole camera. This approach helps to preserve the generality of our discussion, particularly with respect to our navigation results. Many references to camera calibration in the literature: Brown (1971) is an early reference on the subject; Weng et al. (1992) present one of the most common models; and Liebe (2002) describes how calibration helps ensure star tracker accuracy.
In addition to modeling sensor response, we must also model measurement availability during the simulation. As targets enter and leave the field of view, we must add and drop the corresponding rows from the measurement sensitivity matrix, H. Measurement availability is governed by the following heuristics: 1. If the Sun lies within a fixed exclusion angle (35 • ) of the sensor boresight, then no measurements are available. 2. Phobos and Deimos are only available when their apparent visual magnitude is brighter than a threshold value (from our reference star tracker design, we use a conservative value of 5.75). 3. Attitude measurements are available when three or more stars are visible in the field of view (and not blocked by Mars). 4. Phobos and Deimos are not visible during eclipses and transits (i.e., when they pass in front or behind Mars). 5. The Mars position can be estimated as long as the part of the illuminated limb is within the FOV. 6. Phase and terminator angles can only be measured when the entire disk of Mars lies within the FOV.
These rules are intended to be conservative-specialized image processing can likely be used to extend the availability of some types of measurements based on mission needs. For convenience of discussion, images that contain the entire planetary disk are termed Case-I, whereas images that contain only a portion of Mars are termed Case-II (see Fig. 3).
In the remainder of this section, we describe the specific models used to predict these observations and introduce state-dependent error models that capture the variation of measurement noise with imaging conditions. Measurement updates in the estimator are incorporated using concatenation of the available measurements: and the measurement sensitivity matrices: Each component of the sensitivity matrix has the form:

Attitude measurements
Although the inertial orientation of each sensor can be expressed as: the actual orientation of the spacecraft body frame is based on a reference attitude derived from the spacecraft's state estimate. Thus, we write We can further expand the sensor-to-body transforms in terms of the nominal mounting orientationC S B and the mounting error, ψ: and thus We introduce the pseudomeasurement Δ defined to represent the (first-order) error rotation vector between the ideal and measured attitudes.
Equating (25) and (26) gives: We can show that by transforming φ into the S-frame, we can reorder the transforms and remove the effect of the nominalC S B rotation.
The ideal Δ vector can be extracted from the RHS matrix elements. Gaussian attitude measurement noise, ν Δ , can be added to give the noisy attitude error pseudomeasurement Similarly, the estimated pseudomeasurement is extracted from: The observations models in our estimator use the nonlinear rotations in (28) and (30). In order to evaluate the partial derivatives, we make first-order approximations to the RHS. Thus we can show Hence and

Direction measurements
Direction vector observations of Solar system bodies is a common type of OpNav measurement. In this scenario, Phobos and Deimos observations are processed using this model. These rely on the displacement vector s S between the spacecraft and the target.
We can calculate the true vector from: where the target position, ρ P , can be predicted from ephemeris calculations. In this study, we use JPL's SPICE toolbox (Acton 1996) for all of our ephemeris calculations. Now the predicted value observations must be calculated differently. The estimator can predict the P-frame ephemerides accurately, but transforming these into the S-frame relies on the state estimates: We assume that we can estimate the detector coordinates of the center of the body, correcting as necessary for partial illumination of the body. A simple pinhole model for the camera gives: Pittelkau's (2003) relationship for the measurement sensitivity matrix for focal-plane measurements allows us to calculate sensitivity with respect to the state components. The general form of these matrices are calculated from: where, after dropping the reference-frame subscript, Thus, the nonzero partial derivatives are This model is most suitable when the targets are close to being point sources. In our baseline scenario, Phobos and Deimos never exceed an apparent angular radius of two pixels allowing the images of the moons to be localized and identified in the same manner as an imaged star. Other orbital trajectories might bring the spacecraft much closer to the moons, but only in the very final approach will the moons have appreciable size.
Algorithms that track stars through multiple images (e.g., see Samaan et al. 2005) could be adapted to track the moon locations in the detector images and distinguish them from background stars. When the line of sight to the target center of mass diverges from the center of illumination, corrections to our simple observation models may be necessary to preserve measurement accuracy. If the apparent size of the target remains small-e.g., a handful of pixels across-a simple correction for the asymmetric illumination might be sufficient, but very large targets must be treated as extended objects.

Position measurement
When estimating the direction to an extended target such as Mars, we must explicitly account for the effects of illumination and distortion within the field of view.
Although forming these pseudomeasurements requires additional processing, the apparent size of the target provides some additional information (i.e., range to the target) that is not present in a pure direction measurement.
"Appendix D" outlines our approach to accurately estimating ρ, ρ c and θ c from the detector-plane observations. The range, r , can be calculated from (8), and the full pseudomeasurement for s S is calculated using (96). We recognize that our approach to determining s S from the detector image is not particularly novel-it merely serves as a useful placeholder for the low-level image processing. Examples abound in the literature of different techniques including (Christian 2015(Christian , 2016Bu et al. 2017).
The Mars observation vector s S can be written in terms of our filter state. The true, noisy observations are given by and the modeled measurement by The measurement sensitivity matrices are then

Phase angle measurements
The illuminated fraction of Mars gives a fairly direct estimate (Enright 2010) of the angle χ between the Sun, Mars, and the observer. Working in the plane defined by the Sun, the spacecraft and Mars, we can show that if b is the semi-minor axis of the terminator curve (see Fig. 4), then The true value can be obtained from vector calculations The measurement sensitivity can be found by taking the (implicit) partial derivatives of both sides of this expression and rearranging to get:

Terminator angle measurements
The terminator angle refers to the orientation of the terminator ellipse on the imaged disk of Mars in the sensor field of view (see Fig. 4). Although this ellipse is a construct of the planet being projected to the sensor, we define it in vector space for simpler relation with other geometries. The major axis of this ellipse is perpendicular to the plane that contains the spacecraft, the Sun, and Mars. Expressed as a unit vector, the major axis must point in the direction: The observed terminator angle is caused by the projection of this vector into Frame-T (correcting for projective distortion). This can be expressed simply in terms of the components of κ T : where The partial derivatives needed to form H are of the form and using (12), the corresponding partial derivatives of κ T are then The rotation C T S is state dependent, so we must also evaluate the corresponding partial derivatives of ρ c and θ c . We can show where Z P is the sensor boresight vector, rotated into the P-frame (i.e., the third row of C S P ). Thus, we can write See (91) for the identity used to evaluate the latter terms in (60) and (61).
The θ c relationships can be derived in a similar fashion and thus ∂θ c

Measurement error models
During the approach scenario, the appearance of the Mars disk grows and shifts in the FOV. It is reasonable to assume that the measurement errors present in the derived measurements will change as well. This section presents a series of error models designed to capture most of these effects. We employ a parametric model for the measurements of r S and an empirical model for errors in χ and β.
Mathematically, we consider all measurements to be of the form: where Thus, we wish to approximate the measurement covariance terms, R (x). Each type of measurement has its own error characteristics.

Attitude covariance
In general, star tracker attitude accuracy depends on the number and geometric distribution of the stars in view. However, for simplicity and generality, we adopt a simple and popular The errors in the cross-boresight direction, σ bs , are identical, while the roll error, σ roll , is usually substantially higher.

Direction covariance
The measurement accuracy for the moon direction measurements depends on the star tracker centroiding accuracy. Our reference trajectory does not approach either moon very closely, and as a result, the apparent radius of Phobos is only about two pixels at closest approach (Deimos is significantly smaller). Thus, measuring moon and star directions requires very similar image processing.
Centroid noise generally increases for dim stars. Figure 5 shows how the ST-16RT centroid error varies with brightness. This combines on-orbit sensitivity data with laboratory temporal noise measurements. These predictions are extrapolated beyond magnitude 4.0 and will be quite conservative for the dimmest detectable sources (i.e., high visual magnitudes).
The photometric brightness of planetary bodies is frequently approximated with the following relationship where V (1, 0) is the (tabulated) visual magnitude at opposition and a distance of 1 au. For this expression, the distances between Mars and the spacecraft, r , and Mars and the Sun, r Sun are also measured in au. In this paper, we have used V Phobos (1, 0) = 11.8 and V Deimos (1, 0) = 12.89. These were taken from the NASA-Goddard Mars Factsheet (Williams 2017). The phase function P (χ) is often approximated as the diffuse reflection from a Lambertian sphere (Lester et al. 1979) Using (69), (70), and the error relation of Fig. 5, we are able to approximate values of σ moon , expressed in pixels. The detector-plane covariance is then:

Position covariance
We have used the model developed by Hikes et al. (2017), for the Mars position covariance. This model provides a convenient parametric expression in terms of a few key quantities derived from geometric and detector properties. Adapting their model to the notational conventions of this paper gives the following expression: The quantities R and D are given by and where N limb is the number of points extracted from the limb curve, and ξ is the half-width angle of the illuminated limb, measured from the center of the planet. The latter quantity is calculated from the observation geometry at each time step; N limb is calculated to maintain a constant linear sample density in the FOV The covariance calculated in (72) is expressed in terms of a reference frame aligned with cusp and limb directions. We can transform this covariance into the S-frame so that we can use it in our filter calculations.

Empirical covariance models
The model presented by Hikes et al. assumes a fit to the illuminated limb and gives an estimate of the covariance of s S . Our measurements of phase and terminator angles rely on fitting curves to the terminator as well as the limb. To produce these secondary error models, we generated a large number of simulated images, analyzed them using image processing routines derived by Enright (2010), and calculated error statistics from the resulting estimates. We then formed a series of simple interpolants that can be used to approximate the measurement errors under different operational conditions. The models capture interesting error behaviors without the need to simulate individual images during estimator execution.
Close examination of the algorithm performance suggested that under most circumstances, variances of χ and β were generally uncorrelated with other errors and depended primarily on the range from Mars, r .
We recognize that limiting the error interpolants to depend only on r will fail to capture some distinctive sensor behaviors. The prior tests of this algorithm (Enright 2010) showed significant increases in error when the planet is completely illuminated (i.e., χ = 0) or barely illuminated (i.e., χ = π/2). Although these effects are generally important, they are largely irrelevant in our particular scenario. At the start of the simulation, enough of the Mars disk is visible to avoid the increased error regime, and as the spacecraft nears the χ = 0 point, Mars appears larger than the FOV and must be treated as a Case-II observation. For different sensor and trajectory combinations in which the whole Mars disk would still lie within the FOV when χ = 0, we would still recommend treating this situation as a Case-II observation due to the elevated errors measuring phase and terminator angles.

Noise estimates
Both the empirical and parametric error models depend on the limb and terminator noise parameter, σ pix . This quantity captures the variability in precisely locating the corresponding curves on the detector. In order to determine representative values for this noise parameter, we examined moon images from the ST-16RT as well as imagery from the Cassini and Voyager-1 and Voyager-2 missions. The selected images captures differences in phase angle, terminator angle, and apparent target size. Using Enright (2010) pixel-level routines to extract the boundary curves, we equated σ pix with the standard deviation in the residual errors after fitting curves to the limb and terminator. Depending on the image characteristic and threshold settings used for image processing, these trials predicted σ pix values in the range of 0.25 − 2.0 pixels. For the remainder of this paper, we have adopted the conservative value of σ pix = 2.0, although we are confident that algorithmic tuning could improve on this value.

Sample error predictions
Figures 6 and 7 show the principal 1 − σ errors at various ranges during approach. We have included data from both the parametric model and the empirical trials. Both sets of curves are virtually identical. Most error components get smaller as the spacecraft approaches Mars with the exception of the limb-ward error-when Mars grows very large in the FOV, it becomes difficult to distinguish between small changes in radius and the limb-ward position. Figure 8 shows the resulting errors in both χ and β. Both quantities exhibit similar behaviors, getting smaller as the spacecraft nears Mars, but the effect is particularly noticeable in β.
The error models presented in this section capture many of the important variations of sensor performance seen throughout the simulated scenario. We acknowledge that the quantitative error predictions are dependent on the specific image processing algorithms used as well as the preprocessing logic employed to identify perimeter pixels and discriminate between limb and terminator curves. Nevertheless, we feel that the quantitative results are credible and reasonably conservative predictions of the overall performance that we can expect.

Filter formulation
The navigation filter studied in this paper is based on a continuous-discrete extended Kalman filter (EKF) framework. In this section, we describe the important elements of the filter

1-Error (km)
Sun-ward Cusp-ward Range Sun-ward (Hikes, et al.) Cusp-ward (Hikes, et al.) Range (Hikes, et al.) formulation and highlight the ways in which it diverges from standard implementations. Generally, we integrate continuous dynamics in order to propagate state and covariance and make discrete measurement updates.

Filter measurements
Measurement updates for our navigation are applied using the standard equations for Kalman gain and state and covariance update:

1-Error (Radian)
The nominal set of filter measurements looks like The exact number of measurements available during these updates (and hence the specific form of h) depends on which targets actually lie within the sensor's FOV. Each step involves calculating the Kalman gain, the measurement sensitivity matrix, and applying updates to both state and covariance. Measurements from multiple sensors can be processed together.

Propagation
Our estimator uses a simplified version of the system dynamics described in Sect. 2.2. The attitude error dynamics uses the same first-order process detailed in (17), but the velocity dynamics omits the perturbation terms, instead using only the simplified dynamics of (14). To partially compensate for this modeling error, we introduce a process noise covariance term Q v , in the propagation of state covariance. Thus the estimator dynamics equations can be summarized as:ẋ These equations are used directly to propagate the state estimate. The state covariance is propagated using the standard EKF linearization: where the state Jacobian is given by and the combined process noise covariance is: We note in passing that care must be taken when numerically integrating (81) and (82). Many elements of P become quite small-particularly the attitude terms-and this may require some variable scaling (or careful choice of integration tolerances).

Experiments and simulations
In the preceding sections, we have derived a general framework for state determination during a hyperbolic Mars approach using a single star sensor pointing toward Mars. Here we evaluate the performance of this system and assess whether there are any fundamental limitations to using nanosatellite-class star trackers for this scenario. We use only a single sensor due to the fact that the only instance that is devoid of stars occurs when the planet completely fills the camera field of view. Otherwise, our Mars-pointing rule ensures that most of the field of view remains unobstructed, and that at least three stars are present at all stages of the approach to allow for conventional star tracking. Through these trials, we wish to establish the general effectiveness of nanosatellite optical navigation, the relative importance of the different measurements that we make, and the utility of our analysis framework to evaluate changes to the mission and satellite design.
Evaluating system performance relies on establishing relevant quantitative metrics. The estimated navigation solution can serve a variety of mission activities (e.g, flyby, orbit insertion, coordinated observations), so the most important metric will depend on the particular mission under study. We propose three metrics that will be useful in a variety of situations: -Position and Velocity Error. These are the magnitudes of the spacecraft position and velocity errors. They can be evaluated at periapse, or anywhere along the spacecraft trajectory. -Periapse Radius Error. This is the accuracy of the periapse radius estimate (as derived from the osculating orbital elements, see "Appendix C"). -Peritime Error. This is the accuracy of determining when the spacecraft reaches periapse.
It is most important when planning orbital insertion maneuvers. Like the periapse radius, it is calculated from osculating elements.
A series of Monte Carlo trials allow us to calculate mean and standard errors for our chosen performance metrics. For each trial, we vary the initial state estimates according to the distribution described by our initial state covariance. Varying the actual initial statei.e., through x 0 rather than x 0 -could also be illuminating, but we note that the spacecraft trajectory during a hyperbolic approach is quite sensitive to initial velocity. Small velocity changes can dramatically alter the approach trajectory, which in turn will affect measurement availability. Performance over very different scenarios then becomes hard to assess with simple statistics. For the following experiments, we run N trial = 5000 trials with each configuration. Thus, the standard error in the performance statistics is about 1%. This should be adequate for most types of preliminary analysis. Table 2 shows the basic parameters that govern the simulation. Our reference scenario is based on a hypothetical mission arriving at Mars in January 2019. Our sensor design is based on the demonstrated performance of the Sinclair Interplanetary ST-16RT star tracker. We make a few simplifying assumption (e.g., a circular field of view), but otherwise use fairly conservative estimates of the sensor precision. From our experience with ST-16RT flight data, sensor's attitude performance is not adversely affected by the presence of the bright Earth limb within the field of view.

Experiment 1: basic performance evaluation
The first set of tests evaluates the nominal behavior of the navigation filter and establish baseline expectations of performance. We enable all of the identified measurements in this trial and maintain a Mars-pointing orientation throughout the trajectory. The spacecraft is equipped with a single optical sensor, and the pointing rule ensures that a significant portion of the field of view remains unobstructed, permitting conventional star tracking at as well as the navigational measurements.
Most of the measurements are available throughout the approach. Figure 9 shows a schedule of the availability over time. The serendipitous timing of the scenario allows good observations of both Phobos and Deimos for most of the trajectory. Our simulations assume the moons to be point sources on the detector, and this is justified as the apparent size of the moons in the detector during the approach reaches a maximum of 2 pixels for Phobos and 0.4 pixels for Deimos. Occultations and transits are few and the moons orbit close enough to Mars that they spend most of their time in the sensor FOV. This ensures the near constant availability of these secondary measurements. In the last few hours before periapse, the moons are no longer reliably found in the field of view, and disk of Mars becomes too large to permit terminator and phase angle measurements. Figure 10 shows the position error over the course of the approach as well as the 3 − σ bound calculated from the diagonal elements of the state covariance. Although initially large, as the spacecraft nears Mars, the solution improves dramatically. The error bound drops below 20 km about 10 hr before periapse, and below 5 km in the final 2 h.
A spacecraft approaching Mars will rely on knowledge of its own trajectory to trigger mission-critical activities (e.g., insertion burns) as the spacecraft approaches periapse. To that end, we can assess the estimator solution and evaluate how well state estimate can predict the satellite's approach to periapse. Figure 11 shows the error in the estimator's prediction of the peritime. At the very start of the approach, the error is quite large-about 17 min.-but Periapse Radius (km) Estimate Truth this converges rapidly as the overall state estimate converges. In this trial, the error remains below 50 s in the final 24 h and lies between 0 and 5 s in the 9 h before periapse. The estimated periapse radius (Fig. 12) offers an alternate perspective on estimator performance. The convergence of the radius estimate is fairly slow, but still falls below 3 km about 10 h prior to periapse. Because the true periapse radius is essentially an osculating quantity, we see significant motion during early periods of the scenario when these third-body gravitational effects are most significant. Because the filter converges after the period of significant perturbation, it is difficult to assess the tracking performance of the estimator.
We can compare the typical trajectories discussed above with results derived from a set of Monte Carlo trials. The first row of Table 3 shows mean and standard deviation values for our baseline case. The statistics show that the selected trajectories are typical of the estimator's behavior. Where mean values are negligible, we have calculated the corresponding standard deviation. Although our key figures of merit are derived from the trajectory state-i.e., position and velocity-the estimator's ability to measure attitude and sensor misalignment is also important. The overall sensor alignment is accurate to about 15 μrad about three axes. In the trials, we did note a fair amount of variability in this quantity, and we attribute this to poor separability between time-varying attitude error and the constant mounting bias. Furthermore the net error, i.e., ψ + φ S , was more stable from trial to trial than the estimates of ψ alone.

Experiment 2: relative measurement contributions
The first experiment presented results from our baseline configuration using all available measurements during filter updates. Many of the individual measurements require their own specialized processing code-particularly those targets whose appearance changes significantly during the approach. Thus, it is important to understand the contributions that each makes to the overall system effectiveness. Measurements that offer little improvement in the navigation solution can be omitted from the overall system, saving both development effort and online processing time.
The first set of Experiment-2 rows in Table 3 considers estimator performance in the absence of a single type of measurement. Attitude and Mars position are considered essential and are not omitted from any of the trials. Estimator performance is only marginally affected by the loss of the terminator angle, but the moon directions and phase angles all make valuable contributions to accuracy. Estimator performance degrades considerably when any of these are unavailable, particularly for the position and velocity error magnitudes. The results suggest that the Position + Terminator results are worse than just the Position alone, but we suspect that this is merely an artifact of the Monte Carlo trials.
A second set of trials considers small sets of observations. The navigation solution is still viable when using only Mars position measurements albeit with significantly increased error. Adding in the Mars Phase or the moons restores virtually all of the full-measurement performance case.
In examining the performance data, we notice some differences from metric to metric. Position and velocity error are significantly affected by the estimator configuration while changes in peritime and periapse radius are relatively small. The misalignment error relies primarily on the attitude information and consequently is significantly affected by the loss of terminator angle. The low error value when the moons are excluded seems anomalous and may be an artifact of the Monte Carlo testing. The behavior of this quantity is rather erratic and examining some of the individual trials suggests that poor separability between the attitude error and sensor misalignment makes it difficult to estimate ψ reliably.

Experiment 3: pointing rules and sensor placement
In this last experiment, we examine how our optical navigation framework can be used as a design tool to evaluate different spacecraft configurations and operational strategies. We emphasize that our aim is to illustrate the utility of our optical navigation framework rather than present optimal or general results.
The first design trial considered the addition of a second star tracker mounted such that its boresight was aligned with the orbit normal. This second sensor does not make any Mars observations but provides additional attitude information. The position and velocity performance improved relative to the baseline Experiment-1 trial. Although a second sensor might not be justified based on performance alone, fusing attitude data from multiple sensors offers overall robustness to the system that are not captured by our selected metrics.
In the second design scenario, we consider the role of spacecraft attitude rules on estimator performance. The baseline pointing rule is designed to place the limb of Mars just off-center in the field of view. This rule is quite effective and as we saw in Experiment-1, many observations are available through periapse passage. However, in a real mission, other operational activities-e.g., orbit insertion, flyby observations-may impose their own attitude constraints as the spacecraft reaches periapse. To capture one possible constraint, we introduce a velocity-aligned pointing rule that kept the spacecraft z-axis aligned with the velocity vector. This might represent an appropriate attitude for orbit insertion. Both of these rules are defined in Sect. 2.1. Figure 13 shows the measurement availability for the velocity-aligned mission scenario. In contrast to Fig. 9, we see that Mars observations are only available up until about 4 h before periapse. Although the scenario begins with Mars in the field of view, it eventually drifts out of sight as the position and velocity vectors diverge. The estimator is left without any position-related updates during the terminal phases of the scenario leading to the dramatic increase in error detailed in Table 3. The degraded accuracy may be sufficient for some missions but represents about a 40-fold increase in position error.
The relatively poor navigation performance stems from the limited availability of Mars measurements during the final hours before periapse. By changing the star tracker mounting on the spacecraft bus, we can restore some of the performance lost in the naive case. In the nominal approach, Mars appears slightly off-center in the −Y direction (as seen in Frame-S), and then gradually drifts out of view (moving toward −Y ). Merely by canting the sensor so that Mars starts as far toward +Y as possible, we can prolong the amount of time that the planet remains in view. A mounting angle of 8.8 • will place the center of Mars 6.5 • off boresight at the beginning of the scenario. Figure 14 shows the resulting change in availability due to this fixed mounting angle. Instead of losing the Mars observations 4 h before periapse, the planet now stays in view until about t − t p = 1.2 hours. The estimator performance is worse (4.4 km in position), than the Mars-pointing trials but shows significant improvement over the original velocity-aligned case. We can investigate this approach further, increasing the sensor FOV to 10.5 deg and the cant-angle to 11.8 deg. We justify this change by considering that our baseline star tracker, the ST-16, has a rectangular, 7.5 × 10.5 deg, FOV. With this new design, the spacecraft can observe Mars up until about 50 min. before periapse. With this availability improvement, the position and velocity errors drop accordingly.
Although this specific result represents an interesting performance compromise, we feel that the more important message is that our estimator framework is responsive to this sort of trade-off analysis. If needed, we could consider alternate design and operational options (e.g., different fields of view, hybrid pointing rules) and evaluate the resulting impact on performance.

Conclusion
In this paper, we have developed an optical navigation framework that provides position and velocity estimates during a hyperbolic Mars approach. The estimator uses a Cartesian EKF to fuse a variety of measurements derived from observations of Mars and its moons. Viewing conditions vary greatly during approach and the error models capture corresponding variations in measurement noise caused by variations in target shape, illumination, and size. The baseline sensor sensitivity and accuracy was modeled on the performance of the Sinclair Interplanetary nanosatellite star tracker.
The estimator performance is quite encouraging, particularly when the satellite can maneuver to keep Mars in view throughout the whole approach. Position errors (1−σ ) are typically a few hundred meters and velocity errors are about 15 cm/s or better. The navigation solution's estimate of periapse arrival is extremely accurate and should be sufficient for direct entry or orbit insertion operations. Periapse radius estimates are also very good, and the estimator can successfully track unmodeled gravitational perturbations from the Sun and Jupiter.
All of the sensor measurements make useful contributions to the solution, but some are more valuable than others. Many of the more complex derived measurements such as phase and terminator angles are generally not worth the increased computational requirements. Simple filter formulations using only Mars position and moon directions baseline configuration. Our experiments also indicate that observations made during the final hours have a huge impact in performance. Operational constraints that repoint the spacecraft near periapse can significantly degrade performance but relatively simple changes in sensor mounting can mitigate these effects.
The trials examined in this study have reinforced the importance of measurements made during the final hours before periapse. Several unmodeled effects will affect the appearance of Mars in this flight regime. These include the oblateness of Mars and the distortion of the limb boundary due to the atmosphere. We are currently examining how these additional factors may affect system performance. We expect that improved image processing algorithms will be necessary in order to maintain the assumed performance. Furthermore, this expected increase in computational requirements will likely make some of the marginal value measurements, such as terminator angle, even less attractive.
The overall performance of our estimator is encouraging, but without a specific target mission, it is difficult to categorically say that the solution accuracy is good enough. Recent Mars missions such as the Mars Science Laboratory (MSL) can provide some insight. For example, the MSL navigation requirements (Martin-Mur et al. 2011) specified that the navigation system must provide 2 km error (3 − σ ) at atmospheric entry with a measurement cutoff of 6 h prior to arrival. Our navigation solution does not meet this specification; we can achieve the accuracy, but not with the specified cutoff. More extensive optical redesign (e.g., narrowing the field of view) might allow the navigation filter to achieve a low error solution further from periapse, minimizing the need for angled sensor mounting, and increasing the operational utility of the navigation solution. We recognize that there is room to further improve the navigational performance, but nevertheless feel that the demonstrated accuracy would be sufficient for other missions.
That nanosatellite-class star trackers can be used as effective guidance sensors for Mars approach carries significance beyond this specific scenario. If a small modest-accuracy star tracker can provide autonomous guidance using only Mars position and radius, then many other missions (e.g., lunar or asteroid missions) can potentially benefit as well. Thus, we see our estimator as a potentially enabling technology for a new class of innovative missions. Open Access This article is distributed under the terms of the Creative Commons Attribution 4.0 International License (http://creativecommons.org/licenses/by/4.0/), which permits unrestricted use, distribution, and reproduction in any medium, provided you give appropriate credit to the original author(s) and the source, provide a link to the Creative Commons license, and indicate if changes were made.

A Mathematical notation
In this section, we detail a few of the notational conventions used throughout this manuscript. For the reader's reference, Table 4 contains the convention for variable diacritics used in this paper. Throughout our development, we frequently use the unary cross-matrix operator. If a = a 1 a 2 a 3 T , then Finally, we denote a principal axis rotation about axisj as C j . The sense of these rotations provides the transform into the rotated frame. For example, if Frame-B is obtained from Frame-A by an x-axis rotation through an angle θ , then:

B Small rotations
In our representation of spacecraft state, we use error rotation vectors to represent small angle rotations in attitude error, φ and mounting error, ψ. A common, first-order approximation is to treat this a series of Euler angle rotations (these commute for small angles): In our testing, this representation led to problems with the conditioning of the state covariance matrix. Using a second-order rotational approximation gave better numerical behavior. Rather than a series of Euler angles, we now assume that φ represents an axis-angle rotation of the form φ = φa where a is a unit vector representing the axis of rotation. We know C (φ) = cos φI + (1 − cos φ) aa T − sin φa × , which to second-order is When we have measurements of the form we can compute the derivative using the rotational approximation. This gives:

C Basic orbit determination
Calculating the osculating orbital elements from the position and velocity is a convenient means of assessing the performance of the estimator. We use the common following common relations to calculate Keplerian elements from the position, r P , and velocity, r P . The vis-viva equation give a relation for the semimajor axis, a: The angular momentum of the orbit is: This allows us to find the semilatus rectum, p and the eccentricity, e.

D Planetary observations
Care must be taken when analyzing extended images of Mars to account for perspective projection. Even if we neglect planetary oblateness and assume that the target is spherical, the projection of the planetary disk onto the detector plane will not form a circular arc.
Other authors have shown that the resulting curve is an ellipse (Christian 2015) and that a correction is necessary to accurately calculate the direction vector to the planet, s S , from the measured ellipse center (Bu et al. 2017). Projection effects can be removed by fitting the horizon curve in vector space [e.g., see work by Christian (2017)], or by making appropriate transformations to the detector coordinates. In order to maintain a consistent approach across all of our different pseudo-measurements, we have chosen the latter approach and sketch out the necessary transformations here.
To re-project coordinates extracted from a detector image, we first express s S in spherical coordinates: We can relate this vector to the corresponding detector coordinates and then re-project these points into pseudo-Cartesian coordinates where radial distances from the origin are true arclength distances: This is only a partial solution to the projection problem-although radial distances from the boresight represent true arclengths, azimuthal distances are distorted. To permit a better fit to the target planet's limb, we shift the origin of our new coordinates to (ρ c , θ c ). This is equivalent to an azimuthal equidistant map projection. The new arclength, λ, and azimuthal, α, coordinates are then found from: cos λ = cos ρ cos ρ c + sin ρ sin ρ c cos (θ − θ c ) sin α = sin ρ sin (θ − θ c ) sin λ cos α = cos ρ − cos λ cos ρ c sin λ sin ρ c .
Similarly, the corresponding rectangular coordinates are: Provided our estimated center was close to the true center, circle fits to the z coordinates will find the true center and angular radius of the planet. Converting back to boresight coordinates relies on applying the spherical geometry relations one more time: cos ρ = cos λ cos ρ c + sin λ sin ρ c cos α sin (θ − θ c ) = sin λ sin α sin ρ cos (θ − θ c ) = cos λ − cos ρ cos ρ c sin ρ sin ρ c .