Examples with FieldCartesianOrbit org.orekit.orbits.FieldCartesianOrbit used on opensource projects

Example 6 with FieldCartesianOrbit

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

the class AbstractForceModelTest method checkStateJacobianVsFiniteDifferences.

protected void checkStateJacobianVsFiniteDifferences(final SpacecraftState state0, final ForceModel forceModel, final AttitudeProvider provider, final double dP, final double checkTolerance, final boolean print) throws OrekitException {
    double[][] finiteDifferencesJacobian = Differentiation.differentiate(state -> forceModel.acceleration(state, forceModel.getParameters()).toArray(), 3, provider, OrbitType.CARTESIAN, PositionAngle.MEAN, dP, 5).value(state0);
    DSFactory factory = new DSFactory(6, 1);
    Field<DerivativeStructure> field = factory.getDerivativeField();
    final FieldAbsoluteDate<DerivativeStructure> fDate = new FieldAbsoluteDate<>(field, state0.getDate());
    final Vector3D p = state0.getPVCoordinates().getPosition();
    final Vector3D v = state0.getPVCoordinates().getVelocity();
    final Vector3D a = state0.getPVCoordinates().getAcceleration();
    final TimeStampedFieldPVCoordinates<DerivativeStructure> fPVA = new TimeStampedFieldPVCoordinates<>(fDate, new FieldVector3D<>(factory.variable(0, p.getX()), factory.variable(1, p.getY()), factory.variable(2, p.getZ())), new FieldVector3D<>(factory.variable(3, v.getX()), factory.variable(4, v.getY()), factory.variable(5, v.getZ())), new FieldVector3D<>(factory.constant(a.getX()), factory.constant(a.getY()), factory.constant(a.getZ())));
    final TimeStampedFieldAngularCoordinates<DerivativeStructure> fAC = new TimeStampedFieldAngularCoordinates<>(fDate, new FieldRotation<>(field, state0.getAttitude().getRotation()), new FieldVector3D<>(field, state0.getAttitude().getSpin()), new FieldVector3D<>(field, state0.getAttitude().getRotationAcceleration()));
    final FieldSpacecraftState<DerivativeStructure> fState = new FieldSpacecraftState<>(new FieldCartesianOrbit<>(fPVA, state0.getFrame(), state0.getMu()), new FieldAttitude<>(state0.getFrame(), fAC), field.getZero().add(state0.getMass()));
    FieldVector3D<DerivativeStructure> dsJacobian = forceModel.acceleration(fState, forceModel.getParameters(fState.getDate().getField()));
    Vector3D dFdPXRef = new Vector3D(finiteDifferencesJacobian[0][0], finiteDifferencesJacobian[1][0], finiteDifferencesJacobian[2][0]);
    Vector3D dFdPXRes = new Vector3D(dsJacobian.getX().getPartialDerivative(1, 0, 0, 0, 0, 0), dsJacobian.getY().getPartialDerivative(1, 0, 0, 0, 0, 0), dsJacobian.getZ().getPartialDerivative(1, 0, 0, 0, 0, 0));
    Vector3D dFdPYRef = new Vector3D(finiteDifferencesJacobian[0][1], finiteDifferencesJacobian[1][1], finiteDifferencesJacobian[2][1]);
    Vector3D dFdPYRes = new Vector3D(dsJacobian.getX().getPartialDerivative(0, 1, 0, 0, 0, 0), dsJacobian.getY().getPartialDerivative(0, 1, 0, 0, 0, 0), dsJacobian.getZ().getPartialDerivative(0, 1, 0, 0, 0, 0));
    Vector3D dFdPZRef = new Vector3D(finiteDifferencesJacobian[0][2], finiteDifferencesJacobian[1][2], finiteDifferencesJacobian[2][2]);
    Vector3D dFdPZRes = new Vector3D(dsJacobian.getX().getPartialDerivative(0, 0, 1, 0, 0, 0), dsJacobian.getY().getPartialDerivative(0, 0, 1, 0, 0, 0), dsJacobian.getZ().getPartialDerivative(0, 0, 1, 0, 0, 0));
    Vector3D dFdVXRef = new Vector3D(finiteDifferencesJacobian[0][3], finiteDifferencesJacobian[1][3], finiteDifferencesJacobian[2][3]);
    Vector3D dFdVXRes = new Vector3D(dsJacobian.getX().getPartialDerivative(0, 0, 0, 1, 0, 0), dsJacobian.getY().getPartialDerivative(0, 0, 0, 1, 0, 0), dsJacobian.getZ().getPartialDerivative(0, 0, 0, 1, 0, 0));
    Vector3D dFdVYRef = new Vector3D(finiteDifferencesJacobian[0][4], finiteDifferencesJacobian[1][4], finiteDifferencesJacobian[2][4]);
    Vector3D dFdVYRes = new Vector3D(dsJacobian.getX().getPartialDerivative(0, 0, 0, 0, 1, 0), dsJacobian.getY().getPartialDerivative(0, 0, 0, 0, 1, 0), dsJacobian.getZ().getPartialDerivative(0, 0, 0, 0, 1, 0));
    Vector3D dFdVZRef = new Vector3D(finiteDifferencesJacobian[0][5], finiteDifferencesJacobian[1][5], finiteDifferencesJacobian[2][5]);
    Vector3D dFdVZRes = new Vector3D(dsJacobian.getX().getPartialDerivative(0, 0, 0, 0, 0, 1), dsJacobian.getY().getPartialDerivative(0, 0, 0, 0, 0, 1), dsJacobian.getZ().getPartialDerivative(0, 0, 0, 0, 0, 1));
    if (print) {
        System.out.println("dF/dPX ref: " + dFdPXRef.getX() + " " + dFdPXRef.getY() + " " + dFdPXRef.getZ());
        System.out.println("dF/dPX res: " + dFdPXRes.getX() + " " + dFdPXRes.getY() + " " + dFdPXRes.getZ());
        System.out.println("dF/dPY ref: " + dFdPYRef.getX() + " " + dFdPYRef.getY() + " " + dFdPYRef.getZ());
        System.out.println("dF/dPY res: " + dFdPYRes.getX() + " " + dFdPYRes.getY() + " " + dFdPYRes.getZ());
        System.out.println("dF/dPZ ref: " + dFdPZRef.getX() + " " + dFdPZRef.getY() + " " + dFdPZRef.getZ());
        System.out.println("dF/dPZ res: " + dFdPZRes.getX() + " " + dFdPZRes.getY() + " " + dFdPZRes.getZ());
        System.out.println("dF/dPX ref norm: " + dFdPXRef.getNorm() + ", abs error: " + Vector3D.distance(dFdPXRef, dFdPXRes) + ", rel error: " + (Vector3D.distance(dFdPXRef, dFdPXRes) / dFdPXRef.getNorm()));
        System.out.println("dF/dPY ref norm: " + dFdPYRef.getNorm() + ", abs error: " + Vector3D.distance(dFdPYRef, dFdPYRes) + ", rel error: " + (Vector3D.distance(dFdPYRef, dFdPYRes) / dFdPYRef.getNorm()));
        System.out.println("dF/dPZ ref norm: " + dFdPZRef.getNorm() + ", abs error: " + Vector3D.distance(dFdPZRef, dFdPZRes) + ", rel error: " + (Vector3D.distance(dFdPZRef, dFdPZRes) / dFdPZRef.getNorm()));
        System.out.println("dF/dVX ref norm: " + dFdVXRef.getNorm() + ", abs error: " + Vector3D.distance(dFdVXRef, dFdVXRes) + ", rel error: " + (Vector3D.distance(dFdVXRef, dFdVXRes) / dFdVXRef.getNorm()));
        System.out.println("dF/dVY ref norm: " + dFdVYRef.getNorm() + ", abs error: " + Vector3D.distance(dFdVYRef, dFdVYRes) + ", rel error: " + (Vector3D.distance(dFdVYRef, dFdVYRes) / dFdVYRef.getNorm()));
        System.out.println("dF/dVZ ref norm: " + dFdVZRef.getNorm() + ", abs error: " + Vector3D.distance(dFdVZRef, dFdVZRes) + ", rel error: " + (Vector3D.distance(dFdVZRef, dFdVZRes) / dFdVZRef.getNorm()));
    }
    checkdFdP(dFdPXRef, dFdPXRes, checkTolerance);
    checkdFdP(dFdPYRef, dFdPYRes, checkTolerance);
    checkdFdP(dFdPZRef, dFdPZRes, checkTolerance);
    checkdFdP(dFdVXRef, dFdVXRes, checkTolerance);
    checkdFdP(dFdVYRef, dFdVYRes, checkTolerance);
    checkdFdP(dFdVZRef, dFdVZRes, checkTolerance);
}

