1   /* Copyright 2002-2025 CS GROUP
2    * Licensed to CS GROUP (CS) under one or more
3    * contributor license agreements.  See the NOTICE file distributed with
4    * this work for additional information regarding copyright ownership.
5    * CS licenses this file to You under the Apache License, Version 2.0
6    * (the "License"); you may not use this file except in compliance with
7    * the License.  You may obtain a copy of the License at
8    *
9    *   http://www.apache.org/licenses/LICENSE-2.0
10   *
11   * Unless required by applicable law or agreed to in writing, software
12   * distributed under the License is distributed on an "AS IS" BASIS,
13   * WITHOUT WARRANTIES OR CONDITIONS OF ANY KIND, either express or implied.
14   * See the License for the specific language governing permissions and
15   * limitations under the License.
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.CalculusFieldElement;
24  import org.hipparchus.util.FastMath;
25  import org.hipparchus.util.MathUtils;
26  import org.hipparchus.util.Precision;
27  import org.junit.jupiter.api.Assertions;
28  import org.junit.jupiter.api.BeforeEach;
29  import org.junit.jupiter.api.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.BistaticRange;
37  import org.orekit.estimation.measurements.BistaticRangeMeasurementCreator;
38  import org.orekit.estimation.measurements.BistaticRangeRate;
39  import org.orekit.estimation.measurements.BistaticRangeRateMeasurementCreator;
40  import org.orekit.estimation.measurements.EstimatedMeasurement;
41  import org.orekit.estimation.measurements.EstimatedMeasurementBase;
42  import org.orekit.estimation.measurements.EstimationModifier;
43  import org.orekit.estimation.measurements.GroundStation;
44  import org.orekit.estimation.measurements.ObservedMeasurement;
45  import org.orekit.estimation.measurements.Range;
46  import org.orekit.estimation.measurements.RangeRate;
47  import org.orekit.estimation.measurements.RangeRateMeasurementCreator;
48  import org.orekit.estimation.measurements.TDOA;
49  import org.orekit.estimation.measurements.TDOAMeasurementCreator;
50  import org.orekit.estimation.measurements.TurnAroundRange;
51  import org.orekit.estimation.measurements.TurnAroundRangeMeasurementCreator;
52  import org.orekit.estimation.measurements.TwoWayRangeMeasurementCreator;
53  import org.orekit.estimation.measurements.gnss.Phase;
54  import org.orekit.estimation.measurements.gnss.PhaseMeasurementCreator;
55  import org.orekit.frames.TopocentricFrame;
56  import org.orekit.gnss.PredefinedGnssSignal;
57  import org.orekit.models.earth.ionosphere.IonosphericDelayModel;
58  import org.orekit.models.earth.ionosphere.IonosphericModel;
59  import org.orekit.models.earth.ionosphere.KlobucharIonoModel;
60  import org.orekit.orbits.OrbitType;
61  import org.orekit.orbits.PositionAngleType;
62  import org.orekit.propagation.FieldSpacecraftState;
63  import org.orekit.propagation.Propagator;
64  import org.orekit.propagation.SpacecraftState;
65  import org.orekit.propagation.conversion.NumericalPropagatorBuilder;
66  import org.orekit.time.AbsoluteDate;
67  import org.orekit.time.FieldAbsoluteDate;
68  import org.orekit.utils.ParameterDriver;
69  
70  public class IonoModifierTest {
71  
72      /** ionospheric model. */
73      private KlobucharIonoModel model;
74  
75      /** frequency [Hz]. */
76      private double frequency;
77  
78      @BeforeEach
79      public void setUp() throws Exception {
80          // Navigation message data
81          // .3820D-07   .1490D-07  -.1790D-06   .0000D-00          ION ALPHA
82          // .1430D+06   .0000D+00  -.3280D+06   .1130D+06          ION BETA
83          model = new KlobucharIonoModel(new double[]{.3820e-07, .1490e-07, -.1790e-06, 0},
84                                         new double[]{.1430e+06, 0, -.3280e+06, .1130e+06});
85          // GPS L1 in HZ
86          frequency = PredefinedGnssSignal.G01.getFrequency();
87      }
88  
89      @Test
90      public void testPhaseIonoModifier() {
91  
92          Context context = EstimationTestUtils.eccentricContext("regular-data:potential:tides");
93  
94          final NumericalPropagatorBuilder propagatorBuilder =
95                          context.createBuilder(OrbitType.KEPLERIAN, PositionAngleType.TRUE, true,
96                                                1.0e-6, 60.0, 0.001);
97  
98          // create perfect range measurements
99          for (final GroundStation station : context.stations) {
100             station.getClockOffsetDriver().setSelected(true);
101             station.getEastOffsetDriver().setSelected(true);
102             station.getNorthOffsetDriver().setSelected(true);
103             station.getZenithOffsetDriver().setSelected(true);
104         }
105         final Propagator propagator = EstimationTestUtils.createPropagator(context.initialOrbit,
106                                                                            propagatorBuilder);
107         final double groundClockOffset =  12.0e-6;
108         for (final GroundStation station : context.stations) {
109             station.getClockOffsetDriver().setValue(groundClockOffset);
110         }
111         final double satClockOffset    = 345.0e-6;
112         final List<ObservedMeasurement<?>> measurements =
113                         EstimationTestUtils.createMeasurements(propagator,
114                                                                new PhaseMeasurementCreator(context,
115                                                                                            PredefinedGnssSignal.G01, 0,
116                                                                                            satClockOffset),
117                                                                1.0, 3.0, 300.0);
118         propagator.clearStepHandlers();
119 
120 
121         final PhaseIonosphericDelayModifier modifier = new PhaseIonosphericDelayModifier(model, frequency);
122 
123         for (final ObservedMeasurement<?> measurement : measurements) {
124             final AbsoluteDate date = measurement.getDate();
125 
126             final SpacecraftState refstate = propagator.propagate(date);
127 
128             Phase phase = (Phase) measurement;
129             EstimatedMeasurementBase<Phase> evalNoMod = phase.estimateWithoutDerivatives(12, 17, new SpacecraftState[] { refstate });
130             Assertions.assertEquals(12, evalNoMod.getIteration());
131             Assertions.assertEquals(17, evalNoMod.getCount());
132 
133 
134             // add modifier
135             phase.addModifier(modifier);
136             boolean found = false;
137             for (final EstimationModifier<Phase> existing : phase.getModifiers()) {
138                 found = found || existing == modifier;
139             }
140             Assertions.assertTrue(found);
141             //
142             EstimatedMeasurement<Phase> eval = phase.estimate(0, 0,  new SpacecraftState[] { refstate });
143             Assertions.assertEquals(evalNoMod.getStatus(), eval.getStatus());
144             eval.setStatus(EstimatedMeasurement.Status.REJECTED);
145             Assertions.assertEquals(EstimatedMeasurement.Status.REJECTED, eval.getStatus());
146             eval.setStatus(evalNoMod.getStatus());
147 
148             try {
149                 eval.getParameterDerivatives(new ParameterDriver("extra", 0, 1, -1, +1), new AbsoluteDate());
150                 Assertions.fail("an exception should have been thrown");
151             } catch (OrekitIllegalArgumentException oiae) {
152                 Assertions.assertEquals(OrekitMessages.UNSUPPORTED_PARAMETER_NAME, oiae.getSpecifier());
153             }
154 
155             final double diffMeters = (eval.getEstimatedValue()[0] - evalNoMod.getEstimatedValue()[0]) * phase.getWavelength();
156             Assertions.assertTrue(diffMeters < 0);
157             Assertions.assertEquals(0.0, diffMeters, 30.0);
158 
159             Assertions.assertEquals(1,
160                                     eval.getAppliedEffects().entrySet().stream().
161                                     filter(e -> e.getKey().getEffectName().equals("ionosphere")).count());
162         }
163     }
164 
165     @Test
166     public void testPhaseEstimatedIonoModifier() {
167 
168         Context context = EstimationTestUtils.eccentricContext("regular-data:potential:tides");
169 
170         final NumericalPropagatorBuilder propagatorBuilder =
171                         context.createBuilder(OrbitType.KEPLERIAN, PositionAngleType.TRUE, true,
172                                               1.0e-6, 60.0, 0.001);
173 
174         // create perfect range measurements
175         for (final GroundStation station : context.stations) {
176             station.getClockOffsetDriver().setSelected(true);
177             station.getEastOffsetDriver().setSelected(true);
178             station.getNorthOffsetDriver().setSelected(true);
179             station.getZenithOffsetDriver().setSelected(true);
180         }
181         final Propagator propagator = EstimationTestUtils.createPropagator(context.initialOrbit,
182                                                                            propagatorBuilder);
183         final double groundClockOffset =  12.0e-6;
184         for (final GroundStation station : context.stations) {
185             station.getClockOffsetDriver().setValue(groundClockOffset);
186         }
187         final double satClockOffset    = 345.0e-6;
188         final List<ObservedMeasurement<?>> measurements =
189                         EstimationTestUtils.createMeasurements(propagator,
190                                                                new PhaseMeasurementCreator(context,
191                                                                                            PredefinedGnssSignal.G01, 0,
192                                                                                            satClockOffset),
193                                                                1.0, 3.0, 300.0);
194         propagator.clearStepHandlers();
195 
196 
197         final IonosphericModel mockModel = new MockIonosphericModel(12.0);
198         mockModel.getParametersDrivers().get(0).setSelected(true);
199         final PhaseIonosphericDelayModifier modifier = new PhaseIonosphericDelayModifier(mockModel, frequency);
200 
201         for (final ObservedMeasurement<?> measurement : measurements) {
202             final AbsoluteDate date = measurement.getDate();
203 
204             final SpacecraftState refstate = propagator.propagate(date);
205 
206             Phase phase = (Phase) measurement;
207             EstimatedMeasurementBase<Phase> evalNoMod = phase.estimateWithoutDerivatives(12, 17,new SpacecraftState[] { refstate });
208             Assertions.assertEquals(12, evalNoMod.getIteration());
209             Assertions.assertEquals(17, evalNoMod.getCount());
210 
211             // add modifier
212             phase.addModifier(modifier);
213             boolean found = false;
214             for (final EstimationModifier<Phase> existing : phase.getModifiers()) {
215                 found = found || existing == modifier;
216             }
217             Assertions.assertTrue(found);
218             //
219             EstimatedMeasurement<Phase> eval = phase.estimate(0, 0,  new SpacecraftState[] { refstate });
220             Assertions.assertEquals(evalNoMod.getStatus(), eval.getStatus());
221             eval.setStatus(EstimatedMeasurement.Status.REJECTED);
222             Assertions.assertEquals(EstimatedMeasurement.Status.REJECTED, eval.getStatus());
223             eval.setStatus(evalNoMod.getStatus());
224 
225             try {
226                 eval.getParameterDerivatives(new ParameterDriver("extra", 0, 1, -1, +1), new AbsoluteDate());
227                 Assertions.fail("an exception should have been thrown");
228             } catch (OrekitIllegalArgumentException oiae) {
229                 Assertions.assertEquals(OrekitMessages.UNSUPPORTED_PARAMETER_NAME, oiae.getSpecifier());
230             }
231 
232             final double diffMeters = (eval.getEstimatedValue()[0] - evalNoMod.getEstimatedValue()[0]) * phase.getWavelength();
233             Assertions.assertEquals(-12.0, diffMeters, 0.1);
234 
235         }
236     }
237 
238     @Test
239     public void testRangeIonoModifier() {
240 
241         Context context = EstimationTestUtils.eccentricContext("regular-data:potential:tides");
242 
243         final NumericalPropagatorBuilder propagatorBuilder =
244                         context.createBuilder(OrbitType.KEPLERIAN, PositionAngleType.TRUE, true,
245                                               1.0e-6, 60.0, 0.001);
246 
247         // create perfect range measurements
248         for (final GroundStation station : context.stations) {
249             station.getClockOffsetDriver().setSelected(true);
250             station.getEastOffsetDriver().setSelected(true);
251             station.getNorthOffsetDriver().setSelected(true);
252             station.getZenithOffsetDriver().setSelected(true);
253         }
254         final Propagator propagator = EstimationTestUtils.createPropagator(context.initialOrbit,
255                                                                            propagatorBuilder);
256         final List<ObservedMeasurement<?>> measurements =
257                         EstimationTestUtils.createMeasurements(propagator,
258                                                                new TwoWayRangeMeasurementCreator(context),
259                                                                1.0, 3.0, 300.0);
260         propagator.clearStepHandlers();
261 
262 
263         final RangeIonosphericDelayModifier modifier = new RangeIonosphericDelayModifier(model, frequency);
264 
265         for (final ObservedMeasurement<?> measurement : measurements) {
266             final AbsoluteDate date = measurement.getDate();
267 
268             final SpacecraftState refstate = propagator.propagate(date);
269 
270             Range range = (Range) measurement;
271             EstimatedMeasurementBase<Range> evalNoMod = range.estimateWithoutDerivatives(12, 17, new SpacecraftState[] { refstate });
272             Assertions.assertEquals(12, evalNoMod.getIteration());
273             Assertions.assertEquals(17, evalNoMod.getCount());
274 
275             // add modifier
276             range.addModifier(modifier);
277             boolean found = false;
278             for (final EstimationModifier<Range> existing : range.getModifiers()) {
279                 found = found || existing == modifier;
280             }
281             Assertions.assertTrue(found);
282             //
283             EstimatedMeasurement<Range> eval = range.estimate(0, 0,  new SpacecraftState[] { refstate });
284             Assertions.assertEquals(evalNoMod.getStatus(), eval.getStatus());
285             eval.setStatus(EstimatedMeasurement.Status.REJECTED);
286             Assertions.assertEquals(EstimatedMeasurement.Status.REJECTED, eval.getStatus());
287             eval.setStatus(evalNoMod.getStatus());
288 
289             try {
290                 eval.getParameterDerivatives(new ParameterDriver("extra", 0, 1, -1, +1), new AbsoluteDate());
291                 Assertions.fail("an exception should have been thrown");
292             } catch (OrekitIllegalArgumentException oiae) {
293                 Assertions.assertEquals(OrekitMessages.UNSUPPORTED_PARAMETER_NAME, oiae.getSpecifier());
294             }
295 
296             final double diffMeters = eval.getEstimatedValue()[0] - evalNoMod.getEstimatedValue()[0];
297             // TODO: check threshold
298             Assertions.assertEquals(0.0, diffMeters, 30.0);
299 
300         }
301     }
302 
303     @Test
304     public void testRangeRateIonoModifier() {
305 
306         Context context = EstimationTestUtils.eccentricContext("regular-data:potential:tides");
307 
308         final NumericalPropagatorBuilder propagatorBuilder =
309                         context.createBuilder(OrbitType.KEPLERIAN, PositionAngleType.TRUE, true,
310                                               1.0e-6, 60.0, 0.001);
311 
312         // create perfect range measurements
313         for (final GroundStation station : context.stations) {
314             station.getClockOffsetDriver().setSelected(true);
315             station.getEastOffsetDriver().setSelected(true);
316             station.getNorthOffsetDriver().setSelected(true);
317             station.getZenithOffsetDriver().setSelected(true);
318         }
319         final Propagator propagator = EstimationTestUtils.createPropagator(context.initialOrbit,
320                                                                            propagatorBuilder);
321         final double satClkDrift = 3.2e-10;
322         final List<ObservedMeasurement<?>> measurements =
323                         EstimationTestUtils.createMeasurements(propagator,
324                                                                new RangeRateMeasurementCreator(context, false, satClkDrift),
325                                                                1.0, 3.0, 300.0);
326         propagator.clearStepHandlers();
327 
328         final RangeRateIonosphericDelayModifier modifier = new RangeRateIonosphericDelayModifier(model, frequency, true);
329 
330         for (final ObservedMeasurement<?> measurement : measurements) {
331             final AbsoluteDate date = measurement.getDate();
332 
333             final SpacecraftState refstate = propagator.propagate(date);
334 
335             RangeRate rangeRate = (RangeRate) measurement;
336             EstimatedMeasurementBase<RangeRate> evalNoMod = rangeRate.estimateWithoutDerivatives(new SpacecraftState[] { refstate });
337 
338             // add modifier
339             rangeRate.addModifier(modifier);
340 
341             //
342             EstimatedMeasurementBase<RangeRate> eval = rangeRate.estimateWithoutDerivatives(new SpacecraftState[] { refstate });
343 
344             final double diffMetersSec = eval.getEstimatedValue()[0] - evalNoMod.getEstimatedValue()[0];
345             // TODO: check threshold
346             Assertions.assertEquals(0.0, diffMetersSec, 0.015);
347 
348         }
349     }
350 
351     @Test
352     public void testTurnAroundRangeIonoModifier() {
353 
354         Context context = EstimationTestUtils.eccentricContext("regular-data:potential:tides");
355 
356         final NumericalPropagatorBuilder propagatorBuilder =
357                         context.createBuilder(OrbitType.KEPLERIAN, PositionAngleType.TRUE, true,
358                                               1.0e-6, 60.0, 0.001);
359 
360         // Create perfect turn-around measurements
361         for (Map.Entry<GroundStation, GroundStation> entry : context.TARstations.entrySet()) {
362             final GroundStation    primaryStation = entry.getKey();
363             final GroundStation    secondaryStation  = entry.getValue();
364             primaryStation.getClockOffsetDriver().setSelected(true);
365             primaryStation.getEastOffsetDriver().setSelected(true);
366             primaryStation.getNorthOffsetDriver().setSelected(true);
367             primaryStation.getZenithOffsetDriver().setSelected(true);
368             secondaryStation.getClockOffsetDriver().setSelected(false);
369             secondaryStation.getEastOffsetDriver().setSelected(true);
370             secondaryStation.getNorthOffsetDriver().setSelected(true);
371             secondaryStation.getZenithOffsetDriver().setSelected(true);
372         }
373         final Propagator propagator = EstimationTestUtils.createPropagator(context.initialOrbit,
374                                                                            propagatorBuilder);
375         final List<ObservedMeasurement<?>> measurements =
376                         EstimationTestUtils.createMeasurements(propagator,
377                                                                new TurnAroundRangeMeasurementCreator(context),
378                                                                1.0, 3.0, 300.0);
379         propagator.clearStepHandlers();
380 
381 
382         final TurnAroundRangeIonosphericDelayModifier modifier = new TurnAroundRangeIonosphericDelayModifier(model, frequency);
383 
384         for (final ObservedMeasurement<?> measurement : measurements) {
385             final AbsoluteDate date = measurement.getDate();
386 
387             final SpacecraftState refstate = propagator.propagate(date);
388 
389             TurnAroundRange turnAroundRange = (TurnAroundRange) measurement;
390             EstimatedMeasurementBase<TurnAroundRange> evalNoMod = turnAroundRange.estimateWithoutDerivatives(12, 17, new SpacecraftState[] { refstate });
391             Assertions.assertEquals(12, evalNoMod.getIteration());
392             Assertions.assertEquals(17, evalNoMod.getCount());
393 
394             // Add modifier
395             turnAroundRange.addModifier(modifier);
396             boolean found = false;
397             for (final EstimationModifier<TurnAroundRange> existing : turnAroundRange.getModifiers()) {
398                 found = found || existing == modifier;
399             }
400             Assertions.assertTrue(found);
401             //
402             EstimatedMeasurement<TurnAroundRange> eval = turnAroundRange.estimate(12, 17, new SpacecraftState[] { refstate });
403             Assertions.assertEquals(evalNoMod.getStatus(), eval.getStatus());
404             eval.setStatus(EstimatedMeasurement.Status.REJECTED);
405             Assertions.assertEquals(EstimatedMeasurement.Status.REJECTED, eval.getStatus());
406             eval.setStatus(evalNoMod.getStatus());
407 
408             try {
409                 eval.getParameterDerivatives(new ParameterDriver("extra", 0, 1, -1, +1));
410                 Assertions.fail("an exception should have been thrown");
411             } catch (OrekitIllegalArgumentException oiae) {
412                 Assertions.assertEquals(OrekitMessages.UNSUPPORTED_PARAMETER_NAME, oiae.getSpecifier());
413             }
414 
415             final double diffMeters = eval.getEstimatedValue()[0] - evalNoMod.getEstimatedValue()[0];
416             // TODO: check threshold
417             Assertions.assertEquals(0.0, diffMeters, 30.0);
418 
419         }
420     }
421 
422     @Test
423     public void testBistaticRangeIonoModifier() {
424 
425         Context context = EstimationTestUtils.eccentricContext("regular-data:potential:tides");
426 
427         final NumericalPropagatorBuilder propagatorBuilder =
428                         context.createBuilder(OrbitType.KEPLERIAN, PositionAngleType.TRUE, true,
429                                               1.0e-6, 60.0, 0.001);
430         final Propagator propagator = EstimationTestUtils.createPropagator(context.initialOrbit,
431                                                                            propagatorBuilder);
432         // create perfect range measurements
433         final GroundStation emitter = context.BRRstations.getKey();
434         emitter.getClockOffsetDriver().setSelected(true);
435         emitter.getEastOffsetDriver().setSelected(true);
436         emitter.getNorthOffsetDriver().setSelected(true);
437         emitter.getZenithOffsetDriver().setSelected(true);
438         final GroundStation receiver = context.BRRstations.getValue();
439         receiver.getClockOffsetDriver().setSelected(true);
440         receiver.getEastOffsetDriver().setSelected(true);
441         receiver.getNorthOffsetDriver().setSelected(true);
442         receiver.getZenithOffsetDriver().setSelected(true);
443         final List<ObservedMeasurement<?>> measurements =
444                         EstimationTestUtils.createMeasurements(propagator,
445                                                                new BistaticRangeMeasurementCreator(context),
446                                                                1.0, 3.0, 300.0);
447         propagator.clearStepHandlers();
448 
449         final BistaticRangeIonosphericDelayModifier modifier =
450                         new BistaticRangeIonosphericDelayModifier(model, frequency);
451 
452         for (final ObservedMeasurement<?> measurement : measurements) {
453             BistaticRange biRange = (BistaticRange) measurement;
454             final SpacecraftState refstate = propagator.propagate(biRange.getDate());
455 
456             // Estimate without modifier
457             EstimatedMeasurementBase<BistaticRange> evalNoMod = biRange.estimateWithoutDerivatives(new SpacecraftState[] { refstate });
458 
459             // add modifier
460             biRange.addModifier(modifier);
461 
462             // Estimate with modifier
463             EstimatedMeasurementBase<BistaticRange> eval = biRange.estimateWithoutDerivatives(new SpacecraftState[] { refstate });
464 
465             final double diffMeters = eval.getEstimatedValue()[0] - evalNoMod.getEstimatedValue()[0];
466             // TODO: check threshold
467             Assertions.assertTrue(diffMeters < 12.0);
468             Assertions.assertTrue(diffMeters >  4.0);
469         }
470     }
471 
472     @Test
473     public void testBistaticRangeRateIonoModifier() {
474 
475         Context context = EstimationTestUtils.eccentricContext("regular-data:potential:tides");
476 
477         final NumericalPropagatorBuilder propagatorBuilder =
478                         context.createBuilder(OrbitType.KEPLERIAN, PositionAngleType.TRUE, true,
479                                               1.0e-6, 60.0, 0.001);
480         final Propagator propagator = EstimationTestUtils.createPropagator(context.initialOrbit,
481                                                                            propagatorBuilder);
482         // create perfect range rate measurements
483         final GroundStation emitter = context.BRRstations.getKey();
484         emitter.getEastOffsetDriver().setSelected(true);
485         emitter.getNorthOffsetDriver().setSelected(true);
486         emitter.getZenithOffsetDriver().setSelected(true);
487         final GroundStation receiver = context.BRRstations.getValue();
488         receiver.getClockOffsetDriver().setSelected(true);
489         receiver.getEastOffsetDriver().setSelected(true);
490         receiver.getNorthOffsetDriver().setSelected(true);
491         receiver.getZenithOffsetDriver().setSelected(true);
492         final List<ObservedMeasurement<?>> measurements =
493                         EstimationTestUtils.createMeasurements(propagator,
494                                                                new BistaticRangeRateMeasurementCreator(context),
495                                                                1.0, 3.0, 300.0);
496         propagator.clearStepHandlers();
497 
498         final BistaticRangeRateIonosphericDelayModifier modifier =
499                         new BistaticRangeRateIonosphericDelayModifier(model, frequency);
500 
501         for (final ObservedMeasurement<?> measurement : measurements) {
502             BistaticRangeRate biRangeRate = (BistaticRangeRate) measurement;
503             final SpacecraftState refstate = propagator.propagate(biRangeRate.getDate());
504 
505             // Estimate without modifier
506             EstimatedMeasurementBase<BistaticRangeRate> evalNoMod = biRangeRate.estimateWithoutDerivatives(new SpacecraftState[] { refstate });
507 
508             // add modifier
509             biRangeRate.addModifier(modifier);
510 
511             // Estimate with modifier
512             EstimatedMeasurementBase<BistaticRangeRate> eval = biRangeRate.estimateWithoutDerivatives(new SpacecraftState[] { refstate });
513 
514             final double diffMetersSec = eval.getEstimatedValue()[0] - evalNoMod.getEstimatedValue()[0];
515             // TODO: check threshold
516             final double epsilon = 1e-5;
517             Assertions.assertTrue(Precision.compareTo(diffMetersSec,  0.002, epsilon) < 0);
518             Assertions.assertTrue(Precision.compareTo(diffMetersSec, -0.010, epsilon) > 0);
519 
520         }
521     }
522 
523     @Test
524     public void testTDOAIonoModifier() {
525 
526         Context context = EstimationTestUtils.eccentricContext("regular-data:potential:tides");
527 
528         final NumericalPropagatorBuilder propagatorBuilder =
529                         context.createBuilder(OrbitType.KEPLERIAN, PositionAngleType.TRUE, true,
530                                               1.0e-6, 60.0, 0.001);
531         final Propagator propagator = EstimationTestUtils.createPropagator(context.initialOrbit,
532                                                                            propagatorBuilder);
533         // create perfect range measurements
534         final GroundStation emitter = context.TDOAstations.getKey();
535         emitter.getClockOffsetDriver().setSelected(true);
536         emitter.getEastOffsetDriver().setSelected(true);
537         emitter.getNorthOffsetDriver().setSelected(true);
538         emitter.getZenithOffsetDriver().setSelected(true);
539         final GroundStation receiver = context.TDOAstations.getValue();
540         receiver.getClockOffsetDriver().setSelected(true);
541         receiver.getEastOffsetDriver().setSelected(true);
542         receiver.getNorthOffsetDriver().setSelected(true);
543         receiver.getZenithOffsetDriver().setSelected(true);
544         final List<ObservedMeasurement<?>> measurements =
545                         EstimationTestUtils.createMeasurements(propagator,
546                                                                new TDOAMeasurementCreator(context),
547                                                                1.0, 3.0, 300.0);
548         propagator.clearStepHandlers();
549 
550         final TDOAIonosphericDelayModifier modifier =
551                         new TDOAIonosphericDelayModifier(model, frequency);
552 
553         for (final ObservedMeasurement<?> measurement : measurements) {
554             TDOA tdoa = (TDOA) measurement;
555             final SpacecraftState refState = propagator.propagate(tdoa.getDate());
556 
557             // Estimate without modifier
558             EstimatedMeasurementBase<TDOA> evalNoMod = tdoa.estimateWithoutDerivatives(new SpacecraftState[] { refState });
559 
560             // add modifier
561             tdoa.addModifier(modifier);
562 
563             // Estimate with modifier
564             EstimatedMeasurement<TDOA> eval = tdoa.estimate(0, 0, new SpacecraftState[] { refState });
565 
566             final double diffSec = eval.getEstimatedValue()[0] - evalNoMod.getEstimatedValue()[0];
567 
568             // TODO: check threshold
569             final double epsilon = 1.e-11;
570             Assertions.assertTrue(Precision.compareTo(diffSec, 2.90e-9, epsilon) < 0);
571             Assertions.assertTrue(Precision.compareTo(diffSec, 0.85e-9, epsilon) > 0);
572         }
573     }
574 
575     @Test
576     public void testAngularIonoModifier() {
577 
578         Context context = EstimationTestUtils.eccentricContext("regular-data:potential:tides");
579 
580         final NumericalPropagatorBuilder propagatorBuilder =
581                         context.createBuilder(OrbitType.KEPLERIAN, PositionAngleType.TRUE, true,
582                                               1.0e-6, 60.0, 0.001);
583 
584         // create perfect range measurements
585         for (final GroundStation station : context.stations) {
586             station.getClockOffsetDriver().setSelected(true);
587             station.getEastOffsetDriver().setSelected(true);
588             station.getNorthOffsetDriver().setSelected(true);
589             station.getZenithOffsetDriver().setSelected(true);
590         }
591         final Propagator propagator = EstimationTestUtils.createPropagator(context.initialOrbit,
592                                                                            propagatorBuilder);
593         final List<ObservedMeasurement<?>> measurements =
594                         EstimationTestUtils.createMeasurements(propagator,
595                                                                new AngularAzElMeasurementCreator(context),
596                                                                1.0, 3.0, 300.0);
597         propagator.clearStepHandlers();
598 
599 
600         final AngularIonosphericDelayModifier modifier = new AngularIonosphericDelayModifier(model, frequency);
601 
602         for (final ObservedMeasurement<?> measurement : measurements) {
603             final AbsoluteDate date = measurement.getDate();
604 
605             final SpacecraftState refstate = propagator.propagate(date);
606 
607             AngularAzEl angular = (AngularAzEl) measurement;
608             EstimatedMeasurementBase<AngularAzEl> evalNoMod = angular.estimateWithoutDerivatives(new SpacecraftState[] { refstate });
609 
610             // add modifier
611             angular.addModifier(modifier);
612             //
613             EstimatedMeasurementBase<AngularAzEl> eval = angular.estimateWithoutDerivatives(new SpacecraftState[] { refstate });
614 
615             final double diffAz = MathUtils.normalizeAngle(eval.getEstimatedValue()[0], evalNoMod.getEstimatedValue()[0]) - evalNoMod.getEstimatedValue()[0];
616             final double diffEl = MathUtils.normalizeAngle(eval.getEstimatedValue()[1], evalNoMod.getEstimatedValue()[1]) - evalNoMod.getEstimatedValue()[1];
617             // TODO: check threshold
618             Assertions.assertEquals(0.0, diffAz, 5.0e-5);
619             Assertions.assertEquals(0.0, diffEl, 5.0e-6);
620         }
621     }
622 
623     @Test
624     public void testKlobucharIonoModel() {
625         Context context = EstimationTestUtils.eccentricContext("regular-data:potential:tides");
626 
627         final NumericalPropagatorBuilder propagatorBuilder =
628                         context.createBuilder(OrbitType.KEPLERIAN, PositionAngleType.TRUE, true,
629                                               1.0e-6, 60.0, 0.001);
630 
631         // create perfect range measurements
632         for (final GroundStation station : context.stations) {
633             station.getClockOffsetDriver().setSelected(true);
634             station.getEastOffsetDriver().setSelected(true);
635             station.getNorthOffsetDriver().setSelected(true);
636             station.getZenithOffsetDriver().setSelected(true);
637         }
638         final Propagator propagator = EstimationTestUtils.createPropagator(context.initialOrbit,
639                                                                            propagatorBuilder);
640         final List<ObservedMeasurement<?>> measurements =
641                         EstimationTestUtils.createMeasurements(propagator,
642                                                                new TwoWayRangeMeasurementCreator(context),
643                                                                1.0, 3.0, 300.0);
644         propagator.clearStepHandlers();
645 
646         for (final ObservedMeasurement<?> measurement : measurements) {
647             // parameter corresponding to station position offset
648             final GroundStation   station = ((Range) measurement).getStation();
649             final AbsoluteDate    date    = measurement.getDate();
650             final SpacecraftState state   = propagator.propagate(date);
651 
652             double delayMeters = model.pathDelay(state, station.getBaseFrame(), frequency, model.getParameters());
653 
654             final double epsilon = 1e-6;
655             Assertions.assertTrue(Precision.compareTo(delayMeters, 15., epsilon) < 0);
656             Assertions.assertTrue(Precision.compareTo(delayMeters, 0., epsilon) > 0);
657         }
658 
659     }
660 
661     private static class MockIonosphericModel implements IonosphericModel, IonosphericDelayModel {
662 
663         /** Driver for the ionospheric delay.*/
664         private final ParameterDriver ionoDelay;
665 
666         /** Constructor.
667          * @param delay initial ionospheric delay
668          */
669         public MockIonosphericModel(final double delay) {
670             ionoDelay = new ParameterDriver("ionospheric delay",
671                                             delay, FastMath.scalb(1.0, 0), 0.0, Double.POSITIVE_INFINITY);
672         }
673 
674         @Override
675         public double pathDelay(final SpacecraftState state, final TopocentricFrame baseFrame,
676                                 final double frequency, double[] parameters) {
677             return parameters[0];
678         }
679 
680         @Override
681         public double pathDelay(final SpacecraftState state,
682                                 final TopocentricFrame baseFrame, final AbsoluteDate receptionDate,
683                                 final double frequency, double[] parameters) {
684             return parameters[0];
685         }
686 
687         @Override
688         public <T extends CalculusFieldElement<T>> T pathDelay(final FieldSpacecraftState<T> state, final TopocentricFrame baseFrame,
689                                                            final double frequency, final  T[] parameters) {
690             return parameters[0];
691         }
692 
693         @Override
694         public <T extends CalculusFieldElement<T>> T pathDelay(final FieldSpacecraftState<T> state,
695                                                                final TopocentricFrame baseFrame,
696                                                                final FieldAbsoluteDate<T> receptionDate,
697                                                                final double frequency, final  T[] parameters) {
698             return parameters[0];
699         }
700 
701         @Override
702         public List<ParameterDriver> getParametersDrivers() {
703             return Collections.singletonList(ionoDelay);
704         }
705 
706     }
707 }
708 
709