Search in sources :

Example 51 with Transform

use of org.orekit.frames.Transform in project Orekit by CS-SI.

the class AngularRaDecMeasurementCreator method handleStep.

public void handleStep(final SpacecraftState currentState, final boolean isLast) throws OrekitException {
    for (final GroundStation station : context.stations) {
        final AbsoluteDate date = currentState.getDate();
        final Frame inertial = currentState.getFrame();
        final Vector3D position = currentState.getPVCoordinates().getPosition();
        if (station.getBaseFrame().getElevation(position, inertial, date) > FastMath.toRadians(30.0)) {
            final UnivariateSolver solver = new BracketingNthOrderBrentSolver(1.0e-12, 5);
            final double downLinkDelay = solver.solve(1000, new UnivariateFunction() {

                public double value(final double x) throws OrekitExceptionWrapper {
                    try {
                        final Transform t = station.getOffsetToInertial(inertial, date.shiftedBy(x));
                        final double d = Vector3D.distance(position, t.transformPosition(Vector3D.ZERO));
                        return d - x * Constants.SPEED_OF_LIGHT;
                    } catch (OrekitException oe) {
                        throw new OrekitExceptionWrapper(oe);
                    }
                }
            }, -1.0, 1.0);
            // Satellite position at signal departure
            final Vector3D satelliteAtDeparture = currentState.shiftedBy(-downLinkDelay).getPVCoordinates().getPosition();
            // Initialize measurement
            final double[] angular = new double[2];
            final double[] sigma = { 1.0, 1.0 };
            final double[] baseweight = { 10.0, 10.0 };
            // Inertial frame used
            final Frame inertialFrame = context.initialOrbit.getFrame();
            // Station position at signal arrival
            // Set the reference date of offset drivers arbitrarily to avoid exception
            station.getPrimeMeridianOffsetDriver().setReferenceDate(date);
            station.getPolarOffsetXDriver().setReferenceDate(date);
            station.getPolarOffsetYDriver().setReferenceDate(date);
            final Transform offsetToInertialArrival = station.getOffsetToInertial(inertialFrame, date);
            final Vector3D stationPArrival = offsetToInertialArrival.transformVector(Vector3D.ZERO);
            // Vector station position at signal arrival - satellite at signal departure
            // In inertial reference frame
            final Vector3D staSat = satelliteAtDeparture.subtract(stationPArrival);
            // Compute measurement
            // Right ascension
            final double baseRightAscension = staSat.getAlpha();
            angular[0] = MathUtils.normalizeAngle(baseRightAscension, 0.0);
            // Declination
            angular[1] = staSat.getDelta();
            addMeasurement(new AngularRaDec(station, inertialFrame, date, angular, sigma, baseweight));
        }
    }
}
Also used : Frame(org.orekit.frames.Frame) OrekitExceptionWrapper(org.orekit.errors.OrekitExceptionWrapper) UnivariateFunction(org.hipparchus.analysis.UnivariateFunction) UnivariateSolver(org.hipparchus.analysis.solvers.UnivariateSolver) AbsoluteDate(org.orekit.time.AbsoluteDate) BracketingNthOrderBrentSolver(org.hipparchus.analysis.solvers.BracketingNthOrderBrentSolver) Vector3D(org.hipparchus.geometry.euclidean.threed.Vector3D) OrekitException(org.orekit.errors.OrekitException) Transform(org.orekit.frames.Transform)

Example 52 with Transform

use of org.orekit.frames.Transform in project Orekit by CS-SI.

the class DragForceTest method accelerationDerivatives.

