Use of radio science and navigation data for the construction of planetary ephemerides

Keywords: Radio science, navigation, deep space tracking, planetary ephemerides, Saturn, Mercury, P9, Kuiper belt, gravity, Cassini, BepiColombo

The central theme of the thesis concerns the exploitation of radio tracking measurements for the development of planetary ephemerides, in particular, applied on two research topics: 1) the analysis of navigation data of the Cassini mission to enhance the ephemeris of Saturn and to increase our knowledge of the outer Solar system; 2) the simulation of BepiColombo measurements collected during the orbital phase at Mercury, for assessing their potential contribution on the INPOP planetary ephemerides. The first research aims at reprocessing Cassini radio tracking data by exploiting the current knowledge of the Saturnian system developed throughout the mission, i.e. the availability of accurate satellite ephemerides and precise gravity solutions for Saturn, Titan and the other major moons. This allows the production of more precise normal points, which are able to constrain the orbit of the planet at meterslevel for 13 years (almost half of its revolution) and to provide invaluable insights on the mass of the Kuiper belt. The results show a reduction of a factor four on normal points uncertainties with respect to previous analyses, providing tighter constraints on the acceptance regions of the Planet 9 as well.

The second research topic focuses on the end-to-end simulations of BepiColombo MORE experiment, aimed at assessing the achievable accuracy on the MPO orbit reconstruction. A thorough study of MORE capabilities, through a realistic covariance analysis, is in fact fundamental to evaluate BepiColombo contribution on planetary ephemerides, since the uncertainties of Mercury's normal points are deduced from the mapped covariance of the spacecraft state.

v

Résumé

Le thème central de la thèse concerne l'utilisation des données de radio tracking pour le développement d'éphémérides planétaires, en particulier, dans deux cas: 1) l'analyse de données de navigation de la mission NASA Cassini pour améliorer les éphémérides de Saturne et augmenter notre connaissance du système solaire externe; 2) simulation des données radio de la mission ESA BepiColombo, qui seront collectées durant la phase orbital à Mercure, pour évaluer leur contribution sur le développement des éphémérides planétaire de l'Intégrateur Numérique Planétaire de l'Observatoire de Paris (INPOP). Le premier sujet de recherche essaie de traiter les données de navigation de la sonde Cassini autour de Saturne en utilisant la connaissance mise à jour du système Saturnien: éphémérides précises pour les lunes du système et meilleure caractérisation de la gravité de Titan et des autres lunes principales. Ça permis la création des points normaux plus précis, capables de contraindre l'orbite de Saturne pour 13 ans (la moitié de la période de révolution autour du Soleil) au niveau des mètres et de donner précieux informations sur le système solaire externe, en particulier sur la masse de la Kuiper belt et sur la possible position de la neuvième planète, P9. Les nouvelles données montrent une réduction de l'incertitude d'un facteur 5 en respect aux analyses précédentes. La deuxième partie de la thèse se concentre sur la production des simulations réalistes des données radio que le Mercury Orbiter Radio-science Experiment (MORE) de la sonde BepiColombo mesurera durant la phase scientifique de sa mission autour de Mercure. Des points normaux sont après produits avec une incertitude déduite de la matrice de covariance de l'état de la sonde, estimé en utilisant ces données simulées. Ces points sont donc traité par le weighted-least square estimateur d'INPOP pour quantifier l'impact que les données de BepiColombo auront sur le développement des éphémérides planétaire, en particulier pour contraindre l'orbite de Mercure et des paramètres relativistes.

Mots clés.

Radio science, navigation, deep space tracking, éphémérides planétaire, Saturne, Mercure, P9, Kuiper belt, gravité, Cassini, BepiColombo vii

Introduction

Deep space tracking has nowadays reached incredible levels of accuracy, with radiometric systems able to position a spacecraft with respect to a planet with accuracies of few meters. These achievements resulted in a huge improvement for planetary ephemerides construction, opening a whole new field of possibilities for Solar system dynamics and fundamental physics studies. Recent decades have marked a relevant step forward in our understanding of the outer Solar system thanks to the continuous discoveries of new TNOs [START_REF] Jewitt | Discovery of the candidate kuiper belt object 1992 qb1[END_REF][START_REF] Brown | Discovery of a candidate inner oort cloud planetoid[END_REF][START_REF] Trujillo | A sedna-like body with a perihelion of 80 astronomical units[END_REF][START_REF] Becker | Raw 2018 invited talks[END_REF]. An extensive review of the current knowledge of the trans-Neptunian Solar system is given in [START_REF] Prialnik | The Trans-Neptunian Solar System[END_REF]. Although these new discoveries have provided fundamental insights into the complex dynamics of TNOs, new questions regarding the forces sculpting their convoluted orbits arise, in particular for the anomalous clustering of the argument of perihelion observed on a series of long-term stable TNOs. In 2016, [START_REF] Batygin | Evidence for a Distant Giant Planet in the Solar System[END_REF] proposed the presence of a ninth planet, beyond the orbit of Neptune, that is able to explain the observed anomalies [START_REF] Batygin | The planet nine hypothesis[END_REF] offers a thorough overview of the P9 hypothesis). Since then, many attempts to locate the elusive planet have followed [START_REF] Fienga | Constraints on the location of a possible 9th planet derived from the cassini data[END_REF][START_REF] Folkner | Sensitivity of Saturn's orbit to a hypothetical distant planet[END_REF]Holman and Payne, 2016a,b;[START_REF] Fienga | New constraints on the location of p9 obtained with the inpop19a planetary ephemeris[END_REF]. However, as [START_REF] Pitjeva | Mass of the kuiper belt[END_REF] point out, a better knowledge of the masses involved is mandatory for disentangling the potential gravitational signal of P9: in particular, the mass of TNOs located in between the 2:1 and 3:2 mean motion resonances with Neptune, forming the so-called Kuiper belt. Hopefully, modern ephemerides, especially the Saturn orbit inferred from Cassini data, can help us to constrain the cumulative mass of these objects. The spacecraft Cassini completed its mission by plunging into Saturn's atmosphere on 15 September 2017. Nevertheless, there is still much to be done with its incredible legacy. Among the vast amount of scientific data gathered by the spacecraft during almost two decades of mission, the radiometric measurements collected for navigation and radio-science purposes represent a valuable tool to precisely locate Saturn within the Solar system. A good estimate of the orbit of the spacecraft relative to Saturn, by means of ground-referenced measurements, allows us to constrain the position of the planet with respect to Earth with metre-level accuracy. An approach based on the use of normal points has been proven [START_REF] Standish | The observational basis for JPL's DE 200, the planetary ephemerides of the Astronomical Almanac[END_REF][START_REF] Moyer | Formulation for Observed and Computed Values of Deep Space Network Data Types for Navigation[END_REF]. Normal points are derived measurements of the signal propagation round-trip light-time in between a ground station and the spacecraft, which is computed using the estimated trajectory. The process of reconstructing the spacecraft orbit, referred to as OD, relies on radio signals exchanged between a ground station and the spacecraft, aimed at measuring the relative distance and radial velocity. More precisely, range data measure the distance as the delay in time due to the signal propagation, whilst range-rate observables measure the Doppler shift on the signal reference frequency due to the relative motion.

In the dissertation, after a brief introduction to planetary ephemerides, we describe the method we follow to produce new normal points for the positioning of the barycentre of Saturn's system, based on a re-analysis of Cassini navigation data. We then present the updated planetary ephemeris INPOP19a we have produced upon these new measurements, processed with the WLS filter of INPOP. We provide the updates on the outer Solar system brought by INPOP19a, in particular a new estimate of the Kuiper belt mass and new constraints on the location of P9. The thesis examines also the case of BepiColombo, the ESA-JAXA mission to Mercury launched in October 2018. The mission boasts a cutting-edge radio system, the most advanced ever flown on an interplanetary probe, specifically developed for the radio science experiment MORE. The radio tracking data of MORE, powered by the KaT instrument, will enable to study the gravity field of the innermost planet of the Solar system, and to exploit its characteristic environment to test fundamental physics. Through a complete end-to-end simulation of the measurements gathered during the orbital phase around Mercury, we provide a precise quantification of the expected accuracies achievable on the positioning of the spacecraft and how planetary ephemerides will benefit from it. The thesis is structured as follows:

• In Chapter 1, a brief introduction to planetary ephemerides and the mathematical description of the Solar System is presented. A presentation of the INPOP planetary ephemerides is also given, with a focus on the INPOP17a release;

• In Chapter 2, we provide the fundamentals of orbit determination with a description of the most common observations used for deep space missions.

We then present our analysis of Cassini radio tracking data in detail, and we introduce the resulting INPOP19a planetary ephemerides. Lastly, we show the updates on the outer Solar System inferred from INPOP19a, in particular the last estimate of the Kuiper belt mass and the new constraints on the location of P9;

• In Chapter 3, we provide an overview of the ESA-JAXA mission BepiColombo to Mercury and its experiments, followed by a brief report on the current knowledge of the planet. We then present the results of the end-to-end simulations of the orbital phase of MORE, testing a new technique (TDMC) for improving its performances and assessing its support to the other instruments, in particular the laser altimeter BELA;

• In the conclusive Chapter (4), we provide a summary of the results achieved in the thesis work, with a brief discussion on their implications.

Chapter 1

Planetary ephemerides 1.1 Introduction

Since the first Babylonian tablets of the 1st millenium BC, the Tables of Toledo of the XII century AD, to modern ephemerides, the need for describing the motion of celestial bodies has always been a constant presence throughout human history, driven by a deep fascination for the sky above our heads.

The home of the ancient gods or the sphere of fixed star, different representations of the celestial domain coexisted, not always peacefully, down the ages. Amid religious and cultural clashes, the astronomical approach has always played a pivotal role, in its attempt to observe and predict the apparent motions of the Moon, the Sun and the planets, and to understand and describe with scientific method the mechanisms that rule the universe.

Although the first testimonies of an heliocentric representation for the Solar System date back to the III century BC, with the studies of the Hellenistic astronomer Aristarchus of Samos, Ancient Roman astronomer Ptolemy's treatise Almagest, supported by Aristotelianism, canonised the geocentric theory and established it as the dominant vision in the medieval European schools. It is only in the XVI century, with the paradigm shift of the Copernican revolution, that the heliocentric model became the reference representation.

In this extremely slow, awareness process, the breakthrough started with the publication of Copernicus' De revolutionibus orbium coelestium (1543). Copernicus' work triggered a flourishing of astronomical studies, including Tycho Brahe's revolutionary observation techniques (1576-1601), Kepler's laws of planetary motion (1609-1619) and Galileo's pamphlet Sidereus Nuncius (1610), culminated with the publication of the Philosophiae Naturalis Principia Mathematica (1687) and the formulation of the universal law of gravitation by Isaac Newton. At that point, for the first time, we were able to understand and predict the motion of celestial bodies with 1. Planetary ephemerides mathematical rigorousness, finding a long waited deterministic order in the universe.

Numerical planetary ephemerides

Almost three centuries after the Newton's theory of gravitation, the launch of the first interplanetary missions in the 1960's represented a turning point in our understanding of the Solar System and our capabilities to construct planetary ephemerides. The introduction of the radiometric measurements from deep space probes, lowered the uncertainties by several order of magnitudes with respect to previous datasets, which were limited to ground-based angular observations. This observational revolution forced a transformation on our approach to planetary ephemerides. The previously accepted analytical methods failed to reproduce the enhanced accuracies of the new measurements, and were replaced by numerical integrators. Since 1960's, JPL has been publishing and constantly updating planetary ephemerides for the Solar System. The fourth iteration, DE96 [START_REF] Standish | JPL Development Ephemeris Number 96[END_REF], was published in 1975 and represented a significant improvement over the previous version, DE69 [START_REF] O'handley | JPL Development Ephemeris Number 69[END_REF], which dated back to 1969, with better processing of the data and refined equations of motions. In 1977, JPL published the first numerical planetary ephemerides, DE102 [START_REF] Newhall | De 102: a numerically integrated ephemeris of the moon and planets spanning forty-four centuries[END_REF], covering the entire available observational set, with an integration time spanning from 1411 BC to 3002 AD.

The numerical ephemeris representation is given in terms of interpolating Chebyshev polynomials. They have been chosen for their stability during evaluation, while providing a near-minimax interpolation. An extensive presentation is provided in [START_REF] Rivlin | The Chebyshev Polynomials[END_REF].

Mathematical description

The bodies integrated in DE102 included the eight planets, plus Pluto and the Moon, with a model for the Lunar librations. The gravitational interactions between the planets were computed by considering their point-mass accelerations, in the isotropic PPN n-body metric [START_REF] Bertotti | Experimental Gravitation[END_REF].

Point-mass Newtonian and relativistic acceleration

With these assumptions, the point-mass Newtonian gravity plus the point-mass relativistic perturbative acceleration for body i, from each other body j, is given by [START_REF] Moyer | Mathematical formulation of the Double-Precision Orbit Determination Program (DPODP)[END_REF]:

ri = j =i µ j (r j -r i ) r 3 ij    1 - 2(β + γ) c 2 k =i µ k r ik - 2β -1 c 2 k =j µ k r jk + γ v i c 2 + (1 + γ) v j c 2 - 2(1 + γ) c 2 ṙi • ṙj - 3 2c 2 (r i -r j ) • ṙj r ij 2 + 1 2c 2 (r j -r i ) • rj    + 1 c 2 j =i µ j r 3 ij (r i -r j ) • [2(1 + γ) ṙi -(1 + 2γ) ṙj ] ( ṙi -ṙj ) + 3 + 4γ 2c 2 j =i µ j rj r ij + 5 m=1 µ a (r a -r i ) r 3 ia
(1.1) where r i , ṙi and ri are the Solar System barycentric position, velocity and acceleration vectors of body i; r ij = |r j -r i |; v i = | ṙi |; µ j = GM j , with G being the universal gravitational constant and M j the mass of body j; β is the PPN parameter measuring the non-linearity in the superposition of gravity; γ is the PPN parameter that measures the space curvature by unit rest mass; and c is the speed of light. The terms with index a refer to the gravitational acceleration from the asteroids: Ceres, Pallas, Vesta, Iris and Bamberga.

Acceleration due to extended body gravity

In addition, for sake of the long-term accuracy of the integrated Lunar orbit, the authors considered the Earth's and Moon's extended gravity field. Extended gravity, which accounts for the deviations of a body from spherical symmetry, can be expressed using the spherical harmonic expansion theory (see [START_REF] Kaula | Theory of Satellite Geodesy -Applications of Satellites to Geodesy[END_REF] for a concise but exhaustive treatment). Then, for a generic body, the gravitational potential can be written as

U (r, λ, φ) = - GM r ∞ l=0 R ref r l l m=0 P lm (sin φ) C lm cos mλ + S lm sin mλ (1.2)
where r is the radial distance, φ the latitude and λ the longitude (see Fig. 1.2); C lm and S lm represent the unnormalised coefficients of degree l and order m at the reference radius R ref , and P lm (sin φ) is the associated Legendre function, which may be defined as follows where P l (x) are the Legendre polynomials of degree l, named after the French mathematician Adrien-Marie Legendre, and defined as

P l0 (x) = 1 2 l l! d l dx l x 2 -1 l (1.3) P lm (x) = 1 -x 2 m 2 d m dx m P l (x) , m = 0, 1, . . . , l (1.4)
P 0 (x) = 1 P 1 (x) = x P 2 (x) = 3 2 x 2 - 1 2 . . . P n+1 (x) = (2n + 1)P n (x) -nP n-1 (x) (n + 1) (1.5)
The coefficients C lm and S lm are defined as

C l0 = 1 M R ref V ρ(r )r P l (sin φ)dV , l ≥ 0 (1.6) C lm = 2 M R ref (l -m)! (l + m)! V ρ(r )r P l (sin φ) cos mλdV , l, m ≥ 0, m ≤ l (1.7) S lm = 2 M R ref (l -m)! (l + m)! V ρ(r )r P l (sin φ) sin mλdV , l, m ≥ 0, m ≤ l (1.8)
It should be mentioned that, because of sin mλ = 0, for m = 0, then S l0 = 0, for l = 0, 1, 2, . . .

(1.9)

The terms of degree l and order m = 0 do not depend on the longitude λ, and are referred as zonal harmonic functions. Their latitude dependence is defined by the Legendre polynomials P l (sin φ), which have exactly l different roots in the interval I = [-1, +1], we can thus distinguish exactly l + 1 latitude zones on the unit sphere, at the borders of which the Legendre polynomial changes sign. Let us define the zonal harmonic coefficient of degree l as

J l = -C l0 (1.10)
In Fig. 1.1 (a) is given an example of zonal harmonic function, with degree l = 6. The terms with l = m, which are called sectorial, have an associated Legendre function P ll (sin φ), which, according to Eq. (1.3), neither depend on latitude nor longitude. They are either multiplied by sin lλ or cos lλ, therefore their sign is constant on latitude but it changes at l equally spaced longitudes (see Fig. The contribution to the inertial acceleration of an extended body arising from the interaction of its own figure with an external point mass can be expressed in a body-fixed frame ξηζ as the one represented in Fig. 1.2; by dividing the zonal and tesseral contributions, the acceleration due to extended body gravity is [START_REF] Moyer | Mathematical formulation of the Double-Precision Orbit Determination Program (DPODP)[END_REF])

    ξ η ζ    = - µ r 2                  lmax l=2 J l R ref r l     (l + 1)P l (sin φ) 0 -cos φP l (sin φ)     zonal harmonics contribution + imax l=1 R ref r l l m=1     -(l + 1)P lm (sin φ)[C lm cos mλ + S lm sin mλ] m sec φP lm (sin φ)[-C lm sin mλ + S lm cos mλ] cos φP lm (sin φ)[C lm cos mλ + S lm sin mλ]     tesseral harmonics contribution                  (1.11)
where r is the separation between the of two bodies centre-of-mass, l max is the maximum degree considered in the expansion, and the primes denote differentiation with respect to the argument sin φ. The accelerations are then transformed into the inertial frame by application of the appropriate rotation matrix.

Tidal accelerations

Because of the accuracy of LLR measurements, a detailed model of the Earth-Moon tidal interaction is required. The tides raised from the Sun and the Moon upon the Earth affect the Moon orbit. The deformation of the Earth is characterised by the Love numbers k 2j , where the orders j = 0, 1, and 2 correspond to tides with longperiod, diurnal, and semi-diurnal periods, respectively. The model accounts for tidal dissipation by introducing a time-delay representation, which includes a response time for the Earth deformation [START_REF] Folkner | The Planetary and Lunar Ephemerides DE430 and DE431[END_REF]. The delay τ varies for longperiod, diurnal, and semi-diurnal responses, therefore different τ j are considered for each order j. The long period zonal tides do not depend on the Earth's rotation, so τ 0 = 0. The acceleration of the Moon due to the distorted Earth depends on the position of the Moon with respect to the Earth r and on the modified position vector r * j for the tide-raising body t (the Moon or the Sun), evaluated at an earlier time t -τ j , which is given for each order j by

r * j = R z (-θτ j )r t (t -τ j ) (1.12)
where R z (-θτ j ) represents a right-handed rotation of the vector r t (t -τ j ) by an angle θτ j about the Earth's rotation axis. Therefore, the resulting acceleration on the Moon can be expressed in cylindrical coordinates, with the z axis perpendicular to the Earth's equator, so that r = ρ + z and the time-delayed position of the tide-raising body is r * j = ρ * j + z * j . For each tide-raising body t, the expression is provided by [START_REF] Folkner | The Planetary and Lunar Ephemerides DE430 and DE431[END_REF] as

