Examples with FieldVector3D org.hipparchus.geometry.euclidean.threed.FieldVector3D used on opensource projects

Example 96 with FieldVector3D

use of org.hipparchus.geometry.euclidean.threed.FieldVector3D in project Orekit by CS-SI.

the class AttitudesSequenceTest method doTestDayNightSwitchField.

private <T extends RealFieldElement<T>> void doTestDayNightSwitchField(final Field<T> field) throws OrekitException {
    // Initial state definition : date, orbit
    final FieldAbsoluteDate<T> initialDate = new FieldAbsoluteDate<>(field, 2004, 01, 01, 23, 30, 00.000, TimeScalesFactory.getUTC());
    final FieldVector3D<T> position = new FieldVector3D<>(field, new Vector3D(-6142438.668, 3492467.560, -25767.25680));
    final FieldVector3D<T> velocity = new FieldVector3D<>(field, new Vector3D(505.8479685, 942.7809215, 7435.922231));
    final FieldOrbit<T> initialOrbit = new FieldKeplerianOrbit<>(new FieldPVCoordinates<>(position, velocity), FramesFactory.getEME2000(), initialDate, Constants.EIGEN5C_EARTH_MU);
    // Attitudes sequence definition
    EventsLogger logger = new EventsLogger();
    final AttitudesSequence attitudesSequence = new AttitudesSequence();
    final AttitudeProvider dayObservationLaw = new LofOffset(initialOrbit.getFrame(), LOFType.VVLH, RotationOrder.XYZ, FastMath.toRadians(20), FastMath.toRadians(40), 0);
    final AttitudeProvider nightRestingLaw = new LofOffset(initialOrbit.getFrame(), LOFType.VVLH);
    final PVCoordinatesProvider sun = CelestialBodyFactory.getSun();
    final PVCoordinatesProvider earth = CelestialBodyFactory.getEarth();
    final EclipseDetector ed = new EclipseDetector(sun, 696000000., earth, Constants.WGS84_EARTH_EQUATORIAL_RADIUS).withHandler(new ContinueOnEvent<EclipseDetector>() {

        private static final long serialVersionUID = 1L;

        int count = 0;

        public EventHandler.Action eventOccurred(final SpacecraftState s, final EclipseDetector d, final boolean increasing) {
            setInEclipse(s.getDate(), !increasing);
            if (count++ == 7) {
                return Action.STOP;
            } else {
                switch(count % 3) {
                    case 0:
                        return Action.CONTINUE;
                    case 1:
                        return Action.RESET_DERIVATIVES;
                    default:
                        return Action.RESET_STATE;
                }
            }
        }
    });
    final EventDetector monitored = logger.monitorDetector(ed);
    final Handler dayToNightHandler = new Handler(dayObservationLaw, nightRestingLaw);
    final Handler nightToDayHandler = new Handler(nightRestingLaw, dayObservationLaw);
    attitudesSequence.addSwitchingCondition(dayObservationLaw, nightRestingLaw, monitored, false, true, 300.0, AngularDerivativesFilter.USE_RRA, dayToNightHandler);
    attitudesSequence.addSwitchingCondition(nightRestingLaw, dayObservationLaw, monitored, true, false, 300.0, AngularDerivativesFilter.USE_RRA, nightToDayHandler);
    FieldSpacecraftState<T> initialState = new FieldSpacecraftState<>(initialOrbit);
    initialState = initialState.addAdditionalState("fortyTwo", field.getZero().add(42.0));
    if (ed.g(initialState.toSpacecraftState()) >= 0) {
        // initial position is in daytime
        setInEclipse(initialDate.toAbsoluteDate(), false);
        attitudesSequence.resetActiveProvider(dayObservationLaw);
    } else {
        // initial position is in nighttime
        setInEclipse(initialDate.toAbsoluteDate(), true);
        attitudesSequence.resetActiveProvider(nightRestingLaw);
    }
    // Propagator : consider the analytical Eckstein-Hechler model
    final FieldPropagator<T> propagator = new FieldEcksteinHechlerPropagator<T>(initialOrbit, attitudesSequence, 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);
    // Register the switching events to the propagator
    attitudesSequence.registerSwitchEvents(field, propagator);
    propagator.setMasterMode(field.getZero().add(60.0), new FieldOrekitFixedStepHandler<T>() {

        public void handleStep(FieldSpacecraftState<T> currentState, boolean isLast) throws OrekitException {
            // the Earth position in spacecraft frame should be along spacecraft Z axis
            // during night time and away from it during day time due to roll and pitch offsets
            final FieldVector3D<T> earth = currentState.toTransform().transformPosition(Vector3D.ZERO);
            final T pointingOffset = FieldVector3D.angle(earth, Vector3D.PLUS_K);
            // the g function is the eclipse indicator, its an angle between Sun and Earth limb,
            // positive when Sun is outside of Earth limb, negative when Sun is hidden by Earth limb
            final double eclipseAngle = ed.g(currentState.toSpacecraftState());
            if (currentState.getDate().durationFrom(lastChange).getReal() > 300) {
                if (inEclipse) {
                    Assert.assertTrue(eclipseAngle <= 0);
                    Assert.assertEquals(0.0, pointingOffset.getReal(), 1.0e-6);
                } else {
                    Assert.assertTrue(eclipseAngle >= 0);
                    Assert.assertEquals(0.767215, pointingOffset.getReal(), 1.0e-6);
                }
            } else {
                // we are in transition
                Assert.assertTrue(pointingOffset.getReal() + " " + (0.767215 - pointingOffset.getReal()), pointingOffset.getReal() <= 0.7672155);
            }
        }
    });
    // Propagate from the initial date for the fixed duration
    propagator.propagate(initialDate.shiftedBy(12600.));
    // as we have 2 switch events (even if they share the same underlying event detector),
    // and these events are triggered at both eclipse entry and exit, we get 8
    // raw events on 2 orbits
    Assert.assertEquals(8, logger.getLoggedEvents().size());
    // we have 4 attitudes switch on 2 orbits, 2 of each type
    Assert.assertEquals(2, dayToNightHandler.dates.size());
    Assert.assertEquals(2, nightToDayHandler.dates.size());
}

