Search in sources :

Example 16 with ParameterDriver

use of org.orekit.utils.ParameterDriver in project Orekit by CS-SI.

the class DSConverter method getState.

/**
 * Get the state with the number of parameters consistent with force model.
 * @param forceModel force model
 * @return state with the number of parameters consistent with force model
 */
public FieldSpacecraftState<DerivativeStructure> getState(final ForceModel forceModel) {
    // count the required number of parameters
    int nbParams = 0;
    for (final ParameterDriver driver : forceModel.getParametersDrivers()) {
        if (driver.isSelected()) {
            ++nbParams;
        }
    }
    // fill in intermediate slots
    while (dsStates.size() < nbParams + 1) {
        dsStates.add(null);
    }
    if (dsStates.get(nbParams) == null) {
        // it is the first time we need this number of parameters
        // we need to create the state
        final DSFactory factory = new DSFactory(freeStateParameters + nbParams, 1);
        final FieldSpacecraftState<DerivativeStructure> s0 = dsStates.get(0);
        // orbit
        final FieldPVCoordinates<DerivativeStructure> pv0 = s0.getPVCoordinates();
        final FieldOrbit<DerivativeStructure> dsOrbit = new FieldCartesianOrbit<>(new TimeStampedFieldPVCoordinates<>(s0.getDate().toAbsoluteDate(), extend(pv0.getPosition(), factory), extend(pv0.getVelocity(), factory), extend(pv0.getAcceleration(), factory)), s0.getFrame(), s0.getMu());
        // attitude
        final FieldAngularCoordinates<DerivativeStructure> ac0 = s0.getAttitude().getOrientation();
        final FieldAttitude<DerivativeStructure> dsAttitude = new FieldAttitude<>(s0.getAttitude().getReferenceFrame(), new TimeStampedFieldAngularCoordinates<>(dsOrbit.getDate(), extend(ac0.getRotation(), factory), extend(ac0.getRotationRate(), factory), extend(ac0.getRotationAcceleration(), factory)));
        // mass
        final DerivativeStructure dsM = extend(s0.getMass(), factory);
        dsStates.set(nbParams, new FieldSpacecraftState<>(dsOrbit, dsAttitude, dsM));
    }
    return dsStates.get(nbParams);
}
Also used : DerivativeStructure(org.hipparchus.analysis.differentiation.DerivativeStructure) DSFactory(org.hipparchus.analysis.differentiation.DSFactory) ParameterDriver(org.orekit.utils.ParameterDriver) FieldCartesianOrbit(org.orekit.orbits.FieldCartesianOrbit) FieldAttitude(org.orekit.attitudes.FieldAttitude)

Example 17 with ParameterDriver

use of org.orekit.utils.ParameterDriver in project Orekit by CS-SI.

the class DSConverter method getParameters.

/**
 * Get the force model parameters.
 * @param state state as returned by {@link #getState(ForceModel)}
 * @param forceModel force model associated with the parameters
 * @return force model parameters
 * @since 9.0
 */
public DerivativeStructure[] getParameters(final FieldSpacecraftState<DerivativeStructure> state, final ForceModel forceModel) {
    final DSFactory factory = state.getMass().getFactory();
    final ParameterDriver[] drivers = forceModel.getParametersDrivers();
    final DerivativeStructure[] parameters = new DerivativeStructure[drivers.length];
    int index = freeStateParameters;
    for (int i = 0; i < drivers.length; ++i) {
        parameters[i] = drivers[i].isSelected() ? factory.variable(index++, drivers[i].getValue()) : factory.constant(drivers[i].getValue());
    }
    return parameters;
}
Also used : DerivativeStructure(org.hipparchus.analysis.differentiation.DerivativeStructure) DSFactory(org.hipparchus.analysis.differentiation.DSFactory) ParameterDriver(org.orekit.utils.ParameterDriver)

Example 18 with ParameterDriver

use of org.orekit.utils.ParameterDriver in project Orekit by CS-SI.

the class PartialDerivativesEquations method computeDerivatives.

/**
 * {@inheritDoc}
 */
