I succedded to run the code correcly on my environement, and now I’m trying to use this code from my own data.
Here the steps I followed :
Simulation of an a trajectory from Keplerian parameters & Keplerian propagator (got 720 * TimestampedPVCoordinates)
Exctration of some points of this trajectory to fit with the field of view of my Radar (got 131 * TimestampedPVCoordinates)
Noizing of those points from Sigmas of my Radar passing threw azimuth/elevation (got 131 * TimestampedPVCoordinates)
IOD (Gibbs, Herrick-Gibbs and others) from 3 points (first, middle and last) of those Noised measurements (got KeplerianOrbit)
Now I joining the code of GorgiAstro, just before the definition of the propagator. I use the DormandPrince853IntegratorBuilder and I have some troubles with the convergence of the BatchLS.
First of all, I have sometime an error message saying that it’s impossible to compute the excentric anonaly from the true anomaly (50 interations). Sometimes I have another message saying that the Jaccobian matrix is singular for my orbit. Sometimes the BatchLS is not able to converge after 1000 iterations, and finally sometimes a result is computed, but I’m not really satified with them (less good than the IOD estimation).
I was wondering if the problem is not comming from the fact that my estimator should use a KeplerianPropagator in order to fit with the initial Simulated orbit, But the object BatchLSEstimator is not able to take as an input a KeplerianPropagatorBuilder, it just needs a NumericalPropagatorBuilder.
I’m just a begginer with OD problems and orekit so it’s a complete mess right now, I don’t really know what to do. Is the Levenberg-Marquart method is not supposed to converge all the time ?
Can’t I just put a KeplerianPropagatorBuilder in input of the BatchLSEstimator ?
Beware that if you directly uses position at time tk to set up a measurement at time tk, you miss both the uplink time and the downlink time, i.e. the time it takes for the signal to go from radar to satellite and back from satellite to Earth. Both times are used in Orekit Range measurement, so it implies what you simulate is not the same as what Orekit computes and you will not be able to get accurate results. If you also takes other modifiers in the orbit determination (like tropospheric delays), they will also not be present in your simulated data.
If you want to simulate radar measurements (which are basically two-way range) and also consider modifiers like tropospheric delays, then you should better use the measurement generation feature. The same is true for azimuth/elevation (but here you will miss only the downlink time as the measurement does not care about what the full signal path is, only the last path from satellite to ground is relevant in this measurement).
If you configure your NumericalPropagatorBuilder to use only central attraction and not any perturbation (i.e. not full Earth potential, no solara radiation pressure, no drag, no 3rd body, not anything), then it will be equivalent to a Keplerian motion.
Levenberg-Marquardt is a very robust method. Of course if there is some inconsistency in the data, it may fail to converge, but for regular solvable problems, it should solve them.
Even if you computed your measurement in a naive way (forgetting about uplink/downlink and delays), you should have converged. You should not get an accurate orbit, but you should get something sensible. So I guess there is something wrong elsewhere. Are you sure your simulated measurements are consistent with the propagation? Is the force model correct (i.e. Keplerian only) in numerical propagator? Are units corrects (remember Orekit uses SI units everywhere, i.e. meters for Range and radians for azimuth/elevation)? Are dates correct? Are the frames OK (terrestrial frame for the radar, inertial frame for the spacecraft)?
Note yet because the KeplerianPropagator that is built by KeplerianPropagatorBuilder does not compute partial derivatives. There is ongoing work (Bryan as an intern for this just now) to allow using analytical propagators in OD, but this is still Work In Progress, it will be published after the internship, perhaps near the end of the year.
The initial 6 months internship has been changed to a 4 months internship. Therefore, the internship objective has been changed. Now it focuses on performing Orbit Determination using TLEPropagator (which works very well). The last step to do in order to have a generalized analytical orbit determination instead of a TLE orbit determination is a generalization of the interfaces and methods developed for the TLE orbit determination. I don’t think it is a difficult task but my intern, Thomas, will not have enough time to do that properly before the end of its internship (18th of September). Indeed, I would prefer more results with the TLE batch least squares orbit determination and also I would like a Kalman Filter orbit determination using TLE propagator.
However, before publishing the TLE orbit determination into Orekit, the generalization of the API will be perform.
For this filtering step, I did not used the Class Range. In fact, I just convert my TimeStampedPVCoordinates into Azimuth, Elevation & Range in my Groundsation Frame. After that, a simple filter is applied on those parameters (if cond OK, keep this TimeStampedPVCoordinates).
Nevertheless, I created Range objects to add them to the estimator. If I understand well what you said, I should take into account the downlink, i.e., simulate a time at which the radar takes the measurement (real time tk of the position + time needed to the information to go from satellite to the radar, or Range / Speed of light modulo the Tropsheric delays). I can imagine that for LEO satellites the delay is really small ?
About the units I only uses radians and meters, the dates are AbsoluteDate objects, generate by Orekit. And finally the Frames are the following, also I think everything is ok here:
Finally, about the KeplerianPropagator as an input of BatchLSEstimator, it’s not at all a problem is I can use a simple NumericalPropagator as your suggested.
Hi guys, I finally found why my process did not succed to converge.
In fact I forgot to add AngularAzEl Measurements to the estimator, and with only Ranges but not Angles, the Batch LS was not able to find out the solution (Quite obvious now).
Also, I’m trying to reuse the code of @yzokras to generate properly Radar Measurements in order to solve the small error due to downlink and uplink shifts.
