Incorrect covariance interpolation from RTN frame

Orekit Python Wrapper usage
1

I was recently implementing covariance interpolation from an OEM, referencing the solution in Covariance Matrix Interpolation from OEM . The end goal of this is getting Pc parity between a SpaceTrack CDM (Ephem/Ephem) and the actual ephems that were used in that CDM.

The issue I am seeing, is that the interpolated covariance (when the input OEM covariance is in RTN) is coming out different than the CDM and also does not match with the non-interpolated values on either side of the TCA from the ephem. This seems to be only for the Radial component. My process is as follows:

I loaded in the two OEM ephems (OEM Parser → propagator), and used an ExtremumApproachDetector to get the TCA and miss distance/velocity - this matched with the CDM. I then set up a covariance interpolator similar to the link above and use it to build an Ephemeris() object. Then finally extract the interpolated covariance through the .getCovariance() method. Code as follows:

# Covariance interpolator
covariance_kep_hermite = StateCovarianceKeplerianHermiteInterpolator(
    OrbitHermiteInterpolator(<EME2000>),
    LOFType.QSW
)

# New propagator
ephem_propagator = Ephemeris(
    <ArrayList(SpacecraftState)>,  # Tabulated states from the OEM
    SpacecraftStateInterpolator(5, <EME2000 Frame>),  # States are in J2000
    <ArrayList(StateCovariance)>,  # RTN covariance built with LOFType.QSW, while J2000 covariance build as such with <EME2000>/OrbitType/PositionAngleType
    covariance_kep_hermite
)

# Extract covariance at TCA
tca_java = AbsoluteDate('2026-01-23T00:58:50.750', TimeScalesFactory.getUTC())
interp_cov = ephem_propagator.getCovariance(tca_java)

# Convert interpolated covariance back into RTN/RIC to compare against CDM & OEM
interp_orbit = ephem_propagator.propagateOrbit(tca_java)
interp_cov_ric = interp_cov.changeCovarianceFrame(interp_orbit, LOFType.QSW)

This is where things get weird. My first ephem (testsat-1.oem) has states and covariance both in J2000. The interpolated covariance here matches almost exactly with the CDM. However my second ephem (testsat-2.oem) had states in J2000 and covariance in RTN. Here I can see the diagonal RIC covariance values on either side of the interpolated point and the Radial value is nowhere near the interpolated one. Note that both ephems have 1 minute steps. Here are the results:

# interp_cov_ric is from the previous block. cdm_cov is the CR_R, CT_T, CN_N values. Units are the covariance units of m^2

# testsat-1.oem -> states & covar in J2000
interp_cov_ric   = [3999.9587648430606, 852573168.0795757, 2012.599168118526]
cdm_cov          = [3999.885341602069,  852565808.4392604, 2007.709022392526]

# testsat-2.oem -> states in J2000, covar in RTN
interp_cov_ric   = [13007.472576136963, 13170728.940087868, 27.1341201964874]
cdm_cov          = [6224.041894496001,  13170283.90046027,  26.96531832610106]

# testsat-2.oem RTN diagonal covar at step BEFORE interpolated point (...00:58:42):
oem_cov_rtn_pre  = [6231.1060419000005, 13161702.127,       26.870291514]
# testsat-2.oem RTN diagonal covar at step AFTER interpolated point (...00:59:42):
oem_cov_rtn_post = [6188.1730422,       13220598.384,       27.384216506999998]

Is there something wrong in my implementation that is approximately doubling the Radial covariance? The interpolator should be using the two points (as default) on either side of the interpolated one, but its value is nowhere near.

Attaching abridged ephems for reference (they produce the same results)

testsat-1.oem (7.0 KB)

testsat-2.oem (6.1 KB)

2

Hi,

You could try to replace the StateCovarianceKeplerianHermiteInterpolator with a StateCovarianceBlender. The former is known to artificially inflate the interpolated covariance, and also to produce non-positive definite matrices.

See this note for more context.

3

Hello @aerow610,

Thank you for reporting this, it does seem strange that only the radial value is wrong. I’m investigating the issue. This might be an issue common to StateCovarianceKeplerianHermiteInterpolator and StateCovarianceBlender.

Cheers,
Vincent

4

I spent more time struggling against WSL than investigating but i got your answer :slight_smile: !

You are most probably loading your testsat-2 covariances with LOFType.QSW which produce your results.

Use LOFType.QSW_INERTIAL and this shall fix your problem.

Cheers,
Vincent

5

Hi @Vincent, thanks for looking into this - the QSW_INERTIAL did the trick!

Note that it worked only when originally building the StateCovariance() from the OEM tabulated values. If I left that as QSW and put QSW_INERTIAL in the StateCovarianceBlender instead, it did not work.

I am curious though, why did that make the difference?

I also tried out a couple other things yesterday which had the same issues in radial.

  • Usage of the StateCovarianceBlender produced the same results as StateCovarianceKeplerianHermiteInterpolator, but good to note that the former should be used going forward.
  • In the initial building of the ArrayList(StateCovariance), I created an intermediate covariance with LOFType.QSW, and then used the .changeCovarianceFrame() method to convert it to J2000. Interestingly was getting the same (wrong) results

Final question going forward: if I have this type of ephem with RTN covariance, should I build everything covariance-related using the LOFType.QSW_INERTIAL, or build it with <EME20000 Frame> / <OrbitType> / <PositionAngleType>? Specifically when building StateCovariance() & StateCovarianceBlender(), and also when executing the interpolated covariance’s .changeCovarianceFrame() method to get back to RIC

6

Hello @aerow610,

You are welcome ! It’s true that i could have provided a more complete answer, sorry about that :downcast_face_with_sweat:.

CCSDS CDM Blue book vs Pink book

It is to be expected that it worked only when you built the StateCovariance from the OEM tabulated values with this frame. Intuitively speaking, local orbital frame are not seen as pseudo-inertial. However they are indeed fixed at epoch when computing a state covariance.

This is not explicitely stated in Annex A of CCSDS ODM Blue book (see this) but this is properly explained in Annex F of CDM pink book version(see this) which also provide a link to the list of possible frames (see this) (please don’t mind the ODM / CDM mix here).

Blending vs Hermite interpolation

As explained by @afvy and in the technical note he linked (many thanks by the way), using Hermite interpolation may produce non-positive definite matrices. Given your OEM time step, it will be better to use blending to ensure this does not happen because the difference in interpolated data can be considered negligible for such case (you can look at Figure 13 inside the technical note) .

As long as you load your data properly, you can express it however you like. However you should keep in mind that going from orbital elements (Circular, Keplerian or Equinoctial) to Cartesian/LOF will perform a linearization (which you cannot completely avoid depending on your use case). My advice would be to stay consistent throughout your code so that your future self and other people will understand what you do later on :sweat_smile:.

Cheers,
Vincent