Search in sources :

Example 11 with OrbitType

use of org.orekit.orbits.OrbitType in project Orekit by CS-SI.

the class KalmanOrbitDeterminationTest method testLageos2.

@Test
public // Orbit determination for Lageos2 based on SLR (range) measurements
void testLageos2() throws URISyntaxException, IllegalArgumentException, IOException, OrekitException, ParseException {
    // Print results on console
    final boolean print = false;
    // input in tutorial resources directory/output
    final String inputPath = KalmanOrbitDeterminationTest.class.getClassLoader().getResource("orbit-determination/Lageos2/od_test_Lageos2.in").toURI().getPath();
    final File input = new File(inputPath);
    // configure Orekit data acces
    Utils.setDataRoot("orbit-determination/Lageos2:potential/icgem-format");
    GravityFieldFactory.addPotentialCoefficientsReader(new ICGEMFormatReader("eigen-6s-truncated", true));
    // Choice of an orbit type to use
    // Default for test is Cartesian
    final OrbitType orbitType = OrbitType.CARTESIAN;
    // Initial orbital Cartesian covariance matrix
    // These covariances are derived from the deltas between initial and reference orbits
    // So in a way they are "perfect"...
    // Cartesian covariance matrix initialization
    final RealMatrix cartesianOrbitalP = MatrixUtils.createRealDiagonalMatrix(new double[] { 1e4, 4e3, 1, 5e-3, 6e-5, 1e-4 });
    // Orbital Cartesian process noise matrix (Q)
    final RealMatrix cartesianOrbitalQ = MatrixUtils.createRealDiagonalMatrix(new double[] { 1.e-4, 1.e-4, 1.e-4, 1.e-10, 1.e-10, 1.e-10 });
    // Initial measurement covariance matrix and process noise matrix
    final RealMatrix measurementP = MatrixUtils.createRealDiagonalMatrix(new double[] { 1., 1., 1., 1. });
    final RealMatrix measurementQ = MatrixUtils.createRealDiagonalMatrix(new double[] { 1e-6, 1e-6, 1e-6, 1e-6 });
    // Kalman orbit determination run.
    ResultKalman kalmanLageos2 = run(input, orbitType, print, cartesianOrbitalP, cartesianOrbitalQ, null, null, measurementP, measurementQ);
    // Definition of the accuracy for the test
    final double distanceAccuracy = 0.86;
    final double velocityAccuracy = 4.12e-3;
    // Tests
    // Note: The reference initial orbit is the same as in the batch LS tests
    // -----
    // Number of measurements processed
    final int numberOfMeas = 258;
    Assert.assertEquals(numberOfMeas, kalmanLageos2.getNumberOfMeasurements());
    // Estimated position and velocity
    final Vector3D estimatedPos = kalmanLageos2.getEstimatedPV().getPosition();
    final Vector3D estimatedVel = kalmanLageos2.getEstimatedPV().getVelocity();
    // Reference position and velocity at initial date (same as in batch LS test)
    final Vector3D refPos0 = new Vector3D(-5532131.956902, 10025696.592156, -3578940.040009);
    final Vector3D refVel0 = new Vector3D(-3871.275109, -607.880985, 4280.972530);
    // Run the reference until Kalman last date
    final Orbit refOrbit = runReference(input, orbitType, refPos0, refVel0, null, kalmanLageos2.getEstimatedPV().getDate());
    final Vector3D refPos = refOrbit.getPVCoordinates().getPosition();
    final Vector3D refVel = refOrbit.getPVCoordinates().getVelocity();
    // Check distances
    final double dP = Vector3D.distance(refPos, estimatedPos);
    final double dV = Vector3D.distance(refVel, estimatedVel);
    Assert.assertEquals(0.0, dP, distanceAccuracy);
    Assert.assertEquals(0.0, dV, velocityAccuracy);
    // Print orbit deltas
    if (print) {
        System.out.println("Test performances:");
        System.out.format("\t%-30s\n", "ΔEstimated / Reference");
        System.out.format(Locale.US, "\t%-10s %20.6f\n", "ΔP [m]", dP);
        System.out.format(Locale.US, "\t%-10s %20.6f\n", "ΔV [m/s]", dV);
    }
    // Test on measurements parameters
    final List<DelegatingDriver> list = new ArrayList<DelegatingDriver>();
    list.addAll(kalmanLageos2.measurementsParameters.getDrivers());
    sortParametersChanges(list);
    // Batch LS values
    // final double[] stationOffSet = { 1.659203,  0.861250,  -0.885352 };
    // final double rangeBias = -0.286275;
    final double[] stationOffSet = { 0.298867, -0.137456, 0.013315 };
    final double rangeBias = 0.002390;
    Assert.assertEquals(stationOffSet[0], list.get(0).getValue(), distanceAccuracy);
    Assert.assertEquals(stationOffSet[1], list.get(1).getValue(), distanceAccuracy);
    Assert.assertEquals(stationOffSet[2], list.get(2).getValue(), distanceAccuracy);
    Assert.assertEquals(rangeBias, list.get(3).getValue(), distanceAccuracy);
    // test on statistic for the range residuals
    final long nbRange = 258;
    // Batch LS values
    // final double[] RefStatRange = { -2.431135, 2.218644, 0.038483, 0.982017 };
    final double[] RefStatRange = { -23.561314, 20.436464, 0.964164, 5.687187 };
    Assert.assertEquals(nbRange, kalmanLageos2.getRangeStat().getN());
    Assert.assertEquals(RefStatRange[0], kalmanLageos2.getRangeStat().getMin(), distanceAccuracy);
    Assert.assertEquals(RefStatRange[1], kalmanLageos2.getRangeStat().getMax(), distanceAccuracy);
    Assert.assertEquals(RefStatRange[2], kalmanLageos2.getRangeStat().getMean(), distanceAccuracy);
    Assert.assertEquals(RefStatRange[3], kalmanLageos2.getRangeStat().getStandardDeviation(), distanceAccuracy);
}
Also used : ICGEMFormatReader(org.orekit.forces.gravity.potential.ICGEMFormatReader) EquinoctialOrbit(org.orekit.orbits.EquinoctialOrbit) CartesianOrbit(org.orekit.orbits.CartesianOrbit) KeplerianOrbit(org.orekit.orbits.KeplerianOrbit) Orbit(org.orekit.orbits.Orbit) CircularOrbit(org.orekit.orbits.CircularOrbit) ArrayList(java.util.ArrayList) GeodeticPoint(org.orekit.bodies.GeodeticPoint) RealMatrix(org.hipparchus.linear.RealMatrix) Vector3D(org.hipparchus.geometry.euclidean.threed.Vector3D) OrbitType(org.orekit.orbits.OrbitType) DelegatingDriver(org.orekit.utils.ParameterDriversList.DelegatingDriver) File(java.io.File) Test(org.junit.Test)

