Search in sources :

Example 46 with Context

use of org.orekit.estimation.Context in project Orekit by CS-SI.

the class ModelTest method testPerfectValue.

@Test
public void testPerfectValue() throws OrekitException {
    final Context context = EstimationTestUtils.eccentricContext("regular-data:potential:tides");
    final NumericalPropagatorBuilder propagatorBuilder = context.createBuilder(OrbitType.KEPLERIAN, PositionAngle.TRUE, true, 1.0e-6, 60.0, 0.001);
    final NumericalPropagatorBuilder[] builders = { propagatorBuilder };
    // create perfect PV measurements
    final Propagator propagator = EstimationTestUtils.createPropagator(context.initialOrbit, propagatorBuilder);
    final List<ObservedMeasurement<?>> measurements = EstimationTestUtils.createMeasurements(propagator, new PVMeasurementCreator(), 0.0, 1.0, 300.0);
    final ParameterDriversList estimatedMeasurementsParameters = new ParameterDriversList();
    for (ObservedMeasurement<?> measurement : measurements) {
        for (final ParameterDriver driver : measurement.getParametersDrivers()) {
            if (driver.isSelected()) {
                estimatedMeasurementsParameters.add(driver);
            }
        }
    }
    // create model
    final ModelObserver modelObserver = new ModelObserver() {

        /**
         * {@inheritDoc}
         */
        @Override
        public void modelCalled(final Orbit[] newOrbits, final Map<ObservedMeasurement<?>, EstimatedMeasurement<?>> newEvaluations) {
            Assert.assertEquals(1, newOrbits.length);
            Assert.assertEquals(0, context.initialOrbit.getDate().durationFrom(newOrbits[0].getDate()), 1.0e-15);
            Assert.assertEquals(0, Vector3D.distance(context.initialOrbit.getPVCoordinates().getPosition(), newOrbits[0].getPVCoordinates().getPosition()), 1.0e-15);
            Assert.assertEquals(measurements.size(), newEvaluations.size());
        }
    };
    final Model model = new Model(builders, measurements, estimatedMeasurementsParameters, modelObserver);
    model.setIterationsCounter(new Incrementor(100));
    model.setEvaluationsCounter(new Incrementor(100));
    // Test forward propagation flag to true
    assertEquals(true, model.isForwardPropagation());
    // evaluate model on perfect start point
    final double[] normalizedProp = propagatorBuilder.getSelectedNormalizedParameters();
    final double[] normalized = new double[normalizedProp.length + estimatedMeasurementsParameters.getNbParams()];
    System.arraycopy(normalizedProp, 0, normalized, 0, normalizedProp.length);
    int i = normalizedProp.length;
    for (final ParameterDriver driver : estimatedMeasurementsParameters.getDrivers()) {
        normalized[i++] = driver.getNormalizedValue();
    }
    Pair<RealVector, RealMatrix> value = model.value(new ArrayRealVector(normalized));
    int index = 0;
    for (ObservedMeasurement<?> measurement : measurements) {
        for (int k = 0; k < measurement.getDimension(); ++k) {
            // the value is already a weighted residual
            Assert.assertEquals(0.0, value.getFirst().getEntry(index++), 1.6e-7);
        }
    }
    Assert.assertEquals(index, value.getFirst().getDimension());
}
Also used : Context(org.orekit.estimation.Context) ArrayRealVector(org.hipparchus.linear.ArrayRealVector) Incrementor(org.hipparchus.util.Incrementor) ParameterDriver(org.orekit.utils.ParameterDriver) RealMatrix(org.hipparchus.linear.RealMatrix) ParameterDriversList(org.orekit.utils.ParameterDriversList) NumericalPropagatorBuilder(org.orekit.propagation.conversion.NumericalPropagatorBuilder) ArrayRealVector(org.hipparchus.linear.ArrayRealVector) RealVector(org.hipparchus.linear.RealVector) Propagator(org.orekit.propagation.Propagator) Map(java.util.Map) ObservedMeasurement(org.orekit.estimation.measurements.ObservedMeasurement) PVMeasurementCreator(org.orekit.estimation.measurements.PVMeasurementCreator) Test(org.junit.Test)

