Search in sources :

Example 11 with LevenbergMarquardtOptimizer

use of org.hipparchus.optim.nonlinear.vector.leastsquares.LevenbergMarquardtOptimizer in project Orekit by CS-SI.

the class OrbitDeterminationTest method createEstimator.

/**
 * Set up estimator.
 * @param parser input file parser
 * @param propagatorBuilder propagator builder
 * @return estimator
 * @throws NoSuchElementException if input parameters are missing
 * @throws OrekitException if some propagator parameters cannot be retrieved
 */
private BatchLSEstimator createEstimator(final KeyValueFileParser<ParameterKey> parser, final NumericalPropagatorBuilder propagatorBuilder) throws NoSuchElementException, OrekitException {
    final boolean optimizerIsLevenbergMarquardt;
    if (!parser.containsKey(ParameterKey.ESTIMATOR_OPTIMIZATION_ENGINE)) {
        optimizerIsLevenbergMarquardt = true;
    } else {
        final String engine = parser.getString(ParameterKey.ESTIMATOR_OPTIMIZATION_ENGINE);
        optimizerIsLevenbergMarquardt = engine.toLowerCase().contains("levenberg");
    }
    final LeastSquaresOptimizer optimizer;
    if (optimizerIsLevenbergMarquardt) {
        // we want to use a Levenberg-Marquardt optimization engine
        final double initialStepBoundFactor;
        if (!parser.containsKey(ParameterKey.ESTIMATOR_LEVENBERG_MARQUARDT_INITIAL_STEP_BOUND_FACTOR)) {
            initialStepBoundFactor = 100.0;
        } else {
            initialStepBoundFactor = parser.getDouble(ParameterKey.ESTIMATOR_LEVENBERG_MARQUARDT_INITIAL_STEP_BOUND_FACTOR);
        }
        optimizer = new LevenbergMarquardtOptimizer().withInitialStepBoundFactor(initialStepBoundFactor);
    } else {
        // we want to use a Gauss-Newton optimization engine
        optimizer = new GaussNewtonOptimizer(Decomposition.QR);
    }
    final double convergenceThreshold;
    if (!parser.containsKey(ParameterKey.ESTIMATOR_NORMALIZED_PARAMETERS_CONVERGENCE_THRESHOLD)) {
        convergenceThreshold = 1.0e-3;
    } else {
        convergenceThreshold = parser.getDouble(ParameterKey.ESTIMATOR_NORMALIZED_PARAMETERS_CONVERGENCE_THRESHOLD);
    }
    final int maxIterations;
    if (!parser.containsKey(ParameterKey.ESTIMATOR_MAX_ITERATIONS)) {
        maxIterations = 10;
    } else {
        maxIterations = parser.getInt(ParameterKey.ESTIMATOR_MAX_ITERATIONS);
    }
    final int maxEvaluations;
    if (!parser.containsKey(ParameterKey.ESTIMATOR_MAX_EVALUATIONS)) {
        maxEvaluations = 20;
    } else {
        maxEvaluations = parser.getInt(ParameterKey.ESTIMATOR_MAX_EVALUATIONS);
    }
    final BatchLSEstimator estimator = new BatchLSEstimator(optimizer, propagatorBuilder);
    estimator.setParametersConvergenceThreshold(convergenceThreshold);
    estimator.setMaxIterations(maxIterations);
    estimator.setMaxEvaluations(maxEvaluations);
    return estimator;
}
Also used : LevenbergMarquardtOptimizer(org.hipparchus.optim.nonlinear.vector.leastsquares.LevenbergMarquardtOptimizer) GaussNewtonOptimizer(org.hipparchus.optim.nonlinear.vector.leastsquares.GaussNewtonOptimizer) LeastSquaresOptimizer(org.hipparchus.optim.nonlinear.vector.leastsquares.LeastSquaresOptimizer) GeodeticPoint(org.orekit.bodies.GeodeticPoint)

Example 12 with LevenbergMarquardtOptimizer

use of org.hipparchus.optim.nonlinear.vector.leastsquares.LevenbergMarquardtOptimizer 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 13 with LevenbergMarquardtOptimizer

use of org.hipparchus.optim.nonlinear.vector.leastsquares.LevenbergMarquardtOptimizer in project Orekit by CS-SI.

the class BiasTest method testEstimateBias.

