Search in sources :

Example 66 with Orbit

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

the class SolarRadiationPressureTest method testGlobalStateJacobianIsotropicSingle.

@Test
public void testGlobalStateJacobianIsotropicSingle() throws OrekitException {
    // initialization
    AbsoluteDate date = new AbsoluteDate(new DateComponents(2003, 03, 01), new TimeComponents(13, 59, 27.816), TimeScalesFactory.getUTC());
    double i = FastMath.toRadians(98.7);
    double omega = FastMath.toRadians(93.0);
    double OMEGA = FastMath.toRadians(15.0 * 22.5);
    Orbit orbit = new KeplerianOrbit(7201009.7124401, 1e-3, i, omega, OMEGA, 0, PositionAngle.MEAN, FramesFactory.getEME2000(), date, Constants.EIGEN5C_EARTH_MU);
    OrbitType integrationType = OrbitType.CARTESIAN;
    double[][] tolerances = NumericalPropagator.tolerances(0.01, orbit, integrationType);
    NumericalPropagator propagator = new NumericalPropagator(new DormandPrince853Integrator(1.0e-3, 120, tolerances[0], tolerances[1]));
    propagator.setOrbitType(integrationType);
    SolarRadiationPressure forceModel = new SolarRadiationPressure(CelestialBodyFactory.getSun(), Constants.WGS84_EARTH_EQUATORIAL_RADIUS, new IsotropicRadiationSingleCoefficient(2.5, 0.7));
    propagator.addForceModel(forceModel);
    SpacecraftState state0 = new SpacecraftState(orbit);
    checkStateJacobian(propagator, state0, date.shiftedBy(3.5 * 3600.0), 1e3, tolerances[0], 2.0e-5);
}
Also used : EquinoctialOrbit(org.orekit.orbits.EquinoctialOrbit) CartesianOrbit(org.orekit.orbits.CartesianOrbit) FieldKeplerianOrbit(org.orekit.orbits.FieldKeplerianOrbit) KeplerianOrbit(org.orekit.orbits.KeplerianOrbit) Orbit(org.orekit.orbits.Orbit) DateComponents(org.orekit.time.DateComponents) TimeComponents(org.orekit.time.TimeComponents) FieldAbsoluteDate(org.orekit.time.FieldAbsoluteDate) AbsoluteDate(org.orekit.time.AbsoluteDate) SpacecraftState(org.orekit.propagation.SpacecraftState) FieldSpacecraftState(org.orekit.propagation.FieldSpacecraftState) NumericalPropagator(org.orekit.propagation.numerical.NumericalPropagator) FieldNumericalPropagator(org.orekit.propagation.numerical.FieldNumericalPropagator) FieldKeplerianOrbit(org.orekit.orbits.FieldKeplerianOrbit) KeplerianOrbit(org.orekit.orbits.KeplerianOrbit) OrbitType(org.orekit.orbits.OrbitType) DormandPrince853Integrator(org.hipparchus.ode.nonstiff.DormandPrince853Integrator) AbstractLegacyForceModelTest(org.orekit.forces.AbstractLegacyForceModelTest) Test(org.junit.Test)

Example 67 with Orbit

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

the class SolarRadiationPressureTest method testGlobalStateJacobianBox.

