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 org.hamcrest.MatcherAssert;
20  import org.hipparchus.CalculusFieldElement;
21  import org.hipparchus.geometry.euclidean.threed.FieldRotation;
22  import org.hipparchus.geometry.euclidean.threed.FieldVector3D;
23  import org.hipparchus.geometry.euclidean.threed.Rotation;
24  import org.hipparchus.geometry.euclidean.threed.RotationConvention;
25  import org.hipparchus.geometry.euclidean.threed.Vector3D;
26  import org.hipparchus.util.Binary64;
27  import org.hipparchus.util.Binary64Field;
28  import org.hipparchus.util.FastMath;
29  import org.junit.jupiter.api.Assertions;
30  import org.junit.jupiter.api.BeforeEach;
31  import org.junit.jupiter.api.Test;
32  import org.orekit.OrekitMatchers;
33  import org.orekit.Utils;
34  import org.orekit.data.DataContext;
35  import org.orekit.time.AbsoluteDate;
36  import org.orekit.time.DateComponents;
37  import org.orekit.time.FieldAbsoluteDate;
38  import org.orekit.time.TimeComponents;
39  import org.orekit.time.TimeScalesFactory;
40  import org.orekit.utils.FieldPVCoordinates;
41  import org.orekit.utils.IERSConventions;
42  import org.orekit.utils.PVCoordinates;
43  
44  
45  public class GTODProviderTest {
46  
47      @Test
48      public void testAASReferenceLEO() {
49  
50          // this reference test has been extracted from the following paper:
51          // Implementation Issues Surrounding the New IAU Reference Systems for Astrodynamics
52          // David A. Vallado, John H. Seago, P. Kenneth Seidelmann
53          // http://www.centerforspace.com/downloads/files/pubs/AAS-06-134.pdf
54          Utils.setLoaders(IERSConventions.IERS_1996,
55                           Utils.buildEOPList(IERSConventions.IERS_1996, ITRFVersion.ITRF_2008, new double[][] {
56                               { 53098, -0.4399619, 0.0015563, -0.140682, 0.333309, -0.052195, -0.003875, Double.NaN, Double.NaN },
57                               { 53099, -0.4399619, 0.0015563, -0.140682, 0.333309, -0.052195, -0.003875, Double.NaN, Double.NaN },
58                               { 53100, -0.4399619, 0.0015563, -0.140682, 0.333309, -0.052195, -0.003875, Double.NaN, Double.NaN },
59                               { 53101, -0.4399619, 0.0015563, -0.140682, 0.333309, -0.052195, -0.003875, Double.NaN, Double.NaN },
60                               { 53102, -0.4399619, 0.0015563, -0.140682, 0.333309, -0.052195, -0.003875, Double.NaN, Double.NaN },
61                               { 53103, -0.4399619, 0.0015563, -0.140682, 0.333309, -0.052195, -0.003875, Double.NaN, Double.NaN },
62                               { 53104, -0.4399619, 0.0015563, -0.140682, 0.333309, -0.052195, -0.003875, Double.NaN, Double.NaN },
63                               { 53105, -0.4399619, 0.0015563, -0.140682, 0.333309, -0.052195, -0.003875, Double.NaN, Double.NaN }
64                           }));
65          AbsoluteDate t0 = new AbsoluteDate(new DateComponents(2004, 4, 6),
66                                             new TimeComponents(7, 51, 28.386009),
67                                             TimeScalesFactory.getUTC());
68  
69          // PEF iau76
70          PVCoordinates pvPEF =
71             new PVCoordinates(new Vector3D(-1033475.0313, 7901305.5856, 6380344.5328),
72                               new Vector3D(-3225.632747, -2872.442511, 5531.931288));
73  
74          // it seems the induced effect of pole nutation correction δΔψ on the equation of the equinoxes
75          // was not taken into account in the reference paper, so we fix it here for the test
76          final double dDeltaPsi =
77                  FramesFactory.getEOPHistory(IERSConventions.IERS_1996, true).getEquinoxNutationCorrection(t0)[0];
78          final double epsilonA = IERSConventions.IERS_1996.getMeanObliquityFunction().value(t0);
79          final Transform fix =
80                  new Transform(t0, new Rotation(Vector3D.PLUS_K,
81                                                 dDeltaPsi * FastMath.cos(epsilonA),
82                                                 RotationConvention.FRAME_TRANSFORM));
83  
84          // TOD iau76
85          PVCoordinates pvTOD =
86              new PVCoordinates(new Vector3D(5094514.7804, 6127366.4612, 6380344.5328),
87                                new Vector3D(-4746.088567, 786.077222, 5531.931288));
88  
89          Transform t = FramesFactory.getTOD(IERSConventions.IERS_1996, true).
90                  getTransformTo(FramesFactory.getGTOD(IERSConventions.IERS_1996, true), t0);
91          checkPV(fix.transformPVCoordinates(pvPEF), t.transformPVCoordinates(pvTOD), 0.00942, 3.12e-5);
92          StaticTransform st = FramesFactory.getTOD(IERSConventions.IERS_1996, true).
93                  getStaticTransformTo(FramesFactory.getGTOD(IERSConventions.IERS_1996, true), t0);
94          MatcherAssert.assertThat(
95                  st.getTranslation(),
96                  OrekitMatchers.vectorCloseTo(t.getTranslation(), 0));
97          MatcherAssert.assertThat(
98                  Rotation.distance(st.getRotation(), t.getRotation()),
99                  OrekitMatchers.closeTo(0, 0));
100 
101         // if we forget to apply nutation corrections, results are much worse, which is expected
102         t = FramesFactory.getTOD(false).getTransformTo(FramesFactory.getGTOD(false), t0);
103         checkPV(fix.transformPVCoordinates(pvPEF), t.transformPVCoordinates(pvTOD), 257.49, 0.13955);
104 
105     }
106 
107     @Test
108     public void testAASReferenceGEO() {
109 
110         // this reference test has been extracted from the following paper:
111         // Implementation Issues Surrounding the New IAU Reference Systems for Astrodynamics
112         // David A. Vallado, John H. Seago, P. Kenneth Seidelmann
113         // http://www.centerforspace.com/downloads/files/pubs/AAS-06-134.pdf
114         Utils.setLoaders(IERSConventions.IERS_1996,
115                          Utils.buildEOPList(IERSConventions.IERS_1996, ITRFVersion.ITRF_2008, new double[][] {
116                              { 53153, -0.4709050,  0.0000000, -0.083853,  0.467217, -0.053614, -0.004494, Double.NaN, Double.NaN },
117                              { 53154, -0.4709050,  0.0000000, -0.083853,  0.467217, -0.053614, -0.004494, Double.NaN, Double.NaN },
118                              { 53155, -0.4709050,  0.0000000, -0.083853,  0.467217, -0.053614, -0.004494, Double.NaN, Double.NaN },
119                              { 53156, -0.4709050,  0.0000000, -0.083853,  0.467217, -0.053614, -0.004494, Double.NaN, Double.NaN },
120                              { 53157, -0.4709050,  0.0000000, -0.083853,  0.467217, -0.053614, -0.004494, Double.NaN, Double.NaN },
121                              { 53158, -0.4709050,  0.0000000, -0.083853,  0.467217, -0.053614, -0.004494, Double.NaN, Double.NaN },
122                              { 53159, -0.4709050,  0.0000000, -0.083853,  0.467217, -0.053614, -0.004494, Double.NaN, Double.NaN },
123                              { 53160, -0.4709050,  0.0000000, -0.083853,  0.467217, -0.053614, -0.004494, Double.NaN, Double.NaN }
124                          }));
125         AbsoluteDate t0 = new AbsoluteDate(new DateComponents(2004, 6, 1),
126                                            TimeComponents.H00,
127                                            TimeScalesFactory.getUTC());
128 
129         Transform t = FramesFactory.getTOD(IERSConventions.IERS_1996, true).
130                 getTransformTo(FramesFactory.getGTOD(IERSConventions.IERS_1996, true), t0);
131         // TOD iau76
132         PVCoordinates pvTOD =
133             new PVCoordinates(new Vector3D(-40577427.7501, -11500096.1306, 10293.2583),
134                               new Vector3D(837.552338, -2957.524176, -0.928772));
135 
136         // PEF iau76
137         PVCoordinates pvPEF =
138             new PVCoordinates(new Vector3D(24796919.2956, -34115870.9001, 10293.2583),
139                               new Vector3D(-0.979178, -1.476540, -0.928772));
140 
141         // it seems the induced effect of pole nutation correction δΔψ on the equation of the equinoxes
142         // was not taken into account in the reference paper, so we fix it here for the test
143         final double dDeltaPsi =
144                 FramesFactory.getEOPHistory(IERSConventions.IERS_1996, true).getEquinoxNutationCorrection(t0)[0];
145         final double epsilonA = IERSConventions.IERS_1996.getMeanObliquityFunction().value(t0);
146         final Transform fix =
147                 new Transform(t0, new Rotation(Vector3D.PLUS_K,
148                                                dDeltaPsi * FastMath.cos(epsilonA),
149                                                RotationConvention.FRAME_TRANSFORM));
150 
151         checkPV(fix.transformPVCoordinates(pvPEF), t.transformPVCoordinates(pvTOD), 0.0503, 3.59e-4);
152 
153         // if we forget to apply nutation corrections, results are much worse, which is expected
154         t = FramesFactory.getTOD(false).getTransformTo(FramesFactory.getGTOD(false), t0);
155         checkPV(fix.transformPVCoordinates(pvPEF), t.transformPVCoordinates(pvTOD), 1458.27, 3.847e-4);
156 
157     }
158 
159     @Test
160     public void testAASReferenceGEOField() {
161 
162         // this reference test has been extracted from the following paper:
163         // Implementation Issues Surrounding the New IAU Reference Systems for Astrodynamics
164         // David A. Vallado, John H. Seago, P. Kenneth Seidelmann
165         // http://www.centerforspace.com/downloads/files/pubs/AAS-06-134.pdf
166         Utils.setLoaders(IERSConventions.IERS_1996,
167                          Utils.buildEOPList(IERSConventions.IERS_1996, ITRFVersion.ITRF_2008, new double[][] {
168                              { 53153, -0.4709050,  0.0000000, -0.083853,  0.467217, -0.053614, -0.004494, Double.NaN, Double.NaN },
169                              { 53154, -0.4709050,  0.0000000, -0.083853,  0.467217, -0.053614, -0.004494, Double.NaN, Double.NaN },
170                              { 53155, -0.4709050,  0.0000000, -0.083853,  0.467217, -0.053614, -0.004494, Double.NaN, Double.NaN },
171                              { 53156, -0.4709050,  0.0000000, -0.083853,  0.467217, -0.053614, -0.004494, Double.NaN, Double.NaN },
172                              { 53157, -0.4709050,  0.0000000, -0.083853,  0.467217, -0.053614, -0.004494, Double.NaN, Double.NaN },
173                              { 53158, -0.4709050,  0.0000000, -0.083853,  0.467217, -0.053614, -0.004494, Double.NaN, Double.NaN },
174                              { 53159, -0.4709050,  0.0000000, -0.083853,  0.467217, -0.053614, -0.004494, Double.NaN, Double.NaN },
175                              { 53160, -0.4709050,  0.0000000, -0.083853,  0.467217, -0.053614, -0.004494, Double.NaN, Double.NaN }
176                          }));
177         FieldAbsoluteDate<Binary64> t0 = new FieldAbsoluteDate<>(Binary64Field.getInstance(),
178                                                                   new DateComponents(2004, 6, 1),
179                                                                   TimeComponents.H00,
180                                                                   TimeScalesFactory.getUTC());
181 
182         FieldTransform<Binary64> t = FramesFactory.getTOD(IERSConventions.IERS_1996, true).
183                 getTransformTo(FramesFactory.getGTOD(IERSConventions.IERS_1996, true), t0);
184         // TOD iau76
185         PVCoordinates pvTOD =
186             new PVCoordinates(new Vector3D(-40577427.7501, -11500096.1306, 10293.2583),
187                               new Vector3D(837.552338, -2957.524176, -0.928772));
188 
189         // PEF iau76
190         PVCoordinates pvPEF =
191             new PVCoordinates(new Vector3D(24796919.2956, -34115870.9001, 10293.2583),
192                               new Vector3D(-0.979178, -1.476540, -0.928772));
193 
194         // it seems the induced effect of pole nutation correction δΔψ on the equation of the equinoxes
195         // was not taken into account in the reference paper, so we fix it here for the test
196         final Binary64 dDeltaPsi =
197                 FramesFactory.getEOPHistory(IERSConventions.IERS_1996, true).getEquinoxNutationCorrection(t0)[0];
198         final Binary64 epsilonA = IERSConventions.IERS_1996.getMeanObliquityFunction().value(t0);
199         final FieldTransform<Binary64> fix =
200                 new FieldTransform<>(t0, new FieldRotation<>(FieldVector3D.getPlusK(Binary64Field.getInstance()),
201                                                              dDeltaPsi.multiply(epsilonA.cos()),
202                                                              RotationConvention.FRAME_TRANSFORM));
203 
204         checkPV(fix.transformPVCoordinates(pvPEF), t.transformPVCoordinates(pvTOD), 0.0503, 3.59e-4);
205 
206         // if we forget to apply nutation corrections, results are much worse, which is expected
207         t = FramesFactory.getTOD(false).getTransformTo(FramesFactory.getGTOD(false), t0);
208         checkPV(fix.transformPVCoordinates(pvPEF), t.transformPVCoordinates(pvTOD), 1458.27, 3.847e-4);
209 
210     }
211 
212     @Test
213     void testGetKinematicTransform() {
214         // GIVEN
215         final GTODProvider provider = new GTODProvider(IERSConventions.IERS_2010,
216                 FramesFactory.getEOPHistory(IERSConventions.IERS_2010, true),
217                 DataContext.getDefault().getTimeScales());
218         final AbsoluteDate date = AbsoluteDate.ARBITRARY_EPOCH;
219         // WHEN
220         final KinematicTransform kinematicTransform = provider.getKinematicTransform(date);
221         // THEN
222         final Transform transform = provider.getTransform(date);
223         Assertions.assertEquals(date, kinematicTransform.getDate());
224         Assertions.assertEquals(transform.getCartesian().getPosition(), kinematicTransform.getTranslation());
225         Assertions.assertEquals(transform.getCartesian().getVelocity(), kinematicTransform.getVelocity());
226         Assertions.assertEquals(0., Rotation.distance(transform.getRotation(), kinematicTransform.getRotation()));
227         Assertions.assertEquals(transform.getRotationRate(), kinematicTransform.getRotationRate());
228     }
229 
230     @Test
231     void testGetStaticTransform() {
232         // GIVEN
233         final GTODProvider provider = new GTODProvider(IERSConventions.IERS_2010,
234                 FramesFactory.getEOPHistory(IERSConventions.IERS_2010, true),
235                 DataContext.getDefault().getTimeScales());
236         final AbsoluteDate date = AbsoluteDate.ARBITRARY_EPOCH;
237         // WHEN
238         final StaticTransform staticTransform = provider.getStaticTransform(date);
239         // THEN
240         final Transform transform = provider.getTransform(date);
241         Assertions.assertEquals(date, staticTransform.getDate());
242         Assertions.assertEquals(transform.getCartesian().getPosition(), staticTransform.getTranslation());
243         Assertions.assertEquals(0., Rotation.distance(transform.getRotation(), staticTransform.getRotation()));
244     }
245 
246     @BeforeEach
247     public void setUp() {
248         Utils.setDataRoot("compressed-data");
249     }
250 
251     private void checkPV(PVCoordinates reference, PVCoordinates result,
252                          double expectedPositionError, double expectedVelocityError) {
253 
254         Vector3D dP = result.getPosition().subtract(reference.getPosition());
255         Vector3D dV = result.getVelocity().subtract(reference.getVelocity());
256         Assertions.assertEquals(expectedPositionError, dP.getNorm(), 0.01 * expectedPositionError);
257         Assertions.assertEquals(expectedVelocityError, dV.getNorm(), 0.01 * expectedVelocityError);
258     }
259 
260     private <T extends CalculusFieldElement<T>> void checkPV(FieldPVCoordinates<T> reference,
261                                                          FieldPVCoordinates<T> result,
262                                                          double expectedPositionError,
263                                                          double expectedVelocityError) {
264 
265         FieldVector3D<T> dP = result.getPosition().subtract(reference.getPosition());
266         FieldVector3D<T> dV = result.getVelocity().subtract(reference.getVelocity());
267         Assertions.assertEquals(expectedPositionError, dP.getNorm().getReal(), 0.01 * expectedPositionError);
268         Assertions.assertEquals(expectedVelocityError, dV.getNorm().getReal(), 0.01 * expectedVelocityError);
269     }
270 
271 }