@SuppressWarnings("unchecked")
@Test
public void testEstimateBias() throws OrekitException {
    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);
    // create perfect range measurements
    final Propagator propagator = EstimationTestUtils.createPropagator(context.initialOrbit, propagatorBuilder);
    final List<ObservedMeasurement<?>> measurements = EstimationTestUtils.createMeasurements(propagator, new RangeMeasurementCreator(context), 1.0, 3.0, 300.0);
    // create range biases: one bias for each station
    final RandomGenerator random = new Well19937a(0x0c4b69da5d64b35al);
    final Bias<?>[] stationsRangeBiases = new Bias<?>[context.stations.size()];
    final double[] realStationsBiases = new double[context.stations.size()];
    for (int i = 0; i < context.stations.size(); ++i) {
        final TopocentricFrame base = context.stations.get(i).getBaseFrame();
        stationsRangeBiases[i] = new Bias<Range>(new String[] { base.getName() + " range bias" }, new double[] { 0.0 }, new double[] { 1.0 }, new double[] { Double.NEGATIVE_INFINITY }, new double[] { Double.POSITIVE_INFINITY });
        realStationsBiases[i] = 2 * random.nextDouble() - 1;
    }
    // create orbit estimator
    final BatchLSEstimator estimator = new BatchLSEstimator(new LevenbergMarquardtOptimizer(), propagatorBuilder);
    // add the measurements, with both spacecraft and stations biases
    for (final ObservedMeasurement<?> measurement : measurements) {
        final Range range = (Range) measurement;
        for (int i = 0; i < context.stations.size(); ++i) {
            if (range.getStation() == context.stations.get(i)) {
                double biasedRange = range.getObservedValue()[0] + realStationsBiases[i];
                final Range m = new Range(range.getStation(), range.getDate(), biasedRange, range.getTheoreticalStandardDeviation()[0], range.getBaseWeight()[0]);
                m.addModifier((Bias<Range>) stationsRangeBiases[i]);
                estimator.addMeasurement(m);
            }
        }
    }
    estimator.setParametersConvergenceThreshold(1.0e-3);
    estimator.setMaxIterations(10);
    estimator.setMaxEvaluations(20);
    // we want to estimate the biases
    for (Bias<?> bias : stationsRangeBiases) {
        for (final ParameterDriver driver : bias.getParametersDrivers()) {
            driver.setSelected(true);
        }
    }
    EstimationTestUtils.checkFit(context, estimator, 2, 3, 0.0, 7.2e-7, 0.0, 2.1e-6, 0.0, 3.7e-7, 0.0, 1.7e-10);
    for (int i = 0; i < stationsRangeBiases.length; ++i) {
        Assert.assertEquals(realStationsBiases[i], stationsRangeBiases[i].getParametersDrivers().get(0).getValue(), 3.3e-6);
    }
}
Also used : Context(org.orekit.estimation.Context) Bias(org.orekit.estimation.measurements.modifiers.Bias) TopocentricFrame(org.orekit.frames.TopocentricFrame) Well19937a(org.hipparchus.random.Well19937a) Range(org.orekit.estimation.measurements.Range) ParameterDriver(org.orekit.utils.ParameterDriver) RandomGenerator(org.hipparchus.random.RandomGenerator) 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) RangeMeasurementCreator(org.orekit.estimation.measurements.RangeMeasurementCreator) ObservedMeasurement(org.orekit.estimation.measurements.ObservedMeasurement) Test(org.junit.Test)

Example 14 with LevenbergMarquardtOptimizer

use of org.hipparchus.optim.nonlinear.vector.leastsquares.LevenbergMarquardtOptimizer in project Orekit by CS-SI.

the class BatchLSEstimatorTest method testMultiSat.