r = 3 2 M + M M µ t R 5 r 5 k 20 r * 5 0 (2z * 2 0 z + ρ * 2 0 ρ) - 5((zz * 0 ) 2 + 1 2 (ρρ * 0 ) 2 )r r 2 + r * 2 0 r + k 21 r * 5 1 2 (ρ • ρ * 1 )z * 1 + zz * 1 ρ * 1 - 10zz * 1 (ρ • ρ * 1 )r r 2 + k 22 r * 5 2 2 (ρ • ρ * 2 )ρ * 2 -ρ * 2 ρ - 5((ρ • ρ * 2 ) 2 -1 2 (ρρ * 2 ) 2 )r r 2 (1.13)
where M and M are the Earth's and Moon's mass, respectively, and R is the mean Earth's radius.

INPOP planetary ephemerides

INPOP are the planetary ephemerides developed at the IMCCE -Observatoire de Paris. They provide a numerical representation of the orbits of the eight planets of the Solar System, Pluto and the Moon, plus the orbits of 300 asteroids of the main belt. The ephemerides are constructed by fitting more than 55 000 observations from interplanetary missions and ground-based measurements.

The project INPOP started in 2003 from the evolution of analytical ephemerides developed at IMCCE in the early 1980's, from the growing needs of higher accuracy for reproducing the observational data, and for long-term astronomical and paleoclimate analyses of the Earth and Mars. The first public release, INPOP06 [START_REF] Fienga | Inpop06: a new numerical planetary ephemeris[END_REF], is based on a test solution, INPOP05, which was an attempt to replicate DE405 [START_REF] Standish | JPL planetary and lunar ephemerides[END_REF]. The two ephemerides are in fact very close (see Table 1.1) and, although there are few differences, in particular on the SSB computation and the asteroids representation, they mostly share the same dynamical model.

In JPL DE ephemerides the Sun is not integrated, but derived by assuming that SSB remains at the origin of the inertial reference frame throughout the integration i µ * i r i = 0 (1.14) with µ * i , considering terms up to order 1/c 2 , being

µ * i = µ i   1 + 1 2c 2 v 2 i - 1 2c 2 j =i µ j r ij   (1.15)
By derivation of Eq. (1.14) we obtain

μ * i = µ i 2c 2   j =i µ j (r j -r i ) • ( ṙj -ṙi ) r 3 i   (1.16)
The SSB at the origin of the axes is assured by solving the system of differential equations

   i µ * i r i = 0 i (µ * i ṙi + μ * i r i ) = 0
(1.17)

Differently from DE, in INPOP, Eq. (1.17) are only solved at the initial of the integration, at J2000; once the frame is centred on the SSB, the equations of motion of planets and Sun are integrated in this fixed reference frame. By integrating the We consider here a simplified model, composed of the Sun, the planets from Mercury to Pluto and the Moon (no asteroids), all considered as point-mass bodies. At the time origin of integration (J2000), the origin O of the reference frame is at the barycenter G of the system (Eqs. [START_REF] Moyer | Mathematical formulation of the Double-Precision Orbit Determination Program (DPODP)[END_REF]4,5)). The curves show the drift of the barycenter G in the reference frame, that is

OG = µ * i -1 i µ * i r i .
Time interval is in years from J2000 and the coordinates of G are in mm.

with, by derivation of (2) and up to order 1/c 2 ,

μ * i = µ i 2c 2         j i µ j r j -r i • ṙ j + ṙi r 3 i j         • (4)
In INPOP, we have taken the term μ * i r i into account. Moreover, the Sun is treated as the other planets without assuming a fixed SSB. Indeed, because of the approximations that are still made in the computation of µ * i and μ * i , there remains a small drift of the SSB in the fixed reference frame centered on the SSB at J2000. The determination of the SSB at the origin of time (J2000) is obtained by solving the equations

i µ * i r i = 0 i µ * i ṙi + μ * i r i = 0 (5)
where µ * i and μ * i are given by Eqs. ( 2) and (4). Because µ * i and μ * i depend on the planet velocities, an iterative process is needed. Contrarily to the JPL method, these equations are solved only at the initial step of the planetary integration, at J2000. Once the frame is centered on the SSB defined by Eqs. (5) at J2000, the equations of motion of planets and Sun are integrated in this fixed reference frame. Because of the approximations in 1/c 2 , the positions and velocities of the SSB deduced at t still has a very small displacement that can be neglected (Fig. 3).

If, as in the JPL model (Standish 2004), the μ * i term is neglected in the second equation of ( 5), a more important drift appears in the SSB motion (Fig. 4). Although this does not have any impact on the precision of the solutions, we have preferred Sun along with the other planets, a small drift of SSB can be detected (see Fig. 1.3), but, as long as the μ * i term is included, it remains negligible. In DE102 Earth-Moon tidal model, only tides raised from the Moon upon the Earth were considered [START_REF] Newhall | De 102: a numerically integrated ephemeris of the moon and planets spanning forty-four centuries[END_REF]. In INPOP, as in DE releases since DE403, the contribution from the Sun is included as well, according to Eq. (1.13). Also the asteroids representation in INPOP slightly differs from DE one. In DE405, the gravitational perturbation from 300 asteroids is considered, with the contribution from Ceres, Pallas, Vesta affecting all the planets and that from the other 297 objects interesting only Mars and the Earth-Moon system. In INPOP06, the same 300 bodies are considered, however, their gravitational acceleration acts on all the planets. The mass of the asteroids Ceres, Pallas, Vesta, Iris and Bamberga, together with the densities of the three taxonomic classes C, S and M , are fitted to observations in INPOP06, whilst in DE405 only that of Ceres, Pallas and Vesta. The differences in the modelling of the asteroids perturbation is the principal cause for the larger divergence between INPOP06 and DE405 for the inner planets [START_REF] Fienga | Analytical and numerical studies of asteroid perturbations on solar system planet dynamics[END_REF]. Moreover, the gravitational potential of an asteroid ring is included in INPOP06 for the integration of the planets. The ring is modelled as circular, with radius r ring , mass M ring , and centred in SSB. The induced gravitational acceleration for body i 1. Planetary ephemerides Table 1.2. Maximum differences between INPOP17a, INPOP13c and DE436 in the time interval from 1980 to 2020. The differences are reported in RA, Dec and range, in a geocentric frame [START_REF] Viswanathan | Inpop17a planetary ephemerides[END_REF].

Planet or body

Geocentric differences (1980 -2020): where r * (θ) is the position vector of a point of the ring, and (u, v) represents an orthogonal basis in the ring plane. An asteroids ring perturbation on the planets has been later included in DE414 [START_REF] Konopliv | A global solution for the mars static and seasonal gravity, mars orientation, phobos and deimos masses, and mars ephemeris[END_REF] and since EPM2000 [START_REF] Pitjeva | Modern numerical ephemerides of planets and the importance of ranging observations for their creation[END_REF].

INPOP17a-INPOP13c INPOP17a-
After INPOP06, a series of releases has been published over the years, with incremental updates and adjustments to the dynamical model and observational set: INPOP08, including a new filtering process with a priori information [START_REF] Fienga | INPOP08, a 4-D planetary ephemeris: from asteroid and time-scale computations to ESA Mars Express and Venus Express contributions[END_REF]; INPOP10a, featuring a new model for the asteroids [START_REF] Fienga | The inpop10a planetary ephemeris and its applications in fundamental physics[END_REF]; INPOP10e, introducing the new Gaia DPAC observations release [START_REF] Fienga | GAIA DPAC INPOP final release: INPOP10e[END_REF]; INPOP13c, including the new NASA MESSENGER range data for Mercury [START_REF] Fienga | Numerical estimation of the sensitivity of INPOP planetary ephemerides to general relativity parameters[END_REF]; INPOP15a, with new normal points for Saturn from NASA Cassini mission [START_REF] Fienga | Tests of GR with IN-POP15a planetary ephemerides: Estimations of possible supplementary advances of perihelia for Mercury and Saturn[END_REF]; and, at last, INPOP17a [START_REF] Viswanathan | Inpop17a planetary ephemerides[END_REF], whose novelties will be analysed in detail in the next section.

INPOP17a

INPOP17a is the INPOP planetary ephemerides release published in 2016 [START_REF] Viswanathan | Inpop17a planetary ephemerides[END_REF]. It introduces profound innovations over the previous versions IN-POP13c and INPOP15a. In particular, the inclusion of the new measurements for Mars from ESA MEX and NASA MO missions, and for Saturn from the Doppler-only analysis of Cassini radio tracking data from [START_REF] Hees | Constraints on modified Newtonian dynamics theories from radio tracking data of the Cassini spacecraft[END_REF], represents a fundamental addition to INPOP dataset, providing increase consistency between the DE and INPOP ephemerides, especially for Mars (see Table 1.2). A list of the observables set included in the fit of INPOP17a is presented in Table 1.3, along with the statistics of the post-fit residuals (mean and SD) and a comparison with the values of INPOP13c. The mass of 28 additional asteroids are estimated in INPOP17a with respect to INPOP13c, for a total of 168 objects. In both the releases, no asteroids ring is included in the integration model. The results are shown in Fig. 1.4, with a comparison with the values obtained from DE430. Two categories of asteroids are identified, the great perturbers, that are able to produce a perturbation on the orbit of Mars larger than 7 m over 40 years in geocentric distance, and the small perturbers, that induce smaller perturbations. For the great perturbers, the estimated masses are consistent with those from [START_REF] Carry | Density of asteroids[END_REF], computed with different techniques, while for small perturbers, an over-estimations of the masses can be noted in INPOP17a solution, compared to [START_REF] Carry | Density of asteroids[END_REF] and DE430 results. The discrepancy can be explained by the number of asteroid masses adjusted in INPOP17a (168) with respect to that of DE430 (343). INPOP17a introduces also significant improvements over Jupiter and Mars orbits, which are now closer to DE430, thanks to a better modelling of the effects of the solar conjunctions [START_REF] Viswanathan | Inpop17a planetary ephemerides[END_REF]. INPOP19a: updates on the outer Solar System from Cassini

The orbit determination problem

Orbit determination is the science that exploits tracking data, afflicted by measurement errors, in order to determine the best1 estimate of the state of a spacecraft, whose initial state is unknown, using an approximate mathematical model [START_REF] Tapley | Statistical Orbit Determination[END_REF]. Precisely, it is an iterative process that compares a batch of observables to their respective computed values by means of a proper dynamic model (Fig. 2.1). The residuals of this comparison are thus minimised through a performance index in order to get the optimal solution. It is then necessary to iterate the procedure to reach the convergence of the solution, due to the nonlinear nature of the problem.

The classical example of orbit determination problem considers the bodies moving under the influence of a central force; in this simplified case, the minimal set of parameters composing the dynamic system state are the position and velocity vectors of the spacecraft. Let this general state vector at a time t be denoted by X(t). The knowledge of this quantity at some initial time t 0 , X 0 , would allow the trajectory to be known at any time by propagating the equations of motion, with the intrinsic limitations of the mathematical model in use. Generally X 0 is undetermined or only known with a certain approximation, therefore, the aim of the batch estimator is to process a batch of observations in order to estimate the state vector at the initial time. The linearisation procedure

As previously mentioned, the problem is described by nonlinear expressions; both the dynamic equations and the observations dependencies from the state, are generally nonlinear. One may therefore express these relationship as

Ẋ = F (X, t) (2.1) Y i = G(X i , t i ) + ε i , with i = 1, . . . , l (2.2)
where X is the n-dimensional state vector and Y is a o-dimensional set of l observables, where o < n and m = o × l n; ε denotes the systematic errors afflicting the observables. If a reference trajectory X * is available, and thus a related set of computed observables Y * , then the state and the observations vectors can be written as

x(t) = X(t) -X * (t) (2.3) ẋ(t) = Ẋ(t) - Ẋ * (t) (2.4) y(t) = Y (t) -Y * (t) (2.5)
with x(t) and y(t) being the state and observations deviations from a reference trajectory. Therefore, while the reference trajectory stays reasonably close to the real one, the deviation vectors can be expanded in Taylor's series as follows

Ẋ(t) = F (X * , t) + ∂F (t) ∂X(t) * [X(t) -X * (t)] + O F [X(t) -X * (t)] (2.6) Y i = G(X * i , t i ) + ∂G(t) ∂X(t) * [X(t) -X * (t)] + O G [X(t) -X * (t)] + ε i (2.7)
where [ ] * indicates that the partial derivative matrix is evaluated with respect to X * . Neglecting higher order terms (O F ,O G ) and using Eq. (2.3-2.5), Eq. (2.6,2.7) can be written as

ẋ(t) = A(t)x(t) (2.8
) The state transition matrix has the following properties:

y i = Hi x i + ε i (2.
1. Φ(t i , t i ) = I , 2. Φ(t i , t j ) = Φ(t i , t k )Φ(t k , t j ) , 3. Φ(t i , t j ) = Φ -1 (t j , t i ) .
Deriving Eq. ( 2.11) and substituting it in Eq. ( 2.8) we obtain

Φ(t, t i ) = A(t)Φ(t, t i ) (2.
12)

The state transition matrix can now be used to refer any state and any observation vectors to the initial ones:

x i = Φ(t i , t 0 )x 0 (2.13) y i = H i x 0 + ε i (2.14) with H i = Hi Φ(t i , t 0 ) (2.15)
If the dynamic model used to integrate the trajectory and to compute the observables perfectly reflects reality and the measurement noise is white and unbiased, then the residuals ε show a purely Gaussian distribution. In fact, recalling the central limit theorem, the sum of a large set of independent random variables tends to a normal distribution, even if the original variables themselves are not normally distributed.

However, approximations and errors in the model introduce signatures in the residuals; it is thus necessary to adjust the model parameters via a procedure known as WLS estimation, in which the optimal solution is defined to be the set of parameter values that minimises the weighted sum of squares of the residuals vector.

The weighted least squares estimator

The least squares solution, first proposed by Gauss in 1809, selects the estimate of x as that value that minimises the sum of squares of the calculated observation residuals, in particular it is selected to minimise the following performance index .16) with ε = (y -Hx) (2.17) therefore the minimum of J and thus the optimal solution x, are obtained for iterative procedure [START_REF] Thornton | Radiometric tracking techniques for deep-space navigation[END_REF].

J(x) = 1 2 ε T ε , ( 2
∂J ∂x = 0 ---→ x = (H T H) -1 H T y (2.18)
H T H is known as the normal matrix, and the solution only exists if it is nonsingular, i.e. if det(H T H) = 0.

The weighted solution introduces a weighting matrix W that accounts for the observables accuracy. Indeed various type of measurements are used in the OD process, and, moreover, those are often collected by different ground station-stations, and are thus characterised by noises with different SDs. It is therefore useful to assign different weights, namely reliability, to those measurements.

In absence of correlation within observations at different times W is simply a diagonal matrix

W =        w 1 0 • • • 0 0 w 2 • • • 0 . . . . . . . . . . . . 0 0 • • • w l       
where generally w i is proportional to the SD of the i-th observable.

The solution becomes the estimate of the x that minimises the weighted sum of squares of residuals

J w (x) = 1 2 ε T W ε, (2.19) hence for ∂J w ∂x = 0 --→ x = H T W H -1 H T W y (2.20)
It possible to prove that if the noise is white and with a Gaussian distribution, and the weights assigned are the inverse of the observations variances, then

P = H T W H -1 (2.21)
is the covariance matrix associated to the estimate. Therefore the diagonal elements represent the variances of the state vector estimate and the off-diagonal terms the covariances. This result is of paramount importance, since it enables the access to the uncertainties related to the estimated values, indeed, an estimate without its related uncertainty is quite useless, since it does not give any information about the reliability of the result.

In the definition of P , it is possible to identify the H T W H matrix, that in Eq. ( 2.19) has taken the place of the normal matrix, H T H, in Eq. (2.18).

A priori information

Often the set of observations comes along with some a priori information that helps the research of the true state. In general this is information gathered by previous tracking campaigns, or ephemeris provided by thirds, and it is always accompanied by its related uncertainty. Indeed a weight to the a priori information has to be assigned as well, based on its uncertainty. The estimation process now provides for the minimisation of an augmented cost function

J ap (x) = 1 2 ε T W ε weighted solution + 1 2 (x -x) T W (x -x) a priori information , (2.22) then for ∂J ap ∂x = 0 --→ x = H T W H + W -1 H T W y + W x (2.23)
where x is the a priori, unbiased value of x and W is the a priori weighting matrix.

The normal matrix has now took the form, H T W H + W . One must pay particular attention to the assignment of these weights since, they act as constraints on the solution, therefore a misplaced trust, i.e. high weight on a wrong a priori information, can prevent the estimator from converging to the true state, or from converging at all. However, the benefits of a priori information are indisputable; in fact, in addition to speed up the convergence process, they often resolve potential singularity problems as well, related to the nonlinearity of OD, that would result in multiple solutions.

The multi-arc technique

The aforementioned classic solutions may fail when the dynamic to estimate is particularly complex or very long periods of tracking are available, forcing the filter to process a great number of observables in a single batch. In presence of large set of data, the estimation process could get complicated, in particular, the normal matrix, H T W H+ W , could become ill-conditioned, making its inversion impractical. Moreover, if the integration time span is too extended, errors and inaccuracies in the dynamic model may lead to the divergence of the integrated trajectory from the true one.

A possible solution to these problems is represented by the "multi-arc" approach [START_REF] Milani | Multi-arc strategy[END_REF]. This method consists in subdivide the trajectory in shorter arcs, in order to limit the number of observables to process in a single batch, and to contain integration divergencies at a level compatible with the observations noise. The solve-for state vector become then an augmented state, divided in two categories: local parameters, consisting of those quantities pertinent to the single arc, such as the spacecraft state vector, periodic accelerations, manoeuvres etc.; and global parameters, comprising those constant values common to every arc, such as planet's GM, spherical harmonics, tidal coefficients, etc. Hence, it can be expressed as

x = l g (2.24)
with l and g being respectively the local and global parameters vectors. The local term is then further divided into vectors l i , with i referring to the i-th arc. Adopting this formulation, the mapping matrix acquires the following structure

H =        ∂y 1 ∂l 1 0 • • • 0 ∂y 1 ∂g 0 ∂y 2 ∂l 2 • • • 0 ∂y 2 ∂g . . . . . . . . . . . . . . . 0 0 • • • ∂y n ∂ln ∂y n ∂g        (2.25)
The solution can then be computed with Eq. ( 2.22), as previously seen with the WLS with a priori solution.

However, one should note that the arc length choice is of paramount importance for the correct operation of this technique. Indeed, if the arc is too short, the solution may become unstable, but, if it is too long, the introduced advantage may be negligible.

Propagation of the covariance matrix

If a state estimate xj is available at time t j , obtained by using Eq. ( 2.22) as solution of an OD process, we may propagate it to any later time by using the state transition matrix and Eq. ( 2.11):

xk = Φ(t k , t j )x j (2.26)
Equivalently, we may propagate its associated covariance matrix, as follows [START_REF] Tapley | Statistical Orbit Determination[END_REF])

Pk = cov(x k ) ≡ E (x k -x k ) (x k -x k ) T (2.27)
then, by substituting Eq. ( 2.26) in Eq. ( 2.27), we obtain

