1   /* Copyright 2002-2022 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  
20  import org.hipparchus.geometry.euclidean.threed.Vector3D;
21  import org.hipparchus.util.FastMath;
22  import org.junit.Assert;
23  import org.junit.Before;
24  import org.junit.Test;
25  import org.orekit.Utils;
26  import org.orekit.time.AbsoluteDate;
27  import org.orekit.time.DateComponents;
28  import org.orekit.time.TimeComponents;
29  import org.orekit.time.TimeScalesFactory;
30  import org.orekit.utils.PVCoordinates;
31  
32  
33  public class EME2000ProviderTest {
34  
35      @Test
36      public void testAASReferenceLEO() {
37  
38          // this reference test has been extracted from the following paper:
39          // Implementation Issues Surrounding the New IAU Reference Systems for Astrodynamics
40          // David A. Vallado, John H. Seago, P. Kenneth Seidelmann
41          // http://www.centerforspace.com/downloads/files/pubs/AAS-06-134.pdf
42          AbsoluteDate t0 = new AbsoluteDate(new DateComponents(2004, 04, 06),
43                                             new TimeComponents(07, 51, 28.386009),
44                                             TimeScalesFactory.getUTC());
45  
46          Transform t = FramesFactory.getGCRF().getTransformTo(FramesFactory.getEME2000(), t0);
47  
48          PVCoordinates pvGcrfIau2000A =
49              new PVCoordinates(new Vector3D(5102508.9579, 6123011.4038, 6378136.9252),
50                                new Vector3D(-4743.220156, 790.536497, 5533.755728));
51          PVCoordinates pvEME2000EqA =
52              new PVCoordinates(new Vector3D(5102509.0383, 6123011.9758, 6378136.3118),
53                                new Vector3D(-4743.219766, 790.536344, 5533.756084));
54          checkPV(pvEME2000EqA, t.transformPVCoordinates(pvGcrfIau2000A), 1.1e-4, 2.6e-7);
55  
56          PVCoordinates pvGcrfIau2000B =
57              new PVCoordinates(new Vector3D(5102508.9579, 6123011.4012, 6378136.9277),
58                                new Vector3D(-4743.220156, 790.536495, 5533.755729));
59          PVCoordinates pvEME2000EqB =
60              new PVCoordinates(new Vector3D(5102509.0383, 6123011.9733, 6378136.3142),
61                                new Vector3D(-4743.219766, 790.536342, 5533.756085));
62          checkPV(pvEME2000EqB, t.transformPVCoordinates(pvGcrfIau2000B), 7.4e-5, 2.6e-7);
63  
64      }
65  
66      @Test
67      public void testAASReferenceGEO() {
68  
69          // this reference test has been extracted from the following paper:
70          // Implementation Issues Surrounding the New IAU Reference Systems for Astrodynamics
71          // David A. Vallado, John H. Seago, P. Kenneth Seidelmann
72          // http://www.centerforspace.com/downloads/files/pubs/AAS-06-134.pdf
73          AbsoluteDate t0 = new AbsoluteDate(new DateComponents(2004, 06, 01),
74                                             TimeComponents.H00,
75                                             TimeScalesFactory.getUTC());
76  
77          Transform t = FramesFactory.getGCRF().getTransformTo(FramesFactory.getEME2000(), t0);
78  
79          PVCoordinates pvGCRFiau2000A =
80              new PVCoordinates(new Vector3D(-40588150.3617, -11462167.0397, 27143.1974),
81                                new Vector3D(834.787458, -2958.305691, -1.172993));
82          PVCoordinates pvEME2000EqA =
83              new PVCoordinates(new Vector3D(-40588149.5482, -11462169.9118, 27146.8462),
84                                new Vector3D(834.787667, -2958.305632, -1.172963));
85          checkPV(pvEME2000EqA, t.transformPVCoordinates(pvGCRFiau2000A), 5.8e-5, 6.4e-7);
86  
87          PVCoordinates pvGCRFiau2000B =
88              new PVCoordinates(new Vector3D(-40588150.3617, -11462167.0397, 27143.2125),
89                                new Vector3D(834.787458, -2958.305691, -1.172999));
90          PVCoordinates pvEME2000EqB =
91              new PVCoordinates(new Vector3D(-40588149.5481, -11462169.9118, 27146.8613),
92                                new Vector3D(834.787667, -2958.305632, -1.172968));
93          checkPV(pvEME2000EqB, t.transformPVCoordinates(pvGCRFiau2000B), 1.1e-4, 5.5e-7);
94  
95      }
96  
97      @Test
98      public void testSofaBp00() {
99  
100         // the reference value has been computed using the March 2012 version of the SOFA library
101         // http://www.iausofa.org/2012_0301_C.html, with the following code
102         //
103         //        double utc1, utc2, tai1, tai2, tt1, tt2, rb[3][3], rp[3][3], rbp[3][3];
104         //
105         //        // 2004-02-14:00:00:00Z, MJD = 53049, UT1-UTC = -0.4093509
106         //        utc1  = DJM0 + 53049.0;
107         //        utc2  = 0.0;
108         //        iauUtctai(utc1, utc2, &tai1, &tai2);
109         //        iauTaitt(tai1, tai2, &tt1, &tt2);
110         //
111         //        iauBp00(tt1, tt2, rb, rp, rbp);
112         //
113         //        printf("iauBp00(%.20g, %.20g, rb, rp, rbp)\n"
114         //               "  rb  --> %.20g %.20g %.20g\n          %.20g %.20g %.20g\n          %.20g %.20g %.20g\n"
115         //               "  rp  --> %.20g %.20g %.20g\n          %.20g %.20g %.20g\n          %.20g %.20g %.20g\n"
116         //               "  rbp --> %.20g %.20g %.20g\n          %.20g %.20g %.20g\n          %.20g %.20g %.20g\n",
117         //               tt1, tt2,
118         //               rb[0][0],  rb[0][1],  rb[0][2],
119         //               rb[1][0],  rb[1][1],  rb[1][2],
120         //               rb[2][0],  rb[2][1],  rb[2][2],
121         //               rp[0][0],  rp[0][1],  rp[0][2],
122         //               rp[1][0],  rp[1][1],  rp[1][2],
123         //               rp[2][0],  rp[2][1],  rp[2][2],
124         //               rbp[0][0], rbp[0][1], rbp[0][2],
125         //               rbp[1][0], rbp[1][1], rbp[1][2],
126         //               rbp[2][0], rbp[2][1], rbp[2][2]);
127         //
128         // the output of this test reads:
129         //        iauBp00(2453049.5, 0.00074287037037037029902, rb, rp, rbp)
130         //        rb  --> 0.99999999999999422684 -7.0782797441991980175e-08 8.0562171469761337802e-08
131         //                7.0782794778573375197e-08 0.99999999999999689138 3.3060414542221364117e-08
132         //                -8.0562173809869716745e-08 -3.3060408839805516801e-08 0.99999999999999622524
133         //        rp  --> 0.99999949573309343531 -0.00092105778423522924759 -0.00040023257863225548568
134         //                0.00092105778625203805956 0.99999957582617116092 -1.7927962069881782439e-07
135         //                0.00040023257399096032498 -1.8935780260465051583e-07 0.99999991990692216337
136         //        rbp --> 0.99999949570013624278 -0.00092112855376512230675 -0.00040015204695196122638
137         //                0.00092112856903123034591 0.99999957576097886491 -1.4614501776464880046e-07
138         //                0.00040015201181019732432 -2.2244653837776004327e-07 0.99999991993915571253
139 
140         AbsoluteDate date = new AbsoluteDate(2004, 2, 14, TimeScalesFactory.getUTC());
141         Frame eme2000 = FramesFactory.getFrame(Predefined.EME2000);
142         Frame gcrf    = FramesFactory.getFrame(Predefined.GCRF);
143         checkRotation(new double[][] {
144             {  0.99999999999999422684,    -7.0782797441991980175e-08, 8.0562171469761337802e-08 },
145             {  7.0782794778573375197e-08,  0.99999999999999689138,    3.3060414542221364117e-08 },
146             { -8.0562173809869716745e-08, -3.3060408839805516801e-08, 0.99999999999999622524    }
147 
148         }, gcrf.getTransformTo(eme2000, date), 2.5e-16);
149 
150     }
151 
152     @Before
153     public void setUp() {
154         Utils.setDataRoot("compressed-data");
155     }
156 
157     private void checkPV(PVCoordinates reference,
158                          PVCoordinates result, double positionThreshold,
159                          double velocityThreshold) {
160 
161         Vector3D dP = result.getPosition().subtract(reference.getPosition());
162         Vector3D dV = result.getVelocity().subtract(reference.getVelocity());
163         Assert.assertEquals(0, dP.getNorm(), positionThreshold);
164         Assert.assertEquals(0, dV.getNorm(), velocityThreshold);
165     }
166 
167     private void checkRotation(double[][] reference, Transform t, double epsilon) {
168         double[][] mat = t.getRotation().getMatrix();
169         for (int i = 0; i < 3; ++i) {
170             for (int j = 0; j < 3; ++j) {
171                 Assert.assertEquals(reference[i][j], mat[i][j], epsilon * FastMath.abs(reference[i][j]));
172 
173             }
174         }
175     }
176 
177 }