Example 12 with OrbitType

use of org.orekit.orbits.OrbitType in project Orekit by CS-SI.

the class AbstractForceModelTest method integrateShiftedState.

private double[] integrateShiftedState(final NumericalPropagator propagator, final SpacecraftState state0, final AbsoluteDate targetDate, final int index, final double h) throws OrekitException {
    OrbitType orbitType = propagator.getOrbitType();
    PositionAngle angleType = propagator.getPositionAngleType();
    double[] a = new double[6];
    double[] aDot = new double[6];
    orbitType.mapOrbitToArray(state0.getOrbit(), angleType, a, aDot);
    a[index] += h;
    SpacecraftState shiftedState = new SpacecraftState(orbitType.mapArrayToOrbit(a, aDot, angleType, state0.getDate(), state0.getMu(), state0.getFrame()), state0.getAttitude(), state0.getMass());
    propagator.setInitialState(shiftedState);
    SpacecraftState integratedState = propagator.propagate(targetDate);
    orbitType.mapOrbitToArray(integratedState.getOrbit(), angleType, a, null);
    return a;
}
Also used : FieldSpacecraftState(org.orekit.propagation.FieldSpacecraftState) SpacecraftState(org.orekit.propagation.SpacecraftState) PositionAngle(org.orekit.orbits.PositionAngle) OrbitType(org.orekit.orbits.OrbitType)