Example 47 with Context

use of org.orekit.estimation.Context in project Orekit by CS-SI.

the class AngularAzElTest method testStateDerivatives.

@Test
public void testStateDerivatives() throws OrekitException {
    Context context = EstimationTestUtils.geoStationnaryContext("regular-data:potential:tides");
    final NumericalPropagatorBuilder propagatorBuilder = context.createBuilder(OrbitType.EQUINOCTIAL, PositionAngle.TRUE, false, 1.0e-6, 60.0, 0.001);
    // create perfect azimuth-elevation measurements
    final Propagator propagator = EstimationTestUtils.createPropagator(context.initialOrbit, propagatorBuilder);
    final List<ObservedMeasurement<?>> measurements = EstimationTestUtils.createMeasurements(propagator, new AngularAzElMeasurementCreator(context), 0.25, 3.0, 600.0);
    propagator.setSlaveMode();
    // Compute measurements.
    double[] AzerrorsP = new double[3 * measurements.size()];
    double[] AzerrorsV = new double[3 * measurements.size()];
    double[] ElerrorsP = new double[3 * measurements.size()];
    double[] ElerrorsV = new double[3 * measurements.size()];
    int AzindexP = 0;
    int AzindexV = 0;
    int ElindexP = 0;
    int ElindexV = 0;
    for (final ObservedMeasurement<?> measurement : measurements) {
        // parameter corresponding to station position offset
        final GroundStation stationParameter = ((AngularAzEl) measurement).getStation();
        // We intentionally propagate to a date which is close to the
        // real spacecraft state but is *not* the accurate date, by
        // compensating only part of the downlink delay. This is done
        // in order to validate the partial derivatives with respect
        // to velocity. If we had chosen the proper state date, the
        // angular would have depended only on the current position but
        // not on the current velocity.
        final AbsoluteDate datemeas = measurement.getDate();
        SpacecraftState state = propagator.propagate(datemeas);
        final Vector3D stationP = stationParameter.getOffsetToInertial(state.getFrame(), datemeas).transformPosition(Vector3D.ZERO);
        final double meanDelay = AbstractMeasurement.signalTimeOfFlight(state.getPVCoordinates(), stationP, datemeas);
        final AbsoluteDate date = measurement.getDate().shiftedBy(-0.75 * meanDelay);
        state = propagator.propagate(date);
        final EstimatedMeasurement<?> estimated = measurement.estimate(0, 0, new SpacecraftState[] { state });
        Assert.assertEquals(2, estimated.getParticipants().length);
        final double[][] jacobian = estimated.getStateDerivatives(0);
        // compute a reference value using finite differences
        final double[][] finiteDifferencesJacobian = Differentiation.differentiate(new StateFunction() {

            public double[] value(final SpacecraftState state) throws OrekitException {
                return measurement.estimate(0, 0, new SpacecraftState[] { state }).getEstimatedValue();
            }
        }, measurement.getDimension(), propagator.getAttitudeProvider(), OrbitType.CARTESIAN, PositionAngle.TRUE, 250.0, 4).value(state);
        Assert.assertEquals(finiteDifferencesJacobian.length, jacobian.length);
        Assert.assertEquals(finiteDifferencesJacobian[0].length, jacobian[0].length);
        final double smallest = FastMath.ulp((double) 1.0);
        for (int i = 0; i < jacobian.length; ++i) {
            for (int j = 0; j < jacobian[i].length; ++j) {
                double relativeError = FastMath.abs((finiteDifferencesJacobian[i][j] - jacobian[i][j]) / finiteDifferencesJacobian[i][j]);
                if ((FastMath.sqrt(finiteDifferencesJacobian[i][j]) < smallest) && (FastMath.sqrt(jacobian[i][j]) < smallest)) {
                    relativeError = 0.0;
                }
                if (j < 3) {
                    if (i == 0) {
                        AzerrorsP[AzindexP++] = relativeError;
                    } else {
                        ElerrorsP[ElindexP++] = relativeError;
                    }
                } else {
                    if (i == 0) {
                        AzerrorsV[AzindexV++] = relativeError;
                    } else {
                        ElerrorsV[ElindexV++] = relativeError;
                    }
                }
            }
        }
    }
    // median errors on Azimuth
    Assert.assertEquals(0.0, new Median().evaluate(AzerrorsP), 1.1e-10);
    Assert.assertEquals(0.0, new Median().evaluate(AzerrorsV), 5.7e-5);
    // median errors on Elevation
    Assert.assertEquals(0.0, new Median().evaluate(ElerrorsP), 3.5e-11);
    Assert.assertEquals(0.0, new Median().evaluate(ElerrorsV), 1.4e-5);
}
Also used : Context(org.orekit.estimation.Context) Median(org.hipparchus.stat.descriptive.rank.Median) AbsoluteDate(org.orekit.time.AbsoluteDate) SpacecraftState(org.orekit.propagation.SpacecraftState) Vector3D(org.hipparchus.geometry.euclidean.threed.Vector3D) NumericalPropagatorBuilder(org.orekit.propagation.conversion.NumericalPropagatorBuilder) Propagator(org.orekit.propagation.Propagator) StateFunction(org.orekit.utils.StateFunction) Test(org.junit.Test)

