Examples with ForceModel org.orekit.forces.ForceModel used on opensource projects

Example 26 with ForceModel

use of org.orekit.forces.ForceModel in project Orekit by CS-SI.

the class NumericalPropagatorTest method testForceModelInitialized.

@Test
public void testForceModelInitialized() throws OrekitException {
    // setup
    // mutable holders
    SpacecraftState[] actualState = new SpacecraftState[1];
    AbsoluteDate[] actualDate = new AbsoluteDate[1];
    ForceModel force = new ForceModelAdapter() {

        @Override
        public void init(SpacecraftState initialState, AbsoluteDate target) {
            actualState[0] = initialState;
            actualDate[0] = target;
        }
    };
    // action
    propagator.setOrbitType(OrbitType.CARTESIAN);
    propagator.addForceModel(force);
    AbsoluteDate target = initDate.shiftedBy(60);
    propagator.propagate(target);
    // verify
    Assert.assertThat(actualDate[0], CoreMatchers.is(target));
    Assert.assertThat(actualState[0].getDate().durationFrom(initDate), CoreMatchers.is(0.0));
    Assert.assertThat(actualState[0].getPVCoordinates(), OrekitMatchers.pvIs(initialState.getPVCoordinates()));
}

Example 27 with ForceModel

use of org.orekit.forces.ForceModel in project Orekit by CS-SI.

the class PartialDerivativesTest method doTestParametersDerivatives.

