Covariance propagation for GPS uncertainty

Orekit Python Wrapper usage
1

Hello all,

I’m trying to investigate the propagation of GPS orbit uncertainty over time using a covariance analysis. This is my first time working on this kind of problem, so I’m not entirely sure about both the implementation and the underlying physics.

It seems that my code is working, using the StateCovarianceMatrixProvider and the associated State Transition Matrix.

I’m using the CovariancePropagation tutorial to do so and the well described “Orbit and Covariance Interpolation/Blending with Orekit” Document.

The initial sigmas are :

SIGMA_R = 2.5 m
SIGMA_S = 2.5 m
SIGMA_W = 2.5 m
SIGMA_VR = 0.5 m/s
SIGMA_VS = 0.5 m/s
SIGMA_VW = 0.5 m/s

I’m basically using a NumericalPropagator to do so, with a SSO @ 400 km.
I’m using the following handler to retrieve covariance during the propagation (in the RTN frame)

Covariance handler 
@JImplements(OrekitFixedStepHandler)
class Prop_Handler:
    def __init__(self, cov_provider, harvester):
        self.cov_provider = cov_provider
        self.harvester    = harvester
        self.data = {
            "times_s":   [],
            "states_pv": [],
            "cov_qsw":   [],
            "stm":       [],
        }

    @JOverride
    def init(self, s0, t, step):
        pass

    @JOverride
    def handleStep(self, state):
        self.data["times_s"].append(state.getDate().durationFrom(EPOCH))

        pv = state.getPVCoordinates()
        r, v = pv.getPosition(), pv.getVelocity()
        self.data["states_pv"].append([
            r.getX(), r.getY(), r.getZ(),
            v.getX(), v.getY(), v.getZ(),
        ])

        sc_eci = self.cov_provider.getStateCovariance(state)
        sc_qsw = sc_eci.changeCovarianceFrame(state.getOrbit(), LOFType.QSW)
        self.data["cov_qsw"].append(mat6(sc_qsw.getMatrix()))

        stm_j = self.harvester.getStateTransitionMatrix(state)
        self.data["stm"].append(mat6(stm_j))

    @JOverride
    def finish(self, finalState):
        pass

And here is the code :

Covariance propagation code 
STM_NAME = "stm"
COV_NAME = "covariance"

harvester = propagator.setupMatricesComputation(STM_NAME, None, None)

P0_diag = [SIGMA_R**2, SIGMA_S**2, SIGMA_W**2,SIGMA_VR**2, SIGMA_VS**2, SIGMA_VW**2,]
P0_java = MatrixUtils.createRealDiagonalMatrix(P0_diag)
cov_init = StateCovariance(P0_java, EPOCH, gcrf,OrbitType.CARTESIAN, PositionAngleType.TRUE)


cov_provider = StateCovarianceMatrixProvider(COV_NAME, STM_NAME, harvester, cov_init)
propagator.addAdditionalDataProvider(cov_provider)

times_s_list      = []
states_pv_list    = []
cov_eci_list      = []   
cov_qsw_list      = []   
stm_list          = []  


T_PROP  = N_ORBITS * T_ORB
handler = Prop_Handler(cov_provider, harvester)
propagator.getMultiplexer().add(DT_OUTPUT, handler)
final_state = propagator.propagate(EPOCH.shiftedBy(T_PROP))

times_s   = np.array(handler.data["times_s"])   
states_pv = np.array(handler.data["states_pv"]) 
cov_qsw   = np.array(handler.data["cov_qsw"])    
stm_arr   = np.array(handler.data["stm"])         
N         = len(times_s)

sigma_qsw = np.sqrt(np.array([np.diag(cov_qsw[k])[:3] for k in range(N)]))
sc_final    = cov_provider.getStateCovariance(final_state)
P_qsw_final = mat6(sc_final.changeCovarianceFrame(final_state.getOrbit(), LOFType.QSW).getMatrix())

Here are the outputs :

fig1_sigma_qsw
fig1_sigma_qsw1800×750 86.8 KB

fig2_ellipsoids_qsw
fig2_ellipsoids_qsw2100×900 99.7 KB

Questions :

  1. Is it normal to observe a divergence in the along-track sigma? If not, could there be an issue with the code?
  2. How can the apparent periodic behavior of the normal and radial sigmas be explained? (When zooming in, it seems they are actually growing rather than purely periodic.)
  3. Can I do the same by changing the initial sigmas to (200 m, 1000 m, 200 m) to investigate a more general TLE uncertainty?

Sorry if my questions are basic, and thank you for any help you can provide!

Cheers,

2

Those are beautiful results and exactly what you should expect to see!

The way I think about the increase in uncertainty (in intrack direction) is that it’s due to the initial uncertainty in velocity. The position and velocity are correlated, so as you propagate forward in time that velocity uncertainty contributes to the position uncertainty. The reason that intrack/crosstrack/radial all behave differently is due to the orbit dynamics; delta-v in the intrack direction has a larger influence than a delta-v in the cross-track direction.

You should definitely play around with the initial covariance. Keep in mind the initial velocity uncertainty has an (arguably larger) influence on the growth of the position covariance over time.

1 Like