1   /* Copyright 2002-2026 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;
18  
19  import java.util.HashMap;
20  import java.util.List;
21  import java.util.Locale;
22  import java.util.Map;
23  
24  import org.hipparchus.analysis.UnivariateVectorFunction;
25  import org.hipparchus.analysis.differentiation.FiniteDifferencesDifferentiator;
26  import org.hipparchus.analysis.differentiation.Gradient;
27  import org.hipparchus.analysis.differentiation.GradientField;
28  import org.hipparchus.analysis.differentiation.UnivariateDifferentiableVectorFunction;
29  import org.hipparchus.geometry.euclidean.threed.Rotation;
30  import org.hipparchus.geometry.euclidean.threed.Vector3D;
31  import org.hipparchus.linear.RealMatrix;
32  import org.hipparchus.optim.nonlinear.vector.leastsquares.LevenbergMarquardtOptimizer;
33  import org.hipparchus.random.RandomGenerator;
34  import org.hipparchus.random.Well19937a;
35  import org.hipparchus.util.FastMath;
36  import org.hipparchus.util.Precision;
37  import org.junit.jupiter.api.Assertions;
38  import org.junit.jupiter.api.Test;
39  import org.orekit.Utils;
40  import org.orekit.bodies.BodyShape;
41  import org.orekit.bodies.GeodeticPoint;
42  import org.orekit.bodies.OneAxisEllipsoid;
43  import org.orekit.errors.OrekitException;
44  import org.orekit.errors.OrekitMessages;
45  import org.orekit.estimation.Context;
46  import org.orekit.estimation.EstimationTestUtils;
47  import org.orekit.estimation.leastsquares.BatchLSEstimator;
48  import org.orekit.frames.FieldTransform;
49  import org.orekit.frames.Frame;
50  import org.orekit.frames.FramesFactory;
51  import org.orekit.frames.StaticTransform;
52  import org.orekit.frames.TopocentricFrame;
53  import org.orekit.frames.Transform;
54  import org.orekit.orbits.OrbitType;
55  import org.orekit.orbits.PositionAngleType;
56  import org.orekit.propagation.Propagator;
57  import org.orekit.propagation.conversion.NumericalPropagatorBuilder;
58  import org.orekit.time.AbsoluteDate;
59  import org.orekit.time.FieldAbsoluteDate;
60  import org.orekit.time.clocks.QuadraticClockModel;
61  import org.orekit.utils.Constants;
62  import org.orekit.utils.FieldPVCoordinates;
63  import org.orekit.utils.FieldPVCoordinatesProvider;
64  import org.orekit.utils.IERSConventions;
65  import org.orekit.utils.PVCoordinates;
66  import org.orekit.utils.PVCoordinatesProvider;
67  import org.orekit.utils.ParameterDriver;
68  
69  class GroundStationTest {
70  
71      @Test
72      void getPVCoordinatesProvider() {
73          // GIVEN
74          EstimationTestUtils.eccentricContext("regular-data:potential:tides");
75          final Frame ecef = FramesFactory.getGTOD(true);
76          final GeodeticPoint geodeticPoint = new GeodeticPoint(1., 2., 3.);
77          final TopocentricFrame topocentricFrame = new TopocentricFrame(new OneAxisEllipsoid(Constants.WGS84_EARTH_EQUATORIAL_RADIUS,
78                  Constants.WGS84_EARTH_FLATTENING, ecef), geodeticPoint, "");
79          final GroundStation station = new GroundStation(topocentricFrame);
80          final AbsoluteDate date = AbsoluteDate.ARBITRARY_EPOCH;
81          for (final ParameterDriver driver: selectAllDrivers(station)) {
82              driver.setReferenceDate(date);
83          }
84          station.getPolarOffsetYDriver().setValue(0.1);
85          station.getEastOffsetDriver().setValue(1);
86          station.getNorthOffsetDriver().setValue(-1);
87          station.getZenithOffsetDriver().setValue(2);
88          final PVCoordinatesProvider pvCoordinatesProvider = station.getPVCoordinatesProvider();
89          // WHEN
90          final PVCoordinates actualPV = pvCoordinatesProvider.getPVCoordinates(date, ecef);
91          // THEN
92          final Vector3D expectedPosition = station.getOffsetToInertial(ecef, date, true).transformPosition(Vector3D.ZERO);
93          Assertions.assertArrayEquals(expectedPosition.toArray(), actualPV.getPosition().toArray(), 1e-8);
94          Assertions.assertEquals(0., actualPV.getVelocity().getNorm2());
95      }
96  
97      @Test
98      void getPVCoordinatesProviderGetPosition() {
99          // GIVEN
100         EstimationTestUtils.eccentricContext("regular-data:potential:tides");
101         final Frame ecef = FramesFactory.getGTOD(true);
102         final GeodeticPoint geodeticPoint = new GeodeticPoint(1., 2., 3.);
103         final TopocentricFrame topocentricFrame = new TopocentricFrame(new OneAxisEllipsoid(Constants.WGS84_EARTH_EQUATORIAL_RADIUS,
104                 Constants.WGS84_EARTH_FLATTENING, ecef), geodeticPoint, "");
105         final GroundStation station = new GroundStation(topocentricFrame);
106         final AbsoluteDate date = AbsoluteDate.ARBITRARY_EPOCH;
107         for (final ParameterDriver driver: selectAllDrivers(station)) {
108             driver.setReferenceDate(date);
109         }
110         station.getClockBiasDriver().setValue(0.1);
111         station.getPolarOffsetYDriver().setValue(2);
112         station.getEastOffsetDriver().setValue(1);
113         station.getNorthOffsetDriver().setValue(-1);
114         station.getZenithOffsetDriver().setValue(2);
115         final PVCoordinatesProvider pvCoordinatesProvider = station.getPVCoordinatesProvider();
116         final Frame eci = FramesFactory.getEME2000();
117         // WHEN
118         final Vector3D actualPosition = pvCoordinatesProvider.getPosition(date, eci);
119         // THEN
120         final PVCoordinates pvCoordinates = pvCoordinatesProvider.getPVCoordinates(date, eci);
121         Assertions.assertEquals(pvCoordinates.getPosition(), actualPosition);
122     }
123 
124     @Test
125     void getPVCoordinatesProviderGetVelocity() {
126         // GIVEN
127         EstimationTestUtils.eccentricContext("regular-data:potential:tides");
128         final Frame ecef = FramesFactory.getGTOD(true);
129         final GeodeticPoint geodeticPoint = new GeodeticPoint(1., 2., 3.);
130         final TopocentricFrame topocentricFrame = new TopocentricFrame(new OneAxisEllipsoid(Constants.WGS84_EARTH_EQUATORIAL_RADIUS,
131                 Constants.WGS84_EARTH_FLATTENING, ecef), geodeticPoint, "");
132         final GroundStation station = new GroundStation(topocentricFrame);
133         final AbsoluteDate date = AbsoluteDate.ARBITRARY_EPOCH;
134         for (final ParameterDriver driver: selectAllDrivers(station)) {
135             driver.setReferenceDate(date);
136         }
137         station.getClockBiasDriver().setValue(0.1);
138         station.getPolarOffsetYDriver().setValue(2);
139         station.getEastOffsetDriver().setValue(1);
140         station.getNorthOffsetDriver().setValue(-1);
141         station.getZenithOffsetDriver().setValue(2);
142         final PVCoordinatesProvider pvCoordinatesProvider = station.getPVCoordinatesProvider();
143         final Frame eci = FramesFactory.getEME2000();
144         // WHEN
145         final Vector3D actualVelocity = pvCoordinatesProvider.getVelocity(date, eci);
146         // THEN
147         final PVCoordinates pvCoordinates = pvCoordinatesProvider.getPVCoordinates(date, eci);
148         Assertions.assertEquals(pvCoordinates.getVelocity(), actualVelocity);
149     }
150 
151     @Test
152     void getFieldPVCoordinatesProvider() {
153         // GIVEN
154         EstimationTestUtils.eccentricContext("regular-data:potential:tides");
155         final Frame ecef = FramesFactory.getGTOD(true);
156         final GeodeticPoint geodeticPoint = new GeodeticPoint(1., 2., 3.);
157         final TopocentricFrame topocentricFrame = new TopocentricFrame(new OneAxisEllipsoid(Constants.WGS84_EARTH_EQUATORIAL_RADIUS,
158                 Constants.WGS84_EARTH_FLATTENING, ecef), geodeticPoint, "");
159         final GroundStation station = new GroundStation(topocentricFrame);
160         final AbsoluteDate date = AbsoluteDate.ARBITRARY_EPOCH;
161         for (final ParameterDriver driver: selectAllDrivers(station)) {
162             driver.setReferenceDate(date);
163         }
164         station.getClockBiasDriver().setValue(-0.1);
165         station.getPolarOffsetXDriver().setValue(3);
166         station.getEastOffsetDriver().setValue(1);
167         station.getNorthOffsetDriver().setValue(-1);
168         station.getZenithOffsetDriver().setValue(2);
169         final int freeParams = 0;
170         final FieldPVCoordinatesProvider<Gradient> pvCoordinatesProvider = station.getFieldPVCoordinatesProvider(freeParams,
171                 new HashMap<>());
172         final Frame eci = FramesFactory.getEME2000();
173         final FieldAbsoluteDate<Gradient> fieldDate = new FieldAbsoluteDate<>(GradientField.getField(freeParams), date);
174         // WHEN
175         final FieldPVCoordinates<Gradient> pvCoordinates = pvCoordinatesProvider.getPVCoordinates(fieldDate, eci);
176         // THEN
177         Assertions.assertEquals(pvCoordinates.getPosition(), pvCoordinatesProvider.getPosition(fieldDate, eci));
178         final PVCoordinates nonFieldPV = station.getPVCoordinatesProvider().getPVCoordinates(date, eci);
179         Assertions.assertArrayEquals(pvCoordinates.getVelocity().toVector3D().toArray(), nonFieldPV.getVelocity().toArray(), 1e-12);
180     }
181 
182     @Test
183     void testEstimateClockOffset() {
184 
185         Context context = EstimationTestUtils.eccentricContext("regular-data:potential:tides");
186 
187         final NumericalPropagatorBuilder propagatorBuilder =
188                         context.createNumerical(OrbitType.KEPLERIAN, PositionAngleType.TRUE, true,
189                                               1.0e-6, 60.0, 0.001);
190 
191         // create perfect range measurements
192         final Propagator propagator = EstimationTestUtils.createPropagator(context.initialOrbit,
193                                                                            propagatorBuilder);
194         final List<ObservedMeasurement<?>> measurements =
195                         EstimationTestUtils.createMeasurements(propagator,
196                                                                new TwoWayRangeMeasurementCreator(context),
197                                                                1.0, 3.0, 300.0);
198 
199         // change one station clock
200         final TopocentricFrame base  = context.stations.get(0).getBaseFrame();
201         final BodyShape parent       = base.getParentShape();
202         final double deltaClock      = 0.00084532;
203         final String changedSuffix   = "-changed";
204         final QuadraticClockModel blankClock = new QuadraticClockModel(context.initialOrbit.getDate(), 0.0, 0.0, 0.0);
205         final GroundStation changed  = new GroundStation(new TopocentricFrame(parent, base.getPoint(),
206                                                                               base.getName() + changedSuffix), context.ut1.getEOPHistory(),
207                                                          blankClock, context.stations.get(0).getDisplacements());
208 
209         // create orbit estimator
210         final BatchLSEstimator estimator = new BatchLSEstimator(new LevenbergMarquardtOptimizer(),
211                                                                 propagatorBuilder);
212         for (final ObservedMeasurement<?> measurement : measurements) {
213             final Range range = (Range) measurement;
214             final String name = range.getObserver().getName() + changedSuffix;
215                 if (changed.getBaseFrame().getName().equals(name)) {
216                     estimator.addMeasurement(new Range(changed, range.isTwoWay(),
217                                                        range.getDate().shiftedBy(deltaClock),
218                                                        range.getObservedValue()[0],
219                                                        range.getTheoreticalStandardDeviation()[0],
220                                                        range.getBaseWeight()[0],
221                                                        range.getSatellites().get(0)));
222                 } else {
223                     estimator.addMeasurement(range);
224                 }
225         }
226         estimator.setParametersConvergenceThreshold(1.0e-3);
227         estimator.setMaxIterations(100);
228         estimator.setMaxEvaluations(200);
229 
230         // we want to estimate station clock offset
231         changed.getClockBiasDriver().setSelected(true);
232         changed.getEastOffsetDriver().setSelected(false);
233         changed.getNorthOffsetDriver().setSelected(false);
234         changed.getZenithOffsetDriver().setSelected(false);
235 
236         EstimationTestUtils.checkFit(false, context, estimator, 2, 3,
237                                      0.0, 9.5e-7,
238                                      0.0, 2.4e-6,
239                                      0.0, 1.8e-7,
240                                      0.0, 8e-11);
241         Assertions.assertEquals(deltaClock, changed.getClockBiasDriver().getValue(), 9.6e-11);
242 
243         RealMatrix normalizedCovariances = estimator.getOptimum().getCovariances(1.0e-10);
244         RealMatrix physicalCovariances   = estimator.getPhysicalCovariances(1.0e-10);
245         Assertions.assertFalse(changed.isSpaceBased());
246         Assertions.assertEquals(7,        normalizedCovariances.getRowDimension());
247         Assertions.assertEquals(7,        normalizedCovariances.getColumnDimension());
248         Assertions.assertEquals(7,        physicalCovariances.getRowDimension());
249         Assertions.assertEquals(7,        physicalCovariances.getColumnDimension());
250         Assertions.assertEquals(4.185e-9, physicalCovariances.getEntry(6, 6), 3.0e-13);
251 
252     }
253 
254     @Test
255     void testClockOffsetValues() {
256 
257         Context context = EstimationTestUtils.eccentricContext("regular-data:potential:tides");
258 
259         // change one station clock
260         final TopocentricFrame base  = context.stations.get(0).getBaseFrame();
261         final BodyShape parent       = base.getParentShape();
262         final String changedSuffix   = "-changed";
263         final QuadraticClockModel quadraticClock = new QuadraticClockModel(context.initialOrbit.getDate(), 3.0e-9, 2.0e-9, 1.0e-9);
264         final GroundStation changed  = new GroundStation(new TopocentricFrame(parent, base.getPoint(),
265                                                                               base.getName() + changedSuffix), context.ut1.getEOPHistory(),
266                                                          quadraticClock, context.stations.get(0).getDisplacements());
267 
268         final AbsoluteDate shiftedDate = context.initialOrbit.getDate().shiftedBy(60.0); // one minute
269         final double bias  = changed.getOffsetValue(shiftedDate);
270         final double drift = changed.getOffsetRate(shiftedDate);
271 
272         Assertions.assertEquals(3.723e-6, bias, 1e-10);
273         Assertions.assertEquals(1.22e-7, drift, 1e-10);
274 
275     }
276 
277     @Test
278     void testEstimateStationPosition() {
279 
280         Context context = EstimationTestUtils.eccentricContext("regular-data:potential:tides");
281 
282         final NumericalPropagatorBuilder propagatorBuilder =
283                         context.createNumerical(OrbitType.KEPLERIAN, PositionAngleType.TRUE, true,
284                                               1.0e-6, 60.0, 0.001);
285 
286         // create perfect range measurements
287         final Propagator propagator = EstimationTestUtils.createPropagator(context.initialOrbit,
288                                                                            propagatorBuilder);
289         final List<ObservedMeasurement<?>> measurements =
290                         EstimationTestUtils.createMeasurements(propagator,
291                                                                new TwoWayRangeMeasurementCreator(context),
292                                                                1.0, 3.0, 300.0);
293 
294         // move one station
295         final RandomGenerator random = new Well19937a(0x4adbecfc743bda60L);
296         final TopocentricFrame base = context.stations.get(0).getBaseFrame();
297         final BodyShape parent = base.getParentShape();
298         final Vector3D baseOrigin = parent.transform(base.getPoint());
299         final Vector3D deltaTopo = new Vector3D(2 * random.nextDouble() - 1,
300                                                 2 * random.nextDouble() - 1,
301                                                 2 * random.nextDouble() - 1);
302         final StaticTransform topoToParent = base.getStaticTransformTo(parent.getBodyFrame(), (AbsoluteDate) null);
303         final Vector3D deltaParent   = topoToParent.transformVector(deltaTopo);
304         final String movedSuffix     = "-moved";
305         final GroundStation moved = new GroundStation(new TopocentricFrame(parent,
306                                                                            parent.transform(baseOrigin.subtract(deltaParent),
307                                                                                             parent.getBodyFrame(),
308                                                                                             null),
309                                                                            base.getName() + movedSuffix), context.ut1.getEOPHistory(),
310                                                       context.stations.get(0).getDisplacements());
311 
312         // create orbit estimator
313         final BatchLSEstimator estimator = new BatchLSEstimator(new LevenbergMarquardtOptimizer(),
314                                                                 propagatorBuilder);
315         for (final ObservedMeasurement<?> measurement : measurements) {
316             final Range range = (Range) measurement;
317             final String name = range.getObserver().getName() + movedSuffix;
318                 if (moved.getBaseFrame().getName().equals(name)) {
319                     estimator.addMeasurement(new Range(moved, range.isTwoWay(), range.getDate(),
320                                                        range.getObservedValue()[0],
321                                                        range.getTheoreticalStandardDeviation()[0],
322                                                        range.getBaseWeight()[0],
323                                                        range.getSatellites().get(0)));
324                 } else {
325                     estimator.addMeasurement(range);
326                 }
327         }
328         estimator.setParametersConvergenceThreshold(1.0e-3);
329         estimator.setMaxIterations(100);
330         estimator.setMaxEvaluations(200);
331 
332         // we want to estimate station offsets
333         moved.getClockBiasDriver().setSelected(false);
334         moved.getEastOffsetDriver().setSelected(true);
335         moved.getNorthOffsetDriver().setSelected(true);
336         moved.getZenithOffsetDriver().setSelected(true);
337 
338         EstimationTestUtils.checkFit(false, context, estimator, 2, 3,
339                                      0.0, 6.7e-7,
340                                      0.0, 1.8e-6,
341                                      0.0, 1.3e-6,
342                                      0.0, 6.1e-10);
343         Assertions.assertEquals(deltaTopo.getX(), moved.getEastOffsetDriver().getValue(),   9.4e-7);
344         Assertions.assertEquals(deltaTopo.getY(), moved.getNorthOffsetDriver().getValue(),  6.2e-7);
345         Assertions.assertEquals(deltaTopo.getZ(), moved.getZenithOffsetDriver().getValue(), 2.6e-7);
346 
347         GeodeticPoint result = moved.getOffsetGeodeticPoint((AbsoluteDate) null);
348 
349         GeodeticPoint reference = context.stations.get(0).getBaseFrame().getPoint();
350         Assertions.assertEquals(reference.getLatitude(),  result.getLatitude(),  1e-12);
351         Assertions.assertEquals(reference.getLongitude(), result.getLongitude(), 1e-12);
352         Assertions.assertEquals(reference.getAltitude(),  result.getAltitude(),  2.7e-7);
353 
354         RealMatrix normalizedCovariances = estimator.getOptimum().getCovariances(1.0e-10);
355         RealMatrix physicalCovariances   = estimator.getPhysicalCovariances(1.0e-10);
356         Assertions.assertEquals(9,       normalizedCovariances.getRowDimension());
357         Assertions.assertEquals(9,       normalizedCovariances.getColumnDimension());
358         Assertions.assertEquals(9,       physicalCovariances.getRowDimension());
359         Assertions.assertEquals(9,       physicalCovariances.getColumnDimension());
360         Assertions.assertEquals(0.55431, physicalCovariances.getEntry(6, 6), 1.0e-5);
361         Assertions.assertEquals(0.22694, physicalCovariances.getEntry(7, 7), 1.0e-5);
362         Assertions.assertEquals(0.13106, physicalCovariances.getEntry(8, 8), 1.0e-5);
363     }
364 
365     @Test
366     void testEstimateEOP() {
367 
368         Context linearEOPContext = EstimationTestUtils.eccentricContext("linear-EOP:regular-data/de431-ephemerides:potential:tides");
369 
370         final AbsoluteDate refDate = new AbsoluteDate(2000, 2, 24, linearEOPContext.utc);
371         final double dut10 = 0.3079738;
372         final double lod   = 0.0011000;
373         final double xp0   = 68450.0e-6;
374         final double xpDot =   -50.0e-6;
375         final double yp0   =    60.0e-6;
376         final double ypDot =     2.0e-6;
377         for (double dt = -2 * Constants.JULIAN_DAY; dt < 2 * Constants.JULIAN_DAY; dt += 300.0) {
378             AbsoluteDate date = refDate.shiftedBy(dt);
379             Assertions.assertEquals(dut10 - dt * lod / Constants.JULIAN_DAY,
380                                 linearEOPContext.ut1.getEOPHistory().getUT1MinusUTC(date),
381                                 1.0e-15);
382             Assertions.assertEquals(lod,
383                                 linearEOPContext.ut1.getEOPHistory().getLOD(date),
384                                 1.0e-15);
385             Assertions.assertEquals((xp0 + xpDot * dt / Constants.JULIAN_DAY) * Constants.ARC_SECONDS_TO_RADIANS,
386                                 linearEOPContext.ut1.getEOPHistory().getPoleCorrection(date).getXp(),
387                                 1.0e-15);
388             Assertions.assertEquals((yp0 + ypDot * dt / Constants.JULIAN_DAY) * Constants.ARC_SECONDS_TO_RADIANS,
389                                 linearEOPContext.ut1.getEOPHistory().getPoleCorrection(date).getYp(),
390                                 1.0e-15);
391         }
392         final NumericalPropagatorBuilder linearPropagatorBuilder =
393                         linearEOPContext.createNumerical(OrbitType.KEPLERIAN, PositionAngleType.TRUE, true,
394                                               1.0e-6, 60.0, 0.001);
395 
396         // create perfect range measurements
397         final Propagator propagator = EstimationTestUtils.createPropagator(linearEOPContext.initialOrbit,
398                                                                            linearPropagatorBuilder);
399         final List<ObservedMeasurement<?>> linearMeasurements =
400                         EstimationTestUtils.createMeasurements(propagator,
401                                                                new TwoWayRangeMeasurementCreator(linearEOPContext),
402                                                                1.0, 5.0, 60.0);
403 
404         Utils.clearFactories();
405         Context zeroEOPContext = EstimationTestUtils.eccentricContext("zero-EOP:regular-data/de431-ephemerides:potential:potential:tides");
406         for (double dt = -2 * Constants.JULIAN_DAY; dt < 2 * Constants.JULIAN_DAY; dt += 300.0) {
407             AbsoluteDate date = refDate.shiftedBy(dt);
408             Assertions.assertEquals(0.0,
409                                 zeroEOPContext.ut1.getEOPHistory().getUT1MinusUTC(date),
410                                 1.0e-15);
411             Assertions.assertEquals(0.0,
412                                 zeroEOPContext.ut1.getEOPHistory().getLOD(date),
413                                 1.0e-15);
414             Assertions.assertEquals(0.0,
415                                 zeroEOPContext.ut1.getEOPHistory().getPoleCorrection(date).getXp(),
416                                 1.0e-15);
417             Assertions.assertEquals(0.0,
418                                 zeroEOPContext.ut1.getEOPHistory().getPoleCorrection(date).getYp(),
419                                 1.0e-15);
420         }
421 
422         // create orbit estimator
423         final NumericalPropagatorBuilder zeroPropagatorBuilder =
424                         linearEOPContext.createNumerical(OrbitType.KEPLERIAN, PositionAngleType.TRUE, true,
425                                               1.0e-6, 60.0, 0.001);
426         final BatchLSEstimator estimator = new BatchLSEstimator(new LevenbergMarquardtOptimizer(),
427                                                                 zeroPropagatorBuilder);
428         for (final ObservedMeasurement<?> linearMeasurement : linearMeasurements) {
429             Range linearRange = (Range) linearMeasurement;
430             for (final GroundStation station : zeroEOPContext.stations) {
431                 if (station.getBaseFrame().getName().equals(linearRange.getObserver().getName())) {
432                     Range zeroRange = new Range(station, linearRange.isTwoWay(),
433                                                 linearRange.getDate(),
434                                                 linearRange.getObservedValue()[0],
435                                                 linearRange.getTheoreticalStandardDeviation()[0],
436                                                 linearRange.getBaseWeight()[0],
437                                                 linearRange.getSatellites().get(0));
438                     estimator.addMeasurement(zeroRange);
439                 }
440             }
441         }
442         estimator.setParametersConvergenceThreshold(1.0e-3);
443         estimator.setMaxIterations(100);
444         estimator.setMaxEvaluations(200);
445 
446         // we want to estimate pole and prime meridian
447         GroundStation station = zeroEOPContext.stations.get(0);
448         station.getPrimeMeridianOffsetDriver().setReferenceDate(refDate);
449         station.getPrimeMeridianOffsetDriver().setSelected(true);
450         station.getPrimeMeridianDriftDriver().setSelected(true);
451         station.getPolarOffsetXDriver().setReferenceDate(refDate);
452         station.getPolarOffsetXDriver().setSelected(true);
453         station.getPolarDriftXDriver().setSelected(true);
454         station.getPolarOffsetYDriver().setReferenceDate(refDate);
455         station.getPolarOffsetYDriver().setSelected(true);
456         station.getPolarDriftYDriver().setSelected(true);
457 
458         // just for the fun and to speed up test, we will use orbit determination, *without* estimating orbit
459         for (final ParameterDriver driver : zeroPropagatorBuilder.getOrbitalParametersDrivers().getDrivers()) {
460             driver.setSelected(false);
461         }
462 
463         estimator.estimate();
464 
465         final double computedDut1  = station.getPrimeMeridianOffsetDriver().getValue() / EstimatedEarthFrameProvider.EARTH_ANGULAR_VELOCITY;
466         final double computedLOD   = station.getPrimeMeridianDriftDriver().getValue() * (-Constants.JULIAN_DAY / EstimatedEarthFrameProvider.EARTH_ANGULAR_VELOCITY);
467         final double computedXp    = station.getPolarOffsetXDriver().getValue() / Constants.ARC_SECONDS_TO_RADIANS;
468         final double computedXpDot = station.getPolarDriftXDriver().getValue()  / Constants.ARC_SECONDS_TO_RADIANS * Constants.JULIAN_DAY;
469         final double computedYp    = station.getPolarOffsetYDriver().getValue() / Constants.ARC_SECONDS_TO_RADIANS;
470         final double computedYpDot = station.getPolarDriftYDriver().getValue()  / Constants.ARC_SECONDS_TO_RADIANS * Constants.JULIAN_DAY;
471         Assertions.assertEquals(0.0, FastMath.abs(dut10 - computedDut1),  1.5e-9);
472         Assertions.assertEquals(0.0, FastMath.abs(lod - computedLOD),     9.6e-10);
473         Assertions.assertEquals(0.0, FastMath.abs(xp0 - computedXp),      1.3e-8);
474         Assertions.assertEquals(0.0, FastMath.abs(xpDot - computedXpDot), 8.1e-9);
475         Assertions.assertEquals(0.0, FastMath.abs(yp0 - computedYp),      5.6e-9);
476         Assertions.assertEquals(0.0, FastMath.abs(ypDot - computedYpDot), 5.5e-9);
477 
478         // thresholds to use if orbit is estimated
479         // (i.e. when commenting out the loop above that sets orbital parameters drivers to "not selected")
480 //         Assertions.assertEquals(dut10, computedDut1,  6.6e-3);
481 //         Assertions.assertEquals(lod,   computedLOD,   1.1e-9);
482 //         Assertions.assertEquals(xp0,   computedXp,    3.3e-8);
483 //         Assertions.assertEquals(xpDot, computedXpDot, 2.2e-8);
484 //         Assertions.assertEquals(yp0,   computedYp,    3.3e-8);
485 //         Assertions.assertEquals(ypDot, computedYpDot, 3.8e-8);
486 
487     }
488 
489     @Test
490     void testClockOffsetCartesianDerivativesOctantPxPyPz() {
491         double relativeTolerancePositionValue      =  2.3e-15;
492         double relativeTolerancePositionDerivative =  2.5e-10;
493         double relativeToleranceVelocityValue      =  3.0e-15;
494         double relativeToleranceVelocityDerivative =  1.7e-10;
495         doTestCartesianDerivatives(FastMath.toRadians(35), FastMath.toRadians(20), 1200.0, 100.0,
496                                    relativeTolerancePositionValue, relativeTolerancePositionDerivative,
497                                    relativeToleranceVelocityValue, relativeToleranceVelocityDerivative,
498                                    ".*-clock-.*");
499     }
500 
501     @Test
502     void testClockOffsetAngularDerivativesOctantPxPyPz() {
503         double toleranceRotationValue              =  1.9e-15;
504         double toleranceRotationDerivative         =  4.9e-15;
505         double toleranceRotationRateValue          =  1.1e-19;
506         double toleranceRotationRateDerivative     =  6.3e-19;
507         doTestAngularDerivatives(FastMath.toRadians(35), FastMath.toRadians(20), 1200.0, 100.0,
508                                  toleranceRotationValue,     toleranceRotationDerivative,
509                                  toleranceRotationRateValue, toleranceRotationRateDerivative,
510                                    ".*-clock-.*");
511     }
512 
513     @Test
514     void testClockOffsetCartesianDerivativesOctantPxPyMz() {
515         double relativeTolerancePositionValue      =  1.4e-15;
516         double relativeTolerancePositionDerivative =  1.7e-10;
517         double relativeToleranceVelocityValue      =  2.5e-15;
518         double relativeToleranceVelocityDerivative =  1.8e-10;
519         doTestCartesianDerivatives(FastMath.toRadians(-35), FastMath.toRadians(20), 1200.0, 100.0,
520                                    relativeTolerancePositionValue, relativeTolerancePositionDerivative,
521                                    relativeToleranceVelocityValue, relativeToleranceVelocityDerivative,
522                                    ".*-clock-.*");
523     }
524 
525     @Test
526     void testClockOffsetAngularDerivativesOctantPxPyMz() {
527         double toleranceRotationValue              =  1.7e-15;
528         double toleranceRotationDerivative         =  5.0e-15;
529         double toleranceRotationRateValue          =  1.1e-19;
530         double toleranceRotationRateDerivative     =  6.3e-19;
531         doTestAngularDerivatives(FastMath.toRadians(-35), FastMath.toRadians(20), 1200.0, 100.0,
532                                  toleranceRotationValue,     toleranceRotationDerivative,
533                                  toleranceRotationRateValue, toleranceRotationRateDerivative,
534                                    ".*-clock-.*");
535     }
536 
537     @Test
538     void testClockOffsetCartesianDerivativesOctantPxMyPz() {
539         double relativeTolerancePositionValue      =  1.7e-15;
540         double relativeTolerancePositionDerivative =  2.6e-10;
541         double relativeToleranceVelocityValue      =  2.8e-15;
542         double relativeToleranceVelocityDerivative =  1.8e-10;
543         doTestCartesianDerivatives(FastMath.toRadians(35), FastMath.toRadians(-20), 1200.0, 100.0,
544                                    relativeTolerancePositionValue, relativeTolerancePositionDerivative,
545                                    relativeToleranceVelocityValue, relativeToleranceVelocityDerivative,
546                                    ".*-clock-.*");
547     }
548 
549     @Test
550     void testClockOffsetAngularDerivativesOctantPxMyPz() {
551         double toleranceRotationValue              =  1.6e-15;
552         double toleranceRotationDerivative         =  4.9e-15;
553         double toleranceRotationRateValue          =  1.1e-19;
554         double toleranceRotationRateDerivative     =  6.3e-19;
555         doTestAngularDerivatives(FastMath.toRadians(35), FastMath.toRadians(-20), 1200.0, 100.0,
556                                  toleranceRotationValue,     toleranceRotationDerivative,
557                                  toleranceRotationRateValue, toleranceRotationRateDerivative,
558                                    ".*-clock-.*");
559     }
560 
561     @Test
562     void testClockOffsetCartesianDerivativesOctantPxMyMz() {
563         double relativeTolerancePositionValue      =  1.5e-15;
564         double relativeTolerancePositionDerivative =  1.6e-10;
565         double relativeToleranceVelocityValue      =  2.3e-15;
566         double relativeToleranceVelocityDerivative =  1.8e-10;
567         doTestCartesianDerivatives(FastMath.toRadians(-35), FastMath.toRadians(-20), 1200.0, 100.0,
568                                    relativeTolerancePositionValue, relativeTolerancePositionDerivative,
569                                    relativeToleranceVelocityValue, relativeToleranceVelocityDerivative,
570                                    ".*-clock-.*");
571     }
572 
573     @Test
574     void testClockOffsetAngularDerivativesOctantPxMyMz() {
575         double toleranceRotationValue              =  1.5e-15;
576         double toleranceRotationDerivative         =  5.0e-15;
577         double toleranceRotationRateValue          =  1.1e-19;
578         double toleranceRotationRateDerivative     =  6.3e-19;
579         doTestAngularDerivatives(FastMath.toRadians(-35), FastMath.toRadians(-20), 1200.0, 100.0,
580                                  toleranceRotationValue,     toleranceRotationDerivative,
581                                  toleranceRotationRateValue, toleranceRotationRateDerivative,
582                                    ".*-clock-.*");
583     }
584 
585     @Test
586     void testClockOffsetCartesianDerivativesOctantMxPyPz() {
587         double relativeTolerancePositionValue      =  1.9e-15;
588         double relativeTolerancePositionDerivative =  2.6e-10;
589         double relativeToleranceVelocityValue      =  2.6e-15;
590         double relativeToleranceVelocityDerivative =  2.0e-10;
591         doTestCartesianDerivatives(FastMath.toRadians(150), FastMath.toRadians(20), 1200.0, 100.0,
592                                    relativeTolerancePositionValue, relativeTolerancePositionDerivative,
593                                    relativeToleranceVelocityValue, relativeToleranceVelocityDerivative,
594                                    ".*-clock-.*");
595     }
596 
597     @Test
598     void testClockOffsetAngularDerivativesOctantMxPyPz() {
599         double toleranceRotationValue              =  1.4e-15;
600         double toleranceRotationDerivative         =  5.2e-15;
601         double toleranceRotationRateValue          =  1.1e-19;
602         double toleranceRotationRateDerivative     =  6.3e-19;
603         doTestAngularDerivatives(FastMath.toRadians(150), FastMath.toRadians(20), 1200.0, 100.0,
604                                  toleranceRotationValue,     toleranceRotationDerivative,
605                                  toleranceRotationRateValue, toleranceRotationRateDerivative,
606                                    ".*-clock-.*");
607     }
608 
609     @Test
610     void testClockOffsetCartesianDerivativesOctantMxPyMz() {
611         double relativeTolerancePositionValue      =  1.9e-15;
612         double relativeTolerancePositionDerivative =  2.6e-10;
613         double relativeToleranceVelocityValue      =  2.6e-15;
614         double relativeToleranceVelocityDerivative =  2.0e-10;
615         doTestCartesianDerivatives(FastMath.toRadians(150), FastMath.toRadians(20), 1200.0, 100.0,
616                                    relativeTolerancePositionValue, relativeTolerancePositionDerivative,
617                                    relativeToleranceVelocityValue, relativeToleranceVelocityDerivative,
618                                    ".*-clock-.*");
619     }
620 
621     @Test
622     void testClockOffsetAngularDerivativesOctantMxPyMz() {
623         double toleranceRotationValue              =  1.4e-15;
624         double toleranceRotationDerivative         =  5.2e-15;
625         double toleranceRotationRateValue          =  1.1e-19;
626         double toleranceRotationRateDerivative     =  6.3e-19;
627         doTestAngularDerivatives(FastMath.toRadians(150), FastMath.toRadians(20), 1200.0, 100.0,
628                                  toleranceRotationValue,     toleranceRotationDerivative,
629                                  toleranceRotationRateValue, toleranceRotationRateDerivative,
630                                    ".*-clock-.*");
631     }
632 
633     @Test
634     void testClockOffsetCartesianDerivativesOctantMxMyPz() {
635         double relativeTolerancePositionValue      =  1.5e-15;
636         double relativeTolerancePositionDerivative =  1.7e-10;
637         double relativeToleranceVelocityValue      =  2.9e-15;
638         double relativeToleranceVelocityDerivative =  1.9e-10;
639         doTestCartesianDerivatives(FastMath.toRadians(-150), FastMath.toRadians(-20), 1200.0, 100.0,
640                                    relativeTolerancePositionValue, relativeTolerancePositionDerivative,
641                                    relativeToleranceVelocityValue, relativeToleranceVelocityDerivative,
642                                    ".*-clock-.*");
643     }
644 
645     @Test
646     void testClockOffsetAngularDerivativesOctantMxMyPz() {
647         double toleranceRotationValue              =  1.6e-15;
648         double toleranceRotationDerivative         =  4.9e-15;
649         double toleranceRotationRateValue          =  1.1e-19;
650         double toleranceRotationRateDerivative     =  6.3e-19;
651         doTestAngularDerivatives(FastMath.toRadians(-150), FastMath.toRadians(-20), 1200.0, 100.0,
652                                  toleranceRotationValue,     toleranceRotationDerivative,
653                                  toleranceRotationRateValue, toleranceRotationRateDerivative,
654                                    ".*-clock-.*");
655     }
656 
657     @Test
658     void testClockOffsetCartesianDerivativesOctantMxMyMz() {
659         double relativeTolerancePositionValue      =  1.5e-15;
660         double relativeTolerancePositionDerivative =  1.7e-10;
661         double relativeToleranceVelocityValue      =  2.9e-15;
662         double relativeToleranceVelocityDerivative =  1.9e-10;
663         doTestCartesianDerivatives(FastMath.toRadians(-150), FastMath.toRadians(-20), 1200.0, 100.0,
664                                    relativeTolerancePositionValue, relativeTolerancePositionDerivative,
665                                    relativeToleranceVelocityValue, relativeToleranceVelocityDerivative,
666                                    ".*-clock-.*");
667     }
668 
669     @Test
670     void testClockOffsetAngularDerivativesOctantMxMyMz() {
671         double toleranceRotationValue              =  1.6e-15;
672         double toleranceRotationDerivative         =  4.9e-15;
673         double toleranceRotationRateValue          =  1.1e-19;
674         double toleranceRotationRateDerivative     =  6.3e-19;
675         doTestAngularDerivatives(FastMath.toRadians(-150), FastMath.toRadians(-20), 1200.0, 100.0,
676                                  toleranceRotationValue,     toleranceRotationDerivative,
677                                  toleranceRotationRateValue, toleranceRotationRateDerivative,
678                                    ".*-clock-.*");
679     }
680 
681     @Test
682     void testClockOffsetCartesianDerivativesNearPole() {
683         double relativeTolerancePositionValue      =  1.2e-15;
684         double relativeTolerancePositionDerivative =  1.6e-04;
685         double relativeToleranceVelocityValue      =  1.0e-13;
686         double relativeToleranceVelocityDerivative =  4.3e-09;
687         doTestCartesianDerivatives(FastMath.toRadians(89.99995), FastMath.toRadians(90), 1200.0, 100.0,
688                                    relativeTolerancePositionValue, relativeTolerancePositionDerivative,
689                                    relativeToleranceVelocityValue, relativeToleranceVelocityDerivative,
690                                    ".*-clock-.*");
691     }
692 
693     @Test
694     void testClockOffsetAngularDerivativesNearPole() {
695         double toleranceRotationValue              =  1.5e-15;
696         double toleranceRotationDerivative         =  5.7e-15;
697         double toleranceRotationRateValue          =  1.1e-19;
698         double toleranceRotationRateDerivative     =  6.3e-19;
699         doTestAngularDerivatives(FastMath.toRadians(89.99995), FastMath.toRadians(90), 1200.0, 100.0,
700                                  toleranceRotationValue,     toleranceRotationDerivative,
701                                  toleranceRotationRateValue, toleranceRotationRateDerivative,
702                                    ".*-clock-.*");
703     }
704 
705     @Test
706     void testStationOffsetCartesianDerivativesOctantPxPyPz() {
707         double relativeTolerancePositionValue      =  2.3e-15;
708         double relativeTolerancePositionDerivative =  1.1e-10;
709         double relativeToleranceVelocityValue      =  3.3e-15;
710         double relativeToleranceVelocityDerivative =  1.2e-10;
711         doTestCartesianDerivatives(FastMath.toRadians(35), FastMath.toRadians(20), 1200.0, 100.0,
712                                    relativeTolerancePositionValue, relativeTolerancePositionDerivative,
713                                    relativeToleranceVelocityValue, relativeToleranceVelocityDerivative,
714                                    ".*-offset-East", ".*-offset-North", ".*-offset-Zenith");
715     }
716 
717     @Test
718     void testStationOffsetAngularDerivativesOctantPxPyPz() {
719         double toleranceRotationValue              =  1.9e-15;
720         double toleranceRotationDerivative         =  3.1e-18;
721         double toleranceRotationRateValue          =  1.5e-19;
722         double toleranceRotationRateDerivative     =  Precision.SAFE_MIN;
723         doTestAngularDerivatives(FastMath.toRadians(35), FastMath.toRadians(20), 1200.0, 100.0,
724                                  toleranceRotationValue,     toleranceRotationDerivative,
725                                  toleranceRotationRateValue, toleranceRotationRateDerivative,
726                                  ".*-offset-East", ".*-offset-North", ".*-offset-Zenith");
727     }
728 
729     @Test
730     void testStationOffsetCartesianDerivativesOctantPxPyMz() {
731         double relativeTolerancePositionValue      =  1.4e-15;
732         double relativeTolerancePositionDerivative =  7.3e-11;
733         double relativeToleranceVelocityValue      =  2.8e-15;
734         double relativeToleranceVelocityDerivative =  1.3e-10;
735         doTestCartesianDerivatives(FastMath.toRadians(-35), FastMath.toRadians(20), 1200.0, 100.0,
736                                    relativeTolerancePositionValue, relativeTolerancePositionDerivative,
737                                    relativeToleranceVelocityValue, relativeToleranceVelocityDerivative,
738                                    ".*-offset-East", ".*-offset-North", ".*-offset-Zenith");
739     }
740 
741     @Test
742     void testStationOffsetAngularDerivativesOctantPxPyMz() {
743         double toleranceRotationValue            =  1.7e-15;
744         double toleranceRotationDerivative       =  3.6e-18;
745         double toleranceRotationRateValue        =  1.5e-19;
746         double toleranceRotationRateDerivative   =  Precision.SAFE_MIN;
747         doTestAngularDerivatives(FastMath.toRadians(-35), FastMath.toRadians(20), 1200.0, 100.0,
748                                  toleranceRotationValue,     toleranceRotationDerivative,
749                                  toleranceRotationRateValue, toleranceRotationRateDerivative,
750                                  ".*-offset-East", ".*-offset-North", ".*-offset-Zenith");
751     }
752 
753     @Test
754     void testStationOffsetCartesianDerivativesOctantPxMyPz() {
755         double relativeTolerancePositionValue      =  1.7e-15;
756         double relativeTolerancePositionDerivative =  9.0e-11;
757         double relativeToleranceVelocityValue      =  2.9e-15;
758         double relativeToleranceVelocityDerivative =  8.8e-11;
759         doTestCartesianDerivatives(FastMath.toRadians(35), FastMath.toRadians(-20), 1200.0, 100.0,
760                                    relativeTolerancePositionValue, relativeTolerancePositionDerivative,
761                                    relativeToleranceVelocityValue, relativeToleranceVelocityDerivative,
762                                    ".*-offset-East", ".*-offset-North", ".*-offset-Zenith");
763     }
764 
765     @Test
766     void testStationOffsetAngularDerivativesOctantPxMyPz() {
767         double toleranceRotationValue              =  1.6e-15;
768         double toleranceRotationDerivative         =  3.6e-18;
769         double toleranceRotationRateValue          =  1.5e-19;
770         double toleranceRotationRateDerivative     =  Precision.SAFE_MIN;
771         doTestAngularDerivatives(FastMath.toRadians(35), FastMath.toRadians(-20), 1200.0, 100.0,
772                                  toleranceRotationValue,     toleranceRotationDerivative,
773                                  toleranceRotationRateValue, toleranceRotationRateDerivative,
774                                  ".*-offset-East", ".*-offset-North", ".*-offset-Zenith");
775     }
776 
777     @Test
778     void testStationOffsetCartesianDerivativesOctantPxMyMz() {
779         double relativeTolerancePositionValue      =  1.5e-15;
780         double relativeTolerancePositionDerivative =  4.6e-11;
781         double relativeToleranceVelocityValue      =  2.7e-15;
782         double relativeToleranceVelocityDerivative =  6.8e-11;
783         doTestCartesianDerivatives(FastMath.toRadians(-35), FastMath.toRadians(-20), 1200.0, 100.0,
784                                    relativeTolerancePositionValue, relativeTolerancePositionDerivative,
785                                    relativeToleranceVelocityValue, relativeToleranceVelocityDerivative,
786                                    ".*-offset-East", ".*-offset-North", ".*-offset-Zenith");
787     }
788 
789     @Test
790     void testStationOffsetAngularDerivativesOctantPxMyMz() {
791         double toleranceRotationValue            =  1.6e-15;
792         double toleranceRotationDerivative       =  2.6e-18;
793         double toleranceRotationRateValue        =  1.5e-19;
794         double toleranceRotationRateDerivative   =  Precision.SAFE_MIN;
795         doTestAngularDerivatives(FastMath.toRadians(-35), FastMath.toRadians(-20), 1200.0, 100.0,
796                                  toleranceRotationValue,     toleranceRotationDerivative,
797                                  toleranceRotationRateValue, toleranceRotationRateDerivative,
798                                  ".*-offset-East", ".*-offset-North", ".*-offset-Zenith");
799     }
800 
801     @Test
802     void testStationOffsetCartesianDerivativesOctantMxPyPz() {
803         double relativeTolerancePositionValue      =  1.6e-15;
804         double relativeTolerancePositionDerivative =  9.9e-11;
805         double relativeToleranceVelocityValue      =  2.6e-15;
806         double relativeToleranceVelocityDerivative =  9.5e-11;
807         doTestCartesianDerivatives(FastMath.toRadians(150), FastMath.toRadians(20), 1200.0, 100.0,
808                                    relativeTolerancePositionValue, relativeTolerancePositionDerivative,
809                                    relativeToleranceVelocityValue, relativeToleranceVelocityDerivative,
810                                    ".*-offset-East", ".*-offset-North", ".*-offset-Zenith");
811     }
812 
813     @Test
814     void testStationOffsetAngularDerivativesOctantMxPyPz() {
815         double toleranceRotationValue            =  1.5e-15;
816         double toleranceRotationDerivative       =  3.1e-18;
817         double toleranceRotationRateValue        =  1.5e-19;
818         double toleranceRotationRateDerivative   =  Precision.SAFE_MIN;
819         doTestAngularDerivatives(FastMath.toRadians(150), FastMath.toRadians(20), 1200.0, 100.0,
820                                  toleranceRotationValue,     toleranceRotationDerivative,
821                                  toleranceRotationRateValue, toleranceRotationRateDerivative,
822                                  ".*-offset-East", ".*-offset-North", ".*-offset-Zenith");
823     }
824 
825     @Test
826     void testStationOffsetCartesianDerivativesOctantMxPyMz() {
827         double relativeTolerancePositionValue      =  1.6e-15;
828         double relativeTolerancePositionDerivative =  9.9e-11;
829         double relativeToleranceVelocityValue      =  2.6e-15;
830         double relativeToleranceVelocityDerivative =  9.5e-11;
831         doTestCartesianDerivatives(FastMath.toRadians(150), FastMath.toRadians(20), 1200.0, 100.0,
832                                    relativeTolerancePositionValue, relativeTolerancePositionDerivative,
833                                    relativeToleranceVelocityValue, relativeToleranceVelocityDerivative,
834                                    ".*-offset-East", ".*-offset-North", ".*-offset-Zenith");
835     }
836 
837     @Test
838     void testStationOffsetAngularDerivativesOctantMxPyMz() {
839         double toleranceRotationValue            =  1.5e-15;
840         double toleranceRotationDerivative       =  3.1e-18;
841         double toleranceRotationRateValue        =  1.5e-19;
842         double toleranceRotationRateDerivative   =  Precision.SAFE_MIN;
843         doTestAngularDerivatives(FastMath.toRadians(150), FastMath.toRadians(20), 1200.0, 100.0,
844                                  toleranceRotationValue,     toleranceRotationDerivative,
845                                  toleranceRotationRateValue, toleranceRotationRateDerivative,
846                                  ".*-offset-East", ".*-offset-North", ".*-offset-Zenith");
847     }
848 
849     @Test
850     void testStationOffsetCartesianDerivativesOctantMxMyPz() {
851         double relativeTolerancePositionValue      =  1.5e-15;
852         double relativeTolerancePositionDerivative =  6.2e-11;
853         double relativeToleranceVelocityValue      =  3.2e-15;
854         double relativeToleranceVelocityDerivative =  1.1e-10;
855         doTestCartesianDerivatives(FastMath.toRadians(-150), FastMath.toRadians(-20), 1200.0, 100.0,
856                                    relativeTolerancePositionValue, relativeTolerancePositionDerivative,
857                                    relativeToleranceVelocityValue, relativeToleranceVelocityDerivative,
858                                    ".*-offset-East", ".*-offset-North", ".*-offset-Zenith");
859     }
860 
861     @Test
862     void testStationOffsetAngularDerivativesOctantMxMyPz() {
863         double toleranceRotationValue            =  1.6e-15;
864         double toleranceRotationDerivative       =  3.4e-18;
865         double toleranceRotationRateValue        =  1.5e-19;
866         double toleranceRotationRateDerivative   =  Precision.SAFE_MIN;
867         doTestAngularDerivatives(FastMath.toRadians(-150), FastMath.toRadians(-20), 1200.0, 100.0,
868                                  toleranceRotationValue,     toleranceRotationDerivative,
869                                  toleranceRotationRateValue, toleranceRotationRateDerivative,
870                                  ".*-offset-East", ".*-offset-North", ".*-offset-Zenith");
871     }
872 
873     @Test
874     void testStationOffsetCartesianDerivativesOctantMxMyMz() {
875         double relativeTolerancePositionValue      =  1.5e-15;
876         double relativeTolerancePositionDerivative =  6.2e-11;
877         double relativeToleranceVelocityValue      =  3.2e-15;
878         double relativeToleranceVelocityDerivative =  1.1e-10;
879         doTestCartesianDerivatives(FastMath.toRadians(-150), FastMath.toRadians(-20), 1200.0, 100.0,
880                                    relativeTolerancePositionValue, relativeTolerancePositionDerivative,
881                                    relativeToleranceVelocityValue, relativeToleranceVelocityDerivative,
882                                    ".*-offset-East", ".*-offset-North", ".*-offset-Zenith");
883     }
884 
885     @Test
886     void testStationOffsetAngularDerivativesOctantMxMyMz() {
887         double toleranceRotationValue            =  1.6e-15;
888         double toleranceRotationDerivative       =  3.4e-18;
889         double toleranceRotationRateValue        =  1.5e-19;
890         double toleranceRotationRateDerivative   =  Precision.SAFE_MIN;
891         doTestAngularDerivatives(FastMath.toRadians(-150), FastMath.toRadians(-20), 1200.0, 100.0,
892                                  toleranceRotationValue,     toleranceRotationDerivative,
893                                  toleranceRotationRateValue, toleranceRotationRateDerivative,
894                                  ".*-offset-East", ".*-offset-North", ".*-offset-Zenith");
895     }
896 
897     @Test
898     void testStationOffsetCartesianDerivativesNearPole() {
899         double relativeTolerancePositionValue      =  2.4e-15;
900         double relativeTolerancePositionDerivative =  2e-9;
901         double relativeToleranceVelocityValue      =  7.5e-14;
902         double relativeToleranceVelocityDerivative =  3.9e-10;
903         doTestCartesianDerivatives(FastMath.toRadians(89.99995), FastMath.toRadians(90), 1200.0, 100.0,
904                                    relativeTolerancePositionValue, relativeTolerancePositionDerivative,
905                                    relativeToleranceVelocityValue, relativeToleranceVelocityDerivative,
906                                    ".*-offset-East", ".*-offset-North", ".*-offset-Zenith");
907     }
908 
909     @Test
910     void testStationOffsetAngularDerivativesNearPole() {
911         double toleranceRotationValue              =  3.9e-15;
912         double toleranceRotationDerivative         =  8.0e-02; // near pole, the East and North directions are singular
913         double toleranceRotationRateValue          =  1.5e-19;
914         double toleranceRotationRateDerivative     =  Precision.SAFE_MIN;
915         doTestAngularDerivatives(FastMath.toRadians(89.99995), FastMath.toRadians(90), 1200.0, 100.0,
916                                  toleranceRotationValue,     toleranceRotationDerivative,
917                                  toleranceRotationRateValue, toleranceRotationRateDerivative,
918                                  ".*-offset-East", ".*-offset-North", ".*-offset-Zenith");
919     }
920 
921     @Test
922     void testPolarMotionCartesianDerivativesOctantPxPyPz() {
923         double relativeTolerancePositionValue      =  2.3e-15;
924         double relativeTolerancePositionDerivative =  7.3e-09;
925         double relativeToleranceVelocityValue      =  3.3e-15;
926         double relativeToleranceVelocityDerivative =  4.8e-09;
927         doTestCartesianDerivatives(FastMath.toRadians(35), FastMath.toRadians(20), 1200.0, 1.0,
928                                    relativeTolerancePositionValue, relativeTolerancePositionDerivative,
929                                    relativeToleranceVelocityValue, relativeToleranceVelocityDerivative,
930                                    "polar-offset-X", "polar-drift-X",
931                                    "polar-offset-Y", "polar-drift-Y");
932     }
933 
934     @Test
935     void testPolarMotionAngularDerivativesOctantPxPyPz() {
936         double toleranceRotationValue              =  1.9e-15;
937         double toleranceRotationDerivative         =  6.6e-09;
938         double toleranceRotationRateValue          =  1.5e-19;
939         double toleranceRotationRateDerivative     =  2.2e-13;
940         doTestAngularDerivatives(FastMath.toRadians(35), FastMath.toRadians(20), 1200.0, 0.2,
941                                  toleranceRotationValue,     toleranceRotationDerivative,
942                                  toleranceRotationRateValue, toleranceRotationRateDerivative,
943                                  "polar-offset-X", "polar-drift-X",
944                                  "polar-offset-Y", "polar-drift-Y");
945     }
946 
947     @Test
948     void testPolarMotionCartesianDerivativesOctantPxPyMz() {
949         double relativeTolerancePositionValue      =  1.4e-15;
950         double relativeTolerancePositionDerivative =  4.1e-09;
951         double relativeToleranceVelocityValue      =  2.8e-15;
952         double relativeToleranceVelocityDerivative =  4.7e-09;
953         doTestCartesianDerivatives(FastMath.toRadians(-35), FastMath.toRadians(20), 1200.0, 1.0,
954                                    relativeTolerancePositionValue, relativeTolerancePositionDerivative,
955                                    relativeToleranceVelocityValue, relativeToleranceVelocityDerivative,
956                                    "polar-offset-X", "polar-drift-X",
957                                    "polar-offset-Y", "polar-drift-Y");
958     }
959 
960     @Test
961     void testPolarMotionAngularDerivativesOctantPxPyMz() {
962         double toleranceRotationValue              =  1.7e-15;
963         double toleranceRotationDerivative         =  6.9e-09;
964         double toleranceRotationRateValue          =  1.5e-19;
965         double toleranceRotationRateDerivative     =  2.2e-13;
966         doTestAngularDerivatives(FastMath.toRadians(-35), FastMath.toRadians(20), 1200.0, 0.2,
967                                  toleranceRotationValue,     toleranceRotationDerivative,
968                                  toleranceRotationRateValue, toleranceRotationRateDerivative,
969                                  "polar-offset-X", "polar-drift-X",
970                                  "polar-offset-Y", "polar-drift-Y");
971     }
972 
973     @Test
974     void testPolarMotionCartesianDerivativesOctantPxMyPz() {
975         double relativeTolerancePositionValue      =  1.7e-15;
976         double relativeTolerancePositionDerivative =  7.8e-09;
977         double relativeToleranceVelocityValue      =  2.9e-15;
978         double relativeToleranceVelocityDerivative =  4.4e-09;
979         doTestCartesianDerivatives(FastMath.toRadians(35), FastMath.toRadians(-20), 1200.0, 1.0,
980                                    relativeTolerancePositionValue, relativeTolerancePositionDerivative,
981                                    relativeToleranceVelocityValue, relativeToleranceVelocityDerivative,
982                                    "polar-offset-X", "polar-drift-X",
983                                    "polar-offset-Y", "polar-drift-Y");
984     }
985 
986     @Test
987     void testPolarMotionAngularDerivativesOctantPxMyPz() {
988         double toleranceRotationValue              =  1.6e-15;
989         double toleranceRotationDerivative         =  7.5e-09;
990         double toleranceRotationRateValue          =  1.5e-19;
991         double toleranceRotationRateDerivative     =  2.2e-13;
992         doTestAngularDerivatives(FastMath.toRadians(35), FastMath.toRadians(-20), 1200.0, 0.2,
993                                  toleranceRotationValue,     toleranceRotationDerivative,
994                                  toleranceRotationRateValue, toleranceRotationRateDerivative,
995                                  "polar-offset-X", "polar-drift-X",
996                                  "polar-offset-Y", "polar-drift-Y");
997     }
998 
999     @Test
1000     void testPolarMotionCartesianDerivativesOctantPxMyMz() {
1001         double relativeTolerancePositionValue      =  1.5e-15;
1002         double relativeTolerancePositionDerivative =  4.1e-09;
1003         double relativeToleranceVelocityValue      =  2.7e-15;
1004         double relativeToleranceVelocityDerivative =  4.4e-09;
1005         doTestCartesianDerivatives(FastMath.toRadians(-35), FastMath.toRadians(-20), 1200.0, 1.0,
1006                                    relativeTolerancePositionValue, relativeTolerancePositionDerivative,
1007                                    relativeToleranceVelocityValue, relativeToleranceVelocityDerivative,
1008                                    "polar-offset-X", "polar-drift-X",
1009                                    "polar-offset-Y", "polar-drift-Y");
1010     }
1011 
1012     @Test
1013     void testPolarMotionAngularDerivativesOctantPxMyMz() {
1014         double toleranceRotationValue            =  1.6e-15;
1015         double toleranceRotationDerivative       =  7.3e-09;
1016         double toleranceRotationRateValue        =  1.5e-19;
1017         double toleranceRotationRateDerivative   =  2.2e-13;
1018         doTestAngularDerivatives(FastMath.toRadians(-35), FastMath.toRadians(-20), 1200.0, 0.2,
1019                                  toleranceRotationValue,     toleranceRotationDerivative,
1020                                  toleranceRotationRateValue, toleranceRotationRateDerivative,
1021                                  "polar-offset-X", "polar-drift-X",
1022                                  "polar-offset-Y", "polar-drift-Y");
1023     }
1024 
1025     @Test
1026     void testPolarMotionCartesianDerivativesOctantMxPyPz() {
1027         double relativeTolerancePositionValue      =  1.6e-15;
1028         double relativeTolerancePositionDerivative =  8.4e-09;
1029         double relativeToleranceVelocityValue      =  2.6e-15;
1030         double relativeToleranceVelocityDerivative =  5.0e-09;
1031         doTestCartesianDerivatives(FastMath.toRadians(150), FastMath.toRadians(20), 1200.0, 1.0,
1032                                    relativeTolerancePositionValue, relativeTolerancePositionDerivative,
1033                                    relativeToleranceVelocityValue, relativeToleranceVelocityDerivative,
1034                                    "polar-offset-X", "polar-drift-X",
1035                                    "polar-offset-Y", "polar-drift-Y");
1036     }
1037 
1038     @Test
1039     void testPolarMotionAngularDerivativesOctantMxPyPz() {
1040         double toleranceRotationValue              =  1.4e-15;
1041         double toleranceRotationDerivative         =  6.3e-09;
1042         double toleranceRotationRateValue          =  1.5e-19;
1043         double toleranceRotationRateDerivative     =  2.2e-13;
1044         doTestAngularDerivatives(FastMath.toRadians(150), FastMath.toRadians(20), 1200.0, 0.2,
1045                                  toleranceRotationValue,     toleranceRotationDerivative,
1046                                  toleranceRotationRateValue, toleranceRotationRateDerivative,
1047                                  "polar-offset-X", "polar-drift-X",
1048                                  "polar-offset-Y", "polar-drift-Y");
1049     }
1050 
1051     @Test
1052     void testPolarMotionCartesianDerivativesOctantMxPyMz() {
1053         double relativeTolerancePositionValue      =  1.6e-15;
1054         double relativeTolerancePositionDerivative =  8.4e-09;
1055         double relativeToleranceVelocityValue      =  2.6e-15;
1056         double relativeToleranceVelocityDerivative =  5.0e-09;
1057         doTestCartesianDerivatives(FastMath.toRadians(150), FastMath.toRadians(20), 1200.0, 1.0,
1058                                    relativeTolerancePositionValue, relativeTolerancePositionDerivative,
1059                                    relativeToleranceVelocityValue, relativeToleranceVelocityDerivative,
1060                                    "polar-offset-X", "polar-drift-X",
1061                                    "polar-offset-Y", "polar-drift-Y");
1062     }
1063 
1064     @Test
1065     void testPolarMotionAngularDerivativesOctantMxPyMz() {
1066         double toleranceRotationValue              =  1.4e-15;
1067         double toleranceRotationDerivative         =  6.3e-09;
1068         double toleranceRotationRateValue          =  1.5e-19;
1069         double toleranceRotationRateDerivative     =  2.2e-13;
1070         doTestAngularDerivatives(FastMath.toRadians(150), FastMath.toRadians(20), 1200.0, 0.2,
1071                                  toleranceRotationValue,     toleranceRotationDerivative,
1072                                  toleranceRotationRateValue, toleranceRotationRateDerivative,
1073                                  "polar-offset-X", "polar-drift-X",
1074                                  "polar-offset-Y", "polar-drift-Y");
1075     }
1076 
1077     @Test
1078     void testPolarMotionCartesianDerivativesOctantMxMyPz() {
1079         double relativeTolerancePositionValue      =  1.5e-15;
1080         double relativeTolerancePositionDerivative =  5.0e-09;
1081         double relativeToleranceVelocityValue      =  3.2e-15;
1082         double relativeToleranceVelocityDerivative =  5.3e-09;
1083         doTestCartesianDerivatives(FastMath.toRadians(-150), FastMath.toRadians(-20), 1200.0, 1.0,
1084                                    relativeTolerancePositionValue, relativeTolerancePositionDerivative,
1085                                    relativeToleranceVelocityValue, relativeToleranceVelocityDerivative,
1086                                    "polar-offset-X", "polar-drift-X",
1087                                    "polar-offset-Y", "polar-drift-Y");
1088     }
1089 
1090     @Test
1091     void testPolarMotionAngularDerivativesOctantMxMyPz() {
1092         double toleranceRotationValue              =  1.6e-15;
1093         double toleranceRotationDerivative         =  7.7e-09;
1094         double toleranceRotationRateValue          =  1.5e-19;
1095         double toleranceRotationRateDerivative     =  2.2e-13;
1096         doTestAngularDerivatives(FastMath.toRadians(-150), FastMath.toRadians(-20), 1200.0, 0.2,
1097                                  toleranceRotationValue,     toleranceRotationDerivative,
1098                                  toleranceRotationRateValue, toleranceRotationRateDerivative,
1099                                  "polar-offset-X", "polar-drift-X",
1100                                  "polar-offset-Y", "polar-drift-Y");
1101     }
1102 
1103     @Test
1104     void testPolarMotionCartesianDerivativesOctantMxMyMz() {
1105         double relativeTolerancePositionValue      =  1.5e-15;
1106         double relativeTolerancePositionDerivative =  5.0e-09;
1107         double relativeToleranceVelocityValue      =  3.2e-15;
1108         double relativeToleranceVelocityDerivative =  5.3e-09;
1109         doTestCartesianDerivatives(FastMath.toRadians(-150), FastMath.toRadians(-20), 1200.0, 1.0,
1110                                    relativeTolerancePositionValue, relativeTolerancePositionDerivative,
1111                                    relativeToleranceVelocityValue, relativeToleranceVelocityDerivative,
1112                                    "polar-offset-X", "polar-drift-X",
1113                                    "polar-offset-Y", "polar-drift-Y");
1114     }
1115 
1116     @Test
1117     void testPolarMotionAngularDerivativesOctantMxMyMz() {
1118         double toleranceRotationValue              =  1.6e-15;
1119         double toleranceRotationDerivative         =  7.7e-09;
1120         double toleranceRotationRateValue          =  1.5e-19;
1121         double toleranceRotationRateDerivative     =  2.2e-13;
1122         doTestAngularDerivatives(FastMath.toRadians(-150), FastMath.toRadians(-20), 1200.0, 0.2,
1123                                  toleranceRotationValue,     toleranceRotationDerivative,
1124                                  toleranceRotationRateValue, toleranceRotationRateDerivative,
1125                                  "polar-offset-X", "polar-drift-X",
1126                                  "polar-offset-Y", "polar-drift-Y");
1127     }
1128 
1129     @Test
1130     void testPolarMotionCartesianDerivativesNearPole() {
1131         double relativeTolerancePositionValue      =  1.2e-15;
1132         double relativeTolerancePositionDerivative =  5.7e-09;
1133         double relativeToleranceVelocityValue      =  9.7e-13;
1134         double relativeToleranceVelocityDerivative =  1.2e-09;
1135         doTestCartesianDerivatives(FastMath.toRadians(89.99995), FastMath.toRadians(90), 1200.0, 1.0,
1136                                    relativeTolerancePositionValue, relativeTolerancePositionDerivative,
1137                                    relativeToleranceVelocityValue, relativeToleranceVelocityDerivative,
1138                                    "polar-offset-X", "polar-drift-X",
1139                                    "polar-offset-Y", "polar-drift-Y");
1140     }
1141 
1142     @Test
1143     void testPolarMotionAngularDerivativesNearPole() {
1144         double toleranceRotationValue              =  1.5e-15;
1145         double toleranceRotationDerivative         =  6.5e-09;
1146         double toleranceRotationRateValue          =  1.5e-19;
1147         double toleranceRotationRateDerivative     =  2.2e-13;
1148         doTestAngularDerivatives(FastMath.toRadians(89.99995), FastMath.toRadians(90), 1200.0, 0.2,
1149                                  toleranceRotationValue,     toleranceRotationDerivative,
1150                                  toleranceRotationRateValue, toleranceRotationRateDerivative,
1151                                  "polar-offset-X", "polar-drift-X",
1152                                  "polar-offset-Y", "polar-drift-Y");
1153     }
1154 
1155     @Test
1156     void testPrimeMeridianCartesianDerivativesOctantPxPyPz() {
1157         double relativeTolerancePositionValue      =  2.3e-15;
1158         double relativeTolerancePositionDerivative =  5.7e-09;
1159         double relativeToleranceVelocityValue      =  3.3e-15;
1160         double relativeToleranceVelocityDerivative =  2.7e-09;
1161         doTestCartesianDerivatives(FastMath.toRadians(35), FastMath.toRadians(20), 1200.0, 1.0,
1162                                    relativeTolerancePositionValue, relativeTolerancePositionDerivative,
1163                                    relativeToleranceVelocityValue, relativeToleranceVelocityDerivative,
1164                                    "prime-meridian-offset", "prime-meridian-drift");
1165     }
1166 
1167     @Test
1168     void testPrimeMeridianAngularDerivativesOctantPxPyPz() {
1169         double toleranceRotationValue              =  1.9e-15;
1170         double toleranceRotationDerivative         =  4.9e-09;
1171         double toleranceRotationRateValue          =  1.5e-19;
1172         double toleranceRotationRateDerivative     =  2.0e-13;
1173         doTestAngularDerivatives(FastMath.toRadians(35), FastMath.toRadians(20), 1200.0, 0.2,
1174                                  toleranceRotationValue,     toleranceRotationDerivative,
1175                                  toleranceRotationRateValue, toleranceRotationRateDerivative,
1176                                  "prime-meridian-offset", "prime-meridian-drift");
1177     }
1178 
1179     @Test
1180     void testPrimeMeridianCartesianDerivativesOctantPxPyMz() {
1181         double relativeTolerancePositionValue      =  1.4e-15;
1182         double relativeTolerancePositionDerivative =  3.1e-09;
1183         double relativeToleranceVelocityValue      =  2.8e-15;
1184         double relativeToleranceVelocityDerivative =  3.3e-09;
1185         doTestCartesianDerivatives(FastMath.toRadians(-35), FastMath.toRadians(20), 1200.0, 1.0,
1186                                    relativeTolerancePositionValue, relativeTolerancePositionDerivative,
1187                                    relativeToleranceVelocityValue, relativeToleranceVelocityDerivative,
1188                                    "prime-meridian-offset", "prime-meridian-drift");
1189     }
1190 
1191     @Test
1192     void testPrimeMeridianAngularDerivativesOctantPxPyMz() {
1193         double toleranceRotationValue            =  1.7e-15;
1194         double toleranceRotationDerivative       =  5.7e-09;
1195         double toleranceRotationRateValue        =  1.5e-19;
1196         double toleranceRotationRateDerivative   =  2.0e-13;
1197         doTestAngularDerivatives(FastMath.toRadians(-35), FastMath.toRadians(20), 1200.0, 0.2,
1198                                  toleranceRotationValue,     toleranceRotationDerivative,
1199                                  toleranceRotationRateValue, toleranceRotationRateDerivative,
1200                                  "prime-meridian-offset", "prime-meridian-drift");
1201     }
1202 
1203     @Test
1204     void testPrimeMeridianCartesianDerivativesOctantPxMyPz() {
1205         double relativeTolerancePositionValue      =  1.7e-15;
1206         double relativeTolerancePositionDerivative =  5.3e-09;
1207         double relativeToleranceVelocityValue      =  2.9e-15;
1208         double relativeToleranceVelocityDerivative =  3.4e-09;
1209         doTestCartesianDerivatives(FastMath.toRadians(35), FastMath.toRadians(-20), 1200.0, 1.0,
1210                                    relativeTolerancePositionValue, relativeTolerancePositionDerivative,
1211                                    relativeToleranceVelocityValue, relativeToleranceVelocityDerivative,
1212                                    "prime-meridian-offset", "prime-meridian-drift");
1213     }
1214 
1215     @Test
1216     void testPrimeMeridianAngularDerivativesOctantPxMyPz() {
1217         double toleranceRotationValue              =  1.6e-15;
1218         double toleranceRotationDerivative         =  6.4e-09;
1219         double toleranceRotationRateValue          =  1.5e-19;
1220         double toleranceRotationRateDerivative     =  2.0e-13;
1221         doTestAngularDerivatives(FastMath.toRadians(35), FastMath.toRadians(-20), 1200.0, 0.2,
1222                                  toleranceRotationValue,     toleranceRotationDerivative,
1223                                  toleranceRotationRateValue, toleranceRotationRateDerivative,
1224                                  "prime-meridian-offset", "prime-meridian-drift");
1225     }
1226 
1227     @Test
1228     void testPrimeMeridianCartesianDerivativesOctantPxMyMz() {
1229         double relativeTolerancePositionValue      =  1.5e-15;
1230         double relativeTolerancePositionDerivative =  3.1e-09;
1231         double relativeToleranceVelocityValue      =  2.7e-15;
1232         double relativeToleranceVelocityDerivative =  3.0e-09;
1233         doTestCartesianDerivatives(FastMath.toRadians(-35), FastMath.toRadians(-20), 1200.0, 1.0,
1234                                    relativeTolerancePositionValue, relativeTolerancePositionDerivative,
1235                                    relativeToleranceVelocityValue, relativeToleranceVelocityDerivative,
1236                                    "prime-meridian-offset", "prime-meridian-drift");
1237     }
1238 
1239     @Test
1240     void testPrimeMeridianAngularDerivativesOctantPxMyMz() {
1241         double toleranceRotationValue              =  1.5e-15;
1242         double toleranceRotationDerivative         =  6.3e-09;
1243         double toleranceRotationRateValue          =  1.5e-19;
1244         double toleranceRotationRateDerivative     =  2.0e-13;
1245         doTestAngularDerivatives(FastMath.toRadians(-35), FastMath.toRadians(-20), 1200.0, 0.2,
1246                                  toleranceRotationValue,     toleranceRotationDerivative,
1247                                  toleranceRotationRateValue, toleranceRotationRateDerivative,
1248                                  "prime-meridian-offset", "prime-meridian-drift");
1249     }
1250 
1251     @Test
1252     void testPrimeMeridianCartesianDerivativesOctantMxPyPz() {
1253         double relativeTolerancePositionValue      =  1.6e-15;
1254         double relativeTolerancePositionDerivative =  5.5e-09;
1255         double relativeToleranceVelocityValue      =  2.6e-15;
1256         double relativeToleranceVelocityDerivative =  2.8e-09;
1257         doTestCartesianDerivatives(FastMath.toRadians(150), FastMath.toRadians(20), 1200.0, 1.0,
1258                                    relativeTolerancePositionValue, relativeTolerancePositionDerivative,
1259                                    relativeToleranceVelocityValue, relativeToleranceVelocityDerivative,
1260                                    "prime-meridian-offset", "prime-meridian-drift");
1261     }
1262 
1263     @Test
1264     void testPrimeMeridianAngularDerivativesOctantMxPyPz() {
1265         double toleranceRotationValue            =  1.5e-15;
1266         double toleranceRotationDerivative       =  6.8e-09;
1267         double toleranceRotationRateValue        =  1.5e-19;
1268         double toleranceRotationRateDerivative   =  2.0e-13;
1269         doTestAngularDerivatives(FastMath.toRadians(150), FastMath.toRadians(20), 1200.0, 0.2,
1270                                  toleranceRotationValue,     toleranceRotationDerivative,
1271                                  toleranceRotationRateValue, toleranceRotationRateDerivative,
1272                                  "prime-meridian-offset", "prime-meridian-drift");
1273     }
1274 
1275     @Test
1276     void testPrimeMeridianCartesianDerivativesOctantMxPyMz() {
1277         double relativeTolerancePositionValue      =  1.6e-15;
1278         double relativeTolerancePositionDerivative =  5.5e-09;
1279         double relativeToleranceVelocityValue      =  2.6e-15;
1280         double relativeToleranceVelocityDerivative =  2.8e-09;
1281         doTestCartesianDerivatives(FastMath.toRadians(150), FastMath.toRadians(20), 1200.0, 1.0,
1282                                    relativeTolerancePositionValue, relativeTolerancePositionDerivative,
1283                                    relativeToleranceVelocityValue, relativeToleranceVelocityDerivative,
1284                                    "prime-meridian-offset", "prime-meridian-drift");
1285     }
1286 
1287     @Test
1288     void testPrimeMeridianAngularDerivativesOctantMxPyMz() {
1289         double toleranceRotationValue            =  1.5e-15;
1290         double toleranceRotationDerivative       =  6.8e-09;
1291         double toleranceRotationRateValue        =  1.5e-19;
1292         double toleranceRotationRateDerivative   =  2.0e-13;
1293         doTestAngularDerivatives(FastMath.toRadians(150), FastMath.toRadians(20), 1200.0, 0.2,
1294                                  toleranceRotationValue,     toleranceRotationDerivative,
1295                                  toleranceRotationRateValue, toleranceRotationRateDerivative,
1296                                  "prime-meridian-offset", "prime-meridian-drift");
1297     }
1298 
1299     @Test
1300     void testPrimeMeridianCartesianDerivativesOctantMxMyPz() {
1301         double relativeTolerancePositionValue      =  1.5e-15;
1302         double relativeTolerancePositionDerivative =  3.1e-09;
1303         double relativeToleranceVelocityValue      =  3.2e-15;
1304         double relativeToleranceVelocityDerivative =  3.2e-09;
1305         doTestCartesianDerivatives(FastMath.toRadians(-150), FastMath.toRadians(-20), 1200.0, 1.0,
1306                                    relativeTolerancePositionValue, relativeTolerancePositionDerivative,
1307                                    relativeToleranceVelocityValue, relativeToleranceVelocityDerivative,
1308                                    "prime-meridian-offset", "prime-meridian-drift");
1309     }
1310 
1311     @Test
1312     void testPrimeMeridianAngularDerivativesOctantMxMyPz() {
1313         double toleranceRotationValue              =  1.6e-15;
1314         double toleranceRotationDerivative         =  6.2e-09;
1315         double toleranceRotationRateValue          =  1.5e-19;
1316         double toleranceRotationRateDerivative     =  2.0e-13;
1317         doTestAngularDerivatives(FastMath.toRadians(-150), FastMath.toRadians(-20), 1200.0, 0.2,
1318                                  toleranceRotationValue,     toleranceRotationDerivative,
1319                                  toleranceRotationRateValue, toleranceRotationRateDerivative,
1320                                  "prime-meridian-offset", "prime-meridian-drift");
1321     }
1322 
1323     @Test
1324     void testPrimeMeridianCartesianDerivativesOctantMxMyMz() {
1325         double relativeTolerancePositionValue      =  1.5e-15;
1326         double relativeTolerancePositionDerivative =  3.1e-09;
1327         double relativeToleranceVelocityValue      =  3.2e-15;
1328         double relativeToleranceVelocityDerivative =  3.2e-09;
1329         doTestCartesianDerivatives(FastMath.toRadians(-150), FastMath.toRadians(-20), 1200.0, 1.0,
1330                                    relativeTolerancePositionValue, relativeTolerancePositionDerivative,
1331                                    relativeToleranceVelocityValue, relativeToleranceVelocityDerivative,
1332                                    "prime-meridian-offset", "prime-meridian-drift");
1333     }
1334 
1335     @Test
1336     void testPrimeMeridianAngularDerivativesOctantMxMyMz() {
1337         double toleranceRotationValue              =  1.6e-15;
1338         double toleranceRotationDerivative         =  6.2e-09;
1339         double toleranceRotationRateValue          =  1.5e-19;
1340         double toleranceRotationRateDerivative     =  2.0e-13;
1341         doTestAngularDerivatives(FastMath.toRadians(-150), FastMath.toRadians(-20), 1200.0, 0.2,
1342                                  toleranceRotationValue,     toleranceRotationDerivative,
1343                                  toleranceRotationRateValue, toleranceRotationRateDerivative,
1344                                  "prime-meridian-offset", "prime-meridian-drift");
1345     }
1346 
1347     @Test
1348     void testPrimeMeridianCartesianDerivativesNearPole() {
1349         double relativeTolerancePositionValue      =  1.2e-15;
1350         double relativeTolerancePositionDerivative =  1.6e-03;
1351         double relativeToleranceVelocityValue      =  5.8e-14;
1352         double relativeToleranceVelocityDerivative =  2.6e-08;
1353         doTestCartesianDerivatives(FastMath.toRadians(89.99995), FastMath.toRadians(90), 1200.0, 1.0,
1354                                    relativeTolerancePositionValue, relativeTolerancePositionDerivative,
1355                                    relativeToleranceVelocityValue, relativeToleranceVelocityDerivative,
1356                                    "prime-meridian-offset", "prime-meridian-drift");
1357     }
1358 
1359     @Test
1360     void testPrimeMeridianAngularDerivativesNearPole() {
1361         double toleranceRotationValue              =  1.5e-15;
1362         double toleranceRotationDerivative         =  7.3e-09;
1363         double toleranceRotationRateValue          =  1.5e-19;
1364         double toleranceRotationRateDerivative     =  2.0e-13;
1365         doTestAngularDerivatives(FastMath.toRadians(89.99995), FastMath.toRadians(90), 1200.0, 0.2,
1366                                  toleranceRotationValue,     toleranceRotationDerivative,
1367                                  toleranceRotationRateValue, toleranceRotationRateDerivative,
1368                                  "prime-meridian-offset", "prime-meridian-drift");
1369     }
1370 
1371     @Test
1372     void testNoReferenceDateGradient() {
1373         Utils.setDataRoot("regular-data");
1374         final Frame eme2000 = FramesFactory.getEME2000();
1375         final AbsoluteDate date = AbsoluteDate.J2000_EPOCH;
1376         final OneAxisEllipsoid earth =
1377                         new OneAxisEllipsoid(Constants.WGS84_EARTH_EQUATORIAL_RADIUS,
1378                                              Constants.WGS84_EARTH_FLATTENING,
1379                                              FramesFactory.getITRF(IERSConventions.IERS_2010, true));
1380         final GroundStation station = new GroundStation(new TopocentricFrame(earth,
1381                                                                              new GeodeticPoint(0.1, 0.2, 100),
1382                                                                              "dummy"));
1383         try {
1384             station.getOffsetToInertial(eme2000, date, false);
1385             Assertions.fail("an exception should have been thrown");
1386         } catch (OrekitException oe) {
1387             Assertions.assertEquals(OrekitMessages.NO_REFERENCE_DATE_FOR_PARAMETER, oe.getSpecifier());
1388             Assertions.assertEquals("prime-meridian-offset", oe.getParts()[0]);
1389         }
1390 
1391         try {
1392             int freeParameters = 9;
1393             Map<String, Integer> indices = new HashMap<>();
1394             for (final ParameterDriver driver : station.getParametersDrivers()) {
1395                 indices.put(driver.getNameSpan(date), indices.size());
1396             }
1397             station.getOffsetToInertial(eme2000, date, freeParameters, indices);
1398             Assertions.fail("an exception should have been thrown");
1399         } catch (OrekitException oe) {
1400             Assertions.assertEquals(OrekitMessages.NO_REFERENCE_DATE_FOR_PARAMETER, oe.getSpecifier());
1401             Assertions.assertEquals("prime-meridian-offset", oe.getParts()[0]);
1402         }
1403 
1404     }
1405 
1406     private void doTestCartesianDerivatives(double latitude, double longitude, double altitude, double stepFactor,
1407                                             double relativeTolerancePositionValue, double relativeTolerancePositionDerivative,
1408                                             double relativeToleranceVelocityValue, double relativeToleranceVelocityDerivative,
1409                                             String... parameterPattern) {
1410         Utils.setDataRoot("regular-data");
1411         final Frame eme2000 = FramesFactory.getEME2000();
1412         final AbsoluteDate date = AbsoluteDate.J2000_EPOCH;
1413         final AbsoluteDate date0 = date.shiftedBy(50000);
1414         final OneAxisEllipsoid earth =
1415                         new OneAxisEllipsoid(Constants.WGS84_EARTH_EQUATORIAL_RADIUS,
1416                                              Constants.WGS84_EARTH_FLATTENING,
1417                                              FramesFactory.getITRF(IERSConventions.IERS_2010, true));
1418         final GroundStation station = new GroundStation(new TopocentricFrame(earth,
1419                                                                              new GeodeticPoint(latitude, longitude, altitude),
1420                                                                              "dummy"));
1421         final GradientField gradientField = GradientField.getField(parameterPattern.length);
1422         ParameterDriver[] selectedDrivers = new ParameterDriver[parameterPattern.length];
1423         UnivariateDifferentiableVectorFunction[] dFCartesian = new UnivariateDifferentiableVectorFunction[parameterPattern.length];
1424         final ParameterDriver[] allDrivers = selectAllDrivers(station);
1425         for (ParameterDriver driver : allDrivers) {
1426             driver.setReferenceDate(date0);
1427         }
1428         Map<String, Integer> indices = new HashMap<>();
1429         for (int k = 0; k < dFCartesian.length; ++k) {
1430             for (final ParameterDriver allDriver : allDrivers) {
1431                 if (allDriver.getName().matches(parameterPattern[k])) {
1432                     selectedDrivers[k] = allDriver;
1433                     dFCartesian[k] = differentiatedStationPV(station, eme2000, date, selectedDrivers[k], stepFactor);
1434                     indices.put(selectedDrivers[k].getNameSpan(date0), k);
1435                 }
1436             }
1437         }
1438 
1439         RandomGenerator generator = new Well19937a(0x084d58a19c498a54L);
1440 
1441         double maxPositionValueRelativeError      = 0;
1442         double maxPositionDerivativeRelativeError = 0;
1443         double maxVelocityValueRelativeError      = 0;
1444         double maxVelocityDerivativeRelativeError = 0;
1445         for (int i = 0; i < 1000; ++i) {
1446 
1447             // randomly change one parameter
1448             ParameterDriver changed = allDrivers[generator.nextInt(allDrivers.length)];
1449             changed.setNormalizedValue(2 * generator.nextDouble() - 1);
1450 
1451             // transform to check
1452             FieldTransform<Gradient> t = station.getOffsetToInertial(eme2000, date, parameterPattern.length, indices);
1453             FieldPVCoordinates<Gradient> pv = t.transformPVCoordinates(FieldPVCoordinates.getZero(gradientField));
1454             for (int k = 0; k < dFCartesian.length; ++k) {
1455 
1456                 // reference values and derivatives computed using finite differences
1457                 Gradient[] refCartesian = dFCartesian[k].value(Gradient.variable(1, 0, selectedDrivers[k].getValue()));
1458 
1459                 // position
1460                 final Vector3D refP = new Vector3D(refCartesian[0].getValue(),
1461                                                    refCartesian[1].getValue(),
1462                                                    refCartesian[2].getValue());
1463                 final Vector3D resP = new Vector3D(pv.getPosition().getX().getValue(),
1464                                                    pv.getPosition().getY().getValue(),
1465                                                    pv.getPosition().getZ().getValue());
1466                 maxPositionValueRelativeError =
1467                                 FastMath.max(maxPositionValueRelativeError, Vector3D.distance(refP, resP) / refP.getNorm());
1468                 final Vector3D refPD = new Vector3D(refCartesian[0].getPartialDerivative(0),
1469                                                     refCartesian[1].getPartialDerivative(0),
1470                                                     refCartesian[2].getPartialDerivative(0));
1471                 final Vector3D resPD = new Vector3D(pv.getPosition().getX().getPartialDerivative(k),
1472                                                     pv.getPosition().getY().getPartialDerivative(k),
1473                                                     pv.getPosition().getZ().getPartialDerivative(k));
1474                 maxPositionDerivativeRelativeError =
1475                                 FastMath.max(maxPositionDerivativeRelativeError, Vector3D.distance(refPD, resPD) / refPD.getNorm());
1476 
1477                 // velocity
1478                 final Vector3D refV = new Vector3D(refCartesian[3].getValue(),
1479                                                    refCartesian[4].getValue(),
1480                                                    refCartesian[5].getValue());
1481                 final Vector3D resV = new Vector3D(pv.getVelocity().getX().getValue(),
1482                                                    pv.getVelocity().getY().getValue(),
1483                                                    pv.getVelocity().getZ().getValue());
1484                 maxVelocityValueRelativeError =
1485                                 FastMath.max(maxVelocityValueRelativeError, Vector3D.distance(refV, resV) / refV.getNorm());
1486                 final Vector3D refVD = new Vector3D(refCartesian[3].getPartialDerivative(0),
1487                                                     refCartesian[4].getPartialDerivative(0),
1488                                                     refCartesian[5].getPartialDerivative(0));
1489                 final Vector3D resVD = new Vector3D(pv.getVelocity().getX().getPartialDerivative(k),
1490                                                     pv.getVelocity().getY().getPartialDerivative(k),
1491                                                     pv.getVelocity().getZ().getPartialDerivative(k));
1492                 maxVelocityDerivativeRelativeError =
1493                                 FastMath.max(maxVelocityDerivativeRelativeError, Vector3D.distance(refVD, resVD) / refVD.getNorm());
1494 
1495             }
1496 
1497         }
1498 
1499         if (maxPositionValueRelativeError      > relativeTolerancePositionValue      ||
1500             maxPositionDerivativeRelativeError > relativeTolerancePositionDerivative ||
1501             maxVelocityValueRelativeError      > relativeToleranceVelocityValue      ||
1502             maxVelocityDerivativeRelativeError > relativeToleranceVelocityDerivative) {
1503             print("relativeTolerancePositionValue",      maxPositionValueRelativeError);
1504             print("relativeTolerancePositionDerivative", maxPositionDerivativeRelativeError);
1505             print("relativeToleranceVelocityValue",      maxVelocityValueRelativeError);
1506             print("relativeToleranceVelocityDerivative", maxVelocityDerivativeRelativeError);
1507         }
1508         Assertions.assertEquals(0.0, maxPositionValueRelativeError,      relativeTolerancePositionValue);
1509         Assertions.assertEquals(0.0, maxPositionDerivativeRelativeError, relativeTolerancePositionDerivative);
1510         Assertions.assertEquals(0.0, maxVelocityValueRelativeError,      relativeToleranceVelocityValue);
1511         Assertions.assertEquals(0.0, maxVelocityDerivativeRelativeError, relativeToleranceVelocityDerivative);
1512 
1513     }
1514 
1515     private void doTestAngularDerivatives(double latitude, double longitude, double altitude, double stepFactor,
1516                                           double toleranceRotationValue,     double toleranceRotationDerivative,
1517                                           double toleranceRotationRateValue, double toleranceRotationRateDerivative,
1518                                           String... parameterPattern) {
1519         Utils.setDataRoot("regular-data");
1520         final Frame eme2000 = FramesFactory.getEME2000();
1521         final AbsoluteDate date = AbsoluteDate.J2000_EPOCH;
1522         final AbsoluteDate date0 = date.shiftedBy(50000);
1523         final OneAxisEllipsoid earth =
1524                         new OneAxisEllipsoid(Constants.WGS84_EARTH_EQUATORIAL_RADIUS,
1525                                              Constants.WGS84_EARTH_FLATTENING,
1526                                              FramesFactory.getITRF(IERSConventions.IERS_2010, true));
1527         final GroundStation station = new GroundStation(new TopocentricFrame(earth,
1528                                                                              new GeodeticPoint(latitude, longitude, altitude),
1529                                                                              "dummy"));
1530         ParameterDriver[] selectedDrivers = new ParameterDriver[parameterPattern.length];
1531         UnivariateDifferentiableVectorFunction[] dFAngular   = new UnivariateDifferentiableVectorFunction[parameterPattern.length];
1532         final ParameterDriver[] allDrivers = selectAllDrivers(station);
1533         for (ParameterDriver driver : allDrivers) {
1534             driver.setReferenceDate(date0);
1535         }
1536         Map<String, Integer> indices = new HashMap<>();
1537         for (int k = 0; k < dFAngular.length; ++k) {
1538             for (final ParameterDriver allDriver : allDrivers) {
1539                 if (allDriver.getName().matches(parameterPattern[k])) {
1540                     selectedDrivers[k] = allDriver;
1541                     dFAngular[k] = differentiatedTransformAngular(station, eme2000, date, selectedDrivers[k], stepFactor);
1542                     indices.put(selectedDrivers[k].getNameSpan(date0), k);
1543                 }
1544             }
1545         }
1546         RandomGenerator generator = new Well19937a(0xa01a1d8fe5d80af7L);
1547 
1548         double maxRotationValueError          = 0;
1549         double maxRotationDerivativeError     = 0;
1550         double maxRotationRateValueError      = 0;
1551         double maxRotationRateDerivativeError = 0;
1552         for (int i = 0; i < 1000; ++i) {
1553 
1554             // randomly change one parameter
1555             ParameterDriver changed = allDrivers[generator.nextInt(allDrivers.length)];
1556             changed.setNormalizedValue(2 * generator.nextDouble() - 1);
1557 
1558             // transform to check
1559             FieldTransform<Gradient> t = station.getOffsetToInertial(eme2000, date, parameterPattern.length, indices);
1560             for (int k = 0; k < dFAngular.length; ++k) {
1561 
1562                 // reference values and derivatives computed using finite differences
1563                 Gradient[] refAngular = dFAngular[k].value(Gradient.variable(1, 0, selectedDrivers[k].getValue()));
1564 
1565                 // rotation
1566                 final Rotation refQ = new Rotation(refAngular[0].getValue(),
1567                                                    refAngular[1].getValue(),
1568                                                    refAngular[2].getValue(),
1569                                                    refAngular[3].getValue(),
1570                                                    true);
1571                 final Rotation resQ = t.getRotation().toRotation();
1572                 maxRotationValueError      = FastMath.max(maxRotationValueError, Rotation.distance(refQ, resQ));
1573                 double sign = FastMath.copySign(1.0,
1574                                                 refAngular[0].getValue() * t.getRotation().getQ0().getValue() +
1575                                                 refAngular[1].getValue() * t.getRotation().getQ1().getValue() +
1576                                                 refAngular[2].getValue() * t.getRotation().getQ2().getValue() +
1577                                                 refAngular[3].getValue() * t.getRotation().getQ3().getValue());
1578                 maxRotationDerivativeError = FastMath.max(maxRotationDerivativeError,
1579                                                           FastMath.abs(sign * refAngular[0].getPartialDerivative(0) -
1580                                                                        t.getRotation().getQ0().getPartialDerivative(k)));
1581                 maxRotationDerivativeError = FastMath.max(maxRotationDerivativeError,
1582                                                           FastMath.abs(sign * refAngular[1].getPartialDerivative(0) -
1583                                                                        t.getRotation().getQ1().getPartialDerivative(k)));
1584                 maxRotationDerivativeError = FastMath.max(maxRotationDerivativeError,
1585                                                           FastMath.abs(sign * refAngular[2].getPartialDerivative(0) -
1586                                                                        t.getRotation().getQ2().getPartialDerivative(k)));
1587                 maxRotationDerivativeError = FastMath.max(maxRotationDerivativeError,
1588                                                           FastMath.abs(sign * refAngular[3].getPartialDerivative(0) -
1589                                                                        t.getRotation().getQ3().getPartialDerivative(k)));
1590 
1591                 // rotation rate
1592                 final Vector3D refRate  = new Vector3D(refAngular[4].getValue(), refAngular[5].getValue(), refAngular[6].getValue());
1593                 final Vector3D resRate  = t.getRotationRate().toVector3D();
1594                 final Vector3D refRateD = new Vector3D(refAngular[4].getPartialDerivative(0),
1595                                                        refAngular[5].getPartialDerivative(0),
1596                                                        refAngular[6].getPartialDerivative(0));
1597                 final Vector3D resRateD = new Vector3D(t.getRotationRate().getX().getPartialDerivative(k),
1598                                                        t.getRotationRate().getY().getPartialDerivative(k),
1599                                                        t.getRotationRate().getZ().getPartialDerivative(k));
1600                 maxRotationRateValueError      = FastMath.max(maxRotationRateValueError, Vector3D.distance(refRate, resRate));
1601                 maxRotationRateDerivativeError = FastMath.max(maxRotationRateDerivativeError, Vector3D.distance(refRateD, resRateD));
1602 
1603             }
1604 
1605         }
1606 
1607         if (maxRotationValueError          > toleranceRotationValue          ||
1608             maxRotationDerivativeError     > toleranceRotationDerivative     ||
1609             maxRotationRateValueError      > toleranceRotationRateValue      ||
1610             maxRotationRateDerivativeError > toleranceRotationRateDerivative) {
1611             print("toleranceRotationValue",          maxRotationValueError);
1612             print("toleranceRotationDerivative",     maxRotationDerivativeError);
1613             print("toleranceRotationRateValue",      maxRotationRateValueError);
1614             print("toleranceRotationRateDerivative", maxRotationRateDerivativeError);
1615         }
1616         Assertions.assertEquals(0.0, maxRotationValueError,           toleranceRotationValue);
1617         Assertions.assertEquals(0.0, maxRotationDerivativeError,      toleranceRotationDerivative);
1618         Assertions.assertEquals(0.0, maxRotationRateValueError,       toleranceRotationRateValue);
1619         Assertions.assertEquals(0.0, maxRotationRateDerivativeError,  toleranceRotationRateDerivative);
1620 
1621     }
1622     private void print(String name, double v) {
1623         if (v < Precision.SAFE_MIN) {
1624             System.out.format(Locale.US, "            double %-35s =  Precision.SAFE_MIN;%n", name);
1625         } else {
1626             double s = FastMath.pow(10.0, 1 - FastMath.floor(FastMath.log(v) / FastMath.log(10.0)));
1627             System.out.format(Locale.US, "            double %-35s = %8.1e;%n",
1628                               name, FastMath.ceil(s * v) / s);
1629         }
1630     }
1631 
1632     private UnivariateDifferentiableVectorFunction differentiatedStationPV(final GroundStation station,
1633                                                                            final Frame eme2000,
1634                                                                            final AbsoluteDate date,
1635                                                                            final ParameterDriver driver,
1636                                                                            final double stepFactor) {
1637 
1638         final FiniteDifferencesDifferentiator differentiator =
1639                         new FiniteDifferencesDifferentiator(5, stepFactor * driver.getScale());
1640 
1641         return differentiator.differentiate((UnivariateVectorFunction) x -> {
1642             final double[] result = new double[6];
1643             try {
1644                 final double previouspI = driver.getValue(date);
1645                 driver.setValue(x, new AbsoluteDate());
1646                 Transform t = station.getOffsetToInertial(eme2000, date, false);
1647                 driver.setValue(previouspI, date);
1648                 PVCoordinates stationPV = t.transformPVCoordinates(PVCoordinates.ZERO);
1649                 result[ 0] = stationPV.getPosition().getX();
1650                 result[ 1] = stationPV.getPosition().getY();
1651                 result[ 2] = stationPV.getPosition().getZ();
1652                 result[ 3] = stationPV.getVelocity().getX();
1653                 result[ 4] = stationPV.getVelocity().getY();
1654                 result[ 5] = stationPV.getVelocity().getZ();
1655             } catch (OrekitException oe) {
1656                 Assertions.fail(oe.getLocalizedMessage());
1657             }
1658             return result;
1659         });
1660     }
1661 
1662     private UnivariateDifferentiableVectorFunction differentiatedTransformAngular(final GroundStation station,
1663                                                                                   final Frame eme2000,
1664                                                                                   final AbsoluteDate date,
1665                                                                                   final ParameterDriver driver,
1666                                                                                   final double stepFactor) {
1667 
1668         final FiniteDifferencesDifferentiator differentiator =
1669                         new FiniteDifferencesDifferentiator(5, stepFactor * driver.getScale());
1670 
1671         return differentiator.differentiate(new UnivariateVectorFunction() {
1672             private double previous0 = Double.NaN;
1673             private double previous1 = Double.NaN;
1674             private double previous2 = Double.NaN;
1675             private double previous3 = Double.NaN;
1676             @Override
1677             public double[] value(double x) {
1678                 final double[] result = new double[7];
1679                 try {
1680                     final double previouspI = driver.getValue(date);
1681                     driver.setValue(x, date);
1682                     Transform t = station.getOffsetToInertial(eme2000, date, false);
1683                     driver.setValue(previouspI, date);
1684                     final double sign;
1685                     if (Double.isNaN(previous0)) {
1686                         sign = +1;
1687                     } else {
1688                         sign = FastMath.copySign(1.0,
1689                                                  previous0 * t.getRotation().getQ0() +
1690                                                  previous1 * t.getRotation().getQ1() +
1691                                                  previous2 * t.getRotation().getQ2() +
1692                                                  previous3 * t.getRotation().getQ3());
1693                     }
1694                     previous0 = sign * t.getRotation().getQ0();
1695                     previous1 = sign * t.getRotation().getQ1();
1696                     previous2 = sign * t.getRotation().getQ2();
1697                     previous3 = sign * t.getRotation().getQ3();
1698                     result[0] = previous0;
1699                     result[1] = previous1;
1700                     result[2] = previous2;
1701                     result[3] = previous3;
1702                     result[4] = t.getRotationRate().getX();
1703                     result[5] = t.getRotationRate().getY();
1704                     result[6] = t.getRotationRate().getZ();
1705                 } catch (OrekitException oe) {
1706                     Assertions.fail(oe.getLocalizedMessage());
1707                 }
1708                 return result;
1709             }
1710         });
1711     }
1712 
1713     private ParameterDriver[] selectAllDrivers(final GroundStation station) {
1714         return new ParameterDriver[] {
1715             station.getPrimeMeridianOffsetDriver(),
1716             station.getPrimeMeridianDriftDriver(),
1717             station.getPolarOffsetXDriver(),
1718             station.getPolarDriftXDriver(),
1719             station.getPolarOffsetYDriver(),
1720             station.getPolarDriftYDriver(),
1721             station.getClockBiasDriver(),
1722             station.getEastOffsetDriver(),
1723             station.getNorthOffsetDriver(),
1724             station.getZenithOffsetDriver()
1725         };
1726     }
1727 
1728 }
1729