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.attitudes;
18  
19  
20  import java.util.ArrayList;
21  import java.util.List;
22  
23  import org.hipparchus.Field;
24  import org.hipparchus.CalculusFieldElement;
25  import org.hipparchus.geometry.euclidean.threed.Rotation;
26  import org.hipparchus.geometry.euclidean.threed.Vector3D;
27  import org.hipparchus.util.Decimal64;
28  import org.hipparchus.util.Decimal64Field;
29  import org.hipparchus.util.FastMath;
30  import org.junit.After;
31  import org.junit.Assert;
32  import org.junit.Before;
33  import org.junit.Test;
34  import org.orekit.Utils;
35  import org.orekit.bodies.OneAxisEllipsoid;
36  import org.orekit.errors.OrekitException;
37  import org.orekit.frames.Frame;
38  import org.orekit.frames.FramesFactory;
39  import org.orekit.orbits.CircularOrbit;
40  import org.orekit.orbits.FieldOrbit;
41  import org.orekit.orbits.Orbit;
42  import org.orekit.orbits.PositionAngle;
43  import org.orekit.propagation.FieldSpacecraftState;
44  import org.orekit.propagation.Propagator;
45  import org.orekit.propagation.SpacecraftState;
46  import org.orekit.propagation.analytical.KeplerianPropagator;
47  import org.orekit.propagation.sampling.OrekitFixedStepHandler;
48  import org.orekit.time.AbsoluteDate;
49  import org.orekit.time.DateComponents;
50  import org.orekit.time.FieldAbsoluteDate;
51  import org.orekit.time.TimeComponents;
52  import org.orekit.time.TimeScalesFactory;
53  import org.orekit.utils.AngularDerivativesFilter;
54  import org.orekit.utils.FieldPVCoordinates;
55  import org.orekit.utils.IERSConventions;
56  import org.orekit.utils.PVCoordinates;
57  import org.orekit.utils.TimeStampedAngularCoordinates;
58  import org.orekit.utils.TimeStampedFieldPVCoordinates;
59  import org.orekit.utils.TimeStampedPVCoordinates;
60  
61  public class TabulatedProviderTest {
62  
63      // Computation date
64      private AbsoluteDate date;
65  
66      // Reference frame = ITRF
67      private Frame itrf;
68  
69      // Satellite position
70      CircularOrbit circOrbit;
71  
72      // Earth shape
73      OneAxisEllipsoid earthShape;
74  
75      @Test
76      public void testDifferentFrames() {
77          double             samplingRate      = 10.0;
78          int                n                 = 8;
79          AttitudeProvider   referenceProvider = new NadirPointing(circOrbit.getFrame(), earthShape);
80          List<TimeStampedAngularCoordinates> sample = createSample(samplingRate, referenceProvider);
81          TabulatedProvider  provider          = new TabulatedProvider(circOrbit.getFrame(), sample, n,
82                                                                       AngularDerivativesFilter.USE_R);
83          Attitude attE = provider.getAttitude((date, frame) -> new TimeStampedPVCoordinates(date, PVCoordinates.ZERO),
84                                               date,
85                                               circOrbit.getFrame());
86          Assert.assertEquals(circOrbit.getFrame().getName(), attE.getReferenceFrame().getName());
87          Frame gcrf = FramesFactory.getGCRF();
88          Attitude attG = provider.getAttitude((date, frame) -> new TimeStampedPVCoordinates(date, PVCoordinates.ZERO),
89                                               date,
90                                               gcrf);
91          Assert.assertEquals(gcrf.getName(), attG.getReferenceFrame().getName());
92  
93          Assert.assertEquals(1.12e-7,
94                              Rotation.distance(attE.getRotation(), attG.getRotation()), 1.0e-9);
95          Assert.assertEquals(circOrbit.getFrame().getTransformTo(gcrf, date).getRotation().getAngle(),
96                              Rotation.distance(attE.getRotation(), attG.getRotation()),
97                              1.0e-14);
98  
99          FieldAttitude<Decimal64> attG64 =
100                         provider.getAttitude((date, frame) -> new TimeStampedFieldPVCoordinates<>(date,
101                                                                         FieldPVCoordinates.getZero(Decimal64Field.getInstance())),
102                                              new FieldAbsoluteDate<>(Decimal64Field.getInstance(), date),
103                                              gcrf);
104         Assert.assertEquals(gcrf.getName(), attG64.getReferenceFrame().getName());
105 
106     }
107 
108     @Test
109     public void testWithoutRate() {
110         double             samplingRate      = 10.0;
111         double             checkingRate      = 1.0;
112         int                n                 = 8;
113         AttitudeProvider   referenceProvider = new NadirPointing(circOrbit.getFrame(), earthShape);
114         List<TimeStampedAngularCoordinates> sample = createSample(samplingRate, referenceProvider);
115         final double       margin            = samplingRate * n / 2;
116         final AbsoluteDate start             = sample.get(0).getDate().shiftedBy(margin);
117         final AbsoluteDate end               = sample.get(sample.size() - 1).getDate().shiftedBy(-margin);
118         TabulatedProvider  provider          = new TabulatedProvider(circOrbit.getFrame(), sample, n,
119                                                                      AngularDerivativesFilter.USE_R);
120         Assert.assertEquals(0.0, checkError(start, end, checkingRate, referenceProvider, provider), 2.2e-14);
121     }
122 
123     @Test
124     public void testWithRate() {
125         double             samplingRate      = 10.0;
126         double             checkingRate      = 1.0;
127         int                n                 = 8;
128         AttitudeProvider   referenceProvider = new NadirPointing(circOrbit.getFrame(), earthShape);
129         List<TimeStampedAngularCoordinates> sample = createSample(samplingRate, referenceProvider);
130         final double       margin            = samplingRate * n / 2;
131         final AbsoluteDate start             = sample.get(0).getDate().shiftedBy(margin);
132         final AbsoluteDate end               = sample.get(sample.size() - 1).getDate().shiftedBy(-margin);
133         TabulatedProvider  provider          = new TabulatedProvider(circOrbit.getFrame(), sample, n,
134                                                                      AngularDerivativesFilter.USE_RR);
135         Assert.assertEquals(0.0, checkError(start, end, checkingRate, referenceProvider, provider), 1.3e-11);
136     }
137 
138     @Test
139     public void testWithAcceleration() {
140         double             samplingRate      = 10.0;
141         double             checkingRate      = 1.0;
142         int                n                 = 8;
143         AttitudeProvider   referenceProvider = new NadirPointing(circOrbit.getFrame(), earthShape);
144         List<TimeStampedAngularCoordinates> sample = createSample(samplingRate, referenceProvider);
145         final double       margin            = samplingRate * n / 2;
146         final AbsoluteDate start             = sample.get(0).getDate().shiftedBy(margin);
147         final AbsoluteDate end               = sample.get(sample.size() - 1).getDate().shiftedBy(-margin);
148         TabulatedProvider  provider          = new TabulatedProvider(circOrbit.getFrame(), sample, n,
149                                                                      AngularDerivativesFilter.USE_RRA);
150         Assert.assertEquals(0.0, checkError(start, end, checkingRate, referenceProvider, provider), 4.3e-9);
151         checkField(Decimal64Field.getInstance(), provider, circOrbit, circOrbit.getDate(), circOrbit.getFrame());
152     }
153 
154     private List<TimeStampedAngularCoordinates> createSample(double samplingRate, AttitudeProvider referenceProvider) {
155 
156         // reference propagator, using a yaw compensation law
157         final KeplerianPropagator referencePropagator = new KeplerianPropagator(circOrbit);
158         referencePropagator.setAttitudeProvider(referenceProvider);
159 
160         // create sample
161         final List<TimeStampedAngularCoordinates> sample = new ArrayList<TimeStampedAngularCoordinates>();
162         referencePropagator.setStepHandler(samplingRate, currentState -> sample.add(currentState.getAttitude().getOrientation()));
163         referencePropagator.propagate(circOrbit.getDate().shiftedBy(2 * circOrbit.getKeplerianPeriod()));
164 
165         return sample;
166 
167     }
168 
169     private double checkError(final AbsoluteDate start, AbsoluteDate end, double checkingRate,
170                               final AttitudeProvider referenceProvider, TabulatedProvider provider) {
171 
172         // prepare an interpolating provider, using only internal steps
173         // (i.e. ignoring interpolation near boundaries)
174         Propagator interpolatingPropagator = new KeplerianPropagator(circOrbit.shiftedBy(start.durationFrom(circOrbit.getDate())));
175         interpolatingPropagator.setAttitudeProvider(provider);
176 
177         // compute interpolation error on the internal steps .
178         final double[] error = new double[1];
179         interpolatingPropagator.setStepHandler(checkingRate, new OrekitFixedStepHandler() {
180 
181             public void init(SpacecraftState s0, AbsoluteDate t, double step) {
182                 error[0] = 0.0;
183             }
184 
185             public void handleStep(SpacecraftState currentState) {
186                 Attitude interpolated = currentState.getAttitude();
187                 Attitude reference    = referenceProvider.getAttitude(currentState.getOrbit(),
188                                                                       currentState.getDate(),
189                                                                       currentState.getFrame());
190                 double localError = Rotation.distance(interpolated.getRotation(), reference.getRotation());
191                 error[0] = FastMath.max(error[0], localError);
192             }
193 
194         });
195 
196         interpolatingPropagator.propagate(end);
197 
198         return error[0];
199 
200     }
201 
202     private <T extends CalculusFieldElement<T>> void checkField(final Field<T> field, final AttitudeProvider provider,
203                                                             final Orbit orbit, final AbsoluteDate date,
204                                                             final Frame frame) {
205         Attitude attitudeD = provider.getAttitude(orbit, date, frame);
206         final FieldOrbit<T> orbitF = new FieldSpacecraftState<>(field, new SpacecraftState(orbit)).getOrbit();
207         final FieldAbsoluteDate<T> dateF = new FieldAbsoluteDate<>(field, date);
208         FieldAttitude<T> attitudeF = provider.getAttitude(orbitF, dateF, frame);
209         Assert.assertEquals(0.0, Rotation.distance(attitudeD.getRotation(), attitudeF.getRotation().toRotation()), 1.0e-15);
210         Assert.assertEquals(0.0, Vector3D.distance(attitudeD.getSpin(), attitudeF.getSpin().toVector3D()), 1.0e-15);
211         Assert.assertEquals(0.0, Vector3D.distance(attitudeD.getRotationAcceleration(), attitudeF.getRotationAcceleration().toVector3D()), 1.0e-15);
212     }
213 
214     @Before
215     public void setUp() {
216         try {
217             Utils.setDataRoot("regular-data");
218 
219             // Computation date
220             date = new AbsoluteDate(new DateComponents(2008, 04, 07),
221                                     TimeComponents.H00,
222                                     TimeScalesFactory.getUTC());
223 
224             // Body mu
225             final double mu = 3.9860047e14;
226 
227             // Reference frame = ITRF
228             itrf = FramesFactory.getITRF(IERSConventions.IERS_2010, true);
229 
230             //  Satellite position
231             circOrbit =
232                 new CircularOrbit(7178000.0, 0.5e-4, -0.5e-4, FastMath.toRadians(50.), FastMath.toRadians(270.),
233                                        FastMath.toRadians(5.300), PositionAngle.MEAN,
234                                        FramesFactory.getEME2000(), date, mu);
235 
236             // Elliptic earth shape
237             earthShape =
238                 new OneAxisEllipsoid(6378136.460, 1 / 298.257222101, itrf);
239 
240         } catch (OrekitException oe) {
241             Assert.fail(oe.getMessage());
242         }
243 
244     }
245 
246     @After
247     public void tearDown() {
248         date = null;
249         itrf = null;
250         circOrbit = null;
251         earthShape = null;
252     }
253 
254 }
255