1
2
3
4
5
6
7
8
9
10
11
12
13
14
15
16
17 package org.orekit.estimation.measurements.modifiers;
18
19 import java.util.Collections;
20 import java.util.List;
21 import java.util.Map;
22
23 import org.hipparchus.RealFieldElement;
24 import org.hipparchus.util.FastMath;
25 import org.hipparchus.util.MathUtils;
26 import org.hipparchus.util.Precision;
27 import org.junit.Assert;
28 import org.junit.Before;
29 import org.junit.Test;
30 import org.orekit.errors.OrekitIllegalArgumentException;
31 import org.orekit.errors.OrekitMessages;
32 import org.orekit.estimation.Context;
33 import org.orekit.estimation.EstimationTestUtils;
34 import org.orekit.estimation.measurements.AngularAzEl;
35 import org.orekit.estimation.measurements.AngularAzElMeasurementCreator;
36 import org.orekit.estimation.measurements.EstimatedMeasurement;
37 import org.orekit.estimation.measurements.EstimationModifier;
38 import org.orekit.estimation.measurements.GroundStation;
39 import org.orekit.estimation.measurements.ObservedMeasurement;
40 import org.orekit.estimation.measurements.Range;
41 import org.orekit.estimation.measurements.RangeMeasurementCreator;
42 import org.orekit.estimation.measurements.RangeRate;
43 import org.orekit.estimation.measurements.RangeRateMeasurementCreator;
44 import org.orekit.estimation.measurements.TurnAroundRange;
45 import org.orekit.estimation.measurements.TurnAroundRangeMeasurementCreator;
46 import org.orekit.estimation.measurements.gnss.Phase;
47 import org.orekit.estimation.measurements.gnss.PhaseMeasurementCreator;
48 import org.orekit.frames.TopocentricFrame;
49 import org.orekit.gnss.Frequency;
50 import org.orekit.models.earth.ionosphere.IonosphericModel;
51 import org.orekit.models.earth.ionosphere.KlobucharIonoModel;
52 import org.orekit.orbits.OrbitType;
53 import org.orekit.orbits.PositionAngle;
54 import org.orekit.propagation.FieldSpacecraftState;
55 import org.orekit.propagation.Propagator;
56 import org.orekit.propagation.SpacecraftState;
57 import org.orekit.propagation.conversion.NumericalPropagatorBuilder;
58 import org.orekit.time.AbsoluteDate;
59 import org.orekit.utils.ParameterDriver;
60
61 public class IonoModifierTest {
62
63
64 private KlobucharIonoModel model;
65
66
67 private double frequency;
68
69 @Before
70 public void setUp() throws Exception {
71
72
73
74 model = new KlobucharIonoModel(new double[]{.3820e-07, .1490e-07, -.1790e-06, 0},
75 new double[]{.1430e+06, 0, -.3280e+06, .1130e+06});
76
77 frequency = Frequency.G01.getMHzFrequency() * 1.0e6;
78 }
79
80 @Test
81 public void testRangeIonoModifier() {
82
83 Context context = EstimationTestUtils.eccentricContext("regular-data:potential:tides");
84
85 final NumericalPropagatorBuilder propagatorBuilder =
86 context.createBuilder(OrbitType.KEPLERIAN, PositionAngle.TRUE, true,
87 1.0e-6, 60.0, 0.001);
88
89
90 for (final GroundStation station : context.stations) {
91 station.getClockOffsetDriver().setSelected(true);
92 station.getEastOffsetDriver().setSelected(true);
93 station.getNorthOffsetDriver().setSelected(true);
94 station.getZenithOffsetDriver().setSelected(true);
95 }
96 final Propagator propagator = EstimationTestUtils.createPropagator(context.initialOrbit,
97 propagatorBuilder);
98 final List<ObservedMeasurement<?>> measurements =
99 EstimationTestUtils.createMeasurements(propagator,
100 new RangeMeasurementCreator(context),
101 1.0, 3.0, 300.0);
102 propagator.setSlaveMode();
103
104
105 final RangeIonosphericDelayModifier modifier = new RangeIonosphericDelayModifier(model, frequency);
106
107 for (final ObservedMeasurement<?> measurement : measurements) {
108 final AbsoluteDate date = measurement.getDate();
109
110 final SpacecraftState refstate = propagator.propagate(date);
111
112 Range range = (Range) measurement;
113 EstimatedMeasurement<Range> evalNoMod = range.estimate(12, 17, new SpacecraftState[] { refstate });
114 Assert.assertEquals(12, evalNoMod.getIteration());
115 Assert.assertEquals(17, evalNoMod.getCount());
116
117
118 range.addModifier(modifier);
119 boolean found = false;
120 for (final EstimationModifier<Range> existing : range.getModifiers()) {
121 found = found || existing == modifier;
122 }
123 Assert.assertTrue(found);
124
125 EstimatedMeasurement<Range> eval = range.estimate(0, 0, new SpacecraftState[] { refstate });
126 Assert.assertEquals(evalNoMod.getStatus(), eval.getStatus());
127 eval.setStatus(EstimatedMeasurement.Status.REJECTED);
128 Assert.assertEquals(EstimatedMeasurement.Status.REJECTED, eval.getStatus());
129 eval.setStatus(evalNoMod.getStatus());
130
131 try {
132 eval.getParameterDerivatives(new ParameterDriver("extra", 0, 1, -1, +1));
133 Assert.fail("an exception should have been thrown");
134 } catch (OrekitIllegalArgumentException oiae) {
135 Assert.assertEquals(OrekitMessages.UNSUPPORTED_PARAMETER_NAME, oiae.getSpecifier());
136 }
137
138 final double diffMeters = eval.getEstimatedValue()[0] - evalNoMod.getEstimatedValue()[0];
139
140 Assert.assertEquals(0.0, diffMeters, 30.0);
141
142 }
143 }
144
145 @Test
146 public void testPhaseIonoModifier() {
147
148 Context context = EstimationTestUtils.eccentricContext("regular-data:potential:tides");
149
150 final NumericalPropagatorBuilder propagatorBuilder =
151 context.createBuilder(OrbitType.KEPLERIAN, PositionAngle.TRUE, true,
152 1.0e-6, 60.0, 0.001);
153
154
155 for (final GroundStation station : context.stations) {
156 station.getClockOffsetDriver().setSelected(true);
157 station.getEastOffsetDriver().setSelected(true);
158 station.getNorthOffsetDriver().setSelected(true);
159 station.getZenithOffsetDriver().setSelected(true);
160 }
161 final Propagator propagator = EstimationTestUtils.createPropagator(context.initialOrbit,
162 propagatorBuilder);
163 final double groundClockOffset = 12.0e-6;
164 for (final GroundStation station : context.stations) {
165 station.getClockOffsetDriver().setValue(groundClockOffset);
166 }
167 final double satClockOffset = 345.0e-6;
168 final List<ObservedMeasurement<?>> measurements =
169 EstimationTestUtils.createMeasurements(propagator,
170 new PhaseMeasurementCreator(context, Frequency.G01, 0,
171 satClockOffset),
172 1.0, 3.0, 300.0);
173 propagator.setSlaveMode();
174
175
176 final PhaseIonosphericDelayModifier modifier = new PhaseIonosphericDelayModifier(model, frequency);
177
178 for (final ObservedMeasurement<?> measurement : measurements) {
179 final AbsoluteDate date = measurement.getDate();
180
181 final SpacecraftState refstate = propagator.propagate(date);
182
183 Phase phase = (Phase) measurement;
184 EstimatedMeasurement<Phase> evalNoMod = phase.estimate(12, 17, new SpacecraftState[] { refstate });
185 Assert.assertEquals(12, evalNoMod.getIteration());
186 Assert.assertEquals(17, evalNoMod.getCount());
187
188
189 phase.addModifier(modifier);
190 boolean found = false;
191 for (final EstimationModifier<Phase> existing : phase.getModifiers()) {
192 found = found || existing == modifier;
193 }
194 Assert.assertTrue(found);
195
196 EstimatedMeasurement<Phase> eval = phase.estimate(0, 0, new SpacecraftState[] { refstate });
197 Assert.assertEquals(evalNoMod.getStatus(), eval.getStatus());
198 eval.setStatus(EstimatedMeasurement.Status.REJECTED);
199 Assert.assertEquals(EstimatedMeasurement.Status.REJECTED, eval.getStatus());
200 eval.setStatus(evalNoMod.getStatus());
201
202 try {
203 eval.getParameterDerivatives(new ParameterDriver("extra", 0, 1, -1, +1));
204 Assert.fail("an exception should have been thrown");
205 } catch (OrekitIllegalArgumentException oiae) {
206 Assert.assertEquals(OrekitMessages.UNSUPPORTED_PARAMETER_NAME, oiae.getSpecifier());
207 }
208
209 final double diffMeters = (eval.getEstimatedValue()[0] - evalNoMod.getEstimatedValue()[0]) * phase.getWavelength();
210 Assert.assertTrue(diffMeters < 0);
211 Assert.assertEquals(0.0, diffMeters, 30.0);
212
213 }
214 }
215
216 @Test
217 public void testPhaseEstimatedIonoModifier() {
218
219 Context context = EstimationTestUtils.eccentricContext("regular-data:potential:tides");
220
221 final NumericalPropagatorBuilder propagatorBuilder =
222 context.createBuilder(OrbitType.KEPLERIAN, PositionAngle.TRUE, true,
223 1.0e-6, 60.0, 0.001);
224
225
226 for (final GroundStation station : context.stations) {
227 station.getClockOffsetDriver().setSelected(true);
228 station.getEastOffsetDriver().setSelected(true);
229 station.getNorthOffsetDriver().setSelected(true);
230 station.getZenithOffsetDriver().setSelected(true);
231 }
232 final Propagator propagator = EstimationTestUtils.createPropagator(context.initialOrbit,
233 propagatorBuilder);
234 final double groundClockOffset = 12.0e-6;
235 for (final GroundStation station : context.stations) {
236 station.getClockOffsetDriver().setValue(groundClockOffset);
237 }
238 final double satClockOffset = 345.0e-6;
239 final List<ObservedMeasurement<?>> measurements =
240 EstimationTestUtils.createMeasurements(propagator,
241 new PhaseMeasurementCreator(context, Frequency.G01, 0,
242 satClockOffset),
243 1.0, 3.0, 300.0);
244 propagator.setSlaveMode();
245
246
247 final IonosphericModel mockModel = new MockIonosphericModel(12.0);
248 mockModel.getParametersDrivers().get(0).setSelected(true);
249 final PhaseIonosphericDelayModifier modifier = new PhaseIonosphericDelayModifier(mockModel, frequency);
250
251 for (final ObservedMeasurement<?> measurement : measurements) {
252 final AbsoluteDate date = measurement.getDate();
253
254 final SpacecraftState refstate = propagator.propagate(date);
255
256 Phase phase = (Phase) measurement;
257 EstimatedMeasurement<Phase> evalNoMod = phase.estimate(12, 17, new SpacecraftState[] { refstate });
258 Assert.assertEquals(12, evalNoMod.getIteration());
259 Assert.assertEquals(17, evalNoMod.getCount());
260
261
262 phase.addModifier(modifier);
263 boolean found = false;
264 for (final EstimationModifier<Phase> existing : phase.getModifiers()) {
265 found = found || existing == modifier;
266 }
267 Assert.assertTrue(found);
268
269 EstimatedMeasurement<Phase> eval = phase.estimate(0, 0, new SpacecraftState[] { refstate });
270 Assert.assertEquals(evalNoMod.getStatus(), eval.getStatus());
271 eval.setStatus(EstimatedMeasurement.Status.REJECTED);
272 Assert.assertEquals(EstimatedMeasurement.Status.REJECTED, eval.getStatus());
273 eval.setStatus(evalNoMod.getStatus());
274
275 try {
276 eval.getParameterDerivatives(new ParameterDriver("extra", 0, 1, -1, +1));
277 Assert.fail("an exception should have been thrown");
278 } catch (OrekitIllegalArgumentException oiae) {
279 Assert.assertEquals(OrekitMessages.UNSUPPORTED_PARAMETER_NAME, oiae.getSpecifier());
280 }
281
282 final double diffMeters = (eval.getEstimatedValue()[0] - evalNoMod.getEstimatedValue()[0]) * phase.getWavelength();
283 Assert.assertEquals(-12.0, diffMeters, 0.1);
284
285 }
286 }
287
288 @Test
289 public void testTurnAroundRangeIonoModifier() {
290
291 Context context = EstimationTestUtils.eccentricContext("regular-data:potential:tides");
292
293 final NumericalPropagatorBuilder propagatorBuilder =
294 context.createBuilder(OrbitType.KEPLERIAN, PositionAngle.TRUE, true,
295 1.0e-6, 60.0, 0.001);
296
297
298 for (Map.Entry<GroundStation, GroundStation> entry : context.TARstations.entrySet()) {
299 final GroundStation masterStation = entry.getKey();
300 final GroundStation slaveStation = entry.getValue();
301 masterStation.getClockOffsetDriver().setSelected(true);
302 masterStation.getEastOffsetDriver().setSelected(true);
303 masterStation.getNorthOffsetDriver().setSelected(true);
304 masterStation.getZenithOffsetDriver().setSelected(true);
305 slaveStation.getClockOffsetDriver().setSelected(false);
306 slaveStation.getEastOffsetDriver().setSelected(true);
307 slaveStation.getNorthOffsetDriver().setSelected(true);
308 slaveStation.getZenithOffsetDriver().setSelected(true);
309 }
310 final Propagator propagator = EstimationTestUtils.createPropagator(context.initialOrbit,
311 propagatorBuilder);
312 final List<ObservedMeasurement<?>> measurements =
313 EstimationTestUtils.createMeasurements(propagator,
314 new TurnAroundRangeMeasurementCreator(context),
315 1.0, 3.0, 300.0);
316 propagator.setSlaveMode();
317
318
319 final TurnAroundRangeIonosphericDelayModifier modifier = new TurnAroundRangeIonosphericDelayModifier(model, frequency);
320
321 for (final ObservedMeasurement<?> measurement : measurements) {
322 final AbsoluteDate date = measurement.getDate();
323
324 final SpacecraftState refstate = propagator.propagate(date);
325
326 TurnAroundRange turnAroundRange = (TurnAroundRange) measurement;
327 EstimatedMeasurement<TurnAroundRange> evalNoMod = turnAroundRange.estimate(12, 17, new SpacecraftState[] { refstate });
328 Assert.assertEquals(12, evalNoMod.getIteration());
329 Assert.assertEquals(17, evalNoMod.getCount());
330
331
332 turnAroundRange.addModifier(modifier);
333 boolean found = false;
334 for (final EstimationModifier<TurnAroundRange> existing : turnAroundRange.getModifiers()) {
335 found = found || existing == modifier;
336 }
337 Assert.assertTrue(found);
338
339 EstimatedMeasurement<TurnAroundRange> eval = turnAroundRange.estimate(12, 17, new SpacecraftState[] { refstate });
340 Assert.assertEquals(evalNoMod.getStatus(), eval.getStatus());
341 eval.setStatus(EstimatedMeasurement.Status.REJECTED);
342 Assert.assertEquals(EstimatedMeasurement.Status.REJECTED, eval.getStatus());
343 eval.setStatus(evalNoMod.getStatus());
344
345 try {
346 eval.getParameterDerivatives(new ParameterDriver("extra", 0, 1, -1, +1));
347 Assert.fail("an exception should have been thrown");
348 } catch (OrekitIllegalArgumentException oiae) {
349 Assert.assertEquals(OrekitMessages.UNSUPPORTED_PARAMETER_NAME, oiae.getSpecifier());
350 }
351
352 final double diffMeters = eval.getEstimatedValue()[0] - evalNoMod.getEstimatedValue()[0];
353
354 Assert.assertEquals(0.0, diffMeters, 30.0);
355
356 }
357 }
358
359 @Test
360 public void testRangeRateIonoModifier() {
361
362 Context context = EstimationTestUtils.eccentricContext("regular-data:potential:tides");
363
364 final NumericalPropagatorBuilder propagatorBuilder =
365 context.createBuilder(OrbitType.KEPLERIAN, PositionAngle.TRUE, true,
366 1.0e-6, 60.0, 0.001);
367
368
369 for (final GroundStation station : context.stations) {
370 station.getClockOffsetDriver().setSelected(true);
371 station.getEastOffsetDriver().setSelected(true);
372 station.getNorthOffsetDriver().setSelected(true);
373 station.getZenithOffsetDriver().setSelected(true);
374 }
375 final Propagator propagator = EstimationTestUtils.createPropagator(context.initialOrbit,
376 propagatorBuilder);
377 final double satClkDrift = 3.2e-10;
378 final List<ObservedMeasurement<?>> measurements =
379 EstimationTestUtils.createMeasurements(propagator,
380 new RangeRateMeasurementCreator(context, false, satClkDrift),
381 1.0, 3.0, 300.0);
382 propagator.setSlaveMode();
383
384 final RangeRateIonosphericDelayModifier modifier = new RangeRateIonosphericDelayModifier(model, frequency, true);
385
386 for (final ObservedMeasurement<?> measurement : measurements) {
387 final AbsoluteDate date = measurement.getDate();
388
389 final SpacecraftState refstate = propagator.propagate(date);
390
391 RangeRate rangeRate = (RangeRate) measurement;
392 EstimatedMeasurement<RangeRate> evalNoMod = rangeRate.estimate(0, 0, new SpacecraftState[] { refstate });
393
394
395 rangeRate.addModifier(modifier);
396
397
398 EstimatedMeasurement<RangeRate> eval = rangeRate.estimate(0, 0, new SpacecraftState[] { refstate });
399
400 final double diffMetersSec = eval.getEstimatedValue()[0] - evalNoMod.getEstimatedValue()[0];
401
402 Assert.assertEquals(0.0, diffMetersSec, 0.015);
403
404 }
405 }
406
407 @Test
408 public void testAngularIonoModifier() {
409
410 Context context = EstimationTestUtils.eccentricContext("regular-data:potential:tides");
411
412 final NumericalPropagatorBuilder propagatorBuilder =
413 context.createBuilder(OrbitType.KEPLERIAN, PositionAngle.TRUE, true,
414 1.0e-6, 60.0, 0.001);
415
416
417 for (final GroundStation station : context.stations) {
418 station.getClockOffsetDriver().setSelected(true);
419 station.getEastOffsetDriver().setSelected(true);
420 station.getNorthOffsetDriver().setSelected(true);
421 station.getZenithOffsetDriver().setSelected(true);
422 }
423 final Propagator propagator = EstimationTestUtils.createPropagator(context.initialOrbit,
424 propagatorBuilder);
425 final List<ObservedMeasurement<?>> measurements =
426 EstimationTestUtils.createMeasurements(propagator,
427 new AngularAzElMeasurementCreator(context),
428 1.0, 3.0, 300.0);
429 propagator.setSlaveMode();
430
431
432 final AngularIonosphericDelayModifier modifier = new AngularIonosphericDelayModifier(model, frequency);
433
434 for (final ObservedMeasurement<?> measurement : measurements) {
435 final AbsoluteDate date = measurement.getDate();
436
437 final SpacecraftState refstate = propagator.propagate(date);
438
439 AngularAzEl angular = (AngularAzEl) measurement;
440 EstimatedMeasurement<AngularAzEl> evalNoMod = angular.estimate(0, 0, new SpacecraftState[] { refstate });
441
442
443 angular.addModifier(modifier);
444
445 EstimatedMeasurement<AngularAzEl> eval = angular.estimate(0, 0, new SpacecraftState[] { refstate });
446
447 final double diffAz = MathUtils.normalizeAngle(eval.getEstimatedValue()[0], evalNoMod.getEstimatedValue()[0]) - evalNoMod.getEstimatedValue()[0];
448 final double diffEl = MathUtils.normalizeAngle(eval.getEstimatedValue()[1], evalNoMod.getEstimatedValue()[1]) - evalNoMod.getEstimatedValue()[1];
449
450 Assert.assertEquals(0.0, diffAz, 5.0e-5);
451 Assert.assertEquals(0.0, diffEl, 5.0e-6);
452 }
453 }
454
455 @Test
456 public void testKlobucharIonoModel() {
457 Context context = EstimationTestUtils.eccentricContext("regular-data:potential:tides");
458
459 final NumericalPropagatorBuilder propagatorBuilder =
460 context.createBuilder(OrbitType.KEPLERIAN, PositionAngle.TRUE, true,
461 1.0e-6, 60.0, 0.001);
462
463
464 for (final GroundStation station : context.stations) {
465 station.getClockOffsetDriver().setSelected(true);
466 station.getEastOffsetDriver().setSelected(true);
467 station.getNorthOffsetDriver().setSelected(true);
468 station.getZenithOffsetDriver().setSelected(true);
469 }
470 final Propagator propagator = EstimationTestUtils.createPropagator(context.initialOrbit,
471 propagatorBuilder);
472 final List<ObservedMeasurement<?>> measurements =
473 EstimationTestUtils.createMeasurements(propagator,
474 new RangeMeasurementCreator(context),
475 1.0, 3.0, 300.0);
476 propagator.setSlaveMode();
477
478 for (final ObservedMeasurement<?> measurement : measurements) {
479
480 final GroundStation station = ((Range) measurement).getStation();
481 final AbsoluteDate date = ((Range) measurement).getDate();
482 final SpacecraftState state = propagator.propagate(date);
483
484 double delayMeters = model.pathDelay(state, station.getBaseFrame(), frequency, model.getParameters());
485
486 final double epsilon = 1e-6;
487 Assert.assertTrue(Precision.compareTo(delayMeters, 15., epsilon) < 0);
488 Assert.assertTrue(Precision.compareTo(delayMeters, 0., epsilon) > 0);
489 }
490
491 }
492
493 private class MockIonosphericModel implements IonosphericModel {
494
495
496 private static final long serialVersionUID = 5944637011744634693L;
497
498
499 private final ParameterDriver ionoDelay;
500
501
502
503
504 public MockIonosphericModel(final double delay) {
505 ionoDelay = new ParameterDriver("ionospheric delay",
506 delay, FastMath.scalb(1.0, 0), 0.0, Double.POSITIVE_INFINITY);
507 }
508
509 @Override
510 public double pathDelay(final SpacecraftState state, final TopocentricFrame baseFrame,
511 final double frequency, double[] parameters) {
512 return parameters[0];
513 }
514
515 @Override
516 public <T extends RealFieldElement<T>> T pathDelay(final FieldSpacecraftState<T> state, final TopocentricFrame baseFrame,
517 final double frequency, final T[] parameters) {
518 return parameters[0];
519 }
520
521 @Override
522 public List<ParameterDriver> getParametersDrivers() {
523 return Collections.singletonList(ionoDelay);
524 }
525
526 }
527 }
528
529