Example 97 with FieldVector3D

use of org.hipparchus.geometry.euclidean.threed.FieldVector3D in project Orekit by CS-SI.

the class PVCoordinates method toDerivativeStructurePV.

/**
 * Transform the instance to a {@link FieldPVCoordinates}<{@link DerivativeStructure}>.
 * <p>
 * The {@link DerivativeStructure} coordinates correspond to time-derivatives up
 * to the user-specified order. As both the instance components {@link #getPosition() position},
 * {@link #getVelocity() velocity} and {@link #getAcceleration() acceleration} and the
 * {@link DerivativeStructure#getPartialDerivative(int...) derivatives} of the components
 * holds time-derivatives, there are several ways to retrieve these derivatives. If for example
 * the {@code order} is set to 2, then both {@code pv.getPosition().getX().getPartialDerivative(2)},
 * {@code pv.getVelocity().getX().getPartialDerivative(1)} and
 * {@code pv.getAcceleration().getX().getValue()} return the exact same value.
 * </p>
 * <p>
 * If derivation order is 1, the first derivative of acceleration will be computed as a
 * Keplerian-only jerk. If derivation order is 2, the second derivative of velocity (which
 * is also the first derivative of acceleration) will be computed as a Keplerian-only jerk,
 * and the second derivative of acceleration will be computed as a Keplerian-only jounce.
 * </p>
 * @param order derivation order for the vector components (must be either 0, 1 or 2)
 * @return pv coordinates with time-derivatives embedded within the coordinates
 * @exception OrekitException if the user specified order is too large
 * @since 9.2
 */