@Test
public void testGlobalStateJacobianBox() throws OrekitException {
    // initialization
    AbsoluteDate date = new AbsoluteDate(new DateComponents(2003, 03, 01), new TimeComponents(13, 59, 27.816), TimeScalesFactory.getUTC());
    double i = FastMath.toRadians(98.7);
    double omega = FastMath.toRadians(93.0);
    double OMEGA = FastMath.toRadians(15.0 * 22.5);
    Orbit orbit = new KeplerianOrbit(7201009.7124401, 1e-3, i, omega, OMEGA, 0, PositionAngle.MEAN, FramesFactory.getEME2000(), date, Constants.EIGEN5C_EARTH_MU);
    OrbitType integrationType = OrbitType.CARTESIAN;
    double[][] tolerances = NumericalPropagator.tolerances(0.01, orbit, integrationType);
    NumericalPropagator propagator = new NumericalPropagator(new DormandPrince853Integrator(1.0e-3, 120, tolerances[0], tolerances[1]));
    propagator.setOrbitType(integrationType);
    SolarRadiationPressure forceModel = new SolarRadiationPressure(CelestialBodyFactory.getSun(), Constants.WGS84_EARTH_EQUATORIAL_RADIUS, new BoxAndSolarArraySpacecraft(1.5, 2.0, 1.8, CelestialBodyFactory.getSun(), 20.0, Vector3D.PLUS_J, 1.2, 0.7, 0.2));
    propagator.addForceModel(forceModel);
    SpacecraftState state0 = new SpacecraftState(orbit);
    checkStateJacobian(propagator, state0, date.shiftedBy(3.5 * 3600.0), 1e3, tolerances[0], 5.0e-4);
}
Also used : EquinoctialOrbit(org.orekit.orbits.EquinoctialOrbit) CartesianOrbit(org.orekit.orbits.CartesianOrbit) FieldKeplerianOrbit(org.orekit.orbits.FieldKeplerianOrbit) KeplerianOrbit(org.orekit.orbits.KeplerianOrbit) Orbit(org.orekit.orbits.Orbit) DateComponents(org.orekit.time.DateComponents) TimeComponents(org.orekit.time.TimeComponents) FieldAbsoluteDate(org.orekit.time.FieldAbsoluteDate) AbsoluteDate(org.orekit.time.AbsoluteDate) BoxAndSolarArraySpacecraft(org.orekit.forces.BoxAndSolarArraySpacecraft) SpacecraftState(org.orekit.propagation.SpacecraftState) FieldSpacecraftState(org.orekit.propagation.FieldSpacecraftState) NumericalPropagator(org.orekit.propagation.numerical.NumericalPropagator) FieldNumericalPropagator(org.orekit.propagation.numerical.FieldNumericalPropagator) FieldKeplerianOrbit(org.orekit.orbits.FieldKeplerianOrbit) KeplerianOrbit(org.orekit.orbits.KeplerianOrbit) OrbitType(org.orekit.orbits.OrbitType) DormandPrince853Integrator(org.hipparchus.ode.nonstiff.DormandPrince853Integrator) AbstractLegacyForceModelTest(org.orekit.forces.AbstractLegacyForceModelTest) Test(org.junit.Test)

Example 68 with Orbit

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

the class SolarRadiationPressureTest method doTestLocalJacobianIsotropicClassicalVsFiniteDifferences.

private void doTestLocalJacobianIsotropicClassicalVsFiniteDifferences(double deltaT, double dP, double checkTolerance, boolean print) throws OrekitException {
    // initialization
    AbsoluteDate date = new AbsoluteDate(new DateComponents(2003, 03, 01), new TimeComponents(13, 59, 27.816), TimeScalesFactory.getUTC());
    double i = FastMath.toRadians(98.7);
    double omega = FastMath.toRadians(93.0);
    double OMEGA = FastMath.toRadians(15.0 * 22.5);
    Orbit orbit = new KeplerianOrbit(7201009.7124401, 1e-3, i, omega, OMEGA, 0, PositionAngle.MEAN, FramesFactory.getEME2000(), date, Constants.EIGEN5C_EARTH_MU);
    final SolarRadiationPressure forceModel = new SolarRadiationPressure(CelestialBodyFactory.getSun(), Constants.WGS84_EARTH_EQUATORIAL_RADIUS, new IsotropicRadiationClassicalConvention(2.5, 0.7, 0.2));
    checkStateJacobianVsFiniteDifferences(new SpacecraftState(orbit.shiftedBy(deltaT)), forceModel, Propagator.DEFAULT_LAW, dP, checkTolerance, print);
}
Also used : SpacecraftState(org.orekit.propagation.SpacecraftState) FieldSpacecraftState(org.orekit.propagation.FieldSpacecraftState) EquinoctialOrbit(org.orekit.orbits.EquinoctialOrbit) CartesianOrbit(org.orekit.orbits.CartesianOrbit) FieldKeplerianOrbit(org.orekit.orbits.FieldKeplerianOrbit) KeplerianOrbit(org.orekit.orbits.KeplerianOrbit) Orbit(org.orekit.orbits.Orbit) DateComponents(org.orekit.time.DateComponents) FieldKeplerianOrbit(org.orekit.orbits.FieldKeplerianOrbit) KeplerianOrbit(org.orekit.orbits.KeplerianOrbit) TimeComponents(org.orekit.time.TimeComponents) FieldAbsoluteDate(org.orekit.time.FieldAbsoluteDate) AbsoluteDate(org.orekit.time.AbsoluteDate)

