1   /* Copyright 2002-2025 CS GROUP
2    * Licensed to CS GROUP (CS) under one or more
3    * contributor license agreements.  See the NOTICE file distributed with
4    * this work for additional information regarding copyright ownership.
5    * CS licenses this file to You under the Apache License, Version 2.0
6    * (the "License"); you may not use this file except in compliance with
7    * the License.  You may obtain a copy of the License at
8    *
9    *   http://www.apache.org/licenses/LICENSE-2.0
10   *
11   * Unless required by applicable law or agreed to in writing, software
12   * distributed under the License is distributed on an "AS IS" BASIS,
13   * WITHOUT WARRANTIES OR CONDITIONS OF ANY KIND, either express or implied.
14   * See the License for the specific language governing permissions and
15   * limitations under the License.
16   */
17  package org.orekit.frames;
18  
19  import java.util.stream.Collectors;
20  
21  import org.hipparchus.util.FastMath;
22  import org.junit.jupiter.api.AfterEach;
23  import org.junit.jupiter.api.Assertions;
24  import org.junit.jupiter.api.BeforeEach;
25  import org.junit.jupiter.api.Test;
26  import org.orekit.Utils;
27  import org.orekit.time.AbsoluteDate;
28  import org.orekit.time.DateComponents;
29  import org.orekit.time.TimeScale;
30  import org.orekit.time.TimeScalesFactory;
31  import org.orekit.utils.Constants;
32  import org.orekit.utils.IERSConventions;
33  
34  
35  public class PredictedEOPHistoryTest {
36  
37      @Test
38      public void testExtensionDates() {
39  
40          // truncate EOP between 2018 and 2020
41          final int mjdLimit = new DateComponents(2021, 1, 1).getMJD();
42          final EOPHistory truncatedEOP = new EOPHistory(trueEOP.getConventions(),
43                                                         EOPHistory.DEFAULT_INTERPOLATION_DEGREE,
44                                                         trueEOP.getEntries().
45                                                                 stream().
46                                                                 filter(e -> e.getMjd() < mjdLimit).
47                                                                 collect(Collectors.toList()),
48                                                         trueEOP.isSimpleEop(),
49                                                         trueEOP.getTimeScales());
50  
51          // predicted history from truncated data
52          EOPHistory predicted = new PredictedEOPHistory(truncatedEOP,
53                                                         30 * Constants.JULIAN_DAY,
54                                                         new EOPFitter(SingleParameterFitter.createDefaultDut1FitterShortTermPrediction(),
55                                                                       SingleParameterFitter.createDefaultPoleFitterShortTermPrediction(),
56                                                                       SingleParameterFitter.createDefaultPoleFitterShortTermPrediction(),
57                                                                       SingleParameterFitter.createDefaultNutationFitterShortTermPrediction(),
58                                                                       SingleParameterFitter.createDefaultNutationFitterShortTermPrediction()));
59          Assertions.assertEquals(0.0,
60                                  new AbsoluteDate(2018, 1, 1, 0, 0, 0.0, utc).durationFrom(predicted.getStartDate()),
61                                  1.0e-10);
62          Assertions.assertEquals(0.0,
63                                  new AbsoluteDate(2021, 1, 30, 0, 0, 0.0, utc).durationFrom(predicted.getEndDate()),
64                                  1.0e-10);
65      }
66  
67      @Test
68      public void testCommonCoverage() {
69  
70          // truncate EOP between 2018 and 2021
71          final int mjdLimit = new DateComponents(2022, 1, 1).getMJD();
72          final EOPHistory truncatedEOP = new EOPHistory(trueEOP.getConventions(),
73                                                         EOPHistory.DEFAULT_INTERPOLATION_DEGREE,
74                                                         trueEOP.getEntries().
75                                                                 stream().
76                                                                 filter(e -> e.getMjd() < mjdLimit).
77                                                                 collect(Collectors.toList()),
78                                                         trueEOP.isSimpleEop(),
79                                                         trueEOP.getTimeScales());
80  
81          EOPHistory predicted = new PredictedEOPHistory(truncatedEOP,
82                                                         30 * Constants.JULIAN_DAY,
83                                                         new EOPFitter(SingleParameterFitter.createDefaultDut1FitterShortTermPrediction(),
84                                                                       SingleParameterFitter.createDefaultPoleFitterShortTermPrediction(),
85                                                                       SingleParameterFitter.createDefaultPoleFitterShortTermPrediction(),
86                                                                       SingleParameterFitter.createDefaultNutationFitterShortTermPrediction(),
87                                                                       SingleParameterFitter.createDefaultNutationFitterShortTermPrediction()));
88  
89          // check we get the same value as raw EOP (dropping the last interpolated day)
90          double maxErrorUT1  = 0;
91          double maxErrorLOD  = 0;
92          double maxErrorXp   = 0;
93          double maxErrorYp   = 0;
94          double maxErrordPsi = 0;
95          double maxErrordEps = 0;
96          double maxErrordx   = 0;
97          double maxErrordy   = 0;
98          for (double dt = 0;
99               dt < truncatedEOP.getEndDate().durationFrom(truncatedEOP.getStartDate()) - Constants.JULIAN_DAY;
100              dt += 20000.0) {
101             final AbsoluteDate   date      = truncatedEOP.getStartDate().shiftedBy(dt);
102             final PoleCorrection rawPC     = truncatedEOP.getPoleCorrection(date);
103             final PoleCorrection predPC    = predicted.getPoleCorrection(date);
104             final double[]       rawENC    = truncatedEOP.getEquinoxNutationCorrection(date);
105             final double[]       predENC   = predicted.getEquinoxNutationCorrection(date);
106             final double[]       rawNroNC  = truncatedEOP.getNonRotatinOriginNutationCorrection(date);
107             final double[]       predNroNC = predicted.getNonRotatinOriginNutationCorrection(date);
108             maxErrorUT1  = FastMath.max(maxErrorUT1,  FastMath.abs(truncatedEOP.getUT1MinusUTC(date) - predicted.getUT1MinusUTC(date)));
109             maxErrorLOD  = FastMath.max(maxErrorLOD,  FastMath.abs(truncatedEOP.getLOD(date)         - predicted.getLOD(date)));
110             maxErrorXp   = FastMath.max(maxErrorXp,   FastMath.abs(rawPC.getXp()            - predPC.getXp()));
111             maxErrorYp   = FastMath.max(maxErrorYp,   FastMath.abs(rawPC.getYp()            - predPC.getYp()));
112             maxErrordPsi = FastMath.max(maxErrordPsi, FastMath.abs(rawENC[0]                - predENC[0]));
113             maxErrordEps = FastMath.max(maxErrordEps, FastMath.abs(rawENC[1]                - predENC[1]));
114             maxErrordx   = FastMath.max(maxErrordx,   FastMath.abs(rawNroNC[0]              - predNroNC[0]));
115             maxErrordy   = FastMath.max(maxErrordy,   FastMath.abs(rawNroNC[1]              - predNroNC[1]));
116         }
117         Assertions.assertEquals(0.0, maxErrorUT1,  1.0e-15);
118         Assertions.assertEquals(0.0, maxErrorLOD,  1.0e-15);
119         Assertions.assertEquals(0.0, maxErrorXp,   1.0e-15);
120         Assertions.assertEquals(0.0, maxErrorYp,   1.0e-15);
121         Assertions.assertEquals(0.0, maxErrordPsi, 1.0e-15);
122         Assertions.assertEquals(0.0, maxErrordEps, 1.0e-15);
123         Assertions.assertEquals(0.0, maxErrordx,   1.0e-15);
124         Assertions.assertEquals(0.0, maxErrordy,   1.0e-15);
125     }
126 
127     @Test
128     @Deprecated
129     public void testDeprecated() {
130 
131         // truncate EOP between 2018 and 2021
132         final int mjdLimit = new DateComponents(2022, 1, 1).getMJD();
133         final EOPHistory truncatedEOP = new EOPHistory(trueEOP.getConventions(),
134                                                        EOPHistory.DEFAULT_INTERPOLATION_DEGREE,
135                                                        trueEOP.getEntries().
136                                                                stream().
137                                                                filter(e -> e.getMjd() < mjdLimit).
138                                                                collect(Collectors.toList()),
139                                                        trueEOP.isSimpleEop(),
140                                                        trueEOP.getTimeScales());
141 
142         EOPHistory predicted = new PredictedEOPHistory(truncatedEOP,
143                                                        30 * Constants.JULIAN_DAY,
144                                                        new EOPFitter(new SingleParameterFitter(Constants.JULIAN_DAY,
145                                                                                                1.0e-12, 3,
146                                                                                                SingleParameterFitter.SUN_PULSATION,
147                                                                                                2 * SingleParameterFitter.SUN_PULSATION,
148                                                                                                3 * SingleParameterFitter.SUN_PULSATION,
149                                                                                                SingleParameterFitter.MOON_DRACONIC_PULSATION,
150                                                                                                2 * SingleParameterFitter.MOON_DRACONIC_PULSATION,
151                                                                                                3 * SingleParameterFitter.MOON_DRACONIC_PULSATION),
152                                                                      SingleParameterFitter.createDefaultPoleFitterLongTermPrediction(),
153                                                                      SingleParameterFitter.createDefaultPoleFitterLongTermPrediction(),
154                                                                      SingleParameterFitter.createDefaultNutationFitterLongTermPrediction(),
155                                                                      SingleParameterFitter.createDefaultNutationFitterLongTermPrediction()));
156 
157         // check we get the same value as raw EOP (dropping the last interpolated day)
158         double maxErrorUT1  = 0;
159         for (double dt = Constants.JULIAN_DAY; dt < 10 * Constants.JULIAN_DAY; dt += 20000.0) {
160             final AbsoluteDate   date      = truncatedEOP.getEndDate().shiftedBy(dt);
161             maxErrorUT1  = FastMath.max(maxErrorUT1,  FastMath.abs(trueEOP.getUT1MinusUTC(date) - predicted.getUT1MinusUTC(date)));
162         }
163         Assertions.assertEquals(4.563, maxErrorUT1, 0.001);
164 
165     }
166 
167     @Test
168     public void testAccuracyShortTerm() {
169         doTestAccuracy(true, 0.148, 0.580, 1.528, 7.842, 316.165, 6077.182, 40759.430);
170     }
171 
172     @Test
173     public void testAccuracyLongTerm() {
174         doTestAccuracy(false, 1.518, 1.993, 2.298, 3.013, 7.398, 17.582, 27.296);
175     }
176 
177     private void doTestAccuracy(final boolean shortTerm,
178                                 final double expectedMax001, final double expectedMax003,
179                                 final double expectedMax005, final double expectedMax010,
180                                 final double expectedMax030, final double expectedMax060,
181                                 final double expectedMax090) {
182 
183         double maxError001Days = 0;
184         double maxError003Days = 0;
185         double maxError005Days = 0;
186         double maxError010Days = 0;
187         double maxError030Days = 0;
188         double maxError060Days = 0;
189         double maxError090Days = 0;
190         for (int d = 0; d < 365; d += 30) {
191             // set up a prediction
192             final DateComponents dc = new DateComponents(new DateComponents(2021, 1, 1), d);
193             final int mjdLimit = dc.getMJD();
194             final EOPHistory truncatedEOP = new EOPHistory(trueEOP.getConventions(),
195                                                            EOPHistory.DEFAULT_INTERPOLATION_DEGREE,
196                                                            trueEOP.getEntries().
197                                                            stream().
198                                                            filter(e -> e.getMjd() < mjdLimit).
199                                                            collect(Collectors.toList()),
200                                                            trueEOP.isSimpleEop(),
201                                                            trueEOP.getTimeScales());
202             final PredictedEOPHistory predictedEOP = shortTerm ?
203                                                new PredictedEOPHistory(truncatedEOP,
204                                                                        180 * Constants.JULIAN_DAY,
205                                                                        new EOPFitter(SingleParameterFitter.createDefaultDut1FitterShortTermPrediction() ,
206                                                                                      SingleParameterFitter.createDefaultPoleFitterShortTermPrediction(),
207                                                                                      SingleParameterFitter.createDefaultPoleFitterShortTermPrediction(),
208                                                                                      SingleParameterFitter.createDefaultNutationFitterShortTermPrediction(),
209                                                                                      SingleParameterFitter.createDefaultNutationFitterShortTermPrediction())) :
210                                                new PredictedEOPHistory(truncatedEOP,
211                                                                        180 * Constants.JULIAN_DAY,
212                                                                        new EOPFitter(SingleParameterFitter.createDefaultDut1FitterLongTermPrediction() ,
213                                                                                      SingleParameterFitter.createDefaultPoleFitterLongTermPrediction(),
214                                                                                      SingleParameterFitter.createDefaultPoleFitterLongTermPrediction(),
215                                                                                      SingleParameterFitter.createDefaultNutationFitterLongTermPrediction(),
216                                                                                      SingleParameterFitter.createDefaultNutationFitterLongTermPrediction()));
217 
218             // set up two itrf frames, one using true, one using predicted EOP
219             final Frame itrfTrue = FramesFactory.buildUncachedITRF(trueEOP, TimeScalesFactory.getUTC());
220             final Frame itrfPred = FramesFactory.buildUncachedITRF(predictedEOP, TimeScalesFactory.getUTC());
221 
222             final AbsoluteDate t0 = new AbsoluteDate(dc, trueEOP.getTimeScales().getUTC());
223             for (double dt = 0; dt < 180 * Constants.JULIAN_DAY; dt += 43200) {
224                 final AbsoluteDate date = t0.shiftedBy(dt);
225                 final Transform t = itrfTrue.getTransformTo(itrfPred, date);
226                 final double deltaP = t.getRotation().getAngle() * Constants.WGS84_EARTH_EQUATORIAL_RADIUS;
227                 if (dt <= Constants.JULIAN_DAY) {
228                     maxError001Days = FastMath.max(maxError001Days, deltaP);
229                 }
230                 if (dt <= 3 * Constants.JULIAN_DAY) {
231                     maxError003Days = FastMath.max(maxError003Days, deltaP);
232                 }
233                 if (dt <= 5 * Constants.JULIAN_DAY) {
234                     maxError005Days = FastMath.max(maxError005Days, deltaP);
235                 }
236                 if (dt <= 10 * Constants.JULIAN_DAY) {
237                     maxError010Days = FastMath.max(maxError010Days, deltaP);
238                 }
239                 if (dt <= 30 * Constants.JULIAN_DAY) {
240                     maxError030Days = FastMath.max(maxError030Days, deltaP);
241                 }
242                 if (dt <= 60 * Constants.JULIAN_DAY) {
243                     maxError060Days = FastMath.max(maxError060Days, deltaP);
244                 }
245                 if (dt <= 90 * Constants.JULIAN_DAY) {
246                     maxError090Days = FastMath.max(maxError090Days, deltaP);
247                 }
248             }
249 
250         }
251 
252         Assertions.assertEquals(expectedMax001, maxError001Days, 0.001);
253         Assertions.assertEquals(expectedMax003, maxError003Days, 0.001);
254         Assertions.assertEquals(expectedMax005, maxError005Days, 0.001);
255         Assertions.assertEquals(expectedMax010, maxError010Days, 0.001);
256         Assertions.assertEquals(expectedMax030, maxError030Days, 0.001);
257         Assertions.assertEquals(expectedMax060, maxError060Days, 0.001);
258         Assertions.assertEquals(expectedMax090, maxError090Days, 0.001);
259 
260     }
261 
262     @BeforeEach
263     public void setUp() {
264         Utils.setDataRoot("eop-prediction");
265         FramesFactory.clearEOPHistoryLoaders();
266         FramesFactory.addDefaultEOP2000HistoryLoaders("none", "none", "^eopc04_14_IAU2000\\.[0-9][0-9]\\.txt$",
267                                                       "none", "none", "none");
268         utc     = TimeScalesFactory.getUTC();
269         trueEOP = FramesFactory.getEOPHistory(IERSConventions.IERS_2010, true);
270 
271         // raw history between 2018 and 2022
272         Assertions.assertEquals(0.0,
273                                 new AbsoluteDate(2018, 1, 1, 0, 0, 0.0, utc).durationFrom(trueEOP.getStartDate()),
274                                 1.0e-10);
275         Assertions.assertEquals(0.0,
276                                 new AbsoluteDate(2022, 12, 31, 0, 0, 0.0, utc).durationFrom(trueEOP.getEndDate()),
277                                 1.0e-10);
278 
279     }
280 
281     @AfterEach
282     public void tearDown() {
283         utc          = null;
284         trueEOP      = null;
285     }
286 
287     TimeScale  utc;
288     EOPHistory trueEOP;
289 
290 }