Pk = E Φ(t k , t j ) (x j -x j ) (x j -x j ) T Φ T (t k , t j ) (2.28)
Since the state transition matrix is deterministic, we may write

Pk = Φ(t k , t j )P j Φ T (t k , t j ) (2.29)
Equations ( 2.26) and ( 2.29) can be used to map the estimate of the state and its associated covariance matrix to any time.

Observations

The orbit determination problem is based on a set of observations gathered by means of radio signals sent between a spacecraft and a ground station, the most common nowadays comprising range and range-rate measurements. There are multiple configurations that can be adopted in order to collect those measurements:

• One-way link: the signal carrier is directly generated by the onboard ultra stable oscillator and transmitted to the ground station.

• Two-way link: in this configuration the signal is generated by the ground station and transmitted to the spacecraft, where the transponder retransmits it, at a coherent turn-around ratio frequency, to the same ground station (Fig.

2.2a).

• Three-way link: it shares the same architecture of the two-way link, but in this case a separate listen-only antenna receives the downlink signal (Fig. 2.2b).

Communications are based on internationally defined frequency bands (Table 2.1), in order to avoid potential interferences from external sources. Deep space missions mostly rely on X-band links for navigation and radio science; however, the use of Ka-band for the Cassini and Juno missions, and the first results [START_REF] Cappuccio | Report on first inflight data of bepicolombo's mercury orbiter radio-science experiment[END_REF] from BepiColombo data, have proved the great advantages brought by this higher frequency for precision OD and radio science purposes.

Range

In OD range is referred to as the distance along the line-of-sight between the ground station and the spacecraft. It is measured as the delay in time due to the signal propagation at the speed of light c,

ρ = c(t r -t t ) (2.30)
where t r and t t are respectively the receiving and transmitting time.

Since two different clocks are used to register these times, the lack of synchronization between them produces an error in the range measurement. For this reason ρ is usually referred to as pseudorange. This problem can be overtaken by using a two-way link, in which the same clock measures both transmitting and receiving time after a so-called RTLT. Therefore, neglecting the time-delay introduced by the spacecraft signal processing, the measured quantity is the round-trip (or two-way) distance

ρ RT = c(t r -t t ) (2.31)

Sequential ranging

At present, the standard ranging system routinely used for navigation is the Xband sequential ranging. The sequential ranging signal consists of a sequence of sinusoidal tones derived from the station frequency standard, modulated onto the transmitter carrier signal. The spacecraft receiver locks on the uplink signal via a PLL that produces a coherent reference signal used to demodulate the ranging signal. The latter is then modulated onto the downlink carrier coherently with the uplink one but with an offset in frequency, introduced to avoid auto-lock, namely the PLL confusing the downlink with the uplink. At the receiving station, another PLL produces a coherent reference signal from the received one and uses it to demodulate the downlink signal. Range measurements are quantised in steps referred to as RU, whose size depends on the the frequency of the highest component of the code [START_REF] Thornton | Radiometric tracking techniques for deep-space navigation[END_REF]. One RU corresponds to approximately 0.94 ns in round-trip range for a 1 MHz ranging tone, which translates to about 14 cm in one-way distance.

PN ranging

The growing need for more accurate measurements for the challenging navigation and radio science requirements of the next generation of deep-space missions has compelled ESA, NASA and other space agencies, to research new ranging techniques. Classic sequential ranging employs a transparent ranging system, i.e. it uses a non-regenerative (or turnaround) ranging technique where the spacecraft translates the uplink ranging signal to the downlink without code acquisition. In this way the uplink noise is also modulated onto the downlink carrier, incurring a path loss of 1/r 4 , degrading the ranging measurement precision. A solution to this problem is offered by the new PN ranging system. It provides regenerative ranging capabilities, i.e. the on-board transponder demodulates and acquires the ranging code by correlation with a local code replica from the uplink ranging signal, and regenerates the ranging code on the downlink. Regenerative ranging has the potential to increase the ranging SNR by as much as 30 dB in deep space applications.

Range biases

Range measurements are afflicted by potential biases from path delays introduced by the propagation media (e.g. troposphere, ionosphere, and solar plasma), instrumentation at the ground station, and the spacecraft radio system. So far, for sequential ranging, the ground station provides a characterisation of the station delay with metres-level uncertainty before (pre-cal) or after (post-cal) each pass [START_REF] Border | Station Delay Calibration for Ranging Measurements[END_REF]. During the pass, variations of about 10 cm are expected, however they are below the RMS level of the observables, and are typically neglected in the data processing. A common approach for absorbing the station calibration errors is to estimate stochastic range biases, one per pass, in order to avoid the introduction of systematics in the OD analysis.

Range-rate

Range-rate are the time-derivative of a range observables' counterpart, related by the equation

y(t) = 1 c dρ RT dt = d(RTLT) dt (2.32)
where y(t) is called relative frequency shift and is typically measured as the Doppler effect on the signal due to the relative motion between the spacecraft and the ground station.

In first analysis, a signal transmitted from the spacecraft in a one-way link can be written as

s t (t) = {s 0t exp [i (ω 0 t + φ 0 )]} (2.33)
then the received signal at the ground station would be

s r (t) = s 0r exp i ω 0 t + φ 0 + k r G/S (t) -r S/C (t -τ ) (2.34)
where k is the wave-vector, r G/S and r S/C are the position vectors of the ground station and the spacecraft, and τ is the one-way light time between the spacecraft and the station. The term k(r G/S (t) -r S/C (t -τ )) is the equivalent of the number of wavelengths contained in the G/S -S/C distance, multiplied by 2π.

The frequency can now be written as the derivative of the signal phase

f r = 1 2π dΦ r (t) dt = f 0 + k 2π d dt r G/S (t) -r S/C (t -τ ) = = f 0 + f 0 c dρ dt = f 0 1 + 1 c dρ dt (2.35)
hence, it is possible to recognise the frequency shift as

f r -f 0 = f 0 1 + ρ c , (2.36)
where ρ/c is the measured Doppler shift. Therefore, Doppler data are obtained by differencing the received signal with the station frequency reference (Fig. 2.3); this is carried out by using a PLL that replicates a phase-coherent, clean replica of the incoming signal, as happens for the range measurements (Sect. 2.2.1). [START_REF] Thornton | Radiometric tracking techniques for deep-space navigation[END_REF].

Doppler noise sources

The precision of these observables directly affects the level of accuracy achievable in the estimation of the parameters, hence, a precise description of the error budget becomes fundamental in assessing radio science experiments results. The most adopted Doppler noise characterisation in time-domain is the Allan deviation (see [START_REF] Riley | Handbook of Frequency Stability Analysis[END_REF] for an exhaustive definition),

σ(τ ) = 1 2(M -1) M -1 i=1 (y i+1 -y i ) 2
(2.37) where y i is the i-th of M fractional frequency values averaged over the sampling interval τ .

State-of-the-art Doppler systems, such as the one on board the ESA-JAXA BepiColombo mission to Mercury (see Sect. 3.3), rely on multi-frequency links that permit the almost complete cancellation of dispersive noise [START_REF] Bertotti | Doppler tracking of spacecraft with multi-frequency links[END_REF], i.e. solar and ionospheric plasma fluctuations. The noise budget of these links is then dominated by local noises arising in the proximity of the ground station, such as tropospheric scintillation and mechanical noise introduced by the antenna structure [START_REF] Asmar | Spacecraft Doppler tracking: Noise budget and accuracy achievable in precision radio science observations[END_REF]. Deep space antennas present indeed notable sizes, with parabolic dishes ranging from 34 to 70 m in diameter, and are subject to thermal expansions, gravity loadings and wind gusts throughout the tracking passes, introducing relevant disturbances on the measurements; typical values are σ ∼ 1.6×10 -14 or 0.005 mm/s at τ =60 s (Iess et al., 2014a). Tropospheric noise refers to the frequency/phase fluctuations introduced on the signal along its propagation path in the Earth's atmosphere. If its dry component contributing to almost the 90% of the path delay can be effectively calibrated, it is the hardly predictable wet part that afflicts the frequency stability most. This noise is related to the weather conditions and the amount of water-vapour in the air at the ground station, therefore it presents significant seasonal variations, with measured peaks of σ ∼ 10 -13 at τ =60 s on summertime days (Iess et al., 2014a). Even if considerably smaller, instabilities introduced by the ground station FTS represent another important disturbance to the link. Modern Hydrogen-Maser clocks installed at DSN and ESTRACK stations provide σ ∼ 10 -15 at 1 000 s (Iess et al., 2014a), with a dominant white frequency noise behaviour (σ ∝ τ -1/2 ) at lower integration times, 10-1 000 s. Hopefully, their short-term stability can be improved with the adoption of CSOs as the prototype installed at DSS 25, capable of assuring σ ∼ 10 -15 even at τ <100 s [START_REF] Giordano | Tests of Sapphire Crystals Produced with Different Growth Processes for Ultra-stable Microwave Oscillators[END_REF].

Normal Points

Once we dispose of a reconstructed trajectory for the spacecraft, including corrections to the range measurements in the form of range biases or planet's orbit adjustments, we can produce normal points for planetary ephemerides construction. These are derived measurements of the Earth-planet distance, computed combining the information on the position of the spacecraft with respect to the planet itself, derived from OD, with the range observable, properly corrected with the range bias. The normal points are thus virtual measurements, time-tagged in the middle of each tracking pass, providing the signal round-trip light-time between the ground station and the spacecraft, corrected for the Earth's troposphere and ionosphere delays and the relative range bias, properly transformed in time delay. The signal propagation time is computed as

t r -t t = r tr c + ∆ GR , (2.38) with r tr = |r r (t r ) -r t (t t )| (2.39)
where t r and t t are the times at receiver and transmitter, respectively, r r (t r ) and r t (t t ) are the positions of the receiver and transmitter in the inertial frame at the receiving and transmitting time, and ∆ GR is the general relativity correction term [START_REF] Moyer | Formulation for Observed and Computed Values of Deep Space Network Data Types for Navigation[END_REF]. The solution, which requires an iterative process to account for the relative motion of the bodies within the propagation time, is performed using MONTE utilities. Therefore, the normal point is given by ρ = (t r -t t ) uplink + (t r -t t ) dnlink + δ stn + ρ bias (2.40) where (•) uplink and (•) dnlink terms represent the uplink and downlink propagation times, δ stn the possible station delays, and ρ bias the estimated range bias expressed in seconds.

Cassini radio tracking data analysis

The Cassini-Huygens mission

The Cassini-Huygens mission was a joint endeavour of NASA, ESA and ASI that orbited Saturn from July 1 2004, until September 15 2017. Cassini was launched in October 15 1997 from Cape Canaveral, Florida, and it arrived at Saturn after an interplanetary space travel of nearly seven years. During this cruise phase, which included two Venus and one Earth flybys, and a six-months long swing by Jupiter, Cassini registered the first of its sensational measurements. During a solar conjunction, while Cassini was passing behind the solar corona (with SEP≈ 0 • ), a measure of the path-delay on the radio link photons due to the space-time curvature caused by the Sun's mass set a new limit on the testing of general theory of relativity [START_REF] Bertotti | A test of general relativity using radio links with the cassini spacecraft[END_REF]. Since then, in about 13 years touring around the ringed planet and its moons Cassini thoroughly investigated the Saturnian system, registering a countless series of scientific breakthroughs, such as the Huygens probe descent and landing on Saturn's biggest moon, Titan, unveiling its incredible methane cycle [START_REF] Atreya | Titan's methane cycle[END_REF]; the discovery of the Enceladus' hidden ocean and its spectacular southern polar plumes [START_REF] Iess | The gravity field and interior structure of enceladus[END_REF], whose last analyses findings suggest the potential habitability of the icy moon [START_REF] Waite | Cassini finds molecular hydrogen in the enceladus plume: Evidence for hydrothermal processes[END_REF]; the measurements of the mass of the rings, with the implications on their unexpected relative young age [START_REF] Iess | Measurement and implications of saturn's gravity field and ring mass[END_REF]. Cassini, indeed, was a huge orbiting laboratory (see Fig. standard in precision Doppler tracking. In fact in the SCE, Cassini did not directly measure the time-delay of the photons as done in the Viking experiment [START_REF] Reasenberg | Viking relativity experiment -Verification of signal retardation by solar gravity[END_REF], but rather the induced relative frequency shift of the carrier:

∆f f 0 = -d∆t dt ≈ -4 (1 + γ) GM bc 3 db dt (2.41)
where ∆f /f 0 indicates the relative frequency shift, b the impact parameter, and γ is the PPN parameter that expresses how much space curvature g ij is produced by unit rest mass. These Doppler observations were strongly affected by the solar corona plasma scintillation, but thanks to the the effective cancellation of this dispersive noise contribution enabled by the KaT multi-frequency link, [START_REF] Bertotti | A test of general relativity using radio links with the cassini spacecraft[END_REF] were able to measure γ with an unprecedented accuracy of 2.3 × 10 -5 . Unfortunately, right after the experiment, in 2003 an irretrievable failure occurred, making the KaT unusable before Cassini arrival at Saturn. However, this event did not undermine the success of the expedition, whose 4 years primal mission was extended twice, with the Equinox (2008-2010) and Solstice (2010-2017) extensions, plus a four-months final stage, the so-called Grand Finale, which featured 22 deep dives between Saturn's cloud tops and the innermost ring, before the final plunge into the giant planet's atmosphere, on September 15 2017. The 13 years baseline of Cassini radio data permits to track the orbit of Saturn for almost half of its revolution around the Sun, offering a unique insight on its dynamical state and the gravitational accelerations which perturb its orbit. We thus produced new normal points for increasing the accuracy of Saturn's ephemerides, based on a re-analysis of the mission navigation data [START_REF] Di Ruscio | Analysis of cassini radio tracking data for the construction of inpop19a: A new estimate of the kuiper belt mass[END_REF] and on the gravity solutions of Titan flybys [START_REF] Durante | Titan's gravity field and interior structure after Cassini[END_REF] and Grand Finale pericentres [START_REF] Iess | Measurement and implications of saturn's gravity field and ring mass[END_REF].

Analysis of navigation data

