Examples with ICGEMFormatReader org.orekit.forces.gravity.potential.ICGEMFormatReader used on opensource projects

Example 6 with ICGEMFormatReader

use of org.orekit.forces.gravity.potential.ICGEMFormatReader in project Orekit by CS-SI.

the class IntegratedEphemerisTest method setUp.

@Before
public void setUp() {
    Utils.setDataRoot("regular-data:potential/icgem-format");
    GravityFieldFactory.addPotentialCoefficientsReader(new ICGEMFormatReader("eigen-6s-truncated", true));
    // Definition of initial conditions with position and velocity
    Vector3D position = new Vector3D(7.0e6, 1.0e6, 4.0e6);
    Vector3D velocity = new Vector3D(-500.0, 8000.0, 1000.0);
    double mu = 3.9860047e14;
    AbsoluteDate initDate = AbsoluteDate.J2000_EPOCH.shiftedBy(584.);
    initialOrbit = new EquinoctialOrbit(new PVCoordinates(position, velocity), FramesFactory.getEME2000(), initDate, mu);
    // Numerical propagator definition
    double[] absTolerance = { 0.0001, 1.0e-11, 1.0e-11, 1.0e-8, 1.0e-8, 1.0e-8, 0.001 };
    double[] relTolerance = { 1.0e-8, 1.0e-8, 1.0e-8, 1.0e-9, 1.0e-9, 1.0e-9, 1.0e-7 };
    AdaptiveStepsizeIntegrator integrator = new DormandPrince853Integrator(0.001, 500, absTolerance, relTolerance);
    integrator.setInitialStepSize(100);
    numericalPropagator = new NumericalPropagator(integrator);
}

Example 7 with ICGEMFormatReader

use of org.orekit.forces.gravity.potential.ICGEMFormatReader in project Orekit by CS-SI.

the class AttitudesSequenceTest method testResetDuringTransitionForward.

@Test
public void testResetDuringTransitionForward() throws OrekitException {
    // Initial state definition : date, orbit
    final AbsoluteDate initialDate = new AbsoluteDate(2004, 01, 01, 23, 30, 00.000, TimeScalesFactory.getUTC());
    final Vector3D position = new Vector3D(-6142438.668, 3492467.560, -25767.25680);
    final Vector3D velocity = new Vector3D(505.8479685, 942.7809215, 7435.922231);
    final Orbit initialOrbit = new KeplerianOrbit(new PVCoordinates(position, velocity), FramesFactory.getEME2000(), initialDate, Constants.EIGEN5C_EARTH_MU);
    final OneAxisEllipsoid earth = new OneAxisEllipsoid(Constants.WGS84_EARTH_EQUATORIAL_RADIUS, Constants.WGS84_EARTH_FLATTENING, FramesFactory.getITRF(IERSConventions.IERS_2010, true));
    final TopocentricFrame volgograd = new TopocentricFrame(earth, new GeodeticPoint(FastMath.toRadians(48.7), FastMath.toRadians(44.5), 24.0), "Волгоград");
    final AttitudesSequence attitudesSequence = new AttitudesSequence();
    final double transitionTime = 250.0;
    final AttitudeProvider nadirPointing = new NadirPointing(initialOrbit.getFrame(), earth);
    final AttitudeProvider targetPointing = new TargetPointing(initialOrbit.getFrame(), volgograd.getPoint(), earth);
    final ElevationDetector eventDetector = new ElevationDetector(volgograd).withConstantElevation(FastMath.toRadians(5.0)).withHandler(new ContinueOnEvent<>());
    final List<AbsoluteDate> nadirToTarget = new ArrayList<>();
    attitudesSequence.addSwitchingCondition(nadirPointing, targetPointing, eventDetector, true, false, transitionTime, AngularDerivativesFilter.USE_RR, (previous, next, state) -> nadirToTarget.add(state.getDate()));
    final double[][] tolerance = NumericalPropagator.tolerances(10.0, initialOrbit, initialOrbit.getType());
    final AdaptiveStepsizeIntegrator integrator = new DormandPrince853Integrator(0.001, 300.0, tolerance[0], tolerance[1]);
    final NumericalPropagator propagator = new NumericalPropagator(integrator);
    GravityFieldFactory.addPotentialCoefficientsReader(new ICGEMFormatReader("g007_eigen_05c_coef", false));
    propagator.addForceModel(new HolmesFeatherstoneAttractionModel(earth.getBodyFrame(), GravityFieldFactory.getNormalizedProvider(8, 8)));
    propagator.setInitialState(new SpacecraftState(initialOrbit, nadirPointing.getAttitude(initialOrbit, initialOrbit.getDate(), initialOrbit.getFrame())));
    propagator.setAttitudeProvider(attitudesSequence);
    attitudesSequence.registerSwitchEvents(propagator);
    propagator.propagate(initialDate.shiftedBy(6000));
    // check that if we restart a forward propagation from an intermediate state
    // we properly get an interpolated attitude despite we missed the event trigger
    final AbsoluteDate midTransition = nadirToTarget.get(0).shiftedBy(0.5 * transitionTime);
    SpacecraftState state = propagator.propagate(midTransition.shiftedBy(-60), midTransition);
    Rotation nadirR = nadirPointing.getAttitude(state.getOrbit(), state.getDate(), state.getFrame()).getRotation();
    Rotation targetR = targetPointing.getAttitude(state.getOrbit(), state.getDate(), state.getFrame()).getRotation();
    final double reorientationAngle = Rotation.distance(nadirR, targetR);
    Assert.assertEquals(0.5 * reorientationAngle, Rotation.distance(state.getAttitude().getRotation(), nadirR), 0.03 * reorientationAngle);
}