public FieldPVCoordinates<DerivativeStructure> toDerivativeStructurePV(final int order) throws OrekitException {
    final DSFactory factory;
    final DerivativeStructure x0;
    final DerivativeStructure y0;
    final DerivativeStructure z0;
    final DerivativeStructure x1;
    final DerivativeStructure y1;
    final DerivativeStructure z1;
    final DerivativeStructure x2;
    final DerivativeStructure y2;
    final DerivativeStructure z2;
    switch(order) {
        case 0:
            factory = new DSFactory(1, order);
            x0 = factory.build(position.getX());
            y0 = factory.build(position.getY());
            z0 = factory.build(position.getZ());
            x1 = factory.build(velocity.getX());
            y1 = factory.build(velocity.getY());
            z1 = factory.build(velocity.getZ());
            x2 = factory.build(acceleration.getX());
            y2 = factory.build(acceleration.getY());
            z2 = factory.build(acceleration.getZ());
            break;
        case 1:
            {
                factory = new DSFactory(1, order);
                final double r2 = position.getNormSq();
                final double r = FastMath.sqrt(r2);
                final double pvOr2 = Vector3D.dotProduct(position, velocity) / r2;
                final double a = acceleration.getNorm();
                final double aOr = a / r;
                final Vector3D keplerianJerk = new Vector3D(-3 * pvOr2, acceleration, -aOr, velocity);
                x0 = factory.build(position.getX(), velocity.getX());
                y0 = factory.build(position.getY(), velocity.getY());
                z0 = factory.build(position.getZ(), velocity.getZ());
                x1 = factory.build(velocity.getX(), acceleration.getX());
                y1 = factory.build(velocity.getY(), acceleration.getY());
                z1 = factory.build(velocity.getZ(), acceleration.getZ());
                x2 = factory.build(acceleration.getX(), keplerianJerk.getX());
                y2 = factory.build(acceleration.getY(), keplerianJerk.getY());
                z2 = factory.build(acceleration.getZ(), keplerianJerk.getZ());
                break;
            }
        case 2:
            {
                factory = new DSFactory(1, order);
                final double r2 = position.getNormSq();
                final double r = FastMath.sqrt(r2);
                final double pvOr2 = Vector3D.dotProduct(position, velocity) / r2;
                final double a = acceleration.getNorm();
                final double aOr = a / r;
                final Vector3D keplerianJerk = new Vector3D(-3 * pvOr2, acceleration, -aOr, velocity);
                final double v2 = velocity.getNormSq();
                final double pa = Vector3D.dotProduct(position, acceleration);
                final double aj = Vector3D.dotProduct(acceleration, keplerianJerk);
                final Vector3D keplerianJounce = new Vector3D(-3 * (v2 + pa) / r2 + 15 * pvOr2 * pvOr2 - aOr, acceleration, 4 * aOr * pvOr2 - aj / (a * r), velocity);
                x0 = factory.build(position.getX(), velocity.getX(), acceleration.getX());
                y0 = factory.build(position.getY(), velocity.getY(), acceleration.getY());
                z0 = factory.build(position.getZ(), velocity.getZ(), acceleration.getZ());
                x1 = factory.build(velocity.getX(), acceleration.getX(), keplerianJerk.getX());
                y1 = factory.build(velocity.getY(), acceleration.getY(), keplerianJerk.getY());
                z1 = factory.build(velocity.getZ(), acceleration.getZ(), keplerianJerk.getZ());
                x2 = factory.build(acceleration.getX(), keplerianJerk.getX(), keplerianJounce.getX());
                y2 = factory.build(acceleration.getY(), keplerianJerk.getY(), keplerianJounce.getY());
                z2 = factory.build(acceleration.getZ(), keplerianJerk.getZ(), keplerianJounce.getZ());
                break;
            }
        default:
            throw new OrekitException(OrekitMessages.OUT_OF_RANGE_DERIVATION_ORDER, order);
    }
    return new FieldPVCoordinates<>(new FieldVector3D<>(x0, y0, z0), new FieldVector3D<>(x1, y1, z1), new FieldVector3D<>(x2, y2, z2));
}