During the Cassini tour of the Saturnian system, in order to allow the flight operations team to navigate the spacecraft, tracking from the DSN stations gathered approximately six hours per day of two-way X-band Doppler and range measurements. The orbit of Cassini was dictated by a complex optimisation process aimed at maximising the number of encounters with the moons while limiting the propellant consumption. The actual trajectory included about one flyby of Titan per orbit, with additional flybys of the other major satellites, such as Rhea, Dione, and Enceladus. To this end, the trajectory was made possible by the moon's gravity assists and specific OTMs that were performed to direct Cassini towards the next encounter [START_REF] Brown | A full mission summary (1997-2017) of the attitude control performance during orbit trim maneuvers performed by the cassini-huygens spacecraft[END_REF]. Moreover, the need for precise pointing of the HGA, along with the requirements set by other instruments, demanded fine control of the spacecraft attitude. This was achieved using either the eight 1-N RCS thrusters or the reaction wheel assembly. A slight unbalance of the RCS thrusts introduced undesired mm/s-magnitude ∆V s on the spacecraft, referred to as small forces [START_REF] Lee | Cassini spacecraft attitude control system performance and lessons learned, 1997-2017[END_REF]. In addition, spacecraft dynamics were affected by SRP, the anisotropic acceleration produced by the RTGs, and the drag of the upper atmosphere of Titan during the closer flybys of the moon (with an altitude < 1000 km) [START_REF] Pelletier | Atmospheric drag model for cassini orbit determination during low altitude titan flybys[END_REF]. As result, reconstructing the orbit was an extremely complicated task, making de facto Cassini one of the most complex space missions ever navigated. To give some numbers, Cassini performed a total of 162 targeted moon flybys and 360 successful OTMs in 13 years. The reconstruction process is nowadays a wellestablished procedure [START_REF] Roth | Cassini in-flight navigation adaptations[END_REF], thanks to the comprehensive experience of the spacecraft dynamics developed throughout the whole mission, the enhanced precision of Saturn's satellites ephemerides (Jacobson, 2016a;[START_REF] Boone | Resolution of orbit determination prediction instabilities at titan during cassini's solstice mission[END_REF], and current knowledge of the gravity fields of the major bodies, which have been measured on dedicated flybys [START_REF] Iess | Measurement and implications of saturn's gravity field and ring mass[END_REF][START_REF] Durante | Titan's gravity field and interior structure after Cassini[END_REF][START_REF] Iess | The gravity field and interior structure of enceladus[END_REF][START_REF] Tortora | Rhea gravity field and interior modeling from cassini data analysis[END_REF][START_REF] Zannoni | The gravity field and interior structure of dione[END_REF]. We limited our analysis to three different periods of the 13-year mission: from Febru- ary to July 2006 (including T10-T15 flybys)5 , from November 2008 to November 2009 (T47-T62 plus E8-E9), and from February to April 2011 (T74-T75). We chose these three intervals trying to extend the data coverage, while looking for periods of the mission during which Cassini performed consecutive flybys of Titan without encounters with the other moons. In this way, we are able to achieve a more accurate reconstruction of the orbit. We then cut data characterised by SEP angles < 30 • because of the large solar plasma scintillation (Iess et al., 2014a;[START_REF] Asmar | Spacecraft Doppler tracking: Noise budget and accuracy achievable in precision radio science observations[END_REF], which would significantly degrade and potentially bias the estimate of Cassini trajectory. We analysed the data using an integration time of 60 s for Doppler and 300 s for range measurements, ignoring the observables with elevation < 15 • . The post-fit residuals merged from the three analysed periods are shown in Fig. 2.8, while the residuals of a typical reconstructed arc are depicted in Fig. 2.6, providing, in this case, the best noise levels of the measurements, obtained near a solar opposition (SEP > 140 • ). In this favourable condition Doppler noise amounts to 0.6 mHz (equivalent to ∼0.02 mm/s), and for range data a RMS of 2.6 DSN RU; while for the whole analysed dataset, the average noise is ∼1.3 mHz for range-rate and ∼2.9 RU for range data. In order to reconstruct the orbit, we divided the trajectory in arcs that span two con-secutive flybys. This approach was also followed by the navigation team [START_REF] Bellerose | Cassini navigation: The road to consistent sub-kilometer accuracy satellite encounters[END_REF]. The extension of these arcs is mainly limited by the dynamical model accuracy in predicting the outbound trajectory of the satellite encounter. Although, a recent work from [START_REF] Lari | Chaotic orbit determination in the context of the juice mission[END_REF] shows that, even with a theoretically perfect model, an intrinsic limit in obtaining a unique trajectory exists in presence of multiple flybys, due to the chaotic nature of the problem. The resulting arc subdivision provides for an overlap between two consecutive flybys (e.g. if the i-th arc includes T52 and T53 as in Fig. 2.6, the (i+1)-th arc includes T53 and T54); this allows for a more robust estimate of the orbit and permits, for each flyby, to choose the arc that provides the most accurate solution.

In our OD solution, using navigation reconstructions [START_REF] Antreasian | Cassini orbit determination performance during the first eight orbits of the Saturn satellite tour[END_REF][START_REF] Antreasian | Cassini orbit determination performance during saturn satellite tour: august 2005 -january 2006[END_REF][START_REF] Antreasian | Cassini orbit determination results: January 2006 -end of prime mission[END_REF][START_REF] Pelletier | Cassini orbit determination performance[END_REF][START_REF] Bellerose | Cassini navigation: The road to consistent sub-kilometer accuracy satellite encounters[END_REF] as a priori information, we then solved for the spacecraft initial conditions, corrections to OTMs direction and amplitude, small forces ∆V components, RTG acceleration, a scale factor for SRP, and stochastic accelerations. The latter are used for compensating remaining mismodellings in the spacecraft dynamics; in particular, constant accelerations to the level of 5 × 10 -13 km/s 2 for each spacecraft-fixed frame axis are estimated and updated every eight hours. Moreover, during closer Titan flybys, we account for drag coefficients and corrections to the predicted RCS thrusts, which are used to counter the atmospheric torques and maintain the desired attitude. A representation of the magnitude of the main accelerations acting on the spacecraft during a typical low-altitude flyby of Titan is given in Fig. 2.7. The plot shows the relevance of non-gravitational accelerations in the evolution of the trajectory and how their correct modelling plays a decisive role in obtaining a good reconstruction of the orbit. Because of the limited knowledge of such forces, in some passes additional stochastic accelerations to the level of 5 × 10 -9 km/s 2 along the spacecraft Z-axis6 and 5 × 10 -11 km/s 2 on X and Y were estimated and updated every five minutes. For the gravity fields of Titan and Saturn, which still introduce a major effect on the spacecraft dynamics, we have taken the estimated values from [START_REF] Durante | Titan's gravity field and interior structure after Cassini[END_REF] and [START_REF] Iess | Measurement and implications of saturn's gravity field and ring mass[END_REF] solutions, without producing any further adjustments in our fit. The orbits and point-mass gravity accelerations from the other Solar System planets were computed using INPOP17a planetary ephemerides [START_REF] Viswanathan | Inpop17a planetary ephemerides[END_REF], while gravity interactions with the other satellites of Saturn were derived from JPL SAT389 (Jacobson, 2016a) and SAT393 [START_REF] Jacobson | Satellite ephemeris file release: Sat393[END_REF]. For the ephemerides of Titan we used the solution provided by the navigation team [START_REF] Antreasian | Cassini orbit determination performance during the first eight orbits of the Saturn satellite tour[END_REF][START_REF] Antreasian | Cassini orbit determination performance during saturn satellite tour: august 2005 -january 2006[END_REF][START_REF] Antreasian | Cassini orbit determination results: January 2006 -end of prime mission[END_REF][START_REF] Pelletier | Cassini orbit determination performance[END_REF][START_REF] Bellerose | Cassini navigation: The road to consistent sub-kilometer accuracy satellite encounters[END_REF]. The reconstruction of Titan orbit by the navigation team is sufficiently accurate (Boone and Bellerose, 2017) to obtain a good fit of the data (see Fig. 2.8) except for the arcs of 2006, which were reconstructed in an early phase of the mission when limited knowledge of the Saturnian system afflicted the navigation reconstruction. For those arcs requiring further adjustments, we estimated minor corrections to Titan state, of the order of tens of metres. An assessment of the impact of the uncertainty of Titan ephemeris on orbit reconstruction is provided in [START_REF] Boone | Resolution of orbit determination prediction instabilities at titan during cassini's solstice mission[END_REF]. Although we used range measurements in our analysis, these were largely de-weighted, since they are affected by range biases (see Sect. 2.2.1). Water-vapour radiometer (only for gravity-dedicated passes) and GPS calibrations are available for media delays. However, to compensate for the remaining calibration error, we estimated a common bias on the range observables per tracking pass. We modelled the biases as stochastic parameters, updated every 24 hours, with large a priori uncertainty (500 RU) in absorbing both the calibration residual and the planetary ephemerides mismodelling. Errors in the relative location of Saturn and Earth cause an erroneous positioning of Cassini with respect to the ground station, and thus a miscalculation of the range computed observable. In our fit we did not correct for this term, but we absorbed it with the range biases. The estimated values are reported in Fig. 2.9. A significant signature with an annual frequency and an amplitude of about 50 m is present as a result of the use of INPOP17a. Thanks to the enhanced accuracy of the new points, it is now possible to improve this ephemeris. A fundamental advantage of our analysis with respect to [START_REF] Hees | Constraints on modified Newtonian dynamics theories from radio tracking data of the Cassini spacecraft[END_REF] derives from the choice of processing longer arcs that include moon flybys and OTMs. In this way, it is possible to better constrain the position of Cassini with a continuous orbit and, therefore, to produce more accurate normal points for Saturn ephemeris. The analysis has been performed using JPL's MONTE [START_REF] Evans | Monte: the next generation of mission design and navigation software[END_REF] for both integrating the equations of motion and generating the computed observables with its observation model [START_REF] Moyer | Formulation for Observed and Computed Values of Deep Space Network Data Types for Navigation[END_REF].

Titan gravity flybys

One of the pivotal objectives of the Cassini-Huygens mission was the study of Saturn's biggest moon, Titan. In the 13 years spent touring the Saturnian system, Cassini performed 162 targeted flybys of Saturn's moons, 127 of which of Titan. Among these, nine were dedicated to measuring the gravity field of the moon, plus an additional flyby, T110, primarily devoted to imaging the moon's north polar lakes (see Table 2.2).

The analysis of the ten passes carried out by [START_REF] Durante | Titan's gravity field and interior structure after Cassini[END_REF] provides an improvement on Titan's gravity field determination with respect to the previous analyses of [START_REF] Iess | Gravity field, shape, and moment of inertia of titan[END_REF] and [START_REF] Iess | The tides of titan[END_REF], respectively limited to the first Along with the enhanced solution for Titan's gravity field, Durante et al. ( 2019) obtained an extremely accurate reconstruction of Cassini trajectory arcs around the CA with the moon, thanks to the specific precautions took by the Flight Operations Team in accordance with the Radio Science Team requirements. In fact, for few days before and after the encounters with Titan, no manoeuvres were executed to preserve the dynamical stability of the platform, and to avoid disturbances on the Doppler measurements collected with the HGA constantly pointed towards Earth. The altitude of these flybys ranges between 2 397 and 3 651 km, at such distances drag from the moon's atmosphere (see Fig. 2.7) is negligible, facilitating the gravity reconstruction.

In collaboration with the authors of [START_REF] Durante | Titan's gravity field and interior structure after Cassini[END_REF] we produce new normal points based on their orbit reconstruction. We compute one normal point per tracking pass, for a total of 42 new points for Saturn's ephemerides. The accuracy of these derived measurements is assessed by propagating the complete state covariance matrix of Cassini. The average uncertainty is 7.4 m, with the CA points, constrained by Titan's gravity, exhibiting the highest accuracies.

Grand Finale pericenters

As a result of Cassini's own discoveries on the potential habitability of the oceanharboring Saturn's moons Enceladus and Titan, the mission Planners decided to destroy the spacecraft by plunging it into Saturn's atmosphere, to prevent any risk of contamination. In fact, the high-temperature provoked by the meteoric entry ensured that no contamination from the orbiter could ever infect the interested moons.

Therefore, for the final stage of the mission, the Grand Finale, in April 2017, Cassini was inserted into a series of extremely eccentric, closer orbits, skimming the upper layers of Saturn's atmosphere. The spacecraft flew close to Saturn as it never did before, passing in between the rings and the planet, at distances of 2 600 to 3 900 km from the top clouds. Five (Rev 273, 274, 278, 280 and 284)7 pericentres of these 22 orbits were devoted to measuring the gravity field and rings mass of Saturn.

In these occasions, two-way X-band range-rate data with an additional Ka-band downlink, were gathered using the HGA firmly pointed towards Earth. From these measurements Iess et al. ( 2019), using a similar approach to the one adopted for Titan gravity flybys analysis [START_REF] Durante | Titan's gravity field and interior structure after Cassini[END_REF], produced a solution based on a multi-arc fit, providing an estimate of the zonal harmonic coefficients of the planet up to degree 12 (see Table 2.4).

Exploiting the exquisite orbit reconstruction obtained within these passes, characterised by a metre-level uncertainty with respect to the planet at the CA, we have produced nine normal points with an unprecedented accuracy of ∼ 3 m. These measurements provide a unique perspective on the Saturnian system due to the peculiar geometry of the Grand Finale proximal orbits. In fact, Cassini in this case directly swung by the planet, differently from what was done for the rest of the mission, during which it only flew by the moons. In this way, the derived normal points provide direct constraints on Saturn theoretically independent from the thorny information of the satellite ephemerides, which have driven the uncertainty on the flybys reconstruction [START_REF] Boone | Resolution of orbit determination prediction instabilities at titan during cassini's solstice mission[END_REF]. Moreover, the inclusion of these observations allowed us to extend Cassini dataset until almost the end of the mission, in July 2017.

INPOP19a

INPOP19a is the latest release of the INPOP planetary ephemerides. It includes the orbits of the eight planets of the Solar System, Pluto and the Moon, plus the Earth and Moon librations. Their orbits are constructed solving for their initial conditions at J2000 with a WLS filter and propagated using an Adams-Cowell integrator with fixed step (see Sect. 1.3). This latest version is a significant step forward with respect to its predecessor, INPOP17a, in many different aspects. First of all, IN-POP19a benefits from several improvements in the field of solar plasma calibrations [START_REF] Fienga | INPOP19a planetary ephemerides[END_REF], which are crucial for the correct processing of ranging data in order to avoid the introduction of systematics in the solution.

Another important innovation is the determination of the main belt asteroids masses estimated with a least squares fit with a priori information from the asteroids spectra (Fienga et al., 2019a,b). The analysis combines space and ground-based observations of the physical properties of the asteroids with the INPOP planetary ephemerides dataset for estimating their mass in a BVLS filter along with the planets orbits. The use of a BVLS limits the solution to given intervals defined by the asteroids taxonomic complexes C, S and X from the MP 3 C catalog [START_REF] Delbo | Mp3c catalog[END_REF]. As a result, INPOP19a provides constraints on the masses of 343 asteroids of the main belt, with respect to the 168 of INPOP17a (see [START_REF] Viswanathan | Inpop17a planetary ephemerides[END_REF] and Sect. 1.3.1). The choice of these objects is based on [START_REF] Kuchynka | A new approach to determining asteroid masses from planetary range measurements[END_REF] work.

Another important feature introduced with the INPOP19a release is the inclusion of the observations from the ESA mission Gaia data release 2 (Gaia Collaboration et al., 2018). The orbits of 14 099 asteroids are thus integrated by fitting this data, using the correlation matrix provided by the Gaia DPAC. For convenience, the orbits of these asteroids are not integrated together with the planets, but we iterate the process to account for the updated asteroids orbits and assure the coherence of the system, particularly for Mars. New adjustments of the Moon orbit and rotation [START_REF] Fienga | INPOP19a planetary ephemerides[END_REF] from LLR data should be also noted. At last, a significant improvement over previous releases derives from the enhanced dataset of planetary observations upon which the INPOP19a ephemerides are built.

Dataset

The INPOP19a full dataset is presented in Table 2.5, where the different type of measurements used for the construction of the ephemerides are shown. These include both ground-based optical observations, both VLBI and ranging data gathered by a dozen of interplanetary missions, providing almost a century of data on the dynamics of the Solar System. Ranging data, both raw (e.g. from the martians missions) and derived measurements in the form of range biases (e.g. Cassini and Juno), represent the core of the database, with their invaluable measurement of absolute distance. They contribute for the 65.5% of the entire dataset, while ground-based optical observations for the 34.3% and VLBI 0.2%. The major novelty comes from the use of the normal points for the Saturn system barycentre deduced from the analysis of Cassini radio tracking data described in Sect. 2.4.2, and from OD solutions of the gravity Grand Finale passes and Titan flybys (see Sect. 2.4.3 and 2.4.4). In addition, new normal points for the Jovian barycentre obtained with the gravity science experiment of the Juno mission [START_REF] Iess | Measurement of jupiter's asymmetric gravity field[END_REF] are now available, increasing of orders of magnitude the accuracy of Jupiter ephemerides in INPOP19a. In fact, these points, produced from the first 9 gravity-dedicated perijoves exploiting the Ka-band capabilities of the Juno radiosystem, feature an unprecedented metres-level accuracy.

Beside these points, the dataset includes ten VLBI measurements and 165 normal points from Cassini mission deduced by the JPL Doppler-only analysis of navigation data described in [START_REF] Hees | Constraints on modified Newtonian dynamics theories from radio tracking data of the Cassini spacecraft[END_REF]. The observations span from 2004 to 2014, with an estimated RMS of 25 m. The ∼4 times larger uncertainty of JPL solution is explained by the choice of [START_REF] Hees | Constraints on modified Newtonian dynamics theories from radio tracking data of the Cassini spacecraft[END_REF] to analyse shorter trajectory arcs, cut at each OTM or targeted encounter, and thus unconstrained by the moons flybys.

Table 2.5 offers also a comparison of the WRMS of the post-fit residuals obtained INPOP19a augmented dataset, especially the new Saturn's and Jupiter's normal points, demands for a substantial update of INPOP dynamical model. In fact, in order to fit the whole database with the 13 years time frame of Cassini measurements, we introduce a new modelling of the outer Solar System, in particular of the objects beyond the orbit of Neptune.

Update of INPOP dynamical model

In order to produce a precise representation of the Solar System with a new release of planetary ephemerides, including the extremely accurate measurements gathered by the latest interplanetary missions, an update of INPOP model for the trans-Neptunian Solar System was required. For this reason, INPOP19a includes a detailed modelling of TNOs distribution in the outer space, in order to simulate their gravity perturbations on the planets. The setup consists of nine binary TNOs included in INPOP dynamical model; these objects, listed in Table 2.6, are the same considered for the construction of the EPM ephemerides produced by [START_REF] Pitjeva | Mass of the kuiper belt[END_REF]. Along with the main belt asteroids, the orbits of these TNOs are integrated together with those of the planets; while their masses are not adjusted in the fit, but fixed to their a priori value. They, in fact, all have at least one natural satellite, and their masses are thus very accurately measured by studying their moon dynamics.

In addition to this, INPOP19a updated dynamical model considers the average gravity perturbation produced by the TNOs enclosed in the 3:2 and 2:1 mean motion resonances with Neptune, forming the Kuiper belt. This is done by including the three circular, not inclined rings located at 39.4, 44.0, and 47.5 AU, shown in Fig. 2.10. Each of them contributes for one-sixth, two-thirds, and one-sixth of the total mass, respectively, following the same scheme adopted by [START_REF] Pitjeva | Mass of the kuiper belt[END_REF]. The central ring at 44.0 AU has more mass, as it represents the sum of the two populations of objects, the resonant and classical KBOs, with semi-major axes between 39.4 and 47.5 AU.

Results

The results of the new modelling of the trans-Neptunian Solar System on Cassini data is clearly visible in Fig. 2 For INPOP17a solution, the detailed modelling of TNOs perturbations on the outer planets was not required, since BepiColombo points were not available and the time span of Cassini data was limited to 10 years only (from 2004 to 2014), and with significant lower accuracies (∼25 m for JPL data, [START_REF] Hees | Constraints on modified Newtonian dynamics theories from radio tracking data of the Cassini spacecraft[END_REF]). The addition of the more accurate Grand Finale points and the reanalysis of navigation data have changed the scenario, extending the data sample to over 13 years, almost the full duration of Cassini mission. In this situation, the limited model of INPOP17a is not able to reproduce the data anymore, showing strong signatures in the post-fit residuals, including a relevant bias on Grand Finale points. Such a trend disappear when the binary TNOs and, especially, the TNO ring are included.

Figure 2.11 also offers a qualitative idea of the sensitivity of the different datasets to the trans-Neptunian Solar System modelling. In particular, the top left plot shows the scarce sensitivity of JPL points; a proper fit of these data is achievable without introducing any TNOs perturbation, as obtained with INPOP17a solution. An evidence of this aspect is also given in Table 2.7, where the limited variation registered on the WRMS of these data from INPOP17a and INPOP19a (from 1.1 to 0.99) is reported.

On the contrary, Grand Finale points exhibit the highest sensitivity because of their unique time frame and peculiar geometry perspective. The WRMS in this 2.5) obtained with different INPOP solutions. In particular, the plots show the residuals obtained with INPOP17a, where no modelling of TNOs is included (top left); a solution accounting for the individual perturbations of the 9 binary TNOs (top right); a solution including the TNO ring but none of the binary TNOs (bottom left); and INPOP19a, which includes both the ring model and the binary TNOs (bottom right).

Results

51

case shows a significant improvement, from a value of ∼30 for INPOP17a and close to unity for INPOP19a. Most of the enhancement derives from the introduction of the TNO ring in the dynamical model and the estimate of its mass; however, the contribution of the binary TNOs results to be significant for the early navigation data (Feb.-Jul. 2006) and for the Grand Finale measurements, providing a reduction of the residuals' WRMS from 2.58 to 1.14.

Finally, it is worth noting how the normal points deduced from the Titan gravity flybys (see Sec. 2.4.3) show a dispersion of the residuals larger than their expected accuracy, even when considering the INPOP19a solution (see Fig. 2.11,bottom right). This can be explained by the fact that they are significantly less than the navigation data (42 vs 572), sharing overall the same accuracy. Moreover, the navigation points are mostly concentrated in a limited interval of time (2008)(2009), which leads the filter to produce a better fit of these data to the detriment of the Titan flybys gravity residuals, which are distributed on a larger interval.

However, the contribution of this dataset to our estimate of the Kuiper belt mass is limited. We performed a test by removing these data from the fit, obtaining a solution that is compatible with the original result.

The estimate of the Kuiper belt mass

As mentioned in Sect 2.5.1, while the masses of the nine binary TNOs are fixed, the mass of the three rings modelling the Kuiper belt is estimated in INPOP19a fit, together with the initial conditions of the orbits of the planets and the other estimated parameters. In fact, although the high correlation (∼ 98.5%) between the different masses of each ring leads to an impossible estimation of the separate contributions, their cumulative mass is well determined. With the INPOP19a extended dataset (Table 2.5) we obtain for the Kuiper belt a mass of M ring = (0.061 ± 0.001)M .

(2.42)

Two types of analyses have hitherto been published: those based on KBOs direct observations and those deduced from the KBOs perturbations on planetary ephemerides. In [START_REF] Bernstein | The Size Distribution of Trans-Neptunian Bodies[END_REF], based on HST observations, the authors deduced a distribution of sizes and infer surface density values, finding a mass of the Kuiper belt of about 0.010 M when considering only KBOs with inclination smaller than 5 • . If we add their estimations of excited KBOs (objects with inclination greater than 5 • ), the total mass deduced from [START_REF] Bernstein | The Size Distribution of Trans-Neptunian Bodies[END_REF] becomes 0.018 M . The value is low compared to our estimate, but also in comparison with the [START_REF] Gladman | The Structure of the Kuiper Belt: Size Distribution and Radial Extent[END_REF] results. In the latter work, the authors used a different assumption for the size distribution, leading to an estimation of the mass of the total Kuiper belt in between 0.04 and 0.1 M . This interval of masses nicely frames our result, which was obtained completely independently.

