Doubts on covariance propagation with stepHandler

[MaximeJ]

Hello @toquika,

@toquika

Can you confirm me that with the MatricesHarvester.getStateTransitionMatrix what I retrieve is actually the error state transition matrix (which is nothing else than the Jacobian matrix dY/dY0)?

Indeed what you get is the Jacobian matrix dY/dY0 that is called error state transition matrix \Phi in Vallado’s book.

@toquika

With this implementation, the results are really bad. To get reasonable results, I had to change the handleStep method as follows:

Which propagator are you using ? The harvester.setReferenceState(currentState); does not re-initialize the internal matrices, as stated in the Javadoc:

    /** Set up reference state.
     * <p>
     * This method is called whenever the global propagation reference state changes.
     * This corresponds to the start of propagation in batch least squares orbit determination
     * or at prediction step for each measurement in Kalman filtering. Its goal is to allow
     * the harvester to compute some internal data. Analytical models like TLE use it to
     * compute analytical derivatives, semi-analytical models like DSST use it to compute
     * short periodic terms, numerical models do not use it at all.
     * </p>
     * @param reference reference state to set
     */

If you want to re-initialize from current state I think you need to instantiate a new harvester.

@toquika

In this case, it seems to me that I am always computing the STM based on the initial reference conditions “Y0”. Since the STM comes from a linearization process, isn’t it wrong to always refer to that initial point? That is why in my previous version I thought it was correct to set a different reference state of the harvester at each iteration, and then calculate the new covariance matrix starting from the last one (and not always from the initial covariance matrix).

State transition matrices are “transitive” so this would do the same, meaning:
With: \Phi(t_i, t_j)=(\frac{\partial{y(t_i)}}{\partial{y(t_j)}})_{6 \times 6} (see for example, Montenbrück and Gill, Satellite Orbits, equation (7.1))
\Phi(t_2, t_0) = \Phi(t_2, t_1).\Phi(t_1, t_0) (see Vallado, p.780, after equation (10-27))
Let P_0, P_1, P_2 be your covariances at times t_0, t_1, t_2
P_1 = \Phi(t_1,t_0).P_0.\Phi^T(t_1,t_0)
P_2 = \Phi(t_2,t_0).P_0.\Phi^T(t_2,t_0) (this is what is done in the tutorial)
And:
P_2 = \Phi(t_2,t_1).\Phi(t_1,t_0).P_0.\Phi^T(t_1,t_0).\Phi^T(t_2,t_1)=\Phi(t_2,t_1).P_1.\Phi^T(t_2,t_1) (this is what you would like to do)
You see here that, if you don’t change your propagation model, starting from initial state at t_0 with P_0 or from an intermediate state at t_1 with P_1 will give you the same covariance P_2 at t_2.

Side note: there is a merge request on the forum written by @bcazabonne that intends to add an additional state for covariance propagation.

Hope this helps you,
Maxime