@Test
public void testMultiSat() throws OrekitException {
    Context context = EstimationTestUtils.eccentricContext("regular-data:potential:tides");
    final NumericalPropagatorBuilder propagatorBuilder1 = context.createBuilder(OrbitType.KEPLERIAN, PositionAngle.TRUE, true, 1.0e-6, 60.0, 1.0);
    final NumericalPropagatorBuilder propagatorBuilder2 = context.createBuilder(OrbitType.KEPLERIAN, PositionAngle.TRUE, true, 1.0e-6, 60.0, 1.0);
    // Create perfect inter-satellites range measurements
    final TimeStampedPVCoordinates original = context.initialOrbit.getPVCoordinates();
    final Orbit closeOrbit = new CartesianOrbit(new TimeStampedPVCoordinates(context.initialOrbit.getDate(), original.getPosition().add(new Vector3D(1000, 2000, 3000)), original.getVelocity().add(new Vector3D(-0.03, 0.01, 0.02))), context.initialOrbit.getFrame(), context.initialOrbit.getMu());
    final Propagator closePropagator = EstimationTestUtils.createPropagator(closeOrbit, propagatorBuilder2);
    closePropagator.setEphemerisMode();
    closePropagator.propagate(context.initialOrbit.getDate().shiftedBy(3.5 * closeOrbit.getKeplerianPeriod()));
    final BoundedPropagator ephemeris = closePropagator.getGeneratedEphemeris();
    Propagator propagator1 = EstimationTestUtils.createPropagator(context.initialOrbit, propagatorBuilder1);
    final List<ObservedMeasurement<?>> r12 = EstimationTestUtils.createMeasurements(propagator1, new InterSatellitesRangeMeasurementCreator(ephemeris), 1.0, 3.0, 300.0);
    // create perfect range measurements for first satellite
    propagator1 = EstimationTestUtils.createPropagator(context.initialOrbit, propagatorBuilder1);
    final List<ObservedMeasurement<?>> r1 = EstimationTestUtils.createMeasurements(propagator1, new RangeMeasurementCreator(context), 1.0, 3.0, 300.0);
    // create orbit estimator
    final BatchLSEstimator estimator = new BatchLSEstimator(new LevenbergMarquardtOptimizer(), propagatorBuilder1, propagatorBuilder2);
    for (final ObservedMeasurement<?> interSat : r12) {
        estimator.addMeasurement(interSat);
    }
    for (final ObservedMeasurement<?> range : r1) {
        estimator.addMeasurement(range);
    }
    estimator.setParametersConvergenceThreshold(1.0e-2);
    estimator.setMaxIterations(10);
    estimator.setMaxEvaluations(20);
    estimator.setObserver(new BatchLSObserver() {

        int lastIter = 0;

        int lastEval = 0;

        /**
         * {@inheritDoc}
         */
        @Override
        public void evaluationPerformed(int iterationsCount, int evaluationscount, Orbit[] orbits, ParameterDriversList estimatedOrbitalParameters, ParameterDriversList estimatedPropagatorParameters, ParameterDriversList estimatedMeasurementsParameters, EstimationsProvider evaluationsProvider, Evaluation lspEvaluation) throws OrekitException {
            if (iterationsCount == lastIter) {
                Assert.assertEquals(lastEval + 1, evaluationscount);
            } else {
                Assert.assertEquals(lastIter + 1, iterationsCount);
            }
            lastIter = iterationsCount;
            lastEval = evaluationscount;
            Assert.assertEquals(r12.size() + r1.size(), evaluationsProvider.getNumber());
            try {
                evaluationsProvider.getEstimatedMeasurement(-1);
                Assert.fail("an exception should have been thrown");
            } catch (OrekitException oe) {
                Assert.assertEquals(LocalizedCoreFormats.OUT_OF_RANGE_SIMPLE, oe.getSpecifier());
            }
            try {
                evaluationsProvider.getEstimatedMeasurement(r12.size() + r1.size());
                Assert.fail("an exception should have been thrown");
            } catch (OrekitException oe) {
                Assert.assertEquals(LocalizedCoreFormats.OUT_OF_RANGE_SIMPLE, oe.getSpecifier());
            }
            AbsoluteDate previous = AbsoluteDate.PAST_INFINITY;
            for (int i = 0; i < evaluationsProvider.getNumber(); ++i) {
                AbsoluteDate current = evaluationsProvider.getEstimatedMeasurement(i).getDate();
                Assert.assertTrue(current.compareTo(previous) >= 0);
                previous = current;
            }
        }
    });
    List<DelegatingDriver> parameters = estimator.getOrbitalParametersDrivers(true).getDrivers();
    ParameterDriver a0Driver = parameters.get(0);
    Assert.assertEquals("a[0]", a0Driver.getName());
    a0Driver.setValue(a0Driver.getValue() + 1.2);
    a0Driver.setReferenceDate(AbsoluteDate.GALILEO_EPOCH);
    ParameterDriver a1Driver = parameters.get(6);
    Assert.assertEquals("a[1]", a1Driver.getName());
    a1Driver.setValue(a1Driver.getValue() - 5.4);
    a1Driver.setReferenceDate(AbsoluteDate.GALILEO_EPOCH);
    final Orbit before = new KeplerianOrbit(parameters.get(6).getValue(), parameters.get(7).getValue(), parameters.get(8).getValue(), parameters.get(9).getValue(), parameters.get(10).getValue(), parameters.get(11).getValue(), PositionAngle.TRUE, closeOrbit.getFrame(), closeOrbit.getDate(), closeOrbit.getMu());
    Assert.assertEquals(4.7246, Vector3D.distance(closeOrbit.getPVCoordinates().getPosition(), before.getPVCoordinates().getPosition()), 1.0e-3);
    Assert.assertEquals(0.0010514, Vector3D.distance(closeOrbit.getPVCoordinates().getVelocity(), before.getPVCoordinates().getVelocity()), 1.0e-6);
    EstimationTestUtils.checkFit(context, estimator, 2, 3, 0.0, 2.3e-06, 0.0, 6.6e-06, 0.0, 6.2e-07, 0.0, 2.8e-10);
    final Orbit determined = new KeplerianOrbit(parameters.get(6).getValue(), parameters.get(7).getValue(), parameters.get(8).getValue(), parameters.get(9).getValue(), parameters.get(10).getValue(), parameters.get(11).getValue(), PositionAngle.TRUE, closeOrbit.getFrame(), closeOrbit.getDate(), closeOrbit.getMu());
    Assert.assertEquals(0.0, Vector3D.distance(closeOrbit.getPVCoordinates().getPosition(), determined.getPVCoordinates().getPosition()), 1.6e-6);
    Assert.assertEquals(0.0, Vector3D.distance(closeOrbit.getPVCoordinates().getVelocity(), determined.getPVCoordinates().getVelocity()), 1.6e-9);
    // got a default one
    for (final ParameterDriver driver : estimator.getOrbitalParametersDrivers(true).getDrivers()) {
        if (driver.getName().startsWith("a[")) {
            // user-specified reference date
            Assert.assertEquals(0, driver.getReferenceDate().durationFrom(AbsoluteDate.GALILEO_EPOCH), 1.0e-15);
        } else {
            // default reference date
            Assert.assertEquals(0, driver.getReferenceDate().durationFrom(propagatorBuilder1.getInitialOrbitDate()), 1.0e-15);
        }
    }
}
Also used : CartesianOrbit(org.orekit.orbits.CartesianOrbit) TimeStampedPVCoordinates(org.orekit.utils.TimeStampedPVCoordinates) AbsoluteDate(org.orekit.time.AbsoluteDate) Vector3D(org.hipparchus.geometry.euclidean.threed.Vector3D) ParameterDriversList(org.orekit.utils.ParameterDriversList) BoundedPropagator(org.orekit.propagation.BoundedPropagator) Propagator(org.orekit.propagation.Propagator) OrekitException(org.orekit.errors.OrekitException) KeplerianOrbit(org.orekit.orbits.KeplerianOrbit) DelegatingDriver(org.orekit.utils.ParameterDriversList.DelegatingDriver) BoundedPropagator(org.orekit.propagation.BoundedPropagator) ObservedMeasurement(org.orekit.estimation.measurements.ObservedMeasurement) EstimationsProvider(org.orekit.estimation.measurements.EstimationsProvider) Context(org.orekit.estimation.Context) Evaluation(org.hipparchus.optim.nonlinear.vector.leastsquares.LeastSquaresProblem.Evaluation) Orbit(org.orekit.orbits.Orbit) CartesianOrbit(org.orekit.orbits.CartesianOrbit) KeplerianOrbit(org.orekit.orbits.KeplerianOrbit) ParameterDriver(org.orekit.utils.ParameterDriver) LevenbergMarquardtOptimizer(org.hipparchus.optim.nonlinear.vector.leastsquares.LevenbergMarquardtOptimizer) NumericalPropagatorBuilder(org.orekit.propagation.conversion.NumericalPropagatorBuilder) InterSatellitesRangeMeasurementCreator(org.orekit.estimation.measurements.InterSatellitesRangeMeasurementCreator) RangeMeasurementCreator(org.orekit.estimation.measurements.RangeMeasurementCreator) InterSatellitesRangeMeasurementCreator(org.orekit.estimation.measurements.InterSatellitesRangeMeasurementCreator) Test(org.junit.Test)