Example 7 with FieldCartesianOrbit

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

the class FieldNumericalPropagatorTest method createHyperbolicOrbit.

private static <T extends RealFieldElement<T>> FieldCartesianOrbit<T> createHyperbolicOrbit(Field<T> field) throws OrekitException {
    T zero = field.getZero();
    final FieldAbsoluteDate<T> date = new FieldAbsoluteDate<>(field, "2003-05-01T00:00:20.000", TimeScalesFactory.getUTC());
    final FieldVector3D<T> position = new FieldVector3D<>(zero.add(224267911.905821), zero.add(290251613.109399), zero.add(45534292.777492));
    final FieldVector3D<T> velocity = new FieldVector3D<>(zero.add(-1494.068165293), zero.add(1124.771027677), zero.add(526.915286134));
    final TimeStampedFieldPVCoordinates<T> pv = new TimeStampedFieldPVCoordinates<>(date, position, velocity, FieldVector3D.getZero(field));
    final Frame frame = FramesFactory.getEME2000();
    final double mu = Constants.EIGEN5C_EARTH_MU;
    return new FieldCartesianOrbit<>(pv, frame, mu);
}

Example 8 with FieldCartesianOrbit

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

the class FieldNumericalPropagatorTest method doTestEventDetectionBug.

private <T extends RealFieldElement<T>> void doTestEventDetectionBug(final Field<T> field) throws OrekitException {
    T zero = field.getZero();
    TimeScale utc = TimeScalesFactory.getUTC();
    FieldAbsoluteDate<T> initialDate = new FieldAbsoluteDate<>(field, 2005, 1, 1, 0, 0, 0.0, utc);
    T duration = zero.add(100000.0);
    FieldAbsoluteDate<T> endDate = new FieldAbsoluteDate<>(initialDate, duration);
    // Initialization of the frame EME2000
    Frame EME2000 = FramesFactory.getEME2000();
    // Initial orbit
    double a = 35786000. + 6378137.0;
    double e = 0.70;
    double rApogee = a * (1 + e);
    double vApogee = FastMath.sqrt(mu * (1 - e) / (a * (1 + e)));
    FieldOrbit<T> geo = new FieldCartesianOrbit<>(new FieldPVCoordinates<>(new FieldVector3D<>(zero.add(rApogee), zero, zero), new FieldVector3D<>(zero, zero.add(vApogee), zero)), EME2000, initialDate, mu);
    duration = geo.getKeplerianPeriod();
    endDate = new FieldAbsoluteDate<>(initialDate, duration);
    // Numerical Integration
    final double minStep = 0.001;
    final double maxStep = 1000;
    final double initStep = 60;
    final OrbitType type = OrbitType.EQUINOCTIAL;
    final double[] absTolerance = { 0.001, 1.0e-9, 1.0e-9, 1.0e-6, 1.0e-6, 1.0e-6, 0.001 };
    final double[] relTolerance = { 1.0e-7, 1.0e-4, 1.0e-4, 1.0e-7, 1.0e-7, 1.0e-7, 1.0e-7 };
    AdaptiveStepsizeFieldIntegrator<T> integrator = new DormandPrince853FieldIntegrator<>(field, minStep, maxStep, absTolerance, relTolerance);
    integrator.setInitialStepSize(zero.add(initStep));
    // Numerical propagator based on the integrator
    FieldNumericalPropagator<T> propagator = new FieldNumericalPropagator<>(field, integrator);
    propagator.setOrbitType(type);
    T mass = field.getZero().add(1000.0);
    FieldSpacecraftState<T> initialState = new FieldSpacecraftState<>(geo, mass);
    propagator.setInitialState(initialState);
    propagator.setOrbitType(OrbitType.CARTESIAN);
    // Set the events Detectors
    FieldApsideDetector<T> event1 = new FieldApsideDetector<>(geo);
    propagator.addEventDetector(event1);
    // Set the propagation mode
    propagator.setSlaveMode();
    // Propagate
    FieldSpacecraftState<T> finalState = propagator.propagate(endDate);
    // we should stop long before endDate
    Assert.assertTrue(endDate.durationFrom(finalState.getDate()).getReal() > 40000.0);
}