private void doTestParametersDerivatives(String parameterName, double tolerance, OrbitType... orbitTypes) throws OrekitException {
    OneAxisEllipsoid earth = new OneAxisEllipsoid(Constants.WGS84_EARTH_EQUATORIAL_RADIUS, Constants.WGS84_EARTH_FLATTENING, FramesFactory.getITRF(IERSConventions.IERS_2010, true));
    ForceModel drag = new DragForce(new HarrisPriester(CelestialBodyFactory.getSun(), earth), new IsotropicDrag(2.5, 1.2));
    NormalizedSphericalHarmonicsProvider provider = GravityFieldFactory.getNormalizedProvider(5, 5);
    ForceModel gravityField = new HolmesFeatherstoneAttractionModel(FramesFactory.getITRF(IERSConventions.IERS_2010, true), provider);
    Orbit baseOrbit = new KeplerianOrbit(7000000.0, 0.01, 0.1, 0.7, 0, 1.2, PositionAngle.TRUE, FramesFactory.getEME2000(), AbsoluteDate.J2000_EPOCH, provider.getMu());
    double dt = 900;
    double dP = 1.0;
    for (OrbitType orbitType : orbitTypes) {
        final Orbit initialOrbit = orbitType.convertType(baseOrbit);
        for (PositionAngle angleType : PositionAngle.values()) {
            NumericalPropagator propagator = setUpPropagator(initialOrbit, dP, orbitType, angleType, drag, gravityField);
            propagator.setMu(provider.getMu());
            for (final ForceModel forceModel : propagator.getAllForceModels()) {
                for (final ParameterDriver driver : forceModel.getParametersDrivers()) {
                    driver.setValue(driver.getReferenceValue());
                    driver.setSelected(driver.getName().equals(parameterName));
                }
            }
            PartialDerivativesEquations partials = new PartialDerivativesEquations("partials", propagator);
            final SpacecraftState initialState = partials.setInitialJacobians(new SpacecraftState(initialOrbit));
            propagator.setInitialState(initialState);
            final JacobiansMapper mapper = partials.getMapper();
            PickUpHandler pickUp = new PickUpHandler(mapper, null);
            propagator.setMasterMode(pickUp);
            propagator.propagate(initialState.getDate().shiftedBy(dt));
            double[][] dYdP = pickUp.getdYdP();
            // compute reference Jacobian using finite differences
            double[][] dYdPRef = new double[6][1];
            NumericalPropagator propagator2 = setUpPropagator(initialOrbit, dP, orbitType, angleType, drag, gravityField);
            propagator2.setMu(provider.getMu());
            ParameterDriversList bound = new ParameterDriversList();
            for (final ForceModel forceModel : propagator2.getAllForceModels()) {
                for (final ParameterDriver driver : forceModel.getParametersDrivers()) {
                    if (driver.getName().equals(parameterName)) {
                        driver.setSelected(true);
                        bound.add(driver);
                    } else {
                        driver.setSelected(false);
                    }
                }
            }
            ParameterDriver selected = bound.getDrivers().get(0);
            double p0 = selected.getReferenceValue();
            double h = selected.getScale();
            selected.setValue(p0 - 4 * h);
            propagator2.resetInitialState(arrayToState(stateToArray(initialState, orbitType, angleType, true), orbitType, angleType, initialState.getFrame(), initialState.getDate(), // the mu may have been reset above
            propagator2.getMu(), initialState.getAttitude()));
            SpacecraftState sM4h = propagator2.propagate(initialOrbit.getDate().shiftedBy(dt));
            selected.setValue(p0 - 3 * h);
            propagator2.resetInitialState(arrayToState(stateToArray(initialState, orbitType, angleType, true), orbitType, angleType, initialState.getFrame(), initialState.getDate(), // the mu may have been reset above
            propagator2.getMu(), initialState.getAttitude()));
            SpacecraftState sM3h = propagator2.propagate(initialOrbit.getDate().shiftedBy(dt));
            selected.setValue(p0 - 2 * h);
            propagator2.resetInitialState(arrayToState(stateToArray(initialState, orbitType, angleType, true), orbitType, angleType, initialState.getFrame(), initialState.getDate(), // the mu may have been reset above
            propagator2.getMu(), initialState.getAttitude()));
            SpacecraftState sM2h = propagator2.propagate(initialOrbit.getDate().shiftedBy(dt));
            selected.setValue(p0 - 1 * h);
            propagator2.resetInitialState(arrayToState(stateToArray(initialState, orbitType, angleType, true), orbitType, angleType, initialState.getFrame(), initialState.getDate(), // the mu may have been reset above
            propagator2.getMu(), initialState.getAttitude()));
            SpacecraftState sM1h = propagator2.propagate(initialOrbit.getDate().shiftedBy(dt));
            selected.setValue(p0 + 1 * h);
            propagator2.resetInitialState(arrayToState(stateToArray(initialState, orbitType, angleType, true), orbitType, angleType, initialState.getFrame(), initialState.getDate(), // the mu may have been reset above
            propagator2.getMu(), initialState.getAttitude()));
            SpacecraftState sP1h = propagator2.propagate(initialOrbit.getDate().shiftedBy(dt));
            selected.setValue(p0 + 2 * h);
            propagator2.resetInitialState(arrayToState(stateToArray(initialState, orbitType, angleType, true), orbitType, angleType, initialState.getFrame(), initialState.getDate(), // the mu may have been reset above
            propagator2.getMu(), initialState.getAttitude()));
            SpacecraftState sP2h = propagator2.propagate(initialOrbit.getDate().shiftedBy(dt));
            selected.setValue(p0 + 3 * h);
            propagator2.resetInitialState(arrayToState(stateToArray(initialState, orbitType, angleType, true), orbitType, angleType, initialState.getFrame(), initialState.getDate(), // the mu may have been reset above
            propagator2.getMu(), initialState.getAttitude()));
            SpacecraftState sP3h = propagator2.propagate(initialOrbit.getDate().shiftedBy(dt));
            selected.setValue(p0 + 4 * h);
            propagator2.resetInitialState(arrayToState(stateToArray(initialState, orbitType, angleType, true), orbitType, angleType, initialState.getFrame(), initialState.getDate(), // the mu may have been reset above
            propagator2.getMu(), initialState.getAttitude()));
            SpacecraftState sP4h = propagator2.propagate(initialOrbit.getDate().shiftedBy(dt));
            fillJacobianColumn(dYdPRef, 0, orbitType, angleType, h, sM4h, sM3h, sM2h, sM1h, sP1h, sP2h, sP3h, sP4h);
            for (int i = 0; i < 6; ++i) {
                Assert.assertEquals(dYdPRef[i][0], dYdP[i][0], FastMath.abs(dYdPRef[i][0] * tolerance));
            }
        }
    }
}

Example 28 with ForceModel

use of org.orekit.forces.ForceModel in project Orekit by CS-SI.

the class PartialDerivativesTest method setUpPropagator.

private NumericalPropagator setUpPropagator(Orbit orbit, double dP, OrbitType orbitType, PositionAngle angleType, ForceModel... models) throws OrekitException {
    final double minStep = 0.001;
    final double maxStep = 1000;
    double[][] tol = NumericalPropagator.tolerances(dP, orbit, orbitType);
    NumericalPropagator propagator = new NumericalPropagator(new DormandPrince853Integrator(minStep, maxStep, tol[0], tol[1]));
    propagator.setOrbitType(orbitType);
    propagator.setPositionAngleType(angleType);
    for (ForceModel model : models) {
        propagator.addForceModel(model);
    }
    return propagator;
}