Example 13 with OrbitType

use of org.orekit.orbits.OrbitType in project Orekit by CS-SI.

the class FieldNumericalPropagatorTest method doTestPropagationTypesElliptical.

private <T extends RealFieldElement<T>> void doTestPropagationTypesElliptical(Field<T> field) throws OrekitException {
    T zero = field.getZero();
    // setup
    final FieldAbsoluteDate<T> initDate = FieldAbsoluteDate.getJ2000Epoch(field);
    FieldSpacecraftState<T> initialState;
    FieldNumericalPropagator<T> propagator;
    final FieldVector3D<T> position = new FieldVector3D<>(zero.add(7.0e6), zero.add(1.0e6), zero.add(4.0e6));
    final FieldVector3D<T> velocity = new FieldVector3D<>(zero.add(-500.0), zero.add(8000.0), zero.add(1000.0));
    final FieldOrbit<T> orbit = new FieldEquinoctialOrbit<>(new FieldPVCoordinates<>(position, velocity), FramesFactory.getEME2000(), initDate, mu);
    initialState = new FieldSpacecraftState<>(orbit);
    OrbitType type = OrbitType.EQUINOCTIAL;
    double[][] tolerance = NumericalPropagator.tolerances(0.001, orbit.toOrbit(), type);
    AdaptiveStepsizeFieldIntegrator<T> integrator = new DormandPrince853FieldIntegrator<>(field, 0.001, 200, tolerance[0], tolerance[1]);
    integrator.setInitialStepSize(zero.add(60));
    propagator = new FieldNumericalPropagator<>(field, integrator);
    propagator.setOrbitType(type);
    propagator.setInitialState(initialState);
    ForceModel gravityField = new HolmesFeatherstoneAttractionModel(FramesFactory.getITRF(IERSConventions.IERS_2010, true), GravityFieldFactory.getNormalizedProvider(5, 5));
    propagator.addForceModel(gravityField);
    // Propagation of the initial at t + dt
    final FieldPVCoordinates<T> pv = initialState.getPVCoordinates();
    final T dP = zero.add(0.001);
    final T dV = pv.getPosition().getNormSq().multiply(pv.getVelocity().getNorm()).reciprocal().multiply(dP.multiply(initialState.getMu()));
    final FieldPVCoordinates<T> pvcM = propagateInType(initialState, dP, OrbitType.CARTESIAN, PositionAngle.MEAN, propagator);
    final FieldPVCoordinates<T> pviM = propagateInType(initialState, dP, OrbitType.CIRCULAR, PositionAngle.MEAN, propagator);
    final FieldPVCoordinates<T> pveM = propagateInType(initialState, dP, OrbitType.EQUINOCTIAL, PositionAngle.MEAN, propagator);
    final FieldPVCoordinates<T> pvkM = propagateInType(initialState, dP, OrbitType.KEPLERIAN, PositionAngle.MEAN, propagator);
    final FieldPVCoordinates<T> pvcE = propagateInType(initialState, dP, OrbitType.CARTESIAN, PositionAngle.ECCENTRIC, propagator);
    final FieldPVCoordinates<T> pviE = propagateInType(initialState, dP, OrbitType.CIRCULAR, PositionAngle.ECCENTRIC, propagator);
    final FieldPVCoordinates<T> pveE = propagateInType(initialState, dP, OrbitType.EQUINOCTIAL, PositionAngle.ECCENTRIC, propagator);
    final FieldPVCoordinates<T> pvkE = propagateInType(initialState, dP, OrbitType.KEPLERIAN, PositionAngle.ECCENTRIC, propagator);
    final FieldPVCoordinates<T> pvcT = propagateInType(initialState, dP, OrbitType.CARTESIAN, PositionAngle.TRUE, propagator);
    final FieldPVCoordinates<T> pviT = propagateInType(initialState, dP, OrbitType.CIRCULAR, PositionAngle.TRUE, propagator);
    final FieldPVCoordinates<T> pveT = propagateInType(initialState, dP, OrbitType.EQUINOCTIAL, PositionAngle.TRUE, propagator);
    final FieldPVCoordinates<T> pvkT = propagateInType(initialState, dP, OrbitType.KEPLERIAN, PositionAngle.TRUE, propagator);
    Assert.assertEquals(0, pvcM.getPosition().subtract(pveT.getPosition()).getNorm().getReal() / dP.getReal(), 3.0);
    Assert.assertEquals(0, pvcM.getVelocity().subtract(pveT.getVelocity()).getNorm().getReal() / dV.getReal(), 2.0);
    Assert.assertEquals(0, pviM.getPosition().subtract(pveT.getPosition()).getNorm().getReal() / dP.getReal(), 0.6);
    Assert.assertEquals(0, pviM.getVelocity().subtract(pveT.getVelocity()).getNorm().getReal() / dV.getReal(), 0.4);
    Assert.assertEquals(0, pvkM.getPosition().subtract(pveT.getPosition()).getNorm().getReal() / dP.getReal(), 0.5);
    Assert.assertEquals(0, pvkM.getVelocity().subtract(pveT.getVelocity()).getNorm().getReal() / dV.getReal(), 0.3);
    Assert.assertEquals(0, pveM.getPosition().subtract(pveT.getPosition()).getNorm().getReal() / dP.getReal(), 0.2);
    Assert.assertEquals(0, pveM.getVelocity().subtract(pveT.getVelocity()).getNorm().getReal() / dV.getReal(), 0.2);
    Assert.assertEquals(0, pvcE.getPosition().subtract(pveT.getPosition()).getNorm().getReal() / dP.getReal(), 3.0);
    Assert.assertEquals(0, pvcE.getVelocity().subtract(pveT.getVelocity()).getNorm().getReal() / dV.getReal(), 2.0);
    Assert.assertEquals(0, pviE.getPosition().subtract(pveT.getPosition()).getNorm().getReal() / dP.getReal(), 0.03);
    Assert.assertEquals(0, pviE.getVelocity().subtract(pveT.getVelocity()).getNorm().getReal() / dV.getReal(), 0.04);
    Assert.assertEquals(0, pvkE.getPosition().subtract(pveT.getPosition()).getNorm().getReal() / dP.getReal(), 0.4);
    Assert.assertEquals(0, pvkE.getVelocity().subtract(pveT.getVelocity()).getNorm().getReal() / dV.getReal(), 0.3);
    Assert.assertEquals(0, pveE.getPosition().subtract(pveT.getPosition()).getNorm().getReal() / dP.getReal(), 0.2);
    Assert.assertEquals(0, pveE.getVelocity().subtract(pveT.getVelocity()).getNorm().getReal() / dV.getReal(), 0.07);
    Assert.assertEquals(0, pvcT.getPosition().subtract(pveT.getPosition()).getNorm().getReal() / dP.getReal(), 3.0);
    Assert.assertEquals(0, pvcT.getVelocity().subtract(pveT.getVelocity()).getNorm().getReal() / dV.getReal(), 2.0);
    Assert.assertEquals(0, pviT.getPosition().subtract(pveT.getPosition()).getNorm().getReal() / dP.getReal(), 0.3);
    Assert.assertEquals(0, pviT.getVelocity().subtract(pveT.getVelocity()).getNorm().getReal() / dV.getReal(), 0.2);
    Assert.assertEquals(0, pvkT.getPosition().subtract(pveT.getPosition()).getNorm().getReal() / dP.getReal(), 0.4);
    Assert.assertEquals(0, pvkT.getVelocity().subtract(pveT.getVelocity()).getNorm().getReal() / dV.getReal(), 0.2);
}
Also used : DormandPrince853FieldIntegrator(org.hipparchus.ode.nonstiff.DormandPrince853FieldIntegrator) ForceModel(org.orekit.forces.ForceModel) FieldVector3D(org.hipparchus.geometry.euclidean.threed.FieldVector3D) FieldEquinoctialOrbit(org.orekit.orbits.FieldEquinoctialOrbit) OrbitType(org.orekit.orbits.OrbitType) HolmesFeatherstoneAttractionModel(org.orekit.forces.gravity.HolmesFeatherstoneAttractionModel)