public double[] computeDerivatives(final SpacecraftState s, final double[] pDot) throws OrekitException {
    // initialize acceleration Jacobians to zero
    final int paramDim = selected.getNbParams();
    final int dim = 3;
    final double[][] dAccdParam = new double[dim][paramDim];
    final double[][] dAccdPos = new double[dim][dim];
    final double[][] dAccdVel = new double[dim][dim];
    final DSConverter fullConverter = new DSConverter(s, 6, propagator.getAttitudeProvider());
    final DSConverter posOnlyConverter = new DSConverter(s, 3, propagator.getAttitudeProvider());
    // compute acceleration Jacobians, finishing with the largest force: Newtonian attraction
    for (final ForceModel forceModel : propagator.getAllForceModels()) {
        final DSConverter converter = forceModel.dependsOnPositionOnly() ? posOnlyConverter : fullConverter;
        final FieldSpacecraftState<DerivativeStructure> dsState = converter.getState(forceModel);
        final DerivativeStructure[] parameters = converter.getParameters(dsState, forceModel);
        final FieldVector3D<DerivativeStructure> acceleration = forceModel.acceleration(dsState, parameters);
        final double[] derivativesX = acceleration.getX().getAllDerivatives();
        final double[] derivativesY = acceleration.getY().getAllDerivatives();
        final double[] derivativesZ = acceleration.getZ().getAllDerivatives();
        // update Jacobians with respect to state
        addToRow(derivativesX, 0, converter.getFreeStateParameters(), dAccdPos, dAccdVel);
        addToRow(derivativesY, 1, converter.getFreeStateParameters(), dAccdPos, dAccdVel);
        addToRow(derivativesZ, 2, converter.getFreeStateParameters(), dAccdPos, dAccdVel);
        int index = converter.getFreeStateParameters();
        for (ParameterDriver driver : forceModel.getParametersDrivers()) {
            if (driver.isSelected()) {
                final int parameterIndex = map.get(driver);
                ++index;
                dAccdParam[0][parameterIndex] += derivativesX[index];
                dAccdParam[1][parameterIndex] += derivativesY[index];
                dAccdParam[2][parameterIndex] += derivativesZ[index];
            }
        }
    }
    // the variational equations of the complete state Jacobian matrix have the following form:
    // [        |        ]   [                 |                  ]   [     |     ]
    // [  Adot  |  Bdot  ]   [  dVel/dPos = 0  |  dVel/dVel = Id  ]   [  A  |  B  ]
    // [        |        ]   [                 |                  ]   [     |     ]
    // ---------+---------   ------------------+------------------- * ------+------
    // [        |        ]   [                 |                  ]   [     |     ]
    // [  Cdot  |  Ddot  ] = [    dAcc/dPos    |     dAcc/dVel    ]   [  C  |  D  ]
    // [        |        ]   [                 |                  ]   [     |     ]
    // The A, B, C and D sub-matrices and their derivatives (Adot ...) are 3x3 matrices
    // The expanded multiplication above can be rewritten to take into account
    // the fixed values found in the sub-matrices in the left factor. This leads to:
    // [ Adot ] = [ C ]
    // [ Bdot ] = [ D ]
    // [ Cdot ] = [ dAcc/dPos ] * [ A ] + [ dAcc/dVel ] * [ C ]
    // [ Ddot ] = [ dAcc/dPos ] * [ B ] + [ dAcc/dVel ] * [ D ]
    // The following loops compute these expressions taking care of the mapping of the
    // (A, B, C, D) matrices into the single dimension array p and of the mapping of the
    // (Adot, Bdot, Cdot, Ddot) matrices into the single dimension array pDot.
    // copy C and E into Adot and Bdot
    final int stateDim = 6;
    final double[] p = s.getAdditionalState(getName());
    System.arraycopy(p, dim * stateDim, pDot, 0, dim * stateDim);
    // compute Cdot and Ddot
    for (int i = 0; i < dim; ++i) {
        final double[] dAdPi = dAccdPos[i];
        final double[] dAdVi = dAccdVel[i];
        for (int j = 0; j < stateDim; ++j) {
            pDot[(dim + i) * stateDim + j] = dAdPi[0] * p[j] + dAdPi[1] * p[j + stateDim] + dAdPi[2] * p[j + 2 * stateDim] + dAdVi[0] * p[j + 3 * stateDim] + dAdVi[1] * p[j + 4 * stateDim] + dAdVi[2] * p[j + 5 * stateDim];
        }
    }
    for (int k = 0; k < paramDim; ++k) {
        // the variational equations of the parameters Jacobian matrix are computed
        // one column at a time, they have the following form:
        // [      ]   [                 |                  ]   [   ]   [                  ]
        // [ Edot ]   [  dVel/dPos = 0  |  dVel/dVel = Id  ]   [ E ]   [  dVel/dParam = 0 ]
        // [      ]   [                 |                  ]   [   ]   [                  ]
        // --------   ------------------+------------------- * ----- + --------------------
        // [      ]   [                 |                  ]   [   ]   [                  ]
        // [ Fdot ] = [    dAcc/dPos    |     dAcc/dVel    ]   [ F ]   [    dAcc/dParam   ]
        // [      ]   [                 |                  ]   [   ]   [                  ]
        // The E and F sub-columns and their derivatives (Edot, Fdot) are 3 elements columns.
        // The expanded multiplication and addition above can be rewritten to take into
        // account the fixed values found in the sub-matrices in the left factor. This leads to:
        // [ Edot ] = [ F ]
        // [ Fdot ] = [ dAcc/dPos ] * [ E ] + [ dAcc/dVel ] * [ F ] + [ dAcc/dParam ]
        // The following loops compute these expressions taking care of the mapping of the
        // (E, F) columns into the single dimension array p and of the mapping of the
        // (Edot, Fdot) columns into the single dimension array pDot.
        // copy F into Edot
        final int columnTop = stateDim * stateDim + k;
        pDot[columnTop] = p[columnTop + 3 * paramDim];
        pDot[columnTop + paramDim] = p[columnTop + 4 * paramDim];
        pDot[columnTop + 2 * paramDim] = p[columnTop + 5 * paramDim];
        // compute Fdot
        for (int i = 0; i < dim; ++i) {
            final double[] dAdPi = dAccdPos[i];
            final double[] dAdVi = dAccdVel[i];
            pDot[columnTop + (dim + i) * paramDim] = dAccdParam[i][k] + dAdPi[0] * p[columnTop] + dAdPi[1] * p[columnTop + paramDim] + dAdPi[2] * p[columnTop + 2 * paramDim] + dAdVi[0] * p[columnTop + 3 * paramDim] + dAdVi[1] * p[columnTop + 4 * paramDim] + dAdVi[2] * p[columnTop + 5 * paramDim];
        }
    }
    // these equations have no effect on the main state itself
    return null;
}
Also used : ForceModel(org.orekit.forces.ForceModel) DerivativeStructure(org.hipparchus.analysis.differentiation.DerivativeStructure) ParameterDriver(org.orekit.utils.ParameterDriver)