Example 15 with LevenbergMarquardtOptimizer

use of org.hipparchus.optim.nonlinear.vector.leastsquares.LevenbergMarquardtOptimizer in project Orekit by CS-SI.

the class BatchLSEstimatorTest method testKeplerRangeWithOnBoardAntennaOffset.

/**
 * Perfect range measurements with a biased start and an on-board antenna range offset
 * @throws OrekitException
 */
@Test
public void testKeplerRangeWithOnBoardAntennaOffset() throws OrekitException {
    Context context = EstimationTestUtils.eccentricContext("regular-data:potential:tides");
    final NumericalPropagatorBuilder propagatorBuilder = context.createBuilder(OrbitType.KEPLERIAN, PositionAngle.TRUE, true, 1.0e-6, 60.0, 1.0);
    propagatorBuilder.setAttitudeProvider(new LofOffset(propagatorBuilder.getFrame(), LOFType.LVLH));
    final Vector3D antennaPhaseCenter = new Vector3D(-1.2, 2.3, -0.7);
    // create perfect range measurements with antenna offset
    final Propagator propagator = EstimationTestUtils.createPropagator(context.initialOrbit, propagatorBuilder);
    final List<ObservedMeasurement<?>> measurements = EstimationTestUtils.createMeasurements(propagator, new RangeMeasurementCreator(context, antennaPhaseCenter), 1.0, 3.0, 300.0);
    // create orbit estimator
    final BatchLSEstimator estimator = new BatchLSEstimator(new LevenbergMarquardtOptimizer(), propagatorBuilder);
    final OnBoardAntennaRangeModifier obaModifier = new OnBoardAntennaRangeModifier(antennaPhaseCenter);
    for (final ObservedMeasurement<?> range : measurements) {
        ((Range) range).addModifier(obaModifier);
        estimator.addMeasurement(range);
    }
    estimator.setParametersConvergenceThreshold(1.0e-2);
    estimator.setMaxIterations(10);
    estimator.setMaxEvaluations(20);
    estimator.setObserver(new BatchLSObserver() {

        int lastIter = 0;

        int lastEval = 0;

        /**
         * {@inheritDoc}
         */
        @Override
        public void evaluationPerformed(int iterationsCount, int evaluationscount, Orbit[] orbits, ParameterDriversList estimatedOrbitalParameters, ParameterDriversList estimatedPropagatorParameters, ParameterDriversList estimatedMeasurementsParameters, EstimationsProvider evaluationsProvider, Evaluation lspEvaluation) throws OrekitException {
            if (iterationsCount == lastIter) {
                Assert.assertEquals(lastEval + 1, evaluationscount);
            } else {
                Assert.assertEquals(lastIter + 1, iterationsCount);
            }
            lastIter = iterationsCount;
            lastEval = evaluationscount;
            Assert.assertEquals(measurements.size(), evaluationsProvider.getNumber());
            try {
                evaluationsProvider.getEstimatedMeasurement(-1);
                Assert.fail("an exception should have been thrown");
            } catch (OrekitException oe) {
                Assert.assertEquals(LocalizedCoreFormats.OUT_OF_RANGE_SIMPLE, oe.getSpecifier());
            }
            try {
                evaluationsProvider.getEstimatedMeasurement(measurements.size());
                Assert.fail("an exception should have been thrown");
            } catch (OrekitException oe) {
                Assert.assertEquals(LocalizedCoreFormats.OUT_OF_RANGE_SIMPLE, oe.getSpecifier());
            }
            AbsoluteDate previous = AbsoluteDate.PAST_INFINITY;
            for (int i = 0; i < evaluationsProvider.getNumber(); ++i) {
                AbsoluteDate current = evaluationsProvider.getEstimatedMeasurement(i).getDate();
                Assert.assertTrue(current.compareTo(previous) >= 0);
                previous = current;
            }
        }
    });
    ParameterDriver aDriver = estimator.getOrbitalParametersDrivers(true).getDrivers().get(0);
    Assert.assertEquals("a", aDriver.getName());
    aDriver.setValue(aDriver.getValue() + 1.2);
    aDriver.setReferenceDate(AbsoluteDate.GALILEO_EPOCH);
    EstimationTestUtils.checkFit(context, estimator, 2, 3, 0.0, 2.0e-5, 0.0, 5.2e-5, 0.0, 2.7e-5, 0.0, 1.1e-8);
    // got a default one
    for (final ParameterDriver driver : estimator.getOrbitalParametersDrivers(true).getDrivers()) {
        if ("a".equals(driver.getName())) {
            // user-specified reference date
            Assert.assertEquals(0, driver.getReferenceDate().durationFrom(AbsoluteDate.GALILEO_EPOCH), 1.0e-15);
        } else {
            // default reference date
            Assert.assertEquals(0, driver.getReferenceDate().durationFrom(propagatorBuilder.getInitialOrbitDate()), 1.0e-15);
        }
    }
}
Also used : Context(org.orekit.estimation.Context) Evaluation(org.hipparchus.optim.nonlinear.vector.leastsquares.LeastSquaresProblem.Evaluation) Orbit(org.orekit.orbits.Orbit) CartesianOrbit(org.orekit.orbits.CartesianOrbit) KeplerianOrbit(org.orekit.orbits.KeplerianOrbit) Range(org.orekit.estimation.measurements.Range) ParameterDriver(org.orekit.utils.ParameterDriver) AbsoluteDate(org.orekit.time.AbsoluteDate) LevenbergMarquardtOptimizer(org.hipparchus.optim.nonlinear.vector.leastsquares.LevenbergMarquardtOptimizer) OnBoardAntennaRangeModifier(org.orekit.estimation.measurements.modifiers.OnBoardAntennaRangeModifier) Vector3D(org.hipparchus.geometry.euclidean.threed.Vector3D) ParameterDriversList(org.orekit.utils.ParameterDriversList) NumericalPropagatorBuilder(org.orekit.propagation.conversion.NumericalPropagatorBuilder) BoundedPropagator(org.orekit.propagation.BoundedPropagator) Propagator(org.orekit.propagation.Propagator) OrekitException(org.orekit.errors.OrekitException) RangeMeasurementCreator(org.orekit.estimation.measurements.RangeMeasurementCreator) InterSatellitesRangeMeasurementCreator(org.orekit.estimation.measurements.InterSatellitesRangeMeasurementCreator) LofOffset(org.orekit.attitudes.LofOffset) ObservedMeasurement(org.orekit.estimation.measurements.ObservedMeasurement) EstimationsProvider(org.orekit.estimation.measurements.EstimationsProvider) Test(org.junit.Test)

