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.hamcrest.MatcherAssert;
20  import org.hipparchus.CalculusFieldElement;
21  import org.hipparchus.Field;
22  import org.hipparchus.analysis.UnivariateFunction;
23  import org.hipparchus.analysis.differentiation.DSFactory;
24  import org.hipparchus.analysis.differentiation.FiniteDifferencesDifferentiator;
25  import org.hipparchus.analysis.differentiation.UnivariateDifferentiableFunction;
26  import org.hipparchus.complex.Complex;
27  import org.hipparchus.complex.ComplexField;
28  import org.hipparchus.geometry.euclidean.threed.FieldVector3D;
29  import org.hipparchus.geometry.euclidean.threed.Vector3D;
30  import org.hipparchus.linear.FieldMatrixPreservingVisitor;
31  import org.hipparchus.linear.MatrixUtils;
32  import org.hipparchus.util.Binary64;
33  import org.hipparchus.util.Binary64Field;
34  import org.hipparchus.util.FastMath;
35  import org.hipparchus.util.MathArrays;
36  import org.hipparchus.util.MathUtils;
37  import org.junit.jupiter.api.Assertions;
38  import org.junit.jupiter.api.BeforeEach;
39  import org.junit.jupiter.api.Test;
40  import org.junit.jupiter.params.ParameterizedTest;
41  import org.junit.jupiter.params.provider.EnumSource;
42  import org.orekit.Utils;
43  import org.orekit.errors.OrekitIllegalArgumentException;
44  import org.orekit.errors.OrekitMessages;
45  import org.orekit.frames.Frame;
46  import org.orekit.frames.FramesFactory;
47  import org.orekit.frames.Transform;
48  import org.orekit.time.AbsoluteDate;
49  import org.orekit.time.FieldAbsoluteDate;
50  import org.orekit.time.TimeScalesFactory;
51  import org.orekit.utils.Constants;
52  import org.orekit.utils.FieldPVCoordinates;
53  import org.orekit.utils.PVCoordinates;
54  import org.orekit.utils.TimeStampedFieldPVCoordinates;
55  import org.orekit.utils.TimeStampedPVCoordinates;
56  
57  import java.util.function.Function;
58  
59  import static org.orekit.OrekitMatchers.relativelyCloseTo;
60  
61  
62  class FieldCircularOrbitTest {
63  
64      // Body mu
65      private double mu;
66  
67      @BeforeEach
68      public void setUp() {
69  
70          Utils.setDataRoot("regular-data");
71  
72          // Body mu
73          mu = 3.9860047e14;
74  
75      }
76  
77      @ParameterizedTest
78      @EnumSource(PositionAngleType.class)
79      void testWithCachedPositionAngleType(final PositionAngleType positionAngleType) {
80          // GIVEN
81          final Vector3D position = new Vector3D(-29536113.0, 30329259.0, -100125.0);
82          final Vector3D velocity = new Vector3D(-2194.0, -2141.0, -8.0);
83          final PVCoordinates pvCoordinates = new PVCoordinates(position, velocity);
84          final double muEarth = 3.9860047e14;
85          final CartesianOrbit cartesianOrbit = new CartesianOrbit(pvCoordinates, FramesFactory.getEME2000(), AbsoluteDate.ARBITRARY_EPOCH, muEarth);
86          final Binary64Field field = Binary64Field.getInstance();
87          final FieldCircularOrbit<Binary64> fieldOrbit = new FieldCircularOrbit<>(field, cartesianOrbit);
88          // WHEN
89          final FieldCircularOrbit<Binary64> orbit = fieldOrbit.withCachedPositionAngleType(positionAngleType);
90          // THEN
91          Assertions.assertEquals(fieldOrbit.getFrame(), orbit.getFrame());
92          Assertions.assertEquals(fieldOrbit.getDate(), orbit.getDate());
93          Assertions.assertEquals(fieldOrbit.getMu(), orbit.getMu());
94          final Vector3D relativePosition = fieldOrbit.getPosition(orbit.getFrame()).subtract(
95                  orbit.getPosition()).toVector3D();
96          Assertions.assertEquals(0., relativePosition.getNorm(), 1e-6);
97          Assertions.assertEquals(fieldOrbit.hasNonKeplerianAcceleration(),
98                  orbit.hasNonKeplerianAcceleration());
99      }
100 
101     @Test
102     void test1721() {
103         // GIVEN
104         final PVCoordinates pvCoordinates = new PVCoordinates(new Vector3D(1, 2, 3),
105                 Vector3D.MINUS_K.scalarMultiply(0.1), new Vector3D(0, 1));
106         final CircularOrbit circularOrbit = new CircularOrbit(new TimeStampedPVCoordinates(AbsoluteDate.ARBITRARY_EPOCH,
107                 pvCoordinates), FramesFactory.getEME2000(), 1.);
108         final Binary64Field field = Binary64Field.getInstance();
109         final FieldCircularOrbit<Binary64> fieldCircularOrbit = new FieldCircularOrbit<>(field, circularOrbit);
110         // WHEN
111         final FieldVector3D<Binary64> velocity = fieldCircularOrbit.getVelocity();
112         // THEN
113         final FieldVector3D<Binary64> difference = velocity.subtract(pvCoordinates.getVelocity());
114         Assertions.assertEquals(0., difference.getNorm().getReal(), 1e-10);
115     }
116 
117     @Test
118     void testInFrameNonKeplerian() {
119         testTemplateInFrame(Vector3D.MINUS_J);
120     }
121 
122     @Test
123     void testInFrameKeplerian() {
124         testTemplateInFrame(Vector3D.ZERO);
125     }
126 
127     private void testTemplateInFrame(final Vector3D acceleration) {
128         // GIVEN
129         final Vector3D position = new Vector3D(-29536113.0, 30329259.0, -100125.0);
130         final Vector3D velocity = new Vector3D(-2194.0, -2141.0, -8.0);
131         final PVCoordinates pvCoordinates = new PVCoordinates(position, velocity, acceleration);
132         final double muEarth = 3.9860047e14;
133         final CartesianOrbit cartesianOrbit = new CartesianOrbit(pvCoordinates, FramesFactory.getEME2000(),
134                 AbsoluteDate.ARBITRARY_EPOCH, muEarth);
135         final KeplerianOrbit keplerianOrbit = new KeplerianOrbit(cartesianOrbit);
136         final FieldCircularOrbit<Binary64> fieldOrbit = new FieldCircularOrbit<>(Binary64Field.getInstance(), keplerianOrbit);
137         // WHEN
138         final FieldCircularOrbit<Binary64> fieldOrbitWithOtherFrame = fieldOrbit.inFrame(FramesFactory.getGCRF());
139         // THEN
140         Assertions.assertNotEquals(fieldOrbit.getFrame(), fieldOrbitWithOtherFrame.getFrame());
141         Assertions.assertEquals(fieldOrbit.getDate(), fieldOrbitWithOtherFrame.getDate());
142         Assertions.assertEquals(fieldOrbit.getMu(), fieldOrbitWithOtherFrame.getMu());
143         final FieldVector3D<Binary64> relativePosition = fieldOrbit.getPosition(fieldOrbitWithOtherFrame.getFrame()).subtract(
144                 fieldOrbitWithOtherFrame.getPosition());
145         Assertions.assertEquals(0., relativePosition.getNorm().getReal(), 1e-6);
146         Assertions.assertEquals(fieldOrbit.hasNonKeplerianAcceleration(),
147                 fieldOrbitWithOtherFrame.hasNonKeplerianAcceleration());
148     }
149 
150     @Test
151     void testCircularToEquinoc() {
152         doTestCircularToEquinoctialEll(Binary64Field.getInstance());
153     }
154 
155     @Test
156     void testCircToEquinoc() {
157         doTestCircularToEquinoctialCirc(Binary64Field.getInstance());
158     }
159 
160     @Test
161     void testAnomalyCirc() {
162         doTestAnomalyCirc(Binary64Field.getInstance());
163     }
164 
165     @Test
166     void testAnomalyEll() {
167         doTestAnomalyEll(Binary64Field.getInstance());
168     }
169 
170     @Test
171     void testCircToCart() {
172         doTestCircularToCartesian(Binary64Field.getInstance());
173     }
174 
175     @Test
176     void testCircToKepl() {
177         doTestCircularToKeplerian(Binary64Field.getInstance());
178     }
179 
180     @Test
181     void testGeometryCirc() {
182         doTestGeometryCirc(Binary64Field.getInstance());
183     }
184 
185     @Test
186     void testGeometryEll() {
187         doTestGeometryEll(Binary64Field.getInstance());
188     }
189 
190     @Test
191     void testJacobianFinited() {
192         doTestJacobianFinitedifferences(Binary64Field.getInstance());
193     }
194 
195     @Test
196     void testJacoabianReference() {
197         doTestJacobianReference(Binary64Field.getInstance());
198     }
199 
200     @Test
201     void testNumericalIssue25() {
202         doTestNumericalIssue25(Binary64Field.getInstance());
203     }
204 
205     @Test
206     void testPerfectlyEquatorial() {
207         doTestPerfectlyEquatorial(Binary64Field.getInstance());
208     }
209 
210     @Test
211     void testPositionVelocityNormsCirc() {
212         doTestPositionVelocityNormsCirc(Binary64Field.getInstance());
213     }
214 
215     @Test
216     void testPositionVelocity() {
217         doTestPositionVelocityNormsEll(Binary64Field.getInstance());
218     }
219 
220     @Test
221     void testSymmetryCir() {
222         doTestSymmetryCir(Binary64Field.getInstance());
223     }
224 
225     @Test
226     void testSymmetryEll() {
227         doTestSymmetryEll(Binary64Field.getInstance());
228     }
229 
230     @Test
231     void testErrors()  {
232         Assertions.assertThrows(IllegalArgumentException.class,
233                                 () -> doTestNonInertialFrame(Binary64Field.getInstance()));
234     }
235 
236     @Test
237     void testHyperbolic1() {
238         doTestHyperbolic1(Binary64Field.getInstance());
239     }
240 
241     @Test
242     void testHyperbolic2() {
243         doTestHyperbolic2(Binary64Field.getInstance());
244     }
245 
246     @Test
247     void testToOrbitWithoutDerivatives() {
248         doTestToOrbitWithoutDerivatives(Binary64Field.getInstance());
249     }
250 
251     @Test
252     void testToOrbitWithDerivatives() {
253         doTestToOrbitWithDerivatives(Binary64Field.getInstance());
254     }
255 
256     @Test
257     void testDerivativesConversionSymmetry() {
258         doTestDerivativesConversionSymmetry(Binary64Field.getInstance());
259     }
260 
261     @Test
262     void testToString() {
263         doTestToString(Binary64Field.getInstance());
264     }
265 
266     @Test
267     void testNonKeplerianDerivatives() {
268         doTestNonKeplerianDerivatives(Binary64Field.getInstance());
269     }
270 
271     @Test
272     void testPositionAngleDerivatives() {
273         doTestPositionAngleDerivatives(Binary64Field.getInstance());
274     }
275 
276     @Test
277     void testEquatorialRetrograde() {
278         doTestEquatorialRetrograde(Binary64Field.getInstance());
279     }
280 
281     @Test
282     void testCopyNonKeplerianAcceleration() {
283         doTestCopyNonKeplerianAcceleration(Binary64Field.getInstance());
284     }
285 
286     @Test
287     void testNormalize() {
288         doTestNormalize(Binary64Field.getInstance());
289     }
290 
291     @Test
292     void testWithKeplerianRates() {
293         // GIVEN
294         final ComplexField field = ComplexField.getInstance();
295         final CircularOrbit expectedOrbit = createOrbitTestFromCircularOrbit(true);
296         final FieldCircularOrbit<Complex> fieldOrbit = new FieldCircularOrbit<>(field, expectedOrbit);
297         // WHEN
298         final FieldCircularOrbit<Complex> actualFieldOrbit = fieldOrbit.withKeplerianRates();
299         // THEN
300         Assertions.assertFalse(actualFieldOrbit.hasNonKeplerianRates());
301         Assertions.assertEquals(fieldOrbit.getMu(), actualFieldOrbit.getMu());
302         Assertions.assertEquals(fieldOrbit.getDate(), actualFieldOrbit.getDate());
303         Assertions.assertEquals(fieldOrbit.getFrame(), actualFieldOrbit.getFrame());
304         Assertions.assertEquals(fieldOrbit.getPosition(), actualFieldOrbit.getPosition());
305     }
306 
307     @Test
308     void testFromCircularOrbitWithDerivatives() {
309         // GIVEN
310         final ComplexField field = ComplexField.getInstance();
311         final CircularOrbit expectedOrbit = createOrbitTestFromCircularOrbit(true);
312         // WHEN
313         final FieldCircularOrbit<Complex> fieldOrbit = new FieldCircularOrbit<>(field, expectedOrbit);
314         // THEN
315         compareFieldOrbitToOrbit(expectedOrbit, fieldOrbit);
316     }
317 
318     @Test
319     void testFromCircularOrbitWithoutDerivatives() {
320         // GIVEN
321         final ComplexField field = ComplexField.getInstance();
322         final CircularOrbit expectedOrbit = createOrbitTestFromCircularOrbit(false);
323         // WHEN
324         final FieldCircularOrbit<Complex> fieldOrbit = new FieldCircularOrbit<>(field, expectedOrbit);
325         // THEN
326         compareFieldOrbitToOrbit(expectedOrbit, fieldOrbit);
327     }
328 
329     private CircularOrbit createOrbitTestFromCircularOrbit(final boolean withDerivatives) {
330         final PositionAngleType positionAngleType = PositionAngleType.TRUE;
331         final AbsoluteDate date = AbsoluteDate.ARBITRARY_EPOCH;
332         final Frame frame = FramesFactory.getEME2000();
333         final double a = 10000.e3;
334         final double ex = 1e-4;
335         final double ey = 1e-3;
336         final double i = 1.;
337         final double raan = -2.;
338         final double alphaV = 0.5;
339         if (withDerivatives) {
340             final double derivative = 0.;
341             return new CircularOrbit(a, ex, ey, i, raan, alphaV, derivative, derivative, derivative, derivative,
342                     derivative, derivative, positionAngleType, frame, date, mu);
343         } else {
344             return new CircularOrbit(a, ex, ey, i, raan, alphaV, positionAngleType, frame, date, mu);
345         }
346     }
347 
348     private <T extends CalculusFieldElement<T>> void compareFieldOrbitToOrbit(final CircularOrbit orbit,
349                                                                               final FieldCircularOrbit<T> fieldOrbit) {
350         Assertions.assertEquals(orbit.getFrame(), fieldOrbit.getFrame());
351         Assertions.assertEquals(orbit.getMu(), fieldOrbit.getMu().getReal());
352         Assertions.assertEquals(orbit.getDate(), fieldOrbit.getDate().toAbsoluteDate());
353         Assertions.assertEquals(orbit.getA(), fieldOrbit.getA().getReal());
354         Assertions.assertEquals(orbit.getCircularEx(), fieldOrbit.getCircularEx().getReal());
355         Assertions.assertEquals(orbit.getCircularEy(), fieldOrbit.getCircularEy().getReal());
356         Assertions.assertEquals(orbit.getRightAscensionOfAscendingNode(),
357                 fieldOrbit.getRightAscensionOfAscendingNode().getReal());
358         Assertions.assertEquals(orbit.getI(), fieldOrbit.getI().getReal());
359         Assertions.assertEquals(orbit.getAlphaV(), fieldOrbit.getAlphaV().getReal());
360         Assertions.assertEquals(orbit.hasNonKeplerianAcceleration(), fieldOrbit.hasNonKeplerianAcceleration());
361         Assertions.assertEquals(orbit.getADot(), fieldOrbit.getADot().getReal());
362         Assertions.assertEquals(orbit.getCircularExDot(), fieldOrbit.getCircularExDot().getReal());
363         Assertions.assertEquals(orbit.getCircularEyDot(), fieldOrbit.getCircularEyDot().getReal());
364         Assertions.assertEquals(orbit.getRightAscensionOfAscendingNodeDot(),
365                 fieldOrbit.getRightAscensionOfAscendingNodeDot().getReal());
366         Assertions.assertEquals(orbit.getIDot(), fieldOrbit.getIDot().getReal());
367         Assertions.assertEquals(orbit.getAlphaVDot(), fieldOrbit.getAlphaVDot().getReal());
368     }
369 
370     @Test
371     void testCoverageCachedPositionAngleTypeWithRates() {
372         // GIVEN
373         final double semiMajorAxis = 1e4;
374         final double eccentricity = 0.;
375         final double expectedAlpha = 0.;
376         final double expectedAlphaDot = 0.;
377         final AbsoluteDate date = AbsoluteDate.ARBITRARY_EPOCH;
378         final Binary64Field field = Binary64Field.getInstance();
379         final Binary64 zero = field.getZero();
380         // WHEN & THEN
381         for (final PositionAngleType inputPositionAngleType: PositionAngleType.values()) {
382             for (final PositionAngleType cachedPositionAngleType: PositionAngleType.values()) {
383                 final FieldCircularOrbit<Binary64> fieldOrbit = new FieldCircularOrbit<>(
384                         zero.newInstance(semiMajorAxis), zero.newInstance(eccentricity), zero, zero, zero,
385                         zero.newInstance(expectedAlpha), zero, zero, zero, zero, zero,
386                         zero.newInstance(expectedAlphaDot), inputPositionAngleType, cachedPositionAngleType,
387                         FramesFactory.getGCRF(), new FieldAbsoluteDate<>(field, date), zero.newInstance(mu));
388                 Assertions.assertEquals(cachedPositionAngleType, fieldOrbit.getCachedPositionAngleType());
389                 Assertions.assertEquals(expectedAlpha, fieldOrbit.getAlphaV().getReal());
390                 Assertions.assertEquals(expectedAlpha, fieldOrbit.getAlphaM().getReal());
391                 Assertions.assertEquals(expectedAlpha, fieldOrbit.getAlphaE().getReal());
392                 Assertions.assertEquals(expectedAlphaDot, fieldOrbit.getAlphaVDot().getReal());
393                 Assertions.assertEquals(expectedAlphaDot, fieldOrbit.getAlphaMDot().getReal());
394                 Assertions.assertEquals(expectedAlphaDot, fieldOrbit.getAlphaEDot().getReal());
395             }
396         }
397     }
398 
399     @Test
400     void testGetAlphaVersusDouble() {
401         // GIVEN
402         final double semiMajorAxis = 1e7;
403         final double eccentricity = 1e-2;
404         final double expectedAlpha = 4.;
405         final AbsoluteDate date = AbsoluteDate.ARBITRARY_EPOCH;
406         final Binary64Field field = Binary64Field.getInstance();
407         final Binary64 zero = field.getZero();
408         // WHEN & THEN
409         for (final PositionAngleType inputPositionAngleType: PositionAngleType.values()) {
410             for (final PositionAngleType cachedPositionAngleType: PositionAngleType.values()) {
411                 final FieldCircularOrbit<Binary64> fieldOrbit = new FieldCircularOrbit<>(
412                         zero.newInstance(semiMajorAxis), zero.newInstance(eccentricity), zero, zero, zero,
413                         zero.newInstance(expectedAlpha), inputPositionAngleType, cachedPositionAngleType,
414                         FramesFactory.getGCRF(), new FieldAbsoluteDate<>(field, date), zero.newInstance(mu));
415                 final CircularOrbit circularOrbit = fieldOrbit.toOrbit();
416                 Assertions.assertEquals(circularOrbit.getAlphaV(), fieldOrbit.getAlphaV().getReal());
417                 Assertions.assertEquals(circularOrbit.getAlphaM(), fieldOrbit.getAlphaM().getReal());
418                 Assertions.assertEquals(circularOrbit.getAlphaE(), fieldOrbit.getAlphaE().getReal());
419             }
420         }
421     }
422 
423     private <T extends CalculusFieldElement<T>> void doTestCircularToEquinoctialEll(Field<T> field) {
424 
425         T zero =  field.getZero();
426         FieldAbsoluteDate<T> date = new FieldAbsoluteDate<>(field);
427 
428         T ix = zero.add(1.200e-04);
429         T iy = zero.add(-1.16e-04);
430         T i  = ix.multiply(ix).add(iy.multiply(iy)).divide(4).sqrt().asin().multiply(2);
431 
432 
433         T raan = iy.atan2(ix);
434 
435         // elliptic orbit
436         FieldCircularOrbit<T> circ =
437             new FieldCircularOrbit<>(zero.add(42166712.0), zero.add(0.5), zero.add(-0.5), i, raan,
438                                      zero.add(5.300).subtract(raan), PositionAngleType.MEAN,
439                                      FramesFactory.getEME2000(), date, zero.add(mu));
440         FieldVector3D<T> pos = circ.getPosition();
441         FieldVector3D<T> vit = circ.getVelocity();
442 
443         FieldPVCoordinates<T> pvCoordinates = new FieldPVCoordinates<>(pos, vit);
444 
445 
446 
447         FieldEquinoctialOrbit<T> param = new FieldEquinoctialOrbit<>(pvCoordinates, FramesFactory.getEME2000(), date, zero.add(mu));
448 
449         Assertions.assertEquals(param.getA().getReal(),  circ.getA().getReal(), Utils.epsilonTest * circ.getA().getReal());
450         Assertions.assertEquals(param.getEquinoctialEx().getReal(), circ.getEquinoctialEx().getReal(), Utils.epsilonE * FastMath.abs(circ.getE().getReal()));
451         Assertions.assertEquals(param.getEquinoctialEy().getReal(), circ.getEquinoctialEy().getReal(), Utils.epsilonE * FastMath.abs(circ.getE().getReal()));
452         Assertions.assertEquals(param.getHx().getReal(), circ.getHx().getReal(), Utils.epsilonAngle * FastMath.abs(circ.getI().getReal()));
453         Assertions.assertEquals(param.getHy().getReal(), circ.getHy().getReal(), Utils.epsilonAngle * FastMath.abs(circ.getI().getReal()));
454 
455 
456         Assertions.assertEquals(MathUtils.normalizeAngle(param.getLv().getReal(), circ.getLv().getReal()),
457                                 circ.getLv().getReal(),
458                                 Utils.epsilonAngle * FastMath.abs(circ.getLv().getReal()));
459 
460         Assertions.assertFalse(circ.hasNonKeplerianAcceleration());
461     }
462 
463     private <T extends CalculusFieldElement<T>> void doTestToOrbitWithoutDerivatives(Field<T> field) {
464         T zero =  field.getZero();
465         FieldAbsoluteDate<T> date = new FieldAbsoluteDate<>(field);
466 
467         FieldVector3D<T> position = new FieldVector3D<>(zero.add(7.0e6), zero.add(1.0e6), zero.add(4.0e6));
468         FieldVector3D<T> velocity = new FieldVector3D<>(zero.add(-500.0), zero.add(8000.0), zero.add(1000.0));
469         FieldPVCoordinates<T> pvCoordinates = new FieldPVCoordinates<>(position, velocity);
470         FieldCircularOrbit<T>  fieldOrbit = new FieldCircularOrbit<>(pvCoordinates, FramesFactory.getEME2000(), date, zero.add(mu));
471         CircularOrbit orbit = fieldOrbit.toOrbit();
472         Assertions.assertFalse(orbit.hasNonKeplerianAcceleration());
473         MatcherAssert.assertThat(orbit.getA(),                             relativelyCloseTo(fieldOrbit.getA().getReal(),                             0));
474         MatcherAssert.assertThat(orbit.getCircularEx(),                    relativelyCloseTo(fieldOrbit.getCircularEx().getReal(),                    0));
475         MatcherAssert.assertThat(orbit.getCircularEy(),                    relativelyCloseTo(fieldOrbit.getCircularEy().getReal(),                    0));
476         MatcherAssert.assertThat(orbit.getI(),                             relativelyCloseTo(fieldOrbit.getI().getReal(),                             0));
477         MatcherAssert.assertThat(orbit.getRightAscensionOfAscendingNode(), relativelyCloseTo(fieldOrbit.getRightAscensionOfAscendingNode().getReal(), 0));
478         MatcherAssert.assertThat(orbit.getAlphaV(),                        relativelyCloseTo(fieldOrbit.getAlphaV().getReal(),                        0));
479     }
480 
481     private <T extends CalculusFieldElement<T>> void doTestToOrbitWithDerivatives(Field<T> field) {
482         T zero =  field.getZero();
483         FieldAbsoluteDate<T> date = new FieldAbsoluteDate<>(field);
484 
485         FieldVector3D<T> position = new FieldVector3D<>(zero.add(7.0e6), zero.add(1.0e6), zero.add(4.0e6));
486         FieldVector3D<T> velocity = new FieldVector3D<>(zero.add(-500.0), zero.add(8000.0), zero.add(1000.0));
487         T r2 = position.getNormSq();
488         T r = r2.sqrt();
489         final FieldVector3D<T> keplerianAcceleration = new FieldVector3D<>(r.multiply(r2).reciprocal().multiply(zero.add(mu).negate()),
490                 position);
491         final FieldVector3D<T> nonKeplerianAcceleration = keplerianAcceleration.scalarMultiply(1.1);
492         final FieldPVCoordinates<T> pvCoordinates = new FieldPVCoordinates<>(position, velocity, nonKeplerianAcceleration);
493         FieldCircularOrbit<T>  fieldOrbit = new FieldCircularOrbit<>(pvCoordinates, FramesFactory.getEME2000(), date, zero.add(mu));
494         CircularOrbit orbit = fieldOrbit.toOrbit();
495         Assertions.assertTrue(orbit.hasNonKeplerianAcceleration());
496         MatcherAssert.assertThat(orbit.getA(),                                 relativelyCloseTo(fieldOrbit.getA().getReal(),                                 0));
497         MatcherAssert.assertThat(orbit.getCircularEx(),                        relativelyCloseTo(fieldOrbit.getCircularEx().getReal(),                        0));
498         MatcherAssert.assertThat(orbit.getCircularEy(),                        relativelyCloseTo(fieldOrbit.getCircularEy().getReal(),                        0));
499         MatcherAssert.assertThat(orbit.getI(),                                 relativelyCloseTo(fieldOrbit.getI().getReal(),                                 0));
500         MatcherAssert.assertThat(orbit.getRightAscensionOfAscendingNode(),     relativelyCloseTo(fieldOrbit.getRightAscensionOfAscendingNode().getReal(),     0));
501         MatcherAssert.assertThat(orbit.getAlphaV(),                            relativelyCloseTo(fieldOrbit.getAlphaV().getReal(),                            0));
502         MatcherAssert.assertThat(orbit.getADot(),                              relativelyCloseTo(fieldOrbit.getADot().getReal(),                              0));
503         MatcherAssert.assertThat(orbit.getEquinoctialExDot(),                  relativelyCloseTo(fieldOrbit.getEquinoctialExDot().getReal(),                  0));
504         MatcherAssert.assertThat(orbit.getEquinoctialEyDot(),                  relativelyCloseTo(fieldOrbit.getEquinoctialEyDot().getReal(),                  0));
505         MatcherAssert.assertThat(orbit.getHxDot(),                             relativelyCloseTo(fieldOrbit.getHxDot().getReal(),                             0));
506         MatcherAssert.assertThat(orbit.getHyDot(),                             relativelyCloseTo(fieldOrbit.getHyDot().getReal(),                             0));
507         MatcherAssert.assertThat(orbit.getLvDot(),                             relativelyCloseTo(fieldOrbit.getLvDot().getReal(),                             0));
508         MatcherAssert.assertThat(orbit.getLEDot(),                             relativelyCloseTo(fieldOrbit.getLEDot().getReal(),                             0));
509         MatcherAssert.assertThat(orbit.getLMDot(),                             relativelyCloseTo(fieldOrbit.getLMDot().getReal(),                             0));
510         MatcherAssert.assertThat(orbit.getEDot(),                              relativelyCloseTo(fieldOrbit.getEDot().getReal(),                              0));
511         MatcherAssert.assertThat(orbit.getIDot(),                              relativelyCloseTo(fieldOrbit.getIDot().getReal(),                              0));
512         MatcherAssert.assertThat(orbit.getCircularExDot(),                     relativelyCloseTo(fieldOrbit.getCircularExDot().getReal(),                     0));
513         MatcherAssert.assertThat(orbit.getCircularEyDot(),                     relativelyCloseTo(fieldOrbit.getCircularEyDot().getReal(),                     0));
514         MatcherAssert.assertThat(orbit.getRightAscensionOfAscendingNodeDot(),  relativelyCloseTo(fieldOrbit.getRightAscensionOfAscendingNodeDot().getReal(),  0));
515         MatcherAssert.assertThat(orbit.getAlphaVDot(),                         relativelyCloseTo(fieldOrbit.getAlphaVDot().getReal(),                         0));
516         MatcherAssert.assertThat(orbit.getAlphaEDot(),                         relativelyCloseTo(fieldOrbit.getAlphaEDot().getReal(),                         0));
517         MatcherAssert.assertThat(orbit.getAlphaMDot(),                         relativelyCloseTo(fieldOrbit.getAlphaMDot().getReal(),                         0));
518         MatcherAssert.assertThat(orbit.getAlphaDot(PositionAngleType.TRUE),        relativelyCloseTo(fieldOrbit.getAlphaDot(PositionAngleType.TRUE).getReal(),        0));
519         MatcherAssert.assertThat(orbit.getAlphaDot(PositionAngleType.ECCENTRIC),   relativelyCloseTo(fieldOrbit.getAlphaDot(PositionAngleType.ECCENTRIC).getReal(),   0));
520         MatcherAssert.assertThat(orbit.getAlphaDot(PositionAngleType.MEAN),        relativelyCloseTo(fieldOrbit.getAlphaDot(PositionAngleType.MEAN).getReal(),        0));
521     }
522 
523     private <T extends CalculusFieldElement<T>> void doTestCircularToEquinoctialCirc(Field<T> field) {
524 
525         T zero =  field.getZero();
526         FieldAbsoluteDate<T> date = new FieldAbsoluteDate<>(field);
527 
528         T ix = zero.add(1.200e-04);
529         T iy = zero.add(-1.16e-04);
530         T i  = ix.multiply(ix).add(iy.multiply(iy)).divide(4).sqrt().asin().multiply(2);
531         T raan = iy.atan2(ix);
532 
533         // circular orbit
534         FieldEquinoctialOrbit<T> circCir =
535             new FieldEquinoctialOrbit<>(zero.add(42166712.0), zero.add(0.1e-10), zero.add(-0.1e-10), i, raan,
536                                         raan.negate().add(5.300), PositionAngleType.MEAN,
537                                         FramesFactory.getEME2000(), date, zero.add(mu));
538         FieldVector3D<T> posCir = circCir.getPosition();
539         FieldVector3D<T> vitCir = circCir.getVelocity();
540 
541         FieldPVCoordinates<T> pvCoordinates = new FieldPVCoordinates<>( posCir, vitCir);
542 
543         FieldEquinoctialOrbit<T> paramCir = new FieldEquinoctialOrbit<>(pvCoordinates, FramesFactory.getEME2000(), date, zero.add(mu));
544         Assertions.assertEquals(paramCir.getA().getReal(), circCir.getA().getReal(), Utils.epsilonTest * circCir.getA().getReal());
545         Assertions.assertEquals(paramCir.getEquinoctialEx().getReal(), circCir.getEquinoctialEx().getReal(), Utils.epsilonEcir * FastMath.abs(circCir.getE().getReal()));
546         Assertions.assertEquals(paramCir.getEquinoctialEy().getReal(), circCir.getEquinoctialEy().getReal(), Utils.epsilonEcir * FastMath.abs(circCir.getE().getReal()));
547         Assertions.assertEquals(paramCir.getHx().getReal(), circCir.getHx().getReal(), Utils.epsilonAngle * FastMath.abs(circCir.getI().getReal()));
548         Assertions.assertEquals(paramCir.getHy().getReal(), circCir.getHy().getReal(), Utils.epsilonAngle * FastMath.abs(circCir.getI().getReal()));
549         Assertions.assertEquals(MathUtils.normalizeAngle(paramCir.getLv().getReal(), circCir.getLv().getReal()), circCir.getLv().getReal(), Utils.epsilonAngle * FastMath.abs(circCir.getLv().getReal()));
550 
551     }
552 
553     private <T extends CalculusFieldElement<T>> void doTestCircularToCartesian(Field<T> field) {
554         T zero =  field.getZero();
555         FieldAbsoluteDate<T> date = new FieldAbsoluteDate<>(field);
556 
557         T ix = zero.add(1.200e-04);
558         T iy = zero.add(-1.16e-04);
559         T i  = ix.multiply(ix).add(iy.multiply(iy)).divide(4).sqrt().asin().multiply(2);
560         T raan = iy.atan2(ix);
561         T cosRaan = raan.cos();
562         T sinRaan = raan.sin();
563         T exTilde = zero.add(-7.900e-6);
564         T eyTilde = zero.add(1.100e-4);
565         T ex = exTilde.multiply(cosRaan).add(eyTilde.multiply(sinRaan));
566         T ey = eyTilde.multiply(cosRaan).subtract(exTilde.multiply(sinRaan));
567 
568         FieldCircularOrbit<T> circ=
569             new FieldCircularOrbit<>(zero.add(42166712.0), ex, ey, i, raan,
570                                      raan.negate().add(5.300), PositionAngleType.MEAN,
571                                      FramesFactory.getEME2000(), date, zero.add(mu));
572         FieldVector3D<T> pos = circ.getPosition();
573         FieldVector3D<T> vel = circ.getVelocity();
574 
575         // check 1/a = 2/r  - V2/mu
576         T r = pos.getNorm();
577         T v = vel.getNorm();
578         Assertions.assertEquals(2 / r.getReal() - v.getReal() * v.getReal() / mu, 1 / circ.getA().getReal(), 1.0e-7);
579 
580         Assertions.assertEquals( 0.233745668678733e+08, pos.getX().getReal(), Utils.epsilonTest * r.getReal());
581         Assertions.assertEquals(-0.350998914352669e+08, pos.getY().getReal(), Utils.epsilonTest * r.getReal());
582         Assertions.assertEquals(-0.150053723123334e+04, pos.getZ().getReal(), Utils.epsilonTest * r.getReal());
583 
584         Assertions.assertEquals(2558.7096558809967, vel.getX().getReal(), Utils.epsilonTest * v.getReal());
585         Assertions.assertEquals(1704.1586039092576, vel.getY().getReal(), Utils.epsilonTest * v.getReal());
586         Assertions.assertEquals(   0.5013093577879, vel.getZ().getReal(), Utils.epsilonTest * v.getReal());
587 
588     }
589 
590     private <T extends CalculusFieldElement<T>> void doTestCircularToKeplerian(Field<T> field) {
591         T zero =  field.getZero();
592         FieldAbsoluteDate<T> date = new FieldAbsoluteDate<>(field);
593 
594         T ix   = zero.add(1.20e-4);
595         T iy   = zero.add(-1.16e-4);
596         T i    = ix.multiply(ix).add(iy.multiply(iy)).divide(4).sqrt().asin().multiply(2);
597         T raan = iy.atan2(ix);
598         T cosRaan = raan.cos();
599         T sinRaan = raan.sin();
600         T exTilde = zero.add(-7.900e-6);
601         T eyTilde = zero.add(1.100e-4);
602         T ex = exTilde.multiply(cosRaan).add(eyTilde.multiply(sinRaan));
603         T ey = eyTilde.multiply(cosRaan).subtract(exTilde.multiply(sinRaan));
604 
605         FieldCircularOrbit<T> circ=
606             new FieldCircularOrbit<>(zero.add(42166712.0), ex, ey, i, raan,
607                                      raan.negate().add(5.300), PositionAngleType.MEAN,
608                                      FramesFactory.getEME2000(), date, zero.add(mu));
609         FieldKeplerianOrbit<T> kep = new FieldKeplerianOrbit<>(circ);
610 
611         Assertions.assertEquals(42166712.000, circ.getA().getReal(), Utils.epsilonTest * kep.getA().getReal());
612         Assertions.assertEquals(0.110283316961361e-03, kep.getE().getReal(), Utils.epsilonE * FastMath.abs(kep.getE().getReal()));
613         Assertions.assertEquals(0.166901168553917e-03, kep.getI().getReal(),
614                      Utils.epsilonAngle * FastMath.abs(kep.getI().getReal()));
615         Assertions.assertEquals(MathUtils.normalizeAngle(-3.87224326008837, kep.getPerigeeArgument().getReal()),
616                      kep.getPerigeeArgument().getReal(),
617                      Utils.epsilonTest * 6 * FastMath.abs(kep.getPerigeeArgument().getReal())); //numerical propagation we changed to 6 time the precision used
618         Assertions.assertEquals(MathUtils.normalizeAngle(5.51473467358854, kep.getRightAscensionOfAscendingNode().getReal()),
619                      kep.getRightAscensionOfAscendingNode().getReal(),
620                      Utils.epsilonTest * FastMath.abs(kep.getRightAscensionOfAscendingNode().getReal()));
621 
622         Assertions.assertEquals(MathUtils.normalizeAngle(zero.add(3.65750858649982), kep.getMeanAnomaly()).getReal(),
623                      kep.getMeanAnomaly().getReal(),
624                      Utils.epsilonTest * 5 * FastMath.abs(kep.getMeanAnomaly().getReal())); //numerical propagation we changed to 6 time the precision used
625 
626     }
627 
628     private <T extends CalculusFieldElement<T>> void doTestHyperbolic1(Field<T> field) {
629         T zero =  field.getZero();
630         FieldAbsoluteDate<T> date = new FieldAbsoluteDate<>(field);
631         try {
632             new FieldCircularOrbit<>(zero.add(-42166712.0), zero.add(1.9), zero.add(0.5), zero.add(0.01), zero.add(-0.02), zero.add( 5.300),
633                                      PositionAngleType.MEAN,  FramesFactory.getEME2000(), date, zero.add(mu));
634             Assertions.fail("an exception should have been thrown");
635         } catch (OrekitIllegalArgumentException oe) {
636             Assertions.assertEquals(OrekitMessages.HYPERBOLIC_ORBIT_NOT_HANDLED_AS, oe.getSpecifier());
637         }
638     }
639 
640     private <T extends CalculusFieldElement<T>>  void doTestHyperbolic2(Field<T> field) {
641         T zero =  field.getZero();
642        FieldAbsoluteDate<T> date = new FieldAbsoluteDate<>(field);
643         FieldOrbit<T> orbit = new FieldKeplerianOrbit<>(zero.add(-42166712.0), zero.add(1.9), zero.add(0.5), zero.add(0.01), zero.add(-0.02), zero.add( 5.300),
644                                                         PositionAngleType.MEAN,  FramesFactory.getEME2000(), date, zero.add(mu));
645         try {
646             new FieldCircularOrbit<>(orbit.getPVCoordinates(), orbit.getFrame(), orbit.getMu());
647             Assertions.fail("an exception should have been thrown");
648         } catch (OrekitIllegalArgumentException oe) {
649             Assertions.assertEquals(OrekitMessages.HYPERBOLIC_ORBIT_NOT_HANDLED_AS, oe.getSpecifier());
650         }
651     }
652 
653     private <T extends CalculusFieldElement<T>> void doTestAnomalyEll(Field<T> field) {
654         T zero =  field.getZero();
655         T one =  field.getOne();
656         FieldAbsoluteDate<T> date = new FieldAbsoluteDate<>(field);
657         // elliptic orbit
658         FieldVector3D<T> position = new FieldVector3D<>(zero.add(7.0e6), zero.add(1.0e6), zero.add(4.0e6));
659         FieldVector3D<T> velocity = new FieldVector3D<>(zero.add(-500.0), zero.add(8000.0), zero.add(1000.0));
660 
661         FieldPVCoordinates<T> pvCoordinates = new FieldPVCoordinates<>( position, velocity);
662 
663         FieldCircularOrbit<T>  p   = new FieldCircularOrbit<>(pvCoordinates, FramesFactory.getEME2000(), date, zero.add(mu));
664         FieldKeplerianOrbit<T> kep = new FieldKeplerianOrbit<>(p);
665 
666         T e       = p.getE();
667         T eRatio  = one.subtract(e).divide(one.add(e)).sqrt();
668         T raan    = kep.getRightAscensionOfAscendingNode();
669         T paPraan = kep.getPerigeeArgument().add(raan);
670 
671         T lv = zero.add(1.1);
672         // formulations for elliptic case
673         T lE = lv.subtract(paPraan).divide(2).tan().multiply(eRatio).atan().multiply(2).add(paPraan);
674         T lM = lE.subtract(e.multiply(lE.subtract(paPraan).sin()));
675 
676         p = new FieldCircularOrbit<>(p.getA() , p.getCircularEx(), p.getCircularEy(),
677                                      p.getRightAscensionOfAscendingNode(),
678                                      p.getAlphaV(), lv.subtract(raan), PositionAngleType.TRUE, p.getFrame(), date, zero.add(mu));
679         Assertions.assertEquals(p.getAlphaV().getReal() + raan.getReal(), lv.getReal(), Utils.epsilonAngle * FastMath.abs(lv.getReal()));
680         Assertions.assertEquals(p.getAlphaE().getReal() + raan.getReal(), lE.getReal(), Utils.epsilonAngle * FastMath.abs(lE.getReal()));
681         Assertions.assertEquals(p.getAlphaM().getReal() + raan.getReal(), lM.getReal(), Utils.epsilonAngle * FastMath.abs(lM.getReal()));
682         p = new FieldCircularOrbit<>(p.getA() , p.getCircularEx(), p.getCircularEy(),
683                                      p.getRightAscensionOfAscendingNode(),
684                                      p.getAlphaV(), zero, PositionAngleType.TRUE, p.getFrame(), date, zero.add(mu));
685 
686 
687         p = new FieldCircularOrbit<>(p.getA() , p.getCircularEx(), p.getCircularEy(),
688                                      p.getRightAscensionOfAscendingNode(),
689                                      p.getAlphaV(), lE.subtract(raan), PositionAngleType.ECCENTRIC, p.getFrame(), date, zero.add(mu));
690         Assertions.assertEquals(p.getAlphaV().getReal() + raan.getReal(), lv.getReal(), Utils.epsilonAngle * FastMath.abs(lv.getReal()));
691         Assertions.assertEquals(p.getAlphaE().getReal() + raan.getReal(), lE.getReal(), Utils.epsilonAngle * FastMath.abs(lE.getReal()));
692         Assertions.assertEquals(p.getAlphaM().getReal() + raan.getReal(), lM.getReal(), Utils.epsilonAngle * FastMath.abs(lM.getReal()));
693         p = new FieldCircularOrbit<>(p.getA() , p.getCircularEx(), p.getCircularEy(),
694                                      p.getRightAscensionOfAscendingNode(),
695                                      p.getAlphaV(), zero, PositionAngleType.TRUE, p.getFrame(), date, zero.add(mu));
696 
697         p = new FieldCircularOrbit<>(p.getA() , p.getCircularEx(), p.getCircularEy(),
698                                      p.getRightAscensionOfAscendingNode(),
699                                      p.getAlphaV(), lM.subtract(raan), PositionAngleType.MEAN, p.getFrame(), date, zero.add(mu));
700         Assertions.assertEquals(p.getAlphaV().getReal() + raan.getReal(), lv.getReal(), Utils.epsilonAngle * FastMath.abs(lv.getReal()));
701         Assertions.assertEquals(p.getAlphaE().getReal() + raan.getReal(), lE.getReal(), Utils.epsilonAngle * FastMath.abs(lE.getReal()));
702         Assertions.assertEquals(p.getAlphaM().getReal() + raan.getReal(), lM.getReal(), Utils.epsilonAngle * FastMath.abs(lM.getReal()));
703 
704     }
705 
706     private <T extends CalculusFieldElement<T>> void doTestAnomalyCirc(Field<T> field) {
707         T zero =  field.getZero();
708         FieldAbsoluteDate<T> date = new FieldAbsoluteDate<>(field);
709 
710         FieldVector3D<T> position = new FieldVector3D<>(zero.add(7.0e6), zero.add(1.0e6), zero.add(4.0e6));
711         FieldVector3D<T> velocity = new FieldVector3D<>(zero.add(-500.0), zero.add(8000.0), zero.add(1000.0));
712         FieldPVCoordinates<T> pvCoordinates = new FieldPVCoordinates<>( position, velocity);
713         FieldCircularOrbit<T>  p = new FieldCircularOrbit<>(pvCoordinates, FramesFactory.getEME2000(), date, zero.add(mu));
714         T raan = p.getRightAscensionOfAscendingNode();
715 
716         // circular orbit
717         p = new FieldCircularOrbit<>(p.getA() , zero, zero, p.getRightAscensionOfAscendingNode(),
718                                      p.getAlphaV(), p.getAlphaV(), PositionAngleType.TRUE, p.getFrame(), date, zero.add(mu));
719 
720         T lv = zero.add(1.1);
721         T lE = lv;
722         T lM = lE;
723 
724         p = new FieldCircularOrbit<>(p.getA() , p.getCircularEx(), p.getCircularEy(),
725                                      p.getRightAscensionOfAscendingNode(),
726                                      p.getAlphaV(), lv.subtract(raan), PositionAngleType.TRUE, p.getFrame(), date, zero.add(mu));
727         Assertions.assertEquals(p.getAlphaV().getReal() + raan.getReal(), lv.getReal(), Utils.epsilonAngle * FastMath.abs(lv.getReal()));
728         Assertions.assertEquals(p.getAlphaE().getReal() + raan.getReal(), lE.getReal(), Utils.epsilonAngle * FastMath.abs(lE.getReal()));
729         Assertions.assertEquals(p.getAlphaM().getReal() + raan.getReal(), lM.getReal(), Utils.epsilonAngle * FastMath.abs(lM.getReal()));
730         p = new FieldCircularOrbit<>(p.getA() , p.getCircularEx(), p.getCircularEy(),
731                                      p.getRightAscensionOfAscendingNode(),
732                                      p.getAlphaV(), zero, PositionAngleType.TRUE, p.getFrame(), date, zero.add(mu));
733 
734         p = new FieldCircularOrbit<>(p.getA() , p.getCircularEx(), p.getCircularEy(),
735                                      p.getRightAscensionOfAscendingNode(),
736                                      p.getAlphaV(), lE.subtract(raan), PositionAngleType.ECCENTRIC, p.getFrame(), date, zero.add(mu));
737 
738         Assertions.assertEquals(p.getAlphaV().getReal() + raan.getReal(), lv.getReal(), Utils.epsilonAngle * FastMath.abs(lv.getReal()));
739         Assertions.assertEquals(p.getAlphaE().getReal() + raan.getReal(), lE.getReal(), Utils.epsilonAngle * FastMath.abs(lE.getReal()));
740         Assertions.assertEquals(p.getAlphaM().getReal() + raan.getReal(), lM.getReal(), Utils.epsilonAngle * FastMath.abs(lM.getReal()));
741         p = new FieldCircularOrbit<>(p.getA() , p.getCircularEx(), p.getCircularEy(),
742                                      p.getRightAscensionOfAscendingNode(),
743                                      p.getAlphaV(), zero, PositionAngleType.TRUE, p.getFrame(), date, zero.add(mu));
744 
745         p = new FieldCircularOrbit<>(p.getA() , p.getCircularEx(), p.getCircularEy(),
746                                      p.getRightAscensionOfAscendingNode(),
747                                      p.getAlphaV(), lM.subtract(raan), PositionAngleType.MEAN, p.getFrame(), date, zero.add(mu));
748         Assertions.assertEquals(p.getAlphaV().getReal() + raan.getReal(), lv.getReal(), Utils.epsilonAngle * FastMath.abs(lv.getReal()));
749         Assertions.assertEquals(p.getAlphaE().getReal() + raan.getReal(), lE.getReal(), Utils.epsilonAngle * FastMath.abs(lE.getReal()));
750         Assertions.assertEquals(p.getAlphaM().getReal() + raan.getReal(), lM.getReal(), Utils.epsilonAngle * FastMath.abs(lM.getReal()));
751 
752     }
753 
754     private <T extends CalculusFieldElement<T>> void doTestPositionVelocityNormsEll(Field<T> field) {
755         T zero =  field.getZero();
756         T one =  field.getOne();
757         FieldAbsoluteDate<T> date = new FieldAbsoluteDate<>(field);
758 
759         // elliptic and non equatorial (i retrograde) orbit
760         T hx =  zero.add(1.2);
761         T hy =  zero.add(2.1);
762         T i  = hx.multiply(hx).add(hy.multiply(hy)).sqrt().atan().multiply(2);
763         T raan = hy.atan2(hx);
764         FieldCircularOrbit<T> p =
765             new FieldCircularOrbit<>(zero.add(42166712.0), zero.add(0.5), zero.add(-0.5), i, raan,
766                                      raan.negate().add(0.67), PositionAngleType.TRUE,
767                                      FramesFactory.getEME2000(), date, zero.add(mu));
768 
769         T ex = p.getEquinoctialEx();
770         T ey = p.getEquinoctialEy();
771         T lv = p.getLv();
772         T ksi     = ex.multiply(lv.cos()).add(1).add(ey.multiply(lv.sin()));
773         T nu      = ex.multiply(lv.sin()).subtract(ey.multiply(lv.cos()));
774         T epsilon = one.subtract(ex.multiply(ex)).subtract(ey.multiply(ey)).sqrt();
775 
776         T a  = p.getA();
777         T na = a.reciprocal().multiply(mu).sqrt();
778 
779         Assertions.assertEquals(a.getReal() * epsilon.getReal() * epsilon.getReal() / ksi.getReal(),
780                      p.getPosition().getNorm().getReal(),
781                      Utils.epsilonTest * FastMath.abs(p.getPosition().getNorm().getReal()));
782         Assertions.assertEquals(na.getReal() * FastMath.sqrt(ksi.getReal() * ksi.getReal() + nu.getReal() * nu.getReal()) / epsilon.getReal(),
783                      p.getVelocity().getNorm().getReal(),
784                      Utils.epsilonTest * FastMath.abs(p.getVelocity().getNorm().getReal()));
785 
786     }
787 
788     private <T extends CalculusFieldElement<T>> void doTestNumericalIssue25(Field<T> field) {
789 
790         T zero =  field.getZero();
791         FieldVector3D<T> position = new FieldVector3D<>(zero.add(3782116.14107698), zero.add(416663.11924914), zero.add(5875541.62103057));
792         FieldVector3D<T> velocity = new FieldVector3D<>(zero.add(-6349.7848910501), zero.add(288.4061811651), zero.add(4066.9366759691));
793         FieldCircularOrbit<T> orbit = new FieldCircularOrbit<>(new FieldPVCoordinates<>(position, velocity),
794                                                                FramesFactory.getEME2000(),
795                                                                new FieldAbsoluteDate<>(field, "2004-01-01T23:00:00.000",
796                                                                                        TimeScalesFactory.getUTC()),
797                                                                zero.add(3.986004415E14));
798         Assertions.assertEquals(0.0, orbit.getE().getReal(), 2.0e-14);
799     }
800 
801     private <T extends CalculusFieldElement<T>> void doTestPerfectlyEquatorial(Field<T> field) {
802         T zero =  field.getZero();
803         FieldVector3D<T> position = new FieldVector3D<>(zero.add(-7293947.695148368), zero.add( 5122184.668436634), zero.add(0.0));
804         FieldVector3D<T> velocity = new FieldVector3D<>(zero.add(-3890.4029433398), zero.add( -5369.811285264604), zero.add(0.0));
805         FieldCircularOrbit<T> orbit = new FieldCircularOrbit<>(new FieldPVCoordinates<>(position, velocity),
806                                                                FramesFactory.getEME2000(),
807                                                                new FieldAbsoluteDate<>(field, "2004-01-01T23:00:00.000",
808                                                                                        TimeScalesFactory.getUTC()),
809                                                                zero.add(3.986004415E14));
810         Assertions.assertEquals(0.0, orbit.getI().getReal(), 2.0e-14);
811         Assertions.assertEquals(0.0, orbit.getRightAscensionOfAscendingNode().getReal(), 2.0e-14);
812     }
813 
814     private <T extends CalculusFieldElement<T>> void doTestPositionVelocityNormsCirc(Field<T> field) {
815         T zero =  field.getZero();
816         T one =  field.getOne();
817         FieldAbsoluteDate<T> date = new FieldAbsoluteDate<>(field);
818 
819         // elliptic and non equatorial (i retrograde) orbit
820         T hx =  zero.add(0.1e-8);
821         T hy =  zero.add(0.1e-8);
822         T i  = hx.multiply(hx).add(hy.multiply(hy)).sqrt().atan().multiply(2);
823         T raan = hy.atan2(hx);
824         FieldCircularOrbit<T> pCirEqua =
825             new FieldCircularOrbit<>(zero.add(42166712.0), zero.add(0.1e-8), zero.add(0.1e-8), i, raan,
826                                      raan.negate().add(0.67), PositionAngleType.TRUE,
827                                      FramesFactory.getEME2000(), date, zero.add(mu));
828 
829         T ex = pCirEqua.getEquinoctialEx();
830         T ey = pCirEqua.getEquinoctialEy();
831         T lv = pCirEqua.getLv();
832         T ksi     = ex.multiply(lv.cos()).add(1).add(ey.multiply(lv.sin()));
833         T nu      = ex.multiply(lv.sin()).subtract(ey.multiply(lv.cos()));
834         T epsilon = one.subtract(ex.multiply(ex)).subtract(ey.multiply(ey)).sqrt();
835 
836         T a  = pCirEqua.getA();
837         T na = a.reciprocal().multiply(mu).sqrt();
838 
839         Assertions.assertEquals(a.getReal() * epsilon.getReal() * epsilon.getReal() / ksi.getReal(),
840                             pCirEqua.getPosition().getNorm().getReal(),
841                             Utils.epsilonTest * FastMath.abs(pCirEqua.getPosition().getNorm().getReal()));
842         Assertions.assertEquals(na.getReal() * FastMath.sqrt(ksi.getReal() * ksi.getReal() + nu.getReal() * nu.getReal()) / epsilon.getReal(),
843                             pCirEqua.getVelocity().getNorm().getReal(),
844                             Utils.epsilonTest * FastMath.abs(pCirEqua.getVelocity().getNorm().getReal()));
845     }
846 
847     private <T extends CalculusFieldElement<T>> void doTestGeometryEll(Field<T> field) {
848 
849         T zero =  field.getZero();
850         FieldAbsoluteDate<T> date = new FieldAbsoluteDate<>(field);
851         // elliptic and non equatorial (i retrograde) orbit
852         T hx =  zero.add(1.2);
853         T hy =  zero.add(2.1);
854         T i  = hx.multiply(hx).add(hy.multiply(hy)).sqrt().atan().multiply(2);
855         T raan = hy.atan2(hx);
856         FieldCircularOrbit<T> p =
857             new FieldCircularOrbit<>(zero.add(42166712.0), zero.add(0.5), zero.add(-0.5), i, raan,
858                                      raan.negate().add(0.67), PositionAngleType.TRUE,
859                                      FramesFactory.getEME2000(), date, zero.add(mu));
860 
861         FieldVector3D<T> position = p.getPosition();
862         FieldVector3D<T> velocity = p.getVelocity();
863         FieldVector3D<T> momentum = p.getPVCoordinates().getMomentum().normalize();
864 
865         T apogeeRadius  = p.getA().multiply( p.getE().add(1));
866         T perigeeRadius = p.getA().multiply(p.getE().negate().add(1));
867 
868         for (T alphaV = zero; alphaV.getReal() <= 2 * FastMath.PI; alphaV=alphaV.add(zero.add(2).multiply(FastMath.PI/100.))) {
869             p = new FieldCircularOrbit<>(p.getA() , p.getCircularEx(), p.getCircularEy(), p.getI(),
870                                          p.getRightAscensionOfAscendingNode(),
871                                          alphaV, PositionAngleType.TRUE, p.getFrame(), date, zero.add(mu));
872             position = p.getPosition();
873             // test if the norm of the position is in the range [perigee radius, apogee radius]
874             // Warning: these tests are without absolute value by choice
875             Assertions.assertTrue((position.getNorm().getReal() - apogeeRadius.getReal())  <= (  apogeeRadius.getReal() * Utils.epsilonTest));
876             Assertions.assertTrue((position.getNorm().getReal() - perigeeRadius.getReal()) >= (- perigeeRadius.getReal() * Utils.epsilonTest));
877 
878             position= position.normalize();
879             velocity = p.getVelocity();
880             velocity= velocity.normalize();
881 
882             // at this stage of computation, all the vectors (position, velocity and momemtum) are normalized here
883 
884             // test of orthogonality between position and momentum
885             Assertions.assertTrue(FastMath.abs(Vector3D.dotProduct(position.toVector3D(), momentum.toVector3D())) < Utils.epsilonTest);
886             // test of orthogonality between velocity and momentum
887             Assertions.assertTrue(FastMath.abs(Vector3D.dotProduct(velocity.toVector3D(), momentum.toVector3D())) < Utils.epsilonTest);
888         }
889 
890     }
891 
892     private <T extends CalculusFieldElement<T>> void doTestGeometryCirc(Field<T> field) {
893 
894         T zero =  field.getZero();
895         FieldAbsoluteDate<T> date = new FieldAbsoluteDate<>(field);
896 
897         //  circular and equatorial orbit
898         T hx =  zero.add(0.1e-8);
899         T hy =  zero.add(0.1e-8);
900         T i  =  hx.multiply(hx).add(hy.multiply(hy)).sqrt().atan().multiply(2);
901         T raan = hy.atan2(hx);
902         FieldCircularOrbit<T> pCirEqua =
903             new FieldCircularOrbit<>(zero.add(42166712.0), zero.add(0.1e-8), zero.add(0.1e-8), i, raan,
904                                     raan.negate().add(0.67), PositionAngleType.TRUE,
905                                     FramesFactory.getEME2000(), date, zero.add(mu));
906 
907         FieldVector3D<T> position = pCirEqua.getPosition();
908         FieldVector3D<T> velocity = pCirEqua.getVelocity();
909         FieldVector3D<T> momentum = pCirEqua.getPVCoordinates().getMomentum().normalize();
910 
911         T apogeeRadius  = pCirEqua.getA().multiply( pCirEqua.getE().add(1));
912         T perigeeRadius = pCirEqua.getA().multiply(pCirEqua.getE().negate().add(1));
913         // test if apogee equals perigee
914         Assertions.assertEquals(perigeeRadius.getReal(), apogeeRadius.getReal(), 1.e+4 * Utils.epsilonTest * apogeeRadius.getReal());
915 
916         for (T alphaV = zero; alphaV.getReal() <= 2 * FastMath.PI; alphaV = alphaV.add(zero.add(2 * FastMath.PI/100.))) {
917             pCirEqua = new FieldCircularOrbit<>(pCirEqua.getA() , pCirEqua.getCircularEx(), pCirEqua.getCircularEy(), pCirEqua.getI(),
918                                                pCirEqua.getRightAscensionOfAscendingNode(),
919                                                alphaV, PositionAngleType.TRUE, pCirEqua.getFrame(), date, zero.add(mu));
920             position = pCirEqua.getPosition();
921 
922             // test if the norm pf the position is in the range [perigee radius, apogee radius]
923             Assertions.assertTrue((position.getNorm().getReal() - apogeeRadius.getReal())  <= (  apogeeRadius.getReal() * Utils.epsilonTest));
924             Assertions.assertTrue((position.getNorm().getReal() - perigeeRadius.getReal()) >= (- perigeeRadius.getReal() * Utils.epsilonTest));
925 
926             position= position.normalize();
927             velocity = pCirEqua.getVelocity();
928             velocity= velocity.normalize();
929 
930             // at this stage of computation, all the vectors (position, velocity and momemtum) are normalized here
931 
932             // test of orthogonality between position and momentum
933             Assertions.assertTrue(FastMath.abs(Vector3D.dotProduct(position.toVector3D(), momentum.toVector3D())) < Utils.epsilonTest);
934             // test of orthogonality between velocity and momentum
935             Assertions.assertTrue(FastMath.abs(Vector3D.dotProduct(velocity.toVector3D(), momentum.toVector3D())) < Utils.epsilonTest);
936         }
937     }
938 
939     private <T extends CalculusFieldElement<T>> void doTestSymmetryEll(Field<T> field) {
940 
941         T zero =  field.getZero();
942         FieldAbsoluteDate<T> date = new FieldAbsoluteDate<>(field);
943 
944         // elliptic and non equatorial orbit
945         FieldVector3D<T> position = new FieldVector3D<>(zero.add(4512.9), zero.add(18260.), zero.add(-5127.));
946         FieldVector3D<T> velocity = new FieldVector3D<>(zero.add(134664.6), zero.add(90066.8), zero.add(72047.6));
947         FieldPVCoordinates<T> pvCoordinates = new FieldPVCoordinates<>(position, velocity);
948 
949         FieldCircularOrbit<T> p = new FieldCircularOrbit<>(pvCoordinates, FramesFactory.getEME2000(), date, zero.add(mu));
950 
951         FieldVector3D<T> positionOffset = p.getPosition();
952         FieldVector3D<T> velocityOffset = p.getVelocity();
953 
954         positionOffset = positionOffset.subtract(position);
955         velocityOffset = velocityOffset.subtract(velocity);
956 
957         Assertions.assertEquals(0.0, positionOffset.getNorm().getReal(), position.getNorm().getReal() * Utils.epsilonTest);
958         Assertions.assertEquals(0.0, velocityOffset.getNorm().getReal(), velocity.getNorm().getReal() * Utils.epsilonTest);
959 
960     }
961 
962     private <T extends CalculusFieldElement<T>> void doTestSymmetryCir(Field<T> field) {
963         T zero =  field.getZero();
964         FieldAbsoluteDate<T> date = new FieldAbsoluteDate<>(field);
965 
966         // circular and equatorial orbit
967         FieldVector3D<T> position = new FieldVector3D<>(zero.add(33051.2), zero.add(26184.9), zero.add(-1.3E-5));
968         FieldVector3D<T> velocity = new FieldVector3D<>(zero.add(-60376.2), zero.add(76208.), zero.add(2.7E-4));
969         FieldPVCoordinates<T> pvCoordinates = new FieldPVCoordinates<>(position, velocity);
970 
971         FieldCircularOrbit<T> p = new FieldCircularOrbit<>(pvCoordinates, FramesFactory.getEME2000(), date, zero.add(mu));
972 
973         FieldVector3D<T> positionOffset = p.getPosition().subtract(position);
974         FieldVector3D<T> velocityOffset = p.getVelocity().subtract(velocity);
975 
976         Assertions.assertEquals(0.0, positionOffset.getNorm().getReal(), position.getNorm().getReal() * Utils.epsilonTest);
977         Assertions.assertEquals(0.0, velocityOffset.getNorm().getReal(), velocity.getNorm().getReal() * Utils.epsilonTest);
978 
979     }
980 
981     private <T extends CalculusFieldElement<T>> void doTestNonInertialFrame(Field<T> field) throws IllegalArgumentException {
982         T zero =  field.getZero();
983         FieldAbsoluteDate<T> date = new FieldAbsoluteDate<>(field);
984 
985         FieldVector3D<T> position = new FieldVector3D<>(zero.add(33051.2), zero.add(26184.9), zero.add(-1.3E-5));
986         FieldVector3D<T> velocity = new FieldVector3D<>(zero.add(-60376.2), zero.add(76208.), zero.add(2.7E-4));
987         FieldPVCoordinates<T> pvCoordinates = new FieldPVCoordinates<>( position, velocity);
988         new FieldCircularOrbit<>(pvCoordinates,
989                                  new Frame(FramesFactory.getEME2000(), Transform.IDENTITY, "non-inertial", false),
990                                  date, zero.add(mu));
991     }
992 
993     private <T extends CalculusFieldElement<T>> void doTestJacobianReference(Field<T> field) {
994         T zero =  field.getZero();
995         FieldAbsoluteDate<T> dateTca = new FieldAbsoluteDate<>(field, 2000, 4, 1, 0, 0, 0.000, TimeScalesFactory.getUTC());
996         double mu =  3.986004415e+14;
997         FieldCircularOrbit<T> orbCir = new FieldCircularOrbit<>(zero.add(7000000.0), zero.add(0.01), zero.add(-0.02), zero.add(1.2), zero.add(2.1),
998                                                                 zero.add(0.7), PositionAngleType.MEAN,
999                                                                 FramesFactory.getEME2000(), dateTca, zero.add(mu));
1000 
1001         // the following reference values have been computed using the free software
1002         // version 6.2 of the MSLIB fortran library by the following program:
1003         //        program cir_jacobian
1004         //
1005         //        use mslib
1006         //        implicit none
1007         //
1008         //        integer, parameter :: nb = 11
1009         //        integer :: i,j
1010         //        type(tm_code_retour)      ::  code_retour
1011         //
1012         //        real(pm_reel), parameter :: mu= 3.986004415e+14_pm_reel
1013         //        real(pm_reel),dimension(3)::vit_car,pos_car
1014         //        type(tm_orb_cir)::cir
1015         //        real(pm_reel), dimension(6,6)::jacob
1016         //        real(pm_reel)::norme
1017         //
1018         //
1019         //        cir%a=7000000_pm_reel
1020         //        cir%ex=0.01_pm_reel
1021         //        cir%ey=-0.02_pm_reel
1022         //        cir%i=1.2_pm_reel
1023         //        cir%gom=2.1_pm_reel
1024         //        cir%pso_M=0.7_pm_reel
1025         //
1026         //        call mv_cir_car(mu,cir,pos_car,vit_car,code_retour)
1027         //        write(*,*)code_retour%valeur
1028         //        write(*,1000)pos_car,vit_car
1029         //
1030         //
1031         //        call mu_norme(pos_car,norme,code_retour)
1032         //        write(*,*)norme
1033         //
1034         //        call mv_car_cir (mu, pos_car, vit_car, cir, code_retour, jacob)
1035         //        write(*,*)code_retour%valeur
1036         //
1037         //        write(*,*)"circular = ", cir%a, cir%ex, cir%ey, cir%i, cir%gom, cir%pso_M
1038         //
1039         //        do i = 1,6
1040         //           write(*,*) " ",(jacob(i,j),j=1,6)
1041         //        end do
1042         //
1043         //        1000 format (6(f24.15,1x))
1044         //        end program cir_jacobian
1045         FieldVector3D<T> pRef = new FieldVector3D<>(zero.add(-4106905.105389204807580), zero.add( 3603162.539798960555345), zero.add(4439730.167038885876536));
1046         FieldVector3D<T> vRef = new FieldVector3D<>(zero.add(740.132407342422994), zero.add(-5308.773280141396754), zero.add( 5250.338353483879473));
1047         double[][] jRefR = {
1048             { -1.1535467596325562      ,        1.0120556393573172,        1.2470306024626943,        181.96913090864561,       -1305.2162699469984,        1290.8494448855752 },
1049             { -5.07367368325471104E-008, -1.27870567070456834E-008,  1.31544531338558113E-007, -3.09332106417043592E-005, -9.60781276304445404E-005,  1.91506964883791605E-004 },
1050             { -6.59428471712402018E-008,  1.24561703203882533E-007, -1.41907027322388158E-008,  7.63442601186485441E-005, -1.77446722746170009E-004,  5.99464401287846734E-005 },
1051             {  7.55079920652274275E-008,  4.41606835295069131E-008,  3.40079310688458225E-008,  7.89724635377817962E-005,  4.61868720707717372E-005,  3.55682891687782599E-005 },
1052             { -9.20788748896973282E-008, -5.38521280004949642E-008, -4.14712660805579618E-008,  7.78626692360739821E-005,  4.55378113077967091E-005,  3.50684505810897702E-005 },
1053             {  1.85082436324531617E-008,  1.20506219457886855E-007, -8.31277842285972640E-008,  1.27364008345789645E-004, -1.54770720974742483E-004, -1.78589436862677754E-004 }
1054         };
1055 
1056         T[][] jRef = MathArrays.buildArray(field, 6, 6);
1057 
1058         for (int ii = 0; ii<6 ; ii++){
1059             for (int jj = 0; jj<6 ; jj++){
1060                 jRef[ii][jj] = zero.add(jRefR[ii][jj]);
1061             }
1062         }
1063 
1064         FieldPVCoordinates<T> pv = orbCir.getPVCoordinates();
1065         Assertions.assertEquals(0, pv.getPosition().subtract(pRef).getNorm().getReal(), 3.0e-16 * pRef.getNorm().getReal());
1066         Assertions.assertEquals(0, pv.getVelocity().subtract(vRef).getNorm().getReal(), 2.0e-12 * vRef.getNorm().getReal());
1067 
1068         T[][] jacobian = MathArrays.buildArray(field, 6, 6);
1069 
1070         orbCir.getJacobianWrtCartesian(PositionAngleType.MEAN, jacobian);
1071 
1072         for (int i = 0; i < jacobian.length; i++) {
1073             T[] row    = jacobian[i];
1074             T[] rowRef = jRef[i];
1075             for (int j = 0; j < row.length; j++) {
1076                 Assertions.assertEquals(0, (row[j].getReal() - rowRef[j].getReal()) / rowRef[j].getReal(), 1e-14);
1077             }
1078         }
1079 
1080     }
1081 
1082     private <T extends CalculusFieldElement<T>> void doTestJacobianFinitedifferences(Field<T> field) {
1083         T zero =  field.getZero();
1084 
1085         FieldAbsoluteDate<T> dateTca = new FieldAbsoluteDate<>(field, 2000, 4, 1, 0, 0, 0.000, TimeScalesFactory.getUTC());
1086         double mu =  3.986004415e+14;
1087         FieldCircularOrbit<T> orbCir = new FieldCircularOrbit<>(zero.add(7000000.0), zero.add(0.01), zero.add(-0.02), zero.add(1.2), zero.add(2.1),
1088                                                                 zero.add(0.7), PositionAngleType.MEAN,
1089                                                                 FramesFactory.getEME2000(), dateTca, zero.add(mu));
1090 
1091         for (PositionAngleType type : PositionAngleType.values()) {
1092             T hP = zero.add(2.0);
1093             T[][] finiteDiffJacobian = finiteDifferencesJacobian(type, orbCir, hP);
1094             T[][] jacobian = MathArrays.buildArray(field, 6, 6);
1095             orbCir.getJacobianWrtCartesian(type, jacobian);
1096 
1097             for (int i = 0; i < jacobian.length; i++) {
1098                 T[] row    = jacobian[i];
1099                 T[] rowRef = finiteDiffJacobian[i];
1100 
1101                 for (int j = 0; j < row.length; j++) {
1102                     Assertions.assertEquals(0, (row[j].getReal() - rowRef[j].getReal()) / rowRef[j].getReal(), 8.0e-9);
1103                 }
1104 
1105            }
1106 
1107             T[][] invJacobian = MathArrays.buildArray(field, 6, 6);
1108             orbCir.getJacobianWrtParameters(type, invJacobian);
1109             MatrixUtils.createFieldMatrix(jacobian).
1110                             multiply(MatrixUtils.createFieldMatrix(invJacobian)).
1111             walkInRowOrder(new FieldMatrixPreservingVisitor<T>() {
1112                 public void start(int rows, int columns,
1113                                   int startRow, int endRow, int startColumn, int endColumn) {
1114                 }
1115 
1116                 public void visit(int row, int column, T value) {
1117                     Assertions.assertEquals(row == column ? 1.0 : 0.0, value.getReal(), 3.0e-9);
1118                 }
1119 
1120                 public T end() {
1121                     return null;
1122                 }
1123             });
1124 
1125         }
1126 
1127     }
1128 
1129     private <T extends CalculusFieldElement<T>> T[][] finiteDifferencesJacobian(PositionAngleType type, FieldCircularOrbit<T> orbit, T hP)
1130         {
1131         Field<T> field = hP.getField();
1132         T[][] jacobian = MathArrays.buildArray(field, 6, 6);
1133         for (int i = 0; i < 6; ++i) {
1134             fillColumn(type, i, orbit, hP, jacobian);
1135         }
1136         return jacobian;
1137     }
1138 
1139     private <T extends CalculusFieldElement<T>> void fillColumn(PositionAngleType type, int i, FieldCircularOrbit<T> orbit, T hP, T[][] jacobian) {
1140 
1141         T zero = hP.getField().getZero();
1142         // at constant energy (i.e. constant semi major axis), we have dV = -mu dP / (V * r^2)
1143         // we use this to compute a velocity step size from the position step size
1144         FieldVector3D<T> p = orbit.getPosition();
1145         FieldVector3D<T> v = orbit.getVelocity();
1146         T hV =  hP.multiply(orbit.getMu()).divide(v.getNorm().multiply(p.getNormSq()));
1147 
1148         T h;
1149         FieldVector3D<T> dP = new FieldVector3D<>(zero, zero, zero);
1150         FieldVector3D<T> dV = new FieldVector3D<>(zero, zero, zero);
1151         switch (i) {
1152         case 0:
1153             h = hP;
1154             dP = new FieldVector3D<>(hP, zero, zero);
1155             break;
1156         case 1:
1157             h = hP;
1158             dP = new FieldVector3D<>(zero, hP, zero);
1159             break;
1160         case 2:
1161             h = hP;
1162             dP = new FieldVector3D<>(zero, zero, hP);
1163             break;
1164         case 3:
1165             h = hV;
1166             dV = new FieldVector3D<>(hV, zero, zero);
1167             break;
1168         case 4:
1169             h = hV;
1170             dV = new FieldVector3D<>(zero, hV, zero);
1171             break;
1172         default:
1173             h = hV;
1174             dV = new FieldVector3D<>(zero, zero, hV);
1175             break;
1176         }
1177 
1178         FieldCircularOrbit<T> oM4h = new FieldCircularOrbit<>(new FieldPVCoordinates<>(new FieldVector3D<>(1, p, -4, dP), new FieldVector3D<>(1, v, -4, dV)),
1179                                                               orbit.getFrame(), orbit.getDate(), orbit.getMu());
1180         FieldCircularOrbit<T> oM3h = new FieldCircularOrbit<>(new FieldPVCoordinates<>(new FieldVector3D<>(1, p, -3, dP), new FieldVector3D<>(1, v, -3, dV)),
1181                                                               orbit.getFrame(), orbit.getDate(), orbit.getMu());
1182         FieldCircularOrbit<T> oM2h = new FieldCircularOrbit<>(new FieldPVCoordinates<>(new FieldVector3D<>(1, p, -2, dP), new FieldVector3D<>(1, v, -2, dV)),
1183                                                               orbit.getFrame(), orbit.getDate(), orbit.getMu());
1184         FieldCircularOrbit<T> oM1h = new FieldCircularOrbit<>(new FieldPVCoordinates<>(new FieldVector3D<>(1, p, -1, dP), new FieldVector3D<>(1, v, -1, dV)),
1185                                                               orbit.getFrame(), orbit.getDate(), orbit.getMu());
1186         FieldCircularOrbit<T> oP1h = new FieldCircularOrbit<>(new FieldPVCoordinates<>(new FieldVector3D<>(1, p, +1, dP), new FieldVector3D<>(1, v, +1, dV)),
1187                                                               orbit.getFrame(), orbit.getDate(), orbit.getMu());
1188         FieldCircularOrbit<T> oP2h = new FieldCircularOrbit<>(new FieldPVCoordinates<>(new FieldVector3D<>(1, p, +2, dP), new FieldVector3D<>(1, v, +2, dV)),
1189                                                               orbit.getFrame(), orbit.getDate(), orbit.getMu());
1190         FieldCircularOrbit<T> oP3h = new FieldCircularOrbit<>(new FieldPVCoordinates<>(new FieldVector3D<>(1, p, +3, dP), new FieldVector3D<>(1, v, +3, dV)),
1191                                                               orbit.getFrame(), orbit.getDate(), orbit.getMu());
1192         FieldCircularOrbit<T> oP4h = new FieldCircularOrbit<>(new FieldPVCoordinates<>(new FieldVector3D<>(1, p, +4, dP), new FieldVector3D<>(1, v, +4, dV)),
1193                                                               orbit.getFrame(), orbit.getDate(), orbit.getMu());
1194 
1195 
1196         jacobian[0][i] =(zero.add( -3).multiply(oP4h.getA()                            .subtract( oM4h.getA()))  .add(
1197                          zero.add( 32).multiply(oP3h.getA()                            .subtract( oM3h.getA()))) .subtract(
1198                          zero.add(168).multiply(oP2h.getA()                            .subtract( oM2h.getA()))) .add(
1199                          zero.add(672).multiply(oP1h.getA()                            .subtract( oM1h.getA())))).divide(h.multiply(840));
1200         jacobian[1][i] =(zero.add( -3).multiply(oP4h.getCircularEx()                   .subtract( oM4h.getCircularEx()))  .add(
1201                          zero.add( 32).multiply(oP3h.getCircularEx()                   .subtract( oM3h.getCircularEx()))) .subtract(
1202                          zero.add(168).multiply(oP2h.getCircularEx()                   .subtract( oM2h.getCircularEx()))) .add(
1203                          zero.add(672).multiply(oP1h.getCircularEx()                   .subtract( oM1h.getCircularEx())))).divide(h.multiply(840));
1204         jacobian[2][i] =(zero.add( -3).multiply(oP4h.getCircularEy()                   .subtract( oM4h.getCircularEy()))  .add(
1205                          zero.add( 32).multiply(oP3h.getCircularEy()                   .subtract( oM3h.getCircularEy()))) .subtract(
1206                          zero.add(168).multiply(oP2h.getCircularEy()                   .subtract( oM2h.getCircularEy()))) .add(
1207                          zero.add(672).multiply(oP1h.getCircularEy()                   .subtract( oM1h.getCircularEy())))).divide(h.multiply(840));
1208         jacobian[3][i] =(zero.add( -3).multiply(oP4h.getI()                            .subtract( oM4h.getI()))   .add(
1209                          zero.add( 32).multiply(oP3h.getI()                            .subtract( oM3h.getI())) ) .subtract(
1210                          zero.add(168).multiply(oP2h.getI()                            .subtract( oM2h.getI())) ) .add(
1211                          zero.add(672).multiply(oP1h.getI()                            .subtract( oM1h.getI())))).divide(h.multiply(840));
1212         jacobian[4][i] =(zero.add( -3).multiply(oP4h.getRightAscensionOfAscendingNode().subtract( oM4h.getRightAscensionOfAscendingNode()))   .add(
1213                          zero.add( 32).multiply(oP3h.getRightAscensionOfAscendingNode().subtract( oM3h.getRightAscensionOfAscendingNode()))) .subtract(
1214                          zero.add(168).multiply(oP2h.getRightAscensionOfAscendingNode().subtract( oM2h.getRightAscensionOfAscendingNode()))) .add(
1215                          zero.add(672).multiply(oP1h.getRightAscensionOfAscendingNode().subtract( oM1h.getRightAscensionOfAscendingNode())))).divide(h.multiply(840));
1216         jacobian[5][i] =(zero.add( -3).multiply(oP4h.getAlpha(type)                    .subtract( oM4h.getAlpha(type)))   .add(
1217                          zero.add( 32).multiply(oP3h.getAlpha(type)                    .subtract( oM3h.getAlpha(type)))) .subtract(
1218                          zero.add(168).multiply(oP2h.getAlpha(type)                    .subtract( oM2h.getAlpha(type)))) .add(
1219                          zero.add(672).multiply(oP1h.getAlpha(type)                    .subtract( oM1h.getAlpha(type))))).divide(h.multiply(840));
1220 
1221     }
1222 
1223     private <T extends CalculusFieldElement<T>> void doTestNonKeplerianDerivatives(Field<T> field) {
1224         final T zero = field.getZero();
1225 
1226         final FieldAbsoluteDate<T> date         = new FieldAbsoluteDate<>(field, "2003-05-01T00:00:20.000", TimeScalesFactory.getUTC());
1227         final FieldVector3D<T>     position     = new FieldVector3D<>(field.getZero().add(6896874.444705),  field.getZero().add(1956581.072644),  field.getZero().add(-147476.245054));
1228         final FieldVector3D<T>     velocity     = new FieldVector3D<>(field.getZero().add(166.816407662), field.getZero().add(-1106.783301861), field.getZero().add(-7372.745712770));
1229         final FieldVector3D <T>    acceleration = new FieldVector3D<>(field.getZero().add(-7.466182457944), field.getZero().add(-2.118153357345),  field.getZero().add(0.160004048437));
1230         final TimeStampedFieldPVCoordinates<T> pv = new TimeStampedFieldPVCoordinates<>(date, position, velocity, acceleration);
1231         final Frame frame = FramesFactory.getEME2000();
1232         final T mu   = zero.add(Constants.EIGEN5C_EARTH_MU);
1233         final FieldCircularOrbit<T> orbit = new FieldCircularOrbit<>(pv, frame, mu);
1234 
1235         Assertions.assertEquals(differentiate(pv, frame, mu, FieldCircularOrbit::getA),
1236                             orbit.getADot().getReal(),
1237                             4.3e-8);
1238         Assertions.assertEquals(differentiate(pv, frame, mu, FieldCircularOrbit::getEquinoctialEx),
1239                             orbit.getEquinoctialExDot().getReal(),
1240                             2.1e-15);
1241         Assertions.assertEquals(differentiate(pv, frame, mu, FieldCircularOrbit::getEquinoctialEy),
1242                             orbit.getEquinoctialEyDot().getReal(),
1243                             5.4e-16);
1244         Assertions.assertEquals(differentiate(pv, frame, mu, FieldCircularOrbit::getHx),
1245                             orbit.getHxDot().getReal(),
1246                             1.6e-15);
1247         Assertions.assertEquals(differentiate(pv, frame, mu, FieldCircularOrbit::getHy),
1248                             orbit.getHyDot().getReal(),
1249                             7.3e-17);
1250         Assertions.assertEquals(differentiate(pv, frame, mu, FieldCircularOrbit::getLv),
1251                             orbit.getLvDot().getReal(),
1252                             3.4e-16);
1253         Assertions.assertEquals(differentiate(pv, frame, mu, FieldCircularOrbit::getLE),
1254                             orbit.getLEDot().getReal(),
1255                             3.5e-15);
1256         Assertions.assertEquals(differentiate(pv, frame, mu, FieldCircularOrbit::getLM),
1257                             orbit.getLMDot().getReal(),
1258                             5.3e-15);
1259         Assertions.assertEquals(differentiate(pv, frame, mu, FieldCircularOrbit::getE),
1260                             orbit.getEDot().getReal(),
1261                             6.8e-16);
1262         Assertions.assertEquals(differentiate(pv, frame, mu, FieldCircularOrbit::getI),
1263                             orbit.getIDot().getReal(),
1264                             5.7e-16);
1265         Assertions.assertEquals(differentiate(pv, frame, mu, FieldCircularOrbit::getCircularEx),
1266                             orbit.getCircularExDot().getReal(),
1267                             2.2e-15);
1268         Assertions.assertEquals(differentiate(pv, frame, mu, FieldCircularOrbit::getCircularEy),
1269                             orbit.getCircularEyDot().getReal(),
1270                             5.3e-17);
1271         Assertions.assertEquals(differentiate(pv, frame, mu, FieldCircularOrbit::getAlphaV),
1272                             orbit.getAlphaVDot().getReal(),
1273                             4.3e-15);
1274         Assertions.assertEquals(differentiate(pv, frame, mu, FieldCircularOrbit::getAlphaE),
1275                             orbit.getAlphaEDot().getReal(),
1276                             1.2e-15);
1277         Assertions.assertEquals(differentiate(pv, frame, mu, FieldCircularOrbit::getAlphaM),
1278                             orbit.getAlphaMDot().getReal(),
1279                             3.7e-15);
1280         Assertions.assertEquals(differentiate(pv, frame, mu, shifted -> shifted.getAlpha(PositionAngleType.TRUE)),
1281                             orbit.getAlphaDot(PositionAngleType.TRUE).getReal(),
1282                             4.3e-15);
1283         Assertions.assertEquals(differentiate(pv, frame, mu, shifted -> shifted.getAlpha(PositionAngleType.ECCENTRIC)),
1284                             orbit.getAlphaDot(PositionAngleType.ECCENTRIC).getReal(),
1285                             1.2e-15);
1286         Assertions.assertEquals(differentiate(pv, frame, mu, shifted -> shifted.getAlpha(PositionAngleType.MEAN)),
1287                             orbit.getAlphaDot(PositionAngleType.MEAN).getReal(),
1288                             3.7e-15);
1289 
1290     }
1291 
1292     private <T extends CalculusFieldElement<T>, S extends Function<FieldCircularOrbit<T>, T>>
1293     double differentiate(TimeStampedFieldPVCoordinates<T> pv, Frame frame, T mu, S picker) {
1294         final DSFactory factory = new DSFactory(1, 1);
1295         FiniteDifferencesDifferentiator differentiator = new FiniteDifferencesDifferentiator(8, 0.1);
1296         UnivariateDifferentiableFunction diff =
1297                 differentiator.differentiate((UnivariateFunction) dt -> picker.apply(new FieldCircularOrbit<>(pv.shiftedBy(dt),
1298                                                                                                               frame, mu)).getReal());
1299         return diff.value(factory.variable(0, 0.0)).getPartialDerivative(1);
1300      }
1301 
1302     private <T extends CalculusFieldElement<T>> void doTestPositionAngleDerivatives(final Field<T> field) {
1303         final T zero = field.getZero();
1304 
1305         final FieldAbsoluteDate<T> date         = new FieldAbsoluteDate<>(field, "2003-05-01T00:00:20.000", TimeScalesFactory.getUTC());
1306         final FieldVector3D<T>     position     = new FieldVector3D<>(field.getZero().add(6896874.444705),  field.getZero().add(1956581.072644),  field.getZero().add(-147476.245054));
1307         final FieldVector3D<T>     velocity     = new FieldVector3D<>(field.getZero().add(166.816407662), field.getZero().add(-1106.783301861), field.getZero().add(-7372.745712770));
1308         final FieldVector3D <T>    acceleration = new FieldVector3D<>(field.getZero().add(-7.466182457944), field.getZero().add(-2.118153357345),  field.getZero().add(0.160004048437));
1309         final TimeStampedFieldPVCoordinates<T> pv = new TimeStampedFieldPVCoordinates<>(date, position, velocity, acceleration);
1310         final Frame frame = FramesFactory.getEME2000();
1311         final double mu   = Constants.EIGEN5C_EARTH_MU;
1312         final FieldCircularOrbit<T> orbit = new FieldCircularOrbit<>(pv, frame, zero.add(mu));
1313 
1314         for (PositionAngleType type : PositionAngleType.values()) {
1315             final FieldCircularOrbit<T> rebuilt = new FieldCircularOrbit<>(orbit.getA(),
1316                                                                            orbit.getCircularEx(),
1317                                                                            orbit.getCircularEy(),
1318                                                                            orbit.getI(),
1319                                                                            orbit.getRightAscensionOfAscendingNode(),
1320                                                                            orbit.getAlpha(type),
1321                                                                            orbit.getADot(),
1322                                                                            orbit.getCircularExDot(),
1323                                                                            orbit.getCircularEyDot(),
1324                                                                            orbit.getIDot(),
1325                                                                            orbit.getRightAscensionOfAscendingNodeDot(),
1326                                                                            orbit.getAlphaDot(type),
1327                                                                            type, orbit.getFrame(), orbit.getDate(), orbit.getMu());
1328             MatcherAssert.assertThat(rebuilt.getA().getReal(),                                relativelyCloseTo(orbit.getA().getReal(),                                1));
1329             MatcherAssert.assertThat(rebuilt.getCircularEx().getReal(),                       relativelyCloseTo(orbit.getCircularEx().getReal(),                       1));
1330             MatcherAssert.assertThat(rebuilt.getCircularEy().getReal(),                       relativelyCloseTo(orbit.getCircularEy().getReal(),                       1));
1331             MatcherAssert.assertThat(rebuilt.getE().getReal(),                                relativelyCloseTo(orbit.getE().getReal(),                                1));
1332             MatcherAssert.assertThat(rebuilt.getI().getReal(),                                relativelyCloseTo(orbit.getI().getReal(),                                1));
1333             MatcherAssert.assertThat(rebuilt.getRightAscensionOfAscendingNode().getReal(),    relativelyCloseTo(orbit.getRightAscensionOfAscendingNode().getReal(),    1));
1334             MatcherAssert.assertThat(rebuilt.getADot().getReal(),                             relativelyCloseTo(orbit.getADot().getReal(),                             1));
1335             MatcherAssert.assertThat(rebuilt.getCircularExDot().getReal(),                    relativelyCloseTo(orbit.getCircularExDot().getReal(),                    1));
1336             MatcherAssert.assertThat(rebuilt.getCircularEyDot().getReal(),                    relativelyCloseTo(orbit.getCircularEyDot().getReal(),                    1));
1337             MatcherAssert.assertThat(rebuilt.getEDot().getReal(),                             relativelyCloseTo(orbit.getEDot().getReal(),                             1));
1338             MatcherAssert.assertThat(rebuilt.getIDot().getReal(),                             relativelyCloseTo(orbit.getIDot().getReal(),                             1));
1339             MatcherAssert.assertThat(rebuilt.getRightAscensionOfAscendingNodeDot().getReal(), relativelyCloseTo(orbit.getRightAscensionOfAscendingNodeDot().getReal(), 1));
1340             for (PositionAngleType type2 : PositionAngleType.values()) {
1341                 MatcherAssert.assertThat(rebuilt.getAlpha(type2).getReal(),    relativelyCloseTo(orbit.getAlpha(type2).getReal(),    1));
1342                 MatcherAssert.assertThat(rebuilt.getAlphaDot(type2).getReal(), relativelyCloseTo(orbit.getAlphaDot(type2).getReal(), 1));
1343             }
1344         }
1345 
1346     }
1347 
1348     private <T extends CalculusFieldElement<T>> void doTestEquatorialRetrograde(final Field<T> field) {
1349         T zero = field.getZero();
1350         FieldVector3D<T> position = new FieldVector3D<>(field.getZero().add(10000000.0), field.getZero(), field.getZero());
1351         FieldVector3D<T> velocity = new FieldVector3D<>(field.getZero(), field.getZero().add(-6500.0), field.getZero().add(1.0e-10));
1352         T r2 = position.getNormSq();
1353         T r  = r2.sqrt();
1354         FieldVector3D<T> acceleration = new FieldVector3D<>(r.multiply(r2.reciprocal().multiply(-mu)), position,
1355                                                             field.getOne(), new FieldVector3D<>(field.getZero().add(-0.1),
1356                                                                                                 field.getZero().add(0.2),
1357                                                                                                 field.getZero().add(0.3)));
1358         FieldPVCoordinates<T> pvCoordinates = new FieldPVCoordinates<>(position, velocity, acceleration);
1359         FieldCircularOrbit<T> orbit = new FieldCircularOrbit<>(pvCoordinates, FramesFactory.getEME2000(),
1360                                                                FieldAbsoluteDate.getJ2000Epoch(field), zero.add(mu));
1361         Assertions.assertEquals(10637829.465, orbit.getA().getReal(), 1.0e-3);
1362         Assertions.assertEquals(-738.145, orbit.getADot().getReal(), 1.0e-3);
1363         Assertions.assertEquals(0.05995861, orbit.getE().getReal(), 1.0e-8);
1364         Assertions.assertEquals(-6.523e-5, orbit.getEDot().getReal(), 1.0e-8);
1365         Assertions.assertEquals(FastMath.PI, orbit.getI().getReal(), 2.0e-14);
1366         Assertions.assertEquals(-4.615e-5, orbit.getIDot().getReal(), 1.0e-8);
1367         Assertions.assertTrue(Double.isNaN(orbit.getHx().getReal()));
1368         Assertions.assertTrue(Double.isNaN(orbit.getHxDot().getReal()));
1369         Assertions.assertTrue(Double.isNaN(orbit.getHy().getReal()));
1370         Assertions.assertTrue(Double.isNaN(orbit.getHyDot().getReal()));
1371     }
1372 
1373     private <T extends CalculusFieldElement<T>> void doTestDerivativesConversionSymmetry(Field<T> field) {
1374         T zero = field.getZero();
1375         final FieldAbsoluteDate<T> date = new FieldAbsoluteDate<>(field, "2003-05-01T00:01:20.000", TimeScalesFactory.getUTC());
1376         FieldVector3D<T> position     = new FieldVector3D<>(zero.add(6893443.400234382),
1377                                                             zero.add(1886406.1073757345),
1378                                                             zero.add(-589265.1150359757));
1379         FieldVector3D<T> velocity     = new FieldVector3D<>(zero.add(-281.1261461082365),
1380                                                             zero.add(-1231.6165642450928),
1381                                                             zero.add(-7348.756363469432));
1382         FieldVector3D<T> acceleration = new FieldVector3D<>(zero.add(-7.460341170581685),
1383                                                             zero.add(-2.0415957334584527),
1384                                                             zero.add(0.6393322823627762));
1385         FieldPVCoordinates<T> pvCoordinates = new FieldPVCoordinates<>( position, velocity, acceleration);
1386         FieldCircularOrbit<T> orbit = new FieldCircularOrbit<>(pvCoordinates, FramesFactory.getEME2000(),
1387                                                                date, zero.add(Constants.EIGEN5C_EARTH_MU));
1388         Assertions.assertTrue(orbit.hasNonKeplerianAcceleration());
1389         T r2 = position.getNormSq();
1390         T r  = r2.sqrt();
1391         FieldVector3D<T> keplerianAcceleration = new FieldVector3D<>(r.multiply(r2).reciprocal().multiply(orbit.getMu().negate()),
1392                                                                      position);
1393         Assertions.assertEquals(0.0101, FieldVector3D.distance(keplerianAcceleration, acceleration).getReal(), 1.0e-4);
1394 
1395         for (OrbitType type : OrbitType.values()) {
1396             FieldOrbit<T> converted = type.convertType(orbit);
1397             Assertions.assertTrue(converted.hasNonKeplerianAcceleration());
1398             FieldCircularOrbit<T> rebuilt = (FieldCircularOrbit<T>) OrbitType.CIRCULAR.convertType(converted);
1399             Assertions.assertTrue(rebuilt.hasNonKeplerianAcceleration());
1400             Assertions.assertEquals(orbit.getADot().getReal(),                             rebuilt.getADot().getReal(),                             3.0e-13);
1401             Assertions.assertEquals(orbit.getCircularExDot().getReal(),                    rebuilt.getCircularExDot().getReal(),                    1.0e-15);
1402             Assertions.assertEquals(orbit.getCircularEyDot().getReal(),                    rebuilt.getCircularEyDot().getReal(),                    1.0e-15);
1403             Assertions.assertEquals(orbit.getIDot().getReal(),                             rebuilt.getIDot().getReal(),                             1.0e-15);
1404             Assertions.assertEquals(orbit.getRightAscensionOfAscendingNodeDot().getReal(), rebuilt.getRightAscensionOfAscendingNodeDot().getReal(), 1.0e-15);
1405             Assertions.assertEquals(orbit.getAlphaVDot().getReal(),                        rebuilt.getAlphaVDot().getReal(),                        1.0e-15);
1406         }
1407 
1408     }
1409 
1410     private <T extends CalculusFieldElement<T>> void doTestToString(Field<T> field) {
1411         T zero = field.getZero();
1412         FieldVector3D<T> position = new FieldVector3D<>(field.getZero().add(-29536113.0),
1413                                                         field.getZero().add(30329259.0),
1414                                                         field.getZero().add(-100125.0));
1415         FieldVector3D<T> velocity = new FieldVector3D<>(field.getZero().add(-2194.0),
1416                                                         field.getZero().add(-2141.0),
1417                                                         field.getZero().add(-8.0));
1418         FieldPVCoordinates<T> pvCoordinates = new FieldPVCoordinates<>(position, velocity);
1419         FieldCircularOrbit<T> orbit = new FieldCircularOrbit<>(pvCoordinates, FramesFactory.getEME2000(),
1420                                                                FieldAbsoluteDate.getJ2000Epoch(field), zero.add(mu));
1421         Assertions.assertEquals("circular parameters: {a: 4.225517000282565E7, ex: 0.002082917137146049, ey: 5.173980074371024E-4, i: 0.20189257051515358, raan: -87.91788415673473, alphaV: -137.84099636616548;}",
1422                             orbit.toString());
1423     }
1424 
1425     private <T extends CalculusFieldElement<T>> void doTestCopyNonKeplerianAcceleration(Field<T> field) {
1426 
1427         final T zero = field.getZero();
1428         final Frame eme2000     = FramesFactory.getEME2000();
1429 
1430         // Define GEO satellite position
1431         final FieldVector3D<T> position = new FieldVector3D<>(field.getZero().add(42164140),
1432                         field.getZero(),
1433                         field.getZero());
1434         // Build PVCoodrinates starting from its position and computing the corresponding circular velocity
1435         final FieldPVCoordinates<T> pv  =
1436                         new FieldPVCoordinates<>(position,
1437                                         new FieldVector3D<>(field.getZero(),
1438                                                         position.getNorm().reciprocal().multiply(mu).sqrt(),
1439                                                         field.getZero()));
1440         // Build a KeplerianOrbit in eme2000
1441         final FieldOrbit<T> orbit = new FieldCircularOrbit<>(pv, eme2000, FieldAbsoluteDate.getJ2000Epoch(field), zero.add(mu));
1442 
1443         // Build another KeplerianOrbit as a copy of the first one
1444         final FieldOrbit<T> orbitCopy = new FieldCircularOrbit<>(orbit);
1445 
1446         // Shift the orbit of a time-interval
1447         final FieldOrbit<T> shiftedOrbit     = orbit.shiftedBy(10); // This works good
1448         final FieldOrbit<T> shiftedOrbitCopy = orbitCopy.shiftedBy(10); // This does not work
1449 
1450         Assertions.assertEquals(0.0,
1451                             FieldVector3D.distance(shiftedOrbit.getPosition(),
1452                                                    shiftedOrbitCopy.getPosition()).getReal(),
1453                             1.0e-10);
1454         Assertions.assertEquals(0.0,
1455                             FieldVector3D.distance(shiftedOrbit.getVelocity(),
1456                                                    shiftedOrbitCopy.getVelocity()).getReal(),
1457                             1.0e-10);
1458 
1459     }
1460 
1461     private <T extends CalculusFieldElement<T>> void doTestNormalize(Field<T> field) {
1462         final T zero = field.getZero();
1463         FieldCircularOrbit<T> withoutDerivatives =
1464                         new FieldCircularOrbit<>(zero.newInstance(42166712.0), zero.newInstance(0.005),
1465                                                  zero.newInstance(-0.025), zero.newInstance(1.6),
1466                                                  zero.newInstance(1.25), zero.newInstance(0.4), PositionAngleType.MEAN,
1467                                                  FramesFactory.getEME2000(), FieldAbsoluteDate.getJ2000Epoch(field),
1468                                                  zero.newInstance(mu));
1469         FieldCircularOrbit<T> ref =
1470                         new FieldCircularOrbit<>(zero.newInstance(24000000.0), zero.newInstance(-0.012),
1471                                                  zero.newInstance(0.01), zero.newInstance(0.2),
1472                                                  zero.newInstance(-6.28), zero.newInstance(6.28), PositionAngleType.MEAN,
1473                                                  FramesFactory.getEME2000(), FieldAbsoluteDate.getJ2000Epoch(field),
1474                                                  zero.newInstance(mu));
1475 
1476         FieldCircularOrbit<T> normalized1 = (FieldCircularOrbit<T>) OrbitType.CIRCULAR.normalize(withoutDerivatives, ref);
1477         Assertions.assertFalse(normalized1.hasNonKeplerianAcceleration());
1478         Assertions.assertEquals(0.0, normalized1.getA().subtract(withoutDerivatives.getA()).getReal(),          1.0e-6);
1479         Assertions.assertEquals(0.0, normalized1.getCircularEx().subtract(withoutDerivatives.getCircularEx()).getReal(), 1.0e-10);
1480         Assertions.assertEquals(0.0, normalized1.getCircularEy().subtract(withoutDerivatives.getCircularEy()).getReal(), 1.0e-10);
1481         Assertions.assertEquals(0.0, normalized1.getI().subtract(withoutDerivatives.getI()).getReal(),          1.0e-10);
1482         Assertions.assertEquals(-MathUtils.TWO_PI, normalized1.getRightAscensionOfAscendingNode().subtract(withoutDerivatives.getRightAscensionOfAscendingNode()).getReal(), 1.0e-10);
1483         Assertions.assertEquals(+MathUtils.TWO_PI, normalized1.getAlphaV().subtract(withoutDerivatives.getAlphaV()).getReal(), 1.0e-10);
1484 
1485         T[] p    = MathArrays.buildArray(field, 6);
1486         T[] pDot = MathArrays.buildArray(field, 6);
1487         for (int i = 0; i < pDot.length; ++i) {
1488             pDot[i] = zero.newInstance(i);
1489         }
1490         OrbitType.CIRCULAR.mapOrbitToArray(withoutDerivatives, PositionAngleType.TRUE, p, null);
1491         FieldCircularOrbit<T> withDerivatives = (FieldCircularOrbit<T>) OrbitType.CIRCULAR.mapArrayToOrbit(p, pDot,
1492                                                                                                            PositionAngleType.TRUE,
1493                                                                                                            withoutDerivatives.getDate(),
1494                                                                                                            withoutDerivatives.getMu(),
1495                                                                                                            withoutDerivatives.getFrame());
1496         FieldCircularOrbit<T> normalized2 = (FieldCircularOrbit<T>) OrbitType.CIRCULAR.normalize(withDerivatives, ref);
1497         Assertions.assertTrue(normalized2.hasNonKeplerianAcceleration());
1498         Assertions.assertEquals(0.0, normalized2.getA().subtract(withDerivatives.getA()).getReal(),          1.0e-6);
1499         Assertions.assertEquals(0.0, normalized2.getCircularEx().subtract(withDerivatives.getCircularEx()).getReal(), 1.0e-10);
1500         Assertions.assertEquals(0.0, normalized2.getCircularEy().subtract(withDerivatives.getCircularEy()).getReal(), 1.0e-10);
1501         Assertions.assertEquals(0.0, normalized2.getI().subtract(withDerivatives.getI()).getReal(),          1.0e-10);
1502         Assertions.assertEquals(-MathUtils.TWO_PI, normalized2.getRightAscensionOfAscendingNode().subtract(withDerivatives.getRightAscensionOfAscendingNode()).getReal(), 1.0e-10);
1503         Assertions.assertEquals(+MathUtils.TWO_PI, normalized2.getAlphaV().subtract(withDerivatives.getAlphaV()).getReal(), 1.0e-10);
1504         Assertions.assertEquals(0.0, normalized2.getADot().subtract(withDerivatives.getADot()).getReal(),          1.0e-6);
1505         Assertions.assertEquals(0.0, normalized2.getCircularExDot().subtract(withDerivatives.getCircularExDot()).getReal(), 1.0e-10);
1506         Assertions.assertEquals(0.0, normalized2.getCircularEyDot().subtract(withDerivatives.getCircularEyDot()).getReal(), 1.0e-10);
1507         Assertions.assertEquals(0.0, normalized2.getIDot().subtract(withDerivatives.getIDot()).getReal(),          1.0e-10);
1508         Assertions.assertEquals(0.0, normalized2.getRightAscensionOfAscendingNodeDot().subtract(withDerivatives.getRightAscensionOfAscendingNodeDot()).getReal(), 1.0e-10);
1509         Assertions.assertEquals(0.0, normalized2.getAlphaVDot().subtract(withDerivatives.getAlphaVDot()).getReal(), 1.0e-10);
1510 
1511     }
1512 
1513 }