Example 98 with FieldVector3D

use of org.hipparchus.geometry.euclidean.threed.FieldVector3D in project Orekit by CS-SI.

the class FieldEcksteinHechlerPropagatorTest method doAlmostSphericalBody.

private <T extends RealFieldElement<T>> void doAlmostSphericalBody(Field<T> field) throws OrekitException {
    T zero = field.getZero();
    FieldAbsoluteDate<T> date = new FieldAbsoluteDate<>(field);
    // Definition of initial conditions
    // ---------------------------------
    // with e around e = 1.4e-4 and i = 1.7 rad
    FieldVector3D<T> position = new FieldVector3D<>(zero.add(3220103.), zero.add(69623.), zero.add(6449822.));
    FieldVector3D<T> velocity = new FieldVector3D<>(zero.add(6414.7), zero.add(-2006.), zero.add(-3180.));
    FieldAbsoluteDate<T> initDate = date.shiftedBy(584.);
    FieldOrbit<T> initialOrbit = new FieldEquinoctialOrbit<>(new FieldPVCoordinates<>(position, velocity), FramesFactory.getEME2000(), initDate, provider.getMu());
    // Initialisation to simulate a Keplerian extrapolation
    // To be noticed: in order to simulate a Keplerian extrapolation with the
    // analytical
    // extrapolator, one should put the zonal coefficients to 0. But due to
    // numerical pbs
    // one must put a non 0 value.
    UnnormalizedSphericalHarmonicsProvider kepProvider = GravityFieldFactory.getUnnormalizedProvider(6.378137e6, 3.9860047e14, TideSystem.UNKNOWN, new double[][] { { 0 }, { 0 }, { 0.1e-10 }, { 0.1e-13 }, { 0.1e-13 }, { 0.1e-14 }, { 0.1e-14 } }, new double[][] { { 0 }, { 0 }, { 0 }, { 0 }, { 0 }, { 0 }, { 0 } });
    // Extrapolators definitions
    // -------------------------
    FieldEcksteinHechlerPropagator<T> extrapolatorAna = new FieldEcksteinHechlerPropagator<>(initialOrbit, kepProvider);
    FieldKeplerianPropagator<T> extrapolatorKep = new FieldKeplerianPropagator<>(initialOrbit);
    // Extrapolation at a final date different from initial date
    // ---------------------------------------------------------
    // extrapolation duration in seconds
    double delta_t = 100.0;
    FieldAbsoluteDate<T> extrapDate = initDate.shiftedBy(delta_t);
    FieldSpacecraftState<T> finalOrbitAna = extrapolatorAna.propagate(extrapDate);
    FieldSpacecraftState<T> finalOrbitKep = extrapolatorKep.propagate(extrapDate);
    Assert.assertEquals(finalOrbitAna.getDate().durationFrom(extrapDate).getReal(), 0.0, Utils.epsilonTest);
    // comparison of each orbital parameters
    Assert.assertEquals(finalOrbitAna.getA().getReal(), finalOrbitKep.getA().getReal(), 10 * Utils.epsilonTest * finalOrbitKep.getA().getReal());
    Assert.assertEquals(finalOrbitAna.getEquinoctialEx().getReal(), finalOrbitKep.getEquinoctialEx().getReal(), Utils.epsilonE * finalOrbitKep.getE().getReal());
    Assert.assertEquals(finalOrbitAna.getEquinoctialEy().getReal(), finalOrbitKep.getEquinoctialEy().getReal(), Utils.epsilonE * finalOrbitKep.getE().getReal());
    Assert.assertEquals(MathUtils.normalizeAngle(finalOrbitAna.getHx().getReal(), finalOrbitKep.getHx().getReal()), finalOrbitKep.getHx().getReal(), Utils.epsilonAngle * FastMath.abs(finalOrbitKep.getI().getReal()));
    Assert.assertEquals(MathUtils.normalizeAngle(finalOrbitAna.getHy().getReal(), finalOrbitKep.getHy().getReal()), finalOrbitKep.getHy().getReal(), Utils.epsilonAngle * FastMath.abs(finalOrbitKep.getI().getReal()));
    Assert.assertEquals(MathUtils.normalizeAngle(finalOrbitAna.getLv().getReal(), finalOrbitKep.getLv().getReal()), finalOrbitKep.getLv().getReal(), Utils.epsilonAngle * FastMath.abs(finalOrbitKep.getLv().getReal()));
    Assert.assertEquals(MathUtils.normalizeAngle(finalOrbitAna.getLE().getReal(), finalOrbitKep.getLE().getReal()), finalOrbitKep.getLE().getReal(), Utils.epsilonAngle * FastMath.abs(finalOrbitKep.getLE().getReal()));
    Assert.assertEquals(MathUtils.normalizeAngle(finalOrbitAna.getLM().getReal(), finalOrbitKep.getLM().getReal()), finalOrbitKep.getLM().getReal(), Utils.epsilonAngle * FastMath.abs(finalOrbitKep.getLM().getReal()));
}