Example 8 with ICGEMFormatReader

use of org.orekit.forces.gravity.potential.ICGEMFormatReader in project Orekit by CS-SI.

the class AttitudesSequenceTest method testOutOfSyncCalls.

@Test
public void testOutOfSyncCalls() throws OrekitException {
    // Initial state definition : date, orbit
    final AbsoluteDate initialDate = new AbsoluteDate(2004, 01, 01, 23, 30, 00.000, TimeScalesFactory.getUTC());
    final Vector3D position = new Vector3D(-6142438.668, 3492467.560, -25767.25680);
    final Vector3D velocity = new Vector3D(505.8479685, 942.7809215, 7435.922231);
    final Orbit initialOrbit = new KeplerianOrbit(new PVCoordinates(position, velocity), FramesFactory.getEME2000(), initialDate, Constants.EIGEN5C_EARTH_MU);
    final OneAxisEllipsoid earth = new OneAxisEllipsoid(Constants.WGS84_EARTH_EQUATORIAL_RADIUS, Constants.WGS84_EARTH_FLATTENING, FramesFactory.getITRF(IERSConventions.IERS_2010, true));
    final TopocentricFrame volgograd = new TopocentricFrame(earth, new GeodeticPoint(FastMath.toRadians(48.7), FastMath.toRadians(44.5), 24.0), "Волгоград");
    final AttitudesSequence attitudesSequence = new AttitudesSequence();
    final double transitionTime = 250.0;
    final AttitudeProvider nadirPointing = new NadirPointing(initialOrbit.getFrame(), earth);
    final AttitudeProvider targetPointing = new TargetPointing(initialOrbit.getFrame(), volgograd.getPoint(), earth);
    final ElevationDetector eventDetector = new ElevationDetector(volgograd).withConstantElevation(FastMath.toRadians(5.0)).withHandler(new ContinueOnEvent<>());
    final Handler nadirToTarget = new Handler(nadirPointing, targetPointing);
    attitudesSequence.addSwitchingCondition(nadirPointing, targetPointing, eventDetector, true, false, transitionTime, AngularDerivativesFilter.USE_RR, nadirToTarget);
    final Handler targetToNadir = new Handler(targetPointing, nadirPointing);
    attitudesSequence.addSwitchingCondition(targetPointing, nadirPointing, eventDetector, false, true, transitionTime, AngularDerivativesFilter.USE_RR, targetToNadir);
    final double[][] tolerance = NumericalPropagator.tolerances(10.0, initialOrbit, initialOrbit.getType());
    final AdaptiveStepsizeIntegrator integrator = new DormandPrince853Integrator(0.001, 300.0, tolerance[0], tolerance[1]);
    final NumericalPropagator propagator = new NumericalPropagator(integrator);
    GravityFieldFactory.addPotentialCoefficientsReader(new ICGEMFormatReader("g007_eigen_05c_coef", false));
    propagator.addForceModel(new HolmesFeatherstoneAttractionModel(earth.getBodyFrame(), GravityFieldFactory.getNormalizedProvider(8, 8)));
    propagator.setInitialState(new SpacecraftState(initialOrbit, nadirPointing.getAttitude(initialOrbit, initialOrbit.getDate(), initialOrbit.getFrame())));
    propagator.setAttitudeProvider(attitudesSequence);
    attitudesSequence.registerSwitchEvents(propagator);
    propagator.setMasterMode(10, (state, isLast) -> {
        Attitude nadirAttitude = nadirPointing.getAttitude(state.getOrbit(), state.getDate(), state.getFrame());
        Attitude targetAttitude = targetPointing.getAttitude(state.getOrbit(), state.getDate(), state.getFrame());
        Attitude stateAttitude = state.getAttitude();
        if (nadirToTarget.dates.isEmpty() || state.getDate().durationFrom(nadirToTarget.dates.get(0)) < 0) {
            // we are stabilized in nadir pointing, before first switch
            checkEqualAttitudes(nadirAttitude, stateAttitude);
        } else if (state.getDate().durationFrom(nadirToTarget.dates.get(0)) <= transitionTime) {
            // we are in transition from nadir to target
            checkBetweenAttitudes(nadirAttitude, targetAttitude, stateAttitude);
        } else if (targetToNadir.dates.isEmpty() || state.getDate().durationFrom(targetToNadir.dates.get(0)) < 0) {
            // we are stabilized in target pointing between the two switches
            checkEqualAttitudes(targetAttitude, stateAttitude);
        } else if (state.getDate().durationFrom(targetToNadir.dates.get(0)) <= transitionTime) {
            // we are in transition from target to nadir
            checkBetweenAttitudes(targetAttitude, nadirAttitude, stateAttitude);
        } else {
            // we are stabilized back in nadir pointing, after second switch
            checkEqualAttitudes(nadirAttitude, stateAttitude);
        }
    });
    propagator.propagate(initialDate.shiftedBy(6000));
}