Example 69 with Orbit

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

the class Model method value.

/**
 * {@inheritDoc}
 */
@Override
public Pair<RealVector, RealMatrix> value(final RealVector point) throws OrekitExceptionWrapper {
    try {
        // Set up the propagators parallelizer
        final NumericalPropagator[] propagators = createPropagators(point);
        final Orbit[] orbits = new Orbit[propagators.length];
        for (int i = 0; i < propagators.length; ++i) {
            mappers[i] = configureDerivatives(propagators[i]);
            orbits[i] = propagators[i].getInitialState().getOrbit();
        }
        final PropagatorsParallelizer parallelizer = new PropagatorsParallelizer(Arrays.asList(propagators), configureMeasurements(point));
        // Reset value and Jacobian
        evaluations.clear();
        value.set(0.0);
        for (int i = 0; i < jacobian.getRowDimension(); ++i) {
            for (int j = 0; j < jacobian.getColumnDimension(); ++j) {
                jacobian.setEntry(i, j, 0.0);
            }
        }
        // Run the propagation, gathering residuals on the fly
        if (forwardPropagation) {
            // Propagate forward from firstDate
            parallelizer.propagate(firstDate.shiftedBy(-1.0), lastDate.shiftedBy(+1.0));
        } else {
            // Propagate backward from lastDate
            parallelizer.propagate(lastDate.shiftedBy(+1.0), firstDate.shiftedBy(-1.0));
        }
        observer.modelCalled(orbits, evaluations);
        return new Pair<RealVector, RealMatrix>(value, jacobian);
    } catch (OrekitException oe) {
        throw new OrekitExceptionWrapper(oe);
    }
}
Also used : OrekitExceptionWrapper(org.orekit.errors.OrekitExceptionWrapper) Orbit(org.orekit.orbits.Orbit) NumericalPropagator(org.orekit.propagation.numerical.NumericalPropagator) PropagatorsParallelizer(org.orekit.propagation.PropagatorsParallelizer) OrekitException(org.orekit.errors.OrekitException) Pair(org.hipparchus.util.Pair)

Example 70 with Orbit

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

the class Model method getMeasurementMatrix.

/**
 * Get the normalized measurement matrix H.
 * H contains the partial derivatives of the measurement with respect to the state.
 * H is an NxM matrix where N is the size of the measurement vector and M the size of the state vector.
 * @param predictedSpacecraftState the spacecraft state associated with the measurement
 * @return the normalized measurement matrix H
 * @throws OrekitException if Jacobians cannot be computed
 */