@Override
protected FieldVector3D<DerivativeStructure> accelerationDerivatives(final ForceModel forceModel, final AbsoluteDate date, final Frame frame, final FieldVector3D<DerivativeStructure> position, final FieldVector3D<DerivativeStructure> velocity, final FieldRotation<DerivativeStructure> rotation, final DerivativeStructure mass) throws OrekitException {
    try {
        java.lang.reflect.Field atmosphereField = DragForce.class.getDeclaredField("atmosphere");
        atmosphereField.setAccessible(true);
        Atmosphere atmosphere = (Atmosphere) atmosphereField.get(forceModel);
        java.lang.reflect.Field spacecraftField = DragForce.class.getDeclaredField("spacecraft");
        spacecraftField.setAccessible(true);
        DragSensitive spacecraft = (DragSensitive) spacecraftField.get(forceModel);
        // retrieve derivation properties
        final DSFactory factory = mass.getFactory();
        // get atmosphere properties in atmosphere own frame
        final Frame atmFrame = atmosphere.getFrame();
        final Transform toBody = frame.getTransformTo(atmFrame, date);
        final FieldVector3D<DerivativeStructure> posBodyDS = toBody.transformPosition(position);
        final Vector3D posBody = posBodyDS.toVector3D();
        final Vector3D vAtmBody = atmosphere.getVelocity(date, posBody, atmFrame);
        // to the Atmosphere interface
        if (factory.getCompiler().getOrder() > 1) {
            throw new OrekitException(OrekitMessages.OUT_OF_RANGE_DERIVATION_ORDER, factory.getCompiler().getOrder());
        }
        final double delta = 1.0;
        final double x = posBody.getX();
        final double y = posBody.getY();
        final double z = posBody.getZ();
        final double rho0 = atmosphere.getDensity(date, posBody, atmFrame);
        final double dRhodX = (atmosphere.getDensity(date, new Vector3D(x + delta, y, z), atmFrame) - rho0) / delta;
        final double dRhodY = (atmosphere.getDensity(date, new Vector3D(x, y + delta, z), atmFrame) - rho0) / delta;
        final double dRhodZ = (atmosphere.getDensity(date, new Vector3D(x, y, z + delta), atmFrame) - rho0) / delta;
        final double[] dXdQ = posBodyDS.getX().getAllDerivatives();
        final double[] dYdQ = posBodyDS.getY().getAllDerivatives();
        final double[] dZdQ = posBodyDS.getZ().getAllDerivatives();
        final double[] rhoAll = new double[dXdQ.length];
        rhoAll[0] = rho0;
        for (int i = 1; i < rhoAll.length; ++i) {
            rhoAll[i] = dRhodX * dXdQ[i] + dRhodY * dYdQ[i] + dRhodZ * dZdQ[i];
        }
        final DerivativeStructure rho = factory.build(rhoAll);
        // we consider that at first order the atmosphere velocity in atmosphere frame
        // does not depend on local position; however atmosphere velocity in inertial
        // frame DOES depend on position since the transform between the frames depends
        // on it, due to central body rotation rate and velocity composition.
        // So we use the transform to get the correct partial derivatives on vAtm
        final FieldVector3D<DerivativeStructure> vAtmBodyDS = new FieldVector3D<>(factory.constant(vAtmBody.getX()), factory.constant(vAtmBody.getY()), factory.constant(vAtmBody.getZ()));
        final FieldPVCoordinates<DerivativeStructure> pvAtmBody = new FieldPVCoordinates<>(posBodyDS, vAtmBodyDS);
        final FieldPVCoordinates<DerivativeStructure> pvAtm = toBody.getInverse().transformPVCoordinates(pvAtmBody);
        // now we can compute relative velocity, it takes into account partial derivatives with respect to position
        final FieldVector3D<DerivativeStructure> relativeVelocity = pvAtm.getVelocity().subtract(velocity);
        // compute acceleration with all its partial derivatives
        return spacecraft.dragAcceleration(new FieldAbsoluteDate<>(factory.getDerivativeField(), date), frame, position, rotation, mass, rho, relativeVelocity, forceModel.getParameters(factory.getDerivativeField()));
    } catch (IllegalArgumentException | IllegalAccessException | NoSuchFieldException | SecurityException e) {
        return null;
    }
}
Also used : Frame(org.orekit.frames.Frame) DerivativeStructure(org.hipparchus.analysis.differentiation.DerivativeStructure) DSFactory(org.hipparchus.analysis.differentiation.DSFactory) FieldVector3D(org.hipparchus.geometry.euclidean.threed.FieldVector3D) FieldVector3D(org.hipparchus.geometry.euclidean.threed.FieldVector3D) Vector3D(org.hipparchus.geometry.euclidean.threed.Vector3D) SimpleExponentialAtmosphere(org.orekit.forces.drag.atmosphere.SimpleExponentialAtmosphere) Atmosphere(org.orekit.forces.drag.atmosphere.Atmosphere) FieldPVCoordinates(org.orekit.utils.FieldPVCoordinates) OrekitException(org.orekit.errors.OrekitException) Transform(org.orekit.frames.Transform)

Example 53 with Transform

use of org.orekit.frames.Transform in project Orekit by CS-SI.

the class SimpleExponentialAtmosphereTest method testExpAtmosphere.