Analogously to this work, [START_REF] Pitjeva | Mass of the kuiper belt[END_REF] study KBOs perturbations on EPM planetary ephemerides. The authors analyse two different model for the ring: one numerical, that is the one our analysis is based on (exposed in Sect. 2.5.1), and one analytical; both giving consistent results. They find a value of M P18 ring = (0.01108 ± 0.0025)M , which is significantly lower than our result. However, although their analysis shares most of the dataset with INPOP, it does not include Juno data and the Cassini normal points produced within our work, but only those from the JPL data processing [START_REF] Hees | Constraints on modified Newtonian dynamics theories from radio tracking data of the Cassini spacecraft[END_REF]. Therefore, if we limit our data sample to that used by [START_REF] Pitjeva | Mass of the kuiper belt[END_REF], we obtain a mass of M ring = (0.020 ± 0.003)M , (2.43) which is now consistent at 3σ with their value.

It is worth noting that the masses of the major TNOs included in our model (see Table 2.6) are fixed in INPOP adjustment, while 31 TNO masses are fitted itn [START_REF] Pitjeva | Mass of the kuiper belt[END_REF]. It is possible that part of TNO masses in INPOP are absorbed by the TNO ring mass, inducing a slightly bigger value than that obtained by [START_REF] Pitjeva | Mass of the kuiper belt[END_REF]. If we add the masses of the fixed TNOs to the mass of the ring, the differences decrease as we obtain M total = (0.0243 ± 0.003)M , (2.44) which is then compatible at 2σ with the estimated value of [START_REF] Pitjeva | Mass of the kuiper belt[END_REF], M P18 total = (0.0197 ± 0.0035)M .

Although it is useful to assess the two results, the sum of masses (M total = M ring + M TNOs ) does not provide a rigorous comparison of our estimate with the result from [START_REF] Pitjeva | Mass of the kuiper belt[END_REF]. Indeed, the impossibility to fully disentangle the tidal effect on the planets produced by TNOs with different distances and masses lead to an estimate of a smaller mass for a closer distribution and, vice versa, a higher mass for farther located objects. Therefore, since the distribution of the 31 objects included in [START_REF] Pitjeva | Mass of the kuiper belt[END_REF] fit is not reported, this discrepancy remains a strong candidate for explaining the divergence of the two results.

In conclusion, all these masses are consistent with the Kuiper belt mass obtained from simulated populations analysis, such as [START_REF] Levison | Origin of the structure of the Kuiper belt during a dynamical instability in the orbits of Uranus and Neptune[END_REF].

Alternative ring model

Although we chose as reference model for simulating the Kuiper belt gravity perturbation in INPOP19a the three-rings configuration from [START_REF] Pitjeva | Mass of the kuiper belt[END_REF] (see Fig. 2.10), we explored alternative options, looking for the most suitable representation. In particular, an interesting solution tested, it is the one based on the sampling of the real distribution of known KBOs orbits. A view of this alternative configuration is provided in Fig. 2.12, where the green dots show the resulting spatial distribution of the objects constituting the ring.

By adopting this model we obtain an alternative solution to INPOP19a, with an estimated mass of the Kuiper belt of:

M ring = (0.041 ± 0.001)M , (2.45)
which is about 30% less than the INPOP19a reference solution. A possible explanation for this lower value may comes from the different average distance of the objects constituting the ring in this model with respect to [START_REF] Pitjeva | Mass of the kuiper belt[END_REF] configuration. In fact, in the three-rings model most of the mass (two-thirds) is concentrated in the central ring, at 44.0 AU (see Fig. 2.13 (a)), while in the alternative model, a conspicuous number of objects are located around 39.5 AU, resulting in a mass distribution closer to the planets (see Fig. 2.13 (b)). Moreover, as a result of considering real KBOs orbits, the alternative ring is slightly inclined with respect to the ecliptic and eccentric, this also may play a role in the different outcomes.

Constraints on P9 location

We conducted a series of simulations aimed at constraining the possible location of P9 in the outer Solar System by exploiting the cutting-edge precision offered by INPOP19a planetary ephemerides (see Sect. 2.5). In particular we analysed the effect produced by the potential gravity perturbation of P9 on the orbits of Saturn, constrained by 13 years of extremely accurate data from the Cassini mission (see Sect. 2.4). We thus modelled the planet as a fixed point-mass, due to its great distance from SSB, and we scanned the sky by varying its position in RA and Dec (see Fig. 2.14). For each point we obtained a P9-perturbed INPOP solution that we compared with the reference one, INPOP19a. We used two statistical criteria for selecting the compatible zones: the first one based on the χ 2 likelihood of post-fit residuals of the P9-perturbed solution, and the second based on the computation of the Mahalanobis distance between the reference (INPOP19a) and the P9-perturbed solutions state vectors, normalised with respect to INPOP19a propagated covariance matrix.

P9 modelling

We modelled the gravity perturbation of P9 on the planets as a point mass, fixed in the sky with respect to SSB. In fact, with an expected semimajor axis of 600-800 AU, the orbital period of P9 should vary from 8 000 to 22 500 years, making de facto its movement over the timespan of INPOP dataset negligible. Hence, the induced gravitational acceleration on a planet "A" can be written as

ẍP9 A = GM P9 r 3 [3 (x A • u) u -x A ] (2.46)
where u is the versor defined by the RA and Dec of P9 in ICRF, x A is planet "A" position vector, r is P9 distance from SSB and M P9 is its mass.

The Mahalanobis distance method

The first criterium we chose for assessing the divergency of the P9-perturbed solutions with respect to the reference, INPOP19a, is based on the computation of the distance 8 between the reference and the candidate solution state vectors, normalised with respect to the propagated covariance matrix of INPOP19a. This definition recalls that of the Mahalanobis distance, which is a statistical tool for measuring how many SDs away a point is from the mean of the distribution.

Because of the increased sensitivity of Saturn's ephemerides provided by IN-POP19a dataset, we decided to focus on its orbit to detect the impact of P9. In particular, due to the Earth-referenced nature of the range measurements which constraints Saturn orbit, we compute the Mahalanobis distance directly on the Earth-Saturn vectors of the P9-perturbed solution and INPOP19a. However, the state vector X output of the INPOP solution comprises the barycen-8 is a generalised distance in n-dimension, with n being the size of the state vector 56 2. INPOP19a: updates on the outer Solar System from Cassini tric 9 coordinates of the planets

X = [x Mer , . . . , x EMB , . . . p , GM , J 2 , . . . ] ∈ R n (2.47)
where x A is the barycentric state vector of planet "A" and p represents the number of integrated bodies in INPOP; the size of X is then n = 6p + q, with q being the number of the other constant estimated parameters. Then, in order to properly compare the Earth-Saturn vectors we need to transform the state in an EMB-centred frame. The covariance matrix P needs to be transformed as well, since in this new frame, to the uncertainty of each vector we need to add that of the new origin of the frame, i.e. x EMB .

To do so we define the Jacobian

J(X) = ∂f (X) ∂X = ∂ (X -X EMB ) ∂X (2.48) with X EMB = [x EMB , . . . , x EMB p-1 , 0, 0, . . . ] ∈ R n-1 , X ∈ R n-1
where X is a reduced state vector lacking of the EMB state vector x EMB .

The EMB-centred state vector and covariance matrix at time t are thus computed from the propagated values (see Sect. 2.1.3) as

X E (t) = X (t) -XEMB (t) ∈ R n-1 (2.49) P E (t) = J(t) P (t)J(t) T , (n -1) × (n -1) (2.50)
Considering only Saturn state vector and covariance, we get

x E Sat (t) = x Sat (t) -x EMB (t) (2.51) P E Sat (t) = cov x E Sat (t) = cov (x Sat (t)) + cov (x EMB (t)) -2cov (x Sat (t), x EMB (t))
(2.52) where the term cov(x) represents the variance of x, or its covariance matrix, and cov(x A , x EMB ) is the cross-covariance matrix between x Sat and x EMB . Now, for the i-th P9-perturbed solution, at each time t, we can compute the Mahalanobis distance from the reference solution INPOP19a

d i (t) = x E Sat 19a (t) -x E Sat i (t) • P E Sat (t) -1 • x E Sat 19a (t) -x E Sat i (t) T (2.53)
The compatibility of the i-th solution is thus assessed by measuring the percentage of time the Mahalanobis distance d i (t) is within the 3σ confidence interval of INPOP19a covariance matrix, over the time interval covered by Cassini data (2004-2017, see Table 2.5).

Results

57

The likelihood method In order to asses the impact of P9 on a perturbed solution of INPOP, together with the Mahalanobis distance we used a likelihood test based on the χ 2 of post-fit residuals. This is same criterium used in [START_REF] Bernus | Constraint on the yukawa suppression of the newtonian potential from the planetary ephemeris inpop19a[END_REF] for constraining the Compton wavelength in the framework of a Yukawa suppression of the Newtonian potential. The mχ 2 is defined as ε T W ε, where ε is the residuals vector of size m and W the weighting matrix (see Sect. 2.1.1). If the post-fit residuals are unbiased and follow a normal distribution, then mχ 2 follows an m degrees of freedom χ 2 law, and when m → ∞, for the i-th P9-perturbed solution

z i = m 2 χ 2 i -χ 2 19a → N (0, 1) (2.54)
For the sake of efficiency, we focus on the observations datasets which are sensitive to the perturbation of P9, we can thus define the reduced χ2 computed on the sensitive data only. Assuming that for the non-sensitive sets ε i ≈ ε 19a , then

χ 2 i -χ 2 19a = m m χ2 i -χ2 19a (2.55)
where m is the number of sensitive observations. Since the SD of the residuals of INPOP19a is very close to the instrumental uncertainties (WLS ≈ 1, see Table 2.5), one can set χ2 19a = 1 (see [START_REF] Bernus | Constraint on the yukawa suppression of the newtonian potential from the planetary ephemeris inpop19a[END_REF] for the full demonstration); then we can write

z i = m 2 χ2 i -1 → N (0, 1) (2.56)
We can thus define the likelihood L i of the i-th solution as

L i = 1 - 1 √ 2π z i -∞ exp - x 2 2 dx (2.57)
Therefore, by definition L 19a = 0.5, and any P9-perturbed solution with L ≈ 0.5 is as likely as INPOP19a.

Among the entire INPOP19a dataset (see Sect. 2.5 and Table 2.5), Cassini data result to be the most sensitive. We also include in the sensitive dataset, Mars data, which, although their marginal sensitivity to P9, contribute for the 47% of the entire dataset; Juno points; and MESSENGER data.

Compatibility zones

As a result of the 3 156 simulations of P9-perturbed ephemerides we do not find any clear indication for the presence of the ninth planet. The white portions of the maps correspond to the rejected solutions, with a threshold at 3σ (L < 0.003) [START_REF] Fienga | New constraints on the location of p9 obtained with the inpop19a planetary ephemeris[END_REF]. The results from the likelihood criterium, which could potentially identify ameliorations in the post-fit residuals (for L > 0.5) coming from the P9 gravitational acceleration, are provided in Fig. 2.15, and display the absence of positive markers. The RA-Dec maps, obtained for different P9 configurations with distances from SSB varying from 400 to 800 AU and masses of 5 or 10 M , are in fact almost completely ruled out (white areas), meaning that a significant degradation of the residuals is registered (out of 3σ). However, there are zones, especially for the further P9 options (r P9 > 600 AU), which provide acceptable likelihood. Among these, the case of r P9 =600 AU and M P9 =5M (top left conrner) is particularly interesting, since it best agrees with last theoretical predictions from [START_REF] Batygin | The planet nine hypothesis[END_REF].

By looking at the results from the Mahalanobis distance method shown in Fig. 2.16, we find confirmation of compatibility for the same zones identified by the likelihood test, in particular for the r P9 =600 AU and M P9 =5M configuration (top left conrner). For this case, the direct comparison of the two criteria acceptance zones is reported in Fig. 2.17 (a). The Mahalanobis acceptance area, which is indicated by the magenta lines, overlaps the likelihood one in two zones, which we identify as "D1" and "D2". A magnification of these two compatibility zones is given in Fig. 2.17 (b) and (c), while their coordinates are reported in Table 2.8 [START_REF] Fienga | New constraints on the location of p9 obtained with the inpop19a planetary ephemeris[END_REF]Di Ruscio et al., 2020a). The coordinates here presented correspond to the zones for which the likelihood is larger than 0.05 (2σ).

M P9 r P9 Zone RA Dec 5M 600 AU D1 [18°: 25°] [-5 : -3] D2 [198°: 202°] [5°: 13°]
Chapter 3

Simulations of BepiColombo mission 3.1 BepiColombo mission overview

BepiColombo will be the first European spacecraft to visit Mercury, the innermost planet of the Solar System. The mission is a result of the cooperation of ESA and JAXA space agencies as it consists of two separate spacecrafts, that will study Mercury from two different orbits. The MMO, developed by JAXA, will remain in a larger eccentricity orbit, focusing on the study of Mercury's magnetosphere; while the ESA MPO will study the planet from an inner polar orbit, with a periherm at 480 km and an apoherm at 1 500 km. BepiColombo launched in October 20, 2018 after several delays due to a laborious development, and it will arrive at Mercury in late 2025. During its 7 years-long interplanetary voyage, MPO and MMO, will be escorted by the MTM and the MOSIF, along with they compose the MCS (Fig. 3.1). The MTM, in particular, is designated to provide during the cruise phase solarelectric propulsion and all the other services that will not be necessary in Mercury's orbit, and it will be jettisoned right before the orbit insertion manoeuvre. MOSIF, instead, supplies thermal protection to MMO and serves as mechanical and electrical interface between the two orbiters until their separation at Mercury.

BepiColombo aims at investigating one of the hottest planet of the Solar System, while facing extreme environmental conditions, with temperatures exceeding 350°C. Its development is a result of overwhelming technical challenges. In fact, in addition to the extraordinary endurance required to survive the harsh environment of Mercury, the planet is also incredible arduous to reach. To approach its orbit, it is necessary to decelerate the spacecraft against the Sun's gravitational attraction, in contrast with all the other missions bound for the outer Solar System, and to reach its enormous orbital speed of about 47 km/s. Furthermore, due to Mercury's orbital inclination with respect to the ecliptic, the spacecraft orbital plane shall be changed too. These manoeuvres are performed by adopting a low-thrust solarelectric propulsion, an innovative technology that has been successfully tested by SMART-11 , and by exploiting a series of gravity-assists: one of the Earth (occurred on April 10, 2020; see Fig. 3.2), two of Venus and six of Mercury itself. Once arrived at Mercury, BepiColombo will be inserted in a polar orbit with the use of conventional thrusters, adopting though a special capture technique, based on the "weak stability boundary" theory. This should limit the risks related to the orbit insertion, granting a more reliable solution over the canonical single-burn capture. At this time MMO will be delivered into its mission orbit.

BepiColombo nominal mission is expected to last for one Earth-year, with the possibility of another year extension. During this period, daily tracking from ES-TRACK and DSN stations will grant the collection of a considerable amount of radio science data, regarding one of the lesser studied bodies in the Solar System. Such valuable information shall deepen our knowledge of the planet, and of the history and formation of terrestrial planets, including Earth.

BepiColombo experiments

BepiColombo at launch was assembled in the 6-m tall MCS, for a total of 4.1 tons. Most of the weight is due to the sun-shield and the transfer module that will be ISA Italian Spring Accelerometer is a three-axis high sensitivity accelerometer that will support the study of the planet Mercury's gravity field to an unprecedented level of accuracy.

MPO-MAG Mercury Magnetometer consists of magnetometers on board MPO

and MMO and will gather data by combining measurements from both the spacecrafts.

MERTIS Mercury Radiometer and Thermal Infrared Spectrometer will provide detailed information about the mineralogical composition of Mercury's surface layer.

MGNS Mercury Gamma-ray and Neutron Spectrometer main objective is to determine the elemental compositions of distinguishable regions over the entire surface of Mercury.

MIXS Mercury Imaging X-ray Spectrometer has three main scientific objectives:

to produce global elemental abundance maps of key rock-forming elements, to perform high spatial resolution mapping of these elemental abundances and to confirm that the auroral zone is an intense source of continuum and line X-rays.

MORE Mercury Orbiter Radio science Experiment addresses BepiColombo's scientific goals in geodesy, geophysics and fundamental physics and it will help to determine the gravity field of Mercury as well as the size and physical state of its core. MORE will also provide for a series of general relativity tests exploiting the amazing opportunities offered by Mercury's orbit.

PHEBUS Probing of Hermean Exosphere by Ultraviolet Spectroscopy uses UV emission from Mercury's exosphere to characterise its composition, structure and dynamics, and surface-exosphere connections.

SERENA Search for Exospheric Refilling and Emitted Natural Abundances can provide information on the whole surface-exosphere-magnetosphere coupled system and the processes involved, plus the interactions between energetic particles, the solar wind, micrometeorites and the interplanetary medium.

SIMBIO-SYS Spectrometer and Imager for MPO Bepicolombo integrated Ob-

servatory System is an integrated suite for the imaging and spectroscopic investigation of the Hermean surface that aims to examine the surface geology, volcanism, global tectonics, surface age and composition, and geophysics of Mercury.

SIXIS Solar Intensity X-ray and particles Spectrometer whose objective is to perform measurements of X-rays and particles of solar origin at the position of BepiColombo

Mercury

Mercury is the smallest planet in the Solar System and the closest to the Sun (Table 3.2). It is thus the fastest planet to orbit the star, with a mean velocity of approximately 47 km/s; probably for this reason ancient Romans decided to name it after the swift-footed messenger god. Its short distance from the Sun, that during its 88 days-long elliptical orbit oscillates between 46 and 70 million km, makes him one of the hottest body of the Solar System, with scorching surface temperatures that can reach up to 450°C on the sunward side. However, the lack of a substantial atmosphere ables to trap this heat, causes the greatest temperature swings in the Solar System, registering temperatures of minus 170°C at night. Mercury, in fact, possesses only an exosphere, made up of the atoms blasted off the surface by the solar wind and striking micrometeoroids.

The solar flux hitting the planet is indeed so strong to pull out any attempt of atmosphere and to form a tail of neutral particles on the sun-opposing side. Without an atmosphere to protect it from incoming asteroids Mercury's surface is thus, alike the Moon, filled with craters of past impacts, conferring it its peculiar appearance.

Mercury is extremely challenging to reach and study, due to its 7°off-ecliptic, highly energetic orbit and its severe conditions. This makes it the less investigated planet in the inner Solar System. Only two missions in the history of human space exploration visited it, Mariner 10 in 1974 and, more recently, MESSENGER. Both of them were NASA missions that, beside unveiled a series of interesting enigmas about the planet, rose a cascade of new ones. One of the most extraordinary discovery carried out by Mariner 10 was that Mercury actually possesses its own magnetic field [START_REF] Ness | Magnetic field observations near mercury: Preliminary results from mariner 10[END_REF][START_REF] Ness | Magnetic field of mercury confirmed[END_REF], in contrast with the scientific community expectations. In fact magnetic fields generally arise in presence of a fast rotation and a molten core. Mercury instead is a slow rotator, with a rotational period of approximately 59 Earth-days, and its solid core is expected to have cooled off long ago. In addition, the field results to be strongly asymmetric with respect to the equator, being in the northern hemisphere approximately three times stronger than in the southern. Vice versa, all the other planets possessing an own magnetic field, show a quasi-perfect symmetry, with only very little difference between the two hemispheres. A solution to this anomaly, based on data collected by MESSENGER, has been proposed in [START_REF] Cao | A dynamo explanation for mercury's anomalous magnetic field[END_REF] and explains it with convective dynamo models driven by volumetric buoyancies, able to generate the quasi-steady north-south asymmetries seen at Mercury. Mercury, with its large metallic core, estimated to be the about the 80% of planet's diameter [START_REF] Genova | Geodetic evidence that mercury has a solid inner core[END_REF], is incredibly dense, with a density second only to Earth. Moreover, as the its core cools and solidifies, the planet shrinks, in a process that is shaping its surface with huge scarps and cliffs, some hundreds of km long. The young age of some of these features indicates that the planet, as the Earth, might still be tectonically active.