Example 29 with ForceModel

use of org.orekit.forces.ForceModel in project Orekit by CS-SI.

the class PartialDerivativesTest method testPropagationTypesElliptical.

@Test
public void testPropagationTypesElliptical() throws OrekitException {
    NormalizedSphericalHarmonicsProvider provider = GravityFieldFactory.getNormalizedProvider(5, 5);
    ForceModel gravityField = new HolmesFeatherstoneAttractionModel(FramesFactory.getITRF(IERSConventions.IERS_2010, true), provider);
    Orbit initialOrbit = new KeplerianOrbit(8000000.0, 0.01, 0.1, 0.7, 0, 1.2, PositionAngle.TRUE, FramesFactory.getEME2000(), AbsoluteDate.J2000_EPOCH, provider.getMu());
    double dt = 900;
    double dP = 0.001;
    for (OrbitType orbitType : OrbitType.values()) {
        for (PositionAngle angleType : PositionAngle.values()) {
            // compute state Jacobian using PartialDerivatives
            NumericalPropagator propagator = setUpPropagator(initialOrbit, dP, orbitType, angleType, gravityField);
            PartialDerivativesEquations partials = new PartialDerivativesEquations("partials", propagator);
            final SpacecraftState initialState = partials.setInitialJacobians(new SpacecraftState(initialOrbit));
            propagator.setInitialState(initialState);
            final JacobiansMapper mapper = partials.getMapper();
            PickUpHandler pickUp = new PickUpHandler(mapper, null);
            propagator.setMasterMode(pickUp);
            propagator.propagate(initialState.getDate().shiftedBy(dt));
            double[][] dYdY0 = pickUp.getdYdY0();
            // compute reference state Jacobian using finite differences
            double[][] dYdY0Ref = new double[6][6];
            AbstractIntegratedPropagator propagator2 = setUpPropagator(initialOrbit, dP, orbitType, angleType, gravityField);
            double[] steps = NumericalPropagator.tolerances(1000000 * dP, initialOrbit, orbitType)[0];
            for (int i = 0; i < 6; ++i) {
                propagator2.resetInitialState(shiftState(initialState, orbitType, angleType, -4 * steps[i], i));
                SpacecraftState sM4h = propagator2.propagate(initialState.getDate().shiftedBy(dt));
                propagator2.resetInitialState(shiftState(initialState, orbitType, angleType, -3 * steps[i], i));
                SpacecraftState sM3h = propagator2.propagate(initialState.getDate().shiftedBy(dt));
                propagator2.resetInitialState(shiftState(initialState, orbitType, angleType, -2 * steps[i], i));
                SpacecraftState sM2h = propagator2.propagate(initialState.getDate().shiftedBy(dt));
                propagator2.resetInitialState(shiftState(initialState, orbitType, angleType, -1 * steps[i], i));
                SpacecraftState sM1h = propagator2.propagate(initialState.getDate().shiftedBy(dt));
                propagator2.resetInitialState(shiftState(initialState, orbitType, angleType, 1 * steps[i], i));
                SpacecraftState sP1h = propagator2.propagate(initialState.getDate().shiftedBy(dt));
                propagator2.resetInitialState(shiftState(initialState, orbitType, angleType, 2 * steps[i], i));
                SpacecraftState sP2h = propagator2.propagate(initialState.getDate().shiftedBy(dt));
                propagator2.resetInitialState(shiftState(initialState, orbitType, angleType, 3 * steps[i], i));
                SpacecraftState sP3h = propagator2.propagate(initialState.getDate().shiftedBy(dt));
                propagator2.resetInitialState(shiftState(initialState, orbitType, angleType, 4 * steps[i], i));
                SpacecraftState sP4h = propagator2.propagate(initialState.getDate().shiftedBy(dt));
                fillJacobianColumn(dYdY0Ref, i, orbitType, angleType, steps[i], sM4h, sM3h, sM2h, sM1h, sP1h, sP2h, sP3h, sP4h);
            }
            for (int i = 0; i < 6; ++i) {
                for (int j = 0; j < 6; ++j) {
                    double error = FastMath.abs((dYdY0[i][j] - dYdY0Ref[i][j]) / dYdY0Ref[i][j]);
                    Assert.assertEquals(0, error, 6.0e-2);
                }
            }
        }
    }
}