@Test
public void testExpAtmosphere() throws OrekitException {
    Vector3D posInEME2000 = new Vector3D(10000, Vector3D.PLUS_I);
    AbsoluteDate date = AbsoluteDate.J2000_EPOCH;
    Frame itrf = FramesFactory.getITRF(IERSConventions.IERS_2010, true);
    SimpleExponentialAtmosphere atm = new SimpleExponentialAtmosphere(new OneAxisEllipsoid(Utils.ae, 1.0 / 298.257222101, itrf), 0.0004, 42000.0, 7500.0);
    Vector3D vel = atm.getVelocity(date, posInEME2000, FramesFactory.getEME2000());
    Transform toBody = FramesFactory.getEME2000().getTransformTo(itrf, date);
    Vector3D test = Vector3D.crossProduct(toBody.getRotationRate(), posInEME2000);
    test = test.subtract(vel);
    Assert.assertEquals(0, test.getNorm(), 2.9e-5);
}
Also used : Frame(org.orekit.frames.Frame) OneAxisEllipsoid(org.orekit.bodies.OneAxisEllipsoid) Vector3D(org.hipparchus.geometry.euclidean.threed.Vector3D) Transform(org.orekit.frames.Transform) AbsoluteDate(org.orekit.time.AbsoluteDate) SimpleExponentialAtmosphere(org.orekit.forces.drag.atmosphere.SimpleExponentialAtmosphere) Test(org.junit.Test)

Example 54 with Transform

use of org.orekit.frames.Transform in project Orekit by CS-SI.

the class HolmesFeatherstoneAttractionModelTest method accelerationDerivatives.

@Override
protected FieldVector3D<DerivativeStructure> accelerationDerivatives(final ForceModel forceModel, final AbsoluteDate date, final Frame frame, final FieldVector3D<DerivativeStructure> position, final FieldVector3D<DerivativeStructure> velocity, final FieldRotation<DerivativeStructure> rotation, final DerivativeStructure mass) throws OrekitException {
    try {
        java.lang.reflect.Field bodyFrameField = HolmesFeatherstoneAttractionModel.class.getDeclaredField("bodyFrame");
        bodyFrameField.setAccessible(true);
        Frame bodyFrame = (Frame) bodyFrameField.get(forceModel);
        // get the position in body frame
        final Transform fromBodyFrame = bodyFrame.getTransformTo(frame, date);
        final Transform toBodyFrame = fromBodyFrame.getInverse();
        final Vector3D positionBody = toBodyFrame.transformPosition(position.toVector3D());
        // compute gradient and Hessian
        final GradientHessian gh = gradientHessian((HolmesFeatherstoneAttractionModel) forceModel, date, positionBody);
        // gradient of the non-central part of the gravity field
        final double[] gInertial = fromBodyFrame.transformVector(new Vector3D(gh.getGradient())).toArray();
        // Hessian of the non-central part of the gravity field
        final RealMatrix hBody = new Array2DRowRealMatrix(gh.getHessian(), false);
        final RealMatrix rot = new Array2DRowRealMatrix(toBodyFrame.getRotation().getMatrix());
        final RealMatrix hInertial = rot.transpose().multiply(hBody).multiply(rot);
        // distribute all partial derivatives in a compact acceleration vector
        final double[] derivatives = new double[1 + mass.getFreeParameters()];
        final DerivativeStructure[] accDer = new DerivativeStructure[3];
        for (int i = 0; i < 3; ++i) {
            // first element is value of acceleration (i.e. gradient of field)
            derivatives[0] = gInertial[i];
            // next three elements are one row of the Jacobian of acceleration (i.e. Hessian of field)
            derivatives[1] = hInertial.getEntry(i, 0);
            derivatives[2] = hInertial.getEntry(i, 1);
            derivatives[3] = hInertial.getEntry(i, 2);
            // next elements (three or four depending on mass being used or not) are left as 0
            accDer[i] = mass.getFactory().build(derivatives);
        }
        return new FieldVector3D<>(accDer);
    } catch (IllegalArgumentException | IllegalAccessException | NoSuchFieldException | SecurityException e) {
        return null;
    }
}
Also used : Frame(org.orekit.frames.Frame) DerivativeStructure(org.hipparchus.analysis.differentiation.DerivativeStructure) FieldVector3D(org.hipparchus.geometry.euclidean.threed.FieldVector3D) RealMatrix(org.hipparchus.linear.RealMatrix) Array2DRowRealMatrix(org.hipparchus.linear.Array2DRowRealMatrix) Array2DRowRealMatrix(org.hipparchus.linear.Array2DRowRealMatrix) FieldVector3D(org.hipparchus.geometry.euclidean.threed.FieldVector3D) Vector3D(org.hipparchus.geometry.euclidean.threed.Vector3D) Transform(org.orekit.frames.Transform)

Example 55 with Transform

use of org.orekit.frames.Transform in project Orekit by CS-SI.

the class RangeAnalytic method theoreticalEvaluationValidation.

/**
 * Added for validation
 * Compares directly numeric and analytic computations
 * @param iteration
 * @param evaluation
 * @param state
 * @return
 * @throws OrekitException
 */
