signalTimeOfFlight derivative in terms of Keplerian parameters

Orekit usage
1

Hello!

It has been a while since the last time I wrote something on this forum. But I am still using Orekit everyday and I am always glad to interact with you :slight_smile:

I am still working on the orbit determination problem described here. As @MaximeJ told me, since the derivatives of my measurements (in Cartesian coordinates) are constant ([n_1, n_2, n_3, 0, 0, 0], where n_i is a constant value), the lines of the Jacobian are almost identical,
which is also a classic cause of singular problem (non pseudo-inversible matrix) exception. Therefore, I was asking myself if it was possible to compute the derivatives of my measurements in terms of Keplerian parameters instead of the Cartesian ones. So also the last three elements of my derivative array (i.e. the ones linked to the satellite velocity in the Cartesian case) won’t be zero, which I think could help.

The measurement, in my case, is defined as (I am extending the Range class):

tau = Gradient.cast_(self.signalTimeOfFlight(pvaSSB, pvaDS.getPosition(), pvaDS.getDate()).multiply(projection_factor))

where pvaSSB is the pva triplet for the Solar System Barycentre (assumed to be constant in time) and projection_factor accounts for the fact that my Range refers to the projection of the SSB-SAT range along the direction which is given by [n_1, n_2, n_3].

Then, I get the derivative array as:

derivative_array = tau.getGradient()

which is given to the EstimatedMeasurement object as

estimated_measurement.setStateDerivatives(0, derivative_array).

Is there a way to give to estimated_measurement the derivatives with respect to the Keplerian parameters and then to solve a least square problem (I am using a SequentialBatchLSEstimator) in this way? I am asking this since I always get a singular matrix error when running the filter with a SequentialGaussNewtonOptimizer(QRDecomposer, False, evaluation). It works if I make use of a pseudo-inversion by replacing False with True, but I was interested in trying also this alternative way!

Best regards,
Samuele

2

Hi,

If you want to compute first order partial derivatives of anything w.r.t. Keplerian elements (so six independent variables), you can do it with automatic differentiation by first defining a FieldKeplerianOrbit, of Gradient type. You would need to have defined the orbital elements as differentiable parameters, for example fieldSma = Gradient.variable(6, 0, sma), fieldEcc = Gradient.variable(6, 1, ecc), etc. Then you can perform your computation, as long as is had a Field implementation. At the end, call the getGradient method.

As for modifying the least squares, it should be doable by overloading classes and/or methods.
But if I’m not wrong, if your propagator builder uses Keplerian as OrbitType, the Jacobian of measurements should already be w.r.t. to the initial conditions in that coordinate system. Anyway maybe I’ve missed your point.

Best,
Romain.

3

Hi Serrof,

thank you so much for your reply. It has been very useful.

I found out that I was using a NumericalPropagatorBuilder (giving to it a KeplerianOrbit object) instead of a KeplerianPropagatorBuilder. The results I obtain by using the Keplerian one are very different from those I got using the Numerical one… I think that it is because of the

estimated_measurement.setStateDerivatives(0, derivative_array)

since I am still passing the derivative_array in Cartesian coordinates in the extended Range class. Am I right? Or is there an automatic conversion which is performed internally?

And also when I call the estimator.getPhysicalCovariances() method I should get the results in terms of Keplerian coordinates right? The values I get are very low so I assume they are.

Many thanks!

Samuele

4

Hi @samuele,

This should only change the propagation theory, so you’ll get different derivatives but they should still be with respect to Cartesian or Keplerian elements depending on what you chose.
Keplerian propagator is an analytical two-body-theory propagator, with no perturbations whatsoever. So you need to make sure that this is what you need.

There is no automatic conversion that I know of.
All the built-in measurements in Orekit are performing derivatives computation w.r.t Cartesian coordinates internally in method (AbstractMeasurement).theoreticalEvaluation.

Then, for example for the batch least-square estimator, the derivatives of the measurement w.r.t to desired orbit type (say Keplerian in your case) are computed in AbstractBatchLSModel.fetchEvaluatedMeasurements like this:
(Where k is the propagator index)

// partial derivatives of the current Cartesian coordinates with respect to current orbital state
            final double[][] aCY = new double[6][6];
            final Orbit currentOrbit = evaluationStates[k].getOrbit();
            currentOrbit.getJacobianWrtParameters(builders[p].getPositionAngle(), aCY);
            final RealMatrix dCdY = new Array2DRowRealMatrix(aCY, false);

Here dCdY is the jacobian (i.e. partial derivatives matrix) of Cartesian state parameters w.r.t Keplerian.

            // Jacobian of the measurement with respect to current orbital state
            final RealMatrix dMdC = new Array2DRowRealMatrix(evaluation.getStateDerivatives(k), false);

dMdC is the jacobian of the derivatives of the measurement w.r.t Cartesian coordinates (this is the matrix you computed in derivative_array).

            final RealMatrix dMdY = dMdC.multiply(dCdY);

This final line gives you the jacobian of the measurement wrt Keplerian parameters.

Yes, if the Orbit object in the PropagatorBuilder is given in Keplerian formalism, then the estimation process and the covariance matrix will be performed with Keplerian elements.

Hope this helps,
Maxime

5

Hi @MaximeJ,

thank you very much.

I am now using a simple Keplerian propagator to run a simple test (to actually understand if everything is working correctly). I understand that the use of a KeplerianPropagatorBuilder instead of a NumericalPropagatorBuilder should only change the propagation theory, but now I am really running two different simulations that only differ in the choice of the propagator (I really just replace that line changing the Numerical builder to the Keplerian one, that I build starting from a KeplerianOrbit object) and the results are significantly different. The orbit determination error diverges in the NumericalPropagatorBuilder case while it converges in the KeplerianPropagatorBuilder case…

What could be the reason for this? Before reading your answer, I thought that it could be a problem related to the estimated_measurement.setStateDerivatives(0, derivative_array) line. I thought it could be an issue related to the fact that the lines of the Jacobian were almost identical (in Cartesian terms).

And finally, you mentioned the batch least-square estimator, in which the derivatives of the measurements are computed as you reported. Is it the same considering the sequential batch least-square estimator?

Thank you for your support.

Best,
Samuele