Mars moon ephemerides after 14 years of Mars Express data

The Mars Express (MEX) mission has been successfully operated around Mars since 2004. Among many results, MEX has provided some of the most accurate astrometric data of the two Mars moons, Phobos and Deimos. In this work we present new ephemerides of Mars' moons benefitting from all previously published astrometric data to the most recent MEX SRC data. All in all, observations from 1877 until 2018 and including spacecraft measurements from Mariner 9 to MEX were included. Assuming a homogeneous interior, we fitted Phobos' forced libration amplitude simultaneously with the Martian tidal k2/Q ratio and the initial state of the moons. Our solution of the physical libration 1.09 +/- 0.01 degrees deviates notably from the homogeneous solution. But considering the very low error bar, this may essentially suggest the necessity to consider higher order harmonics, with an improved rotation model, in the future. While most data could be successfully fitted, we found a disagreement between the Mars Reconnaissance Orbiter and the Mars Express astrometric data at the kilometer level probably associated with a biased phase correction. The present solution precision is expected at the level of a few hundreds of meters for Phobos and several hundreds of meters for Deimos for the coming years. The real accuracy of our new ephemerides will have to be confirmed by confrontation with independent observational means.


Introduction
Space astrometry of the natural satellites allows us to reach a precision incomparable with respect to the telescopic observations. Regarding Phobos, we can access a short history of these telescopic observations in Pascu et al. (2014). In order to exploit the combination of space data and ground-based ones for computing ephemerides, we received funding in the scope of the European FP7 program and coordinated the European Satellite Partnership for Computing Ephemerides (ESPaCE) (Thuillot et al. 2012) consortium from 2011 to 2015. ESPaCE brought together seven laboratories (IMCCE-Paris Observatory, Royal Observatory of Belgium ROB, Joint Institute for VLBI in Europe JIVE-ERIC, Technical University of Berlin TUB, German Aerospace Center DLR, French Space Center CNES, Technical University of Delft TUD), that focused on different activities related to the positioning of spacecraft and natural satellites such as digitizing of old telescopic observations, the exploitation of spacecraft radio tracking data to reconstruct spacecraft orbits and the development of planetary moons ephemerides. The ESPaCE portal at http://espace.oma.be provides access to further information. In this context, a fruitful collaboration has been established between the ESPaCE consortium and the MaRS radio-science and the HRSC camera teams of the MEX mission, in preparation of the December 2013 close approach of Phobos.
On the 29th of December 2013, Mars Express (MEX) came as close as 45 km to Mars' moon Phobos -closer than ever before. To guarantee the highest precision on Phobos' position, an astrometric observation campaign of Phobos was performed from MEX' on board camera SRC. Several tens of observations were performed and reduced using the UCAC4 star catalogue to provide strong constraints on Phobos' orbit in the ICRF reference frame (Pasewaldt et al., in prep.). New ephemerides of the Martian satellites were developed as part of the consortium. They are presented in this work.
In Section 2, we detail the numerical model we used. In Section 3, we present the whole set of observations that was used in the present work. The next section shows the astrometric results and provides physical parameters of interest like the amplitude of Phobos' forced libration and Mars' tidal k 2 /Q ratio. In section 5 we discuss an issue we encountered when considering Mars Reconnaissance Orbiter (MRO) observations. The last section is devoted to precision of the extrapolation and comparison with former ephemerides.