protected EstimatedMeasurement<Range> theoreticalEvaluationValidation(final int iteration, final int evaluation, final SpacecraftState state) throws OrekitException {
    // Station & DSFactory attributes from parent Range class
    final GroundStation groundStation = getStation();
    // get the number of parameters used for derivation
    int nbParams = 6;
    final Map<String, Integer> indices = new HashMap<>();
    for (ParameterDriver driver : getParametersDrivers()) {
        if (driver.isSelected()) {
            indices.put(driver.getName(), nbParams++);
        }
    }
    final DSFactory dsFactory = new DSFactory(nbParams, 1);
    final Field<DerivativeStructure> field = dsFactory.getDerivativeField();
    final FieldVector3D<DerivativeStructure> zero = FieldVector3D.getZero(field);
    // Range derivatives are computed with respect to spacecraft state in inertial frame
    // and station position in station's offset frame
    // -------
    // 
    // Parameters:
    // - 0..2 - Px, Py, Pz   : Position of the spacecraft in inertial frame
    // - 3..5 - Vx, Vy, Vz   : Velocity of the spacecraft in inertial frame
    // - 6..8 - QTx, QTy, QTz: Position of the station in station's offset frame
    // Coordinates of the spacecraft expressed as a derivative structure
    final TimeStampedFieldPVCoordinates<DerivativeStructure> pvaDS = getCoordinates(state, 0, dsFactory);
    // transform between station and inertial frame, expressed as a derivative structure
    // The components of station's position in offset frame are the 3 last derivative parameters
    final AbsoluteDate downlinkDate = getDate();
    final FieldAbsoluteDate<DerivativeStructure> downlinkDateDS = new FieldAbsoluteDate<>(field, downlinkDate);
    final FieldTransform<DerivativeStructure> offsetToInertialDownlink = groundStation.getOffsetToInertial(state.getFrame(), downlinkDateDS, dsFactory, indices);
    // Station position in inertial frame at end of the downlink leg
    final TimeStampedFieldPVCoordinates<DerivativeStructure> stationDownlink = offsetToInertialDownlink.transformPVCoordinates(new TimeStampedFieldPVCoordinates<>(downlinkDateDS, zero, zero, zero));
    // Compute propagation times
    // (if state has already been set up to pre-compensate propagation delay,
    // we will have offset == downlinkDelay and transitState will be
    // the same as state)
    // Downlink delay
    final DerivativeStructure tauD = signalTimeOfFlight(pvaDS, stationDownlink.getPosition(), downlinkDateDS);
    // Transit state
    final double delta = downlinkDate.durationFrom(state.getDate());
    final DerivativeStructure tauDMDelta = tauD.negate().add(delta);
    final SpacecraftState transitState = state.shiftedBy(tauDMDelta.getValue());
    // Transit state position (re)computed with derivative structures
    final TimeStampedFieldPVCoordinates<DerivativeStructure> transitStateDS = pvaDS.shiftedBy(tauDMDelta);
    // Station at transit state date (derivatives of tauD taken into account)
    final TimeStampedFieldPVCoordinates<DerivativeStructure> stationAtTransitDate = stationDownlink.shiftedBy(tauD.negate());
    // Uplink delay
    final DerivativeStructure tauU = signalTimeOfFlight(stationAtTransitDate, transitStateDS.getPosition(), transitStateDS.getDate());
    // Prepare the evaluation
    final EstimatedMeasurement<Range> estimated = new EstimatedMeasurement<Range>(this, iteration, evaluation, new SpacecraftState[] { transitState }, null);
    // Range value
    final DerivativeStructure tau = tauD.add(tauU);
    final double cOver2 = 0.5 * Constants.SPEED_OF_LIGHT;
    final DerivativeStructure range = tau.multiply(cOver2);
    estimated.setEstimatedValue(range.getValue());
    // Range partial derivatives with respect to state
    final double[] derivatives = range.getAllDerivatives();
    estimated.setStateDerivatives(0, Arrays.copyOfRange(derivatives, 1, 7));
    // (beware element at index 0 is the value, not a derivative)
    for (final ParameterDriver driver : getParametersDrivers()) {
        final Integer index = indices.get(driver.getName());
        if (index != null) {
            estimated.setParameterDerivatives(driver, derivatives[index + 1]);
        }
    }
    // ----------
    // VALIDATION
    // -----------
    // Computation of the value without DS
    // ----------------------------------
    // Time difference between t (date of the measurement) and t' (date tagged in spacecraft state)
    // Station position at signal arrival
    final Transform topoToInertDownlink = groundStation.getOffsetToInertial(state.getFrame(), downlinkDate);
    final PVCoordinates QDownlink = topoToInertDownlink.transformPVCoordinates(PVCoordinates.ZERO);
    // Downlink time of flight from spacecraft to station
    final double td = signalTimeOfFlight(state.getPVCoordinates(), QDownlink.getPosition(), downlinkDate);
    final double dt = delta - td;
    // Transit state position
    final AbsoluteDate transitT = state.getDate().shiftedBy(dt);
    final SpacecraftState transit = state.shiftedBy(dt);
    final Vector3D transitP = transitState.getPVCoordinates().getPosition();
    // Station position at signal departure
    // First guess
    // AbsoluteDate uplinkDate = downlinkDate.shiftedBy(-getObservedValue()[0] / cOver2);
    // final Transform topoToInertUplink =
    // station.getOffsetFrame().getTransformTo(state.getFrame(), uplinkDate);
    // TimeStampedPVCoordinates QUplink = topoToInertUplink.
    // transformPVCoordinates(new TimeStampedPVCoordinates(uplinkDate, PVCoordinates.ZERO));
    // Station position at transit state date
    final Transform topoToInertAtTransitDate = groundStation.getOffsetToInertial(state.getFrame(), transitT);
    TimeStampedPVCoordinates QAtTransitDate = topoToInertAtTransitDate.transformPVCoordinates(new TimeStampedPVCoordinates(transitT, PVCoordinates.ZERO));
    // Uplink time of flight
    final double tu = signalTimeOfFlight(QAtTransitDate, transitP, transitT);
    // Total time of flight
    final double t = td + tu;
    // Real date and position of station at signal departure
    AbsoluteDate uplinkDate = downlinkDate.shiftedBy(-t);
    TimeStampedPVCoordinates QUplink = topoToInertDownlink.shiftedBy(-t).transformPVCoordinates(new TimeStampedPVCoordinates(uplinkDate, PVCoordinates.ZERO));
    // Range value
    double r = t * cOver2;
    double dR = r - range.getValue();
    // td derivatives / state
    // -----------------------
    // Qt = Master station position at tmeas = t = signal arrival at master station
    final Vector3D vel = state.getPVCoordinates().getVelocity();
    final Vector3D Qt_V = QDownlink.getVelocity();
    final Vector3D Ptr = transit.getPVCoordinates().getPosition();
    final Vector3D Ptr_Qt = QDownlink.getPosition().subtract(Ptr);
    final double dDown = Constants.SPEED_OF_LIGHT * Constants.SPEED_OF_LIGHT * td - Vector3D.dotProduct(Ptr_Qt, vel);
    // Derivatives of the downlink time of flight
    final double dtddPx = -Ptr_Qt.getX() / dDown;
    final double dtddPy = -Ptr_Qt.getY() / dDown;
    final double dtddPz = -Ptr_Qt.getZ() / dDown;
    final double dtddVx = dtddPx * dt;
    final double dtddVy = dtddPy * dt;
    final double dtddVz = dtddPz * dt;
    // From the DS
    final double dtddPxDS = tauD.getPartialDerivative(1, 0, 0, 0, 0, 0, 0, 0, 0);
    final double dtddPyDS = tauD.getPartialDerivative(0, 1, 0, 0, 0, 0, 0, 0, 0);
    final double dtddPzDS = tauD.getPartialDerivative(0, 0, 1, 0, 0, 0, 0, 0, 0);
    final double dtddVxDS = tauD.getPartialDerivative(0, 0, 0, 1, 0, 0, 0, 0, 0);
    final double dtddVyDS = tauD.getPartialDerivative(0, 0, 0, 0, 1, 0, 0, 0, 0);
    final double dtddVzDS = tauD.getPartialDerivative(0, 0, 0, 0, 0, 1, 0, 0, 0);
    // Difference
    final double d_dtddPx = dtddPxDS - dtddPx;
    final double d_dtddPy = dtddPyDS - dtddPy;
    final double d_dtddPz = dtddPzDS - dtddPz;
    final double d_dtddVx = dtddVxDS - dtddVx;
    final double d_dtddVy = dtddVyDS - dtddVy;
    final double d_dtddVz = dtddVzDS - dtddVz;
    // tu derivatives / state
    // -----------------------
    final Vector3D Qt2_Ptr = Ptr.subtract(QUplink.getPosition());
    final double dUp = Constants.SPEED_OF_LIGHT * Constants.SPEED_OF_LIGHT * tu - Vector3D.dotProduct(Qt2_Ptr, Qt_V);
    // test
    // // Speed of the station at tmeas-t
    // // Note: Which one to use in the calculation of dUp ???
    // final Vector3D Qt2_V    = QUplink.getVelocity();
    // final double   dUp      = Constants.SPEED_OF_LIGHT * Constants.SPEED_OF_LIGHT * tu -
    // Vector3D.dotProduct(Qt2_Ptr, Qt2_V);
    // test
    // tu derivatives
    final double dtudPx = 1. / dUp * Qt2_Ptr.dotProduct(Vector3D.PLUS_I.add((Qt_V.subtract(vel)).scalarMultiply(dtddPx)));
    final double dtudPy = 1. / dUp * Qt2_Ptr.dotProduct(Vector3D.PLUS_J.add((Qt_V.subtract(vel)).scalarMultiply(dtddPy)));
    final double dtudPz = 1. / dUp * Qt2_Ptr.dotProduct(Vector3D.PLUS_K.add((Qt_V.subtract(vel)).scalarMultiply(dtddPz)));
    final double dtudVx = dtudPx * dt;
    final double dtudVy = dtudPy * dt;
    final double dtudVz = dtudPz * dt;
    // From the DS
    final double dtudPxDS = tauU.getPartialDerivative(1, 0, 0, 0, 0, 0, 0, 0, 0);
    final double dtudPyDS = tauU.getPartialDerivative(0, 1, 0, 0, 0, 0, 0, 0, 0);
    final double dtudPzDS = tauU.getPartialDerivative(0, 0, 1, 0, 0, 0, 0, 0, 0);
    final double dtudVxDS = tauU.getPartialDerivative(0, 0, 0, 1, 0, 0, 0, 0, 0);
    final double dtudVyDS = tauU.getPartialDerivative(0, 0, 0, 0, 1, 0, 0, 0, 0);
    final double dtudVzDS = tauU.getPartialDerivative(0, 0, 0, 0, 0, 1, 0, 0, 0);
    // Difference
    final double d_dtudPx = dtudPxDS - dtudPx;
    final double d_dtudPy = dtudPyDS - dtudPy;
    final double d_dtudPz = dtudPzDS - dtudPz;
    final double d_dtudVx = dtudVxDS - dtudVx;
    final double d_dtudVy = dtudVyDS - dtudVy;
    final double d_dtudVz = dtudVzDS - dtudVz;
    // Range derivatives / state
    // -----------------------
    // R = Range
    double dRdPx = (dtddPx + dtudPx) * cOver2;
    double dRdPy = (dtddPy + dtudPy) * cOver2;
    double dRdPz = (dtddPz + dtudPz) * cOver2;
    double dRdVx = (dtddVx + dtudVx) * cOver2;
    double dRdVy = (dtddVy + dtudVy) * cOver2;
    double dRdVz = (dtddVz + dtudVz) * cOver2;
    // With DS
    double dRdPxDS = range.getPartialDerivative(1, 0, 0, 0, 0, 0, 0, 0, 0);
    double dRdPyDS = range.getPartialDerivative(0, 1, 0, 0, 0, 0, 0, 0, 0);
    double dRdPzDS = range.getPartialDerivative(0, 0, 1, 0, 0, 0, 0, 0, 0);
    double dRdVxDS = range.getPartialDerivative(0, 0, 0, 1, 0, 0, 0, 0, 0);
    double dRdVyDS = range.getPartialDerivative(0, 0, 0, 0, 1, 0, 0, 0, 0);
    double dRdVzDS = range.getPartialDerivative(0, 0, 0, 0, 0, 1, 0, 0, 0);
    // Diff
    final double d_dRdPx = dRdPxDS - dRdPx;
    final double d_dRdPy = dRdPyDS - dRdPy;
    final double d_dRdPz = dRdPzDS - dRdPz;
    final double d_dRdVx = dRdVxDS - dRdVx;
    final double d_dRdVy = dRdVyDS - dRdVy;
    final double d_dRdVz = dRdVzDS - dRdVz;
    // td derivatives / station
    // -----------------------
    final AngularCoordinates ac = topoToInertDownlink.getAngular().revert();
    final Rotation rotTopoToInert = ac.getRotation();
    final Vector3D omega = ac.getRotationRate();
    final Vector3D dtddQI = Ptr_Qt.scalarMultiply(1. / dDown);
    final double dtddQIx = dtddQI.getX();
    final double dtddQIy = dtddQI.getY();
    final double dtddQIz = dtddQI.getZ();
    final Vector3D dtddQ = rotTopoToInert.applyTo(dtddQI);
    // With DS
    double dtddQxDS = tauD.getPartialDerivative(0, 0, 0, 0, 0, 0, 1, 0, 0);
    double dtddQyDS = tauD.getPartialDerivative(0, 0, 0, 0, 0, 0, 0, 1, 0);
    double dtddQzDS = tauD.getPartialDerivative(0, 0, 0, 0, 0, 0, 0, 0, 1);
    // Diff
    final double d_dtddQx = dtddQxDS - dtddQ.getX();
    final double d_dtddQy = dtddQyDS - dtddQ.getY();
    final double d_dtddQz = dtddQzDS - dtddQ.getZ();
    // tu derivatives / station
    // -----------------------
    // Inertial frame
    final double dtudQIx = 1 / dUp * Qt2_Ptr.dotProduct(Vector3D.MINUS_I.add((Qt_V.subtract(vel)).scalarMultiply(dtddQIx)).subtract(Vector3D.PLUS_I.crossProduct(omega).scalarMultiply(t)));
    final double dtudQIy = 1 / dUp * Qt2_Ptr.dotProduct(Vector3D.MINUS_J.add((Qt_V.subtract(vel)).scalarMultiply(dtddQIy)).subtract(Vector3D.PLUS_J.crossProduct(omega).scalarMultiply(t)));
    final double dtudQIz = 1 / dUp * Qt2_Ptr.dotProduct(Vector3D.MINUS_K.add((Qt_V.subtract(vel)).scalarMultiply(dtddQIz)).subtract(Vector3D.PLUS_K.crossProduct(omega).scalarMultiply(t)));
    // // test
    // final double dtudQIx = 1/dUp*Qt2_Ptr
    // //                        .dotProduct(Vector3D.MINUS_I);
    // //                                    .dotProduct((Qt_V.subtract(vel)).scalarMultiply(dtddQIx));
    // .dotProduct(Vector3D.MINUS_I.crossProduct(omega).scalarMultiply(t));
    // final double dtudQIy = 1/dUp*Qt2_Ptr
    // //                        .dotProduct(Vector3D.MINUS_J);
    // //                                    .dotProduct((Qt_V.subtract(vel)).scalarMultiply(dtddQIy));
    // .dotProduct(Vector3D.MINUS_J.crossProduct(omega).scalarMultiply(t));
    // final double dtudQIz = 1/dUp*Qt2_Ptr
    // //                        .dotProduct(Vector3D.MINUS_K);
    // //                                    .dotProduct((Qt_V.subtract(vel)).scalarMultiply(dtddQIz));
    // .dotProduct(Vector3D.MINUS_K.crossProduct(omega).scalarMultiply(t));
    // 
    // double dtu_dQxDS = tauU.getPartialDerivative(0, 0, 0, 0, 0, 0, 1, 0, 0);
    // double dtu_dQyDS = tauU.getPartialDerivative(0, 0, 0, 0, 0, 0, 0, 1, 0);
    // double dtu_dQzDS = tauU.getPartialDerivative(0, 0, 0, 0, 0, 0, 0, 0, 1);
    // final Vector3D dtudQDS = new Vector3D(dtu_dQxDS, dtu_dQyDS, dtu_dQzDS);
    // final Vector3D dtudQIDS = rotTopoToInert.applyInverseTo(dtudQDS);
    // double dtudQIxDS = dtudQIDS.getX();
    // double dtudQIyDS = dtudQIDS.getY();
    // double dtudQIxzS = dtudQIDS.getZ();
    // // test
    // Topocentric frame
    final Vector3D dtudQI = new Vector3D(dtudQIx, dtudQIy, dtudQIz);
    final Vector3D dtudQ = rotTopoToInert.applyTo(dtudQI);
    // With DS
    double dtudQxDS = tauU.getPartialDerivative(0, 0, 0, 0, 0, 0, 1, 0, 0);
    double dtudQyDS = tauU.getPartialDerivative(0, 0, 0, 0, 0, 0, 0, 1, 0);
    double dtudQzDS = tauU.getPartialDerivative(0, 0, 0, 0, 0, 0, 0, 0, 1);
    // Diff
    final double d_dtudQx = dtudQxDS - dtudQ.getX();
    final double d_dtudQy = dtudQyDS - dtudQ.getY();
    final double d_dtudQz = dtudQzDS - dtudQ.getZ();
    // Range derivatives / station
    // -----------------------
    double dRdQx = (dtddQ.getX() + dtudQ.getX()) * cOver2;
    double dRdQy = (dtddQ.getY() + dtudQ.getY()) * cOver2;
    double dRdQz = (dtddQ.getZ() + dtudQ.getZ()) * cOver2;
    // With DS
    double dRdQxDS = range.getPartialDerivative(0, 0, 0, 0, 0, 0, 1, 0, 0);
    double dRdQyDS = range.getPartialDerivative(0, 0, 0, 0, 0, 0, 0, 1, 0);
    double dRdQzDS = range.getPartialDerivative(0, 0, 0, 0, 0, 0, 0, 0, 1);
    // Diff
    final double d_dRdQx = dRdQxDS - dRdQx;
    final double d_dRdQy = dRdQyDS - dRdQy;
    final double d_dRdQz = dRdQzDS - dRdQz;
    // Print results to avoid warning
    final boolean printResults = false;
    if (printResults) {
        System.out.println("dR = " + dR);
        System.out.println("d_dtddPx = " + d_dtddPx);
        System.out.println("d_dtddPy = " + d_dtddPy);
        System.out.println("d_dtddPz = " + d_dtddPz);
        System.out.println("d_dtddVx = " + d_dtddVx);
        System.out.println("d_dtddVy = " + d_dtddVy);
        System.out.println("d_dtddVz = " + d_dtddVz);
        System.out.println("d_dtudPx = " + d_dtudPx);
        System.out.println("d_dtudPy = " + d_dtudPy);
        System.out.println("d_dtudPz = " + d_dtudPz);
        System.out.println("d_dtudVx = " + d_dtudVx);
        System.out.println("d_dtudVy = " + d_dtudVy);
        System.out.println("d_dtudVz = " + d_dtudVz);
        System.out.println("d_dRdPx = " + d_dRdPx);
        System.out.println("d_dRdPy = " + d_dRdPy);
        System.out.println("d_dRdPz = " + d_dRdPz);
        System.out.println("d_dRdVx = " + d_dRdVx);
        System.out.println("d_dRdVy = " + d_dRdVy);
        System.out.println("d_dRdVz = " + d_dRdVz);
        System.out.println("d_dtddQx = " + d_dtddQx);
        System.out.println("d_dtddQy = " + d_dtddQy);
        System.out.println("d_dtddQz = " + d_dtddQz);
        System.out.println("d_dtudQx = " + d_dtudQx);
        System.out.println("d_dtudQy = " + d_dtudQy);
        System.out.println("d_dtudQz = " + d_dtudQz);
        System.out.println("d_dRdQx = " + d_dRdQx);
        System.out.println("d_dRdQy = " + d_dRdQy);
        System.out.println("d_dRdQz = " + d_dRdQz);
    }
    // Dummy return
    return estimated;
}
Also used : HashMap(java.util.HashMap) TimeStampedPVCoordinates(org.orekit.utils.TimeStampedPVCoordinates) PVCoordinates(org.orekit.utils.PVCoordinates) TimeStampedFieldPVCoordinates(org.orekit.utils.TimeStampedFieldPVCoordinates) TimeStampedPVCoordinates(org.orekit.utils.TimeStampedPVCoordinates) FieldAbsoluteDate(org.orekit.time.FieldAbsoluteDate) AbsoluteDate(org.orekit.time.AbsoluteDate) SpacecraftState(org.orekit.propagation.SpacecraftState) Vector3D(org.hipparchus.geometry.euclidean.threed.Vector3D) FieldVector3D(org.hipparchus.geometry.euclidean.threed.FieldVector3D) AngularCoordinates(org.orekit.utils.AngularCoordinates) DerivativeStructure(org.hipparchus.analysis.differentiation.DerivativeStructure) DSFactory(org.hipparchus.analysis.differentiation.DSFactory) ParameterDriver(org.orekit.utils.ParameterDriver) Rotation(org.hipparchus.geometry.euclidean.threed.Rotation) Transform(org.orekit.frames.Transform) FieldTransform(org.orekit.frames.FieldTransform) FieldAbsoluteDate(org.orekit.time.FieldAbsoluteDate)

