Numerical Propegation Producing Funny Results

Orekit usage
1

Hey Everyone,

I’ve been attempting to use orekit to simulate observations of a radar to test some UCT algorithms. Shouldn’t be too hard, right?

I have been using the numerical propegator in the following code adapted from one of the tutorials:

"private static NumericalPropagator createNumProp(final Orbit orbit, final double mass) throws OrekitException {
final double tol = NumericalPropagator.tolerances(1.0, orbit, orbit.getType());
final double minStep = 1.e-3;
final double maxStep = 1.e+3;
AdaptiveStepsizeIntegrator integrator = new DormandPrince853Integrator(minStep, maxStep, tol[0], tol[1]);
integrator.setInitialStepSize(10.);

NumericalPropagator numProp = new NumericalPropagator(integrator);
numProp.setInitialState(new SpacecraftState(orbit, mass));
numProp.setOrbitType(OrbitType.KEPLERIAN);
return numProp;
}

Blockquote

And then I call that method as follows:

Orbit initialOrbit = new KeplerianOrbit(a, e, i,pa,raan, anomaly, PositionAngle.MEAN , inertialFrame, initialDate, mu);

        final double dP       = 0.001;
        final double minStep  = 0.001;
        final double maxStep  = 1000;
        final double initStep = 60;
        final double[][] tolerance = NumericalPropagator.tolerances(dP, initialOrbit, OrbitType.KEPLERIAN);
        AdaptiveStepsizeIntegrator integrator = new DormandPrince853Integrator(minStep, maxStep, tolerance[0], tolerance[1]);
        integrator.setInitialStepSize(initStep);
       double mass = 1000;
       NumericalPropagator numProp = createNumProp(initialOrbit, mass);

Then at specific intervals I use .getX, .getY and .getZ to produce geocentric positions (I think), and save them to a .csv.

Well this isn’t working and I’m not sure why at all. Instead of an orbit, my output produces a line.Here is a screenshot of a matlab plot of the output. (imgur link) If you need any additional information let me know!

My hair is falling out at an alarming rate! Please help me out!

2

What is the code you are using to print positions? I would expect something like the following. (I haven’t tested this code myself.)

Orbit initialOrbit = new KeplerianOrbit(a, e, i,pa,raan, anomaly, 
        PositionAngle.MEAN , inertialFrame, initialDate, mu);
double mass = 1000;
NumericalPropagator numProp = createNumProp(initialOrbit, mass);
Frame ecef = FramesFactory.getITRF(...)
numProp.setMasterMode(60, (s, l) -> {System.out.println(s.getPVCoordinates(ecf));});
numProp.propagate(endDate);
1 Like
3

Hey thanks for the quick reply. So basically my print statement at this point propegates 120 seconds and then saves the state to a string builder. Late on in the code the string builder is used to save to a .csv file. The code to do so is the following:

Vector3D scloc = finalState.getPVCoordinates().getPosition();
Vector3D radarRead = radarPos.subtract(scloc);
double x = scloc.getX();
double y = scloc.getY();
double z = scloc.getZ();
builder.append(finalState.getDate()+“,”);
builder.append(x+“,”);
builder.append(y+“,”);
builder.append(z+“,”);
builder.append(‘\n’);
lengthShifted = lengthShifted +120;
finalState=
numProp.propagate(finalState.getDate().shiftedBy(120));

Basically I want to have an event detector to trigger these “radar” reads with an event detector, so I only get the relevant ones but first I have to troubleshoot.

In your code you used talked about the frame. Here is my code concerning the frame:

Frame earthFrame = FramesFactory.getITRF(IERSConventions.IERS_2010, true);
BodyShape earth = new OneAxisEllipsoid(Constants.WGS84_EARTH_EQUATORIAL_RADIUS,

Constants.WGS84_EARTH_FLATTENING, earthFrame);

I thought maybe it was my state vector, so the first run through I used the PV coordinates from one of the tutorial programs but now I have it changed to the following:

// Initial state definition : date, orbit
AbsoluteDate initialDate = new AbsoluteDate(2004, 01, 01, 23, 30, 00.000, TimeScalesFactory.getUTC());
double mu = 3.986004415e+14; // gravitation coefficient
Frame inertialFrame = FramesFactory.getEME2000(); // inertial frame for orbit definition
double a = 5.216527285956187e+07;
double e = 2.720000000000000e-03;
double i = 0.1;
double pa = 287.2221/180Math.PI;
double raan = 327.2813/180Math.PI;
double anomaly = 105.4665/180*Math.PI;

Is it maybe that I’m propagating piece-wise in sets of 120? Maybe my frame is wrong?

Thanks again for your help Evan!

Edit: Im debating rebuilding orekit in my eclipse

4

This part is probably a bug. finalState.getPVCoordinates().getPosition() returns the position in ECI, but I’m guessing radarPos is in ECF, unless you’ve omitted the part where you recompute the radar’s position in ECI at each time step.

Calling this in a loop should give correct results, but it is not great in terms of numerical accuracy or computational efficiency. Better to use master mode.

I don’t know why you’re only seeing linear motion as NumericalPropagator should add a central attraction model when you set the initial state. Try explicitly adding the a NewtonianAttraction force model using numProp.addForceModel(new NewtonianAttraction(...)).

If that doesn’t work can you post a concise, complete, and working example, preferably as a JUnit test case with assert statements that show how the observed behavior differs from the expected behavior.

Regards,
Evan

5

This is quite embarrassing but the error actually occured outside of orekit. The data import I used for visualisation on Matlab could not handle the amount of data and did Something strange.

Your help was greatly appreciated, and thank you for helping out on such a great project :slight_smile:

Regards,
Paul