Example 9 with ICGEMFormatReader

use of org.orekit.forces.gravity.potential.ICGEMFormatReader in project Orekit by CS-SI.

the class KalmanOrbitDeterminationTest method testW3B.

@Test
public // Orbit determination for range, azimuth elevation measurements
void testW3B() 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/W3B/od_test_W3.in").toURI().getPath();
    final File input = new File(inputPath);
    // Configure Orekit data access
    Utils.setDataRoot("orbit-determination/W3B: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[] { FastMath.pow(2.4e4, 2), FastMath.pow(1.e5, 2), FastMath.pow(4.e4, 2), FastMath.pow(3.5, 2), FastMath.pow(2., 2), FastMath.pow(0.6, 2) });
    // 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 });
    // Propagation covariance and process noise matrices
    final RealMatrix propagationP = MatrixUtils.createRealDiagonalMatrix(new double[] { // Cd
    FastMath.pow(2., 2), // leak-X
    FastMath.pow(5.7e-6, 2), // leak-X
    FastMath.pow(1.1e-11, 2), // leak-Y
    FastMath.pow(7.68e-7, 2), // leak-Y
    FastMath.pow(1.26e-10, 2), // leak-Z
    FastMath.pow(5.56e-6, 2), // leak-Z
    FastMath.pow(2.79e-10, 2) });
    final RealMatrix propagationQ = MatrixUtils.createRealDiagonalMatrix(new double[] { // Cd
    FastMath.pow(1e-3, 2), // Leaks
    0., // Leaks
    0., // Leaks
    0., // Leaks
    0., // Leaks
    0., // Leaks
    0. });
    // Measurement covariance and process noise matrices
    // az/el bias sigma = 0.06deg
    // range bias sigma = 100m
    final double angularVariance = FastMath.pow(FastMath.toRadians(0.06), 2);
    final double rangeVariance = FastMath.pow(500., 2);
    final RealMatrix measurementP = MatrixUtils.createRealDiagonalMatrix(new double[] { angularVariance, angularVariance, rangeVariance, angularVariance, angularVariance, rangeVariance, angularVariance, angularVariance, rangeVariance, angularVariance, angularVariance, rangeVariance, angularVariance, angularVariance, rangeVariance });
    // Process noise sigma: 1e-6 for all
    final double measQ = FastMath.pow(1e-6, 2);
    final RealMatrix measurementQ = MatrixUtils.createRealIdentityMatrix(measurementP.getRowDimension()).scalarMultiply(measQ);
    // Kalman orbit determination run.
    ResultKalman kalmanW3B = run(input, orbitType, print, cartesianOrbitalP, cartesianOrbitalQ, propagationP, propagationQ, measurementP, measurementQ);
    // Tests
    // -----
    // Definition of the accuracy for the test
    final double distanceAccuracy = 0.1;
    // degrees
    final double angleAccuracy = 1e-5;
    // Number of measurements processed
    final int numberOfMeas = 521;
    Assert.assertEquals(numberOfMeas, kalmanW3B.getNumberOfMeasurements());
    // Test on propagator parameters
    // -----------------------------
    // Batch LS result
    // final double dragCoef  = -0.2154;
    final double dragCoef = 0.1931;
    Assert.assertEquals(dragCoef, kalmanW3B.propagatorParameters.getDrivers().get(0).getValue(), 1e-3);
    final Vector3D leakAcceleration0 = new Vector3D(kalmanW3B.propagatorParameters.getDrivers().get(1).getValue(), kalmanW3B.propagatorParameters.getDrivers().get(3).getValue(), kalmanW3B.propagatorParameters.getDrivers().get(5).getValue());
    // Batch LS results
    // Assert.assertEquals(8.002e-6, leakAcceleration0.getNorm(), 1.0e-8);
    Assert.assertEquals(5.994e-6, leakAcceleration0.getNorm(), 1.0e-8);
    final Vector3D leakAcceleration1 = new Vector3D(kalmanW3B.propagatorParameters.getDrivers().get(2).getValue(), kalmanW3B.propagatorParameters.getDrivers().get(4).getValue(), kalmanW3B.propagatorParameters.getDrivers().get(6).getValue());
    // Batch LS results
    // Assert.assertEquals(3.058e-10, leakAcceleration1.getNorm(), 1.0e-12);
    Assert.assertEquals(1.831e-10, leakAcceleration1.getNorm(), 1.0e-12);
    // Test on measurements parameters
    // -------------------------------
    final List<DelegatingDriver> list = new ArrayList<DelegatingDriver>();
    list.addAll(kalmanW3B.measurementsParameters.getDrivers());
    sortParametersChanges(list);
    // Station CastleRock
    // Batch LS results
    // final double[] CastleAzElBias  = { 0.062701342, -0.003613508 };
    // final double   CastleRangeBias = 11274.4677;
    final double[] CastleAzElBias = { 0.062635, -0.003672 };
    final double CastleRangeBias = 11289.3678;
    Assert.assertEquals(CastleAzElBias[0], FastMath.toDegrees(list.get(0).getValue()), angleAccuracy);
    Assert.assertEquals(CastleAzElBias[1], FastMath.toDegrees(list.get(1).getValue()), angleAccuracy);
    Assert.assertEquals(CastleRangeBias, list.get(2).getValue(), distanceAccuracy);
    // Station Fucino
    // Batch LS results
    // final double[] FucAzElBias  = { -0.053526137, 0.075483886 };
    // final double   FucRangeBias = 13467.8256;
    final double[] FucAzElBias = { -0.053298, 0.075589 };
    final double FucRangeBias = 13482.0715;
    Assert.assertEquals(FucAzElBias[0], FastMath.toDegrees(list.get(3).getValue()), angleAccuracy);
    Assert.assertEquals(FucAzElBias[1], FastMath.toDegrees(list.get(4).getValue()), angleAccuracy);
    Assert.assertEquals(FucRangeBias, list.get(5).getValue(), distanceAccuracy);
    // Station Kumsan
    // Batch LS results
    // final double[] KumAzElBias  = { -0.023574208, -0.054520756 };
    // final double   KumRangeBias = 13512.57594;
    final double[] KumAzElBias = { -0.022805, -0.055057 };
    final double KumRangeBias = 13502.7459;
    Assert.assertEquals(KumAzElBias[0], FastMath.toDegrees(list.get(6).getValue()), angleAccuracy);
    Assert.assertEquals(KumAzElBias[1], FastMath.toDegrees(list.get(7).getValue()), angleAccuracy);
    Assert.assertEquals(KumRangeBias, list.get(8).getValue(), distanceAccuracy);
    // Station Pretoria
    // Batch LS results
    // final double[] PreAzElBias = { 0.030201539, 0.009747877 };
    // final double PreRangeBias = 13594.11889;
    final double[] PreAzElBias = { 0.030353, 0.009658 };
    final double PreRangeBias = 13609.2516;
    Assert.assertEquals(PreAzElBias[0], FastMath.toDegrees(list.get(9).getValue()), angleAccuracy);
    Assert.assertEquals(PreAzElBias[1], FastMath.toDegrees(list.get(10).getValue()), angleAccuracy);
    Assert.assertEquals(PreRangeBias, list.get(11).getValue(), distanceAccuracy);
    // Station Uralla
    // Batch LS results
    // final double[] UraAzElBias = { 0.167814449, -0.12305252 };
    // final double UraRangeBias = 13450.26738;
    final double[] UraAzElBias = { 0.167519, -0.122842 };
    final double UraRangeBias = 13441.7019;
    Assert.assertEquals(UraAzElBias[0], FastMath.toDegrees(list.get(12).getValue()), angleAccuracy);
    Assert.assertEquals(UraAzElBias[1], FastMath.toDegrees(list.get(13).getValue()), angleAccuracy);
    Assert.assertEquals(UraRangeBias, list.get(14).getValue(), distanceAccuracy);
    // Test on statistic for the range residuals
    final long nbRange = 182;
    // statistics for the range residual (min, max, mean, std)
    final double[] RefStatRange = { -12.981, 18.046, -1.133, 5.312 };
    Assert.assertEquals(nbRange, kalmanW3B.getRangeStat().getN());
    Assert.assertEquals(RefStatRange[0], kalmanW3B.getRangeStat().getMin(), distanceAccuracy);
    Assert.assertEquals(RefStatRange[1], kalmanW3B.getRangeStat().getMax(), distanceAccuracy);
    Assert.assertEquals(RefStatRange[2], kalmanW3B.getRangeStat().getMean(), distanceAccuracy);
    Assert.assertEquals(RefStatRange[3], kalmanW3B.getRangeStat().getStandardDeviation(), distanceAccuracy);
    // test on statistic for the azimuth residuals
    final long nbAzi = 339;
    // statistics for the azimuth residual (min, max, mean, std)
    final double[] RefStatAzi = { -0.041441, 0.023473, -0.004426, 0.009911 };
    Assert.assertEquals(nbAzi, kalmanW3B.getAzimStat().getN());
    Assert.assertEquals(RefStatAzi[0], kalmanW3B.getAzimStat().getMin(), angleAccuracy);
    Assert.assertEquals(RefStatAzi[1], kalmanW3B.getAzimStat().getMax(), angleAccuracy);
    Assert.assertEquals(RefStatAzi[2], kalmanW3B.getAzimStat().getMean(), angleAccuracy);
    Assert.assertEquals(RefStatAzi[3], kalmanW3B.getAzimStat().getStandardDeviation(), angleAccuracy);
    // test on statistic for the elevation residuals
    final long nbEle = 339;
    final double[] RefStatEle = { -0.025399, 0.043345, 0.001011, 0.010636 };
    Assert.assertEquals(nbEle, kalmanW3B.getElevStat().getN());
    Assert.assertEquals(RefStatEle[0], kalmanW3B.getElevStat().getMin(), angleAccuracy);
    Assert.assertEquals(RefStatEle[1], kalmanW3B.getElevStat().getMax(), angleAccuracy);
    Assert.assertEquals(RefStatEle[2], kalmanW3B.getElevStat().getMean(), angleAccuracy);
    Assert.assertEquals(RefStatEle[3], kalmanW3B.getElevStat().getStandardDeviation(), angleAccuracy);
    RealMatrix covariances = kalmanW3B.getCovariances();
    Assert.assertEquals(28, covariances.getRowDimension());
    Assert.assertEquals(28, covariances.getColumnDimension());
    // drag coefficient variance
    Assert.assertEquals(0.016349, covariances.getEntry(6, 6), 1.0e-5);
    // leak-X constant term variance
    Assert.assertEquals(2.047303E-13, covariances.getEntry(7, 7), 1.0e-16);
    // leak-Y constant term variance
    Assert.assertEquals(5.462497E-13, covariances.getEntry(9, 9), 1.0e-15);
    // leak-Z constant term variance
    Assert.assertEquals(1.717781E-11, covariances.getEntry(11, 11), 1.0e-15);
    // Test on orbital parameters
    // Done at the end to avoid changing the estimated propagation parameters
    // ----------------------------------------------------------------------
    // Estimated position and velocity
    final Vector3D estimatedPos = kalmanW3B.getEstimatedPV().getPosition();
    final Vector3D estimatedVel = kalmanW3B.getEstimatedPV().getVelocity();
    // Reference position and velocity at initial date (same as in batch LS test)
    final Vector3D refPos0 = new Vector3D(-40541446.255, -9905357.41, 206777.413);
    final Vector3D refVel0 = new Vector3D(759.0685, -1476.5156, 54.793);
    // Gather the selected propagation parameters and initialize them to the values found
    // with the batch LS method
    final ParameterDriversList refPropagationParameters = kalmanW3B.propagatorParameters;
    final double dragCoefRef = -0.215433133145843;
    final double[] leakXRef = { +5.69040439901955E-06, 1.09710906802403E-11 };
    final double[] leakYRef = { -7.66440256777678E-07, 1.25467464335066E-10 };
    final double[] leakZRef = { -5.574055079952E-06, 2.78703463746911E-10 };
    for (DelegatingDriver driver : refPropagationParameters.getDrivers()) {
        switch(driver.getName()) {
            case "drag coefficient":
                driver.setValue(dragCoefRef);
                break;
            case "leak-X[0]":
                driver.setValue(leakXRef[0]);
                break;
            case "leak-X[1]":
                driver.setValue(leakXRef[1]);
                break;
            case "leak-Y[0]":
                driver.setValue(leakYRef[0]);
                break;
            case "leak-Y[1]":
                driver.setValue(leakYRef[1]);
                break;
            case "leak-Z[0]":
                driver.setValue(leakZRef[0]);
                break;
            case "leak-Z[1]":
                driver.setValue(leakZRef[1]);
                break;
        }
    }
    // Run the reference until Kalman last date
    final Orbit refOrbit = runReference(input, orbitType, refPos0, refVel0, refPropagationParameters, kalmanW3B.getEstimatedPV().getDate());
    // Test on last orbit
    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);
    // FIXME: debug - Comparison with batch LS is bad
    final double debugDistanceAccuracy = 234.73;
    final double debugVelocityAccuracy = 0.086;
    Assert.assertEquals(0.0, Vector3D.distance(refPos, estimatedPos), debugDistanceAccuracy);
    Assert.assertEquals(0.0, Vector3D.distance(refVel, estimatedVel), debugVelocityAccuracy);
    // 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);
    }
}

