1
2
3
4
5
6
7
8
9
10
11
12
13
14
15
16
17 package org.orekit.estimation.leastsquares;
18
19 import java.util.ArrayList;
20 import java.util.List;
21
22 import org.hipparchus.exception.LocalizedCoreFormats;
23 import org.hipparchus.geometry.euclidean.threed.Vector3D;
24 import org.hipparchus.linear.RealMatrix;
25 import org.hipparchus.optim.nonlinear.vector.leastsquares.LeastSquaresProblem.Evaluation;
26 import org.hipparchus.optim.nonlinear.vector.leastsquares.LevenbergMarquardtOptimizer;
27 import org.junit.Assert;
28 import org.junit.Test;
29 import org.orekit.attitudes.LofOffset;
30 import org.orekit.errors.OrekitException;
31 import org.orekit.errors.OrekitMessages;
32 import org.orekit.estimation.DSSTContext;
33 import org.orekit.estimation.DSSTEstimationTestUtils;
34 import org.orekit.estimation.measurements.EstimationsProvider;
35 import org.orekit.estimation.measurements.GroundStation;
36 import org.orekit.estimation.measurements.ObservedMeasurement;
37 import org.orekit.estimation.measurements.PVMeasurementCreator;
38 import org.orekit.estimation.measurements.Range;
39 import org.orekit.estimation.measurements.RangeMeasurementCreator;
40 import org.orekit.estimation.measurements.RangeRateMeasurementCreator;
41 import org.orekit.estimation.measurements.modifiers.OnBoardAntennaRangeModifier;
42 import org.orekit.frames.LOFType;
43 import org.orekit.orbits.Orbit;
44 import org.orekit.propagation.Propagator;
45 import org.orekit.propagation.conversion.DSSTPropagatorBuilder;
46 import org.orekit.propagation.semianalytical.dsst.DSSTPropagator;
47 import org.orekit.propagation.semianalytical.dsst.forces.DSSTNewtonianAttraction;
48 import org.orekit.time.AbsoluteDate;
49 import org.orekit.utils.ParameterDriver;
50 import org.orekit.utils.ParameterDriversList;
51
52 public class DSSTBatchLSEstimatorTest {
53
54
55
56
57 @Test
58 public void testKeplerPV() {
59
60 DSSTContext context = DSSTEstimationTestUtils.eccentricContext("regular-data:potential:tides");
61
62 final DSSTPropagatorBuilder propagatorBuilder =
63 context.createBuilder(true, 60.0, 600.0, 1.0);
64
65
66 final Propagator propagator = DSSTEstimationTestUtils.createPropagator(context.initialOrbit,
67 propagatorBuilder);
68 final List<ObservedMeasurement<?>> measurements =
69 DSSTEstimationTestUtils.createMeasurements(propagator,
70 new PVMeasurementCreator(),
71 0.0, 1.0, 300.0);
72
73
74 final BatchLSEstimator estimator = new BatchLSEstimator(new LevenbergMarquardtOptimizer(),
75 propagatorBuilder);
76 for (final ObservedMeasurement<?> measurement : measurements) {
77 estimator.addMeasurement(measurement);
78 }
79 estimator.setParametersConvergenceThreshold(1.0e-2);
80 estimator.setMaxIterations(10);
81 estimator.setMaxEvaluations(20);
82
83 DSSTEstimationTestUtils.checkFit(context, estimator, 1, 2,
84 0.0, 4.8e-9,
85 0.0, 2.6e-8,
86 0.0, 8.9e-9,
87 0.0, 4.4e-12);
88
89 RealMatrix normalizedCovariances = estimator.getOptimum().getCovariances(1.0e-10);
90 RealMatrix physicalCovariances = estimator.getPhysicalCovariances(1.0e-10);
91 Assert.assertEquals(6, normalizedCovariances.getRowDimension());
92 Assert.assertEquals(6, normalizedCovariances.getColumnDimension());
93 Assert.assertEquals(6, physicalCovariances.getRowDimension());
94 Assert.assertEquals(6, physicalCovariances.getColumnDimension());
95 Assert.assertEquals(0.00258, physicalCovariances.getEntry(0, 0), 1.0e-5);
96
97 }
98
99
100 @Test
101 public void testKeplerPVBackward() {
102
103 DSSTContext context = DSSTEstimationTestUtils.eccentricContext("regular-data:potential:tides");
104
105 final DSSTPropagatorBuilder propagatorBuilder =
106 context.createBuilder(true, 60.0, 600.0, 1.0);
107
108
109 final Propagator propagator = DSSTEstimationTestUtils.createPropagator(context.initialOrbit,
110 propagatorBuilder);
111 final List<ObservedMeasurement<?>> measurements =
112 DSSTEstimationTestUtils.createMeasurements(propagator,
113 new PVMeasurementCreator(),
114 0.0, -1.0, 300.0);
115
116
117 final BatchLSEstimator estimator = new BatchLSEstimator(new LevenbergMarquardtOptimizer(),
118 propagatorBuilder);
119 for (final ObservedMeasurement<?> measurement : measurements) {
120 estimator.addMeasurement(measurement);
121 }
122 estimator.setParametersConvergenceThreshold(1.0e-2);
123 estimator.setMaxIterations(10);
124 estimator.setMaxEvaluations(20);
125
126 DSSTEstimationTestUtils.checkFit(context, estimator, 1, 3,
127 0.0, 4.8e-9,
128 0.0, 2.7e-8,
129 0.0, 3.9e-9,
130 0.0, 1.9e-12);
131
132 RealMatrix normalizedCovariances = estimator.getOptimum().getCovariances(1.0e-10);
133 RealMatrix physicalCovariances = estimator.getPhysicalCovariances(1.0e-10);
134 Assert.assertEquals(6, normalizedCovariances.getRowDimension());
135 Assert.assertEquals(6, normalizedCovariances.getColumnDimension());
136 Assert.assertEquals(6, physicalCovariances.getRowDimension());
137 Assert.assertEquals(6, physicalCovariances.getColumnDimension());
138 Assert.assertEquals(0.00258, physicalCovariances.getEntry(0, 0), 1.0e-5);
139
140 }
141
142
143
144
145 @Test
146 public void testKeplerRange() {
147
148 DSSTContext context = DSSTEstimationTestUtils.eccentricContext("regular-data:potential:tides");
149
150 final DSSTPropagatorBuilder propagatorBuilder =
151 context.createBuilder(true, 60.0, 600.0, 1.0);
152
153
154 final Propagator propagator = DSSTEstimationTestUtils.createPropagator(context.initialOrbit,
155 propagatorBuilder);
156 final List<ObservedMeasurement<?>> measurements =
157 DSSTEstimationTestUtils.createMeasurements(propagator,
158 new RangeMeasurementCreator(context),
159 1.0, 3.0, 300.0);
160
161
162 final BatchLSEstimator estimator = new BatchLSEstimator(new LevenbergMarquardtOptimizer(),
163 propagatorBuilder);
164 for (final ObservedMeasurement<?> range : measurements) {
165 estimator.addMeasurement(range);
166 }
167 estimator.setParametersConvergenceThreshold(1.0e-2);
168 estimator.setMaxIterations(10);
169 estimator.setMaxEvaluations(20);
170 estimator.setObserver(new BatchLSObserver() {
171 int lastIter = 0;
172 int lastEval = 0;
173
174 @Override
175 public void evaluationPerformed(int iterationsCount, int evaluationscount,
176 Orbit[] orbits,
177 ParameterDriversList estimatedOrbitalParameters,
178 ParameterDriversList estimatedPropagatorParameters,
179 ParameterDriversList estimatedMeasurementsParameters,
180 EstimationsProvider evaluationsProvider, Evaluation lspEvaluation) {
181 if (iterationsCount == lastIter) {
182 Assert.assertEquals(lastEval + 1, evaluationscount);
183 } else {
184 Assert.assertEquals(lastIter + 1, iterationsCount);
185 }
186 lastIter = iterationsCount;
187 lastEval = evaluationscount;
188 Assert.assertEquals(measurements.size(), evaluationsProvider.getNumber());
189 try {
190 evaluationsProvider.getEstimatedMeasurement(-1);
191 Assert.fail("an exception should have been thrown");
192 } catch (OrekitException oe) {
193 Assert.assertEquals(LocalizedCoreFormats.OUT_OF_RANGE_SIMPLE, oe.getSpecifier());
194 }
195 try {
196 evaluationsProvider.getEstimatedMeasurement(measurements.size());
197 Assert.fail("an exception should have been thrown");
198 } catch (OrekitException oe) {
199 Assert.assertEquals(LocalizedCoreFormats.OUT_OF_RANGE_SIMPLE, oe.getSpecifier());
200 }
201 AbsoluteDate previous = AbsoluteDate.PAST_INFINITY;
202 for (int i = 0; i < evaluationsProvider.getNumber(); ++i) {
203 AbsoluteDate current = evaluationsProvider.getEstimatedMeasurement(i).getDate();
204 Assert.assertTrue(current.compareTo(previous) >= 0);
205 previous = current;
206 }
207 }
208 });
209
210 ParameterDriver aDriver = estimator.getOrbitalParametersDrivers(true).getDrivers().get(0);
211 Assert.assertEquals("a", aDriver.getName());
212 aDriver.setValue(aDriver.getValue() + 1.2);
213 aDriver.setReferenceDate(AbsoluteDate.GALILEO_EPOCH);
214
215 DSSTEstimationTestUtils.checkFit(context, estimator, 2, 3,
216 0.0, 3.1e-6,
217 0.0, 5.7e-6,
218 0.0, 1.3e-6,
219 0.0, 5.2e-10);
220
221
222
223 for (final ParameterDriver driver : estimator.getOrbitalParametersDrivers(true).getDrivers()) {
224 if ("a".equals(driver.getName())) {
225
226 Assert.assertEquals(0, driver.getReferenceDate().durationFrom(AbsoluteDate.GALILEO_EPOCH), 1.0e-15);
227 } else {
228
229 Assert.assertEquals(0, driver.getReferenceDate().durationFrom(propagatorBuilder.getInitialOrbitDate()), 1.0e-15);
230 }
231 }
232
233 }
234
235
236
237
238 @Test
239 public void testKeplerRangeWithOnBoardAntennaOffset() {
240
241 DSSTContext context = DSSTEstimationTestUtils.eccentricContext("regular-data:potential:tides");
242
243 final DSSTPropagatorBuilder propagatorBuilder =
244 context.createBuilder(true, 60.0, 600.0, 1.0);
245 propagatorBuilder.setAttitudeProvider(new LofOffset(propagatorBuilder.getFrame(), LOFType.LVLH));
246 final Vector3D antennaPhaseCenter = new Vector3D(-1.2, 2.3, -0.7);
247
248
249 final Propagator propagator = DSSTEstimationTestUtils.createPropagator(context.initialOrbit,
250 propagatorBuilder);
251 final List<ObservedMeasurement<?>> measurements =
252 DSSTEstimationTestUtils.createMeasurements(propagator,
253 new RangeMeasurementCreator(context, antennaPhaseCenter),
254 1.0, 3.0, 300.0);
255
256
257 final BatchLSEstimator estimator = new BatchLSEstimator(new LevenbergMarquardtOptimizer(),
258 propagatorBuilder);
259 final OnBoardAntennaRangeModifier obaModifier = new OnBoardAntennaRangeModifier(antennaPhaseCenter);
260 for (final ObservedMeasurement<?> range : measurements) {
261 ((Range) range).addModifier(obaModifier);
262 estimator.addMeasurement(range);
263 }
264 estimator.setParametersConvergenceThreshold(1.0e-2);
265 estimator.setMaxIterations(10);
266 estimator.setMaxEvaluations(20);
267 estimator.setObserver(new BatchLSObserver() {
268 int lastIter = 0;
269 int lastEval = 0;
270
271 @Override
272 public void evaluationPerformed(int iterationsCount, int evaluationscount,
273 Orbit[] orbits,
274 ParameterDriversList estimatedOrbitalParameters,
275 ParameterDriversList estimatedPropagatorParameters,
276 ParameterDriversList estimatedMeasurementsParameters,
277 EstimationsProvider evaluationsProvider, Evaluation lspEvaluation) {
278 if (iterationsCount == lastIter) {
279 Assert.assertEquals(lastEval + 1, evaluationscount);
280 } else {
281 Assert.assertEquals(lastIter + 1, iterationsCount);
282 }
283 lastIter = iterationsCount;
284 lastEval = evaluationscount;
285 Assert.assertEquals(measurements.size(), evaluationsProvider.getNumber());
286 try {
287 evaluationsProvider.getEstimatedMeasurement(-1);
288 Assert.fail("an exception should have been thrown");
289 } catch (OrekitException oe) {
290 Assert.assertEquals(LocalizedCoreFormats.OUT_OF_RANGE_SIMPLE, oe.getSpecifier());
291 }
292 try {
293 evaluationsProvider.getEstimatedMeasurement(measurements.size());
294 Assert.fail("an exception should have been thrown");
295 } catch (OrekitException oe) {
296 Assert.assertEquals(LocalizedCoreFormats.OUT_OF_RANGE_SIMPLE, oe.getSpecifier());
297 }
298 AbsoluteDate previous = AbsoluteDate.PAST_INFINITY;
299 for (int i = 0; i < evaluationsProvider.getNumber(); ++i) {
300 AbsoluteDate current = evaluationsProvider.getEstimatedMeasurement(i).getDate();
301 Assert.assertTrue(current.compareTo(previous) >= 0);
302 previous = current;
303 }
304 }
305 });
306
307 ParameterDriver aDriver = estimator.getOrbitalParametersDrivers(true).getDrivers().get(0);
308 Assert.assertEquals("a", aDriver.getName());
309 aDriver.setValue(aDriver.getValue() + 1.2);
310 aDriver.setReferenceDate(AbsoluteDate.GALILEO_EPOCH);
311
312 DSSTEstimationTestUtils.checkFit(context, estimator, 2, 3,
313 0.0, 2.3e-5,
314 0.0, 5.9e-5,
315 0.0, 2.7e-5,
316 0.0, 1.1e-8);
317
318
319
320 for (final ParameterDriver driver : estimator.getOrbitalParametersDrivers(true).getDrivers()) {
321 if ("a".equals(driver.getName())) {
322
323 Assert.assertEquals(0, driver.getReferenceDate().durationFrom(AbsoluteDate.GALILEO_EPOCH), 1.0e-15);
324 } else {
325
326 Assert.assertEquals(0, driver.getReferenceDate().durationFrom(propagatorBuilder.getInitialOrbitDate()), 1.0e-15);
327 }
328 }
329
330 }
331
332 @Test
333 public void testWrappedException() {
334
335 DSSTContext context = DSSTEstimationTestUtils.eccentricContext("regular-data:potential:tides");
336
337 final DSSTPropagatorBuilder propagatorBuilder =
338 context.createBuilder(true, 60.0, 600.0, 1.0);
339
340
341 final Propagator propagator = DSSTEstimationTestUtils.createPropagator(context.initialOrbit,
342 propagatorBuilder);
343 final List<ObservedMeasurement<?>> measurements =
344 DSSTEstimationTestUtils.createMeasurements(propagator,
345 new RangeMeasurementCreator(context),
346 1.0, 3.0, 300.0);
347
348
349 final BatchLSEstimator estimator = new BatchLSEstimator(new LevenbergMarquardtOptimizer(),
350 propagatorBuilder);
351 for (final ObservedMeasurement<?> range : measurements) {
352 estimator.addMeasurement(range);
353 }
354 estimator.setParametersConvergenceThreshold(1.0e-2);
355 estimator.setMaxIterations(10);
356 estimator.setMaxEvaluations(20);
357 estimator.setObserver(new BatchLSObserver() {
358
359 @Override
360 public void evaluationPerformed(int iterationsCount, int evaluationscount,
361 Orbit[] orbits,
362 ParameterDriversList estimatedOrbitalParameters,
363 ParameterDriversList estimatedPropagatorParameters,
364 ParameterDriversList estimatedMeasurementsParameters,
365 EstimationsProvider evaluationsProvider, Evaluation lspEvaluation) throws DummyException {
366 throw new DummyException();
367 }
368 });
369
370 try {
371 DSSTEstimationTestUtils.checkFit(context, estimator, 3, 4,
372 0.0, 1.5e-6,
373 0.0, 3.2e-6,
374 0.0, 3.8e-7,
375 0.0, 1.5e-10);
376 Assert.fail("an exception should have been thrown");
377 } catch (DummyException de) {
378
379 }
380
381 }
382
383 private static class DummyException extends OrekitException {
384 private static final long serialVersionUID = 1L;
385 public DummyException() {
386 super(OrekitMessages.INTERNAL_ERROR);
387 }
388 }
389
390
391
392
393 @Test
394 public void testKeplerRangeRate() {
395
396 DSSTContext context = DSSTEstimationTestUtils.eccentricContext("regular-data:potential:tides");
397
398 final DSSTPropagatorBuilder propagatorBuilder =
399 context.createBuilder(true, 60.0, 600.0, 1.0);
400
401
402 final Propagator propagator = DSSTEstimationTestUtils.createPropagator(context.initialOrbit,
403 propagatorBuilder);
404 final double groundClockDrift = 4.8e-9;
405 for (final GroundStation station : context.stations) {
406 station.getClockDriftDriver().setValue(groundClockDrift);
407 }
408 final double satClkDrift = 3.2e-10;
409 final List<ObservedMeasurement<?>> measurements1 =
410 DSSTEstimationTestUtils.createMeasurements(propagator,
411 new RangeRateMeasurementCreator(context, false, satClkDrift),
412 1.0, 3.0, 300.0);
413
414 final List<ObservedMeasurement<?>> measurements = new ArrayList<ObservedMeasurement<?>>();
415 measurements.addAll(measurements1);
416
417
418 final BatchLSEstimator estimator = new BatchLSEstimator(new LevenbergMarquardtOptimizer(),
419 propagatorBuilder);
420 for (final ObservedMeasurement<?> rangerate : measurements) {
421 estimator.addMeasurement(rangerate);
422 }
423 estimator.setParametersConvergenceThreshold(1.0e-3);
424 estimator.setMaxIterations(10);
425 estimator.setMaxEvaluations(20);
426
427 DSSTEstimationTestUtils.checkFit(context, estimator, 1, 2,
428 0.0, 5.4e-7,
429 0.0, 1.2e-6,
430 0.0, 8.3e-4,
431 0.0, 4.5e-7);
432 }
433
434
435
436
437 @Test
438 public void testKeplerRangeAndRangeRate() {
439
440 DSSTContext context = DSSTEstimationTestUtils.eccentricContext("regular-data:potential:tides");
441
442 final DSSTPropagatorBuilder propagatorBuilder =
443 context.createBuilder(true, 60.0, 600.0, 1.0);
444
445
446 final Propagator propagator = DSSTEstimationTestUtils.createPropagator(context.initialOrbit,
447 propagatorBuilder);
448
449 final List<ObservedMeasurement<?>> measurementsRange =
450 DSSTEstimationTestUtils.createMeasurements(propagator,
451 new RangeMeasurementCreator(context),
452 1.0, 3.0, 300.0);
453 final double groundClockDrift = 4.8e-9;
454 for (final GroundStation station : context.stations) {
455 station.getClockDriftDriver().setValue(groundClockDrift);
456 }
457 final double satClkDrift = 3.2e-10;
458 final List<ObservedMeasurement<?>> measurementsRangeRate =
459 DSSTEstimationTestUtils.createMeasurements(propagator,
460 new RangeRateMeasurementCreator(context, false, satClkDrift),
461 1.0, 3.0, 300.0);
462
463
464 final List<ObservedMeasurement<?>> measurements = new ArrayList<ObservedMeasurement<?>>();
465 measurements.addAll(measurementsRange);
466 measurements.addAll(measurementsRangeRate);
467
468
469 final BatchLSEstimator estimator = new BatchLSEstimator(new LevenbergMarquardtOptimizer(),
470 propagatorBuilder);
471 for (final ObservedMeasurement<?> meas : measurements) {
472 estimator.addMeasurement(meas);
473 }
474 estimator.setParametersConvergenceThreshold(1.0e-3);
475 estimator.setMaxIterations(10);
476 estimator.setMaxEvaluations(20);
477
478
479 DSSTEstimationTestUtils.checkFit(context, estimator, 1, 3,
480 0.0, 4.8e-7,
481 0.0, 1.6e-6,
482 0.0, 4.4e-8,
483 0.0, 1.9e-11);
484 }
485
486 @Test
487 public void testIssue359() {
488 DSSTContext context = DSSTEstimationTestUtils.eccentricContext("regular-data:potential:tides");
489
490 final DSSTPropagatorBuilder propagatorBuilder =
491 context.createBuilder(true, 60.0, 600.0, 1.0);
492
493
494
495 propagatorBuilder.getPropagationParametersDrivers().getDrivers().get(0).setSelected(true);
496
497 final DSSTPropagator propagator = (DSSTPropagator) DSSTEstimationTestUtils.createPropagator(context.initialOrbit,
498 propagatorBuilder);
499 final List<ObservedMeasurement<?>> measurements =
500 DSSTEstimationTestUtils.createMeasurements(propagator,
501 new PVMeasurementCreator(),
502 0.0, 1.0, 300.0);
503
504
505 final BatchLSEstimator estimator = new BatchLSEstimator(new LevenbergMarquardtOptimizer(),
506 propagatorBuilder);
507 for (final ObservedMeasurement<?> measurement : measurements) {
508 estimator.addMeasurement(measurement);
509 }
510 ParameterDriversList estimatedParameters = estimator.getPropagatorParametersDrivers(true);
511
512 final String driverName = DSSTNewtonianAttraction.CENTRAL_ATTRACTION_COEFFICIENT;
513 Assert.assertTrue(propagator.getAllForceModels().get(0) instanceof DSSTNewtonianAttraction);
514 Assert.assertTrue(propagatorBuilder.getAllForceModels().get(0) instanceof DSSTNewtonianAttraction);
515 Assert.assertNotNull(estimatedParameters.findByName(driverName));
516 Assert.assertTrue(propagator.getAllForceModels().get(0).getParametersDrivers().get(0).isSelected());
517 Assert.assertTrue(propagatorBuilder.getAllForceModels().get(0).getParametersDrivers().get(0).isSelected());
518 }
519
520 }