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.orbits;
18  
19  import org.hipparchus.Field;
20  import org.hipparchus.analysis.polynomials.SmoothStepFactory;
21  import org.hipparchus.geometry.euclidean.threed.Vector3D;
22  import org.hipparchus.util.Binary64;
23  import org.hipparchus.util.Binary64Field;
24  import org.hipparchus.util.FastMath;
25  import org.junit.jupiter.api.Assertions;
26  import org.junit.jupiter.api.BeforeEach;
27  import org.junit.jupiter.api.DisplayName;
28  import org.junit.jupiter.api.Test;
29  import org.mockito.Mockito;
30  import org.orekit.Utils;
31  import org.orekit.errors.OrekitIllegalArgumentException;
32  import org.orekit.frames.Frame;
33  import org.orekit.frames.FramesFactory;
34  import org.orekit.propagation.analytical.FieldAbstractAnalyticalPropagator;
35  import org.orekit.propagation.analytical.FieldKeplerianPropagator;
36  import org.orekit.time.AbsoluteDate;
37  import org.orekit.time.FieldAbsoluteDate;
38  import org.orekit.time.FieldTimeInterpolator;
39  import org.orekit.utils.Constants;
40  import org.orekit.utils.FieldPVCoordinates;
41  import org.orekit.utils.PVCoordinates;
42  
43  import java.util.ArrayList;
44  import java.util.List;
45  
46  class FieldOrbitBlenderTest {
47  
48      // Constants
49      final double DEFAULT_MU   = Constants.IERS2010_EARTH_MU;
50      final double EARTH_RADIUS = Constants.IERS2010_EARTH_EQUATORIAL_RADIUS;
51  
52      final private Field<Binary64> field = Binary64Field.getInstance();
53      final private Binary64        one   = field.getOne();
54  
55      private FieldOrbit<Binary64> getDefaultOrbitAtDate(final FieldAbsoluteDate<Binary64> date, final Frame inertialFrame) {
56          final Vector3D      position = new Vector3D(EARTH_RADIUS + 4000000, 0, 0);
57          final Vector3D      velocity = new Vector3D(0, FastMath.sqrt(DEFAULT_MU / position.getNorm()), 0);
58          final PVCoordinates pv       = new PVCoordinates(position, velocity);
59  
60          return new FieldCartesianOrbit<>(new FieldPVCoordinates<>(field, pv), inertialFrame, date, one.multiply(DEFAULT_MU));
61      }
62  
63      @BeforeEach
64      public void setUp() {
65          Utils.setDataRoot("regular-data");
66      }
67  
68      @Test
69      @DisplayName("test blending case with Keplerian dynamic (exact results expected)")
70      void testBlendingWithKeplerianDynamic() {
71          // Given
72          final Frame orbitFrame  = FramesFactory.getEME2000();
73          final Frame outputFrame = FramesFactory.getGCRF();
74  
75          final FieldAbsoluteDate<Binary64> previousTabulatedDate = new FieldAbsoluteDate<>(field, new AbsoluteDate());
76          final FieldAbsoluteDate<Binary64> nextTabulatedDate     = previousTabulatedDate.shiftedBy(3600);
77          final FieldAbsoluteDate<Binary64> interpolationDate1    = previousTabulatedDate.shiftedBy(1200);
78          final FieldAbsoluteDate<Binary64> interpolationDate2    = previousTabulatedDate.shiftedBy(1500);
79          final FieldAbsoluteDate<Binary64> interpolationDate3    = previousTabulatedDate.shiftedBy(2000);
80  
81          final FieldOrbit<Binary64> previousTabulatedOrbit =
82                  getDefaultOrbitAtDate(previousTabulatedDate, orbitFrame);
83          final FieldAbstractAnalyticalPropagator<Binary64> propagator =
84                  new FieldKeplerianPropagator<>(previousTabulatedOrbit);
85  
86          final FieldOrbit<Binary64>       nextTabulatedOrbit = propagator.propagate(nextTabulatedDate).getOrbit();
87          final List<FieldOrbit<Binary64>> orbitSample        = new ArrayList<>();
88          orbitSample.add(previousTabulatedOrbit);
89          orbitSample.add(nextTabulatedOrbit);
90  
91          final SmoothStepFactory.FieldSmoothStepFunction<Binary64> blendingFunction = SmoothStepFactory.getQuadratic(field);
92          final FieldOrbitBlender<Binary64> orbitBlender =
93                  new FieldOrbitBlender<>(blendingFunction, propagator, outputFrame);
94  
95          // When & Then
96          final double epsilon = 1e-8; // 10 nm
97  
98          OrbitBlenderTest.assertOrbit(propagator.propagate(interpolationDate1).getOrbit().toOrbit(),
99                                       orbitBlender.interpolate(interpolationDate1, orbitSample).toOrbit(), epsilon);
100         OrbitBlenderTest.assertOrbit(propagator.propagate(interpolationDate2).getOrbit().toOrbit(),
101                                      orbitBlender.interpolate(interpolationDate2, orbitSample).toOrbit(), epsilon);
102         OrbitBlenderTest.assertOrbit(propagator.propagate(interpolationDate3).getOrbit().toOrbit(),
103                                      orbitBlender.interpolate(interpolationDate3, orbitSample).toOrbit(), epsilon);
104 
105         Assertions.assertEquals(outputFrame, orbitBlender.getOutputInertialFrame());
106     }
107 
108     @Test
109     @DisplayName("test specific case (blending at tabulated date)")
110     void testBlendingAtTabulatedDate() {
111         // Given
112         final Frame orbitFrame  = FramesFactory.getEME2000();
113         final Frame outputFrame = FramesFactory.getGCRF();
114 
115         final FieldAbsoluteDate<Binary64> previousTabulatedDate = new FieldAbsoluteDate<>(field, new AbsoluteDate());
116         final FieldAbsoluteDate<Binary64> nextTabulatedDate     = previousTabulatedDate.shiftedBy(3600);
117 
118         final FieldOrbit<Binary64>       previousTabulatedOrbit = getDefaultOrbitAtDate(previousTabulatedDate, orbitFrame);
119         final FieldOrbit<Binary64>       nextTabulatedOrbit     = getDefaultOrbitAtDate(nextTabulatedDate, orbitFrame);
120         final List<FieldOrbit<Binary64>> orbitSample            = new ArrayList<>();
121         orbitSample.add(previousTabulatedOrbit);
122         orbitSample.add(nextTabulatedOrbit);
123 
124         final SmoothStepFactory.FieldSmoothStepFunction<Binary64> blendingFunction = SmoothStepFactory.getQuadratic(field);
125         final FieldTimeInterpolator<FieldOrbit<Binary64>, Binary64> orbitBlender = new FieldOrbitBlender<>(blendingFunction,
126                                                                                                            new FieldKeplerianPropagator<>(
127                                                                                                                    previousTabulatedOrbit),
128                                                                                                            outputFrame);
129 
130         // When
131         final FieldOrbit<Binary64> blendedOrbit = orbitBlender.interpolate(previousTabulatedDate, orbitSample);
132 
133         // Then
134         OrbitBlenderTest.assertOrbit(previousTabulatedOrbit.toOrbit(), blendedOrbit.toOrbit(), 1e-11);
135     }
136 
137     @Test
138     @DisplayName("Test error thrown when using non inertial frame")
139     void testErrorThrownWhenUsingNonInertialFrame() {
140         // Given
141         @SuppressWarnings("unchecked")
142         final SmoothStepFactory.FieldSmoothStepFunction<Binary64> blendingFunctionMock = Mockito.mock(
143                 SmoothStepFactory.FieldSmoothStepFunction.class);
144 
145         @SuppressWarnings("unchecked")
146         final FieldAbstractAnalyticalPropagator<Binary64> propagatorMock =
147                 Mockito.mock(FieldAbstractAnalyticalPropagator.class);
148 
149         final Frame nonInertialFrame = Mockito.mock(Frame.class);
150         Mockito.when(nonInertialFrame.isPseudoInertial()).thenReturn(false);
151 
152         // When & Then
153         Exception thrown = Assertions.assertThrows(OrekitIllegalArgumentException.class,
154                                                    () -> new FieldOrbitBlender<>(blendingFunctionMock,
155                                                                                  propagatorMock,
156                                                                                  nonInertialFrame));
157 
158         Assertions.assertEquals("non pseudo-inertial frame \"null\"", thrown.getMessage());
159     }
160 
161 }