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 org.hipparchus.exception.LocalizedCoreFormats;
20 import org.hipparchus.geometry.euclidean.threed.Vector3D;
21 import org.hipparchus.linear.RealMatrix;
22 import org.hipparchus.optim.nonlinear.vector.leastsquares.LeastSquaresProblem.Evaluation;
23 import org.hipparchus.optim.nonlinear.vector.leastsquares.LevenbergMarquardtOptimizer;
24 import org.junit.jupiter.api.Assertions;
25 import org.junit.jupiter.api.Test;
26 import org.orekit.attitudes.LofOffset;
27 import org.orekit.errors.OrekitException;
28 import org.orekit.errors.OrekitMessages;
29 import org.orekit.estimation.Context;
30 import org.orekit.estimation.EstimationTestUtils;
31 import org.orekit.estimation.Force;
32 import org.orekit.estimation.measurements.AngularAzElMeasurementCreator;
33 import org.orekit.estimation.measurements.EstimationsProvider;
34 import org.orekit.estimation.measurements.GroundStation;
35 import org.orekit.estimation.measurements.InterSatellitesRangeMeasurementCreator;
36 import org.orekit.estimation.measurements.MultiplexedMeasurement;
37 import org.orekit.estimation.measurements.ObservedMeasurement;
38 import org.orekit.estimation.measurements.PVMeasurementCreator;
39 import org.orekit.estimation.measurements.Range;
40 import org.orekit.estimation.measurements.TwoWayRangeMeasurementCreator;
41 import org.orekit.estimation.measurements.RangeRateMeasurementCreator;
42 import org.orekit.estimation.measurements.modifiers.PhaseCentersRangeModifier;
43 import org.orekit.forces.gravity.NewtonianAttraction;
44 import org.orekit.forces.radiation.RadiationSensitive;
45 import org.orekit.frames.LOFType;
46 import org.orekit.gnss.antenna.FrequencyPattern;
47 import org.orekit.orbits.CartesianOrbit;
48 import org.orekit.orbits.KeplerianOrbit;
49 import org.orekit.orbits.Orbit;
50 import org.orekit.orbits.OrbitType;
51 import org.orekit.orbits.PositionAngleType;
52 import org.orekit.propagation.BoundedPropagator;
53 import org.orekit.propagation.EphemerisGenerator;
54 import org.orekit.propagation.Propagator;
55 import org.orekit.propagation.conversion.NumericalPropagatorBuilder;
56 import org.orekit.propagation.numerical.NumericalPropagator;
57 import org.orekit.time.AbsoluteDate;
58 import org.orekit.time.ChronologicalComparator;
59 import org.orekit.utils.ParameterDriver;
60 import org.orekit.utils.ParameterDriversList;
61 import org.orekit.utils.ParameterDriversList.DelegatingDriver;
62 import org.orekit.utils.TimeStampedPVCoordinates;
63
64 import java.util.ArrayList;
65 import java.util.List;
66 import java.util.stream.IntStream;
67
68 class BatchLSEstimatorTest {
69
70
71
72
73 @Test
74 void testKeplerPVMultipleDrag() {
75
76 Context context = EstimationTestUtils.eccentricContext("regular-data:potential:tides");
77
78 final NumericalPropagatorBuilder propagatorBuilder =
79 context.createBuilder(OrbitType.KEPLERIAN, PositionAngleType.TRUE, true,
80 1.0e-6, 60.0, 1.0, Force.DRAG);
81
82 for (ParameterDriver driver:propagatorBuilder.getPropagationParametersDrivers().getDrivers()) {
83 if (driver.getName().equals("drag coefficient")) {
84 driver.setSelected(true);
85 driver.addSpanAtDate(context.initialOrbit.getDate());
86 }
87 }
88
89
90 final Propagator propagator = EstimationTestUtils.createPropagator(context.initialOrbit,
91 propagatorBuilder);
92 final List<ObservedMeasurement<?>> measurements =
93 EstimationTestUtils.createMeasurements(propagator,
94 new PVMeasurementCreator(),
95 -3.0, 3.0, 300.0);
96
97
98
99
100 final BatchLSEstimator estimator = new BatchLSEstimator(new LevenbergMarquardtOptimizer(),
101 propagatorBuilder);
102 measurements.forEach(estimator::addMeasurement);
103 estimator.setParametersConvergenceThreshold(1.0e-2);
104 estimator.setMaxIterations(10);
105 estimator.setMaxEvaluations(20);
106
107 EstimationTestUtils.checkFit(context, estimator, 1, 2,
108 0.0, 7.8e-8,
109 0.0, 6.0e-7,
110 0.0, 3.2e-7,
111 0.0, 1.3e-10);
112
113
114 List<DelegatingDriver> Orbparameters = estimator.getOrbitalParametersDrivers(true).getDrivers();
115 Assertions.assertEquals(context.initialOrbit.getA(), Orbparameters.get(0).getValue() , 1.0e-8);
116 Assertions.assertEquals(context.initialOrbit.getE(), Orbparameters.get(1).getValue() , 1.0e-12);
117 Assertions.assertEquals(context.initialOrbit.getI(), Orbparameters.get(2).getValue() , 1.0e-12);
118
119 RealMatrix jacobian = estimator.getOptimum().getJacobian();
120 Assertions.assertEquals(8, jacobian.getColumnDimension());
121
122 }
123
124
125
126
127 @Test
128 void testKeplerPV() {
129
130 Context context = EstimationTestUtils.eccentricContext("regular-data:potential:tides");
131
132 final NumericalPropagatorBuilder propagatorBuilder =
133 context.createBuilder(OrbitType.KEPLERIAN, PositionAngleType.TRUE, true,
134 1.0e-6, 60.0, 1.0);
135
136
137 final Propagator propagator = EstimationTestUtils.createPropagator(context.initialOrbit,
138 propagatorBuilder);
139 final List<ObservedMeasurement<?>> measurements =
140 EstimationTestUtils.createMeasurements(propagator,
141 new PVMeasurementCreator(),
142 0.0, 1.0, 300.0);
143
144
145 final BatchLSEstimator estimator = new BatchLSEstimator(new LevenbergMarquardtOptimizer(),
146 propagatorBuilder);
147 measurements.forEach(estimator::addMeasurement);
148 estimator.setParametersConvergenceThreshold(1.0e-2);
149 estimator.setMaxIterations(10);
150 estimator.setMaxEvaluations(20);
151
152 EstimationTestUtils.checkFit(context, estimator, 1, 4,
153 0.0, 2.2e-8,
154 0.0, 1.1e-7,
155 0.0, 1.4e-8,
156 0.0, 6.3e-12);
157
158 RealMatrix normalizedCovariances = estimator.getOptimum().getCovariances(1.0e-10);
159 RealMatrix physicalCovariances = estimator.getPhysicalCovariances(1.0e-10);
160 Assertions.assertEquals(6, normalizedCovariances.getRowDimension());
161 Assertions.assertEquals(6, normalizedCovariances.getColumnDimension());
162 Assertions.assertEquals(6, physicalCovariances.getRowDimension());
163 Assertions.assertEquals(6, physicalCovariances.getColumnDimension());
164 Assertions.assertEquals(0.00258, physicalCovariances.getEntry(0, 0), 1.0e-5);
165
166 }
167
168
169 @Test
170 void testKeplerPVBackward() {
171
172 Context context = EstimationTestUtils.eccentricContext("regular-data:potential:tides");
173
174 final NumericalPropagatorBuilder propagatorBuilder =
175 context.createBuilder(OrbitType.KEPLERIAN, PositionAngleType.TRUE, true,
176 1.0e-6, 60.0, 1.0);
177
178
179 final Propagator propagator = EstimationTestUtils.createPropagator(context.initialOrbit,
180 propagatorBuilder);
181 final List<ObservedMeasurement<?>> measurements =
182 EstimationTestUtils.createMeasurements(propagator,
183 new PVMeasurementCreator(),
184 0.0, -1.0, 300.0);
185
186
187 final BatchLSEstimator estimator = new BatchLSEstimator(new LevenbergMarquardtOptimizer(),
188 propagatorBuilder);
189 measurements.forEach(estimator::addMeasurement);
190 estimator.setParametersConvergenceThreshold(1.0e-2);
191 estimator.setMaxIterations(10);
192 estimator.setMaxEvaluations(20);
193
194 EstimationTestUtils.checkFit(context, estimator, 1, 2,
195 0.0, 8.3e-9,
196 0.0, 5.3e-8,
197 0.0, 5.6e-9,
198 0.0, 1.6e-12);
199
200 RealMatrix normalizedCovariances = estimator.getOptimum().getCovariances(1.0e-10);
201 RealMatrix physicalCovariances = estimator.getPhysicalCovariances(1.0e-10);
202 Assertions.assertEquals(6, normalizedCovariances.getRowDimension());
203 Assertions.assertEquals(6, normalizedCovariances.getColumnDimension());
204 Assertions.assertEquals(6, physicalCovariances.getRowDimension());
205 Assertions.assertEquals(6, physicalCovariances.getColumnDimension());
206 Assertions.assertEquals(0.00258, physicalCovariances.getEntry(0, 0), 1.0e-5);
207
208 }
209
210
211
212
213 @Test
214 void testKeplerRange() {
215
216 Context context = EstimationTestUtils.eccentricContext("regular-data:potential:tides");
217
218 final NumericalPropagatorBuilder propagatorBuilder =
219 context.createBuilder(OrbitType.KEPLERIAN, PositionAngleType.TRUE, true,
220 1.0e-6, 60.0, 1.0);
221
222
223 final Propagator propagator = EstimationTestUtils.createPropagator(context.initialOrbit,
224 propagatorBuilder);
225 final List<ObservedMeasurement<?>> measurements =
226 EstimationTestUtils.createMeasurements(propagator,
227 new TwoWayRangeMeasurementCreator(context),
228 1.0, 3.0, 300.0);
229
230
231 final BatchLSEstimator estimator = new BatchLSEstimator(new LevenbergMarquardtOptimizer(),
232 propagatorBuilder);
233 measurements.forEach(estimator::addMeasurement);
234 estimator.setParametersConvergenceThreshold(1.0e-2);
235 estimator.setMaxIterations(10);
236 estimator.setMaxEvaluations(20);
237 estimator.setObserver(new BatchLSObserver() {
238 int lastIter = 0;
239 int lastEval = 0;
240
241 @Override
242 public void evaluationPerformed(int iterationsCount, int evaluationscount,
243 Orbit[] orbits,
244 ParameterDriversList estimatedOrbitalParameters,
245 ParameterDriversList estimatedPropagatorParameters,
246 ParameterDriversList estimatedMeasurementsParameters,
247 EstimationsProvider evaluationsProvider, Evaluation lspEvaluation) {
248 if (iterationsCount == lastIter) {
249 Assertions.assertEquals(lastEval + 1, evaluationscount);
250 } else {
251 Assertions.assertEquals(lastIter + 1, iterationsCount);
252 }
253 lastIter = iterationsCount;
254 lastEval = evaluationscount;
255 Assertions.assertEquals(measurements.size(), evaluationsProvider.getNumber());
256 try {
257 evaluationsProvider.getEstimatedMeasurement(-1);
258 Assertions.fail("an exception should have been thrown");
259 } catch (OrekitException oe) {
260 Assertions.assertEquals(LocalizedCoreFormats.OUT_OF_RANGE_SIMPLE, oe.getSpecifier());
261 }
262 try {
263 evaluationsProvider.getEstimatedMeasurement(measurements.size());
264 Assertions.fail("an exception should have been thrown");
265 } catch (OrekitException oe) {
266 Assertions.assertEquals(LocalizedCoreFormats.OUT_OF_RANGE_SIMPLE, oe.getSpecifier());
267 }
268 AbsoluteDate previous = AbsoluteDate.PAST_INFINITY;
269 for (int i = 0; i < evaluationsProvider.getNumber(); ++i) {
270 AbsoluteDate current = evaluationsProvider.getEstimatedMeasurement(i).getDate();
271 Assertions.assertTrue(current.compareTo(previous) >= 0);
272 previous = current;
273 }
274 }
275 });
276
277 ParameterDriver aDriver = estimator.getOrbitalParametersDrivers(true).getDrivers().get(0);
278 Assertions.assertEquals("a", aDriver.getName());
279 aDriver.setValue(aDriver.getValue() + 1.2);
280 aDriver.setReferenceDate(AbsoluteDate.GALILEO_EPOCH);
281
282 EstimationTestUtils.checkFit(context, estimator, 2, 3,
283 0.0, 1.2e-6,
284 0.0, 2.8e-6,
285 0.0, 7.0e-7,
286 0.0, 3e-10);
287
288
289
290 for (final ParameterDriver driver : estimator.getOrbitalParametersDrivers(true).getDrivers()) {
291 if ("a".equals(driver.getName())) {
292
293 Assertions.assertEquals(0, driver.getReferenceDate().durationFrom(AbsoluteDate.GALILEO_EPOCH), 1.0e-15);
294 } else {
295
296 Assertions.assertEquals(0,
297 driver.getReferenceDate().durationFrom(propagatorBuilder.getInitialOrbitDate()), 1.0e-15);
298 }
299 }
300
301 }
302
303
304
305
306 @Test
307 void testKeplerRangeWithOnBoardAntennaOffset() {
308
309 Context context = EstimationTestUtils.eccentricContext("regular-data:potential:tides");
310
311 final NumericalPropagatorBuilder propagatorBuilder =
312 context.createBuilder(OrbitType.KEPLERIAN, PositionAngleType.TRUE, true,
313 1.0e-6, 60.0, 1.0);
314 propagatorBuilder.setAttitudeProvider(new LofOffset(propagatorBuilder.getFrame(), LOFType.LVLH));
315 final Vector3D antennaPhaseCenter = new Vector3D(-1.2, 2.3, -0.7);
316
317
318 final Propagator propagator = EstimationTestUtils.createPropagator(context.initialOrbit,
319 propagatorBuilder);
320 final List<ObservedMeasurement<?>> measurements =
321 EstimationTestUtils.createMeasurements(propagator,
322 new TwoWayRangeMeasurementCreator(context,
323 Vector3D.ZERO, null,
324 antennaPhaseCenter, null,
325 0),
326 1.0, 3.0, 300.0);
327
328
329 final BatchLSEstimator estimator = new BatchLSEstimator(new LevenbergMarquardtOptimizer(),
330 propagatorBuilder);
331 final PhaseCentersRangeModifier obaModifier = new PhaseCentersRangeModifier(FrequencyPattern.ZERO_CORRECTION,
332 new FrequencyPattern(antennaPhaseCenter,
333 null));
334 measurements.forEach(m -> {
335 ((Range) m).addModifier(obaModifier);
336 estimator.addMeasurement(m);
337 });
338 estimator.setParametersConvergenceThreshold(1.0e-2);
339 estimator.setMaxIterations(10);
340 estimator.setMaxEvaluations(20);
341 estimator.setObserver(new BatchLSObserver() {
342 int lastIter = 0;
343 int lastEval = 0;
344
345 @Override
346 public void evaluationPerformed(int iterationsCount, int evaluationscount,
347 Orbit[] orbits,
348 ParameterDriversList estimatedOrbitalParameters,
349 ParameterDriversList estimatedPropagatorParameters,
350 ParameterDriversList estimatedMeasurementsParameters,
351 EstimationsProvider evaluationsProvider, Evaluation lspEvaluation) {
352 if (iterationsCount == lastIter) {
353 Assertions.assertEquals(lastEval + 1, evaluationscount);
354 } else {
355 Assertions.assertEquals(lastIter + 1, iterationsCount);
356 }
357 lastIter = iterationsCount;
358 lastEval = evaluationscount;
359 Assertions.assertEquals(measurements.size(), evaluationsProvider.getNumber());
360 try {
361 evaluationsProvider.getEstimatedMeasurement(-1);
362 Assertions.fail("an exception should have been thrown");
363 } catch (OrekitException oe) {
364 Assertions.assertEquals(LocalizedCoreFormats.OUT_OF_RANGE_SIMPLE, oe.getSpecifier());
365 }
366 try {
367 evaluationsProvider.getEstimatedMeasurement(measurements.size());
368 Assertions.fail("an exception should have been thrown");
369 } catch (OrekitException oe) {
370 Assertions.assertEquals(LocalizedCoreFormats.OUT_OF_RANGE_SIMPLE, oe.getSpecifier());
371 }
372 AbsoluteDate previous = AbsoluteDate.PAST_INFINITY;
373 for (int i = 0; i < evaluationsProvider.getNumber(); ++i) {
374 AbsoluteDate current = evaluationsProvider.getEstimatedMeasurement(i).getDate();
375 Assertions.assertTrue(current.compareTo(previous) >= 0);
376 previous = current;
377 }
378 }
379 });
380
381 ParameterDriver aDriver = estimator.getOrbitalParametersDrivers(true).getDrivers().get(0);
382 Assertions.assertEquals("a", aDriver.getName());
383 aDriver.setValue(aDriver.getValue() + 1.2);
384 aDriver.setReferenceDate(AbsoluteDate.GALILEO_EPOCH);
385
386 EstimationTestUtils.checkFit(context, estimator, 2, 3,
387 0.0, 2.0e-5,
388 0.0, 5.3e-5,
389 0.0, 2.7e-5,
390 0.0, 1.1e-8);
391
392
393
394 for (final ParameterDriver driver : estimator.getOrbitalParametersDrivers(true).getDrivers()) {
395 if ("a".equals(driver.getName())) {
396
397 Assertions.assertEquals(0, driver.getReferenceDate().durationFrom(AbsoluteDate.GALILEO_EPOCH), 1.0e-15);
398 } else {
399
400 Assertions.assertEquals(0,
401 driver.getReferenceDate().durationFrom(propagatorBuilder.getInitialOrbitDate()), 1.0e-15);
402 }
403 }
404
405 }
406
407 @Test
408 void testMultiSat() {
409
410 Context context = EstimationTestUtils.eccentricContext("regular-data:potential:tides");
411
412 final NumericalPropagatorBuilder propagatorBuilder1 =
413 context.createBuilder(OrbitType.KEPLERIAN, PositionAngleType.TRUE, true,
414 1.0e-6, 60.0, 1.0e-3);
415 final NumericalPropagatorBuilder propagatorBuilder2 =
416 context.createBuilder(OrbitType.KEPLERIAN, PositionAngleType.TRUE, true,
417 1.0e-6, 60.0, 1.0e-3);
418
419
420 final TimeStampedPVCoordinates original = context.initialOrbit.getPVCoordinates();
421 final Orbit closeOrbit = new CartesianOrbit(new TimeStampedPVCoordinates(context.initialOrbit.getDate(),
422 original.getPosition().add(new Vector3D(1000, 2000, 3000)),
423 original.getVelocity().add(new Vector3D(-0.03, 0.01, 0.02))),
424 context.initialOrbit.getFrame(),
425 context.initialOrbit.getMu());
426 final Propagator closePropagator = EstimationTestUtils.createPropagator(closeOrbit,
427 propagatorBuilder2);
428 final EphemerisGenerator generator = closePropagator.getEphemerisGenerator();
429 closePropagator.propagate(context.initialOrbit.getDate().shiftedBy(3.5 * closeOrbit.getKeplerianPeriod()));
430 final BoundedPropagator ephemeris = generator.getGeneratedEphemeris();
431 Propagator propagator1 = EstimationTestUtils.createPropagator(context.initialOrbit,
432 propagatorBuilder1);
433
434 final double localClockOffset = 0.137e-6;
435 final double remoteClockOffset = 469.0e-6;
436 final List<ObservedMeasurement<?>> r12 =
437 EstimationTestUtils.createMeasurements(propagator1,
438 new InterSatellitesRangeMeasurementCreator(ephemeris,
439 localClockOffset,
440 remoteClockOffset),
441 1.0, 3.0, 300.0);
442
443
444 propagator1 = EstimationTestUtils.createPropagator(context.initialOrbit,
445 propagatorBuilder1);
446 final List<ObservedMeasurement<?>> r1 =
447 EstimationTestUtils.createMeasurements(propagator1,
448 new TwoWayRangeMeasurementCreator(context),
449 1.0, 3.0, 300.0);
450
451
452 final BatchLSEstimator estimator = new BatchLSEstimator(new LevenbergMarquardtOptimizer(),
453 propagatorBuilder1,
454 propagatorBuilder2);
455 r1.forEach(estimator::addMeasurement);
456 r12.forEach(estimator::addMeasurement);
457 estimator.setParametersConvergenceThreshold(1.0e-2);
458 estimator.setMaxIterations(10);
459 estimator.setMaxEvaluations(20);
460 estimator.setObserver(new BatchLSObserver() {
461 int lastIter = 0;
462 int lastEval = 0;
463
464 @Override
465 public void evaluationPerformed(int iterationsCount, int evaluationscount,
466 Orbit[] orbits,
467 ParameterDriversList estimatedOrbitalParameters,
468 ParameterDriversList estimatedPropagatorParameters,
469 ParameterDriversList estimatedMeasurementsParameters,
470 EstimationsProvider evaluationsProvider, Evaluation lspEvaluation) {
471 if (iterationsCount == lastIter) {
472 Assertions.assertEquals(lastEval + 1, evaluationscount);
473 } else {
474 Assertions.assertEquals(lastIter + 1, iterationsCount);
475 }
476 lastIter = iterationsCount;
477 lastEval = evaluationscount;
478 Assertions.assertEquals(r12.size() + r1.size(), evaluationsProvider.getNumber());
479 try {
480 evaluationsProvider.getEstimatedMeasurement(-1);
481 Assertions.fail("an exception should have been thrown");
482 } catch (OrekitException oe) {
483 Assertions.assertEquals(LocalizedCoreFormats.OUT_OF_RANGE_SIMPLE, oe.getSpecifier());
484 }
485 try {
486 evaluationsProvider.getEstimatedMeasurement(r12.size() + r1.size());
487 Assertions.fail("an exception should have been thrown");
488 } catch (OrekitException oe) {
489 Assertions.assertEquals(LocalizedCoreFormats.OUT_OF_RANGE_SIMPLE, oe.getSpecifier());
490 }
491 AbsoluteDate previous = AbsoluteDate.PAST_INFINITY;
492 for (int i = 0; i < evaluationsProvider.getNumber(); ++i) {
493 AbsoluteDate current = evaluationsProvider.getEstimatedMeasurement(i).getDate();
494 Assertions.assertTrue(current.compareTo(previous) >= 0);
495 previous = current;
496 }
497 }
498 });
499
500 List<DelegatingDriver> parameters = estimator.getOrbitalParametersDrivers(true).getDrivers();
501 ParameterDriver a0Driver = parameters.get(0);
502 Assertions.assertEquals("a[0]", a0Driver.getName());
503 a0Driver.setValue(a0Driver.getValue() + 1.2);
504 a0Driver.setReferenceDate(AbsoluteDate.GALILEO_EPOCH);
505
506 ParameterDriver a1Driver = parameters.get(6);
507 Assertions.assertEquals("a[1]", a1Driver.getName());
508 a1Driver.setValue(a1Driver.getValue() - 5.4);
509 a1Driver.setReferenceDate(AbsoluteDate.GALILEO_EPOCH);
510
511 final Orbit before = new KeplerianOrbit(parameters.get( 6).getValue(),
512 parameters.get( 7).getValue(),
513 parameters.get( 8).getValue(),
514 parameters.get( 9).getValue(),
515 parameters.get(10).getValue(),
516 parameters.get(11).getValue(),
517 PositionAngleType.TRUE,
518 closeOrbit.getFrame(),
519 closeOrbit.getDate(),
520 closeOrbit.getMu());
521 Assertions.assertEquals(4.7246, Vector3D.distance(closeOrbit.getPosition(),
522 before.getPosition()), 1.0e-3);
523 Assertions.assertEquals(0.0010514, Vector3D.distance(closeOrbit.getPVCoordinates().getVelocity(),
524 before.getPVCoordinates().getVelocity()), 1.0e-6);
525 EstimationTestUtils.checkFit(context, estimator, 3, 4,
526 0.0, 5e-06,
527 0.0, 1.3e-05,
528 0.0, 8.3e-07,
529 0.0, 3.7e-10);
530
531 final Orbit determined = new KeplerianOrbit(parameters.get( 6).getValue(),
532 parameters.get( 7).getValue(),
533 parameters.get( 8).getValue(),
534 parameters.get( 9).getValue(),
535 parameters.get(10).getValue(),
536 parameters.get(11).getValue(),
537 PositionAngleType.TRUE,
538 closeOrbit.getFrame(),
539 closeOrbit.getDate(),
540 closeOrbit.getMu());
541 Assertions.assertEquals(0.0, Vector3D.distance(closeOrbit.getPosition(),
542 determined.getPosition()), 6.2e-6);
543 Assertions.assertEquals(0.0, Vector3D.distance(closeOrbit.getPVCoordinates().getVelocity(),
544 determined.getPVCoordinates().getVelocity()), 1.6e-9);
545
546
547
548 for (final ParameterDriver driver : estimator.getOrbitalParametersDrivers(true).getDrivers()) {
549 if (driver.getName().startsWith("a[")) {
550
551 Assertions.assertEquals(0, driver.getReferenceDate().durationFrom(AbsoluteDate.GALILEO_EPOCH), 1.0e-15);
552 } else {
553
554 Assertions.assertEquals(0,
555 driver.getReferenceDate().durationFrom(propagatorBuilder1.getInitialOrbitDate()), 1.0e-15);
556 }
557 }
558
559 }
560
561
562
563
564
565
566 @Test
567 void testMultiSatWithParameters() {
568
569
570 final boolean muEstimated = true;
571 final boolean crEstimated1 = false;
572 final boolean caEstimated1 = true;
573 final boolean crEstimated2 = true;
574 final boolean caEstimated2 = false;
575
576
577
578 final Context context = EstimationTestUtils.eccentricContext("regular-data:potential:tides");
579 final NumericalPropagatorBuilder propagatorBuilder1 =
580 context.createBuilder(OrbitType.KEPLERIAN, PositionAngleType.TRUE, true,
581 1.0e-6, 60.0, 1.0e-3, Force.POTENTIAL, Force.SOLAR_RADIATION_PRESSURE);
582
583
584 String satName = "sat 1";
585 for (DelegatingDriver driver:propagatorBuilder1.getPropagationParametersDrivers().getDrivers()) {
586 if (driver.getName().equals("central attraction coefficient")) {
587 driver.setSelected(muEstimated);
588 }
589 if (driver.getName().equals(RadiationSensitive.REFLECTION_COEFFICIENT)) {
590 driver.setName(driver.getName() + " " + satName);
591 driver.setSelected(crEstimated1);
592 }
593 if (driver.getName().equals(RadiationSensitive.ABSORPTION_COEFFICIENT)) {
594 driver.setName(driver.getName() + " " + satName);
595 driver.setSelected(caEstimated1);
596 }
597 }
598
599
600 final Context context2 = EstimationTestUtils.eccentricContext("regular-data:potential:tides");
601 final NumericalPropagatorBuilder propagatorBuilder2 =
602 context2.createBuilder(OrbitType.KEPLERIAN, PositionAngleType.TRUE, true,
603 1.0e-6, 60.0, 1.0e-3, Force.POTENTIAL, Force.SOLAR_RADIATION_PRESSURE);
604
605
606 satName = "sat 2";
607 for (ParameterDriver driver:propagatorBuilder2.getPropagationParametersDrivers().getDrivers()) {
608 if (driver.getName().equals("central attraction coefficient")) {
609 driver.setSelected(muEstimated);
610 }
611 if (driver.getName().equals(RadiationSensitive.REFLECTION_COEFFICIENT)) {
612 driver.setName(driver.getName() + " " + satName);
613 driver.setSelected(crEstimated2);
614 }
615 if (driver.getName().equals(RadiationSensitive.ABSORPTION_COEFFICIENT)) {
616 driver.setName(driver.getName() + " " + satName);
617 driver.setSelected(caEstimated2);
618 }
619 }
620
621
622 final TimeStampedPVCoordinates original = context.initialOrbit.getPVCoordinates();
623 final Orbit closeOrbit = new CartesianOrbit(new TimeStampedPVCoordinates(context.initialOrbit.getDate(),
624 original.getPosition().add(new Vector3D(1000, 2000, 3000)),
625 original.getVelocity().add(new Vector3D(-0.03, 0.01, 0.02))),
626 context.initialOrbit.getFrame(),
627 context.initialOrbit.getMu());
628 final Propagator closePropagator = EstimationTestUtils.createPropagator(closeOrbit,
629 propagatorBuilder2);
630 final EphemerisGenerator generator = closePropagator.getEphemerisGenerator();
631 closePropagator.propagate(context.initialOrbit.getDate().shiftedBy(3.5 * closeOrbit.getKeplerianPeriod()));
632 final BoundedPropagator ephemeris = generator.getGeneratedEphemeris();
633 Propagator propagator1 = EstimationTestUtils.createPropagator(context.initialOrbit,
634 propagatorBuilder1);
635 final List<ObservedMeasurement<?>> r12 =
636 EstimationTestUtils.createMeasurements(propagator1,
637 new InterSatellitesRangeMeasurementCreator(ephemeris, 0., 0.),
638 1.0, 3.0, 120.0);
639
640
641 propagator1 = EstimationTestUtils.createPropagator(context.initialOrbit,
642 propagatorBuilder1);
643 final List<ObservedMeasurement<?>> r1 =
644 EstimationTestUtils.createMeasurements(propagator1,
645 new TwoWayRangeMeasurementCreator(context),
646 1.0, 3.0, 120.0);
647
648
649 propagator1 = EstimationTestUtils.createPropagator(context.initialOrbit,
650 propagatorBuilder1);
651 final List<ObservedMeasurement<?>> a1 =
652 EstimationTestUtils.createMeasurements(propagator1,
653 new AngularAzElMeasurementCreator(context),
654 1.0, 3.0, 120.0);
655
656
657 final BatchLSEstimator estimator = new BatchLSEstimator(new LevenbergMarquardtOptimizer(),
658 propagatorBuilder1,
659 propagatorBuilder2);
660 r12.forEach(estimator::addMeasurement);
661 r1.forEach(estimator::addMeasurement);
662 a1.forEach(estimator::addMeasurement);
663 estimator.setParametersConvergenceThreshold(1.0e-3);
664 estimator.setMaxIterations(10);
665 estimator.setMaxEvaluations(20);
666 estimator.setObserver(new BatchLSObserver() {
667 int lastIter = 0;
668 int lastEval = 0;
669
670 @Override
671 public void evaluationPerformed(int iterationsCount, int evaluationscount,
672 Orbit[] orbits,
673 ParameterDriversList estimatedOrbitalParameters,
674 ParameterDriversList estimatedPropagatorParameters,
675 ParameterDriversList estimatedMeasurementsParameters,
676 EstimationsProvider evaluationsProvider, Evaluation lspEvaluation) {
677 if (iterationsCount == lastIter) {
678 Assertions.assertEquals(lastEval + 1, evaluationscount);
679 } else {
680 Assertions.assertEquals(lastIter + 1, iterationsCount);
681 }
682 lastIter = iterationsCount;
683 lastEval = evaluationscount;
684
685 AbsoluteDate previous = AbsoluteDate.PAST_INFINITY;
686 for (int i = 0; i < evaluationsProvider.getNumber(); ++i) {
687 AbsoluteDate current = evaluationsProvider.getEstimatedMeasurement(i).getDate();
688 Assertions.assertTrue(current.compareTo(previous) >= 0);
689 previous = current;
690 }
691 }
692 });
693
694 List<DelegatingDriver> parameters = estimator.getOrbitalParametersDrivers(true).getDrivers();
695 ParameterDriver a0Driver = parameters.get(0);
696 Assertions.assertEquals("a[0]", a0Driver.getName());
697 a0Driver.setValue(a0Driver.getValue() + 1.2);
698 a0Driver.setReferenceDate(AbsoluteDate.GALILEO_EPOCH);
699
700 ParameterDriver a1Driver = parameters.get(6);
701 Assertions.assertEquals("a[1]", a1Driver.getName());
702 a1Driver.setValue(a1Driver.getValue() - 5.4, null);
703 a1Driver.setReferenceDate(AbsoluteDate.GALILEO_EPOCH);
704
705 final Orbit before = new KeplerianOrbit(parameters.get( 6).getValue(),
706 parameters.get( 7).getValue(),
707 parameters.get( 8).getValue(),
708 parameters.get( 9).getValue(),
709 parameters.get(10).getValue(),
710 parameters.get(11).getValue(),
711 PositionAngleType.TRUE,
712 closeOrbit.getFrame(),
713 closeOrbit.getDate(),
714 closeOrbit.getMu());
715 Assertions.assertEquals(4.7246, Vector3D.distance(closeOrbit.getPosition(),
716 before.getPosition()), 1.0e-3);
717 Assertions.assertEquals(0.0010514, Vector3D.distance(closeOrbit.getPVCoordinates().getVelocity(),
718 before.getPVCoordinates().getVelocity()), 1.0e-6);
719 EstimationTestUtils.checkFit(context, estimator, 5, 6,
720 0.0, 5.3e-06,
721 0.0, 1.4e-05,
722 0.0, 8.8e-07,
723 0.0, 3.6e-10);
724
725 final Orbit determined = new KeplerianOrbit(parameters.get( 6).getValue(),
726 parameters.get( 7).getValue(),
727 parameters.get( 8).getValue(),
728 parameters.get( 9).getValue(),
729 parameters.get(10).getValue(),
730 parameters.get(11).getValue(),
731 PositionAngleType.TRUE,
732 closeOrbit.getFrame(),
733 closeOrbit.getDate(),
734 closeOrbit.getMu());
735 Assertions.assertEquals(0.0, Vector3D.distance(closeOrbit.getPosition(),
736 determined.getPosition()), 5.3e-6);
737 Assertions.assertEquals(0.0, Vector3D.distance(closeOrbit.getPVCoordinates().getVelocity(),
738 determined.getPVCoordinates().getVelocity()), 2.9e-9);
739
740
741
742 for (final ParameterDriver driver : estimator.getOrbitalParametersDrivers(true).getDrivers()) {
743 if (driver.getName().startsWith("a[")) {
744
745 Assertions.assertEquals(0, driver.getReferenceDate().durationFrom(AbsoluteDate.GALILEO_EPOCH), 1.0e-15);
746 } else {
747
748 Assertions.assertEquals(0,
749 driver.getReferenceDate().durationFrom(propagatorBuilder1.getInitialOrbitDate()), 1.0e-15);
750 }
751 }
752
753 }
754
755
756
757
758
759 @Test
760 void testIssue617() {
761
762 Context context = EstimationTestUtils.eccentricContext("regular-data:potential:tides");
763
764 final NumericalPropagatorBuilder propagatorBuilder1 =
765 context.createBuilder(OrbitType.KEPLERIAN, PositionAngleType.TRUE, true,
766 1.0e-6, 60.0, 1.0);
767 final NumericalPropagatorBuilder propagatorBuilder2 =
768 context.createBuilder(OrbitType.KEPLERIAN, PositionAngleType.TRUE, true,
769 1.0e-6, 60.0, 1.0);
770
771
772 final TimeStampedPVCoordinates original = context.initialOrbit.getPVCoordinates();
773 final Orbit closeOrbit = new CartesianOrbit(new TimeStampedPVCoordinates(context.initialOrbit.getDate(),
774 original.getPosition().add(new Vector3D(1000, 2000, 3000)),
775 original.getVelocity().add(new Vector3D(-0.03, 0.01, 0.02))),
776 context.initialOrbit.getFrame(),
777 context.initialOrbit.getMu());
778 final Propagator closePropagator = EstimationTestUtils.createPropagator(closeOrbit,
779 propagatorBuilder2);
780 final EphemerisGenerator generator = closePropagator.getEphemerisGenerator();
781 closePropagator.propagate(context.initialOrbit.getDate().shiftedBy(3.5 * closeOrbit.getKeplerianPeriod()));
782 final BoundedPropagator ephemeris = generator.getGeneratedEphemeris();
783 Propagator propagator1 = EstimationTestUtils.createPropagator(context.initialOrbit,
784 propagatorBuilder1);
785
786 final List<ObservedMeasurement<?>> r12 =
787 EstimationTestUtils.createMeasurements(propagator1,
788 new InterSatellitesRangeMeasurementCreator(ephemeris, 0., 0.),
789 1.0, 3.0, 300.0);
790
791
792 propagator1 = EstimationTestUtils.createPropagator(context.initialOrbit,
793 propagatorBuilder1);
794 final List<ObservedMeasurement<?>> r1 =
795 EstimationTestUtils.createMeasurements(propagator1,
796 new TwoWayRangeMeasurementCreator(context),
797 1.0, 3.0, 300.0);
798
799
800 final BatchLSEstimator estimator = new BatchLSEstimator(new LevenbergMarquardtOptimizer(),
801 propagatorBuilder1,
802 propagatorBuilder2);
803
804 final List<ObservedMeasurement<?>> independentMeasurements = new ArrayList<>();
805 independentMeasurements.addAll(r1);
806 independentMeasurements.addAll(r12);
807
808
809
810
811
812 final List<ObservedMeasurement<?>> multiplexed = multiplexMeasurements(independentMeasurements, 60.0);
813
814 for (final ObservedMeasurement<?> measurement : multiplexed) {
815 estimator.addMeasurement(measurement);
816 }
817
818 estimator.setParametersConvergenceThreshold(1.0e-2);
819 estimator.setMaxIterations(10);
820 estimator.setMaxEvaluations(20);
821 estimator.setObserver(new BatchLSObserver() {
822 int lastIter = 0;
823 int lastEval = 0;
824
825 @Override
826 public void evaluationPerformed(int iterationsCount, int evaluationscount,
827 Orbit[] orbits,
828 ParameterDriversList estimatedOrbitalParameters,
829 ParameterDriversList estimatedPropagatorParameters,
830 ParameterDriversList estimatedMeasurementsParameters,
831 EstimationsProvider evaluationsProvider, Evaluation lspEvaluation) {
832 if (iterationsCount == lastIter) {
833 Assertions.assertEquals(lastEval + 1, evaluationscount);
834 } else {
835 Assertions.assertEquals(lastIter + 1, iterationsCount);
836 }
837 lastIter = iterationsCount;
838 lastEval = evaluationscount;
839 Assertions.assertEquals(multiplexed.size(), evaluationsProvider.getNumber());
840 try {
841 evaluationsProvider.getEstimatedMeasurement(-1);
842 Assertions.fail("an exception should have been thrown");
843 } catch (OrekitException oe) {
844 Assertions.assertEquals(LocalizedCoreFormats.OUT_OF_RANGE_SIMPLE, oe.getSpecifier());
845 }
846 try {
847 evaluationsProvider.getEstimatedMeasurement(r12.size() + r1.size());
848 Assertions.fail("an exception should have been thrown");
849 } catch (OrekitException oe) {
850 Assertions.assertEquals(LocalizedCoreFormats.OUT_OF_RANGE_SIMPLE, oe.getSpecifier());
851 }
852 AbsoluteDate previous = AbsoluteDate.PAST_INFINITY;
853 for (int i = 0; i < evaluationsProvider.getNumber(); ++i) {
854 AbsoluteDate current = evaluationsProvider.getEstimatedMeasurement(i).getDate();
855 Assertions.assertTrue(current.compareTo(previous) >= 0);
856 previous = current;
857 }
858 }
859 });
860
861 List<DelegatingDriver> parameters = estimator.getOrbitalParametersDrivers(true).getDrivers();
862 ParameterDriver a0Driver = parameters.get(0);
863 Assertions.assertEquals("a[0]", a0Driver.getName());
864 a0Driver.setValue(a0Driver.getValue() + 1.2);
865 a0Driver.setReferenceDate(AbsoluteDate.GALILEO_EPOCH);
866
867 ParameterDriver a1Driver = parameters.get(6);
868 Assertions.assertEquals("a[1]", a1Driver.getName());
869 a1Driver.setValue(a1Driver.getValue() - 5.4);
870 a1Driver.setReferenceDate(AbsoluteDate.GALILEO_EPOCH);
871
872 final Orbit before = new KeplerianOrbit(parameters.get( 6).getValue(),
873 parameters.get( 7).getValue(),
874 parameters.get( 8).getValue(),
875 parameters.get( 9).getValue(),
876 parameters.get(10).getValue(),
877 parameters.get(11).getValue(),
878 PositionAngleType.TRUE,
879 closeOrbit.getFrame(),
880 closeOrbit.getDate(),
881 closeOrbit.getMu());
882 Assertions.assertEquals(4.7246, Vector3D.distance(closeOrbit.getPosition(),
883 before.getPosition()), 1.0e-3);
884 Assertions.assertEquals(0.0010514, Vector3D.distance(closeOrbit.getPVCoordinates().getVelocity(),
885 before.getPVCoordinates().getVelocity()), 1.0e-6);
886 EstimationTestUtils.checkFit(context, estimator, 2, 3,
887 0.0, 2.9e-06,
888 0.0, 8.1e-06,
889 0.0, 7.1e-07,
890 0.0, 3.2e-10);
891
892 final Orbit determined = new KeplerianOrbit(parameters.get( 6).getValue(),
893 parameters.get( 7).getValue(),
894 parameters.get( 8).getValue(),
895 parameters.get( 9).getValue(),
896 parameters.get(10).getValue(),
897 parameters.get(11).getValue(),
898 PositionAngleType.TRUE,
899 closeOrbit.getFrame(),
900 closeOrbit.getDate(),
901 closeOrbit.getMu());
902 Assertions.assertEquals(0.0, Vector3D.distance(closeOrbit.getPosition(),
903 determined.getPosition()), 4.6e-6);
904 Assertions.assertEquals(0.0, Vector3D.distance(closeOrbit.getPVCoordinates().getVelocity(),
905 determined.getPVCoordinates().getVelocity()), 1.6e-9);
906
907
908
909 for (final ParameterDriver driver : estimator.getOrbitalParametersDrivers(true).getDrivers()) {
910 if (driver.getName().startsWith("a[")) {
911
912 Assertions.assertEquals(0, driver.getReferenceDate().durationFrom(AbsoluteDate.GALILEO_EPOCH), 1.0e-15);
913 } else {
914
915 Assertions.assertEquals(0,
916 driver.getReferenceDate().durationFrom(propagatorBuilder1.getInitialOrbitDate()), 1.0e-15);
917 }
918 }
919
920 }
921
922 @Test
923 void testWrappedException() {
924
925 Context context = EstimationTestUtils.eccentricContext("regular-data:potential:tides");
926
927 final NumericalPropagatorBuilder propagatorBuilder =
928 context.createBuilder(OrbitType.KEPLERIAN, PositionAngleType.TRUE, true,
929 1.0e-6, 60.0, 1.0);
930
931
932 final Propagator propagator = EstimationTestUtils.createPropagator(context.initialOrbit,
933 propagatorBuilder);
934 final List<ObservedMeasurement<?>> measurements =
935 EstimationTestUtils.createMeasurements(propagator,
936 new TwoWayRangeMeasurementCreator(context),
937 1.0, 3.0, 300.0);
938
939
940 final BatchLSEstimator estimator = new BatchLSEstimator(new LevenbergMarquardtOptimizer(),
941 propagatorBuilder);
942 measurements.forEach(estimator::addMeasurement);
943 estimator.setParametersConvergenceThreshold(1.0e-2);
944 estimator.setMaxIterations(10);
945 estimator.setMaxEvaluations(20);
946 estimator.setObserver(new BatchLSObserver() {
947
948 @Override
949 public void evaluationPerformed(int iterationsCount, int evaluationscount,
950 Orbit[] orbits,
951 ParameterDriversList estimatedOrbitalParameters,
952 ParameterDriversList estimatedPropagatorParameters,
953 ParameterDriversList estimatedMeasurementsParameters,
954 EstimationsProvider evaluationsProvider, Evaluation lspEvaluation) {
955 throw new DummyException();
956 }
957 });
958
959 try {
960 EstimationTestUtils.checkFit(context, estimator, 3, 4,
961 0.0, 1.5e-6,
962 0.0, 3.2e-6,
963 0.0, 3.8e-7,
964 0.0, 1.5e-10);
965 Assertions.fail("an exception should have been thrown");
966 } catch (DummyException de) {
967
968 }
969
970 }
971
972 private static class DummyException extends OrekitException {
973 private static final long serialVersionUID = 1L;
974 public DummyException() {
975 super(OrekitMessages.INTERNAL_ERROR);
976 }
977 }
978
979
980
981
982 @Test
983 void testKeplerRangeRate() {
984
985 Context context = EstimationTestUtils.eccentricContext("regular-data:potential:tides");
986
987 final NumericalPropagatorBuilder propagatorBuilder =
988 context.createBuilder(OrbitType.KEPLERIAN, PositionAngleType.TRUE, true,
989 1.0e-6, 60.0, 1.0);
990
991
992 final Propagator propagator = EstimationTestUtils.createPropagator(context.initialOrbit,
993 propagatorBuilder);
994 final double groundClockDrift = 4.8e-9;
995 for (final GroundStation station : context.stations) {
996 station.getClockDriftDriver().setValue(groundClockDrift);
997 }
998 final double satClkDrift = 3.2e-10;
999 final List<ObservedMeasurement<?>> measurements =
1000 EstimationTestUtils.createMeasurements(propagator,
1001 new RangeRateMeasurementCreator(context, false, satClkDrift),
1002 1.0, 3.0, 300.0);
1003
1004
1005
1006 final BatchLSEstimator estimator = new BatchLSEstimator(new LevenbergMarquardtOptimizer(),
1007 propagatorBuilder);
1008 for (final ObservedMeasurement<?> rangerate : measurements) {
1009 estimator.addMeasurement(rangerate);
1010 }
1011 estimator.setParametersConvergenceThreshold(1.0e-3);
1012 estimator.setMaxIterations(10);
1013 estimator.setMaxEvaluations(20);
1014
1015 EstimationTestUtils.checkFit(context, estimator, 1, 2,
1016 0.0, 5.3e-7,
1017 0.0, 1.3e-6,
1018 0.0, 8.4e-4,
1019 0.0, 5.1e-7);
1020 }
1021
1022
1023
1024
1025 @Test
1026 void testKeplerRangeAndRangeRate() {
1027
1028 Context context = EstimationTestUtils.eccentricContext("regular-data:potential:tides");
1029
1030 final NumericalPropagatorBuilder propagatorBuilder =
1031 context.createBuilder(OrbitType.KEPLERIAN, PositionAngleType.TRUE, true,
1032 1.0e-6, 60.0, 1.0);
1033
1034
1035 final Propagator propagator = EstimationTestUtils.createPropagator(context.initialOrbit,
1036 propagatorBuilder);
1037
1038 final List<ObservedMeasurement<?>> measurementsRange =
1039 EstimationTestUtils.createMeasurements(propagator,
1040 new TwoWayRangeMeasurementCreator(context),
1041 1.0, 3.0, 300.0);
1042 final double groundClockDrift = 4.8e-9;
1043 for (final GroundStation station : context.stations) {
1044 station.getClockDriftDriver().setValue(groundClockDrift);
1045 }
1046 final double satClkDrift = 3.2e-10;
1047 final List<ObservedMeasurement<?>> measurementsRangeRate =
1048 EstimationTestUtils.createMeasurements(propagator,
1049 new RangeRateMeasurementCreator(context, false, satClkDrift),
1050 1.0, 3.0, 300.0);
1051
1052
1053 final List<ObservedMeasurement<?>> measurements = new ArrayList<>();
1054 measurements.addAll(measurementsRange);
1055 measurements.addAll(measurementsRangeRate);
1056
1057
1058 final BatchLSEstimator estimator = new BatchLSEstimator(new LevenbergMarquardtOptimizer(),
1059 propagatorBuilder);
1060 for (final ObservedMeasurement<?> meas : measurements) {
1061 estimator.addMeasurement(meas);
1062 }
1063 estimator.setParametersConvergenceThreshold(1.0e-3);
1064 estimator.setMaxIterations(10);
1065 estimator.setMaxEvaluations(20);
1066
1067
1068 EstimationTestUtils.checkFit(context, estimator, 1, 2,
1069 0.0, 4.6e7,
1070 0.0, 1.8e-6,
1071 0.0, 5.8e-7,
1072 0.0, 2.7e-10);
1073 }
1074
1075
1076
1077
1078 @Test
1079 void testIssue359() {
1080
1081 Context context = EstimationTestUtils.eccentricContext("regular-data:potential:tides");
1082
1083 final NumericalPropagatorBuilder propagatorBuilder =
1084 context.createBuilder(OrbitType.KEPLERIAN, PositionAngleType.TRUE, true,
1085 1.0e-6, 60.0, 1.0);
1086
1087
1088
1089 propagatorBuilder.getPropagationParametersDrivers().getDrivers().get(0).setSelected(true);
1090
1091 final NumericalPropagator propagator = (NumericalPropagator) EstimationTestUtils.createPropagator(context.initialOrbit,
1092 propagatorBuilder);
1093 final List<ObservedMeasurement<?>> measurements =
1094 EstimationTestUtils.createMeasurements(propagator,
1095 new PVMeasurementCreator(),
1096 0.0, 1.0, 300.0);
1097
1098
1099 final BatchLSEstimator estimator = new BatchLSEstimator(new LevenbergMarquardtOptimizer(),
1100 propagatorBuilder);
1101
1102 measurements.forEach(estimator::addMeasurement);
1103
1104 ParameterDriversList estimatedParameters = estimator.getPropagatorParametersDrivers(true);
1105
1106 final String driverName = NewtonianAttraction.CENTRAL_ATTRACTION_COEFFICIENT;
1107 Assertions.assertInstanceOf(NewtonianAttraction.class, propagator.getAllForceModels().get(0));
1108 Assertions.assertInstanceOf(NewtonianAttraction.class, propagatorBuilder.getAllForceModels().get(0));
1109 Assertions.assertNotNull(estimatedParameters.findByName(driverName));
1110 Assertions.assertTrue(propagator.getAllForceModels().get(0).getParameterDriver(driverName).isSelected());
1111 Assertions.assertTrue(propagatorBuilder.getAllForceModels().get(0).getParameterDriver(driverName).isSelected());
1112 }
1113
1114 @Test
1115 void testEstimateOnlyOneOrbitalParameter() {
1116 doTestEstimateOnlySomeOrbitalParameters(new boolean[]{ true, false, false, false, false, false });
1117 doTestEstimateOnlySomeOrbitalParameters(new boolean[]{ false, true, false, false, false, false });
1118 doTestEstimateOnlySomeOrbitalParameters(new boolean[]{ false, false, true, false, false, false });
1119 doTestEstimateOnlySomeOrbitalParameters(new boolean[]{ false, false, false, true, false, false });
1120 doTestEstimateOnlySomeOrbitalParameters(new boolean[]{ false, false, false, false, true, false });
1121 doTestEstimateOnlySomeOrbitalParameters(new boolean[]{ false, false, false, false, false, true });
1122 }
1123
1124 @Test
1125 void testEstimateOnlyFewOrbitalParameters() {
1126 doTestEstimateOnlySomeOrbitalParameters(new boolean[]{ false, true, false, true, false, false });
1127 }
1128
1129 private void doTestEstimateOnlySomeOrbitalParameters(boolean[] orbitalParametersEstimated) {
1130
1131 Context context = EstimationTestUtils.eccentricContext("regular-data:potential:tides");
1132
1133 final NumericalPropagatorBuilder propagatorBuilder =
1134 context.createBuilder(OrbitType.KEPLERIAN, PositionAngleType.TRUE, true,
1135 1.0e-6, 60.0, 1.0e-3);
1136
1137
1138 final Propagator propagator = EstimationTestUtils.createPropagator(context.initialOrbit,
1139 propagatorBuilder);
1140 final List<ObservedMeasurement<?>> measurements =
1141 EstimationTestUtils.createMeasurements(propagator,
1142 new PVMeasurementCreator(),
1143 0.0, 1.0, 300.0);
1144
1145 List<ParameterDriversList.DelegatingDriver> orbitalElementsDrivers = propagatorBuilder.getOrbitalParametersDrivers().getDrivers();
1146 IntStream.range(0, orbitalParametersEstimated.length).forEach(i -> {
1147 final ParameterDriver driver = orbitalElementsDrivers.get(i);
1148 if (orbitalParametersEstimated[i]) {
1149 driver.setSelected(true);
1150 driver.setValue(1.0001 * driver.getReferenceValue());
1151 } else {
1152 driver.setSelected(false);
1153 }
1154 });
1155
1156
1157 final BatchLSEstimator estimator = new BatchLSEstimator(new LevenbergMarquardtOptimizer(),
1158 propagatorBuilder);
1159 measurements.forEach(estimator::addMeasurement);
1160 estimator.setParametersConvergenceThreshold(1.0e-2);
1161 estimator.setMaxIterations(10);
1162 estimator.setMaxEvaluations(20);
1163
1164
1165 estimator.estimate();
1166
1167
1168 IntStream.range(0, orbitalParametersEstimated.length).forEach(i -> {
1169 if (orbitalParametersEstimated[i]) {
1170 final ParameterDriver driver = orbitalElementsDrivers.get(i);
1171 Assertions.assertEquals(driver.getReferenceValue(), driver.getValue(),
1172 driver.getReferenceValue() * 1.0e-3);
1173 }
1174 });
1175
1176 }
1177
1178
1179
1180
1181
1182
1183 private List<ObservedMeasurement<?>> multiplexMeasurements(final List<ObservedMeasurement<?>> independentMeasurements,
1184 final double tol) {
1185 final List<ObservedMeasurement<?>> multiplexed = new ArrayList<>();
1186 independentMeasurements.sort(new ChronologicalComparator());
1187 List<ObservedMeasurement<?>> clump = new ArrayList<>();
1188 for (final ObservedMeasurement<?> measurement : independentMeasurements) {
1189 if (!clump.isEmpty() && measurement.getDate().durationFrom(clump.get(0).getDate()) > tol) {
1190
1191
1192 if (clump.size() == 1) {
1193 multiplexed.add(clump.get(0));
1194 } else {
1195 multiplexed.add(new MultiplexedMeasurement(clump));
1196 }
1197
1198
1199 clump = new ArrayList<>();
1200
1201 }
1202 clump.add(measurement);
1203 }
1204
1205 if (clump.size() == 1) {
1206 multiplexed.add(clump.get(0));
1207 } else {
1208 multiplexed.add(new MultiplexedMeasurement(clump));
1209 }
1210 return multiplexed;
1211 }
1212
1213 }
1214
1215