Example 14 with OrbitType

use of org.orekit.orbits.OrbitType in project Orekit by CS-SI.

the class FieldNumericalPropagatorTest method doTestEphemerisDatesBackward.

private <T extends RealFieldElement<T>> void doTestEphemerisDatesBackward(Field<T> field) throws OrekitException {
    T zero = field.getZero();
    // setup
    TimeScale tai = TimeScalesFactory.getTAI();
    FieldAbsoluteDate<T> initialDate = new FieldAbsoluteDate<>(field, "2015-07-05", tai);
    FieldAbsoluteDate<T> startDate = new FieldAbsoluteDate<>(field, "2015-07-03", tai).shiftedBy(-0.1);
    FieldAbsoluteDate<T> endDate = new FieldAbsoluteDate<>(field, "2015-07-04", tai);
    Frame eci = FramesFactory.getGCRF();
    FieldKeplerianOrbit<T> orbit = new FieldKeplerianOrbit<>(zero.add(600e3 + Constants.WGS84_EARTH_EQUATORIAL_RADIUS), zero, zero, zero, zero, zero, PositionAngle.TRUE, eci, initialDate, mu);
    OrbitType type = OrbitType.CARTESIAN;
    double[][] tol = NumericalPropagator.tolerances(1e-3, orbit.toOrbit(), type);
    FieldNumericalPropagator<T> prop = new FieldNumericalPropagator<>(field, new DormandPrince853FieldIntegrator<>(field, 0.1, 500, tol[0], tol[1]));
    prop.setOrbitType(type);
    prop.resetInitialState(new FieldSpacecraftState<>(new FieldCartesianOrbit<>(orbit)));
    // action
    prop.setEphemerisMode();
    prop.propagate(endDate, startDate);
    FieldBoundedPropagator<T> ephemeris = prop.getGeneratedEphemeris();
    // verify
    TimeStampedFieldPVCoordinates<T> actualPV = ephemeris.getPVCoordinates(startDate, eci);
    TimeStampedFieldPVCoordinates<T> expectedPV = orbit.getPVCoordinates(startDate, eci);
    MatcherAssert.assertThat(actualPV.getPosition().toVector3D(), OrekitMatchers.vectorCloseTo(expectedPV.getPosition().toVector3D(), 1.0));
    MatcherAssert.assertThat(actualPV.getVelocity().toVector3D(), OrekitMatchers.vectorCloseTo(expectedPV.getVelocity().toVector3D(), 1.0));
    MatcherAssert.assertThat(ephemeris.getMinDate().durationFrom(startDate).getReal(), OrekitMatchers.closeTo(0, 0));
    MatcherAssert.assertThat(ephemeris.getMaxDate().durationFrom(endDate).getReal(), OrekitMatchers.closeTo(0, 0));
    // test date
    FieldAbsoluteDate<T> date = endDate.shiftedBy(-0.11);
    Assert.assertEquals(ephemeris.propagate(date).getDate().durationFrom(date).getReal(), 0, 0);
}
Also used : Frame(org.orekit.frames.Frame) TimeScale(org.orekit.time.TimeScale) FieldCartesianOrbit(org.orekit.orbits.FieldCartesianOrbit) FieldKeplerianOrbit(org.orekit.orbits.FieldKeplerianOrbit) OrbitType(org.orekit.orbits.OrbitType) FieldAbsoluteDate(org.orekit.time.FieldAbsoluteDate)