Example 48 with Context

use of org.orekit.estimation.Context in project Orekit by CS-SI.

the class AngularAzElTest method testParameterDerivatives.

@Test
public void testParameterDerivatives() throws OrekitException {
    Context context = EstimationTestUtils.geoStationnaryContext("regular-data:potential:tides");
    final NumericalPropagatorBuilder propagatorBuilder = context.createBuilder(OrbitType.EQUINOCTIAL, PositionAngle.TRUE, false, 1.0e-6, 60.0, 0.001);
    // create perfect azimuth-elevation measurements
    for (final GroundStation station : context.stations) {
        station.getEastOffsetDriver().setSelected(true);
        station.getNorthOffsetDriver().setSelected(true);
        station.getZenithOffsetDriver().setSelected(true);
    }
    final Propagator propagator = EstimationTestUtils.createPropagator(context.initialOrbit, propagatorBuilder);
    final List<ObservedMeasurement<?>> measurements = EstimationTestUtils.createMeasurements(propagator, new AngularAzElMeasurementCreator(context), 0.25, 3.0, 600.0);
    propagator.setSlaveMode();
    for (final ObservedMeasurement<?> measurement : measurements) {
        // parameter corresponding to station position offset
        final GroundStation stationParameter = ((AngularAzEl) measurement).getStation();
        // We intentionally propagate to a date which is close to the
        // real spacecraft state but is *not* the accurate date, by
        // compensating only part of the downlink delay. This is done
        // in order to validate the partial derivatives with respect
        // to velocity. If we had chosen the proper state date, the
        // angular would have depended only on the current position but
        // not on the current velocity.
        final AbsoluteDate datemeas = measurement.getDate();
        final SpacecraftState stateini = propagator.propagate(datemeas);
        final Vector3D stationP = stationParameter.getOffsetToInertial(stateini.getFrame(), datemeas).transformPosition(Vector3D.ZERO);
        final double meanDelay = AbstractMeasurement.signalTimeOfFlight(stateini.getPVCoordinates(), stationP, datemeas);
        final AbsoluteDate date = measurement.getDate().shiftedBy(-0.75 * meanDelay);
        final SpacecraftState state = propagator.propagate(date);
        final ParameterDriver[] drivers = new ParameterDriver[] { stationParameter.getEastOffsetDriver(), stationParameter.getNorthOffsetDriver(), stationParameter.getZenithOffsetDriver() };
        for (int i = 0; i < 3; ++i) {
            final double[] gradient = measurement.estimate(0, 0, new SpacecraftState[] { state }).getParameterDerivatives(drivers[i]);
            Assert.assertEquals(2, measurement.getDimension());
            Assert.assertEquals(2, gradient.length);
            for (final int k : new int[] { 0, 1 }) {
                final ParameterFunction dMkdP = Differentiation.differentiate(new ParameterFunction() {

                    /**
                     * {@inheritDoc}
                     */
                    @Override
                    public double value(final ParameterDriver parameterDriver) throws OrekitException {
                        return measurement.estimate(0, 0, new SpacecraftState[] { state }).getEstimatedValue()[k];
                    }
                }, drivers[i], 3, 50.0);
                final double ref = dMkdP.value(drivers[i]);
                if (ref > 1.e-12) {
                    Assert.assertEquals(ref, gradient[k], 3e-10 * FastMath.abs(ref));
                }
            }
        }
    }
}
Also used : Context(org.orekit.estimation.Context) ParameterDriver(org.orekit.utils.ParameterDriver) AbsoluteDate(org.orekit.time.AbsoluteDate) SpacecraftState(org.orekit.propagation.SpacecraftState) Vector3D(org.hipparchus.geometry.euclidean.threed.Vector3D) NumericalPropagatorBuilder(org.orekit.propagation.conversion.NumericalPropagatorBuilder) ParameterFunction(org.orekit.utils.ParameterFunction) Propagator(org.orekit.propagation.Propagator) OrekitException(org.orekit.errors.OrekitException) Test(org.junit.Test)