Example 19 with ParameterDriver

use of org.orekit.utils.ParameterDriver in project Orekit by CS-SI.

the class SolarRadiationPressureTest method testParameterDerivativeIsotropicSingle.

@Test
public void testParameterDerivativeIsotropicSingle() throws OrekitException {
    final Vector3D pos = new Vector3D(6.46885878304673824e+06, -1.88050918456274318e+06, -1.32931592294715829e+04);
    final Vector3D vel = new Vector3D(2.14718074509906819e+03, 7.38239351251748485e+03, -1.14097953925384523e+01);
    final SpacecraftState state = new SpacecraftState(new CartesianOrbit(new PVCoordinates(pos, vel), FramesFactory.getGCRF(), new AbsoluteDate(2003, 3, 5, 0, 24, 0.0, TimeScalesFactory.getTAI()), Constants.EIGEN5C_EARTH_MU));
    RadiationSensitive rs = new IsotropicRadiationSingleCoefficient(2.5, 0.7);
    SolarRadiationPressure forceModel = new SolarRadiationPressure(CelestialBodyFactory.getSun(), Constants.WGS84_EARTH_EQUATORIAL_RADIUS, rs);
    checkParameterDerivative(state, forceModel, RadiationSensitive.REFLECTION_COEFFICIENT, 1.0, 2.0e-15);
    try {
        rs.radiationPressureAcceleration(state.getDate(), state.getFrame(), state.getPVCoordinates().getPosition(), state.getAttitude().getRotation(), state.getMass(), Vector3D.ZERO, new double[2], RadiationSensitive.ABSORPTION_COEFFICIENT);
        Assert.fail("an exception should have been thrown");
    } catch (OrekitException oe) {
        Assert.assertEquals(OrekitMessages.UNSUPPORTED_PARAMETER_NAME, oe.getSpecifier());
    }
    for (ParameterDriver driver : rs.getRadiationParametersDrivers()) {
        Assert.assertEquals(RadiationSensitive.REFLECTION_COEFFICIENT, driver.getName());
    }
}
Also used : SpacecraftState(org.orekit.propagation.SpacecraftState) FieldSpacecraftState(org.orekit.propagation.FieldSpacecraftState) CartesianOrbit(org.orekit.orbits.CartesianOrbit) FieldVector3D(org.hipparchus.geometry.euclidean.threed.FieldVector3D) Vector3D(org.hipparchus.geometry.euclidean.threed.Vector3D) PVCoordinates(org.orekit.utils.PVCoordinates) FieldPVCoordinates(org.orekit.utils.FieldPVCoordinates) OrekitException(org.orekit.errors.OrekitException) ParameterDriver(org.orekit.utils.ParameterDriver) FieldAbsoluteDate(org.orekit.time.FieldAbsoluteDate) AbsoluteDate(org.orekit.time.AbsoluteDate) AbstractLegacyForceModelTest(org.orekit.forces.AbstractLegacyForceModelTest) Test(org.junit.Test)