Mercury's orbit has long been an enigma for the astronomers. The magnitude of its perihelion precession is indeed impossible to explain by means of Newtonian physics alone. In order to correctly model this anomalous rotation, in fact, it is necessary to appeal to general relativity theory. With this in mind, Albert Einstein, in 1916, included the perihelion precession of Mercury's orbit among the three classical tests of general relativity, along with the deflection of light produced by the Sun, and the gravitational redshift effect of light. The classical third bodies effect of the other planets in the Solar System, together with the Sun oblateness, is indeed only responsible for 531.66 ± 0.69 arcsec per Julian century of the 574 .64 ± 0.69 total precession. The remnant 42 .98 ± 0.04 are related to a relativity effect, i.e. the gravitation being mediated by the curvature of spacetime [START_REF] Einstein | Erklarung der Perihelionbewegung der Merkur aus der allgemeinen Relativitatstheorie[END_REF]). Mercury's peculiar orbit makes it an extraordinary laboratory for fundamental physics science, as foreseen by Einstein. BepiColombo, in fact, with its radio science experiment MORE will exploit this environment providing for a series of cuttingedge tests of general relativity theory.

The radio science experiment MORE

The Mercury Orbiter Radio science Experiment is a complex system involving hardware and software solutions, both on board the spacecraft and at the ground station.

The key element is represented by the KaT [START_REF] Iess | Advanced radio science instrumentation for the mission bepicolombo to mercury[END_REF] on board MPO, which enables to establish a triple multi-frequency link configuration (X/X, X/Ka and Ka/Ka), granting the complete cancellation of the dispersive noises contribution from plasma and ionospheric scintillation (see Sec. 2.2.2). MORE KaT, along with the DST [START_REF] De Tiberis | The x/x/ka-band deep space transponder for the bepicolombo mission to mercury[END_REF], indeed constitute the most advanced radio system ever flown in space. It successfully passed the commissioning phase in December 2018, and its performance are now under scrutiny with dedicated campaign tests 2 . MORE aims at estimating a series of physical quantities of geodesy, geophysics and fundamental physics interest, such as Mercury's gravity field, PPN parameters of 2 The first end-to-end test campaign is scheduled for the 4-11 September 2020, and the second for October 2020 (MORE Team, 2011).

Gravity field experiment goals (geodesy)

• Gravity field's spherical harmonic coefficients up 1σ < 10 -9 to degree and order 25 Degree general relativity, tidal Love numbers, etc. (Table 3.3). This is made possible by means of extremely accurate range and range-rate observables (see Sec. 2.2.2). In particular, these parameters will be estimated through a series of focused experiments: the gravimetry experiment, with the goal of determining the gravity field of Mercury and its rotation state, and the relativity experiment, to determine PPN parameters and other quantities, such as the mass and the dynamic oblateness of the Sun [START_REF] Imperi | An analysis of the geodesy and relativity experiments of bepicolombo[END_REF]. The combination of frequent solar conjunctions with the adoption of a multi-frequency radio system, makes MORE a cutting-edge experiment for fundamental and solar corona physics studies.

Gravity experiment

Mercury's gravity field has been studied by both the previous NASA missions that visited the planet, however, Mariner 10 only provided a low accuracy estimate of the degree and order 2. MESSENGER instead, during its 3-year mission, has thoroughly analysed the planet, providing a detailed characterisation of its gravity field and internal structure [START_REF] Genova | Geodetic evidence that mercury has a solid inner core[END_REF]. BepiColombo, with MORE, will deepen our knowledge of Mercury's gravity field, increasing the accuracy measurements thanks to the triple link configuration and the favourable orbit geometry, which will grant a global coverage of the planet.

The different wavelengths of the gravity field, i.e. the different spherical harmonic degrees, permit to probe the various aspects of the planet's structure. In particular, long wavelengths combined with the rotation state can constrain the internal structure, size and physical state of the core, improving the results from MESSEN-GER [START_REF] Genova | Geodetic evidence that mercury has a solid inner core[END_REF]; while the higher frequency terms, in combination with altimetry data from BELA, allow the structure of the crust and its interface with the mantle to be examined. The gravity field experiment will mostly rely on plasma-calibrated range-rate observables, daily acquired from ESTRACK DSA-3 antenna, at Malargüe, in Argentina. Additional support will be provided by DSN stations, in particular from DSS 25 at Goldstone, California. The requirements of these measurements, in terms of noise budget, are shown in Table (3.4). The final level, comprehensive of the whole link variability, is 1.38 × 10 -14 at 1000 s of integration time. In order to achieve this level of accuracy, a partial calibration (∼80%) of the tropospheric scintillation with water-vapour radiometers is required. Moreover, for the purpose of compensating the strong non-gravitational accelerations arising at Hermean environment (e.g. SRP and planetary albedo) MORE will be supported by the measurements of ISA.

We performed a complete end-to-end simulation of the Hermean phase of MORE experiment in order to evaluate the expected results achievable with the nominal and extended missios. In particular, we conducted a covariance analysis to provide a realistic prediction of the formal accuracies on the estimated parameters. The setup used in the simulations includes a 50×50 gravity field of Mercury based on MESSENGER's HgM008 solution [START_REF] Genova | Geodetic evidence that mercury has a solid inner core[END_REF], point mass gravity perturbations from all the other planets, GM and oblateness of the Sun, general relativity effects, SRP and desaturation manoeuvres. For reconstructing the orbits of the planets and the Moon we used INPOP17a planetary ephemerides (see Sec. 1.3.1). Both the dynamical model and the observation model, used for simulating the ob-servables in a relativistic context, are based on JPL's MONTE software [START_REF] Evans | Monte: the next generation of mission design and navigation software[END_REF]. The measurements are then processed using the Sapienza radio science laboratory's minimum variance filter Oracle, with a constrained3 multi-arc approach (see Sec. 2.1.2), consisting in a division of the MPO trajectory in 365 arcs of 24 h each. The estimated parameters include the spacecraft state vector at the central epoch of each arc and the ∆V s of the desaturation manoeuvres, two per arc, as local parameters; Mercury's gravity field up to degree and order 50, the tidal Love number k 2 and the planet's obliquity and libration as global parameters. In addition, calibration coefficients of ISA's low-frequency (planet's orbital period) and high-frequency (MPO's orbital period) error are estimated, respectively as local and global parameters.

Support to BELA for measuring Mercury's h 2

One of the tasks of MORE is to provide reliable and accurate reconstructions of MPO orbit around Mercury for the other instrument teams. The baseline requirement is set to 1 m on the radial direction, however, with such accurate measurements, MORE will be capable of much higher accuracy (see Fig. 3.4). In particular, BELA altimetric measurements strongly rely on the precise reconstruction of the MPO orbit. In the framework of a collaboration with the instrument team, we provided an ensemble of trajectories to be used in the generation of BELA synthetic observables, aimed at assessing the expected accuracy on Mercury's tidal Love number h 2 . Starting from the results of the end-to-end simulation of MORE gravity experiment, we perturbed the six components of the spacecraft state vectors of each arc with 100 error realisations. The errors δl are produced as samples of random variables following a multivariate Gaussian distribution,

f (δl i ) = 1 (2π) 6 detP i exp - 1 2 δl T i P -1 i δl i (3.1)
where P i is the covariance submatrix of the spacecraft state vector of the i-th arc retrieved from OD with Oracle.

The results of the simulations show that, with BELA measurements supported by MORE, Mercury's h 2 can be determined with an absolute accuracy of ±0.012 [START_REF] Thor | Prospects for measuring mercury´s tidal love number h2 with the bepicolombo laser altimeter[END_REF]. and Ka-band (scientific) tracking periods, respectively. The black dashed lines indicate the two reaction wheels desaturation manoeuvres expected on each arc [START_REF] Imperi | An analysis of the geodesy and relativity experiments of bepicolombo[END_REF][START_REF] Thor | Prospects for measuring mercury´s tidal love number h2 with the bepicolombo laser altimeter[END_REF].

Application of TDMC on MORE gravity experiment

As discussed in Sec. 2.2.2, the dominating noise sources in a plasma calibrated Doppler link, such as the one of MORE, are the local noise, which include tropospheric and antenna mechanical disturbances. These noise contributions appear in the link through their respective transfer functions, based on their nature and the specific configuration of the link. The two-way link architecture, for example, introduces a characteristic modulation in the power spectrum of local noises where T is the signal RTLT. The same antenna is in fact used both to transmit and receive the signal, hence, all those noises arising at the ground station will echo in the downlink after a RTLT. Tropospheric and mechanical disturbances introduce a positive correlation +δ(t -T ), while clock jitters are anticorrelated -δ(t -T ).

h 2 (t) = δ(t) ± δ(t -T ) (3.2)
In 2008, [START_REF] Armstrong | Reducing antenna mechanical noise in precision spacecraft tracking[END_REF] showed that this characteristic can be exploited by suitably combining simultaneous two-and three-way measurements, in order to suppress tropospheric and mechanical noise from the two-way station, resulting in a significant improvement of the overall Doppler link stability. By virtue of its capabilities, the authors refer to the technique with the name Time-Delay Mechanical-noise Cancellation, or TDMC. The Doppler frequency fluctuations for the two-way and three-way links can be schematised, neglecting lower order noise sources, as follows

y 2 (t) = M 2 (t) + M 2 (t -T ) + T 2 (t) + T 2 (t -T ) + C 2 (t) -C 2 (t -T ) + y s (t) (3.3) y 3 (t) = M 3 (t) + M 2 (t -T ) + T 3 (t) + T 2 (t -T ) + C 3 (t) -C 2 (t -T ) + y s (t) (3.4)
where M * , T * and C * represent the antenna mechanical, tropospheric and clock noise, respectively, while y s (t) is the two-way Doppler signal corrected for the ground stations motion. In this formulation, the modulation introduced by the transfer function of Eq. (3.2) is clear. It is worth noting that the echos of the noise introduced by transmitting

Simulations of BepiColombo mission

antenna participate in both links. It is thus possible to combine the two signals for removing these terms

E(t) = y 3 (t) + y 3 (t -T ) -y 2 (t -T ) = (3.5) = M 3 (t) + M 3 (t -T ) + T 3 (t) + T 3 (t -T ) + C 3 (t) -2C 2 (t -T ) + y s (t) (3.6)
As a result of the combination, the new observable E(t) resembles in its structure a fictitious two-way link carried out by the listen-only antenna, operating both as transmitter and receiver. The local noises from the two-way station are thus replaced by those of the three-way one. Hence, the choice of a smaller and stiffer listen-only antenna, located in a favourable dry location, allows a potential orderof-magnitude enhancement of the measurements accuracy, resulting from an abrupt reduction of local noises. We should mention that, as downside of the combination, the ground station clock noise is doubled in E(t).

Between many promising candidates to ancillary (three-way) antenna, the most interesting results to be LLAMA, a 12-m astronomic antenna still under construction, located in the arid desert of Puna de Atacama in Argentina, one of the driest place on Earth, at an altitude of 4 820 m. The antenna is provided by Vertex, the same company who built the near APEX, with which it shares the same characteristics. Analyses from [START_REF] Greve | Mechanical measurements of the alma prototype antennas[END_REF] conducted on APEX showed a factor of 10 greater mechanical stability with respect to DSS 25 standards, due to the limited dimension of the dish and the rigidity and low thermal expansion coefficients of the adopted materials. Moreover, a thorough characterisation of Atacama troposphere has been conducted over the years by ESO with dedicated water-vapour radiometers [START_REF] Eso | Apex weather archive query form[END_REF], showing order-of-magnitude improvements over DSN's Goldstone complex [START_REF] Keihm | Water vapor radiometer measurements of the tropospheric delay fluctuations at goldstone over a full year[END_REF] under favourable conditions (austral winter). BepiColombo results to be the most suitable mission for the application of TDMC for two main reasons. In first place, the relative proximity of Mercury to Earth, with a mean distance throughout the year of 0.34 AU, translates into a limited RTLT of ∼600 s, this offers a larger availability of combined data. Secondly, the distinct vocation of the mission to radio science, expressed through a cutting-edge radio system, will grant outstanding accuracy in Doppler measurements. In particular, all the side noise sources have been reduced to state-of-the-art level; therefore, the Doppler variability results to be primally dominated by antenna mechanical and tropospheric noise. These are indeed the ideal conditions for the application of TDMC, since the reduction of these noise thus entails a significant reduction of the whole noise contribution. Therefore, we performed a complete end-to-end simulation of the gravity experiment the oscillations of the gain factors depicted in the rightmost plots (blue lines). However, an overall gain factor of ~2 is observed, producing precious benefits for all those instruments whose measurements rely on a precise knowledge of the SC position, e. g. the BELA altimeter. for assessing the potential impact of TDMC on MORE measurements. In particular, we compared the results of the baseline two-way configuration from ESTRACK DSA-3 at Malargüe, with those achievable by applying the TDMC technique in a three-way configuration involving DSA-3 and LLAMA.

In order to provide reliable results, we produced a realistic simulation of the different noise contributions of the error budget based on [START_REF] Notaro | Feasibility of an innovative technique for noise reduction in spacecraft doppler tracking[END_REF].

As result of the application of TDMC on BepiColombo gravity experiment, a general enhancement of ∼40% on the estimate of gravity coefficients is registered [START_REF] Di Ruscio | Improvements in BepiColombo and JUICE radio science experiments with a multi-station tracking configuration for the reduction of Doppler noise[END_REF]. In particular a gain of 1.6 is observed on the estimate of Mercury's GM and 1.7 for J 2 and k 2 (see Tab. 3.5). Moreover, the same level of amelioration registered on J 2 is measured on the higher degree coefficients of gravity field as well. This means that the adoption of TDMC could provide an higher resolution of MORE gravity field estimate, enabling refined analyses of the planet's crust structures and their correlations with topography. For what concerns Mercury's rotational state, gains factor of ∼1.7 are achieved on the estimate of RA and Dec of the pole and the longitudinal librations, potentially providing tighter constraints on the interior structure of the planet. The enhancements on the positioning of the spacecraft with respect to the planet are depicted in Fig. 3.6. Since the spacecraft state vector is a local parameter, it is more affected by the diurnal and seasonal variations simulated for the tropospheric noise at the ground stations and this is clearly visible in the oscillations of the gain factors depicted in the rightmost plots (blue lines). Nonetheless, an overall gain factor of ∼2 is observed, producing precious benefits for all those instruments whose measurements rely on a precise knowledge of the spacecraft position, e.g. BELA.

Relativity experiment

The MORE relativity experiment relies on the precise estimate of Mercury's orbit by exploiting the extraordinary accuracy of BepiColombo radio tracking measurements. In fact, Mercury's dynamics in its revolution around the Sun can be used to put new limits on Einstein's theory of general relativity and to test alternative theories. MORE will be able to constrain the orbit of the planet to unprecedented levels, resulting in a significant step forward in planetary ephemerides construction. The experiment cannot prescind from an accurate solution for the ephemerides of the Earth as well, because of the Earth-referenced nature of the radiometric measurements involved in the estimation process. With these premises, the OD problem results to be ill-conditioned, due to an ap-proximate rank deficiency. To overcome this problem, [START_REF] Imperi | An analysis of the geodesy and relativity experiments of bepicolombo[END_REF] solve-for a constrained set of parameters of the Earth's and Mercury's vectors, 8 out of 12, blocking the Earth's position x, y and z, and its component of the velocity, v z . However, this approach could potentially introduce systematics in the solution, due to the level of accuracy expected by BepiColombo ranging system will make the ephemerides uncertainties no longer negligible. Differently from sequential ranging, whose precision and reliability are severely affected by path delays from the ground station electronics, the new 24 Mcps PN ranging system of MORE will represent a breakthrough for Solar System dynamics studies. Its ability to provide "absolute" measurements of distance between the Earth and Mercury with a relative accuracy of few centimetres [START_REF] Cappuccio | Report on first inflight data of bepicolombo's mercury orbiter radio-science experiment[END_REF], ables to compete with laser systems, will change the strategy nowadays followed in planetary ephemerides construction. In the current approach, the radio science or navigation teams produce a converged orbital solution and provide derived measurements, i.e. normal points, to the ephemerides developer, which uses them for updating the set of ephemerides. In this way, which is the same followed for the construction of INPOP19a in Sect. 2.4, no iterations are expected between the orbital and ephemerides solutions and, thus, there are intrinsic inconsistencies which could introduce systematic errors. In fact, different softwares and assumptions are used for the orbital fit and the ephemerides adjustment. So far, this problem has been neglected due to the use of the less precise X-band sequential ranging measurements, that are affected by meters-level calibration errors [START_REF] Border | Station Delay Calibration for Ranging Measurements[END_REF] which heavily degrade the observables accuracy. However, with the performances of the innovative ranging system of the BepiColombo radio science experiment, if we neglect the potential inconsistencies between the two models, and, especially, we do not iterate between the OD and ephemerides fit, there is a concrete risk of obtaining a biased solution.

In this thesis, we work to ameliorate the compatibility between MONTE and INPOP, by aligning the different mathematical models, by integrating the INPOP numerical integrator into MONTE, and by developing a series of tools, such as the propagation of INPOP covariance matrix. As future work, thanks to the close collaboration between GéoAzur and Sapienza Radio Science laboratories, a new paradigm can be explored, a new approach that involves the orbital fit and the ephemerides estimation in a closed iterative loop, with INPOP integrator providing the planets orbits and MONTE's dynamical and observations models producing the spacecraft trajectory and the computed observables. The filter Oracle can then process the residuals (observed-computed observables) and compute the best estimate corrections. These are re-initialised in INPOP and MONTE and iterated to reach the convergence. The self-consistency of this approach directly avoids the risk of introducing systematics and bias in the solution. Moreover, this method represents the prelude for a combined multi-spacecraft fit, that would significantly improve the estimate of relativistic parameters. The increased dataset timeline and different geometries involved would indeed reduce the uncertainties and the correlations between the parameters [START_REF] De Marchi | Testing general relativity in the solar system: present and future perspectives[END_REF]). This will be an unprecedented analysis, since no attempts to process both Doppler and range data from different missions in one unique solution are registered to date. It is in fact an extremely complicated task, which requires a thorough knowledge of the deep space orbit determination problem and experience in planetary ephemerides construction. Beside MESSENGER and BepiColombo that offer an almost 20 year-long baseline on the innermost planet of the Solar System (2011-2027), the best candidates to be included in this type of analysis are NASA MRO, which is still operating after its Mars orbit insertion in 2006, and, for a future perspective, VERITAS, a newly accepted NASA Discovery mission that is expected to enter in orbit around Venus in 2028 and orbit the planet for 4 Venus cycles (2.7 years) with a MORE-like radio tracking system.

Chapter 4

Conclusions

The thesis has presented the contribution of radio tracking data and precision orbit determination of modern deep-space missions for the construction of planetary ephemerides. In particular, two applications are given, through the analysis of navigation and radio science data of NASA Cassini mission to Saturn, and through the end-to-end simulations of ESA-JAXA mission BepiColombo to Mercury.