Aggregations

LevenbergMarquardtOptimizer (org.hipparchus.optim.nonlinear.vector.leastsquares.LevenbergMarquardtOptimizer)15 Test (org.junit.Test)13 NumericalPropagatorBuilder (org.orekit.propagation.conversion.NumericalPropagatorBuilder)13 Context (org.orekit.estimation.Context)12 Propagator (org.orekit.propagation.Propagator)12 ObservedMeasurement (org.orekit.estimation.measurements.ObservedMeasurement)11 BoundedPropagator (org.orekit.propagation.BoundedPropagator)9 RangeMeasurementCreator (org.orekit.estimation.measurements.RangeMeasurementCreator)7 ParameterDriver (org.orekit.utils.ParameterDriver)7 InterSatellitesRangeMeasurementCreator (org.orekit.estimation.measurements.InterSatellitesRangeMeasurementCreator)6 CartesianOrbit (org.orekit.orbits.CartesianOrbit)6 KeplerianOrbit (org.orekit.orbits.KeplerianOrbit)6 Orbit (org.orekit.orbits.Orbit)6 AbsoluteDate (org.orekit.time.AbsoluteDate)6 Evaluation (org.hipparchus.optim.nonlinear.vector.leastsquares.LeastSquaresProblem.Evaluation)5 BatchLSEstimator (org.orekit.estimation.leastsquares.BatchLSEstimator)5 EstimationsProvider (org.orekit.estimation.measurements.EstimationsProvider)5 ParameterDriversList (org.orekit.utils.ParameterDriversList)5 Vector3D (org.hipparchus.geometry.euclidean.threed.Vector3D)4 OrekitException (org.orekit.errors.OrekitException)4