Example 49 with Context

use of org.orekit.estimation.Context in project Orekit by CS-SI.

the class GroundStationTest method testEstimateEOP.

@Test
public void testEstimateEOP() throws OrekitException {
    Context linearEOPContext = EstimationTestUtils.eccentricContext("linear-EOP:regular-data/de431-ephemerides:potential:tides");
    final AbsoluteDate refDate = new AbsoluteDate(2000, 2, 24, linearEOPContext.utc);
    final double dut10 = 0.3079738;
    final double lod = 0.0011000;
    final double xp0 = 68450.0e-6;
    final double xpDot = -50.0e-6;
    final double yp0 = 60.0e-6;
    final double ypDot = 2.0e-6;
    for (double dt = -2 * Constants.JULIAN_DAY; dt < 2 * Constants.JULIAN_DAY; dt += 300.0) {
        AbsoluteDate date = refDate.shiftedBy(dt);
        Assert.assertEquals(dut10 - dt * lod / Constants.JULIAN_DAY, linearEOPContext.ut1.getEOPHistory().getUT1MinusUTC(date), 1.0e-15);
        Assert.assertEquals(lod, linearEOPContext.ut1.getEOPHistory().getLOD(date), 1.0e-15);
        Assert.assertEquals((xp0 + xpDot * dt / Constants.JULIAN_DAY) * Constants.ARC_SECONDS_TO_RADIANS, linearEOPContext.ut1.getEOPHistory().getPoleCorrection(date).getXp(), 1.0e-15);
        Assert.assertEquals((yp0 + ypDot * dt / Constants.JULIAN_DAY) * Constants.ARC_SECONDS_TO_RADIANS, linearEOPContext.ut1.getEOPHistory().getPoleCorrection(date).getYp(), 1.0e-15);
    }
    final NumericalPropagatorBuilder linearPropagatorBuilder = linearEOPContext.createBuilder(OrbitType.KEPLERIAN, PositionAngle.TRUE, true, 1.0e-6, 60.0, 0.001);
    // create perfect range measurements
    final Propagator propagator = EstimationTestUtils.createPropagator(linearEOPContext.initialOrbit, linearPropagatorBuilder);
    final List<ObservedMeasurement<?>> linearMeasurements = EstimationTestUtils.createMeasurements(propagator, new RangeMeasurementCreator(linearEOPContext), 1.0, 5.0, 60.0);
    Utils.clearFactories();
    Context zeroEOPContext = EstimationTestUtils.eccentricContext("zero-EOP:regular-data/de431-ephemerides:potential:potential:tides");
    for (double dt = -2 * Constants.JULIAN_DAY; dt < 2 * Constants.JULIAN_DAY; dt += 300.0) {
        AbsoluteDate date = refDate.shiftedBy(dt);
        Assert.assertEquals(0.0, zeroEOPContext.ut1.getEOPHistory().getUT1MinusUTC(date), 1.0e-15);
        Assert.assertEquals(0.0, zeroEOPContext.ut1.getEOPHistory().getLOD(date), 1.0e-15);
        Assert.assertEquals(0.0, zeroEOPContext.ut1.getEOPHistory().getPoleCorrection(date).getXp(), 1.0e-15);
        Assert.assertEquals(0.0, zeroEOPContext.ut1.getEOPHistory().getPoleCorrection(date).getYp(), 1.0e-15);
    }
    // create orbit estimator
    final NumericalPropagatorBuilder zeroPropagatorBuilder = linearEOPContext.createBuilder(OrbitType.KEPLERIAN, PositionAngle.TRUE, true, 1.0e-6, 60.0, 0.001);
    final BatchLSEstimator estimator = new BatchLSEstimator(new LevenbergMarquardtOptimizer(), zeroPropagatorBuilder);
    for (final ObservedMeasurement<?> linearMeasurement : linearMeasurements) {
        Range linearRange = (Range) linearMeasurement;
        for (final GroundStation station : zeroEOPContext.stations) {
            if (station.getBaseFrame().getName().equals(linearRange.getStation().getBaseFrame().getName())) {
                Range zeroRange = new Range(station, linearRange.getDate(), linearRange.getObservedValue()[0], linearRange.getTheoreticalStandardDeviation()[0], linearRange.getBaseWeight()[0]);
                estimator.addMeasurement(zeroRange);
            }
        }
    }
    estimator.setParametersConvergenceThreshold(1.0e-3);
    estimator.setMaxIterations(100);
    estimator.setMaxEvaluations(200);
    // we want to estimate pole and prime meridian
    GroundStation station = zeroEOPContext.stations.get(0);
    station.getPrimeMeridianOffsetDriver().setReferenceDate(refDate);
    station.getPrimeMeridianOffsetDriver().setSelected(true);
    station.getPrimeMeridianDriftDriver().setSelected(true);
    station.getPolarOffsetXDriver().setReferenceDate(refDate);
    station.getPolarOffsetXDriver().setSelected(true);
    station.getPolarDriftXDriver().setSelected(true);
    station.getPolarOffsetYDriver().setReferenceDate(refDate);
    station.getPolarOffsetYDriver().setSelected(true);
    station.getPolarDriftYDriver().setSelected(true);
    // just for the fun and to speed up test, we will use orbit determination, *without* estimating orbit
    for (final ParameterDriver driver : zeroPropagatorBuilder.getOrbitalParametersDrivers().getDrivers()) {
        driver.setSelected(false);
    }
    estimator.estimate();
    final double computedDut1 = station.getPrimeMeridianOffsetDriver().getValue() / EstimatedEarthFrameProvider.EARTH_ANGULAR_VELOCITY;
    final double computedLOD = station.getPrimeMeridianDriftDriver().getValue() * (-Constants.JULIAN_DAY / EstimatedEarthFrameProvider.EARTH_ANGULAR_VELOCITY);
    final double computedXp = station.getPolarOffsetXDriver().getValue() / Constants.ARC_SECONDS_TO_RADIANS;
    final double computedXpDot = station.getPolarDriftXDriver().getValue() / Constants.ARC_SECONDS_TO_RADIANS * Constants.JULIAN_DAY;
    final double computedYp = station.getPolarOffsetYDriver().getValue() / Constants.ARC_SECONDS_TO_RADIANS;
    final double computedYpDot = station.getPolarDriftYDriver().getValue() / Constants.ARC_SECONDS_TO_RADIANS * Constants.JULIAN_DAY;
    Assert.assertEquals(dut10, computedDut1, 4.3e-10);
    Assert.assertEquals(lod, computedLOD, 4.9e-10);
    Assert.assertEquals(xp0, computedXp, 5.6e-9);
    Assert.assertEquals(xpDot, computedXpDot, 7.2e-9);
    Assert.assertEquals(yp0, computedYp, 1.1e-9);
    Assert.assertEquals(ypDot, computedYpDot, 2.8e-11);
// thresholds to use if orbit is estimated
// (i.e. when commenting out the loop above that sets orbital parameters drivers to "not selected")
// Assert.assertEquals(dut10, computedDut1,  6.6e-3);
// Assert.assertEquals(lod,   computedLOD,   1.1e-9);
// Assert.assertEquals(xp0,   computedXp,    3.3e-8);
// Assert.assertEquals(xpDot, computedXpDot, 2.2e-8);
// Assert.assertEquals(yp0,   computedYp,    3.3e-8);
// Assert.assertEquals(ypDot, computedYpDot, 3.8e-8);
}
Also used : Context(org.orekit.estimation.Context) ParameterDriver(org.orekit.utils.ParameterDriver) FieldAbsoluteDate(org.orekit.time.FieldAbsoluteDate) AbsoluteDate(org.orekit.time.AbsoluteDate) BatchLSEstimator(org.orekit.estimation.leastsquares.BatchLSEstimator) LevenbergMarquardtOptimizer(org.hipparchus.optim.nonlinear.vector.leastsquares.LevenbergMarquardtOptimizer) NumericalPropagatorBuilder(org.orekit.propagation.conversion.NumericalPropagatorBuilder) Propagator(org.orekit.propagation.Propagator) Test(org.junit.Test)

