How Can I Interpolate Ephemerides?

Question

I have two ephemerides for orbital elements—that is, each contains $a$ (semi-major axis), $e$ (ecc.), $i$ (incl.), $\Omega$ (long. asc. node), $\omega$ (arg. of peri.), $\tau$ (time of peri.), $T$ (period), and a timestamp $t$.

I want to compute a new ephemeris at a point somewhere between these two ephemerides. For example, the timestamp ($t$) of the new ephemeris can be calculated by a linear interpolation of the start ($t_0$) and end ($t_1$) ephemerides' times as:$$ t = \text{lerp}(t_0,t_1, x) = t_0 \cdot (1-x) + t_1 \cdot x $$

Unfortunately, I can't seem to get an interpolation for the other parameters that doesn't have discontinuities in it. My question: how should I be interpolating ephemerides? Maybe I need to do this through some intermediate parameters, like $\varpi$, $\lambda$, or $\nu$?

---

_{^{For example, if I take $\tau$ from the first ephemeris and lerp everything else, I get a strong discontinuity when mean anomaly is calculated (as $M:=n (t-\tau)=(2 \pi/T) (t - \tau)$). I suspect this is because the mean anomaly calculation is only valid at the ephemeris, and as e.g. $T$ changes due to interpolation, $M$ changes erroneously.}}

_{^{Here is an example for Mercury. The first three ephemerides come from NASA HORIZONS and are correct. Interpolating the first two ephemerides produces $M=328^\circ$. However, two seconds later, the interpolation (which is now interpolating from the last two ephemerides) produces $M=294^\circ$.}}

Answer 1

OK, I see what you did. You looked at the osculating elements for Mercury's orbit around the solar system barycenter, by using Horizons with these parameters:

As you point out in the attached pastebin, the values jump around a lot:

$$SOE
2451544.500000000 = A.D. 2000-Jan-01 00:00:00.0000 TDB 
 EC= 1.970095254610514E-01 QR= 3.160869485749435E-01 IN= 7.013139646115477E+00
 OM= 4.812103959293303E+01 W = 2.595697040970977E+01 Tp=  2451500.244608603884
 N = 3.993474279540898E+00 MA= 1.767327672710030E+02 TA= 1.777641486117292E+02
 A = 3.936372330648511E-01 AD= 4.711875175547588E-01 PR= 9.014706864254217E+01
2451575.500000000 = A.D. 2000-Feb-01 00:00:00.0000 TDB 
 EC= 2.099217560035118E-01 QR= 2.935320561427645E-01 IN= 7.027294624675849E+00
 OM= 4.850357410323192E+01 W = 3.527448059792021E+01 Tp=  2451590.708947353531
 N = 4.355289183393867E+00 MA= 2.937606361003802E+02 TA= 2.698751266828741E+02
 A = 3.715227679957090E-01 AD= 4.495134798486537E-01 PR= 8.265811633648384E+01
2451604.500000000 = A.D. 2000-Mar-01 00:00:00.0000 TDB 
 EC= 2.241380367466583E-01 QR= 3.093951924304613E-01 IN= 6.978935239526586E+00
 OM= 4.831291617128525E+01 W = 3.114644418385429E+01 Tp=  2451590.362934530713
 N = 3.916530358018232E+00 MA= 5.536824608413040E+01 TA= 7.970140671677829E+01
 A = 3.987760801330001E-01 AD= 4.881569678355390E-01 PR= 9.191809256960804E+01
$$EOE

so interpolation doesn't work well.

The problem: Mercury does NOT orbit the solar system barycenter; it orbits the Sun (What point does Earth actually orbit?)

The barycenter is the center of mass of our solar system. If you're far enough away from our solar system, you can treat the solar system as a point mass at that point.

Inside our solar system, that doesn't work as well for most planets (massive Jupiter is an exception). Example: when Sun-Mercury-Jupiter form a straight line (and Jupiter's on the same side of the Sun as Mercury), Jupiter's gravitational on influence on Mercury is stronger than it is on the Sun.

If you run Horizons using the Sun as the center body:

you'll get:

$$SOE
2451544.500000000 = A.D. 2000-Jan-01 00:00:00.0000 TDB 
 EC= 2.056302512089075E-01 QR= 3.074991199665784E-01 IN= 7.005014199657344E+00
 OM= 4.833053756455964E+01 W = 2.912428058698772E+01 Tp=  2451502.287118767854
 N = 4.092345945977128E+00 MA= 1.727497133778637E+02 TA= 1.751155303115542E+02
 A = 3.870982252717257E-01 AD= 4.666973305768729E-01 PR= 8.796910250314700E+01
2451575.500000000 = A.D. 2000-Feb-01 00:00:00.0000 TDB 
 EC= 2.056312618551657E-01 QR= 3.074988256430905E-01 IN= 7.005012539288613E+00
 OM= 4.833050698201237E+01 W = 2.912403060593708E+01 Tp=  2451590.256131388247
 N = 4.092344011669980E+00 MA= 2.996128340775531E+02 TA= 2.767940962437320E+02
 A = 3.870983472501978E-01 AD= 4.666978688573051E-01 PR= 8.796914408304917E+01
2451604.500000000 = A.D. 2000-Mar-01 00:00:00.0000 TDB 
 EC= 2.056306629761345E-01 QR= 3.074990706849075E-01 IN= 7.005011440208267E+00
 OM= 4.833046370516735E+01 W = 2.912368651066544E+01 Tp=  2451590.256076942664
 N = 4.092343747828102E+00 MA= 5.829102946757737E+01 TA= 8.083993808120411E+01
 A = 3.870983638882191E-01 AD= 4.666976570915307E-01 PR= 8.796914975460210E+01
$$EOE

Notice that most of the parameters change only infinitesimally, since Mercury's orbit around the Sun is fairly stable (for a counterexample, see the Moon's orbit around the Earth, which varies quite a bit).

As a note, I also thought it would be nice if the planets orbited the barycenter and did some work on it, but the numbers just don't work out: https://github.com/barrycarter/bcapps/blob/master/ASTRO/playground4.m