Example 99 with FieldVector3D

use of org.hipparchus.geometry.euclidean.threed.FieldVector3D in project Orekit by CS-SI.

the class FieldKeplerianPropagatorTest method testTuple.

@Test
public void testTuple() throws OrekitException {
    AbsoluteDate initDate = AbsoluteDate.J2000_EPOCH.shiftedBy(584.);
    KeplerianOrbit k0 = new KeplerianOrbit(7209668.0, 0.5e-4, 1.7, 2.1, 2.9, 6.2, PositionAngle.TRUE, FramesFactory.getEME2000(), initDate, mu);
    TimeStampedPVCoordinates pv0 = k0.getPVCoordinates();
    TimeStampedPVCoordinates pv1 = new TimeStampedPVCoordinates(pv0.getDate(), pv0.getPosition(), pv0.getVelocity().add(new Vector3D(2.0, pv0.getVelocity().normalize())));
    KeplerianOrbit k1 = new KeplerianOrbit(pv1, k0.getFrame(), k0.getMu());
    FieldOrbit<Tuple> twoOrbits = new FieldKeplerianOrbit<>(new Tuple(k0.getA(), k1.getA()), new Tuple(k0.getE(), k1.getE()), new Tuple(k0.getI(), k1.getI()), new Tuple(k0.getPerigeeArgument(), k1.getPerigeeArgument()), new Tuple(k0.getRightAscensionOfAscendingNode(), k1.getRightAscensionOfAscendingNode()), new Tuple(k0.getMeanAnomaly(), k1.getMeanAnomaly()), PositionAngle.MEAN, FramesFactory.getEME2000(), new FieldAbsoluteDate<>(initDate, new Tuple(0.0, 0.0)), mu);
    Field<Tuple> field = twoOrbits.getDate().getField();
    FieldPropagator<Tuple> propagator = new FieldKeplerianPropagator<>(twoOrbits);
    Min minTangential = new Min();
    Max maxTangential = new Max();
    Min minRadial = new Min();
    Max maxRadial = new Max();
    propagator.setMasterMode(field.getZero().add(10.0), (s, isLast) -> {
        FieldVector3D<Tuple> p = s.getPVCoordinates().getPosition();
        FieldVector3D<Tuple> v = s.getPVCoordinates().getVelocity();
        Vector3D p0 = new Vector3D(p.getX().getComponent(0), p.getY().getComponent(0), p.getZ().getComponent(0));
        Vector3D v0 = new Vector3D(v.getX().getComponent(0), v.getY().getComponent(0), v.getZ().getComponent(0));
        Vector3D t = v0.normalize();
        Vector3D r = Vector3D.crossProduct(v0, Vector3D.crossProduct(p0, v0)).normalize();
        Vector3D p1 = new Vector3D(p.getX().getComponent(1), p.getY().getComponent(1), p.getZ().getComponent(1));
        double dT = Vector3D.dotProduct(p1.subtract(p0), t);
        double dR = Vector3D.dotProduct(p1.subtract(p0), r);
        minTangential.increment(dT);
        maxTangential.increment(dT);
        minRadial.increment(dR);
        maxRadial.increment(dR);
    });
    propagator.propagate(twoOrbits.getDate().shiftedBy(Constants.JULIAN_DAY / 8));
    Assert.assertEquals(-72525.685, minTangential.getResult(), 1.0e-3);
    Assert.assertEquals(926.046, maxTangential.getResult(), 1.0e-3);
    Assert.assertEquals(-92.800, minRadial.getResult(), 1.0e-3);
    Assert.assertEquals(7739.532, maxRadial.getResult(), 1.0e-3);
}