private RealMatrix getMeasurementMatrix(final SpacecraftState predictedSpacecraftState) throws OrekitException {
    // Number of parameters
    final int nbOrb = estimatedOrbitalParameters.getNbParams();
    final int nbPropag = estimatedPropagationParameters.getNbParams();
    final int nbMeas = estimatedMeasurementsParameters.getNbParams();
    // Initialize measurement matrix H: N x M
    // N: Number of measurements in current measurement
    // M: State vector size
    final RealMatrix measurementMatrix = MatrixUtils.createRealMatrix(predictedMeasurement.getEstimatedValue().length, nbOrb + nbPropag + nbMeas);
    // Predicted orbit
    final Orbit predictedOrbit = predictedSpacecraftState.getOrbit();
    // Observed measurement characteristics
    final ObservedMeasurement<?> observedMeasurement = predictedMeasurement.getObservedMeasurement();
    final double[] sigma = observedMeasurement.getTheoreticalStandardDeviation();
    // Measurement matrix's columns related to orbital parameters
    // ----------------------------------------------------------
    // Partial derivatives of the current Cartesian coordinates with respect to current orbital state
    final double[][] aCY = new double[6][6];
    // dC/dY
    predictedOrbit.getJacobianWrtParameters(builder.getPositionAngle(), aCY);
    final RealMatrix dCdY = new Array2DRowRealMatrix(aCY, false);
    // Jacobian of the measurement with respect to current Cartesian coordinates
    final RealMatrix dMdC = new Array2DRowRealMatrix(predictedMeasurement.getStateDerivatives(0), false);
    // Jacobian of the measurement with respect to current orbital state
    final RealMatrix dMdY = dMdC.multiply(dCdY);
    // Fill the normalized measurement matrix's columns related to estimated orbital parameters
    if (nbOrb > 0) {
        for (int i = 0; i < dMdY.getRowDimension(); ++i) {
            int jOrb = 0;
            for (int j = 0; j < dMdY.getColumnDimension(); ++j) {
                final DelegatingDriver delegating = builder.getOrbitalParametersDrivers().getDrivers().get(j);
                if (delegating.isSelected()) {
                    measurementMatrix.setEntry(i, jOrb++, dMdY.getEntry(i, j) / sigma[i] * delegating.getScale());
                }
            }
        }
    }
    // Jacobian of the measurement with respect to propagation parameters
    if (nbPropag > 0) {
        final double[][] aYPp = new double[6][nbPropag];
        jacobiansMapper.getParametersJacobian(predictedSpacecraftState, aYPp);
        final RealMatrix dYdPp = new Array2DRowRealMatrix(aYPp, false);
        final RealMatrix dMdPp = dMdY.multiply(dYdPp);
        for (int i = 0; i < dMdPp.getRowDimension(); ++i) {
            for (int j = 0; j < nbPropag; ++j) {
                final DelegatingDriver delegating = estimatedPropagationParameters.getDrivers().get(j);
                measurementMatrix.setEntry(i, nbOrb + j, dMdPp.getEntry(i, j) / sigma[i] * delegating.getScale());
            }
        }
    }
    // Jacobian of the measurement with respect to measurement parameters
    if (nbMeas > 0) {
        // Gather the measurement parameters linked to current measurement
        for (final ParameterDriver driver : observedMeasurement.getParametersDrivers()) {
            if (driver.isSelected()) {
                // Derivatives of current measurement w/r to selected measurement parameter
                final double[] aMPm = predictedMeasurement.getParameterDerivatives(driver);
                // Check that the measurement parameter is managed by the filter
                if (measurementParameterColumns.get(driver.getName()) != null) {
                    // Column of the driver in the measurement matrix
                    final int driverColumn = measurementParameterColumns.get(driver.getName());
                    // Fill the corresponding indexes of the measurement matrix
                    for (int i = 0; i < aMPm.length; ++i) {
                        measurementMatrix.setEntry(i, driverColumn, aMPm[i] / sigma[i] * driver.getScale());
                    }
                }
            }
        }
    }
    // Return the normalized measurement matrix
    return measurementMatrix;
}
Also used : RealMatrix(org.hipparchus.linear.RealMatrix) Array2DRowRealMatrix(org.hipparchus.linear.Array2DRowRealMatrix) Orbit(org.orekit.orbits.Orbit) Array2DRowRealMatrix(org.hipparchus.linear.Array2DRowRealMatrix) DelegatingDriver(org.orekit.utils.ParameterDriversList.DelegatingDriver) ParameterDriver(org.orekit.utils.ParameterDriver)

Aggregations

Orbit (org.orekit.orbits.Orbit)211 KeplerianOrbit (org.orekit.orbits.KeplerianOrbit)161 Test (org.junit.Test)153 AbsoluteDate (org.orekit.time.AbsoluteDate)153 SpacecraftState (org.orekit.propagation.SpacecraftState)129 Vector3D (org.hipparchus.geometry.euclidean.threed.Vector3D)99 EquinoctialOrbit (org.orekit.orbits.EquinoctialOrbit)94 CartesianOrbit (org.orekit.orbits.CartesianOrbit)88 FieldAbsoluteDate (org.orekit.time.FieldAbsoluteDate)74 CircularOrbit (org.orekit.orbits.CircularOrbit)68 PVCoordinates (org.orekit.utils.PVCoordinates)66 Frame (org.orekit.frames.Frame)51 NumericalPropagator (org.orekit.propagation.numerical.NumericalPropagator)51 DateComponents (org.orekit.time.DateComponents)48 FieldSpacecraftState (org.orekit.propagation.FieldSpacecraftState)46 Propagator (org.orekit.propagation.Propagator)46 TimeComponents (org.orekit.time.TimeComponents)44 OneAxisEllipsoid (org.orekit.bodies.OneAxisEllipsoid)43 AbstractLegacyForceModelTest (org.orekit.forces.AbstractLegacyForceModelTest)41 FieldKeplerianOrbit (org.orekit.orbits.FieldKeplerianOrbit)39