Example 15 with OrbitType

use of org.orekit.orbits.OrbitType in project Orekit by CS-SI.

the class FieldNumericalPropagatorTest method doTestStopEvent.

private <T extends RealFieldElement<T>> void doTestStopEvent(Field<T> field) throws OrekitException {
    T zero = field.getZero();
    // setup
    final FieldAbsoluteDate<T> initDate = FieldAbsoluteDate.getJ2000Epoch(field);
    FieldSpacecraftState<T> initialState;
    FieldNumericalPropagator<T> propagator;
    final FieldVector3D<T> position = new FieldVector3D<>(zero.add(7.0e6), zero.add(1.0e6), zero.add(4.0e6));
    final FieldVector3D<T> velocity = new FieldVector3D<>(zero.add(-500.0), zero.add(8000.0), zero.add(1000.0));
    final FieldOrbit<T> orbit = new FieldEquinoctialOrbit<>(new FieldPVCoordinates<>(position, velocity), FramesFactory.getEME2000(), initDate, mu);
    initialState = new FieldSpacecraftState<>(orbit);
    OrbitType type = OrbitType.EQUINOCTIAL;
    double[][] tolerance = NumericalPropagator.tolerances(0.001, orbit.toOrbit(), type);
    AdaptiveStepsizeFieldIntegrator<T> integrator = new DormandPrince853FieldIntegrator<>(field, 0.001, 200, tolerance[0], tolerance[1]);
    integrator.setInitialStepSize(zero.add(60));
    propagator = new FieldNumericalPropagator<>(field, integrator);
    propagator.setOrbitType(type);
    propagator.setInitialState(initialState);
    final FieldAbsoluteDate<T> stopDate = initDate.shiftedBy(1000);
    CheckingHandler<FieldDateDetector<T>, T> checking = new CheckingHandler<FieldDateDetector<T>, T>(Action.STOP);
    propagator.addEventDetector(new FieldDateDetector<>(stopDate).withHandler(checking));
    Assert.assertEquals(1, propagator.getEventsDetectors().size());
    checking.assertEvent(false);
    final FieldSpacecraftState<T> finalState = propagator.propagate(initDate.shiftedBy(3200));
    checking.assertEvent(true);
    Assert.assertEquals(0, finalState.getDate().durationFrom(stopDate).getReal(), 1.0e-10);
    propagator.clearEventsDetectors();
    Assert.assertEquals(0, propagator.getEventsDetectors().size());
}
Also used : DormandPrince853FieldIntegrator(org.hipparchus.ode.nonstiff.DormandPrince853FieldIntegrator) FieldDateDetector(org.orekit.propagation.events.FieldDateDetector) FieldVector3D(org.hipparchus.geometry.euclidean.threed.FieldVector3D) FieldEquinoctialOrbit(org.orekit.orbits.FieldEquinoctialOrbit) OrbitType(org.orekit.orbits.OrbitType)

