Possible bug in MeasurementHandler

[xche]

Hi,

Recently I met with a problem that might be a bug in Orekit. And I would like to hear your opinions.

The problem occurs when determining orbits of multiple satellites. In the current software, MeasurementHandler.handleStep only checks the currentDate of interpolator 0 (from line 73 to 78):

> // Current state date for interpolator 0
> final AbsoluteDate currentDate = interpolators.get(0).getCurrentState().getDate();
> if (model.isForwardPropagation()  && next.getDate().compareTo(currentDate) > 0 ||
>     !model.isForwardPropagation() && next.getDate().compareTo(currentDate) < 0) {
>     return;
> }

However, when there are more than one satellites, although PropagatorsParallelizer will ensure them to propagate with the same step, there is still a minor time gap (maybe only a few milliseconds) among the satellites. Because they run independently and may not have exactly the same speed and ends strictly at the same moment. So, there is a chance that satellite 0 finishes this step first and triggers MeasurementHandler to process a measurement that has nothing to do with it.

Let me explain this in an example. Say I need to determine two satellite (0 and 1) orbits simultaneously. The next measurement that needs to process in MeasurementHandler.handleStep is a range observation from satellite 1 at time t0. While PropagatorsParallelizer ensures them to propagate simultaneously to the next step t0, the thread of satellite 0 finishes propagation first (a few milliseconds ahead of satellite 1), it will trigger MeasurementHandler.handleStep. And since the currentDate of interpolators.get(0) indeed reaches t0, it will pass the if-clause above and start to process this range measurement that belongs to satellite 1.

The reason I found this is that I sometimes met with a strange situation that the theta in AbstractODEStateInterpolator.getInterpolatedState is even larger than 1. I traced back to this and I think this part of MeasurementHandler causes the problem. Because the situation above happens, when it starts to process this range measurement, actually the satellite 1 hasn’t reached to t0 yet. So the interpolation will process a step that even larger than interpolation step, and cause theta larger than 1. Although in Orekit it allows to accept theta larger than 1, it still causes the loss of accuracy and requires longer time to process.

What I think would be a fix is to check currentDate of all relevant interpolators (in case more than one satellite involved like inter-satellite range) instead of only check interpolator 0. Something like moving the above if-clause to line 85:

> // get the observed measurement
> final ObservedMeasurement<?> observed = next.getMeasurement();
> // estimate the theoretical measurement
> final SpacecraftState[] states = new SpacecraftState[observed.getSatellites().size()];
> for (int i = 0; i < states.length; ++i) {
>     final ObservableSatellite satellite = observed.getSatellites().get(i);
>     final int satInd = satellite.getPropagatorIndex();
>     final AbsoluteDate currentDate = interpolators.get(satInd).getCurrentState().getDate();
>     if ((model.isForwardPropagation()  && (next.getDate().compareTo(currentDate) > 0)) ||
>         (!model.isForwardPropagation() && (next.getDate().compareTo(currentDate) < 0))) {
>         return;
>     }
>     states[i] = interpolators.get(satInd).getInterpolatedState(next.getDate());
> }

I hope I explain this clearly. And I wonder if I get this correct.

Best regards,
Xingchi

[luc]

I am really surprised by this, so it should be investigated in depth.
The interpolators that MeasurementHandler receives are gathered together by PropagatorsParallelizer between lines 175-180, and the selected date corresponds to the earliest ending propagator, so all propagators should already have reached this date. This is roughly depicted in the ASCII-art in the class javadoc.

Could you look in a debugger what is the range of all propagators in the PropagatorsParallelizer.propagate method that calls MeasurementHandler when this problem occurs?

[xche]

Hi Luc,

Sure, here is some snapshots from a test. I have two satellites to be determined. I used a fixed step integrator with step size of 5 s. I set a breakpoint at line 74 in MeasurementHandler.handleStep. From the first epoch 0, the program will stop twice every step size. The first time, the 2 interpolators got:

You can see that the 2nd interpolator hasn’t reach time 10 yet. And I press continue, then the next stop will be:

You can see that this time the 2nd interpolator reaches time 10. And when the propagation continues, the similar behavior also continues:

If I did not mess something up and got errors, I think this test could be repeated by others too.

[luc]

Could you look at the time field in the softCurrentState instead of the globalCurrentState?
The soft state corresponds to the restricted part and this is what should matter.
To make sure, you could also display the output of interpolator.getCurrentState().getDate() as getCurrentState() should use the softCurrentState.

[xche]

Hi Luc,

here are the snapshots of softCurrentState:

For softCurrentState, the first few steps they indeed have the same time, but after I propagate a bit longer, it shows similar behavior as globalCurrentState, like above.

I also display the output of interpolator.getCurrentState().getDate(), it shows the same as softCurrentState.

final AbsoluteDate currentDate0 = interpolators.get(0).getCurrentState().getDate();
final AbsoluteDate currentDate1 = interpolators.get(1).getCurrentState().getDate();

and the next stop:

[luc]

There is a problem here.
I am trying to locate it, however I do not find anywhere in my work environment any class that starts with ModifiedDormand... as shown in your snapshot (the type of the second step interpolator, the first one being an Adams interpolator). Do you use a customized version of Hipparchus here?

[xche]

Yes, sort of. I use a custom integrator. I use a fixed-step multi-step integrator (Adams-Bashforth & Adams-Moulton method), with a single-step integrator (Runge-Kutta) to start up.

Thanks for your hints. I tested several scenarios. I think the problem occurs when switching the integrator from the single-step integrator to multi-step integrator. I will look in more detail to find out the reason. If you have any clues, I’ll be glad to know.

[xche]

Hi Luc,

I think I found the problem. In my custom Adams interpolator, I accidentally messed up the global states with soft states. After fixing this, the problem seems no longer occur.

I think back then I was confused with this global and soft terms. But now it’s a bit clearer. Many thanks for your help!

[luc]

Fine, it’s a relief.

I agree the terms are not properly chosen (I am the one to blame for this). Also we should probably
move the restrictStep method from AbstractODEStateInterpolator abstract class up to the ODEStateInterpolator interface in Hipparchus, as it is really important and needed for many things,
like events handling.

[xche]

Haha, sorry for the false alarm. Probably a javadoc with more details to distinguish the two terms would be helpful (or maybe there is one but I missed it). But anyway, thanks again for your replies and help!