For the Cassini case, we have provided the results of the reanalysis of navigation data based on the current knowledge of the Saturnian system. We have shown how, by exploiting the last gravity field solutions for the major bodies of the system [START_REF] Iess | Measurement and implications of saturn's gravity field and ring mass[END_REF][START_REF] Durante | Titan's gravity field and interior structure after Cassini[END_REF] and updated satellite ephemerides (Jacobson, 2016a), and by using the combined information from range and range-rate data, an amelioration of a factor four over previous analyses [START_REF] Hees | Constraints on modified Newtonian dynamics theories from radio tracking data of the Cassini spacecraft[END_REF] is achievable on the reconstruction accuracy [START_REF] Di Ruscio | Analysis of cassini radio tracking data for the construction of inpop19a: A new estimate of the kuiper belt mass[END_REF]. A detailed description of the process we have followed for the production of the derived measurements, i.e. normal points, used for planetary ephemerides construction is also given. In particular, we have shown how, by using stochastic range biases constant per each tracking pass for correcting potential errors on the range observables, it is possible to obtain a reliable reconstruction of the error contribution from planetary ephemerides. This approach has also the advantage to provide a realistic assessment on the accuracy of the normal points, measurements that ordinarily lack of any information on their uncertainties. As a result of the application of this approach, we have enriched INPOP planetary ephemerides database with 572 new points from the reanalysed navigation data, with an average accuracy of 6 m, 42 points from the gravity dedicated flybys of Titan, and 9 from the Grand Finale pericentres passes devoted to measure Saturn's gravity field and rings mass. Based upon this augmented dataset, which also includes new points from the radio 4. Conclusions science experiment of NASA Juno mission to Jupiter, we have built a new set of planetary ephemerides, INPOP19a. Among the numerous novelties introduced in the INPOP19a release, including a revised model for the Moon librations from LLR data and a fit of Gaia observations for the positioning of 14 099 asteroids of the main belt, particularly important is the new description of the trans-Neptunian Solar System. The latter is an addition required for fitting the enlarged dataset of INPOP19a (see Tab. 2.5); indeed, the enhanced accuracy of the new measurements along with their extended timeframe, demanded for a thorough characterisation of the outer Solar System, in particular a precise modelling of the gravity accelerations from the Kuiper belt and the other TNOs perturbing Saturn's orbit. The new model adds nine binary TNOs (see Tab. 2.6) to the list of integrated bodies, and a representation of the Kuiper belt based on [START_REF] Pitjeva | Mass of the kuiper belt[END_REF] model, consisting of three circular rings located at 39.4, 44.0, and 47.5 AU, respectively (see Fig. 2.10). The cumulative mass of the three rings is adjusted in the fit, offering a new estimate of the Kuiper belt mass. We found a value of M ring = (0.061 ± 0.001)M which is compatible with theoretical predictions [START_REF] Gladman | The Structure of the Kuiper Belt: Size Distribution and Radial Extent[END_REF][START_REF] Levison | Origin of the structure of the Kuiper belt during a dynamical instability in the orbits of Uranus and Neptune[END_REF]. In order to substantiate our result, an alternative model for representing the influence of the Kuiper belt on planetary ephemerides has been explored. This is based on a ring populated by objects whose orbits have been retrieved by sampling the real orbits of known KBOs (see Fig. 2.12). In this case, we have obtained a mass for the Kuiper belt of M ring = (0.041 ± 0.001)M , which is lower with respect to INPOP19a result. The discrepancy can be ascribed to the different spatial distribution of the two ring models. These results highlight how, with the current accuracies of planetary ephemerides, the estimate of the Kuiper belt mass is model-dependent, due to the strong correlation between distance and mass, and that only an independent constraint from a different observations set can solve this problem. As a result of the new dynamical model and increased accuracies, INPOP19a offers a comprehensive, state-of-the-art description of the Solar System, with interesting new perspectives on the trans-Neptunian space. Indeed, the enhanced constraints provided on the planets orbits, in particular on Jupiter and Saturn, enables a series of investigations, from fundamental physics [START_REF] Bernus | Constraint on the yukawa suppression of the newtonian potential from the planetary ephemeris inpop19a[END_REF] to planetary science (Di Ruscio et al., 2020a,b;[START_REF] Fienga | New constraints on the location of p9 obtained with the inpop19a planetary ephemeris[END_REF]. In the dissertation, we have presented the results of an analysis aimed at locating a potential ninth planet in the outer Solar System. We have used two criteria to test the compatibility of such a body on the latest ephemerides solution: the first based on the computation of the Mahalanobis distance between the P9-perturbed solution and INPOP19a, and the second based on a likelihood test of post-fit residuals. After 3 156 simulations producing ephemerides perturbed by P9 in different RA-Dec positions, both the methods have not identified any clear, positive markers for the existence of P9; however we have located two regions of the sky, compatible at 3σ with INPOP19a, which could potentially harbour a 5M planet at a distance of 600 AU (see Table 2.8 and Fig. 2.17).

The thesis offers also an overview of BepiColombo radio science experiment MORE, with the results of the end-to-end simulation we have conducted for the Hermean phase of the mission. We provide the expected levels of uncertainties achievable in the positioning of the MPO, which is crucial for assessing the impact of the mission on planetary ephemerides construction. In the process, we have analysed the application of a new tracking configuration for the reduction of Doppler noise, TDMC, which could provide a potential amelioration of a factor two under favourable conditions.

To conclude, we have shown how planetary ephemerides can benefit from radiometric measurements of modern interplanetary missions, in particular from those with a specific devotion to radio science, featuring cutting-edge radio systems, such as Cassini and BepiColombo.

MGS Mars Global Surveyor

MMO

Figure 1 . 1 .

 11 Figure 1.1. Examples of zonal (l = 6, m = 0) (a), sectorial (l = m = 7) (b), and tesseral harmonic function (l = 13, m = 7) (c). The white areas represent the zones with positive polynomial values, those black with negative values.

  1.1 (b)). Harmonic functions with m = 0 and m = l are called tesseral functions; they divide the sphere into 2m(l -m) different regions, in m(l -m) of which they assume a positive value (see Fig. 1.1 (c)).

Figure 1 . 2 .

 12 Figure 1.2. The ξηζ body-fixed coordinate system in which extended body gravity accelerations are computed. Longitude (λ) and latitude (φ) angles are shown.

Fig. 3 .

 3 Fig.3. These curves show the drift of SSB in INPOP. We consider here a simplified model, composed of the Sun, the planets from Mercury to Pluto and the Moon (no asteroids), all considered as point-mass bodies. At the time origin of integration (J2000), the origin O of the reference frame is at the barycenter G of the system (Eqs.[START_REF] Moyer | Mathematical formulation of the Double-Precision Orbit Determination Program (DPODP)[END_REF] 4, 5)). The curves show the drift of the barycenter G in the reference frame, that is OG = µ *

Figure 1 . 3 .

 13 Figure 1.3. The curves show the drift of SSB in INPOP along the three axes. Plot from Fienga et al. (2008).

Figure 1 . 4 .

 14 Figure 1.4. Estimated masses for the asteroids orbit from INPOP17a and DE430. The x-axis provides the log of masses extracted from Carry (2012), the y-axis the INPOP17a and DE430 estimated masses. Plot from Viswanathan et al. (2017).

20 2 .

 2 INPOP19a: updates on the outer Solar System from Cassini

Figure 2 . 1 .

 21 Figure 2.1. Sketch summarising the OD process flow; the residuals checking entails an

Figure 2 . 2 .

 22 Figure 2.2. Scheme of two-way (a) and three-way (b) tracking configurations.

Figure 2 . 3 .

 23 Figure 2.3. Doppler extraction procedure. The difference between the transmitted and the received carrier frequencies yields the Doppler tone. A cycle counter measures phase change of the Doppler tone, yielding a measure of range change during the count interval[START_REF] Thornton | Radiometric tracking techniques for deep-space navigation[END_REF].

Figure 2 . 4 .

 24 Figure 2.4. Schematic depicting the error contributions producing biases on range measurements. The erroneous positioning of the planet from ephemerides predictions, together with the errors in the calibration of station delays, are the dominant contributions.

Figure 2 . 5 .

 25 Figure 2.5. Diagrams of the Cassini spacecraft and Huygens probe on the Titan IV launch adapter, highlighting the various components and instruments. The mission features a total of 18 experiments, 6 of which on board Huygens probe. [NASA Jet Propulsion Laboratory]

Figure 2 . 6 .

 26 Figure 2.6. Residuals of two-way Doppler data at 60 s count-time (top) and range data at 300 s (bottom) for a typical arc. The different colours and markers indicate the DSN complex and the specific DSS used for each tracking pass. The red vertical lines indicate the two flybys of Titan performed within the analysed period. In the bottom left boxes, the amount, mean, and SD of residuals are reported.

Figure 2 . 7 .

 27 Figure 2.7. Main accelerations (bottom panel) acting on Cassini during a typical lowaltitude, in particular T50, with the CA (red vertical line) occurred at a distance of 960 km from the surface (top panel).

Figure 2 . 8 .Figure 2 . 9 .

 2829 Figure 2.8. Residuals of two-way Doppler (top) and range (bottom) data for the three analysed periods. The different colours and markers indicate the DSN complex and the specific DSS used for each tracking pass. The red vertical lines indicate the targeted moons flybys.

44 2 .

 2 INPOP19a: updates on the outer Solar System from Cassini

Figure 2 .

 2 Figure 2.10. Kuiper belt model from Pitjeva and Pitjev (2018). The green dots show the distribution of the KBOs forming the three rings at 39.4, 44.0 and 47.5 AU. The central ring is more populated since it contributes for the two-thirds of the total mass, with respect to the one-sixth of the other two.

Figure 2 .

 2 Figure 2.11. Saturn post-fit residuals for the four datasets available for the Cassini mission (see Table 2.5) obtained with different INPOP solutions. In particular, the plots show the residuals obtained with INPOP17a, where no modelling of TNOs is included (top left); a solution accounting for the individual perturbations of the 9 binary TNOs (top right); a solution including the TNO ring but none of the binary TNOs (bottom left); and INPOP19a, which includes both the ring model and the binary TNOs (bottom right).

52 2 .

 2 INPOP19a: updates on the outer Solar System from Cassini

Figure 2 .

 2 Figure 2.12. Alternative model for Kuiper belt. The green dots show the distribution of the objects derived by sampling the real orbits of known KBOs.

Figure 2 . 13 .

 213 Figure 2.13. Histograms of the KBOs distributions for INPOP19a reference ring model (a) and the alternative model (b) for the Kuiper belt.

Figure 2 . 14 .

 214 Figure 2.14. Sampling distribution in RA and Dec for the P9 location analysis. Each point represents a test position for the planet, with its related gravity perturbation on the planetary ephemerides computed using Eq. (2.46); e.g. the highlighted point represents the i-th ephemerides solution perturbed by a P9 located in 200°RA and 20°D ec. The red line indicates the ecliptic plane.

Figure 2 .

 2 Figure 2.15. Likelihood maps for P9-perturbed solutions considering r P9 =400, 500, 600, 650, 700, 750 and 800 AU and M P9 = 5 and 10 M . Each RA-Dec point corresponds to a solution perturbed by P9 placed in that position. The colourbar provides the likelihood of the ephemerides including P9 perturbations with respect to INPOP19a.The white portions of the maps correspond to the rejected solutions, with a threshold at 3σ (L < 0.003)[START_REF] Fienga | New constraints on the location of p9 obtained with the inpop19a planetary ephemeris[END_REF].

Figure 2 .

 2 Figure 2.16. Compatibility based on the Mahalanobis distance criterium for Saturn orbit, considering r P9 =400, 500, 600, 650, 700, 750, 800 AU and M P9 =5 and 10 M . Each RA-Dec point corresponds to a solution perturbed by P9 placed in that position. The colourbar gives the compatibility percentage between the propagated orbit of Saturn from INPOP19a and the P9-perturbed ephemerides[START_REF] Fienga | New constraints on the location of p9 obtained with the inpop19a planetary ephemeris[END_REF].

Figure 2 . 17 .

 217 Figure 2.17. LIkelihood test for the r P9 =600 AU and M P9 =5 M configuration, with the superimposition of the Mahalnobis distance compatibility results depicted in magenta lines (a). The lower panels (b) and (c) provide a magnification of the two acceptance zones compatible for both criteria, D1 and D2, respectively[START_REF] Fienga | New constraints on the location of p9 obtained with the inpop19a planetary ephemeris[END_REF].

Figure 3 . 1 .

 31 Figure 3.1. Elements of the BepiColombo Mercury Composite Spacecraft. From left to right: MTM, MPO, MOSIF, and MMO. [ESA/ATG medialab copyright]

Figure 3 . 2 .

 32 Figure 3.2. Image taken by one of the MCAM "selfie" cameras on board BepiColombo, as the spacecraft zoomed past the planet during its Earth flyby. [ESA/BepiColombo/MTM, CC BY-SA 3.0 IGO]

66 3 .

 3 Simulations of BepiColombo mission

Figure 3 . 3 .

 33 Figure 3.3. Comparison of Earth and Mercury relatives sizes. The picture of Mercury is one of the first from the MESSENGER probe, shot in Oct. 6, 2008 during a flyby of the planet. The Earth is a composition of images taken by NASA satellite Terra.[NASA/JHUAPL/CIW]

70 3 .

 3 Simulations of BepiColombo mission

•

  Mercury's obliquity (43m on the surface) 1σ < 3.7 arcsec • Physical librations amplitude in longitude 1σ ≈ 3.2 arcsec• C m /C a 1σ < 0.05 • C/M R 2 b 1σ < 0.003Other goals in geodesy• Topography of the planet to the accuracy of the 1 m laser altimeter (BELA) • Planetary figure: mean, polar and equatorial radius 1 part in 10 7 (combination of MORE and BELA data) Fundamental physics goals • γ, controlling the deflection of light and the relative time delay 2.5 × 10 -6 • β, controlling the relativistic advance of Mercury's perihelion 5 × 10 -6 • η, controlling the gravitational self-energy contribution 2 × 10 -5 • Time variation of the gravitational constant, d(ln G) dt 3 × 10 -13 years -1 • Sun's oblateness, J 2 2 × 10 -9 Orbit determination goals • Spacecraft position in a Mercurycentric frame 1σ < 1 m • Mercury position in the heliocentric frame 1σ < 10 m Additional scientific goals • Measure solar wind velocity in the acceleration region by correlating uplink and downlink plasma contributions • Measure solar wind velocity in the corona by exploiting different refraction of X and Ka signals • Measure the scintillation index for X and Ka band as function of SEP angle • Characterise and compare performances for the radio links X/Ka and Ka/Ka, especially at small impact parameters (1-5 solar radii) a Ratio between mantle and planet inertia b Condensation coefficient

Figure 3 . 4 .

 34 Figure 3.4. MPO reconstructed position formal uncertainties in the Hermean RTN-frame at the centre of each arc (a). In (b) are reported the propagated uncertainties for a typical one-day arc. The light blue and light red zones depict the X-band (navigation) and Ka-band (scientific) tracking periods, respectively. The black dashed lines indicate the two reaction wheels desaturation manoeuvres expected on each arc[START_REF] Imperi | An analysis of the geodesy and relativity experiments of bepicolombo[END_REF][START_REF] Thor | Prospects for measuring mercury´s tidal love number h2 with the bepicolombo laser altimeter[END_REF].

Figure 3 . 5 .

 35 Figure 3.5. Sketch depicting two-way and three-way link configurations for BepiColombo. Local noise contributions arising at transmitting (t -T ) and receiving (t) times are shown for the two stations (Di Ruscio et al., 2019).

Figure 3 .

 3 Figure 3. MPO's position formal accuracies in the Mercurycentered Radial-Transverse-Normal frame, obtained with the nominal two-way tracking configuration and with the use of TDMC. The blue lines on the right show the respective gain factors registered on the three components.