Aggregations

OrbitType (org.orekit.orbits.OrbitType)69 Test (org.junit.Test)39 NumericalPropagator (org.orekit.propagation.numerical.NumericalPropagator)38 SpacecraftState (org.orekit.propagation.SpacecraftState)35 FieldSpacecraftState (org.orekit.propagation.FieldSpacecraftState)31 DormandPrince853Integrator (org.hipparchus.ode.nonstiff.DormandPrince853Integrator)29 KeplerianOrbit (org.orekit.orbits.KeplerianOrbit)28 Orbit (org.orekit.orbits.Orbit)28 DormandPrince853FieldIntegrator (org.hipparchus.ode.nonstiff.DormandPrince853FieldIntegrator)27 FieldAbsoluteDate (org.orekit.time.FieldAbsoluteDate)25 FieldKeplerianOrbit (org.orekit.orbits.FieldKeplerianOrbit)24 Frame (org.orekit.frames.Frame)23 FieldNumericalPropagator (org.orekit.propagation.numerical.FieldNumericalPropagator)23 AbstractLegacyForceModelTest (org.orekit.forces.AbstractLegacyForceModelTest)21 AbsoluteDate (org.orekit.time.AbsoluteDate)17 PVCoordinates (org.orekit.utils.PVCoordinates)17 AdaptiveStepsizeIntegrator (org.hipparchus.ode.nonstiff.AdaptiveStepsizeIntegrator)16 FieldVector3D (org.hipparchus.geometry.euclidean.threed.FieldVector3D)15 CartesianOrbit (org.orekit.orbits.CartesianOrbit)15 PositionAngle (org.orekit.orbits.PositionAngle)15