Example 50 with Context

use of org.orekit.estimation.Context in project Orekit by CS-SI.

the class TurnAroundRangeAnalyticTest method genericTestStateDerivatives.

/**
 * Generic test function for derivatives with respect to state
 * @param isModifier Use of atmospheric modifiers
 * @param isFiniteDifferences Finite differences reference calculation if true, TurnAroundRange class otherwise
 * @param printResults Print the results ?
 * @throws OrekitException
 */
void genericTestStateDerivatives(final boolean isModifier, final boolean isFiniteDifferences, final boolean printResults, final double refErrorsPMedian, final double refErrorsPMean, final double refErrorsPMax, final double refErrorsVMedian, final double refErrorsVMean, final double refErrorsVMax) throws OrekitException {
    Context context = EstimationTestUtils.eccentricContext("regular-data:potential:tides");
    // Context context = EstimationTestUtils.geoStationnaryContext();
    final NumericalPropagatorBuilder propagatorBuilder = context.createBuilder(OrbitType.KEPLERIAN, PositionAngle.TRUE, true, 1.0e-6, 60.0, 0.001);
    // create perfect range2 measurements
    final Propagator propagator = EstimationTestUtils.createPropagator(context.initialOrbit, propagatorBuilder);
    final List<ObservedMeasurement<?>> measurements = EstimationTestUtils.createMeasurements(propagator, new TurnAroundRangeMeasurementCreator(context), 1.0, 3.0, 300.0);
    propagator.setSlaveMode();
    double[] errorsP = new double[3 * measurements.size()];
    double[] errorsV = new double[3 * measurements.size()];
    int indexP = 0;
    int indexV = 0;
    // Print the results ? Header
    if (printResults) {
        System.out.format(Locale.US, "%-15s  %-15s  %-23s  %-23s  " + "%10s  %10s  %10s  " + "%10s  %10s  %10s  " + "%10s  %10s  %10s  " + "%10s  %10s  %10s%n", "Master Station", "Slave Station", "Measurement Date", "State Date", "ΔdPx", "ΔdPy", "ΔdPz", "ΔdVx", "ΔdVy", "ΔdVz", "rel ΔdPx", "rel ΔdPy", "rel ΔdPz", "rel ΔdVx", "rel ΔdVy", "rel ΔdVz");
    }
    // Loop on the measurements
    for (final ObservedMeasurement<?> measurement : measurements) {
        // Add modifiers if test implies it
        final TurnAroundRangeTroposphericDelayModifier modifier = new TurnAroundRangeTroposphericDelayModifier(SaastamoinenModel.getStandardModel());
        if (isModifier) {
            ((TurnAroundRange) measurement).addModifier(modifier);
        }
        // We intentionally propagate to a date which is close to the
        // real spacecraft state but is *not* the accurate date, by
        // compensating only part of the downlink delay. This is done
        // in order to validate the partial derivatives with respect
        // to velocity. If we had chosen the proper state date, the
        // range would have depended only on the current position but
        // not on the current velocity.
        final double meanDelay = measurement.getObservedValue()[0] / Constants.SPEED_OF_LIGHT;
        final AbsoluteDate date = measurement.getDate().shiftedBy(-0.75 * meanDelay);
        final SpacecraftState state = propagator.propagate(date);
        final EstimatedMeasurement<TurnAroundRange> TAR = new TurnAroundRangeAnalytic((TurnAroundRange) measurement).theoreticalEvaluationAnalytic(0, 0, propagator.getInitialState(), state);
        if (isModifier) {
            modifier.modify(TAR);
        }
        final double[][] jacobian = TAR.getStateDerivatives(0);
        // Jacobian reference value
        final double[][] jacobianRef;
        if (isFiniteDifferences) {
            // Compute a reference value using finite differences
            jacobianRef = Differentiation.differentiate(new StateFunction() {

                public double[] value(final SpacecraftState state) throws OrekitException {
                    return measurement.estimate(0, 0, new SpacecraftState[] { state }).getEstimatedValue();
                }
            }, measurement.getDimension(), propagator.getAttitudeProvider(), OrbitType.CARTESIAN, PositionAngle.TRUE, 2.0, 3).value(state);
        } else {
            // Compute a reference value using TurnAroundRange class function
            jacobianRef = ((TurnAroundRange) measurement).theoreticalEvaluation(0, 0, new SpacecraftState[] { state }).getStateDerivatives(0);
        }
        // //Test: Test point by point with the debugger
        // if (!isFiniteDifferences && !isModifier) {
        // final EstimatedMeasurement<TurnAroundRange> test =
        // new TurnAroundRangeAnalytic((TurnAroundRange)measurement).theoreticalEvaluationValidation(0, 0, state);
        // }
        // //Test
        Assert.assertEquals(jacobianRef.length, jacobian.length);
        Assert.assertEquals(jacobianRef[0].length, jacobian[0].length);
        double[][] dJacobian = new double[jacobian.length][jacobian[0].length];
        double[][] dJacobianRelative = new double[jacobian.length][jacobian[0].length];
        for (int i = 0; i < jacobian.length; ++i) {
            for (int j = 0; j < jacobian[i].length; ++j) {
                dJacobian[i][j] = jacobian[i][j] - jacobianRef[i][j];
                dJacobianRelative[i][j] = FastMath.abs(dJacobian[i][j] / jacobianRef[i][j]);
                if (j < 3) {
                    errorsP[indexP++] = dJacobianRelative[i][j];
                } else {
                    errorsV[indexV++] = dJacobianRelative[i][j];
                }
            }
        }
        // Print results on the console ? Print the Jacobian
        if (printResults) {
            String masterStationName = ((TurnAroundRange) measurement).getMasterStation().getBaseFrame().getName();
            String slaveStationName = ((TurnAroundRange) measurement).getSlaveStation().getBaseFrame().getName();
            System.out.format(Locale.US, "%-15s  %-15s  %-23s  %-23s  " + "%10.3e  %10.3e  %10.3e  " + "%10.3e  %10.3e  %10.3e  " + "%10.3e  %10.3e  %10.3e  " + "%10.3e  %10.3e  %10.3e%n", masterStationName, slaveStationName, measurement.getDate(), date, dJacobian[0][0], dJacobian[0][1], dJacobian[0][2], dJacobian[0][3], dJacobian[0][4], dJacobian[0][5], dJacobianRelative[0][0], dJacobianRelative[0][1], dJacobianRelative[0][2], dJacobianRelative[0][3], dJacobianRelative[0][4], dJacobianRelative[0][5]);
        }
    }
    // End loop on the measurements
    // Compute some statistics
    final double errorsPMedian = new Median().evaluate(errorsP);
    final double errorsPMean = new Mean().evaluate(errorsP);
    final double errorsPMax = new Max().evaluate(errorsP);
    final double errorsVMedian = new Median().evaluate(errorsV);
    final double errorsVMean = new Mean().evaluate(errorsV);
    final double errorsVMax = new Max().evaluate(errorsV);
    // Print the results on console ? Final results
    if (printResults) {
        System.out.println();
        System.out.format(Locale.US, "Relative errors dR/dP -> Median: %6.3e / Mean: %6.3e / Max: %6.3e%n", errorsPMedian, errorsPMean, errorsPMax);
        System.out.format(Locale.US, "Relative errors dR/dV -> Median: %6.3e / Mean: %6.3e / Max: %6.3e%n", errorsVMedian, errorsVMean, errorsVMax);
    }
    // Assert the results / max values depend on the test
    Assert.assertEquals(0.0, errorsPMedian, refErrorsPMedian);
    Assert.assertEquals(0.0, errorsPMean, refErrorsPMean);
    Assert.assertEquals(0.0, errorsPMax, refErrorsPMax);
    Assert.assertEquals(0.0, errorsVMedian, refErrorsVMedian);
    Assert.assertEquals(0.0, errorsVMean, refErrorsVMean);
    Assert.assertEquals(0.0, errorsVMax, refErrorsVMax);
}
Also used : Context(org.orekit.estimation.Context) Mean(org.hipparchus.stat.descriptive.moment.Mean) Max(org.hipparchus.stat.descriptive.rank.Max) Median(org.hipparchus.stat.descriptive.rank.Median) AbsoluteDate(org.orekit.time.AbsoluteDate) SpacecraftState(org.orekit.propagation.SpacecraftState) NumericalPropagatorBuilder(org.orekit.propagation.conversion.NumericalPropagatorBuilder) Propagator(org.orekit.propagation.Propagator) TurnAroundRangeTroposphericDelayModifier(org.orekit.estimation.measurements.modifiers.TurnAroundRangeTroposphericDelayModifier) OrekitException(org.orekit.errors.OrekitException) StateFunction(org.orekit.utils.StateFunction)

