use of org.hipparchus.optim.nonlinear.vector.leastsquares.LevenbergMarquardtOptimizer in project Orekit by CS-SI.
the class BatchLSEstimatorTest method testMultiSatWithParameters.
/**
* A modified version of the previous test with a selection of propagation drivers to estimate
* One common (ยต)
* Some specifics for each satellite (Cr and Ca)
*
* @throws OrekitException
*/
@Test
public void testMultiSatWithParameters() throws OrekitException {
// Test: Set the propagator drivers to estimate for each satellite
final boolean muEstimated = true;
final boolean crEstimated1 = true;
final boolean caEstimated1 = true;
final boolean crEstimated2 = true;
final boolean caEstimated2 = false;
// Builder sat 1
final 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, Force.POTENTIAL, Force.SOLAR_RADIATION_PRESSURE);
// Adding selection of parameters
String satName = "sat 1";
for (DelegatingDriver driver : propagatorBuilder1.getPropagationParametersDrivers().getDrivers()) {
if (driver.getName().equals("central attraction coefficient")) {
driver.setSelected(muEstimated);
}
if (driver.getName().equals(RadiationSensitive.REFLECTION_COEFFICIENT)) {
driver.setName(driver.getName() + " " + satName);
driver.setSelected(crEstimated1);
}
if (driver.getName().equals(RadiationSensitive.ABSORPTION_COEFFICIENT)) {
driver.setName(driver.getName() + " " + satName);
driver.setSelected(caEstimated1);
}
}
// Builder for sat 2
final Context context2 = EstimationTestUtils.eccentricContext("regular-data:potential:tides");
final NumericalPropagatorBuilder propagatorBuilder2 = context2.createBuilder(OrbitType.KEPLERIAN, PositionAngle.TRUE, true, 1.0e-6, 60.0, 1.0, Force.POTENTIAL, Force.SOLAR_RADIATION_PRESSURE);
// Adding selection of parameters
satName = "sat 2";
for (ParameterDriver driver : propagatorBuilder2.getPropagationParametersDrivers().getDrivers()) {
if (driver.getName().equals("central attraction coefficient")) {
driver.setSelected(muEstimated);
}
if (driver.getName().equals(RadiationSensitive.REFLECTION_COEFFICIENT)) {
driver.setName(driver.getName() + " " + satName);
driver.setSelected(crEstimated2);
}
if (driver.getName().equals(RadiationSensitive.ABSORPTION_COEFFICIENT)) {
driver.setName(driver.getName() + " " + satName);
driver.setSelected(caEstimated2);
}
}
// 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;
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, 4, 5, 0.0, 6.0e-06, 0.0, 1.7e-05, 0.0, 4.4e-07, 0.0, 1.7e-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()), 5.8e-6);
Assert.assertEquals(0.0, Vector3D.distance(closeOrbit.getPVCoordinates().getVelocity(), determined.getPVCoordinates().getVelocity()), 3.5e-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);
}
}
}
use of org.hipparchus.optim.nonlinear.vector.leastsquares.LevenbergMarquardtOptimizer in project Orekit by CS-SI.
the class BatchLSEstimatorTest method testKeplerPVBackward.
/**
* Test PV measurements generation and backward propagation in least-square orbit determination.
*/
@Test
public void testKeplerPVBackward() 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);
// 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);
// create orbit estimator
final BatchLSEstimator estimator = new BatchLSEstimator(new LevenbergMarquardtOptimizer(), propagatorBuilder);
for (final ObservedMeasurement<?> measurement : measurements) {
estimator.addMeasurement(measurement);
}
estimator.setParametersConvergenceThreshold(1.0e-2);
estimator.setMaxIterations(10);
estimator.setMaxEvaluations(20);
EstimationTestUtils.checkFit(context, estimator, 1, 2, 0.0, 8.3e-9, 0.0, 5.3e-8, 0.0, 5.6e-9, 0.0, 1.6e-12);
RealMatrix normalizedCovariances = estimator.getOptimum().getCovariances(1.0e-10);
RealMatrix physicalCovariances = estimator.getPhysicalCovariances(1.0e-10);
Assert.assertEquals(6, normalizedCovariances.getRowDimension());
Assert.assertEquals(6, normalizedCovariances.getColumnDimension());
Assert.assertEquals(6, physicalCovariances.getRowDimension());
Assert.assertEquals(6, physicalCovariances.getColumnDimension());
Assert.assertEquals(0.00258, physicalCovariances.getEntry(0, 0), 1.0e-5);
}
use of org.hipparchus.optim.nonlinear.vector.leastsquares.LevenbergMarquardtOptimizer in project Orekit by CS-SI.
the class BatchLSEstimatorTest method testWrappedException.
@Test
public void testWrappedException() 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);
// 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 orbit estimator
final BatchLSEstimator estimator = new BatchLSEstimator(new LevenbergMarquardtOptimizer(), propagatorBuilder);
for (final ObservedMeasurement<?> range : measurements) {
estimator.addMeasurement(range);
}
estimator.setParametersConvergenceThreshold(1.0e-2);
estimator.setMaxIterations(10);
estimator.setMaxEvaluations(20);
estimator.setObserver(new BatchLSObserver() {
/**
* {@inheritDoc}
*/
@Override
public void evaluationPerformed(int iterationsCount, int evaluationscount, Orbit[] orbits, ParameterDriversList estimatedOrbitalParameters, ParameterDriversList estimatedPropagatorParameters, ParameterDriversList estimatedMeasurementsParameters, EstimationsProvider evaluationsProvider, Evaluation lspEvaluation) throws DummyException {
throw new DummyException();
}
});
try {
EstimationTestUtils.checkFit(context, estimator, 3, 4, 0.0, 1.5e-6, 0.0, 3.2e-6, 0.0, 3.8e-7, 0.0, 1.5e-10);
Assert.fail("an exception should have been thrown");
} catch (DummyException de) {
// expected
}
}
use of org.hipparchus.optim.nonlinear.vector.leastsquares.LevenbergMarquardtOptimizer in project Orekit by CS-SI.
the class OrbitDetermination 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;
}
use of org.hipparchus.optim.nonlinear.vector.leastsquares.LevenbergMarquardtOptimizer in project Orekit by CS-SI.
the class BatchLSEstimatorTest method testKeplerRange.
/**
* Perfect range measurements with a biased start
* @throws OrekitException
*/
@Test
public void testKeplerRange() 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);
// 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 orbit estimator
final BatchLSEstimator estimator = new BatchLSEstimator(new LevenbergMarquardtOptimizer(), propagatorBuilder);
for (final ObservedMeasurement<?> range : measurements) {
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, 1.1e-6, 0.0, 2.8e-6, 0.0, 4.0e-7, 0.0, 2.2e-10);
// 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);
}
}
}
Aggregations