Example 30 with ForceModel

use of org.orekit.forces.ForceModel in project Orekit by CS-SI.

the class PartialDerivativesTest method testPropagationTypesHyperbolic.

@Test
public void testPropagationTypesHyperbolic() throws OrekitException {
    NormalizedSphericalHarmonicsProvider provider = GravityFieldFactory.getNormalizedProvider(5, 5);
    ForceModel gravityField = new HolmesFeatherstoneAttractionModel(FramesFactory.getITRF(IERSConventions.IERS_2010, true), provider);
    Orbit initialOrbit = new KeplerianOrbit(new PVCoordinates(new Vector3D(-1551946.0, 708899.0, 6788204.0), new Vector3D(-9875.0, -3941.0, -1845.0)), FramesFactory.getEME2000(), AbsoluteDate.J2000_EPOCH, provider.getMu());
    double dt = 900;
    double dP = 0.001;
    for (OrbitType orbitType : new OrbitType[] { OrbitType.KEPLERIAN, OrbitType.CARTESIAN }) {
        for (PositionAngle angleType : PositionAngle.values()) {
            // compute state Jacobian using PartialDerivatives
            NumericalPropagator propagator = setUpPropagator(initialOrbit, dP, orbitType, angleType, gravityField);
            PartialDerivativesEquations partials = new PartialDerivativesEquations("partials", propagator);
            final SpacecraftState initialState = partials.setInitialJacobians(new SpacecraftState(initialOrbit));
            propagator.setInitialState(initialState);
            final JacobiansMapper mapper = partials.getMapper();
            PickUpHandler pickUp = new PickUpHandler(mapper, null);
            propagator.setMasterMode(pickUp);
            propagator.propagate(initialState.getDate().shiftedBy(dt));
            double[][] dYdY0 = pickUp.getdYdY0();
            // compute reference state Jacobian using finite differences
            double[][] dYdY0Ref = new double[6][6];
            AbstractIntegratedPropagator propagator2 = setUpPropagator(initialOrbit, dP, orbitType, angleType, gravityField);
            double[] steps = NumericalPropagator.tolerances(1000000 * dP, initialOrbit, orbitType)[0];
            for (int i = 0; i < 6; ++i) {
                propagator2.resetInitialState(shiftState(initialState, orbitType, angleType, -4 * steps[i], i));
                SpacecraftState sM4h = propagator2.propagate(initialState.getDate().shiftedBy(dt));
                propagator2.resetInitialState(shiftState(initialState, orbitType, angleType, -3 * steps[i], i));
                SpacecraftState sM3h = propagator2.propagate(initialState.getDate().shiftedBy(dt));
                propagator2.resetInitialState(shiftState(initialState, orbitType, angleType, -2 * steps[i], i));
                SpacecraftState sM2h = propagator2.propagate(initialState.getDate().shiftedBy(dt));
                propagator2.resetInitialState(shiftState(initialState, orbitType, angleType, -1 * steps[i], i));
                SpacecraftState sM1h = propagator2.propagate(initialState.getDate().shiftedBy(dt));
                propagator2.resetInitialState(shiftState(initialState, orbitType, angleType, 1 * steps[i], i));
                SpacecraftState sP1h = propagator2.propagate(initialState.getDate().shiftedBy(dt));
                propagator2.resetInitialState(shiftState(initialState, orbitType, angleType, 2 * steps[i], i));
                SpacecraftState sP2h = propagator2.propagate(initialState.getDate().shiftedBy(dt));
                propagator2.resetInitialState(shiftState(initialState, orbitType, angleType, 3 * steps[i], i));
                SpacecraftState sP3h = propagator2.propagate(initialState.getDate().shiftedBy(dt));
                propagator2.resetInitialState(shiftState(initialState, orbitType, angleType, 4 * steps[i], i));
                SpacecraftState sP4h = propagator2.propagate(initialState.getDate().shiftedBy(dt));
                fillJacobianColumn(dYdY0Ref, i, orbitType, angleType, steps[i], sM4h, sM3h, sM2h, sM1h, sP1h, sP2h, sP3h, sP4h);
            }
            for (int i = 0; i < 6; ++i) {
                for (int j = 0; j < 6; ++j) {
                    double error = FastMath.abs((dYdY0[i][j] - dYdY0Ref[i][j]) / dYdY0Ref[i][j]);
                    Assert.assertEquals(0, error, 1.0e-3);
                }
            }
        }
    }
}