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.Arrays;
22  import java.util.List;
23  
24  import org.hipparchus.CalculusFieldElement;
25  import org.hipparchus.Field;
26  import org.hipparchus.geometry.euclidean.threed.Rotation;
27  import org.hipparchus.geometry.euclidean.threed.RotationConvention;
28  import org.hipparchus.geometry.euclidean.threed.RotationOrder;
29  import org.hipparchus.geometry.euclidean.threed.Vector3D;
30  import org.hipparchus.random.RandomGenerator;
31  import org.hipparchus.random.Well19937a;
32  import org.hipparchus.util.Decimal64Field;
33  import org.hipparchus.util.FastMath;
34  import org.hipparchus.util.MathUtils;
35  import org.junit.After;
36  import org.junit.Assert;
37  import org.junit.Before;
38  import org.junit.Test;
39  import org.orekit.Utils;
40  import org.orekit.bodies.OneAxisEllipsoid;
41  import org.orekit.errors.OrekitException;
42  import org.orekit.errors.OrekitMessages;
43  import org.orekit.frames.Frame;
44  import org.orekit.frames.FramesFactory;
45  import org.orekit.frames.LOFType;
46  import org.orekit.orbits.CircularOrbit;
47  import org.orekit.orbits.FieldOrbit;
48  import org.orekit.orbits.Orbit;
49  import org.orekit.orbits.PositionAngle;
50  import org.orekit.propagation.FieldSpacecraftState;
51  import org.orekit.propagation.Propagator;
52  import org.orekit.propagation.SpacecraftState;
53  import org.orekit.propagation.analytical.KeplerianPropagator;
54  import org.orekit.time.AbsoluteDate;
55  import org.orekit.time.DateComponents;
56  import org.orekit.time.FieldAbsoluteDate;
57  import org.orekit.time.TimeComponents;
58  import org.orekit.time.TimeScalesFactory;
59  import org.orekit.utils.AngularDerivativesFilter;
60  import org.orekit.utils.Constants;
61  import org.orekit.utils.IERSConventions;
62  import org.orekit.utils.PVCoordinates;
63  import org.orekit.utils.TimeStampedAngularCoordinates;
64  
65  
66  public class TabulatedLofOffsetTest {
67  
68      // Computation date
69      private AbsoluteDate date;
70  
71      // Body mu
72      private double mu;
73  
74      // Reference frame = ITRF
75      private Frame itrf;
76  
77      // Earth shape
78      OneAxisEllipsoid earth;
79  
80      //  Satellite position
81      CircularOrbit orbit;
82      PVCoordinates pvSatEME2000;
83  
84      @Test
85      public void testConstantOffset() {
86  
87          RandomGenerator random = new Well19937a(0x1199d4bb8f53d2b6l);
88          for (LOFType type : LOFType.values()) {
89              for (int i = 0; i < 100; ++i) {
90                  double a1 = random.nextDouble() * MathUtils.TWO_PI;
91                  double a2 = random.nextDouble() * MathUtils.TWO_PI;
92                  double a3 = random.nextDouble() * MathUtils.TWO_PI;
93                  LofOffset law          = new LofOffset(orbit.getFrame(), type, RotationOrder.XYZ, a1, a2, a3);
94                  Rotation  offsetAtt    = law.getAttitude(orbit, orbit.getDate(), orbit.getFrame()).getRotation();
95                  LofOffset aligned      = new LofOffset(orbit.getFrame(), type);
96                  Rotation  alignedAtt   = aligned.getAttitude(orbit, orbit.getDate(), orbit.getFrame()).getRotation();
97                  Rotation  offsetProper = offsetAtt.compose(alignedAtt.revert(), RotationConvention.VECTOR_OPERATOR);
98                  TabulatedLofOffset tabulated =
99                          new TabulatedLofOffset(orbit.getFrame(), type,
100                                                Arrays.asList(new TimeStampedAngularCoordinates(orbit.getDate().shiftedBy(-10),
101                                                                                                offsetProper,
102                                                                                                Vector3D.ZERO,
103                                                                                                Vector3D.ZERO),
104                                                              new TimeStampedAngularCoordinates(orbit.getDate().shiftedBy(0),
105                                                                                                offsetProper,
106                                                                                                Vector3D.ZERO,
107                                                                                                Vector3D.ZERO),
108                                                              new TimeStampedAngularCoordinates(orbit.getDate().shiftedBy(+10),
109                                                                                                offsetProper,
110                                                                                                Vector3D.ZERO,
111                                                                                                Vector3D.ZERO)),
112                                                2, AngularDerivativesFilter.USE_R);
113                 Rotation rebuilt = tabulated.getAttitude(orbit, orbit.getDate(), orbit.getFrame()).getRotation();
114                 Assert.assertEquals(0.0, Rotation.distance(offsetAtt, rebuilt), 1.2e-15);
115                 Assert.assertEquals(3, tabulated.getTable().size());
116             }
117         }
118     }
119 
120     @Test
121     public void testYawCompensation() {
122 
123         // create a sample from Yaw compensation law
124         final LOFType type = LOFType.VNC;
125         final List<TimeStampedAngularCoordinates> sample = new ArrayList<TimeStampedAngularCoordinates>();
126         final AttitudeProvider yawCompensLaw =
127                 new YawCompensation(orbit.getFrame(), new NadirPointing(orbit.getFrame(), earth));
128         final Propagator originalPropagator = new KeplerianPropagator(orbit);
129         originalPropagator.setAttitudeProvider(yawCompensLaw);
130         originalPropagator.setStepHandler(1.0, currentState -> {
131                 Rotation  offsetAtt    = currentState.getAttitude().getRotation();
132                 LofOffset aligned      = new LofOffset(currentState.getFrame(), type);
133                 Rotation  alignedAtt   = aligned.getAttitude(currentState.getOrbit(), currentState.getDate(),
134                                                              currentState.getFrame()).getRotation();
135                 Rotation  offsetProper = offsetAtt.compose(alignedAtt.revert(), RotationConvention.VECTOR_OPERATOR);
136                 sample.add(new TimeStampedAngularCoordinates(currentState.getDate(),
137                                                              offsetProper, Vector3D.ZERO, Vector3D.ZERO));
138             });
139         final AbsoluteDate endDate = orbit.getDate().shiftedBy(2000);
140         originalPropagator.propagate(endDate);
141         originalPropagator.clearStepHandlers();
142 
143         // use the sample and compare it to original
144         final BoundedAttitudeProvider tabulated = new TabulatedLofOffset(orbit.getFrame(), type, sample,
145                                                                          6, AngularDerivativesFilter.USE_RR);
146         Assert.assertEquals(0., orbit.getDate().durationFrom(tabulated.getMinDate()), Double.MIN_VALUE);
147         Assert.assertEquals(0., endDate.durationFrom(tabulated.getMaxDate()), Double.MIN_VALUE);
148         final Propagator rebuildingPropagator = new KeplerianPropagator(orbit);
149         rebuildingPropagator.setAttitudeProvider(tabulated);
150         rebuildingPropagator.setStepHandler(0.3, currentState -> {
151                 final SpacecraftState rebuilt = originalPropagator.propagate(currentState.getDate());
152                 final Rotation r1 = currentState.getAttitude().getRotation();
153                 final Rotation r2 = rebuilt.getAttitude().getRotation();
154                 Assert.assertEquals(0.0, Rotation.distance(r1, r2), 7.0e-6);
155                 checkField(Decimal64Field.getInstance(), tabulated,
156                            currentState.getOrbit(), currentState.getDate(), currentState.getFrame());
157             });
158         rebuildingPropagator.propagate(orbit.getDate().shiftedBy(50), orbit.getDate().shiftedBy(1950));
159 
160     }
161 
162     private <T extends CalculusFieldElement<T>> void checkField(final Field<T> field, final AttitudeProvider provider,
163                                                                 final Orbit orbit, final AbsoluteDate date,
164                                                                 final Frame frame)
165         {
166         Attitude attitudeD = provider.getAttitude(orbit, date, frame);
167         final FieldOrbit<T> orbitF = new FieldSpacecraftState<>(field, new SpacecraftState(orbit)).getOrbit();
168         final FieldAbsoluteDate<T> dateF = new FieldAbsoluteDate<>(field, date);
169         FieldAttitude<T> attitudeF = provider.getAttitude(orbitF, dateF, frame);
170         Assert.assertEquals(0.0, Rotation.distance(attitudeD.getRotation(), attitudeF.getRotation().toRotation()), 1.0e-15);
171         Assert.assertEquals(0.0, Vector3D.distance(attitudeD.getSpin(), attitudeF.getSpin().toVector3D()), 1.0e-15);
172         Assert.assertEquals(0.0, Vector3D.distance(attitudeD.getRotationAcceleration(), attitudeF.getRotationAcceleration().toVector3D()), 1.0e-15);
173     }
174 
175     @Test
176     public void testNonPseudoInertialFrame() {
177         final Frame itrf = FramesFactory.getITRF(IERSConventions.IERS_2010, true);
178         try {
179             new TabulatedLofOffset(itrf, LOFType.QSW,
180                                    Arrays.asList(new TimeStampedAngularCoordinates(orbit.getDate().shiftedBy(-10),
181                                                                                    Rotation.IDENTITY, Vector3D.ZERO, Vector3D.ZERO),
182                                                  new TimeStampedAngularCoordinates(orbit.getDate().shiftedBy(0),
183                                                                                    Rotation.IDENTITY, Vector3D.ZERO, Vector3D.ZERO),
184                                                  new TimeStampedAngularCoordinates(orbit.getDate().shiftedBy(+10),
185                                                                                    Rotation.IDENTITY, Vector3D.ZERO, Vector3D.ZERO)),
186                                    2, AngularDerivativesFilter.USE_R);
187         } catch (OrekitException oe) {
188             Assert.assertEquals(OrekitMessages.NON_PSEUDO_INERTIAL_FRAME, oe.getSpecifier());
189             Assert.assertEquals(itrf.getName(), oe.getParts()[0]);
190         }
191     }
192 
193     @Before
194     public void setUp() {
195         try {
196 
197             Utils.setDataRoot("regular-data");
198 
199             // Computation date
200             date = new AbsoluteDate(new DateComponents(2008, 04, 07),
201                                     TimeComponents.H00,
202                                     TimeScalesFactory.getUTC());
203 
204             // Body mu
205             mu = 3.9860047e14;
206 
207             // Reference frame = ITRF
208             itrf = FramesFactory.getITRF(IERSConventions.IERS_2010, true);
209 
210             // Elliptic earth shape
211             earth = new OneAxisEllipsoid(Constants.WGS84_EARTH_EQUATORIAL_RADIUS,
212                                          Constants.WGS84_EARTH_FLATTENING,
213                                          itrf);
214 
215             //  Satellite position
216             orbit =
217                 new CircularOrbit(7178000.0, 0.5e-8, -0.5e-8, FastMath.toRadians(50.), FastMath.toRadians(150.),
218                                        FastMath.toRadians(5.300), PositionAngle.MEAN,
219                                        FramesFactory.getEME2000(), date, mu);
220             pvSatEME2000 = orbit.getPVCoordinates();
221 
222 
223         } catch (OrekitException oe) {
224             Assert.fail(oe.getMessage());
225         }
226 
227     }
228 
229     @After
230     public void tearDown() {
231         date = null;
232         itrf = null;
233         earth = null;
234     }
235 
236 }
237