Example 100 with FieldVector3D

use of org.hipparchus.geometry.euclidean.threed.FieldVector3D in project Orekit by CS-SI.

the class FieldKeplerianPropagatorTest method doTestNoDerivatives.

private <T extends RealFieldElement<T>> void doTestNoDerivatives(Field<T> field) throws OrekitException {
    T zero = field.getZero();
    for (OrbitType type : OrbitType.values()) {
        // create an initial orbit with non-Keplerian acceleration
        final FieldAbsoluteDate<T> date = new FieldAbsoluteDate<>(field, 2003, 9, 16, 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(505.848), zero.add(942.781), zero.add(7435.922));
        final FieldVector3D<T> keplerAcceleration = new FieldVector3D<>(position.getNormSq().reciprocal().multiply(-mu), position.normalize());
        final FieldVector3D<T> nonKeplerAcceleration = new FieldVector3D<>(zero.add(0.001), zero.add(0.002), zero.add(0.003));
        final FieldVector3D<T> acceleration = keplerAcceleration.add(nonKeplerAcceleration);
        final TimeStampedFieldPVCoordinates<T> pva = new TimeStampedFieldPVCoordinates<>(date, position, velocity, acceleration);
        final FieldOrbit<T> initial = type.convertType(new FieldCartesianOrbit<>(pva, FramesFactory.getEME2000(), mu));
        Assert.assertEquals(type, initial.getType());
        // the derivatives are available at this stage
        checkDerivatives(initial, true);
        FieldKeplerianPropagator<T> propagator = new FieldKeplerianPropagator<>(initial);
        Assert.assertEquals(type, propagator.getInitialState().getOrbit().getType());
        // non-Keplerian derivatives are explicitly removed when building the Keplerian-only propagator
        checkDerivatives(propagator.getInitialState().getOrbit(), false);
        FieldPVCoordinates<T> initPV = propagator.getInitialState().getOrbit().getPVCoordinates();
        Assert.assertEquals(nonKeplerAcceleration.getNorm().getReal(), FieldVector3D.distance(acceleration, initPV.getAcceleration()).getReal(), 2.0e-15);
        Assert.assertEquals(0.0, FieldVector3D.distance(keplerAcceleration, initPV.getAcceleration()).getReal(), 5.0e-15);
        T dt = initial.getKeplerianPeriod().multiply(0.2);
        FieldOrbit<T> orbit = propagator.propagateOrbit(initial.getDate().shiftedBy(dt));
        Assert.assertEquals(type, orbit.getType());
        // at the end, we don't have non-Keplerian derivatives
        checkDerivatives(orbit, false);
        // using shiftedBy on the initial orbit, non-Keplerian derivatives would have been preserved
        checkDerivatives(initial.shiftedBy(dt), true);
    }
}