Example 10 with ICGEMFormatReader

use of org.orekit.forces.gravity.potential.ICGEMFormatReader in project Orekit by CS-SI.

the class KalmanOrbitDeterminationTest method runReference.

/**
 * Use the physical models in the input file
 * Incorporate the initial reference values
 * And run the propagation until the last measurement to get the reference orbit at the same date
 * as the Kalman filter
 * @param input Input configuration file
 * @param orbitType Orbit type to use (calculation and display)
 * @param refPosition Initial reference position
 * @param refVelocity Initial reference velocity
 * @param refPropagationParameters Reference propagation parameters
 * @param kalmanFinalDate The final date of the Kalman filter
 * @return The reference orbit at the same date as the Kalman filter
 * @throws IOException Input file cannot be opened
 * @throws IllegalArgumentException Issue in key/value reading of input file
 * @throws OrekitException An Orekit exception... should be explicit
 * @throws ParseException Parsing of the input file or measurement file failed
 */
private Orbit runReference(final File input, final OrbitType orbitType, final Vector3D refPosition, final Vector3D refVelocity, final ParameterDriversList refPropagationParameters, final AbsoluteDate kalmanFinalDate) throws IOException, IllegalArgumentException, OrekitException, ParseException {
    // Read input parameters
    KeyValueFileParser<ParameterKey> parser = new KeyValueFileParser<ParameterKey>(ParameterKey.class);
    parser.parseInput(input.getAbsolutePath(), new FileInputStream(input));
    // Gravity field
    GravityFieldFactory.addPotentialCoefficientsReader(new ICGEMFormatReader("eigen-5c.gfc", true));
    final NormalizedSphericalHarmonicsProvider gravityField = createGravityField(parser);
    // Orbit initial guess
    Orbit initialRefOrbit = new CartesianOrbit(new PVCoordinates(refPosition, refVelocity), parser.getInertialFrame(ParameterKey.INERTIAL_FRAME), parser.getDate(ParameterKey.ORBIT_DATE, TimeScalesFactory.getUTC()), gravityField.getMu());
    // Convert to desired orbit type
    initialRefOrbit = orbitType.convertType(initialRefOrbit);
    // IERS conventions
    final IERSConventions conventions;
    if (!parser.containsKey(ParameterKey.IERS_CONVENTIONS)) {
        conventions = IERSConventions.IERS_2010;
    } else {
        conventions = IERSConventions.valueOf("IERS_" + parser.getInt(ParameterKey.IERS_CONVENTIONS));
    }
    // Central body
    final OneAxisEllipsoid body = createBody(parser);
    // Propagator builder
    final NumericalPropagatorBuilder propagatorBuilder = createPropagatorBuilder(parser, conventions, gravityField, body, initialRefOrbit);
    // Force the selected propagation parameters to their reference values
    if (refPropagationParameters != null) {
        for (DelegatingDriver refDriver : refPropagationParameters.getDrivers()) {
            for (DelegatingDriver driver : propagatorBuilder.getPropagationParametersDrivers().getDrivers()) {
                if (driver.getName().equals(refDriver.getName())) {
                    driver.setValue(refDriver.getValue());
                }
            }
        }
    }
    // Build the reference propagator
    final NumericalPropagator propagator = propagatorBuilder.buildPropagator(propagatorBuilder.getSelectedNormalizedParameters());
    // Propagate until last date and return the orbit
    return propagator.propagate(kalmanFinalDate).getOrbit();
}