Example 20 with ParameterDriver

use of org.orekit.utils.ParameterDriver in project Orekit by CS-SI.

the class Range method theoreticalEvaluation.

/**
 * {@inheritDoc}
 */
@Override
protected EstimatedMeasurement<Range> theoreticalEvaluation(final int iteration, final int evaluation, final SpacecraftState[] states) throws OrekitException {
    final SpacecraftState state = states[getPropagatorsIndices().get(0)];
    // Range derivatives are computed with respect to spacecraft state in inertial frame
    // and station parameters
    // ----------------------
    // 
    // Parameters:
    // - 0..2 - Position of the spacecraft in inertial frame
    // - 3..5 - Velocity of the spacecraft in inertial frame
    // - 6..n - station parameters (station offsets, pole, prime meridian...)
    int nbParams = 6;
    final Map<String, Integer> indices = new HashMap<>();
    for (ParameterDriver driver : getParametersDrivers()) {
        if (driver.isSelected()) {
            indices.put(driver.getName(), nbParams++);
        }
    }
    final DSFactory factory = new DSFactory(nbParams, 1);
    final Field<DerivativeStructure> field = factory.getDerivativeField();
    final FieldVector3D<DerivativeStructure> zero = FieldVector3D.getZero(field);
    // Coordinates of the spacecraft expressed as a derivative structure
    final TimeStampedFieldPVCoordinates<DerivativeStructure> pvaDS = getCoordinates(state, 0, factory);
    // 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 = station.getOffsetToInertial(state.getFrame(), downlinkDateDS, factory, 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 delta == tauD 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 deltaMTauD = tauD.negate().add(delta);
    final SpacecraftState transitState = state.shiftedBy(deltaMTauD.getValue());
    // Transit state (re)computed with derivative structures
    final TimeStampedFieldPVCoordinates<DerivativeStructure> transitStateDS = pvaDS.shiftedBy(deltaMTauD);
    // 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());
    final TimeStampedFieldPVCoordinates<DerivativeStructure> stationUplink = stationDownlink.shiftedBy(-tauD.getValue() - tauU.getValue());
    // Prepare the evaluation
    final EstimatedMeasurement<Range> estimated = new EstimatedMeasurement<Range>(this, iteration, evaluation, new SpacecraftState[] { transitState }, new TimeStampedPVCoordinates[] { stationUplink.toTimeStampedPVCoordinates(), transitStateDS.toTimeStampedPVCoordinates(), stationDownlink.toTimeStampedPVCoordinates() });
    // Range value
    final double cOver2 = 0.5 * Constants.SPEED_OF_LIGHT;
    final DerivativeStructure tau = tauD.add(tauU);
    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]);
        }
    }
    return estimated;
}
Also used : HashMap(java.util.HashMap) DerivativeStructure(org.hipparchus.analysis.differentiation.DerivativeStructure) DSFactory(org.hipparchus.analysis.differentiation.DSFactory) ParameterDriver(org.orekit.utils.ParameterDriver) FieldAbsoluteDate(org.orekit.time.FieldAbsoluteDate) AbsoluteDate(org.orekit.time.AbsoluteDate) SpacecraftState(org.orekit.propagation.SpacecraftState) FieldAbsoluteDate(org.orekit.time.FieldAbsoluteDate)

Aggregations

ParameterDriver (org.orekit.utils.ParameterDriver)80 AbsoluteDate (org.orekit.time.AbsoluteDate)33 SpacecraftState (org.orekit.propagation.SpacecraftState)32 NumericalPropagatorBuilder (org.orekit.propagation.conversion.NumericalPropagatorBuilder)27 Test (org.junit.Test)23 Propagator (org.orekit.propagation.Propagator)23 Context (org.orekit.estimation.Context)21 ParameterDriversList (org.orekit.utils.ParameterDriversList)20 OrekitException (org.orekit.errors.OrekitException)19 Vector3D (org.hipparchus.geometry.euclidean.threed.Vector3D)16 FieldAbsoluteDate (org.orekit.time.FieldAbsoluteDate)16 ObservedMeasurement (org.orekit.estimation.measurements.ObservedMeasurement)15 Orbit (org.orekit.orbits.Orbit)15 ArrayList (java.util.ArrayList)14 DerivativeStructure (org.hipparchus.analysis.differentiation.DerivativeStructure)14 ParameterFunction (org.orekit.utils.ParameterFunction)14 DSFactory (org.hipparchus.analysis.differentiation.DSFactory)13 HashMap (java.util.HashMap)11 KeplerianOrbit (org.orekit.orbits.KeplerianOrbit)11 RealMatrix (org.hipparchus.linear.RealMatrix)10