Figure 3 . 6 .

 36 Figure 3.6. MPO position formal uncertainties in RTN-frame obtained with the baseline two-way configuration from Malargüe and with the application of TDMC using LLAMA as three-way antenna. The right-side plots report the gain factor expressed as σ 2way /σ TDMC ratio. Plot from Di Ruscio et al. (2019).

  Mercury Magnetospheric Orbiter MO Mars Odyssey MONTE Mission Analysis, Operations, and Navigation Toolkit Environment MORE Mercury Orbiter Radio Science Experiment MOSIF Magnetospheric Orbiter Sunshield and Interface Structure MPO Mercury Planetary Orbiter MRO Mars Reconnaissance Orbiter MTM Mercury Transfer Module NASA National Aeronautics and Space Administration OD orbit determination OTM orbital trim manoeuvre P9 Planet 9 PLL phase-locked loop PN pseudo-noise PPN parametrised post-Newtonian RA right ascension RCS reaction control system RMS root mean square RTG radioisotope thermoelectric generator RTLT round-trip light-time RTN radial transverse normal RU range unit SCE Solar Conjunction Experiment SD standard deviation

  Maximum differences between DE405 and INPOP06 in the ICRF frame. The comparisons are made over three time intervals centred at the initial conditions, J2000: ±30 yr, ±100 yr, and from -400 to +200 yr (the complete interval of DE405). Tablefrom[START_REF] Fienga | Inpop06: a new numerical planetary ephemeris[END_REF].

	1.3 INPOP planetary ephemerides		13
	A. Fienga et al.: INPOP06. A new numerical planetary ephemeris	317
	ongitude in INPOP06 (step size
	for the Earth-Moon Barycenter,
	given in microarcseconds (µas)
	d over 10 000 years.
	10 000 yr
	ε i (µas) 1.5 × 10 -1 41 288 ε i (µm) 4901 9.3 × 10 -3 5335 7.3 × 10 -3 461 4.2 × 10 -4 34 9.0 × 10 -6 37 5.3 × 10 -6 23 1.7 × 10 -6 29 1.3 × 10 -6 38 1.3 × 10 -6 1.3 2 5 1 3	Planet or body Mercury Venus EMB Mars Jupiter	Differences over the time interval: ±30 yr ±100 yr -400, +200 yr Heliocentric frame (m) 6.2 26 0.5 26 1.6 7.7 58 517 3 31 108
		Saturn	20	36
		Uranus	38	45
		Neptune	35	77
		Pluto	31	119
		Moon *	0.009	0.034	0.540
			Longitude λ (mas)
		Mercury	174	749	6
		Venus	78	534	6
		EMB	32	286	2
	in the longitude of the Moon. f of the difference after one way	Mars Jupiter	354 57	6 600 276	42 1
	line is obtained by least square	Saturn	16	51
		Uranus	4	14
		Neptune	0.7	14
	Few elements differ between	Pluto	0.9	7
		Moon *	60	170	6
	center		Latitude φ (mas)
	Sun is not integrated in the and velocity are determined ids, assuming that the Solar	Mercury Venus	65 30	268 207	2 2
	t the origin of the inertial ref-	EMB	13	113
	c position, u i the barycentric ody i, it is assumed that (see	Mars Jupiter	152 23	2 602 102	16
	(1)	Saturn Uranus	5 0.6	17 5
		Neptune	0.2	5
	,	Pluto	0.4	3
		Moon *	25	73	2
	(2)

1. Planetary ephemerides

Table 1.1. * Geocentric frame

Table 1 .3. List

 1 of observables included in INPOP17a dataset. The statistics (mean and SD) of INPOP17a and INPOP13c residuals are reported as well[START_REF] Viswanathan | Inpop17a planetary ephemerides[END_REF].

	Type	n o	Period	Mean	SD
				17a	13c	17a	13c
			Mercury
	Direct range (m)	462	1971-Apr/1997-Aug -167.975 -101.524 1 605.318	861.494
	Mariner range (m)	2	1974-Mar/1976-Mar	-75.645 -196.405	78.743	19.636
	MESSENGER range (m)	269	2011-Mar/2013-Feb	1.250	4.008	8.603	12.387
			Venus
	Direct range (m)	489	1965-Dec/1990-Jan	497.493	504.569 2 236.789 2 237.636
	VLBI (mas)	46	1990-Sep/2013-Feb	1.195	1.591	2.356	2.575
	VEX range (m)	24 249	2006-Apr/2011-Jun	2.762	2.362	6.576	6.693
			Mars
	Viking range (m)	1 257	1976-Jun/1982-Nov	-8.203	-1.494	60.833	41.189
	VLBI (mas)	194	1989-Feb/2013-Nov	0.124	0.116	0.348	0.333
	Pathfinder ran-ge (m)	90	1997-Jun/Sep	-0.233	19.324	13.268	14.096
	MGS range (m)	2 417	1999-Apr/2006-Sep	4.370	3.890	0.362	3.777
	MO range (m)	21 673	2002-Jan/2014-Jan	2.155	2.470	5.797	7.200
	MEX range (m)	29 203 2005-Mar/2016-May	-0.059	4.669	3.552	23.361
			Jupiter
	Optical RA/ Dec (as)	6 394	1914-May/2008-Jun	-0.005 -0.047	-0.039 -0.048	0.308 0.302	0.297 0.301
	Flybys RA/ Dec (mas)	5	1974-Dec/2001-Jan	2.396 -9.148	2.554 -10.853	2.351 9.895	2.961 11.425
	Flybys range (m)	5	1974-Dec/2001-Jan -644.521 -985.957	1414.411	1775.627
	VLBI (mas)	24	1996-Jul/1997-Dec	-0.702	-0.450	11.376	11.069
			Saturn
	Optical RA/ Dec (as)	7 945 1924-Mar/2008-May	0.023 -0.012	-0.006 -0.012	0.304 0.266	0.293 0.266
	Cassini:
	• range (m)	169 2004-May/2014-May	5.059 -471.270	31.618	340.340
	• VLBI RA/ • Dec (mas)	10	2004-Sep/2009-Apr	0.172 -0.203	0.113 -0.115	0.553 0.272	0.630 0.331
			Uranus
	Optical RA/ Dec (as)	13 008	1914-Jun/2011-Sep	-0.003 -0.027	0.007 -0.006	0.216 0.234	0.205 0.234
	Flybys RA/ Dec (mas)	1	1986-Jan	-0.021 -0.064	-0.021 -0.028	--	--
	Flybys range (m)	1	1986-Jan	-0.080	20.771	-	-
					[Continued on next page]

1. Planetary ephemerides Table 1.3. -

  Continued from previous page.

	Type	n o	Period	Mean		SD
				17a	13c	17a	13c
			Neptune
	Optical RA/ Dec (as) Flybys RA/ Dec (mas) Chapter 2	5 375 1	1913-Dec/2007-Nov 1989-Aug	0.008 -0.011 -0.023 -0.015	0.003 -0.002 -0.005 -0.011	0.259 0.303 --	0.258 0.299 --
	Flybys range (m)	1	1989-Aug	-2.869	51.507	-	-
			Pluto
	Optical RA/ Dec (as)	3 052	1914-Jan/2008-Jun	-0.035 0.010	0.020 0.001	0.566 0.479	0.574 0.525
	HST RA/ Dec (as)	5	1998-Feb	-0.052 -0.001	-0.018 -0.026	0.044 0.048	0.044 0.048
	Occultations RA/ Dec (as)	13	2005-Jun/2009-Sep	-0.009 0.008	-0.100 0.000	0.045 0.026	0.044 0.027

Table 2 .

 2

	S	2 110 -2 120	2 290 -2 300
	X	7 145 -7 190	8 400 -8 450
	Ka	34 200 -34 700	31 800 -32 300

1. International frequencies adopted for deep-space communications.

Band Uplink frequency (MHz) Downlink frequency (Mhz)

Table 2 . 2 .

 22 Summary of gravity dedicated flybys of Titan reporting: the CA date, altitude, latitude and longitude, Titan's mean anomaly M and SEP angle.

	Flyby	Epoch	Alt. (km) Lat. (°) Long. (°) M (°) SEP (°)
	T011	27 Feb 2006	1 812	-0.2	255.6	173	150.2
	T022	28 Dec 2006	1 297	45.4	355.9	197	130.6
	T033	29 Jun 2007	1 933	8.4	63.1	15	46.4
	T045	31 Jul 2008	1 614	-43.5	162.7	346	30.0
	T068	20 May 2010	1 397	-48.9	241.1	82	120.3
	T074	18 Feb 2011	3 651	1.0	113.4	159	131.2
	T089	17 Feb 2013	1 978	21.0	203.1	66	106.0
	T099	6 Mar 2014	1 500	-31.1	181.0	64	111.0
	T110 a 16 Mar 2015	2 274	74.8	263.1	248	108.6
	T122	10 Aug 2016	1 698	12.4	234.4	304	112.7

a flyby dedicated to VIMS, tracking data collected with LGA four and six passes. The extended dataset, in fact, allows the estimation of the coefficients up to degree 5 of the gravitational potential spherical harmonic expansion. The estimated values and the associated formal uncertainties, deduced by

[START_REF] Durante | Titan's gravity field and interior structure after Cassini[END_REF]

by processing the complete dataset of Titan gravity flybys (see Table

2

.2) in a multi-arc fit (see

Sect. 2.1.2)

, are reported in

Table 2.3.

Table 2 . 3 .

 23 Titan's gravity field solution from[START_REF] Durante | Titan's gravity field and interior structure after Cassini[END_REF] including the gravitational parameter GM , the unnormalised potential coefficients with a reference radius of 2575 km, and the tidal Love number k 2 , with their respective 1σ uncertainties.

	Parameter	Value	Parameter	Value
	GM	8978.1383 ± 0.0003 km 3 /s 2
	Unnormalized harmonic coefficients (×10 6 ):
	J 2	33.089 ± 0.609
	C 21	0.513 ± 0.215	S 21	0.612 ± 0.359
	C 22	10.385 ± 0.084	S 22	-0.064 ± 0.066
	J 3	-0.179 ± 0.720
	C 31	1.481 ± 0.254	S 31	0.811 ± 0.402
	C 32	0.183 ± 0.153	S 32	-0.027 ± 0.099
	C 33	-0.222 ± 0.017	S 33	-0.226 ± 0.019
	J 4	-1.077 ± 1.844
	C 41	-0.842 ± 0.299	S 41	0.191 ± 0.717
	C 42	0.183 ± 0.107	S 42	0.198 ± 0.106
	C 43	-0.012 ± 0.039	S 43	-0.062 ± 0.033
	C 44	-0.014 ± 0.003	S 44	-0.012 ± 0.004
	J 5	1.118 ± 2.022
	C 51	0.361 ± 0.406	S 51	0.267 ± 0.604
	C 52	-0.097 ± 0.118	S 52	0.044 ± 0.094
	C 53	-0.016 ± 0.019	S 53	-0.004 ± 0.012
	C 54	0.007 ± 0.004	S 54	-0.002 ± 0.004
	C 55	0.000 ± 0.001	S 55	0.000 ± 0.001
	k 2	0.616 ± 0.067

Table 2 .

 2 4. Saturn's gravity field solution and associated uncertainties at 1σ obtained by[START_REF] Iess | Measurement and implications of saturn's gravity field and ring mass[END_REF] from Grand Finale pericentres. The harmonic coefficients reported are unnormalized, and relative to a reference radius of 60 330 km.

	Coefficients (×10 6 )	Value
	J 2	16 290.573 ± 0.028
	J 3	0.059 ± 0.023
	J 4	-935.314 ± 0.037
	J 5	-0.224 ± 0.054
	J 6	86.340 ± 0.087
	J 7	0.108 ± 0.122
	J 8	-14.624 ± 0.205
	J 9	0.369 ± 0.260
	J 10	4.672 ± 0.420
	J 11	-0.317 ± 0.458
	J 12	-0.997 ± 0.672

Table 2 . 5 .

 25 List of datasets included in INPOP19a fit. Columns 1 and 2 provide information on the mission, type, and number of observations. Column 3 gives the time interval covered, while column 4 lists the related accuracies provided by space agencies or navigation teams. Finally, in the last two columns the WRMS for INPOP19a and IN-POP17a are reported. The measurement sets produced within this work are highlighted in grey.

	Type	n o	Period	Average	WRMS
				accuracy	19a	17a
			Mercury
	Direct range (m)	462	1971-Apr/1997-Aug	900.0	0.95	0.96
	Mariner range (m)	2	1974-Mar/1976-Mar	100.0	0.37	0.78
	MESSENGER range (m)	1 096	2011-Mar/2014-Apr	5.0	0.82	1.29
			Venus
	Direct range (m)	489	1965-Dec/1990-Jan	1 400.0	0.98	0.98
	VLBI (mas)	68	1990-Sep/2013-Feb	2.0	1.13	1.178
	VEX range (m)	24 783	2006-Apr/2011-Jun	7.0	0.93	0.93
			Mars
	VLBI (mas)	194	1989-Feb/2013-Nov	0.3	1.26	1.16
	MGS range (m)	2 459	1999-Apr/2006-Sep	2.0	0.93	1.31
	MRO/MO range (m)	20 985	2002-Feb/2014-Jan	1.2	1.07	1.91
	MEX range (m):
	• INPOP17a interval	29 203 2005-Mar/2016-May	2.0	0.97	1.26
	• INPOP19a interval	30 669 2005-Mar/2017-May	2.0	0.98	3.37
			Jupiter
	Optical RA/Dec (as)	6 416	1924-May/2008-Jun	0.3	1.0	1.0
	Flybys RA/Dec (mas)	5	1974-Dec/2001-Jan	4.0/12.0	0.94/1.0	0.58/0.82
	Flybys range (m)	5	1974-Dec/2001-Jan	2 000.0	0.98	0.71
	VLBI (mas)	24	1996-Jul/1997-Dec	11.0	1.01	1.03
	Juno range (m)	9	2016-Aug/2018-Sep	20.0	0.945	116.0
			Saturn
	Optical RA/Dec (as)	7 826 1924-Mar/2008-May	0.3	0.96/0.87	0.96/0.87
	Cassini:
	• VLBI RA/Dec (mas)	10	2004-Sep/2009-Apr	0.6/0.3	0.97/0.99	0.92/0.91
	• JPL range (m)	165 2004-May/2014-May	25.0	0.99	1.01
	• Navigation + Titan gra-vity flybys range (m)	614	2006-Jan/2016-Aug	6.0	1.01	2.64
	• Grand Finale range (m)	9	2017-May/Jul	3.0	1.14	29.0
			Uranus
	Optical RA/Dec (as)	12 893	1924-Aug/2011-Sep	0.2/0.3	1.09/0.82	1.09/0.82
	Flybys RA/Dec (mas)	1	1986-Jan	50.0	0.12/0.42	0.42/1.23
	Flybys range (m)	1	1986-Jan	50.0	0.92	0.002
			Neptune
	Optical RA/Dec (as)	5 254	1924-Jan/2007-Nov	0.25/0.3 1.008/0.97 1.008/0.97
	Flybys RA/Dec (mas)	1	1989-Aug	15.0	0.11/0.15	1.0/1.57
	Flybys range (m)	1	1989-Aug	2.0	1.14	1.42

updates on the outer Solar System from CassiniTable 2 . 6 .

 26 List of the binary TNOs included in INPOP19a dynamical model, providing their IAU identification number and mass. Values extracted from[START_REF] Pitjeva | Mass of the kuiper belt[END_REF].

	TNO IAU GM × 10 14 mass ×10 -20
	number	(AU 3 d -2 )	(kg)
	136199	251.9160	169.3357
	136108	60.36990	40.58011
	136472	44.98510	30.23858
	90482	9.554830	6.422671
	208996	8.007340	5.382462
	50000	7.235460	4.863611
	174567	3.994670	2.685181
	120347	1.934696	1.300485
	55637	1.880380	1.263975

Table 2 . 7 .

 27 .11. The plot shows the post-fit residuals of the four datasets computed with different solutions of INPOP: the INPOP17a release (top left), where no model for the TNO is included, either individually or with a ring; INPOP19a (bottom right), which features both a TNO ring model and nine binary TNOs (see Sect. 2.5.1); an intermediate solution (top right) which includes Cassini data WRMS values for multiple INPOP solutions. The first columnindicates the specific Cassini dataset, while in the remaining columns, the WRMS for different INPOP fits are reported. TNOs, but without the TNO ring model; and another case, that does not consider the perturbations produced by the binary TNOs, but the effect of the ring. It is evident as the combined inclusion of the individual perturbation of nine binary TNOs, together with the adjustment of the mass of a TNO ring, significantly improves the post-fit residuals. This becomes particularly important if we consider an interval of time spread over several decades, such the one offered by INPOP19a.

	2.6 Results

Table 2 . 8 .

 28 Possible zones for the presence of P9 according to the Mahalanobis distance and likelihood compatibility criteria with respect to INPOP19a planetary ephemerides.

Table 3 .

 3 1.MPO and MMO main characteristics.

		MPO	MMO
	Stabilisation	3-axis stabilised	15-rpm spin stabilised
	Orientation	Nadir pointing	Spin axis at 90°to Sun
	Orbit	Polar orbit, period of 2.3 h 480×1 500 km	Polar orbit, period of 9.3 h 590×11 640 km
	Payload mass	80 kg	45 kg
	Total mass	1 150 kg (in Mercury orbit)	275 kg (in Mercury orbit)
	Communications link X/Ka-band	X-band
	Antenna	1 m high-temperature resistant 0.8 m phased array HGA steerable HGA
	[http://sci.esa.int/bepicolombo/47346-fact-sheet/]

Table 3 . 2 .

 32 Mercury's and Earth's main parameters.

	Parameter	Mercury	Earth	Ratio ( / )
	Mass (×10 24 kg)	0.33011	5.9724	18.0832
	Volume (×10 10 km 3 )	6.083	108.321	17.7936
	Mean radius (km)	2 439.7	6 371.0	2.6110
	Ellipticity (Flattening)	0.00000	0.00335	n/a
	Mean density (kg/m 3 )	5 427	5 514	1.016
	Surface Gravity (m/s 2 )	3.70	9.81	2.64
	GM (×10 5 km/s 2 )	0.22032	3.968	18.0918
	Solar irradiance (W/m 2 )	9 082.7	1 361.0	0.15
	J 2 (×10 6 )	50.3	1 082.63	21.52
	Semimajor axis (×10 6 km)	57.91	149.60	2.58
	Sidereal orbit period (days)	87.969	365.256	4.152
	Mean orbital velocity (km/s)	47.36	29.78	0.63
	Orbit inclination (deg)	7.00	0.00	n/a
	Orbit eccentricity	0.2056	0.0167	0.0812
	Rotation period (Earth-days)	59	1	0.0169
	[https://nssdc.gsfc.nasa.gov/planetary/factsheet/mercuryfact.html]

Table 3 . 3 .

 33 List of MORE's scientific goals

3.3 The radio science experiment MORE 71Table 3 .

 3 4. MORE's range-rate error budget(MORE Team, 2011).

	Source (two-way, @1000 s) Single-link value Plasma-free value
	Space segment		9.40 × 10 -15
	DST + TWTA a in X/X	6.00 × 10 -15
	DST + TWTA in X/Ka	1.00 × 10 -14	1.14 × 10 -15
	KaT in Ka/Ka	1.00 × 10 -15
	RFDA b	1.00 × 10 -15	1.11 × 10 -15
	HGA assembly	7.33 × 10 -15	8.14 × 10 -15
	MPO structure	1.00 × 10 -15	1.11 × 10 -15
	Ground segment		6.37 × 10 -15
	FTS	2.00E-15	2.22E-15
	G/S electronics X/X (tx and rx chains) G/S electronics in X/Ka	1.53 × 10 -15 1.53 × 10 -15	1.54 × 10 -15
	G/S electronics in Ka/Ka	1.53 × 10 -15
	Ground antenna	5.00 × 10 -15	5.55 × 10 -15
	Station location	1.20 × 10 -15	1.33 × 10 -15
	EOP c	6.00 × 10 -16	6.66 × 10 -16
	Earth solid tides	5.00 × 10 -16	5.55 × 10 -16
	Media		7.79 × 10 -15
	Residual plasma	5.00 × 10 -16	5.55 × 10 -16
	Residual troposphere	7.00 × 10 -15	7.77 × 10 -15
	End-to-end		1.38 × 10 -14

a Traveling Wave Tube Amplifier b Radio Frequency Distribution Assembly c Earth Orientation Parameters

Table 3 .5.

 3 Results of the covariance analysis of BepiColombo MORE gravity experiment for the reference two-way and the TDMC-enhanced solutions. In particular the formal uncertainties of Mercury's GM , J 2 , tidal Love number k 2 and rotational state (pole direction in RA and Dec and libration) are reported[START_REF] Di Ruscio | Improvements in BepiColombo and JUICE radio science experiments with a multi-station tracking configuration for the reduction of Doppler noise[END_REF].

	GM	J 2	k 2	RA	Dec	Lib.
	(km 3 /s 2 )			(arcsec) (arcsec) (arcsec)
			Two-way (DSA-3):
	2.1 × 10 -4 7.4 × 10 -11 5.2 × 10 -4	0.53	0.18	0.35
			TDMC (DSA-3 / LLAMA):
	3.3 × 10 -4 4.3 × 10 -11 5.2 × 10 -4	0.29	0.11	0.20

Planetary ephemerides

1.3 INPOP planetary ephemerides

best, i.e. optimal in statistical sense.

pseudo because in the general Lambert's problem two positions are given, while in OD the initial conditions comprise the spacecraft initial position and speed.

INPOP19a: updates on the outer Solar System from Cassini

5 712 kg with fuel, Huygens probe and adapter; 2 125 kg unfueled orbiter alone.

A second test (SCE2), expected in 2003, was canceled by the Cassini Project after the success of the first one, due to technical problems with the reactions wheels.

The flyby labelling provides for the initial of the interested moon (e.g. T for Titan, E for Enceladus) followed by the encounter number.

Z-axis corresponds to the HGA pointing axis (see Fig.2.5).

Data from a sixth pass dedicated to gravity (Rev 275) were lost for configuration problems at the ground station.

SMART-1, acronym for Small Missions for Advance Research and Technology-1, was an ESA mission specifically developed to test solar-electric propulsion and other deep-space technologies, while performing scientific observations of the Moon.

Simulations of BepiColombo mission

Inter-arc constraints on the spacecraft state vectors are given to the level of 1 m.

Acknowledgments

I would like to express my gratitude to my supervisors, Prof. Luciano Iess and Prof. Agnes Fienga, which have patiently followed and mentored me throughout the PhD, providing me invaluable opportunities and precious advice. I would love to thank my colleagues of the Radio Science laboratory, which, beside offering me fruitful suggestions for my work, have livened up my days in the laboratory with their friendship; they have also proven to be great travel companions and "exceptionally experienced canyoneers". I would like to thank Sara and my family as well, which have supported me during these years.

Abbreviations and Acronyms