Example 9 with FieldCartesianOrbit

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

the class FieldNumericalPropagatorTest method createEllipticOrbit.

private static <T extends RealFieldElement<T>> FieldCartesianOrbit<T> createEllipticOrbit(Field<T> field) throws OrekitException {
    T zero = field.getZero();
    final FieldAbsoluteDate<T> date = new FieldAbsoluteDate<>(field, "2003-05-01T00:00:20.000", TimeScalesFactory.getUTC());
    final FieldVector3D<T> position = new FieldVector3D<>(zero.add(6896874.444705), zero.add(1956581.072644), zero.add(-147476.245054));
    final FieldVector3D<T> velocity = new FieldVector3D<>(zero.add(166.816407662), zero.add(-1106.783301861), zero.add(-7372.745712770));
    final TimeStampedFieldPVCoordinates<T> pv = new TimeStampedFieldPVCoordinates<>(date, position, velocity, FieldVector3D.getZero(field));
    final Frame frame = FramesFactory.getEME2000();
    final double mu = Constants.EIGEN5C_EARTH_MU;
    return new FieldCartesianOrbit<>(pv, frame, mu);
}

Example 10 with FieldCartesianOrbit

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

the class FieldApsideDetectorTest method doTestSimple.

private <T extends RealFieldElement<T>> void doTestSimple(Field<T> field) throws OrekitException {
    final T zero = field.getZero();
    final TimeScale utc = TimeScalesFactory.getUTC();
    final FieldVector3D<T> position = new FieldVector3D<>(zero.add(-6142438.668), zero.add(3492467.56), zero.add(-25767.257));
    final FieldVector3D<T> velocity = new FieldVector3D<>(zero.add(506.0), zero.add(943.0), zero.add(7450));
    final FieldAbsoluteDate<T> date = new FieldAbsoluteDate<>(field, 2003, 9, 16, utc);
    final FieldOrbit<T> orbit = new FieldCartesianOrbit<>(new FieldPVCoordinates<>(position, velocity), FramesFactory.getEME2000(), date, Constants.EIGEN5C_EARTH_MU);
    FieldEcksteinHechlerPropagator<T> propagator = new FieldEcksteinHechlerPropagator<>(orbit, Constants.EIGEN5C_EARTH_EQUATORIAL_RADIUS, Constants.EIGEN5C_EARTH_MU, Constants.EIGEN5C_EARTH_C20, Constants.EIGEN5C_EARTH_C30, Constants.EIGEN5C_EARTH_C40, Constants.EIGEN5C_EARTH_C50, Constants.EIGEN5C_EARTH_C60);
    FieldEventDetector<T> detector = new FieldApsideDetector<>(propagator.getInitialState().getOrbit()).withMaxCheck(zero.add(600.0)).withThreshold(zero.add(1.0e-12)).withHandler(new FieldContinueOnEvent<FieldApsideDetector<T>, T>());
    Assert.assertEquals(600.0, detector.getMaxCheckInterval().getReal(), 1.0e-15);
    Assert.assertEquals(1.0e-12, detector.getThreshold().getReal(), 1.0e-15);
    Assert.assertEquals(AbstractDetector.DEFAULT_MAX_ITER, detector.getMaxIterationCount());
    FieldEventsLogger<T> logger = new FieldEventsLogger<>();
    propagator.addEventDetector(logger.monitorDetector(detector));
    propagator.propagate(propagator.getInitialState().getOrbit().getDate().shiftedBy(Constants.JULIAN_DAY));
    Assert.assertEquals(30, logger.getLoggedEvents().size());
    for (FieldLoggedEvent<T> e : logger.getLoggedEvents()) {
        FieldKeplerianOrbit<T> o = (FieldKeplerianOrbit<T>) OrbitType.KEPLERIAN.convertType(e.getState().getOrbit());
        double expected = e.isIncreasing() ? 0.0 : FastMath.PI;
        Assert.assertEquals(expected, MathUtils.normalizeAngle(o.getMeanAnomaly().getReal(), expected), 4.0e-14);
    }
}