1   /* Copyright 2002-2018 CS Systèmes d'Information
2    * Licensed to CS Systèmes d'Information (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.utils;
18  
19  
20  import org.hipparchus.analysis.UnivariateFunction;
21  import org.hipparchus.analysis.UnivariateVectorFunction;
22  import org.hipparchus.analysis.differentiation.DSFactory;
23  import org.hipparchus.analysis.differentiation.DerivativeStructure;
24  import org.hipparchus.analysis.differentiation.FiniteDifferencesDifferentiator;
25  import org.hipparchus.analysis.differentiation.UnivariateDifferentiableFunction;
26  import org.hipparchus.analysis.differentiation.UnivariateDifferentiableVectorFunction;
27  import org.hipparchus.util.Decimal64;
28  import org.hipparchus.util.FastMath;
29  import org.hipparchus.util.MathUtils;
30  import org.junit.Assert;
31  import org.junit.Before;
32  import org.junit.Test;
33  import org.orekit.Utils;
34  import org.orekit.errors.OrekitException;
35  import org.orekit.frames.EOPHistory;
36  import org.orekit.frames.FramesFactory;
37  import org.orekit.frames.ITRFVersion;
38  import org.orekit.time.AbsoluteDate;
39  import org.orekit.time.TimeVectorFunction;
40  import org.orekit.time.DateComponents;
41  import org.orekit.time.FieldAbsoluteDate;
42  import org.orekit.time.TimeScalarFunction;
43  import org.orekit.time.TimeScale;
44  import org.orekit.time.TimeScalesFactory;
45  import org.orekit.time.UT1Scale;
46  
47  
48  public class IERSConventionsTest {
49  
50      @Test
51      public void testConventionsNumber() {
52          Assert.assertEquals(3, IERSConventions.values().length);
53      }
54  
55      @Test
56      public void testIERS1996NutationAngles() throws OrekitException {
57  
58          // the reference value has been computed using the March 2012 version of the SOFA library
59          // http://www.iausofa.org/2012_0301_C.html, with the following code
60          //
61          //        double utc1, utc2, tai1, tai2, tt1, tt2, dpsi, deps, epsa;
62          //
63          //        // 2004-02-14:00:00:00Z, MJD = 53049, UT1-UTC = -0.4093509
64          //        utc1  = DJM0 + 53049.0;
65          //        utc2  = 0.0;
66          //        iauUtctai(utc1, utc2, &tai1, &tai2);
67          //        iauTaitt(tai1, tai2, &tt1, &tt2);
68          //
69          //        iauNut80(tt1, tt2, &dpsi, &deps);
70          //        epsa = iauObl80(tt1, tt2);
71          //
72          //        printf("iauNut80(%.20g, %.20g, dpsi, deps)\n  --> %.20g %.20g\n",
73          //               tt1, tt2, dpsi, deps);
74          //        printf("iau0bl80(%.20g, %.20g)\n  --> %.20g\n",
75          //               tt1, tt2, epsa);
76          //
77          // the output of this test reads:
78          //      iauNut80(2453049.5, 0.00074287037037037029902, dpsi, deps)
79          //        --> -5.3059154211478291722e-05 3.2051803135750973851e-05
80          //      iau0bl80(2453049.5, 0.00074287037037037029902)
81          //        --> 0.40908345528446415917
82  
83          // the thresholds below have been set up large enough to have the test pass with
84          // the default Orekit setting and the reference SOFA setting. There are implementation
85          // differences between the two libraries, as the Delaunay parameters are not the
86          // same (in Orekit, we are compliant with page 23 in chapter 5 of IERS conventions 1996.
87          // If the content of the IERS-conventions/1996/nutation-arguments.txt file by Orekit is
88          // changed to match SOFA setting instead of IERS conventions as follows:
89          //
90          //        # Mean Anomaly of the Moon
91          //        F1 ≡ l = 485866.733″ + 1717915922.633″t + 31.310″t² + 0.064″t³
92          //
93          //       # Mean Anomaly of the Sun
94          //        F2 ≡ l' = 1287099.804″ + 129596581.224″t − 0.577″t² − 0.012t³
95          //
96          //       # L − Ω
97          //        F3 ≡ F = 335778.877″ + 1739527263.137″t − 13.257″t² + 0.011″t³
98          //
99          //       # Mean Elongation of the Moon from the Sun
100         //        F4 ≡ D = 1072261.307″ + 1602961601.328″t − 6.891″t² + 0.019″t³
101         //
102         //       # Mean Longitude of the Ascending Node of the Moon
103         //        F5 ≡ Ω = 450160.280″ − 6962890.539″t + 7.455″t² + 0.008″t³
104         //
105         // then the thresholds for the test can be reduced to 3e-13 for ∆ψ and 7e-14 for ∆ε.
106         // We decided to stick with IERS published reference values for the default Orekit setting.
107         // The differences are nevertheless quite small (4.8e-11 radians is sub-millimeter level
108         // in low Earth orbit).
109         AbsoluteDate date = new AbsoluteDate(2004, 2, 14, TimeScalesFactory.getUTC());
110         double[] angles= IERSConventions.IERS_1996.getNutationFunction().value(date);
111         Assert.assertEquals(-5.3059154211478291722e-05, angles[0], 4.8e-11); // 3e-13 with SOFA values
112         Assert.assertEquals(3.2051803135750973851e-05,  angles[1], 1.3e-11); // 7e-14 with SOFA values
113 
114         double epsilonA = IERSConventions.IERS_1996.getMeanObliquityFunction().value(date);
115         Assert.assertEquals(0.40908345528446415917,     epsilonA, 1.0e-16);
116 
117     }
118 
119     @Test
120     public void testGMST82FieldConsistency() throws OrekitException {
121         final TimeScale ut1 = TimeScalesFactory.getUT1(IERSConventions.IERS_1996, true);
122         checkScalarFunctionConsistency(IERSConventions.IERS_1996.getGMSTFunction(ut1),
123                                        AbsoluteDate.J2000_EPOCH.shiftedBy(-0.4 * Constants.JULIAN_DAY),
124                                        0.8 * Constants.JULIAN_DAY, 600.0, 10.0, 9.0e-16, 8.0e-13);
125     }
126 
127     @Test
128     public void testGMST82Sofa() throws OrekitException {
129 
130         // the reference value has been computed using the March 2012 version of the SOFA library
131         // http://www.iausofa.org/2012_0301_C.html, with the following code
132         //
133         //        double utc1, utc2, ut11, ut12, gmst;
134         //
135         //        // 2004-02-14:00:00:00Z, MJD = 53049, UT1-UTC = -0.4093509
136         //        utc1  = DJM0 + 53049.0;
137         //        utc2  = 0.0;
138         //        iauUtcut1(utc1, utc2, -0.4093509, &ut11, &ut12);
139         //        gmst = iauGmst82(ut11, ut12);
140         //        printf("iaugmst82(%.20g, %.20g)\n  --> %.20g\n", ut11, ut12, gmst);
141         //
142         //        // 2004-02-29:00:00:00Z, MJD = 53064, UT1-UTC = -0.4175723
143         //        utc1 = DJM0 + 53064.0;
144         //        utc2 = 0.0;
145         //        iauUtcut1(utc1, utc2, -0.4175723, &ut11, &ut12);
146         //        gmst = iauGmst82(ut11, ut12);
147         //        printf("iaugmst82(%.20g, %.20g)\n  --> %.20g\n", ut11, ut12, gmst);
148         //
149         // the output of this test reads:
150         //      iaugmst82(2453049.5, -4.7378576388888813016e-06)
151         //        --> 2.5021977627453466653
152         //      iaugmst82(2453064.5, -4.8330127314815448519e-06)
153         //        --> 2.7602390405411441066
154 
155         final TimeScale ut1 = TimeScalesFactory.getUT1(IERSConventions.IERS_1996, true);
156         final TimeScalarFunction gmst82 = IERSConventions.IERS_1996.getGMSTFunction(ut1);
157         AbsoluteDate date = new AbsoluteDate(2004, 2, 14, TimeScalesFactory.getUTC());
158         double gmst = MathUtils.normalizeAngle(gmst82.value(date), 0.0);
159         Assert.assertEquals(2.5021977627453466653, gmst, 2.0e-13);
160         date = new AbsoluteDate(2004, 2, 29, TimeScalesFactory.getUTC());
161         gmst = MathUtils.normalizeAngle(gmst82.value(date), 0.0);
162         Assert.assertEquals(2.7602390405411441066, gmst, 4.0e-13);
163 
164     }
165 
166     @Test
167     public void testGMST00FieldConsistency() throws OrekitException {
168         final TimeScale ut1 = TimeScalesFactory.getUT1(IERSConventions.IERS_2003, true);
169         checkScalarFunctionConsistency(IERSConventions.IERS_2003.getGMSTFunction(ut1),
170                                        AbsoluteDate.J2000_EPOCH.shiftedBy(-0.4 * Constants.JULIAN_DAY),
171                                        0.8 * Constants.JULIAN_DAY, 600.0, 10.0, 2.0e-15, 8.0e-13);
172     }
173 
174     @Test
175     public void testGMST00Sofa() throws OrekitException {
176 
177         // the reference value has been computed using the March 2012 version of the SOFA library
178         // http://www.iausofa.org/2012_0301_C.html, with the following code
179         //
180         //        double utc1, utc2, tai1, tai2, tt1, tt2, ut11, ut12, gmst;
181         //
182         //        // 2004-02-14:00:00:00Z, MJD = 53049, UT1-UTC = -0.4093509
183         //        utc1  = DJM0 + 53049.0;
184         //        utc2  = 0.0;
185         //        iauUtctai(utc1, utc2, &tai1, &tai2);
186         //        iauTaitt(tai1, tai2, &tt1, &tt2);
187         //        iauUtcut1(utc1, utc2, -0.4093509, &ut11, &ut12);
188         //        gmst = iauGmst00(ut11, ut12, tt1, tt2);
189         //        printf("iaugmst00(%.20g, %.20g, %.20g, %.20g)\n  --> %.20g\n",
190         //               ut11, ut12, tt1, tt2, gmst);
191         //
192         //        // 2004-02-29:00:00:00Z, MJD = 53064, UT1-UTC = -0.4175723
193         //        utc1 = DJM0 + 53064.0;
194         //        utc2 = 0.0;
195         //        iauUtctai(utc1, utc2, &tai1, &tai2);
196         //        iauTaitt(tai1, tai2, &tt1, &tt2);
197         //        iauUtcut1(utc1, utc2, -0.4175723, &ut11, &ut12);
198         //        gmst = iauGmst00(ut11, ut12, tt1, tt2);
199         //        printf("iaugmst00(%.20g, %.20g, %.20g, %.20g)\n  --> %.20g\n",
200         //               ut11, ut12, tt1, tt2, gmst);
201         //
202         // the output of this test reads:
203         //      iaugmst00(2453049.5, -4.7378576388888813016e-06, 2453049.5, 0.00074287037037037029902)
204         //        --> 2.5021977786243789765
205         //      iaugmst00(2453064.5, -4.8330127314815448519e-06, 2453064.5, 0.00074287037037037029902)
206         //        --> 2.7602390558728311376
207 
208         final TimeScale ut1 = TimeScalesFactory.getUT1(IERSConventions.IERS_2003, true);
209         final TimeScalarFunction gmst00 = IERSConventions.IERS_2003.getGMSTFunction(ut1);
210         AbsoluteDate date = new AbsoluteDate(2004, 2, 14, TimeScalesFactory.getUTC());
211         double gmst = MathUtils.normalizeAngle(gmst00.value(date), 0.0);
212         Assert.assertEquals(2.5021977786243789765, gmst, 1.0e-15);
213         date = new AbsoluteDate(2004, 2, 29, TimeScalesFactory.getUTC());
214         gmst = MathUtils.normalizeAngle(gmst00.value(date), 0.0);
215         Assert.assertEquals(2.7602390558728311376, gmst, 1.0e-15);
216 
217     }
218 
219     @Test
220     public void testGMST06FieldConsistency() throws OrekitException {
221         final TimeScale ut1 = TimeScalesFactory.getUT1(IERSConventions.IERS_2010, true);
222         checkScalarFunctionConsistency(IERSConventions.IERS_2010.getGMSTFunction(ut1),
223                                        AbsoluteDate.J2000_EPOCH.shiftedBy(-0.4 * Constants.JULIAN_DAY),
224                                        0.8 * Constants.JULIAN_DAY, 600.0, 10.0, 2.0e-15, 8.0e-13);
225     }
226 
227     @Test
228     public void testGMST06Sofa() throws OrekitException {
229 
230         // the reference value has been computed using the March 2012 version of the SOFA library
231         // http://www.iausofa.org/2012_0301_C.html, with the following code
232         //
233         //        double utc1, utc2, tai1, tai2, tt1, tt2, ut11, ut12, gmst;
234         //
235         //        // 2004-02-14:00:00:00Z, MJD = 53049, UT1-UTC = -0.4093509
236         //        utc1  = DJM0 + 53049.0;
237         //        utc2  = 0.0;
238         //        iauUtctai(utc1, utc2, &tai1, &tai2);
239         //        iauTaitt(tai1, tai2, &tt1, &tt2);
240         //        iauUtcut1(utc1, utc2, -0.4093509, &ut11, &ut12);
241         //        gmst = iauGmst06(ut11, ut12, tt1, tt2);
242         //        printf("iaugmst06(%.20g, %.20g, %.20g, %.20g)\n  --> %.20g\n",
243         //               ut11, ut12, tt1, tt2, gmst);
244         //
245         //        // 2004-02-29:00:00:00Z, MJD = 53064, UT1-UTC = -0.4175723
246         //        utc1 = DJM0 + 53064.0;
247         //        utc2 = 0.0;
248         //        iauUtctai(utc1, utc2, &tai1, &tai2);
249         //        iauTaitt(tai1, tai2, &tt1, &tt2);
250         //        iauUtcut1(utc1, utc2, -0.4175723, &ut11, &ut12);
251         //        gmst = iauGmst06(ut11, ut12, tt1, tt2);
252         //        printf("iaugmst06(%.20g, %.20g, %.20g, %.20g)\n  --> %.20g\n",
253         //               ut11, ut12, tt1, tt2, gmst);
254         //
255         // the output of this test reads:
256         //      iaugmst06(2453049.5, -4.7378576388888813016e-06, 2453049.5, 0.00074287037037037029902)
257         //        --> 2.5021977784096232078
258         //      iaugmst06(2453064.5, -4.8330127314815448519e-06, 2453064.5, 0.00074287037037037029902)
259         //        --> 2.7602390556555129741
260 
261         final TimeScale ut1 = TimeScalesFactory.getUT1(IERSConventions.IERS_2010, true);
262         final TimeScalarFunction gmst06 = IERSConventions.IERS_2010.getGMSTFunction(ut1);
263         AbsoluteDate date = new AbsoluteDate(2004, 2, 14, TimeScalesFactory.getUTC());
264         double gmst = MathUtils.normalizeAngle(gmst06.value(date), 0.0);
265         Assert.assertEquals(2.5021977784096232078, gmst, 1.0e-15);
266         date = new AbsoluteDate(2004, 2, 29, TimeScalesFactory.getUTC());
267         gmst = MathUtils.normalizeAngle(gmst06.value(date), 0.0);
268         Assert.assertEquals(2.7602390556555129741, gmst, 3.0e-15);
269 
270     }
271 
272     @Test
273     public void testERA1996FieldConsistency() throws OrekitException {
274         final TimeScale ut1 = TimeScalesFactory.getUT1(IERSConventions.IERS_1996, true);
275         checkScalarFunctionConsistency(IERSConventions.IERS_1996.getEarthOrientationAngleFunction(ut1),
276                                        AbsoluteDate.J2000_EPOCH.shiftedBy(-0.4 * Constants.JULIAN_DAY),
277                                        0.8 * Constants.JULIAN_DAY, 1800.0, 10.0, 2.0e-15, 8.0e-13);
278     }
279 
280     @Test
281     public void testERA2003FieldConsistency() throws OrekitException {
282         final TimeScale ut1 = TimeScalesFactory.getUT1(IERSConventions.IERS_2003, true);
283         checkScalarFunctionConsistency(IERSConventions.IERS_2003.getEarthOrientationAngleFunction(ut1),
284                                        AbsoluteDate.J2000_EPOCH.shiftedBy(-0.4 * Constants.JULIAN_DAY),
285                                        0.8 * Constants.JULIAN_DAY, 1800.0, 10.0, 2.0e-15, 8.0e-13);
286     }
287 
288     @Test
289     public void testERA2010FieldConsistency() throws OrekitException {
290         final TimeScale ut1 = TimeScalesFactory.getUT1(IERSConventions.IERS_2010, true);
291         checkScalarFunctionConsistency(IERSConventions.IERS_2010.getEarthOrientationAngleFunction(ut1),
292                                        AbsoluteDate.J2000_EPOCH.shiftedBy(-0.4 * Constants.JULIAN_DAY),
293                                        0.8 * Constants.JULIAN_DAY, 1800.0, 10.0, 2.0e-15, 8.0e-13);
294     }
295 
296     @Test
297     public void testERASofa() throws OrekitException {
298 
299         // the reference value has been computed using the March 2012 version of the SOFA library
300         // http://www.iausofa.org/2012_0301_C.html, with the following code
301         //
302         //        double tai1, tai2, ut11, ut12, era, taiutc;
303         //        taiutc = 32.0;
304         //
305         //        // 2004-02-14:00:00:00Z, MJD = 53049, UT1-UTC = -0.4093509
306         //        tai1 = DJM0 + 53049.0;
307         //        tai2 = taiutc / DAYSEC;
308         //        iauTaiut1(tai1, tai2, -0.4093509 - taiutc, &ut11, &ut12);
309         //        era = iauEra00(ut11, ut12);
310         //        printf("iauera00(%.20g, %.20g)\n  --> %.20g\n", ut11, ut12, era);
311         //
312         //        // 2004-02-29:00:00:00Z, MJD = 53064, UT1-UTC = -0.4175723
313         //        tai1 = DJM0 + 53064.0;
314         //        tai2 = taiutc / DAYSEC;
315         //        iauTaiut1(tai1, tai2, -0.4175723 - taiutc, &ut11, &ut12);
316         //        era = iauEra00(ut11, ut12);
317         //        printf("iauera00(%.20g, %.20g)\n  --> %.20g\n", ut11, ut12, era);
318         //
319         // the output of this test reads:
320         //      iauera00(2453049.5, -4.7378576388888813016e-06)
321         //        --> 2.5012766511308228701
322         //      iauera00(2453064.5, -4.8330127314815448519e-06)
323         //        --> 2.7593087452455264952
324 
325         final TimeScale ut1 = TimeScalesFactory.getUT1(IERSConventions.IERS_2003, true);
326         final TimeScalarFunction era00 =
327                 IERSConventions.IERS_2003.getEarthOrientationAngleFunction(ut1);
328         AbsoluteDate date = new AbsoluteDate(2004, 2, 14, TimeScalesFactory.getUTC());
329         double era  = MathUtils.normalizeAngle(era00.value(date), 0.0);
330         Assert.assertEquals(2.5012766511308228701, era, 1.0e-15);
331         date = new AbsoluteDate(2004, 2, 29, TimeScalesFactory.getUTC());
332         era  = MathUtils.normalizeAngle(era00.value(date), 0.0);
333         Assert.assertEquals(2.7593087452455264952, era, 1.0e-15);
334 
335     }
336 
337     @Test
338     public void testGST94FieldConsistency() throws OrekitException {
339         EOPHistory eopHistory = FramesFactory.getEOPHistory(IERSConventions.IERS_1996, true);
340         checkScalarFunctionConsistency(IERSConventions.IERS_1996.getGASTFunction(TimeScalesFactory.getUT1(eopHistory),
341                                                                                  eopHistory),
342                                        AbsoluteDate.J2000_EPOCH.shiftedBy(-0.4 * Constants.JULIAN_DAY),
343                                        0.8 * Constants.JULIAN_DAY, 1800.0, 10.0, 9.0e-16, 8.0e-13);
344     }
345 
346     @Test
347     public void testGST94FieldConsistencyNoEOP() throws OrekitException {
348         checkScalarFunctionConsistency(IERSConventions.IERS_1996.getGASTFunction(TimeScalesFactory.getUT1(null), null),
349                                        AbsoluteDate.J2000_EPOCH.shiftedBy(-0.4 * Constants.JULIAN_DAY),
350                                        0.8 * Constants.JULIAN_DAY, 1800.0, 10.0, 9.0e-16, 2.0e-16);
351     }
352 
353     @Test
354     public void testGST94Sofa() throws OrekitException {
355 
356         // the reference value has been computed using the March 2012 version of the SOFA library
357         // http://www.iausofa.org/2012_0301_C.html, with the following code
358         //
359         //        double utc1, utc2, ut11, ut12, tai1, tai2, tt1, tt2, gmst82, eqeq94;
360         //
361         //        // 2004-02-14:00:00:00Z, MJD = 53049, UT1-UTC = -0.4093509
362         //        utc1  = DJM0 + 53049.0;
363         //        utc2  = 0.0;
364         //        iauUtcut1(utc1, utc2, -0.4093509, &ut11, &ut12);
365         //        iauUtctai(utc1, utc2, &tai1, &tai2);
366         //        iauTaitt(tai1, tai2, &tt1, &tt2);
367         //        gmst82 = iauGmst82(ut11, ut12);
368         //        eqeq94 = iauEqeq94(tt1, tt2);
369         //        printf("iauGmst82(%.20g, %.20g)\n  --> %.20g\n", ut11, ut12, gmst82);
370         //        printf("iauEqeq94(%.20g, %.20g)\n  --> %.20g\n", tt1, tt2, eqeq94);
371         //        printf(" gmst82 + eqeq94  --> %.20g\n", gmst82 + eqeq94);
372         //
373         //        // 2004-02-29:00:00:00Z, MJD = 53064, UT1-UTC = -0.4175723
374         //        utc1 = DJM0 + 53064.0;
375         //        utc2 = 0.0;
376         //        iauUtcut1(utc1, utc2, -0.4175723, &ut11, &ut12);
377         //        iauUtctai(utc1, utc2, &tai1, &tai2);
378         //        iauTaitt(tai1, tai2, &tt1, &tt2);
379         //        gmst82 = iauGmst82(ut11, ut12);
380         //        eqeq94 = iauEqeq94(tt1, tt2);
381         //        printf("iauGmst82(%.20g, %.20g)\n  --> %.20g\n", ut11, ut12, gmst82);
382         //        printf("iauEqeq94(%.20g, %.20g)\n  --> %.20g\n", tt1, tt2, eqeq94);
383         //        printf(" gmst82 + eqeq94  --> %.20g\n", gmst82 + eqeq94);
384         //
385         // the output of this test reads:
386         //      iauGmst82(2453049.5, -4.7378576388888813016e-06)
387         //        --> 2.5021977627453466653
388         //      iauEqeq94(2453049.5, 0.00074287037037037029902)
389         //        --> -4.8671604682267536886e-05
390         //       gmst82 + eqeq94  --> 2.5021490911406645274
391         //      iauGmst82(2453064.5, -4.8330127314815448519e-06)
392         //        --> 2.7602390405411441066
393         //      iauEqeq94(2453064.5, 0.00074287037037037029902)
394         //        --> -4.8893054690771762302e-05
395         //       gmst82 + eqeq94  --> 2.7601901474864534158
396 
397         // As can be seen in the code above, we didn't call iauGst94, because as
398         // stated in the function header comment in SOFA source files:
399         //        "accuracy has been compromised for the sake of
400         //         convenience in that UT is used instead of TDB (or TT) to compute
401         //         the equation of the equinoxes."
402         // So we rather performed the date conversion and then called ourselves iauGmst82
403         // with a date in UTC and eqe94 with a date in TT, restoring full accuracy in the
404         // SOFA computation for this test.
405 
406         // The thresholds below have been set up large enough to have the test pass with
407         // the default Orekit setting and the reference SOFA setting. There are implementation
408         // differences between the two libraries, as the Delaunay parameters are not the
409         // same (in Orekit, we are compliant with page 23 in chapter 5 of IERS conventions 1996.
410         // If the content of the IERS-conventions/1996/nutation-arguments.txt file by Orekit is
411         // changed to match SOFA setting instead of IERS conventions as follows:
412         //
413         //        # Mean Anomaly of the Moon
414         //        F1 ≡ l = 485866.733″ + 1717915922.633″t + 31.310″t² + 0.064″t³
415         //
416         //       # Mean Anomaly of the Sun
417         //        F2 ≡ l' = 1287099.804″ + 129596581.224″t − 0.577″t² − 0.012t³
418         //
419         //       # L − Ω
420         //        F3 ≡ F = 335778.877″ + 1739527263.137″t − 13.257″t² + 0.011″t³
421         //
422         //       # Mean Elongation of the Moon from the Sun
423         //        F4 ≡ D = 1072261.307″ + 1602961601.328″t − 6.891″t² + 0.019″t³
424         //
425         //       # Mean Longitude of the Ascending Node of the Moon
426         //        F5 ≡ Ω = 450160.280″ − 6962890.539″t + 7.455″t² + 0.008″t³
427         //
428         // then the thresholds for the test can be reduced to 4.2e-13 for the first test,
429         // and 4.7e-13 for the second test.
430         // We decided to stick with IERS published reference values for the default Orekit setting.
431         // The differences are nevertheless quite small (9e-11 radians is sub-millimeter level
432         // in low Earth orbit).
433         Utils.setLoaders(IERSConventions.IERS_1996,
434                          Utils.buildEOPList(IERSConventions.IERS_1996, ITRFVersion.ITRF_2008, new double[][] {
435                              { 53047, -0.4093509, 0.00, 0.00, 0.00, 0.00, 0.00, 0.00, 0.00 },
436                              { 53048, -0.4093509, 0.00, 0.00, 0.00, 0.00, 0.00, 0.00, 0.00 },
437                              { 53049, -0.4093509, 0.00, 0.00, 0.00, 0.00, 0.00, 0.00, 0.00 },
438                              { 53050, -0.4093509, 0.00, 0.00, 0.00, 0.00, 0.00, 0.00, 0.00 },
439                              { 53051, -0.4093509, 0.00, 0.00, 0.00, 0.00, 0.00, 0.00, 0.00 },
440                              { 53052, -0.4093509, 0.00, 0.00, 0.00, 0.00, 0.00, 0.00, 0.00 },
441                              { 53053, -0.4093509, 0.00, 0.00, 0.00, 0.00, 0.00, 0.00, 0.00 },
442                              { 53054, -0.4093509, 0.00, 0.00, 0.00, 0.00, 0.00, 0.00, 0.00 },
443                              // ...
444                              { 53059, -0.4175723, 0.00, 0.00, 0.00, 0.00, 0.00, 0.00, 0.00 },
445                              { 53060, -0.4175723, 0.00, 0.00, 0.00, 0.00, 0.00, 0.00, 0.00 },
446                              { 53061, -0.4175723, 0.00, 0.00, 0.00, 0.00, 0.00, 0.00, 0.00 },
447                              { 53062, -0.4175723, 0.00, 0.00, 0.00, 0.00, 0.00, 0.00, 0.00 },
448                              { 53063, -0.4175723, 0.00, 0.00, 0.00, 0.00, 0.00, 0.00, 0.00 },
449                              { 53064, -0.4175723, 0.00, 0.00, 0.00, 0.00, 0.00, 0.00, 0.00 },
450                              { 53065, -0.4175723, 0.00, 0.00, 0.00, 0.00, 0.00, 0.00, 0.00 },
451                              { 53066, -0.4175723, 0.00, 0.00, 0.00, 0.00, 0.00, 0.00, 0.00 }
452                          }));
453         EOPHistory eopHistory = FramesFactory.getEOPHistory(IERSConventions.IERS_1996, true);
454         final TimeScalarFunction gst94 =
455                 IERSConventions.IERS_1996.getGASTFunction(TimeScalesFactory.getUT1(eopHistory),
456                                                           eopHistory);
457         AbsoluteDate date = new AbsoluteDate(2004, 2, 14, TimeScalesFactory.getUTC());
458         double gst = MathUtils.normalizeAngle(gst94.value(date), 0.0);
459         Assert.assertEquals(2.5021490911406645274, gst, 4.4e-11); // 4.2e-13 with SOFA values
460         date = new AbsoluteDate(2004, 2, 29, TimeScalesFactory.getUTC());
461         gst = MathUtils.normalizeAngle(gst94.value(date), 0.0);
462         Assert.assertEquals(2.7601901474864534158, gst, 9.0e-11); // 4.7e-13 with SOFA values
463 
464     }
465 
466     @Test
467     public void testGST00AFieldConsistency() throws OrekitException {
468         EOPHistory eopHistory = FramesFactory.getEOPHistory(IERSConventions.IERS_2003, true);
469         checkScalarFunctionConsistency(IERSConventions.IERS_2003.getGASTFunction(TimeScalesFactory.getUT1(eopHistory), eopHistory),
470                                        AbsoluteDate.J2000_EPOCH.shiftedBy(-0.4 * Constants.JULIAN_DAY),
471                                        0.8 * Constants.JULIAN_DAY, 1800.0, 10.0, 2.0e-15, 8.0e-13);
472     }
473 
474     @Test
475     public void testGST00AFieldConsistencyNoEOP() throws OrekitException {
476         checkScalarFunctionConsistency(IERSConventions.IERS_2003.getGASTFunction(TimeScalesFactory.getUT1(null), null),
477                                        AbsoluteDate.J2000_EPOCH.shiftedBy(-0.4 * Constants.JULIAN_DAY),
478                                        0.8 * Constants.JULIAN_DAY, 1800.0, 10.0, 2.0e-15, 2.0e-16);
479     }
480 
481     @Test
482     public void testGST00ASofa() throws OrekitException {
483 
484         // the reference value has been computed using the March 2012 version of the SOFA library
485         // http://www.iausofa.org/2012_0301_C.html, with the following code
486         //
487         //        double utc1, utc2, tai1, tai2, tt1, tt2, ut11, ut12, gst;
488         //
489         //        // 2004-02-14:00:00:00Z, MJD = 53049, UT1-UTC = -0.4093509
490         //        utc1  = DJM0 + 53049.0;
491         //        utc2  = 0.0;
492         //        iauUtctai(utc1, utc2, &tai1, &tai2);
493         //        iauTaitt(tai1, tai2, &tt1, &tt2);
494         //        iauUtcut1(utc1, utc2, -0.4093509, &ut11, &ut12);
495         //        gst = iauGst00a(ut11, ut12, tt1, tt2);
496         //        printf("iaugst00a(%.20g, %.20g, %.20g, %.20g)\n  --> %.20g\n",
497         //               ut11, ut12, tt1, tt2, gst);
498         //
499         //        // 2004-02-29:00:00:00Z, MJD = 53064, UT1-UTC = -0.4175723
500         //        utc1 = DJM0 + 53064.0;
501         //        utc2 = 0.0;
502         //        iauUtctai(utc1, utc2, &tai1, &tai2);
503         //        iauTaitt(tai1, tai2, &tt1, &tt2);
504         //        iauUtcut1(utc1, utc2, -0.4175723, &ut11, &ut12);
505         //        gst = iauGst00a(ut11, ut12, tt1, tt2);
506         //        printf("iaugst00a(%.20g, %.20g, %.20g, %.20g)\n  --> %.20g\n",
507         //               ut11, ut12, tt1, tt2, gst);
508         //
509         // the output of this test reads:
510         //      iaugst00a(2453049.5, -4.7378576388888813016e-06, 2453049.5, 0.00074287037037037029902)
511         //        --> 2.5021491024360624778
512         //      iaugst00a(2453064.5, -4.8330127314815448519e-06, 2453064.5, 0.00074287037037037029902)
513         //        --> 2.7601901615614221619
514 
515         Utils.setLoaders(IERSConventions.IERS_2003,
516                          Utils.buildEOPList(IERSConventions.IERS_2003, ITRFVersion.ITRF_2008, new double[][] {
517                              { 53047, -0.4093509, 0.00, 0.00, 0.00, 0.00, 0.00, 0.00, 0.00 },
518                              { 53048, -0.4093509, 0.00, 0.00, 0.00, 0.00, 0.00, 0.00, 0.00 },
519                              { 53049, -0.4093509, 0.00, 0.00, 0.00, 0.00, 0.00, 0.00, 0.00 },
520                              { 53050, -0.4093509, 0.00, 0.00, 0.00, 0.00, 0.00, 0.00, 0.00 },
521                              { 53051, -0.4093509, 0.00, 0.00, 0.00, 0.00, 0.00, 0.00, 0.00 },
522                              { 53052, -0.4093509, 0.00, 0.00, 0.00, 0.00, 0.00, 0.00, 0.00 },
523                              { 53053, -0.4093509, 0.00, 0.00, 0.00, 0.00, 0.00, 0.00, 0.00 },
524                              { 53054, -0.4093509, 0.00, 0.00, 0.00, 0.00, 0.00, 0.00, 0.00 },
525                              // ...
526                              { 53059, -0.4175723, 0.00, 0.00, 0.00, 0.00, 0.00, 0.00, 0.00 },
527                              { 53060, -0.4175723, 0.00, 0.00, 0.00, 0.00, 0.00, 0.00, 0.00 },
528                              { 53061, -0.4175723, 0.00, 0.00, 0.00, 0.00, 0.00, 0.00, 0.00 },
529                              { 53062, -0.4175723, 0.00, 0.00, 0.00, 0.00, 0.00, 0.00, 0.00 },
530                              { 53063, -0.4175723, 0.00, 0.00, 0.00, 0.00, 0.00, 0.00, 0.00 },
531                              { 53064, -0.4175723, 0.00, 0.00, 0.00, 0.00, 0.00, 0.00, 0.00 },
532                              { 53065, -0.4175723, 0.00, 0.00, 0.00, 0.00, 0.00, 0.00, 0.00 },
533                              { 53066, -0.4175723, 0.00, 0.00, 0.00, 0.00, 0.00, 0.00, 0.00 }
534                          }));
535         EOPHistory eopHistory = FramesFactory.getEOPHistory(IERSConventions.IERS_2003, true);
536         final TimeScalarFunction gst00a =
537                 IERSConventions.IERS_2003.getGASTFunction(TimeScalesFactory.getUT1(eopHistory), eopHistory);
538         AbsoluteDate date = new AbsoluteDate(2004, 2, 14, TimeScalesFactory.getUTC());
539         double gst = MathUtils.normalizeAngle(gst00a.value(date), 0.0);
540         Assert.assertEquals(2.5021491024360624778, gst, 2.0e-13);
541         date = new AbsoluteDate(2004, 2, 29, TimeScalesFactory.getUTC());
542         gst = MathUtils.normalizeAngle(gst00a.value(date), 0.0);
543         Assert.assertEquals(2.7601901615614221619, gst, 3.0e-13);
544 
545     }
546 
547     @Test
548     public void testGST06FieldConsistency() throws OrekitException {
549         EOPHistory eopHistory = FramesFactory.getEOPHistory(IERSConventions.IERS_2010, true);
550         checkScalarFunctionConsistency(IERSConventions.IERS_2010.getGASTFunction(TimeScalesFactory.getUT1(eopHistory), eopHistory),
551                                        AbsoluteDate.J2000_EPOCH.shiftedBy(-0.4 * Constants.JULIAN_DAY),
552                                        0.8 * Constants.JULIAN_DAY, 1800.0, 10.0, 2.0e-15, 8.0e-13);
553     }
554 
555     @Test
556     public void testGST06FieldConsistencyNoEOP() throws OrekitException {
557         checkScalarFunctionConsistency(IERSConventions.IERS_2010.getGASTFunction(TimeScalesFactory.getUT1(null), null),
558                                        AbsoluteDate.J2000_EPOCH.shiftedBy(-0.4 * Constants.JULIAN_DAY),
559                                        0.8 * Constants.JULIAN_DAY, 1800.0, 10.0, 2.0e-15, 2.0e-16);
560     }
561 
562     @Test
563     public void testGST06Sofa() throws OrekitException {
564 
565         // the reference value has been computed using the March 2012 version of the SOFA library
566         // http://www.iausofa.org/2012_0301_C.html, with the following code
567         //
568         //        double utc1, utc2, tai1, tai2, tt1, tt2, ut11, ut12, gst;
569         //
570         //        // 2004-02-14:00:00:00Z, MJD = 53049, UT1-UTC = -0.4093509
571         //        utc1  = DJM0 + 53049.0;
572         //        utc2  = 0.0;
573         //        iauUtctai(utc1, utc2, &tai1, &tai2);
574         //        iauTaitt(tai1, tai2, &tt1, &tt2);
575         //        iauUtcut1(utc1, utc2, -0.4093509, &ut11, &ut12);
576         //        gst = iauGst06a(ut11, ut12, tt1, tt2);
577         //        printf("iaugst06a(%.20g, %.20g, %.20g, %.20g)\n  --> %.20g\n",
578         //               ut11, ut12, tt1, tt2, gst);
579         //
580         //        // 2004-02-29:00:00:00Z, MJD = 53064, UT1-UTC = -0.4175723
581         //        utc1 = DJM0 + 53064.0;
582         //        utc2 = 0.0;
583         //        iauUtctai(utc1, utc2, &tai1, &tai2);
584         //        iauTaitt(tai1, tai2, &tt1, &tt2);
585         //        iauUtcut1(utc1, utc2, -0.4175723, &ut11, &ut12);
586         //        gst = iauGst06a(ut11, ut12, tt1, tt2);
587         //        printf("iaugst06a(%.20g, %.20g, %.20g, %.20g)\n  --> %.20g\n",
588         //               ut11, ut12, tt1, tt2, gst);
589         //
590         // the output of this test reads:
591         //      iaugst06a(2453049.5, -4.7378576388888813016e-06, 2453049.5, 0.00074287037037037029902)
592         //        --> 2.5021491022006503435
593         //      iaugst06a(2453064.5, -4.8330127314815448519e-06, 2453064.5, 0.00074287037037037029902)
594         //        --> 2.7601901613234058885
595 
596         Utils.setLoaders(IERSConventions.IERS_2010,
597                          Utils.buildEOPList(IERSConventions.IERS_2010, ITRFVersion.ITRF_2008, new double[][] {
598                              { 53047, -0.4093509, 0.00, 0.00, 0.00, 0.00, 0.00, 0.00, 0.00 },
599                              { 53048, -0.4093509, 0.00, 0.00, 0.00, 0.00, 0.00, 0.00, 0.00 },
600                              { 53049, -0.4093509, 0.00, 0.00, 0.00, 0.00, 0.00, 0.00, 0.00 },
601                              { 53050, -0.4093509, 0.00, 0.00, 0.00, 0.00, 0.00, 0.00, 0.00 },
602                              { 53051, -0.4093509, 0.00, 0.00, 0.00, 0.00, 0.00, 0.00, 0.00 },
603                              { 53052, -0.4093509, 0.00, 0.00, 0.00, 0.00, 0.00, 0.00, 0.00 },
604                              { 53053, -0.4093509, 0.00, 0.00, 0.00, 0.00, 0.00, 0.00, 0.00 },
605                              { 53054, -0.4093509, 0.00, 0.00, 0.00, 0.00, 0.00, 0.00, 0.00 },
606                              // ...
607                              { 53059, -0.4175723, 0.00, 0.00, 0.00, 0.00, 0.00, 0.00, 0.00 },
608                              { 53060, -0.4175723, 0.00, 0.00, 0.00, 0.00, 0.00, 0.00, 0.00 },
609                              { 53061, -0.4175723, 0.00, 0.00, 0.00, 0.00, 0.00, 0.00, 0.00 },
610                              { 53062, -0.4175723, 0.00, 0.00, 0.00, 0.00, 0.00, 0.00, 0.00 },
611                              { 53063, -0.4175723, 0.00, 0.00, 0.00, 0.00, 0.00, 0.00, 0.00 },
612                              { 53064, -0.4175723, 0.00, 0.00, 0.00, 0.00, 0.00, 0.00, 0.00 },
613                              { 53065, -0.4175723, 0.00, 0.00, 0.00, 0.00, 0.00, 0.00, 0.00 },
614                              { 53066, -0.4175723, 0.00, 0.00, 0.00, 0.00, 0.00, 0.00, 0.00 }
615                          }));
616         EOPHistory eopHistory = FramesFactory.getEOPHistory(IERSConventions.IERS_2010, true);
617         final TimeScalarFunction gst06 =
618                 IERSConventions.IERS_2010.getGASTFunction(TimeScalesFactory.getUT1(eopHistory), eopHistory);
619         AbsoluteDate date = new AbsoluteDate(2004, 2, 14, TimeScalesFactory.getUTC());
620         double gst = MathUtils.normalizeAngle(gst06.value(date), 0.0);
621         Assert.assertEquals(2.5021491022006503435, gst, 1.3e-12);
622         date = new AbsoluteDate(2004, 2, 29, TimeScalesFactory.getUTC());
623         gst = MathUtils.normalizeAngle(gst06.value(date), 0.0);
624         Assert.assertEquals(2.7601901613234058885, gst, 1.2e-12);
625 
626     }
627 
628     @Test
629     public void testGMST2000vs82() throws OrekitException {
630 
631         final TimeScalarFunction gmst82 =
632                 IERSConventions.IERS_1996.getGMSTFunction(TimeScalesFactory.getUT1(IERSConventions.IERS_1996, true));
633         final TimeScalarFunction gmst00 =
634                 IERSConventions.IERS_2003.getGMSTFunction(TimeScalesFactory.getUT1(IERSConventions.IERS_2003, true));
635         final DSFactory factory = new DSFactory(1, 1);
636         final FieldAbsoluteDate<DerivativeStructure> ds2000 =
637                         FieldAbsoluteDate.getJ2000Epoch(factory.getDerivativeField());
638         for (double dt = 0; dt < 10 * Constants.JULIAN_YEAR; dt += 10 * Constants.JULIAN_DAY) {
639             FieldAbsoluteDate<DerivativeStructure> date = ds2000.shiftedBy(factory.variable(0, dt));
640             DerivativeStructure delta = gmst00.value(date).subtract(gmst82.value(date));
641             // GMST82 and GMST2000 are similar to about 15 milli-arcseconds between 2000 and 2010
642             Assert.assertEquals(0.0, MathUtils.normalizeAngle(delta.getValue(), 0.0), 7.1e-8);
643             Assert.assertEquals(0.0, delta.getPartialDerivative(1), 1.0e-15);
644         }
645     }
646 
647     @Test
648     public void testGMST2000vs2006() throws OrekitException {
649 
650         final UT1Scale ut1 = TimeScalesFactory.getUT1(IERSConventions.IERS_2010, true);
651         final TimeScalarFunction gmst00 = IERSConventions.IERS_2003.getGMSTFunction(ut1);
652         final TimeScalarFunction gmst06 = IERSConventions.IERS_2010.getGMSTFunction(ut1);
653         final DSFactory factory = new DSFactory(1, 1);
654         final FieldAbsoluteDate<DerivativeStructure> ds2000 =
655                         FieldAbsoluteDate.getJ2000Epoch(factory.getDerivativeField());
656         for (double dt = 0; dt < 10 * Constants.JULIAN_YEAR; dt += 10 * Constants.JULIAN_DAY) {
657             FieldAbsoluteDate<DerivativeStructure> date = ds2000.shiftedBy(factory.variable(0, dt));
658             DerivativeStructure delta = gmst06.value(date).subtract(gmst00.value(date));
659             // GMST2006 and GMST2000 are similar to about 0.15 milli-arcseconds between 2000 and 2010
660             Assert.assertEquals(0.0, MathUtils.normalizeAngle(delta.getValue(), 0.0), 7e-10);
661             Assert.assertEquals(0.0, delta.getPartialDerivative(1), 3.0e-18);
662         }
663     }
664 
665     @Test
666     public void testIAU1994ResolutionC7Discontinuity() throws OrekitException {
667         TimeVectorFunction nutation = IERSConventions.IERS_1996.getNutationFunction();
668         AbsoluteDate switchDate = new AbsoluteDate(1997, 2, 27, TimeScalesFactory.getUTC());
669         double h = 0.01;
670         for (double dt = -1.0 - h / 2; dt <= 1.0 + h /2; dt += h) {
671             AbsoluteDate d = switchDate.shiftedBy(dt);
672             final double currentCorr = nutation.value(d)[2];
673                 if (dt < 0) {
674                     Assert.assertEquals(0.0, currentCorr, 1.0e-20);
675                 } else {
676                     Assert.assertEquals(-7.87098e-12, currentCorr, 1.0e-15);
677                 }
678         }
679     }
680 
681     @Test
682     public void testTidalCorrection1996() throws OrekitException {
683 
684         // the reference value has been computed using interp_old.f from
685         // ftp://hpiers.obspm.fr/iers/models/ with the following driver program:
686         //
687         //        PROGRAM OR_TEST
688         //        IMPLICIT NONE
689         //        INTEGER J, N
690         //        PARAMETER (N = 4)
691         //        INTEGER MJD(N)
692         //        DOUBLE PRECISION TTTAI, TAIUTC
693         //        PARAMETER (TAIUTC = 32.000D0)
694         //        PARAMETER (TTTAI  = 32.184D0)
695         //        DOUBLE PRECISION RJD(N), X(N), Y(N), T(N), H(5)
696         //        DOUBLE PRECISION RJDINT, XINT, YINT, TINT, CORX, CORY, CORT
697         //        DATA(MJD(J),    T(J),        X(J),      Y(J),
698         //       &     J=1,4)/
699         //       &    52653, -0.2979055D0, -0.120344D0, 0.217095D0,
700         //       &    52654, -0.2984238D0, -0.121680D0, 0.219400D0,
701         //       &    52655, -0.2987682D0, -0.122915D0, 0.221760D0,
702         //       &    52656, -0.2989957D0, -0.124248D0, 0.224294D0/
703         //  C
704         //        DATA(H(J),J=1,5)/0.0D0, 3600.0D0, 7200.0D0, 43200.0D0, 86400.0D0/
705         //  C
706         //        DO 10 J = 1, N
707         //            RJD(J) = MJD(J) + (TTTAI + TAIUTC) / 86400.0D0
708         //   10   CONTINUE
709         //  C
710         //        DO 20 J = 1, 5
711         //            RJDINT = RJD(2) + H(J) / 86400.0D0
712         //            CALL INTERP(RJD,X,Y,T,N,RJDINT,XINT,YINT,TINT)
713         //            WRITE(6, 30) H(J), TINT, XINT, YINT
714         //   20   CONTINUE
715         //   30   FORMAT(F7.1,3(1X, F20.17))
716         //  C
717         //        END PROGRAM
718         //
719         // the output of this test reads:
720         //           0.0 -0.29839889612705234 -0.12191552977388567  0.21921143351558756
721         //        3600.0 -0.29841696994736727 -0.12208276003864491  0.21925416583607277
722         //        7200.0 -0.29843402412052122 -0.12218102086455683  0.21930263880545320
723         //       43200.0 -0.29866785146390035 -0.12250027826538630  0.22103779809979979
724         //       86400.0 -0.29874248853173840 -0.12308592577174847  0.22161565557764881
725 
726         Utils.setLoaders(IERSConventions.IERS_1996,
727                          Utils.buildEOPList(IERSConventions.IERS_1996, ITRFVersion.ITRF_2008, new double[][] {
728                              {  52653,  -0.2979055,   0.0005744,  -0.120344,   0.217095, 0.0, 0.0, 0.0, 0.0 },
729                              {  52654,  -0.2984238,   0.0004224,  -0.121680,   0.219400, 0.0, 0.0, 0.0, 0.0 },
730                              {  52655,  -0.2987682,   0.0002878,  -0.122915,   0.221760, 0.0, 0.0, 0.0, 0.0 },
731                              {  52656,  -0.2989957,   0.0001778,  -0.124248,   0.224294, 0.0, 0.0, 0.0, 0.0 }
732                          }));
733         EOPHistory eopHistory = FramesFactory.getEOPHistory(IERSConventions.IERS_1996, false);
734 
735         final AbsoluteDate t0 = new AbsoluteDate(new DateComponents(DateComponents.MODIFIED_JULIAN_EPOCH, 52654),
736                                            TimeScalesFactory.getUTC());
737         for (double[] row : new double[][] {
738             {     0.0, -0.29839889612705234, -0.12191552977388567,  0.21921143351558756 },
739             {  3600.0, -0.29841696994736727, -0.12208276003864491,  0.21925416583607277 },
740             {  7200.0, -0.29843402412052122, -0.12218102086455683,  0.21930263880545320 },
741             { 43200.0, -0.29866785146390035, -0.12250027826538630,  0.22103779809979979 },
742             { 86400.0, -0.29874248853173840, -0.12308592577174847,  0.22161565557764881 }
743         }) {
744             AbsoluteDate date = t0.shiftedBy(row[0]);
745             Assert.assertEquals(row[1], eopHistory.getUT1MinusUTC(date), 8.8e-11);
746             Assert.assertEquals(row[2] * Constants.ARC_SECONDS_TO_RADIANS,
747                                 eopHistory.getPoleCorrection(date).getXp(),
748                                 3.2e-14);
749             Assert.assertEquals(row[3] * Constants.ARC_SECONDS_TO_RADIANS,
750                                 eopHistory.getPoleCorrection(date).getYp(),
751                                 8.2e-15);
752         }
753 
754     }
755 
756     @Test
757     public void testTidalCorrection2003() throws OrekitException {
758 
759         // the reference value has been computed using the September 2007 version
760         // of interp.f from ftp://hpiers.obspm.fr/iers/models/ adding input/output
761         // parameters for LOD (which was already computed in the underlying routines),
762         // with the following driver program:
763         //
764         //        PROGRAM OR_TEST
765         //        IMPLICIT NONE
766         //        INTEGER J, N
767         //        PARAMETER (N = 4)
768         //        INTEGER MJD(N)
769         //        DOUBLE PRECISION TTTAI, TAIUTC
770         //        PARAMETER (TAIUTC = 32.000D0)
771         //        PARAMETER (TTTAI  = 32.184D0)
772         //        DOUBLE PRECISION RJD(N), UT1(N), LOD(N), X(N), Y(N), H(5)
773         //        DOUBLE PRECISION RJD_INT, UT1_INT, LOD_INT, X_INT, Y_INT
774         //        DATA(MJD(J),   UT1(J),     LOD(J),      X(J),      Y(J),
775         //       &     J=1,4)/
776         //       &    52653, -0.2979055D0, 0.0005744D0, -0.120344D0, 0.217095D0,
777         //       &    52654, -0.2984238D0, 0.0004224D0, -0.121680D0, 0.219400D0,
778         //       &    52655, -0.2987682D0, 0.0002878D0, -0.122915D0, 0.221760D0,
779         //       &    52656, -0.2989957D0, 0.0001778D0, -0.124248D0, 0.224294D0/
780         //  C
781         //        DATA(H(J),J=1,5)/0.0D0, 3600.0D0, 7200.0D0, 43200.0D0, 86400.0D0/
782         //  C
783         //        DO 10 J = 1, N
784         //            RJD(J) = MJD(J) + (TTTAI + TAIUTC) / 86400.0D0
785         //   10   CONTINUE
786         //  C
787         //        DO 20 J = 1, 5
788         //            RJD_INT = RJD(2) + H(J) / 86400.0D0
789         //            CALL INTERP(RJD,X,Y,UT1,LOD,N,
790         //       &                RJD_INT,X_INT,Y_INT,UT1_INT,LOD_INT)
791         //            WRITE(6, 30) H(J),UT1_INT,LOD_INT,X_INT,Y_INT
792         //   20   CONTINUE
793         //   30   FORMAT(F7.1,4(1X, F20.17))
794         //  C
795         //        END PROGRAM
796         //
797         // the output of this test reads:
798         //           0.0 -0.29840026968370659  0.00045312852893139 -0.12196223480123573  0.21922730818562719
799         //        3600.0 -0.29841834564816189  0.00041710864863793 -0.12213345007640604  0.21927433626001305
800         //        7200.0 -0.29843503870494986  0.00039207574457087 -0.12222881007999241  0.21932415788122142
801         //       43200.0 -0.29866930257052676  0.00042895046506082 -0.12247697694276605  0.22105450666130921
802         //       86400.0 -0.29874235341010519  0.00035460263868306 -0.12312252389660779  0.22161364352515728
803 
804         Utils.setLoaders(IERSConventions.IERS_2003,
805                          Utils.buildEOPList(IERSConventions.IERS_2003, ITRFVersion.ITRF_2008, new double[][] {
806                              {  52653,  -0.2979055,   0.0005744,  -0.120344,   0.217095, 0.0, 0.0, 0.0, 0.0 },
807                              {  52654,  -0.2984238,   0.0004224,  -0.121680,   0.219400, 0.0, 0.0, 0.0, 0.0 },
808                              {  52655,  -0.2987682,   0.0002878,  -0.122915,   0.221760, 0.0, 0.0, 0.0, 0.0 },
809                              {  52656,  -0.2989957,   0.0001778,  -0.124248,   0.224294, 0.0, 0.0, 0.0, 0.0 }
810                          }));
811         EOPHistory eopHistory = FramesFactory.getEOPHistory(IERSConventions.IERS_2003, false);
812 
813         final AbsoluteDate t0 = new AbsoluteDate(new DateComponents(DateComponents.MODIFIED_JULIAN_EPOCH, 52654),
814                                                  TimeScalesFactory.getUTC());
815 
816         for (double[] row : new double[][] {
817             {     0.0, -0.29840026968370659,  0.00045312852893139, -0.12196223480123573,  0.21922730818562719 },
818             {  3600.0, -0.29841834564816189,  0.00041710864863793, -0.12213345007640604,  0.21927433626001305 },
819             {  7200.0, -0.29843503870494986,  0.00039207574457087, -0.12222881007999241,  0.21932415788122142 },
820             { 43200.0, -0.29866930257052676,  0.00042895046506082, -0.12247697694276605,  0.22105450666130921 },
821             { 86400.0, -0.29874235341010519,  0.00035460263868306, -0.12312252389660779,  0.22161364352515728 }
822         }) {
823             AbsoluteDate date = t0.shiftedBy(row[0]);
824             Assert.assertEquals(row[1], eopHistory.getUT1MinusUTC(date), 4.3e-8);
825             Assert.assertEquals(row[2], eopHistory.getLOD(date), 1.4e-7);
826             Assert.assertEquals(row[3] * Constants.ARC_SECONDS_TO_RADIANS,
827                                 eopHistory.getPoleCorrection(date).getXp(),
828                                 1.6e-10);
829             Assert.assertEquals(row[4] * Constants.ARC_SECONDS_TO_RADIANS,
830                                 eopHistory.getPoleCorrection(date).getYp(),
831                                 0.7e-10);
832         }
833 
834     }
835 
836     @Test
837     public void testTidalCorrection2010() throws OrekitException {
838 
839         // the reference value has been computed using the September 2007 version
840         // of interp.f from ftp://hpiers.obspm.fr/iers/models/ adding input/output
841         // parameters for LOD (which was already computed in the underlying routines),
842         // with the following driver program:
843         //
844         //        PROGRAM OR_TEST
845         //        IMPLICIT NONE
846         //        INTEGER J, N
847         //        PARAMETER (N = 4)
848         //        INTEGER MJD(N)
849         //        DOUBLE PRECISION TTTAI, TAIUTC
850         //        PARAMETER (TAIUTC = 32.000D0)
851         //        PARAMETER (TTTAI  = 32.184D0)
852         //        DOUBLE PRECISION RJD(N), UT1(N), LOD(N), X(N), Y(N), H(5)
853         //        DOUBLE PRECISION RJD_INT, UT1_INT, LOD_INT, X_INT, Y_INT
854         //        DATA(MJD(J),   UT1(J),     LOD(J),      X(J),      Y(J),
855         //       &     J=1,4)/
856         //       &    52653, -0.2979055D0, 0.0005744D0, -0.120344D0, 0.217095D0,
857         //       &    52654, -0.2984238D0, 0.0004224D0, -0.121680D0, 0.219400D0,
858         //       &    52655, -0.2987682D0, 0.0002878D0, -0.122915D0, 0.221760D0,
859         //       &    52656, -0.2989957D0, 0.0001778D0, -0.124248D0, 0.224294D0/
860         //  C
861         //        DATA(H(J),J=1,5)/0.0D0, 3600.0D0, 7200.0D0, 43200.0D0, 86400.0D0/
862         //  C
863         //        DO 10 J = 1, N
864         //            RJD(J) = MJD(J) + (TTTAI + TAIUTC) / 86400.0D0
865         //   10   CONTINUE
866         //  C
867         //        DO 20 J = 1, 5
868         //            RJD_INT = RJD(2) + H(J) / 86400.0D0
869         //            CALL INTERP(RJD,X,Y,UT1,LOD,N,
870         //       &                RJD_INT,X_INT,Y_INT,UT1_INT,LOD_INT)
871         //            WRITE(6, 30) H(J),UT1_INT,LOD_INT,X_INT,Y_INT
872         //   20   CONTINUE
873         //   30   FORMAT(F7.1,4(1X, F20.17))
874         //  C
875         //        END PROGRAM
876         //
877         // the output of this test reads:
878         //           0.0 -0.29840026968370659  0.00045312852893139 -0.12196223480123573  0.21922730818562719
879         //        3600.0 -0.29841834564816189  0.00041710864863793 -0.12213345007640604  0.21927433626001305
880         //        7200.0 -0.29843503870494986  0.00039207574457087 -0.12222881007999241  0.21932415788122142
881         //       43200.0 -0.29866930257052676  0.00042895046506082 -0.12247697694276605  0.22105450666130921
882         //       86400.0 -0.29874235341010519  0.00035460263868306 -0.12312252389660779  0.22161364352515728
883 
884         Utils.setLoaders(IERSConventions.IERS_2010,
885                          Utils.buildEOPList(IERSConventions.IERS_2010, ITRFVersion.ITRF_2008, new double[][] {
886                              {  52653,  -0.2979055,   0.0005744,  -0.120344,   0.217095, 0.0, 0.0, 0.0, 0.0 },
887                              {  52654,  -0.2984238,   0.0004224,  -0.121680,   0.219400, 0.0, 0.0, 0.0, 0.0 },
888                              {  52655,  -0.2987682,   0.0002878,  -0.122915,   0.221760, 0.0, 0.0, 0.0, 0.0 },
889                              {  52656,  -0.2989957,   0.0001778,  -0.124248,   0.224294, 0.0, 0.0, 0.0, 0.0 }
890                          }));
891         EOPHistory eopHistory = FramesFactory.getEOPHistory(IERSConventions.IERS_2010, false);
892 
893         final AbsoluteDate t0 = new AbsoluteDate(new DateComponents(DateComponents.MODIFIED_JULIAN_EPOCH, 52654),
894                                                  TimeScalesFactory.getUTC());
895 
896         for (double[] row : new double[][] {
897             {     0.0, -0.29840026968370659,  0.00045312852893139, -0.12196223480123573,  0.21922730818562719 },
898             {  3600.0, -0.29841834564816189,  0.00041710864863793, -0.12213345007640604,  0.21927433626001305 },
899             {  7200.0, -0.29843503870494986,  0.00039207574457087, -0.12222881007999241,  0.21932415788122142 },
900             { 43200.0, -0.29866930257052676,  0.00042895046506082, -0.12247697694276605,  0.22105450666130921 },
901             { 86400.0, -0.29874235341010519,  0.00035460263868306, -0.12312252389660779,  0.22161364352515728 }
902         }) {
903             AbsoluteDate date = t0.shiftedBy(row[0]);
904             Assert.assertEquals(row[1], eopHistory.getUT1MinusUTC(date), 2.5e-12);
905             Assert.assertEquals(row[2], eopHistory.getLOD(date), 4.3e-11);
906             Assert.assertEquals(row[3] * Constants.ARC_SECONDS_TO_RADIANS,
907                                 eopHistory.getPoleCorrection(date).getXp(),
908                                 1.6e-10);
909             Assert.assertEquals(row[4] * Constants.ARC_SECONDS_TO_RADIANS,
910                                 eopHistory.getPoleCorrection(date).getYp(),
911                                 0.7e-10);
912         }
913 
914     }
915 
916     @Test
917     public void testGMSTRate1996() throws OrekitException {
918         final UT1Scale ut1 = TimeScalesFactory.getUT1(IERSConventions.IERS_1996, true);
919         checkGMSTRate(IERSConventions.IERS_1996.getGMSTFunction(ut1),
920                       IERSConventions.IERS_1996.getGMSTRateFunction(ut1),
921                       AbsoluteDate.J2000_EPOCH.shiftedBy(-0.4 * Constants.JULIAN_DAY),
922                       0.8 * Constants.JULIAN_DAY, 600.0, 10.0, 8.0e-13);
923     }
924 
925     @Test
926     public void testGMSTRate2003() throws OrekitException {
927         final UT1Scale ut1 = TimeScalesFactory.getUT1(IERSConventions.IERS_2003, true);
928         checkGMSTRate(IERSConventions.IERS_2003.getGMSTFunction(ut1),
929                       IERSConventions.IERS_2003.getGMSTRateFunction(ut1),
930                       AbsoluteDate.J2000_EPOCH.shiftedBy(-0.4 * Constants.JULIAN_DAY),
931                       0.8 * Constants.JULIAN_DAY, 600.0, 10.0, 8.0e-13);
932     }
933 
934     @Test
935     public void testGMSTRate2010() throws OrekitException {
936         final UT1Scale ut1 = TimeScalesFactory.getUT1(IERSConventions.IERS_2010, true);
937         checkGMSTRate(IERSConventions.IERS_2010.getGMSTFunction(ut1),
938                       IERSConventions.IERS_2010.getGMSTRateFunction(ut1),
939                       AbsoluteDate.J2000_EPOCH.shiftedBy(-0.4 * Constants.JULIAN_DAY),
940                       0.8 * Constants.JULIAN_DAY, 600.0, 10.0, 8.0e-13);
941     }
942 
943     private void checkGMSTRate(final TimeScalarFunction gmst, final TimeScalarFunction gmstRate,
944                                final AbsoluteDate date, final double span, final double sampleStep,
945                                final double h, final double tolerance)
946         throws OrekitException {
947         UnivariateDifferentiableFunction differentiated =
948                         new FiniteDifferencesDifferentiator(4, h).differentiate(new UnivariateFunction() {
949                             @Override
950                             public double value(final double dt) {
951                                 return gmst.value(date.shiftedBy(dt));
952                             }
953                         });
954         DSFactory factory = new DSFactory(1, 1);
955         double maxRateError = 0;
956         for (double dt = 0; dt < span; dt += sampleStep) {
957             double rateRef = differentiated.value(factory.variable(0, dt)).getPartialDerivative(1);
958             double rate    = gmstRate.value(date.shiftedBy(dt));
959             maxRateError   = FastMath.max(maxRateError, FastMath.abs(rateRef - rate));
960         }
961         Assert.assertEquals(0, maxRateError, tolerance);
962     }
963 
964     @Test
965     public void testFieldGMSTRate1996() throws OrekitException {
966         final UT1Scale ut1 = TimeScalesFactory.getUT1(IERSConventions.IERS_1996, true);
967         checkFieldGMSTRate(IERSConventions.IERS_1996.getGMSTFunction(ut1),
968                            IERSConventions.IERS_1996.getGMSTRateFunction(ut1),
969                            AbsoluteDate.J2000_EPOCH.shiftedBy(-0.4 * Constants.JULIAN_DAY),
970                            0.8 * Constants.JULIAN_DAY, 600.0, 10.0, Double.MIN_VALUE);
971     }
972 
973     @Test
974     public void testFieldGMSTRate2003() throws OrekitException {
975         final UT1Scale ut1 = TimeScalesFactory.getUT1(IERSConventions.IERS_2003, true);
976         checkFieldGMSTRate(IERSConventions.IERS_2003.getGMSTFunction(ut1),
977                            IERSConventions.IERS_2003.getGMSTRateFunction(ut1),
978                            AbsoluteDate.J2000_EPOCH.shiftedBy(-0.4 * Constants.JULIAN_DAY),
979                            0.8 * Constants.JULIAN_DAY, 600.0, 10.0, Double.MIN_VALUE);
980     }
981 
982     @Test
983     public void testFieldGMSTRate2010() throws OrekitException {
984         final UT1Scale ut1 = TimeScalesFactory.getUT1(IERSConventions.IERS_2010, true);
985         checkFieldGMSTRate(IERSConventions.IERS_2010.getGMSTFunction(ut1),
986                            IERSConventions.IERS_2010.getGMSTRateFunction(ut1),
987                            AbsoluteDate.J2000_EPOCH.shiftedBy(-0.4 * Constants.JULIAN_DAY),
988                            0.8 * Constants.JULIAN_DAY, 600.0, 10.0, Double.MIN_VALUE);
989     }
990 
991     private void checkFieldGMSTRate(final TimeScalarFunction gmst, final TimeScalarFunction gmstRate,
992                                     final AbsoluteDate date, final double span, final double sampleStep,
993                                     final double h, final double tolerance)
994         throws OrekitException {
995         double maxError = 0;
996         for (double dt = 0; dt < span; dt += sampleStep) {
997             double rateDouble    = gmstRate.value(date.shiftedBy(dt));
998             double rateDecimal64 = gmstRate.value(new FieldAbsoluteDate<>(date, new Decimal64(dt))).doubleValue();
999             maxError = FastMath.max(maxError, FastMath.abs(rateDouble - rateDecimal64));
1000         }
1001         Assert.assertEquals(0, maxError, tolerance);
1002     }
1003 
1004     @Test
1005     public void testPrecessionFunction1996FieldConsistency() throws OrekitException {
1006         checkVectorFunctionConsistency(IERSConventions.IERS_1996.getPrecessionFunction(), 3,
1007                                        AbsoluteDate.J2000_EPOCH.shiftedBy(-0.4 * Constants.JULIAN_DAY),
1008                                        0.8 * Constants.JULIAN_DAY, 600.0, 10.0, 2.0e-22, 3.0e-19);
1009     }
1010 
1011     @Test
1012     public void testPrecessionFunction2003FieldConsistency() throws OrekitException {
1013         checkVectorFunctionConsistency(IERSConventions.IERS_2003.getPrecessionFunction(), 3,
1014                                        AbsoluteDate.J2000_EPOCH.shiftedBy(-0.4 * Constants.JULIAN_DAY),
1015                                        0.8 * Constants.JULIAN_DAY, 600.0, 10.0, 6.0e-17, 6.0e-18);
1016     }
1017 
1018     @Test
1019     public void testPrecessionFunction2010FieldConsistency() throws OrekitException {
1020         checkVectorFunctionConsistency(IERSConventions.IERS_2010.getPrecessionFunction(), 3,
1021                                        AbsoluteDate.J2000_EPOCH.shiftedBy(-0.4 * Constants.JULIAN_DAY),
1022                                        0.8 * Constants.JULIAN_DAY, 600.0, 10.0, 6.0e-17, 6.0e-18);
1023     }
1024 
1025     @Test
1026     public void testNutationFunction1996FieldConsistency() throws OrekitException {
1027         checkVectorFunctionConsistency(IERSConventions.IERS_1996.getNutationFunction(), 3,
1028                                        AbsoluteDate.J2000_EPOCH.shiftedBy(-0.4 * Constants.JULIAN_DAY),
1029                                        0.8 * Constants.JULIAN_DAY, 600.0, 10.0, 2.0e-19, 4.0e-21);
1030     }
1031 
1032     @Test
1033     public void testNutationFunction2003FieldConsistency() throws OrekitException {
1034         checkVectorFunctionConsistency(IERSConventions.IERS_2003.getNutationFunction(), 3,
1035                                        AbsoluteDate.J2000_EPOCH.shiftedBy(-0.4 * Constants.JULIAN_DAY),
1036                                        0.8 * Constants.JULIAN_DAY, 600.0, 10.0, 3.0e-19, 6.0e-21);
1037     }
1038 
1039     @Test
1040     public void testNutationFunction2010FieldConsistency() throws OrekitException {
1041         checkVectorFunctionConsistency(IERSConventions.IERS_2010.getNutationFunction(), 3,
1042                                        AbsoluteDate.J2000_EPOCH.shiftedBy(-0.4 * Constants.JULIAN_DAY),
1043                                        0.8 * Constants.JULIAN_DAY, 600.0, 10.0, 5.0e-19, 6.0e-21);
1044     }
1045 
1046     @Test
1047     public void testXYSpXY2Function1996FieldConsistency() throws OrekitException {
1048         checkVectorFunctionConsistency(IERSConventions.IERS_1996.getXYSpXY2Function(), 3,
1049                                        AbsoluteDate.J2000_EPOCH.shiftedBy(-0.4 * Constants.JULIAN_DAY),
1050                                        0.8 * Constants.JULIAN_DAY, 600.0, 10.0, 6.0e-20, 2.0e-21);
1051     }
1052 
1053     @Test
1054     public void testXYSpXY2Function2003FieldConsistency() throws OrekitException {
1055         checkVectorFunctionConsistency(IERSConventions.IERS_2003.getXYSpXY2Function(), 3,
1056                                        AbsoluteDate.J2000_EPOCH.shiftedBy(-0.4 * Constants.JULIAN_DAY),
1057                                        0.8 * Constants.JULIAN_DAY, 600.0, 10.0, 2.0e-19, 3.0e-21);
1058     }
1059 
1060     @Test
1061     public void testXYSpXY2Function2010FieldConsistency() throws OrekitException {
1062         checkVectorFunctionConsistency(IERSConventions.IERS_2010.getXYSpXY2Function(), 3,
1063                                        AbsoluteDate.J2000_EPOCH.shiftedBy(-0.4 * Constants.JULIAN_DAY),
1064                                        0.8 * Constants.JULIAN_DAY, 600.0, 10.0, 2.0e-19, 3.0e-21);
1065     }
1066 
1067     @Test
1068     public void testMeanObliquityFunction1996FieldConsistency() throws OrekitException {
1069         checkScalarFunctionConsistency(IERSConventions.IERS_1996.getMeanObliquityFunction(),
1070                                        AbsoluteDate.J2000_EPOCH.shiftedBy(-0.4 * Constants.JULIAN_DAY),
1071                                        0.8 * Constants.JULIAN_DAY, 600.0, 10.0, 6.0e-17, 6.0e-18);
1072     }
1073 
1074     @Test
1075     public void testMeanObliquityFunction2003FieldConsistency() throws OrekitException {
1076         checkScalarFunctionConsistency(IERSConventions.IERS_2003.getMeanObliquityFunction(),
1077                                        AbsoluteDate.J2000_EPOCH.shiftedBy(-0.4 * Constants.JULIAN_DAY),
1078                                        0.8 * Constants.JULIAN_DAY, 600.0, 10.0, 6.0e-17, 6.0e-18);
1079     }
1080 
1081     @Test
1082     public void testMeanObliquityFunction2010FieldConsistency() throws OrekitException {
1083         checkScalarFunctionConsistency(IERSConventions.IERS_2010.getMeanObliquityFunction(),
1084                                        AbsoluteDate.J2000_EPOCH.shiftedBy(-0.4 * Constants.JULIAN_DAY),
1085                                        0.8 * Constants.JULIAN_DAY, 600.0, 10.0, 6.0e-17, 6.0e-18);
1086     }
1087 
1088     @Test
1089     public void testEOPTidalCorrection1996FieldConsistency() throws OrekitException {
1090         checkVectorFunctionConsistency(IERSConventions.IERS_1996.getEOPTidalCorrection(), 4,
1091                                        AbsoluteDate.J2000_EPOCH.shiftedBy(-0.4 * Constants.JULIAN_DAY),
1092                                        0.8 * Constants.JULIAN_DAY, 600.0, 10.0, 2.0e-17, 6.0e-20);
1093     }
1094 
1095     @Test
1096     public void testEOPTidalCorrection2003FieldConsistency() throws OrekitException {
1097         checkVectorFunctionConsistency(IERSConventions.IERS_2003.getEOPTidalCorrection(), 4,
1098                                        AbsoluteDate.J2000_EPOCH.shiftedBy(-0.4 * Constants.JULIAN_DAY),
1099                                        0.8 * Constants.JULIAN_DAY, 600.0, 10.0, 2.0e-17, 2.0e-19);
1100     }
1101 
1102     @Test
1103     public void testEOPTidalCorrection2010FieldConsistency() throws OrekitException {
1104         checkVectorFunctionConsistency(IERSConventions.IERS_2010.getEOPTidalCorrection(), 4,
1105                                        AbsoluteDate.J2000_EPOCH.shiftedBy(-0.4 * Constants.JULIAN_DAY),
1106                                        0.8 * Constants.JULIAN_DAY, 600.0, 10.0, 2.0e-17, 2.0e-19);
1107     }
1108 
1109     @Test
1110     public void testTideFrequencyDependenceFunction1996FieldConsistency() throws OrekitException {
1111         final UT1Scale ut1 = TimeScalesFactory.getUT1(IERSConventions.IERS_1996, true);
1112         checkVectorFunctionConsistency(IERSConventions.IERS_1996.getTideFrequencyDependenceFunction(ut1), 5,
1113                                        AbsoluteDate.J2000_EPOCH.shiftedBy(-0.4 * Constants.JULIAN_DAY),
1114                                        0.8 * Constants.JULIAN_DAY, 600.0, 10.0, 4.0e-24, 4.0e-22);
1115     }
1116 
1117     @Test
1118     public void testTideFrequencyDependenceFunction2003FieldConsistency() throws OrekitException {
1119         final UT1Scale ut1 = TimeScalesFactory.getUT1(IERSConventions.IERS_2003, true);
1120         checkVectorFunctionConsistency(IERSConventions.IERS_2003.getTideFrequencyDependenceFunction(ut1), 5,
1121                                        AbsoluteDate.J2000_EPOCH.shiftedBy(-0.4 * Constants.JULIAN_DAY),
1122                                        0.8 * Constants.JULIAN_DAY, 600.0, 10.0, 4.0e-24, 4.0e-22);
1123     }
1124 
1125     @Test
1126     public void testTideFrequencyDependenceFunction2010FieldConsistency() throws OrekitException {
1127         final UT1Scale ut1 = TimeScalesFactory.getUT1(IERSConventions.IERS_2010, true);
1128         checkVectorFunctionConsistency(IERSConventions.IERS_2010.getTideFrequencyDependenceFunction(ut1), 5,
1129                                        AbsoluteDate.J2000_EPOCH.shiftedBy(-0.4 * Constants.JULIAN_DAY),
1130                                        0.8 * Constants.JULIAN_DAY, 600.0, 10.0, 4.0e-24, 4.0e-22);
1131     }
1132 
1133     @Test
1134     public void testSolidPoleTide1996FieldConsistency() throws OrekitException {
1135         EOPHistory eopHistory = FramesFactory.getEOPHistory(IERSConventions.IERS_1996, false);
1136         checkVectorFunctionConsistency(IERSConventions.IERS_1996.getSolidPoleTide(eopHistory), 2,
1137                                        AbsoluteDate.J2000_EPOCH.shiftedBy(-0.4 * Constants.JULIAN_DAY),
1138                                        0.8 * Constants.JULIAN_DAY, 600.0, 10.0, 3.0e-25, 4.0e-26);
1139     }
1140 
1141     @Test
1142     public void testSolidPoleTide2003FieldConsistency() throws OrekitException {
1143         EOPHistory eopHistory = FramesFactory.getEOPHistory(IERSConventions.IERS_2003, false);
1144         checkVectorFunctionConsistency(IERSConventions.IERS_2003.getSolidPoleTide(eopHistory), 2,
1145                                        AbsoluteDate.J2000_EPOCH.shiftedBy(-0.4 * Constants.JULIAN_DAY),
1146                                        0.8 * Constants.JULIAN_DAY, 600.0, 10.0, 2.0e-25, 3.0e-26);
1147     }
1148 
1149     @Test
1150     public void testSolidPoleTide2010FieldConsistency() throws OrekitException {
1151         EOPHistory eopHistory = FramesFactory.getEOPHistory(IERSConventions.IERS_2010, false);
1152         checkVectorFunctionConsistency(IERSConventions.IERS_2010.getSolidPoleTide(eopHistory), 2,
1153                                        AbsoluteDate.J2000_EPOCH.shiftedBy(-0.4 * Constants.JULIAN_DAY),
1154                                        0.8 * Constants.JULIAN_DAY, 600.0, 10.0, 2.0e-25, 2.0e-26);
1155     }
1156 
1157     @Test
1158     public void testOceanPoleTide1996FieldConsistency() throws OrekitException {
1159         EOPHistory eopHistory = FramesFactory.getEOPHistory(IERSConventions.IERS_1996, false);
1160         checkVectorFunctionConsistency(IERSConventions.IERS_1996.getOceanPoleTide(eopHistory), 2,
1161                                        AbsoluteDate.J2000_EPOCH.shiftedBy(-0.4 * Constants.JULIAN_DAY),
1162                                        0.8 * Constants.JULIAN_DAY, 600.0, 10.0, Double.MIN_VALUE, Double.MIN_VALUE);
1163     }
1164 
1165     @Test
1166     public void testOceanPoleTide2003FieldConsistency() throws OrekitException {
1167         EOPHistory eopHistory = FramesFactory.getEOPHistory(IERSConventions.IERS_2003, false);
1168         checkVectorFunctionConsistency(IERSConventions.IERS_2003.getOceanPoleTide(eopHistory), 2,
1169                                        AbsoluteDate.J2000_EPOCH.shiftedBy(-0.4 * Constants.JULIAN_DAY),
1170                                        0.8 * Constants.JULIAN_DAY, 600.0, 10.0, Double.MIN_VALUE, Double.MIN_VALUE);
1171     }
1172 
1173     @Test
1174     public void testOceanPoleTide2010FieldConsistency() throws OrekitException {
1175         EOPHistory eopHistory = FramesFactory.getEOPHistory(IERSConventions.IERS_2010, false);
1176         checkVectorFunctionConsistency(IERSConventions.IERS_2010.getOceanPoleTide(eopHistory), 2,
1177                                        AbsoluteDate.J2000_EPOCH.shiftedBy(-0.4 * Constants.JULIAN_DAY),
1178                                        0.8 * Constants.JULIAN_DAY, 600.0, 10.0, 2.0e-26, 2.0e-27);
1179     }
1180 
1181     private void checkScalarFunctionConsistency(final TimeScalarFunction function,
1182                                                 final AbsoluteDate date, final double span, final double sampleStep,
1183                                                 final double h, final double valueTolerance, final double derivativeTolerance) {
1184 
1185         UnivariateDifferentiableFunction differentiated =
1186                 new FiniteDifferencesDifferentiator(4, h).differentiate(new UnivariateFunction() {
1187                     @Override
1188                     public double value(final double dt) {
1189                         return function.value(date.shiftedBy(dt));
1190                     }
1191                 });
1192 
1193         DSFactory factory = new DSFactory(1, 1);
1194         FieldAbsoluteDate<DerivativeStructure> dsDate = new FieldAbsoluteDate<>(date, factory.constant(0.0));
1195         double maxValueError = 0;
1196         double maxDerivativeError = 0;
1197         for (double dt = 0; dt < span; dt += sampleStep) {
1198             DerivativeStructure dsdt = factory.variable(0, dt);
1199             DerivativeStructure yRef = differentiated.value(dsdt);
1200             DerivativeStructure y    = function.value(dsDate.shiftedBy(dsdt));
1201             maxValueError      = FastMath.max(maxValueError,
1202                                               FastMath.abs(yRef.getValue() - y.getValue()));
1203             maxDerivativeError = FastMath.max(maxDerivativeError,
1204                                               FastMath.abs(yRef.getPartialDerivative(1) - y.getPartialDerivative(1)));
1205         }
1206         Assert.assertEquals(0, maxValueError,      valueTolerance);
1207         Assert.assertEquals(0, maxDerivativeError, derivativeTolerance);
1208 
1209     }
1210 
1211     private void checkVectorFunctionConsistency(final TimeVectorFunction function, final int dimension,
1212                                                 final AbsoluteDate date, final double span, final double sampleStep,
1213                                                 final double h, final double valueTolerance, final double derivativeTolerance) {
1214 
1215         UnivariateDifferentiableVectorFunction differentiated =
1216                 new FiniteDifferencesDifferentiator(4, h).differentiate(new UnivariateVectorFunction() {
1217                     @Override
1218                     public double[] value(final double dt) {
1219                         return function.value(date.shiftedBy(dt));
1220                     }
1221                 });
1222 
1223         DSFactory factory = new DSFactory(1, 1);
1224         FieldAbsoluteDate<DerivativeStructure> dsDate = new FieldAbsoluteDate<>(date, factory.constant(0.0));
1225         double maxValueError = 0;
1226         double maxDerivativeError = 0;
1227         for (double dt = 0; dt < span; dt += sampleStep) {
1228             DerivativeStructure dsdt = factory.variable(0, dt);
1229             DerivativeStructure[] yRef = differentiated.value(dsdt);
1230             DerivativeStructure[] y    = function.value(dsDate.shiftedBy(dsdt));
1231             Assert.assertEquals(dimension, yRef.length);
1232             Assert.assertEquals(dimension, y.length);
1233             for (int i = 0; i < dimension; ++i) {
1234                 maxValueError      = FastMath.max(maxValueError,
1235                                                   FastMath.abs(yRef[i].getValue() - y[i].getValue()));
1236                 maxDerivativeError = FastMath.max(maxDerivativeError,
1237                                                   FastMath.abs(yRef[i].getPartialDerivative(1) - y[i].getPartialDerivative(1)));
1238             }
1239         }
1240         Assert.assertEquals(0, maxValueError,      valueTolerance);
1241         Assert.assertEquals(0, maxDerivativeError, derivativeTolerance);
1242 
1243     }
1244 
1245     @Before
1246     public void setUp() {
1247         Utils.setDataRoot("compressed-data");
1248     }
1249 
1250 }