Aggregations

Transform (org.orekit.frames.Transform)75 Vector3D (org.hipparchus.geometry.euclidean.threed.Vector3D)56 AbsoluteDate (org.orekit.time.AbsoluteDate)33 Frame (org.orekit.frames.Frame)28 FieldTransform (org.orekit.frames.FieldTransform)26 SpacecraftState (org.orekit.propagation.SpacecraftState)26 FieldVector3D (org.hipparchus.geometry.euclidean.threed.FieldVector3D)25 PVCoordinates (org.orekit.utils.PVCoordinates)23 TimeStampedPVCoordinates (org.orekit.utils.TimeStampedPVCoordinates)21 Test (org.junit.Test)20 Rotation (org.hipparchus.geometry.euclidean.threed.Rotation)18 FieldAbsoluteDate (org.orekit.time.FieldAbsoluteDate)17 GeodeticPoint (org.orekit.bodies.GeodeticPoint)13 TopocentricFrame (org.orekit.frames.TopocentricFrame)12 OneAxisEllipsoid (org.orekit.bodies.OneAxisEllipsoid)11 OrekitException (org.orekit.errors.OrekitException)11 DerivativeStructure (org.hipparchus.analysis.differentiation.DerivativeStructure)10 FieldPVCoordinates (org.orekit.utils.FieldPVCoordinates)10 FieldRotation (org.hipparchus.geometry.euclidean.threed.FieldRotation)8 CircularOrbit (org.orekit.orbits.CircularOrbit)8