Aggregations

Context (org.orekit.estimation.Context)74 Propagator (org.orekit.propagation.Propagator)67 NumericalPropagatorBuilder (org.orekit.propagation.conversion.NumericalPropagatorBuilder)67 Test (org.junit.Test)60 AbsoluteDate (org.orekit.time.AbsoluteDate)49 ObservedMeasurement (org.orekit.estimation.measurements.ObservedMeasurement)40 SpacecraftState (org.orekit.propagation.SpacecraftState)35 Vector3D (org.hipparchus.geometry.euclidean.threed.Vector3D)28 ParameterDriver (org.orekit.utils.ParameterDriver)21 OrekitException (org.orekit.errors.OrekitException)18 Median (org.hipparchus.stat.descriptive.rank.Median)17 RangeMeasurementCreator (org.orekit.estimation.measurements.RangeMeasurementCreator)17 Orbit (org.orekit.orbits.Orbit)17 ParameterDriversList (org.orekit.utils.ParameterDriversList)16 ArrayList (java.util.ArrayList)14 Max (org.hipparchus.stat.descriptive.rank.Max)14 BoundedPropagator (org.orekit.propagation.BoundedPropagator)13 RealMatrix (org.hipparchus.linear.RealMatrix)12 LevenbergMarquardtOptimizer (org.hipparchus.optim.nonlinear.vector.leastsquares.LevenbergMarquardtOptimizer)12 StateFunction (org.orekit.utils.StateFunction)11