Modeling
We used a very similar approach as described in Lainey et al. (2007). In particular, we used the NOE numerical code to model the orbits of Phobos and Deimos. The Mars rotation model considered was the one of Konopliv et al. (2016) with the Mars A&A proofs: manuscript no.  gravity field MRO120D that we truncated at degree and order 12. Table 1 provides the numerical values of these coefficients. The masses of Mars (42828.3750104 km 3 /sec 2 ), Phobos (7.11 × 10 −4 km 3 /sec 2 ), Deimos (9.46 × 10 −5 km 3 /sec 2 ) and the Mars Love number k 2 (0.169) were borrowed from this gravity solution. The planetary ephemerides were INPOP19a (Fienga et al. 2020) that introduced the perturbations of the Sun, the Moon and all the other planets.
Three amendments on the dynamics were introduced in comparison to Lainey et al. (2007). The first two were the introduction of general relativity and tidal-cross effects (Lainey et al. 2017) that consist in the cross action of the tidal bulges raised on Mars by different tidal raising bodies. While rather negligible at the current level of accuracy for Mars, we added these effects for completeness. Tides on Mars' moons were neglected.
A third but important dynamical effect was the introduction of the forced libration on the rotation of Phobos (Jacobson (2010)). Indeed, such perturbation introduces a secular effect on the periapsis of Phobos (Borderies & Yoder (1990), Lainey et al. (2019)), barely maskable in the fitting procedure. The secular drift, associated with the quadrupole field of Phobos on its periapsis (or any moon in somewhat similar configuration) was given by Borderies & Yoder (1990) and recalled by Jacobson (2010) to be where a, e, n, M and ̟ are the traditional osculating Keplerian semi-major axis, eccentricity, mean motion and mean anomaly while R, J 2 = −c 20 and c 22 denote the mean radius and unnormalized gravity coefficients of Phobos. Here we explicitly see the action of the libration amplitude A on the periapsis. Our A is −A in Jacobson (2010). Of course, this action of the quadrupole field of the satellite is just one effect affecting the pericenter. The precession of the pericenter is mainly due to the flattening of Mars.
Since the gravity field of Phobos is hard to determine directly from the current data (Yang et al. 2019), we relied on  who used Phobos' shape and homogeneous density hypothesis to derive Phobos' gravity field. They obtained for the first gravity field coefficients c 20 = −0.066127 and c 22 = 0.009917, assuming a Phobos' radius of 14.0km.

Observation sets
Most data sets used in this work were already processed in Lainey et al. (2007) and Jacobson (2010), and benefited greatly from the astrometric catalogue of Morley (1989). A large astrometric database is available on the Natural Satellites DataBase (NSDB) server (Arlot & Emelyanov 2009). It turned out that some of the old measurements were extremely biased and therefore were eventually not included in the fit. Moreover, ground observations with residuals larger than 2 arcseconds were removed. A rejection criteria beyond 3σ was applied to ground observations. Last, spacecraft data were weighted in the same way as in Lainey et al. (2019). Here we present two important new data sets, not available to the former works.

