1   /* Copyright 2002-2022 CS GROUP
2    * Licensed to CS GROUP (CS) under one or more
3    * contributor license agreements.  See the NOTICE file distributed with
4    * this work for additional information regarding copyright ownership.
5    * CS licenses this file to You under the Apache License, Version 2.0
6    * (the "License"); you may not use this file except in compliance with
7    * the License.  You may obtain a copy of the License at
8    *
9    *   http://www.apache.org/licenses/LICENSE-2.0
10   *
11   * Unless required by applicable law or agreed to in writing, software
12   * distributed under the License is distributed on an "AS IS" BASIS,
13   * WITHOUT WARRANTIES OR CONDITIONS OF ANY KIND, either express or implied.
14   * See the License for the specific language governing permissions and
15   * limitations under the License.
16   */
17  package org.orekit.propagation.analytical;
18  
19  
20  import java.util.ArrayList;
21  import java.util.List;
22  
23  import org.hamcrest.MatcherAssert;
24  import org.hipparchus.CalculusFieldElement;
25  import org.hipparchus.exception.DummyLocalizable;
26  import org.hipparchus.exception.LocalizedCoreFormats;
27  import org.hipparchus.geometry.euclidean.threed.Vector3D;
28  import org.hipparchus.util.FastMath;
29  import org.hipparchus.util.MathUtils;
30  import org.junit.After;
31  import org.junit.Assert;
32  import org.junit.Before;
33  import org.junit.Test;
34  import org.orekit.OrekitMatchers;
35  import org.orekit.Utils;
36  import org.orekit.attitudes.Attitude;
37  import org.orekit.attitudes.AttitudeProvider;
38  import org.orekit.attitudes.FieldAttitude;
39  import org.orekit.attitudes.LofOffset;
40  import org.orekit.bodies.BodyShape;
41  import org.orekit.bodies.GeodeticPoint;
42  import org.orekit.bodies.OneAxisEllipsoid;
43  import org.orekit.errors.OrekitException;
44  import org.orekit.forces.maneuvers.ImpulseManeuver;
45  import org.orekit.frames.Frame;
46  import org.orekit.frames.FramesFactory;
47  import org.orekit.frames.LOFType;
48  import org.orekit.frames.TopocentricFrame;
49  import org.orekit.orbits.CartesianOrbit;
50  import org.orekit.orbits.CircularOrbit;
51  import org.orekit.orbits.EquinoctialOrbit;
52  import org.orekit.orbits.KeplerianOrbit;
53  import org.orekit.orbits.Orbit;
54  import org.orekit.orbits.OrbitType;
55  import org.orekit.orbits.PositionAngle;
56  import org.orekit.propagation.AdditionalStateProvider;
57  import org.orekit.propagation.BoundedPropagator;
58  import org.orekit.propagation.EphemerisGenerator;
59  import org.orekit.propagation.Propagator;
60  import org.orekit.propagation.SpacecraftState;
61  import org.orekit.propagation.events.AbstractDetector;
62  import org.orekit.propagation.events.AltitudeDetector;
63  import org.orekit.propagation.events.ApsideDetector;
64  import org.orekit.propagation.events.DateDetector;
65  import org.orekit.propagation.events.ElevationDetector;
66  import org.orekit.propagation.events.NodeDetector;
67  import org.orekit.propagation.events.handlers.ContinueOnEvent;
68  import org.orekit.propagation.sampling.OrekitFixedStepHandler;
69  import org.orekit.propagation.sampling.OrekitStepHandler;
70  import org.orekit.propagation.sampling.OrekitStepInterpolator;
71  import org.orekit.time.AbsoluteDate;
72  import org.orekit.time.DateComponents;
73  import org.orekit.time.FieldAbsoluteDate;
74  import org.orekit.time.TimeComponents;
75  import org.orekit.time.TimeScale;
76  import org.orekit.time.TimeScalesFactory;
77  import org.orekit.utils.Constants;
78  import org.orekit.utils.FieldPVCoordinatesProvider;
79  import org.orekit.utils.IERSConventions;
80  import org.orekit.utils.PVCoordinates;
81  import org.orekit.utils.PVCoordinatesProvider;
82  import org.orekit.utils.TimeStampedPVCoordinates;
83  
84  
85  public class KeplerianPropagatorTest {
86  
87      // Body mu
88      private double mu;
89  
90      /**
91       * Check that the date returned by {@link KeplerianPropagator#propagate(AbsoluteDate)}
92       * is the same as the date passed to propagate().
93       */
94      @Test
95      public void testPropagationDate() {
96          // setup
97          AbsoluteDate initDate = AbsoluteDate.J2000_EPOCH;
98          // date s.t. target - date rounds down when represented as a double.
99          AbsoluteDate target =
100                 initDate.shiftedBy(20.0).shiftedBy(FastMath.ulp(20.0) / 4);
101         Orbit ic = new KeplerianOrbit(6378137 + 500e3, 1e-3, 0, 0, 0, 0,
102                 PositionAngle.TRUE, FramesFactory.getGCRF(), initDate, mu);
103         Propagator propagator = new KeplerianPropagator(ic);
104 
105         // action
106         SpacecraftState actual = propagator.propagate(target);
107 
108         // verify
109         Assert.assertEquals(target, actual.getDate());
110     }
111 
112     @Test
113     public void testEphemerisModeWithHandler() {
114         // setup
115         AbsoluteDate initDate = AbsoluteDate.GPS_EPOCH;
116         Orbit ic = new KeplerianOrbit(6378137 + 500e3, 1e-3, 0, 0, 0, 0,
117                 PositionAngle.TRUE, FramesFactory.getGCRF(), initDate, mu);
118         Propagator propagator = new KeplerianPropagator(ic);
119         AbsoluteDate end = initDate.shiftedBy(90 * 60);
120 
121         // action
122         final List<SpacecraftState> states = new ArrayList<>();
123         final EphemerisGenerator generator = propagator.getEphemerisGenerator();
124         propagator.setStepHandler(interpolator -> {
125             states.add(interpolator.getCurrentState());
126             states.add(interpolator.getPreviousState());
127         });
128         propagator.propagate(end);
129         final BoundedPropagator ephemeris = generator.getGeneratedEphemeris();
130 
131         //verify
132         Assert.assertTrue(states.size() > 1); // got some data
133         for (SpacecraftState state : states) {
134             PVCoordinates actual =
135                     ephemeris.propagate(state.getDate()).getPVCoordinates();
136             MatcherAssert.assertThat(actual, OrekitMatchers.pvIs(state.getPVCoordinates()));
137         }
138     }
139 
140     @Test
141     public void testAdditionalState() {
142         AbsoluteDate initDate = AbsoluteDate.GPS_EPOCH;
143         Orbit ic = new KeplerianOrbit(6378137 + 500e3, 1e-3, 0, 0, 0, 0, PositionAngle.TRUE, FramesFactory.getGCRF(), initDate, mu);
144         Propagator propagator = new KeplerianPropagator(ic);
145         SpacecraftState initialState = propagator.getInitialState().addAdditionalState("myState", 4.2);
146         propagator.resetInitialState(initialState);
147         AbsoluteDate end = initDate.shiftedBy(90 * 60);
148         EphemerisGenerator generator = propagator.getEphemerisGenerator();
149         SpacecraftState finalStateKeplerianPropagator = propagator.propagate(end);
150         BoundedPropagator ephemeris = generator.getGeneratedEphemeris();
151         SpacecraftState ephemerisInitialState = ephemeris.getInitialState();
152         SpacecraftState finalStateBoundedPropagator = ephemeris.propagate(end);
153         Assert.assertEquals(4.2, finalStateKeplerianPropagator.getAdditionalState("myState")[0], 1.0e-15);
154         Assert.assertEquals(4.2, ephemerisInitialState.getAdditionalState("myState")[0], 1.0e-15);
155         Assert.assertEquals(4.2, finalStateBoundedPropagator.getAdditionalState("myState")[0], 1.0e-15);
156     }
157 
158     @Test
159     public void sameDateCartesian() {
160 
161         // Definition of initial conditions with position and velocity
162         //------------------------------------------------------------
163         Vector3D position = new Vector3D(7.0e6, 1.0e6, 4.0e6);
164         Vector3D velocity = new Vector3D(-500.0, 8000.0, 1000.0);
165 
166         AbsoluteDate initDate = AbsoluteDate.J2000_EPOCH.shiftedBy(584.);
167         Orbit initialOrbit = new EquinoctialOrbit(new PVCoordinates(position, velocity),
168                                                   FramesFactory.getEME2000(), initDate, mu);
169 
170         // Extrapolator definition
171         // -----------------------
172         KeplerianPropagator extrapolator = new KeplerianPropagator(initialOrbit);
173 
174         // Extrapolation at the initial date
175         // ---------------------------------
176         double delta_t = 0.0; // extrapolation duration in seconds
177         AbsoluteDate extrapDate = initDate.shiftedBy(delta_t);
178 
179         SpacecraftState finalOrbit = extrapolator.propagate(extrapDate);
180 
181         double a = finalOrbit.getA();
182         // another way to compute n
183         double n = FastMath.sqrt(finalOrbit.getMu()/FastMath.pow(a, 3));
184 
185         Assert.assertEquals(n*delta_t,
186                             finalOrbit.getLM() - initialOrbit.getLM(),
187                             Utils.epsilonTest * FastMath.abs(n*delta_t));
188         Assert.assertEquals(MathUtils.normalizeAngle(finalOrbit.getLM(), initialOrbit.getLM()), initialOrbit.getLM(), Utils.epsilonAngle * FastMath.abs(initialOrbit.getLM()));
189 
190         Assert.assertEquals(finalOrbit.getA(), initialOrbit.getA(), Utils.epsilonTest * initialOrbit.getA());
191         Assert.assertEquals(finalOrbit.getE(), initialOrbit.getE(), Utils.epsilonE * initialOrbit.getE());
192         Assert.assertEquals(MathUtils.normalizeAngle(finalOrbit.getI(), initialOrbit.getI()), initialOrbit.getI(), Utils.epsilonAngle * FastMath.abs(initialOrbit.getI()));
193 
194     }
195 
196     @Test
197     public void sameDateKeplerian() {
198         // Definition of initial conditions with Keplerian parameters
199         //-----------------------------------------------------------
200         AbsoluteDate initDate = AbsoluteDate.J2000_EPOCH.shiftedBy(584.);
201         Orbit initialOrbit = new KeplerianOrbit(7209668.0, 0.5e-4, 1.7, 2.1, 2.9,
202                                                 6.2, PositionAngle.TRUE,
203                                                 FramesFactory.getEME2000(), initDate, mu);
204 
205         // Extrapolator definition
206         // -----------------------
207         KeplerianPropagator extrapolator = new KeplerianPropagator(initialOrbit);
208 
209         // Extrapolation at the initial date
210         // ---------------------------------
211         double delta_t = 0.0; // extrapolation duration in seconds
212         AbsoluteDate extrapDate = initDate.shiftedBy(delta_t);
213 
214         SpacecraftState finalOrbit = extrapolator.propagate(extrapDate);
215 
216         double a = finalOrbit.getA();
217         // another way to compute n
218         double n = FastMath.sqrt(finalOrbit.getMu()/FastMath.pow(a, 3));
219 
220         Assert.assertEquals(n*delta_t,
221                      finalOrbit.getLM() - initialOrbit.getLM(),
222                      Utils.epsilonTest * FastMath.max(100., FastMath.abs(n*delta_t)));
223         Assert.assertEquals(MathUtils.normalizeAngle(finalOrbit.getLM(), initialOrbit.getLM()), initialOrbit.getLM(), Utils.epsilonAngle * FastMath.abs(initialOrbit.getLM()));
224 
225         Assert.assertEquals(finalOrbit.getA(), initialOrbit.getA(), Utils.epsilonTest * initialOrbit.getA());
226         Assert.assertEquals(finalOrbit.getE(), initialOrbit.getE(), Utils.epsilonE * initialOrbit.getE());
227         Assert.assertEquals(MathUtils.normalizeAngle(finalOrbit.getI(), initialOrbit.getI()), initialOrbit.getI(), Utils.epsilonAngle * FastMath.abs(initialOrbit.getI()));
228 
229     }
230 
231     @Test
232     public void propagatedCartesian() {
233 
234         // Definition of initial conditions with position and velocity
235         //------------------------------------------------------------
236         Vector3D position = new Vector3D(7.0e6, 1.0e6, 4.0e6);
237         Vector3D velocity = new Vector3D(-500.0, 8000.0, 1000.0);
238         double mu = 3.9860047e14;
239 
240         AbsoluteDate initDate = AbsoluteDate.J2000_EPOCH.shiftedBy(584.);
241         Orbit initialOrbit = new EquinoctialOrbit(new PVCoordinates(position, velocity),
242                                                   FramesFactory.getEME2000(), initDate, mu);
243 
244         // Extrapolator definition
245         // -----------------------
246         KeplerianPropagator extrapolator = new KeplerianPropagator(initialOrbit);
247 
248         // Extrapolation at a final date different from initial date
249         // ---------------------------------------------------------
250         double delta_t = 100000.0; // extrapolation duration in seconds
251         AbsoluteDate extrapDate = initDate.shiftedBy(delta_t);
252 
253         SpacecraftState finalOrbit = extrapolator.propagate(extrapDate);
254 
255 
256         // computation of (M final - M initial) with another method
257         double a = finalOrbit.getA();
258         // another way to compute n
259         double n = FastMath.sqrt(finalOrbit.getMu()/FastMath.pow(a, 3));
260 
261         Assert.assertEquals(n * delta_t,
262                             finalOrbit.getLM() - initialOrbit.getLM(),
263                             Utils.epsilonAngle);
264 
265         // computation of M final orbit
266         double LM = finalOrbit.getLE()
267         - finalOrbit.getEquinoctialEx()*FastMath.sin(finalOrbit.getLE())
268         + finalOrbit.getEquinoctialEy()*FastMath.cos(finalOrbit.getLE());
269 
270         Assert.assertEquals(LM , finalOrbit.getLM() , Utils.epsilonAngle);
271 
272         // test of tan ((LE - Lv)/2) :
273         Assert.assertEquals(FastMath.tan((finalOrbit.getLE() - finalOrbit.getLv())/2.),
274                      tangLEmLv(finalOrbit.getLv(), finalOrbit.getEquinoctialEx(), finalOrbit.getEquinoctialEy()),
275                      Utils.epsilonAngle);
276 
277         // test of evolution of M vs E: LM = LE - ex*sin(LE) + ey*cos(LE)
278         // with ex and ey the same for initial and final orbit
279         double deltaM = finalOrbit.getLM() - initialOrbit.getLM();
280         double deltaE = finalOrbit.getLE() - initialOrbit.getLE();
281         double delta  = finalOrbit.getEquinoctialEx() * (FastMath.sin(finalOrbit.getLE()) - FastMath.sin(initialOrbit.getLE()))
282         - finalOrbit.getEquinoctialEy() * (FastMath.cos(finalOrbit.getLE()) - FastMath.cos(initialOrbit.getLE()));
283 
284         Assert.assertEquals(deltaM, deltaE - delta, Utils.epsilonAngle);
285 
286         // the orbital elements except for Mean/True/Eccentric latitude arguments are the same
287         Assert.assertEquals(finalOrbit.getA(), initialOrbit.getA(), Utils.epsilonTest * initialOrbit.getA());
288         Assert.assertEquals(finalOrbit.getEquinoctialEx(), initialOrbit.getEquinoctialEx(), Utils.epsilonE);
289         Assert.assertEquals(finalOrbit.getEquinoctialEy(), initialOrbit.getEquinoctialEy(), Utils.epsilonE);
290         Assert.assertEquals(finalOrbit.getHx(), initialOrbit.getHx(), Utils.epsilonAngle);
291         Assert.assertEquals(finalOrbit.getHy(), initialOrbit.getHy(), Utils.epsilonAngle);
292 
293         // for final orbit
294         double ex = finalOrbit.getEquinoctialEx();
295         double ey = finalOrbit.getEquinoctialEy();
296         double hx = finalOrbit.getHx();
297         double hy = finalOrbit.getHy();
298         double LE = finalOrbit.getLE();
299 
300         double ex2 = ex*ex;
301         double ey2 = ey*ey;
302         double hx2 = hx*hx;
303         double hy2 = hy*hy;
304         double h2p1 = 1. + hx2 + hy2;
305         double beta = 1. / (1. + FastMath.sqrt(1. - ex2 - ey2));
306 
307         double x3 = -ex + (1.- beta*ey2)*FastMath.cos(LE) + beta*ex*ey*FastMath.sin(LE);
308         double y3 = -ey + (1. -beta*ex2)*FastMath.sin(LE) + beta*ex*ey*FastMath.cos(LE);
309 
310         Vector3D U = new Vector3D((1. + hx2 - hy2)/ h2p1,
311                                   (2.*hx*hy)/h2p1,
312                                   (-2.*hy)/h2p1);
313 
314         Vector3D V = new Vector3D((2.*hx*hy)/ h2p1,
315                                   (1.- hx2+ hy2)/h2p1,
316                                   (2.*hx)/h2p1);
317 
318         Vector3D r = new Vector3D(finalOrbit.getA(), new Vector3D(x3, U, y3, V));
319 
320         Assert.assertEquals(finalOrbit.getPVCoordinates().getPosition().getNorm(), r.getNorm(), Utils.epsilonTest * r.getNorm());
321 
322     }
323 
324     @Test
325     public void propagatedKeplerian() {
326 
327         // Definition of initial conditions with Keplerian parameters
328         //-----------------------------------------------------------
329         AbsoluteDate initDate = AbsoluteDate.J2000_EPOCH.shiftedBy(584.);
330         Orbit initialOrbit = new KeplerianOrbit(7209668.0, 0.5e-4, 1.7, 2.1, 2.9,
331                                                 6.2, PositionAngle.TRUE,
332                                                 FramesFactory.getEME2000(), initDate, mu);
333 
334         // Extrapolator definition
335         // -----------------------
336         KeplerianPropagator extrapolator = new KeplerianPropagator(initialOrbit);
337 
338         // Extrapolation at a final date different from initial date
339         // ---------------------------------------------------------
340         double delta_t = 100000.0; // extrapolation duration in seconds
341         AbsoluteDate extrapDate = initDate.shiftedBy(delta_t);
342 
343         SpacecraftState finalOrbit = extrapolator.propagate(extrapDate);
344         Assert.assertEquals(6092.3362422560844633, finalOrbit.getKeplerianPeriod(), 1.0e-12);
345         Assert.assertEquals(0.001031326088602888358, finalOrbit.getKeplerianMeanMotion(), 1.0e-16);
346 
347         // computation of (M final - M initial) with another method
348         double a = finalOrbit.getA();
349         // another way to compute n
350         double n = FastMath.sqrt(finalOrbit.getMu()/FastMath.pow(a, 3));
351 
352         Assert.assertEquals(n * delta_t,
353                      finalOrbit.getLM() - initialOrbit.getLM(),
354                      Utils.epsilonAngle);
355 
356         // computation of M final orbit
357         double LM = finalOrbit.getLE()
358         - finalOrbit.getEquinoctialEx()*FastMath.sin(finalOrbit.getLE())
359         + finalOrbit.getEquinoctialEy()*FastMath.cos(finalOrbit.getLE());
360 
361         Assert.assertEquals(LM , finalOrbit.getLM() , Utils.epsilonAngle);
362 
363         // test of tan ((LE - Lv)/2) :
364         Assert.assertEquals(FastMath.tan((finalOrbit.getLE() - finalOrbit.getLv())/2.),
365                      tangLEmLv(finalOrbit.getLv(), finalOrbit.getEquinoctialEx(), finalOrbit.getEquinoctialEy()),
366                      Utils.epsilonAngle);
367 
368         // test of evolution of M vs E: LM = LE - ex*sin(LE) + ey*cos(LE)
369         // with ex and ey the same for initial and final orbit
370         double deltaM = finalOrbit.getLM() - initialOrbit.getLM();
371         double deltaE = finalOrbit.getLE() - initialOrbit.getLE();
372         double delta  = finalOrbit.getEquinoctialEx() * (FastMath.sin(finalOrbit.getLE()) - FastMath.sin(initialOrbit.getLE())) - finalOrbit.getEquinoctialEy() * (FastMath.cos(finalOrbit.getLE()) - FastMath.cos(initialOrbit.getLE()));
373 
374         Assert.assertEquals(deltaM, deltaE - delta, Utils.epsilonAngle);
375 
376         // the orbital elements except for Mean/True/Eccentric latitude arguments are the same
377         Assert.assertEquals(finalOrbit.getA(), initialOrbit.getA(), Utils.epsilonTest * initialOrbit.getA());
378         Assert.assertEquals(finalOrbit.getEquinoctialEx(), initialOrbit.getEquinoctialEx(), Utils.epsilonE);
379         Assert.assertEquals(finalOrbit.getEquinoctialEy(), initialOrbit.getEquinoctialEy(), Utils.epsilonE);
380         Assert.assertEquals(finalOrbit.getHx(), initialOrbit.getHx(), Utils.epsilonAngle);
381         Assert.assertEquals(finalOrbit.getHy(), initialOrbit.getHy(), Utils.epsilonAngle);
382 
383         // for final orbit
384         double ex = finalOrbit.getEquinoctialEx();
385         double ey = finalOrbit.getEquinoctialEy();
386         double hx = finalOrbit.getHx();
387         double hy = finalOrbit.getHy();
388         double LE = finalOrbit.getLE();
389 
390         double ex2 = ex*ex;
391         double ey2 = ey*ey;
392         double hx2 = hx*hx;
393         double hy2 = hy*hy;
394         double h2p1 = 1. + hx2 + hy2;
395         double beta = 1. / (1. + FastMath.sqrt(1. - ex2 - ey2));
396 
397         double x3 = -ex + (1.- beta*ey2)*FastMath.cos(LE) + beta*ex*ey*FastMath.sin(LE);
398         double y3 = -ey + (1. -beta*ex2)*FastMath.sin(LE) + beta*ex*ey*FastMath.cos(LE);
399 
400         Vector3D U = new Vector3D((1. + hx2 - hy2)/ h2p1,
401                                   (2.*hx*hy)/h2p1,
402                                   (-2.*hy)/h2p1);
403 
404         Vector3D V = new Vector3D((2.*hx*hy)/ h2p1,
405                                   (1.- hx2+ hy2)/h2p1,
406                                   (2.*hx)/h2p1);
407 
408         Vector3D r = new Vector3D(finalOrbit.getA(), new Vector3D(x3, U, y3, V));
409 
410         Assert.assertEquals(finalOrbit.getPVCoordinates().getPosition().getNorm(), r.getNorm(), Utils.epsilonTest * r.getNorm());
411 
412     }
413 
414     @Test(expected = OrekitException.class)
415     public void wrongAttitude() {
416         KeplerianOrbit orbit =
417             new KeplerianOrbit(1.0e10, 1.0e-4, 1.0e-2, 0, 0, 0, PositionAngle.TRUE,
418                                FramesFactory.getEME2000(), AbsoluteDate.J2000_EPOCH, 3.986004415e14);
419         AttitudeProvider wrongLaw = new AttitudeProvider() {
420             public Attitude getAttitude(PVCoordinatesProvider pvProv, AbsoluteDate date, Frame frame) {
421                 throw new OrekitException(new DummyLocalizable("gasp"), new RuntimeException());
422             }
423             public <T extends CalculusFieldElement<T>> FieldAttitude<T> getAttitude(FieldPVCoordinatesProvider<T> pvProv,
424                                                                                 FieldAbsoluteDate<T> date, Frame frame)
425                 {
426                 throw new OrekitException(new DummyLocalizable("gasp"), new RuntimeException());
427             }
428         };
429         KeplerianPropagator propagator = new KeplerianPropagator(orbit, wrongLaw);
430         propagator.propagate(AbsoluteDate.J2000_EPOCH.shiftedBy(10.0));
431     }
432 
433     @Test(expected = OrekitException.class)
434     public void testStepException() {
435         final KeplerianOrbit orbit =
436             new KeplerianOrbit(7.8e6, 0.032, 0.4, 0.1, 0.2, 0.3, PositionAngle.TRUE,
437                                FramesFactory.getEME2000(), AbsoluteDate.J2000_EPOCH, 3.986004415e14);
438         KeplerianPropagator propagator = new KeplerianPropagator(orbit);
439         propagator.getMultiplexer().add(new OrekitStepHandler() {
440             public void handleStep(OrekitStepInterpolator interpolator) {
441             }
442             public void finish(SpacecraftState finalState) {
443                 throw new OrekitException((Throwable) null, new DummyLocalizable("dummy error"));
444             }
445         });
446 
447         propagator.propagate(orbit.getDate().shiftedBy(-3600));
448 
449     }
450 
451     @Test(expected = OrekitException.class)
452     public void tesWrapedAttitudeException() {
453         final KeplerianOrbit orbit =
454             new KeplerianOrbit(7.8e6, 0.032, 0.4, 0.1, 0.2, 0.3, PositionAngle.TRUE,
455                                FramesFactory.getEME2000(), AbsoluteDate.J2000_EPOCH, 3.986004415e14);
456         KeplerianPropagator propagator = new KeplerianPropagator(orbit,
457                                                                  new AttitudeProvider() {
458                                                                     public Attitude getAttitude(PVCoordinatesProvider pvProv, AbsoluteDate date,
459                                                                                                 Frame frame)
460                                                                         {
461                                                                         throw new OrekitException((Throwable) null,
462                                                                                                   new DummyLocalizable("dummy error"));
463                                                                     }
464                                                                     public <T extends CalculusFieldElement<T>> FieldAttitude<T> getAttitude(FieldPVCoordinatesProvider<T> pvProv,
465                                                                                                                                         FieldAbsoluteDate<T> date, Frame frame)
466                                                                         {
467                                                                         throw new OrekitException((Throwable) null,
468                                                                                                   new DummyLocalizable("dummy error"));
469                                                                     }
470                                                                 });
471         propagator.propagate(orbit.getDate().shiftedBy(10.09));
472     }
473 
474     @Test
475     public void ascendingNode() {
476         final KeplerianOrbit orbit =
477             new KeplerianOrbit(7.8e6, 0.032, 0.4, 0.1, 0.2, 0.3, PositionAngle.TRUE,
478                                FramesFactory.getEME2000(), AbsoluteDate.J2000_EPOCH, 3.986004415e14);
479         KeplerianPropagator propagator = new KeplerianPropagator(orbit);
480         propagator.addEventDetector(new NodeDetector(orbit, FramesFactory.getITRF(IERSConventions.IERS_2010, true)));
481         AbsoluteDate farTarget = AbsoluteDate.J2000_EPOCH.shiftedBy(10000.0);
482         SpacecraftState propagated = propagator.propagate(farTarget);
483         PVCoordinates pv = propagated.getPVCoordinates(FramesFactory.getITRF(IERSConventions.IERS_2010, true));
484         Assert.assertTrue(farTarget.durationFrom(propagated.getDate()) > 3500.0);
485         Assert.assertTrue(farTarget.durationFrom(propagated.getDate()) < 4000.0);
486         Assert.assertEquals(0, pv.getPosition().getZ(), 2.0e-6);
487         Assert.assertTrue(pv.getVelocity().getZ() > 0);
488     }
489 
490     @Test
491     public void stopAtTargetDate() {
492         final KeplerianOrbit orbit =
493             new KeplerianOrbit(7.8e6, 0.032, 0.4, 0.1, 0.2, 0.3, PositionAngle.TRUE,
494                                FramesFactory.getEME2000(), AbsoluteDate.J2000_EPOCH, 3.986004415e14);
495         KeplerianPropagator propagator = new KeplerianPropagator(orbit);
496         Frame itrf = FramesFactory.getITRF(IERSConventions.IERS_2010, true);
497         propagator.addEventDetector(new NodeDetector(orbit, itrf).withHandler(new ContinueOnEvent<NodeDetector>()));
498         AbsoluteDate farTarget = orbit.getDate().shiftedBy(10000.0);
499         SpacecraftState propagated = propagator.propagate(farTarget);
500         Assert.assertEquals(0.0, FastMath.abs(farTarget.durationFrom(propagated.getDate())), 1.0e-3);
501     }
502 
503     @Test
504     public void perigee() {
505         final KeplerianOrbit orbit =
506             new KeplerianOrbit(7.8e6, 0.032, 0.4, 0.1, 0.2, 0.3, PositionAngle.TRUE,
507                                FramesFactory.getEME2000(), AbsoluteDate.J2000_EPOCH, 3.986004415e14);
508         KeplerianPropagator propagator = new KeplerianPropagator(orbit);
509         propagator.addEventDetector(new ApsideDetector(orbit));
510         AbsoluteDate farTarget = AbsoluteDate.J2000_EPOCH.shiftedBy(10000.0);
511         SpacecraftState propagated = propagator.propagate(farTarget);
512         PVCoordinates pv = propagated.getPVCoordinates(FramesFactory.getITRF(IERSConventions.IERS_2010, true));
513         Assert.assertTrue(farTarget.durationFrom(propagated.getDate()) > 3000.0);
514         Assert.assertTrue(farTarget.durationFrom(propagated.getDate()) < 3500.0);
515         Assert.assertEquals(orbit.getA() * (1.0 - orbit.getE()), pv.getPosition().getNorm(), 1.0e-6);
516     }
517 
518     @Test
519     public void altitude() {
520         final KeplerianOrbit orbit =
521             new KeplerianOrbit(7.8e6, 0.032, 0.4, 0.1, 0.2, 0.3, PositionAngle.TRUE,
522                                FramesFactory.getEME2000(), AbsoluteDate.J2000_EPOCH, 3.986004415e14);
523         KeplerianPropagator propagator = new KeplerianPropagator(orbit);
524         BodyShape bodyShape =
525             new OneAxisEllipsoid(6378137.0, 1.0 / 298.257222101, FramesFactory.getITRF(IERSConventions.IERS_2010, true));
526         AltitudeDetector detector =
527             new AltitudeDetector(0.05 * orbit.getKeplerianPeriod(),
528                                  1500000, bodyShape);
529         Assert.assertEquals(1500000, detector.getAltitude(), 1.0e-12);
530         propagator.addEventDetector(detector);
531         AbsoluteDate farTarget = AbsoluteDate.J2000_EPOCH.shiftedBy(10000.0);
532         SpacecraftState propagated = propagator.propagate(farTarget);
533         Assert.assertTrue(farTarget.durationFrom(propagated.getDate()) > 5400.0);
534         Assert.assertTrue(farTarget.durationFrom(propagated.getDate()) < 5500.0);
535         GeodeticPoint gp = bodyShape.transform(propagated.getPVCoordinates().getPosition(),
536                                                propagated.getFrame(), propagated.getDate());
537         Assert.assertEquals(1500000, gp.getAltitude(), 0.1);
538     }
539 
540     @Test
541     public void date() {
542         final KeplerianOrbit orbit =
543             new KeplerianOrbit(7.8e6, 0.032, 0.4, 0.1, 0.2, 0.3, PositionAngle.TRUE,
544                                FramesFactory.getEME2000(), AbsoluteDate.J2000_EPOCH, 3.986004415e14);
545         KeplerianPropagator propagator = new KeplerianPropagator(orbit);
546         final AbsoluteDate stopDate = AbsoluteDate.J2000_EPOCH.shiftedBy(500.0);
547         propagator.addEventDetector(new DateDetector(stopDate));
548         AbsoluteDate farTarget = AbsoluteDate.J2000_EPOCH.shiftedBy(10000.0);
549         SpacecraftState propagated = propagator.propagate(farTarget);
550         Assert.assertEquals(0, stopDate.durationFrom(propagated.getDate()), 1.0e-10);
551     }
552 
553     @Test
554     public void setting() {
555         final KeplerianOrbit orbit =
556             new KeplerianOrbit(7.8e6, 0.032, 0.4, 0.1, 0.2, 0.3, PositionAngle.TRUE,
557                                FramesFactory.getEME2000(), AbsoluteDate.J2000_EPOCH, 3.986004415e14);
558         KeplerianPropagator propagator = new KeplerianPropagator(orbit);
559         final OneAxisEllipsoid earthShape =
560             new OneAxisEllipsoid(6378136.460, 1 / 298.257222101, FramesFactory.getITRF(IERSConventions.IERS_2010, true));
561         final TopocentricFrame topo =
562             new TopocentricFrame(earthShape, new GeodeticPoint(0.389, -2.962, 0), null);
563         propagator.addEventDetector(new ElevationDetector(60, AbstractDetector.DEFAULT_THRESHOLD, topo).withConstantElevation(0.09));
564         AbsoluteDate farTarget = AbsoluteDate.J2000_EPOCH.shiftedBy(10000.0);
565         SpacecraftState propagated = propagator.propagate(farTarget);
566         final double elevation = topo.getElevation(propagated.getPVCoordinates().getPosition(),
567                                                    propagated.getFrame(),
568                                                    propagated.getDate());
569         final double zVelocity = propagated.getPVCoordinates(topo).getVelocity().getZ();
570         Assert.assertTrue(farTarget.durationFrom(propagated.getDate()) > 7800.0);
571         Assert.assertTrue(farTarget.durationFrom(propagated.getDate()) < 7900.0);
572         Assert.assertEquals(0.09, elevation, 1.0e-9);
573         Assert.assertTrue(zVelocity < 0);
574     }
575 
576     @Test
577     public void fixedStep() {
578         final KeplerianOrbit orbit =
579             new KeplerianOrbit(7.8e6, 0.032, 0.4, 0.1, 0.2, 0.3, PositionAngle.TRUE,
580                                FramesFactory.getEME2000(), AbsoluteDate.J2000_EPOCH, 3.986004415e14);
581         KeplerianPropagator propagator = new KeplerianPropagator(orbit);
582         final double step = 100.0;
583         final int[] counter = new int[] {0};  // mutable int
584         propagator.setStepHandler(step, new OrekitFixedStepHandler() {
585             private AbsoluteDate previous;
586             public void handleStep(SpacecraftState currentState) {
587                 if (previous != null) {
588                     Assert.assertEquals(step, currentState.getDate().durationFrom(previous), 1.0e-10);
589                 }
590                 // check state is accurate
591                 PVCoordinates expected = new KeplerianPropagator(orbit)
592                         .propagate(currentState.getDate()).getPVCoordinates();
593                 MatcherAssert.assertThat(
594                         currentState.getPVCoordinates(),
595                         OrekitMatchers.pvIs(expected));
596                 previous = currentState.getDate();
597                 counter[0]++;
598             }
599         });
600         AbsoluteDate farTarget = AbsoluteDate.J2000_EPOCH.shiftedBy(10000.0);
601         propagator.propagate(farTarget);
602         // check the step handler was executed
603         Assert.assertEquals(
604                 counter[0],
605                 (int) (farTarget.durationFrom(orbit.getDate()) / step) + 1);
606     }
607 
608     @Test
609     public void variableStep() {
610         final KeplerianOrbit orbit =
611             new KeplerianOrbit(7.8e6, 0.032, 0.4, 0.1, 0.2, 0.3, PositionAngle.TRUE,
612                                FramesFactory.getEME2000(), AbsoluteDate.J2000_EPOCH, 3.986004415e14);
613         KeplerianPropagator propagator = new KeplerianPropagator(orbit);
614         final double step = orbit.getKeplerianPeriod() / 100;
615         final int[] counter = new int[] {0};  // mutable int
616         propagator.setStepHandler(new OrekitStepHandler() {
617             private AbsoluteDate t = orbit.getDate();
618             @Override
619             public void handleStep(OrekitStepInterpolator interpolator) {
620                 // check the states provided by the interpolator are accurate.
621                 do {
622                     PVCoordinates expected = new KeplerianPropagator(orbit)
623                             .propagate(t).getPVCoordinates();
624                     MatcherAssert.assertThat(
625                             interpolator.getInterpolatedState(t).getPVCoordinates(),
626                             OrekitMatchers.pvIs(expected));
627                     t = t.shiftedBy(step);
628                     counter[0]++;
629                 } while (t.compareTo(interpolator.getCurrentState().getDate()) <= 0);
630             }
631         });
632         AbsoluteDate farTarget = AbsoluteDate.J2000_EPOCH.shiftedBy(10000.0);
633         propagator.propagate(farTarget);
634         // check the step handler was executed
635         Assert.assertEquals(
636                 counter[0],
637                 (int) (farTarget.durationFrom(orbit.getDate()) / step) + 1);
638     }
639 
640     @Test
641     public void ephemeris() {
642         final KeplerianOrbit orbit =
643             new KeplerianOrbit(7.8e6, 0.032, 0.4, 0.1, 0.2, 0.3, PositionAngle.TRUE,
644                                FramesFactory.getEME2000(), AbsoluteDate.J2000_EPOCH, 3.986004415e14);
645         KeplerianPropagator propagator = new KeplerianPropagator(orbit);
646         AbsoluteDate farTarget = AbsoluteDate.J2000_EPOCH.shiftedBy(10000.0);
647         final EphemerisGenerator generator = propagator.getEphemerisGenerator();
648         propagator.propagate(farTarget);
649         BoundedPropagator ephemeris = generator.getGeneratedEphemeris();
650         Assert.assertEquals(0.0, ephemeris.getMinDate().durationFrom(orbit.getDate()), 1.0e10);
651         Assert.assertEquals(0.0, ephemeris.getMaxDate().durationFrom(farTarget), 1.0e10);
652     }
653 
654     @Test
655     public void testIssue14() {
656         AbsoluteDate initialDate = AbsoluteDate.J2000_EPOCH;
657         final KeplerianOrbit initialOrbit =
658             new KeplerianOrbit(7.8e6, 0.032, 0.4, 0.1, 0.2, 0.3, PositionAngle.TRUE,
659                                FramesFactory.getEME2000(), initialDate, 3.986004415e14);
660         KeplerianPropagator propagator = new KeplerianPropagator(initialOrbit);
661 
662         propagator.getEphemerisGenerator();
663         propagator.propagate(initialDate.shiftedBy(initialOrbit.getKeplerianPeriod()));
664         PVCoordinates pv1 = propagator.getPVCoordinates(initialDate, FramesFactory.getEME2000());
665 
666         final EphemerisGenerator generator = propagator.getEphemerisGenerator();
667         propagator.propagate(initialDate.shiftedBy(initialOrbit.getKeplerianPeriod()));
668         PVCoordinates pv2 = generator.getGeneratedEphemeris().getPVCoordinates(initialDate, FramesFactory.getEME2000());
669 
670         Assert.assertEquals(0.0, pv1.getPosition().subtract(pv2.getPosition()).getNorm(), 1.0e-15);
671         Assert.assertEquals(0.0, pv1.getVelocity().subtract(pv2.getVelocity()).getNorm(), 1.0e-15);
672 
673     }
674 
675     @Test
676     public void testIssue107() {
677         final TimeScale utc = TimeScalesFactory.getUTC();
678         final Vector3D position = new Vector3D(-6142438.668, 3492467.56, -25767.257);
679         final Vector3D velocity = new Vector3D(505.848, 942.781, 7435.922);
680         final AbsoluteDate date = new AbsoluteDate(2003, 9, 16, utc);
681         final Orbit orbit = new CircularOrbit(new PVCoordinates(position,  velocity),
682                                               FramesFactory.getEME2000(), date, mu);
683 
684         Propagator propagator = new KeplerianPropagator(orbit) {
685             AbsoluteDate lastDate = AbsoluteDate.PAST_INFINITY;
686 
687             protected SpacecraftState basicPropagate(final AbsoluteDate date) {
688                 if (date.compareTo(lastDate) < 0) {
689                     throw new OrekitException(LocalizedCoreFormats.SIMPLE_MESSAGE,
690                                                    "no backward propagation allowed");
691                 }
692                 lastDate = date;
693                 return super.basicPropagate(date);
694             }
695         };
696 
697         SpacecraftState finalState = propagator.propagate(date.shiftedBy(3600.0));
698         Assert.assertEquals(3600.0, finalState.getDate().durationFrom(date), 1.0e-15);
699 
700     }
701 
702     @Test
703     public void testMu() {
704         final KeplerianOrbit orbit1 =
705                 new KeplerianOrbit(7.8e6, 0.032, 0.4, 0.1, 0.2, 0.3, PositionAngle.TRUE,
706                                    FramesFactory.getEME2000(), AbsoluteDate.J2000_EPOCH,
707                                    Constants.WGS84_EARTH_MU);
708         final KeplerianOrbit orbit2 =
709                 new KeplerianOrbit(7.8e6, 0.032, 0.4, 0.1, 0.2, 0.3, PositionAngle.TRUE,
710                                    FramesFactory.getEME2000(), AbsoluteDate.J2000_EPOCH,
711                                    Constants.EIGEN5C_EARTH_MU);
712         final AbsoluteDate target = orbit1.getDate().shiftedBy(10000.0);
713         PVCoordinates pv1       = new KeplerianPropagator(orbit1).propagate(target).getPVCoordinates();
714         PVCoordinates pv2       = new KeplerianPropagator(orbit2).propagate(target).getPVCoordinates();
715         PVCoordinates pvWithMu1 = new KeplerianPropagator(orbit2, orbit1.getMu()).propagate(target).getPVCoordinates();
716         Assert.assertEquals(0.026054, Vector3D.distance(pv1.getPosition(), pv2.getPosition()),       1.0e-6);
717         Assert.assertEquals(0.0,      Vector3D.distance(pv1.getPosition(), pvWithMu1.getPosition()), 1.0e-15);
718     }
719 
720     @Test
721     public void testResetStateForward() {
722         final Frame eme2000 = FramesFactory.getEME2000();
723         final AbsoluteDate date = new AbsoluteDate(new DateComponents(2008, 6, 23),
724                                                    new TimeComponents(14, 0, 0),
725                                                    TimeScalesFactory.getUTC());
726         final Orbit orbit = new KeplerianOrbit(8000000.0, 0.01, 0.87, 2.44, 0.21, -1.05, PositionAngle.MEAN,
727                                            eme2000,
728                                            date, Constants.EIGEN5C_EARTH_MU);
729         final KeplerianPropagator propagator = new KeplerianPropagator(orbit, new LofOffset(eme2000, LOFType.LVLH));
730 
731         // maneuver along Z in attitude aligned with LVLH will change orbital plane
732         final AbsoluteDate maneuverDate = date.shiftedBy(1000.0);
733         propagator.addEventDetector(new ImpulseManeuver<>(new DateDetector(maneuverDate),
734                                                           new Vector3D(0.0, 0.0, -100.0),
735                                                           350.0));
736 
737         final Vector3D initialNormal = orbit.getPVCoordinates().getMomentum();
738         propagator.setStepHandler(60.0, state -> {
739             final Vector3D currentNormal = state.getPVCoordinates().getMomentum();
740             if (state.getDate().isBefore(maneuverDate)) {
741                 Assert.assertEquals(0.000, Vector3D.angle(initialNormal, currentNormal), 1.0e-3);
742             } else {
743                 Assert.assertEquals(0.014, Vector3D.angle(initialNormal, currentNormal), 1.0e-3);
744             }
745         });
746 
747         propagator.propagate(orbit.getDate().shiftedBy(1500.0));
748 
749     }
750 
751     @Test
752     public void testResetStateBackward() {
753         final Frame eme2000 = FramesFactory.getEME2000();
754         final AbsoluteDate date = new AbsoluteDate(new DateComponents(2008, 6, 23),
755                                                    new TimeComponents(14, 0, 0),
756                                                    TimeScalesFactory.getUTC());
757         final Orbit orbit = new KeplerianOrbit(8000000.0, 0.01, 0.87, 2.44, 0.21, -1.05, PositionAngle.MEAN,
758                                            eme2000,
759                                            date, Constants.EIGEN5C_EARTH_MU);
760         final KeplerianPropagator propagator = new KeplerianPropagator(orbit, new LofOffset(eme2000, LOFType.LVLH));
761 
762         // maneuver along Z in attitude aligned with LVLH will change orbital plane
763         final AbsoluteDate maneuverDate = date.shiftedBy(-1000.0);
764         propagator.addEventDetector(new ImpulseManeuver<>(new DateDetector(maneuverDate),
765                                                           new Vector3D(0.0, 0.0, -100.0),
766                                                           350.0));
767 
768         final Vector3D initialNormal = orbit.getPVCoordinates().getMomentum();
769         propagator.setStepHandler(60.0, state -> {
770             final Vector3D currentNormal = state.getPVCoordinates().getMomentum();
771             if (state.getDate().isAfter(maneuverDate)) {
772                 Assert.assertEquals(0.000, Vector3D.angle(initialNormal, currentNormal), 1.0e-3);
773             } else {
774                 Assert.assertEquals(0.014, Vector3D.angle(initialNormal, currentNormal), 1.0e-3);
775             }
776         });
777 
778         propagator.propagate(orbit.getDate().shiftedBy(-1500.0));
779 
780     }
781 
782     @Test
783     public void testNoDerivatives() {
784         for (OrbitType type : OrbitType.values()) {
785 
786             // create an initial orbit with non-Keplerian acceleration
787             final AbsoluteDate date         = new AbsoluteDate(2003, 9, 16, TimeScalesFactory.getUTC());
788             final Vector3D     position     = new Vector3D(-6142438.668, 3492467.56, -25767.257);
789             final Vector3D     velocity     = new Vector3D(505.848, 942.781, 7435.922);
790             final Vector3D     keplerAcceleration = new Vector3D(-mu / position.getNormSq(), position.normalize());
791             final Vector3D     nonKeplerAcceleration = new Vector3D(0.001, 0.002, 0.003);
792             final Vector3D     acceleration = keplerAcceleration.add(nonKeplerAcceleration);
793             final TimeStampedPVCoordinates pva = new TimeStampedPVCoordinates(date, position, velocity, acceleration);
794             final Orbit initial = type.convertType(new CartesianOrbit(pva, FramesFactory.getEME2000(), mu));
795             Assert.assertEquals(type, initial.getType());
796 
797             // the derivatives are available at this stage
798             checkDerivatives(initial, true);
799 
800             KeplerianPropagator propagator = new KeplerianPropagator(initial);
801             Assert.assertEquals(type, propagator.getInitialState().getOrbit().getType());
802 
803             // non-Keplerian derivatives are explicitly removed when building the Keplerian-only propagator
804             checkDerivatives(propagator.getInitialState().getOrbit(), false);
805 
806             PVCoordinates initPV = propagator.getInitialState().getOrbit().getPVCoordinates();
807             Assert.assertEquals(nonKeplerAcceleration.getNorm(), Vector3D.distance(acceleration, initPV.getAcceleration()), 2.0e-15);
808             Assert.assertEquals(0.0,
809                                 Vector3D.distance(keplerAcceleration, initPV.getAcceleration()),
810                                 4.0e-15);
811 
812             double dt = 0.2 * initial.getKeplerianPeriod();
813             Orbit orbit = propagator.propagateOrbit(initial.getDate().shiftedBy(dt));
814             Assert.assertEquals(type, orbit.getType());
815 
816             // at the end, we don't have non-Keplerian derivatives
817             checkDerivatives(orbit, false);
818 
819             // using shiftedBy on the initial orbit, non-Keplerian derivatives would have been preserved
820             checkDerivatives(initial.shiftedBy(dt), true);
821 
822         }
823     }
824 
825     @Test
826     public void testStackableGenerators() {
827         final Frame eme2000 = FramesFactory.getEME2000();
828         final AbsoluteDate date = new AbsoluteDate(new DateComponents(2008, 6, 23),
829                                                    new TimeComponents(14, 0, 0),
830                                                    TimeScalesFactory.getUTC());
831         final Orbit orbit = new KeplerianOrbit(8000000.0, 0.01, 0.87, 2.44, 0.21, -1.05, PositionAngle.MEAN,
832                                            eme2000,
833                                            date, Constants.EIGEN5C_EARTH_MU);
834         final KeplerianPropagator propagator = new KeplerianPropagator(orbit, new LofOffset(eme2000, LOFType.LVLH));
835 
836         // we have A → B → C → D → E → F but register them in a different order
837         propagator.addAdditionalStateProvider(new DependentGenerator("F", "E"));
838         propagator.addAdditionalStateProvider(new DependentGenerator("B", "A"));
839         propagator.addAdditionalStateProvider(new DependentGenerator("E", "D"));
840         propagator.addAdditionalStateProvider(new DependentGenerator("C", "B"));
841         propagator.addAdditionalStateProvider(new DependentGenerator("A", null));
842         propagator.addAdditionalStateProvider(new DependentGenerator("D", "C"));
843 
844         final SpacecraftState finalState = propagator.propagate(orbit.getDate().shiftedBy(3600.0));
845         Assert.assertEquals(1,   finalState.getAdditionalState("A").length);
846         Assert.assertEquals(1.0, finalState.getAdditionalState("A")[0], 1.0e-15);
847         Assert.assertEquals(1,   finalState.getAdditionalState("B").length);
848         Assert.assertEquals(2.0, finalState.getAdditionalState("B")[0], 1.0e-15);
849         Assert.assertEquals(1,   finalState.getAdditionalState("C").length);
850         Assert.assertEquals(3.0, finalState.getAdditionalState("C")[0], 1.0e-15);
851         Assert.assertEquals(1,   finalState.getAdditionalState("D").length);
852         Assert.assertEquals(4.0, finalState.getAdditionalState("D")[0], 1.0e-15);
853         Assert.assertEquals(1,   finalState.getAdditionalState("E").length);
854         Assert.assertEquals(5.0, finalState.getAdditionalState("E")[0], 1.0e-15);
855         Assert.assertEquals(1,   finalState.getAdditionalState("F").length);
856         Assert.assertEquals(6.0, finalState.getAdditionalState("F")[0], 1.0e-15);
857 
858     }
859 
860     @Test
861     public void testCircularDependency() {
862         final Frame eme2000 = FramesFactory.getEME2000();
863         final AbsoluteDate date = new AbsoluteDate(new DateComponents(2008, 6, 23),
864                                                    new TimeComponents(14, 0, 0),
865                                                    TimeScalesFactory.getUTC());
866         final Orbit orbit = new KeplerianOrbit(8000000.0, 0.01, 0.87, 2.44, 0.21, -1.05, PositionAngle.MEAN,
867                                            eme2000,
868                                            date, Constants.EIGEN5C_EARTH_MU);
869         final KeplerianPropagator propagator = new KeplerianPropagator(orbit, new LofOffset(eme2000, LOFType.LVLH));
870 
871         // here, the dependency creates a loop, which is detected and its adders ignored
872         propagator.addAdditionalStateProvider(new DependentGenerator("F", "E"));
873         propagator.addAdditionalStateProvider(new DependentGenerator("B", "A"));
874         propagator.addAdditionalStateProvider(new DependentGenerator("E", "D"));
875         propagator.addAdditionalStateProvider(new DependentGenerator("C", "B"));
876         propagator.addAdditionalStateProvider(new DependentGenerator("A", null));
877         propagator.addAdditionalStateProvider(new DependentGenerator("D", "F"));
878 
879         final SpacecraftState finalState = propagator.propagate(orbit.getDate().shiftedBy(3600.0));
880         Assert.assertEquals(1,   finalState.getAdditionalState("A").length);
881         Assert.assertEquals(1.0, finalState.getAdditionalState("A")[0], 1.0e-15);
882         Assert.assertEquals(1,   finalState.getAdditionalState("B").length);
883         Assert.assertEquals(2.0, finalState.getAdditionalState("B")[0], 1.0e-15);
884         Assert.assertEquals(1,   finalState.getAdditionalState("C").length);
885         Assert.assertEquals(3.0, finalState.getAdditionalState("C")[0], 1.0e-15);
886         Assert.assertFalse(finalState.hasAdditionalState("D"));
887         Assert.assertFalse(finalState.hasAdditionalState("E"));
888         Assert.assertFalse(finalState.hasAdditionalState("F"));
889 
890     }
891 
892     private void checkDerivatives(final Orbit orbit, final boolean expectedDerivatives) {
893         Assert.assertEquals(expectedDerivatives, orbit.hasDerivatives());
894         Assert.assertNotEquals(expectedDerivatives, Double.isNaN(orbit.getADot()));
895         Assert.assertNotEquals(expectedDerivatives, Double.isNaN(orbit.getEquinoctialExDot()));
896         Assert.assertNotEquals(expectedDerivatives, Double.isNaN(orbit.getEquinoctialEyDot()));
897         Assert.assertNotEquals(expectedDerivatives, Double.isNaN(orbit.getHxDot()));
898         Assert.assertNotEquals(expectedDerivatives, Double.isNaN(orbit.getHyDot()));
899         Assert.assertNotEquals(expectedDerivatives, Double.isNaN(orbit.getLEDot()));
900         Assert.assertNotEquals(expectedDerivatives, Double.isNaN(orbit.getLvDot()));
901         Assert.assertNotEquals(expectedDerivatives, Double.isNaN(orbit.getLMDot()));
902         Assert.assertNotEquals(expectedDerivatives, Double.isNaN(orbit.getEDot()));
903         Assert.assertNotEquals(expectedDerivatives, Double.isNaN(orbit.getIDot()));
904     }
905 
906     private static double tangLEmLv(double Lv, double ex, double ey){
907         // tan ((LE - Lv) /2)) =
908         return (ey*FastMath.cos(Lv) - ex*FastMath.sin(Lv)) /
909         (1 + ex*FastMath.cos(Lv) + ey*FastMath.sin(Lv) + FastMath.sqrt(1 - ex*ex - ey*ey));
910     }
911 
912     @Before
913     public void setUp() {
914         Utils.setDataRoot("regular-data");
915         mu  = 3.9860047e14;
916     }
917 
918     @After
919     public void tearDown() {
920         mu   = Double.NaN;
921     }
922 
923     private static class DependentGenerator implements AdditionalStateProvider {
924 
925         private final String name;
926         private final String dependency;
927 
928         DependentGenerator(final String name, final String dependency) {
929             this.name       = name;
930             this.dependency = dependency;
931         }
932 
933         public String getName() {
934             return name;
935         }
936 
937         public boolean yield(final SpacecraftState state) {
938             return dependency != null && state.getAdditionalStatesValues().getEntry(dependency) == null;
939         }
940 
941         public double[] getAdditionalState(final SpacecraftState state) {
942             return new double[] {
943                 dependency == null ? 1.0 : state.getAdditionalState(dependency)[0] + 1.0
944             };
945         }
946 
947     }
948 
949 }
950