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