Observations from digitized plates
Following the discovery of Phobos and Deimos, and the detection of the long suspected secular accelerations in their longitudes (Sharpless 1945), the United States Naval Observatory (USNO) began a thirty-year  program of photographic observations of the Martian satellites (Pascu 1977(Pascu , 1978(Pascu , 1979(Pascu , 2012. These observations were among the most accurate (Morley 1989) and were used to support all space reconnaissance projects of the Martian system.
Photographic observations were begun in 1967 and continued at every opposition through 1997. They were taken with the USNO 61-inch astrometric reflector in Flagstaff, Arizona, and the USNO 26-inch refractor in Washington, D.C. with the use of special filters. Several Schott 5 inch x 7 inch x 3 mm GG14 (yellow) filters were polished optically flat. In the center of each a small, thin metallic nichrome film having an optical density of about 3.0 was deposited by evaporation. The GG14 was chosen to accommodate both telescopes, and nichrome was chosen because it transmits neutrally in the visual bandwidth. The function of the small nichrome filter was to reduce the intensity of the planetary image to that of the satellites, producing a measurable image of the planetary disk. A number of Kodak emulsions were used, including 103aJ, 103aG, and IIIaJ. More details about the observations and astrometric results are available in Robert et al. (2015).
Four hundred twenty-five plates were selected and transmitted to Royal Observatory of Belgium (ROB) to be digitized (Robert et al. 2011;de Cuyper et al. 2011). Each plate contains two to three exposures shifted in the RA direction. Measured (x, y) plate positions were corrected for instrumental and spherical effects, and the reductions were performed using four or six suitable constant functional models (Robert et al. 2011(Robert et al. , 2014  The observed positions of Mars, Phobos, and Deimos were compared with their theoretical computed positions given by DE430 planetary ephemeris ) and NOE MarsSatV1_0 satellite ephemerides (Lainey et al. 2007). The key point was that the NOE MarsSatV1_0/DE430 astrometric residuals (O-C) for all observations had an rms of 47.8 mas, 60.5 mas and 50.7 mas, for Mars, Phobos and Deimos, respectively. Overall intersatellite RMS (O-C) is 39 mas. These rms correspond to the observation accuracies over thirty years, providing observations that are similar in accuracy to old space data but with a more Gaussian profile and a larger time span.

New MEX/SRC data
MEX is in a highly elliptical and nearly polar orbit about Mars reaching well beyond the almost circular path of Phobos (Jaumann et al. 2007). Due to the asphericity of the planet's gravitational field its longitude of ascending node is drifting westwards and its argument of pericenter is drifting in the direction opposite to the spacecraft motion. At the same time its orbital period is similar to the one of the inner satellite of Mars resulting in several flyby opportunities in consecutive orbits at intervals of about five to six months. Whenever a close encounter between the space probe and this Martian moon occurs high resolution observations of Phobos by the onboard cameras are planned, see section 3.2.1. In the past years new opportunities to observe Phobos and Deimos in conjunction with other so- Table 1. Mars gravity field used in the present work (MRO120D (Konopliv et al. (2016)) truncated at degree and order 12. Coefficients are normalised and Mars' radius is equal to 3396.0 km.  -0.32345607088980E-06 -0.16550385123550E-05 S 12,12 -0.89210620377340E-07 lar system objects or stars (so-called mutual events, see section 3.2.2) were identified and more regularly planned. These observations do not necessarily require the proximity to the body (Ziese and Willner 2018).
The Super Resolution Channel (SRC) is part of the High Resolution Stereo Camera (HRSC) onboard MEX. It fea-tures a compact Maksutov-Cassegrain optical system, a 1K by 1K interline-transfer CCD, and fast read-out electronics ).

Close encounters
During an approach maneuver the camera is pointed to a fixed location on the celestial sphere and is rotated around the boresight axis such that the lower and upper image borders are parallel to the relative velocity vector of Phobos w.r.t. the spacecraft.
In the usual eight-image sequences the first and the last pictures are long-time exposures to detect the faint light of background stars. Exposure times of the images in between are adjusted to the brighter Phobos surface. Pointing corrections in sample and line are derived from star images and then interpolated linearly for the pictures of Phobos .
For the purpose of astrometry, SRC images have been focussed using an image-derived point spread function in a Richardson-Lucy deconvolution (Michael & Neukum 2009). The image positions of stars have been determined by fitting a 2D Gaussian profile to the pixel value distributions (Duxbury 2012, pers. comm.). Pointing corrections have been derived w.r.t. the predictions of reference star positions based on the UCAC4 catalog. The directions to Phobos have been measured by applying the limb-fit approach based on the latest 3D shape model by .
Between May 2013 and March 2014 a larger number of flybys than usual have been performed in order to support the Phobos gravity field experiment. Therefore data from 38 approach maneuvers with distances ranging from 350 to 14 000 km have been evaluated. Summing all up, 340 astrometric observations of Phobos in right ascension and declination coordinates have been provided. Estimated uncertainties vary from hundreds of meters up to one km (Pasewaldt et al., in prep.).
Hitherto, the estimation of uncertainties is still incomplete. This is due to the occasional presence of spacecraft jitter during image acquisition (see Pasewaldt et al. (2015) for details). For sequences with continuous observations of reference stars or other celestial objects this causes no problem. But for several eight-image series with only two star pictures an additional uncertainty is introduced that is only hard to quantify. Therefore, we distinguished our Phobos measurements into those performed during linear variation and those carried out during nonlinear variation in camera pointing.
Besides the measurement of reference objects, the CCD line coordinate of a Phobos position could provide information on the presence or absence of spacecraft oscillations during image recording. Based on the imaging geometry described above this should be almost constant throughout a sequence. From at least three successive Phobos observations we can calculate two line coordinate differences.
While only positive or negative differences have been classified as linear variations, alternating positive and negative differences have been categorized as non-linear variations in pointing. If in addition the pre-fit residuals were deviating by more than three sigma from the arithmetic mean, we marked this as an outlier. For series with only two Phobos measurements we could not determine whether the outlier has been related to a non-linear change in pointing or not. Further on, in some images we fitted the shape model-derived limb to only very short limb point arcs.
This classification is not an equivalent substitute for a proper weight estimate, but it gives at least a suggestion for the following evaluation of our data.

Mutual events observations
Longer imaging sequences of mutual events showing either both Martian moons or one of the moons with Jupiter or Saturn in the background were analyzed by Ziese and Willner (2018). These observations can be obtained more frequently than direct observations during flybys. Moreover, a wide range of MEX orbit positions are possible. However, in contrast to close encounter images, Phobos and Deimos can also cross the image plane diagonally. It also occurs that the moons are only partially visible.
To determine the bodies' locations within the image, simulated images are computed. The simulations are based on the ephemerides model mar097, a rotational model by Stark et al. (see Archinal (2018) and Archinal et al. (2019)) and shape models of Phobos and Deimos derived by  and Thomas et al. (2000), respectively. To achieve a better agreement with the observation, the simulation is convolved with a point spread function describing the image distortion of the SRC. Here a subset of the PSF derived by Duxbury et al. (2011) was applied. Matching the illuminated simulation against the observation provides the bodies' line and sample coordinates with sub-pixel accuracy. For images showing either Phobos and/or Deimos with a far distant object, the derived right ascension and declination of the moon (as seen from MEX) are determined. The exhaustive list of MEX data we used is given in Table 4. In particular, due a current limitation of our software, mutual events evolving both moons simultaneously (i.e. without the presence of Jupiter or Saturn) have not been included in the current solution.

Adjustment results
All data were fitted using the RMS of each data set as a priori uncertainty. Data points higher than 2 arcsec for ground observations were systematically removed. Then, a 3 sigma rejection criteria was applied. For spacecraft observations, we used the published estimated precision. However, to assess a better weight, we rescaled these uncertainties by a scalar chosen to get about 66% of astrometric residuals within 1σ level for each satellite and coordinate.
We provide in Figures 1-3 and Tables 2-6 the astrometric residuals of all ground and space data used in the development of the present ephemerides as well as initial conditions of Phobos and Deimos after fit. The high precision of both the reduction of photographic plates and the new MEX data is confirmed, with typical residuals at the level of 40-50 mas and 300 meters, respectively. The combination of the existence of two Mars moons with the use of inter-satellite fitting often provides both opposite mean values and similar sigma. Such feature breaks up when both moons could not be observed always simultaneously, as indicated by a different number of observation per moon. We note that the MEX residuals may be associated with an error on the astrometric calibration of the field and the error on the spacecraft position. In particular, it was found that the long imaging sequence of mutual events (see subsection 3.2) introduced a significant bias in the fitting procedure. While the jitter effect of the spacecraft could be removed thanks to the mutual event opportunity, the error on MEX' position in space remained. An error of the order of a few hundred of meters on MEX' reconstructed orbit was found, much higher than the initially estimated error at the level of 20-25 meters (Rosenblatt et al. 2008). It appeared that the increase of the error on MEX' position in 3-dimensional space arose from an unfortunate consequence of the new tracking strategy of MEX. Indeed, since 2011 ESA has decreased the amount of tracking data around periapsis, entailing a poorer constraint on MEX' motion.
Thanks to the new sets, we formally gain an order of magnitude on the uncertainty (1σ) of the Martian tidal quality factor Article number, page 4 of 11 V. Lainey et al.: Mars moon ephemerides after 14 years of Mars Express data and the physical libration of Phobos. Assuming k 2 =0.169, we obtained Q = 93.68 ± 0.18. Since the eccentricity of Phobos' orbit is small (less than 2%), the secular variation of its semi-major axis can be easily approximated using (Kaula 1964 where m and M denote Phobos' and Mars' mass, respectively. From (2) we get for the secular acceleration on longitude 1 2 dn dt ≃ + 9n 2 k 2 mR 5 4QMa 5 .
Using the approximation n 2 a 3 ≃ GM, and the value 9378km for Phobos semi-major axis, we obtain for Phobos' tidal acceleration 1.24 × 10 −3 deg/yr 2 . A more accurate estimation can be determined by running the full numerical model with/without tidal dissipation inside Mars. This provides for Phobos' tidal acceleration (1.257 ± 0.003) × 10 −3 deg/yr 2 . In practice, this orbital decay and associated longitude acceleration is what allows us to estimate Mars' Q. It is noteworthy to recall that we do not include higher Love numbers (k 3 , k 4 ...) in our fit. Hence, our estimation of Q related to k 2 should in reality be more considered as a lumped Q parameter. Considering the high precision we now get on this tidal parameter (less than 0.2%), it would probably make sense to add at least k 3 . Nevertheless, it would be practically impossible to solve for two different Q, one associated with each tidal frequency. Phobos' physical libration was simultaneously found to be A = 1.09 ± 0.01 degrees (Table 7). This last value is at 5σ below its theoretical value of 1.14 degrees ) assuming homogeneity. Nevertheless, considering the very small error bar (1%), it is pretty clear that in the harmonics expansion of Phobos a higher order degree should be considered, as well as an improved rotation model like suggested by Rambaux et al. (2012), before drawing any conclusion on Phobos' interior.

Precision versus accuracy
While fitting our model to astrometric data, we discovered a disagreement between MRO 1 and MEX data. In particular, while both sets did provide small residuals independently, they appeared to be strongly deteriorated when considered together for Deimos. Since both data sets roughly cover the same years, an error on the modeling sounds unlikely. Unfortunately, we could not use Phobos data to get a deeper understanding of the possible bias in these sets since MRO data of Phobos are of much less quality. This bias was simultaneously observed at JPL and interpreted as a bias of 55. ′′ 775 and −7. ′′ 321 in right ascension and declination of MEX data, respectively (Robert Jacobson, private communication). Since this bias appeared roughly constant over years, MEX observations were supposed to be offset. However, due to MEX' specific orbit, most MEX observations were performed in similar geometrical condition. As a consequence, a constant bias on the predicted position of Deimos may provide similar astrometric residuals.
A closer look at the MRO data showed that the errors of a few kilometers on Deimos observations were related to a difference lower than 0.2 pixel, only. Such a small difference is related to the scale of the image since Deimos was far away when observed for astrometric purposes by the MRO spacecraft. In particular, the Martian moons are known to have a complex shape, and the fit of their center of figure (assumed to be approximately equal to their center of mass) may be inaccurate when not properly resolved. Last but not least, using independent measurements, Ziese and Willner (2018) showed that our former solution NOE-4-2015-b (which did not use MRO data) was in much better agreement for Deimos than mar097. Hence, we privileged here again the MEX data and removed the MRO data during the fitting procedure.

Comparison of ephemerides and extrapolation precision
We compared our ephemerides with our former solution NOE-4-2015-b as well as the ephemerides mar097 developed at JPL (Jacobson 2010). In particular, NOE-4-2015-b solution used a pretty similar dynamical model and observation sets to the present work. However, it do not benefit from the new MEX observations of Ziese and Willner (2018), as well as the unpublished ones used here. Figure 4 shows differences for both moons in distance in 3-dimensional space. In the case of Phobos the differences at the time close to 2010 are particularly small (at the level of few kilometers) as a consequence of MEX data that have already been available in 2010. Then, differences increase with time due to small modeling/weight differences during the fitting procedure. The new MEX data used in this study shall foster our confidence in the current ephemerides. On the other hand, in the case of Deimos differences increase linearly over time, as a consequence of the different treatment of the MRO data. As already said, such difference was already pointed out by Ziese and Willner (2018). Figure 5 shows the formal uncertainty of our new ephemeris of Phobos and Deimos. Phobos' ephemeris uncertainty is typically few hundred of meters for the present time. Due to its larger distance to MEX, Deimos ephemeris is less well constrained, but expected to be precise within several hundreds of meters. The real accuracy of our new ephemerides will have to be confirmed by confrontation with independent observational means.

Conclusions
We developed new ephemerides of Phobos and Deimos, fitted to the largest time span of observations, including recent unpublished astrometric data obtained by evaluating images taken by the HRSC aboard Mars Express. During the fitting process, we could solve for the physical libration of Phobos, found to be equal to 1.09 ± 0.01 degrees. Considering the very low uncertainty, this may suggest to consider higher order harmonics with an improved rotation model in the future. We confirmed an inconsistency between recent Mars Reconnaissance Orbiter and Mars Express data on Deimos observations. Since Mars Reconnaissance Orbiter's observations were performed from a large distance we imputed the inconsistency to an erroneous center of figure determination. Our ephemerides NOE-4-2020 are available in SPICE format directly from ftp://ftp.imcce.fr/pub/ephem/satel/NOE/MARS/2020/.