1
2
3
4
5
6
7
8
9
10
11
12
13
14
15
16
17 package org.orekit.propagation.events;
18
19 import java.util.Locale;
20 import java.util.function.Function;
21
22 import org.hamcrest.CoreMatchers;
23 import org.hamcrest.MatcherAssert;
24 import org.hipparchus.exception.LocalizedCoreFormats;
25 import org.hipparchus.geometry.euclidean.threed.Vector3D;
26 import org.hipparchus.ode.ODEIntegrator;
27 import org.hipparchus.ode.events.Action;
28 import org.hipparchus.ode.nonstiff.ClassicalRungeKuttaIntegrator;
29 import org.hipparchus.ode.nonstiff.DormandPrince853Integrator;
30 import org.hipparchus.util.FastMath;
31 import org.junit.Assert;
32 import org.junit.Before;
33 import org.junit.Test;
34 import org.orekit.Utils;
35 import org.orekit.errors.OrekitException;
36 import org.orekit.frames.Frame;
37 import org.orekit.frames.FramesFactory;
38 import org.orekit.orbits.CircularOrbit;
39 import org.orekit.orbits.EquinoctialOrbit;
40 import org.orekit.orbits.KeplerianOrbit;
41 import org.orekit.orbits.Orbit;
42 import org.orekit.orbits.OrbitType;
43 import org.orekit.orbits.PositionAngle;
44 import org.orekit.propagation.Propagator;
45 import org.orekit.propagation.SpacecraftState;
46 import org.orekit.propagation.analytical.KeplerianPropagator;
47 import org.orekit.propagation.events.handlers.ContinueOnEvent;
48 import org.orekit.propagation.events.handlers.EventHandler;
49 import org.orekit.propagation.events.handlers.StopOnEvent;
50 import org.orekit.propagation.numerical.NumericalPropagator;
51 import org.orekit.propagation.sampling.OrekitFixedStepHandler;
52 import org.orekit.time.AbsoluteDate;
53 import org.orekit.time.TimeScale;
54 import org.orekit.time.TimeScalesFactory;
55 import org.orekit.utils.Constants;
56 import org.orekit.utils.PVCoordinates;
57 import org.orekit.utils.PVCoordinatesProvider;
58
59 public class EventDetectorTest {
60
61 private double mu;
62
63 @Test
64 public void testEventHandlerInit() {
65 final TimeScale utc = TimeScalesFactory.getUTC();
66 final Vector3D position = new Vector3D(-6142438.668, 3492467.56, -25767.257);
67 final Vector3D velocity = new Vector3D(505.848, 942.781, 7435.922);
68 final AbsoluteDate date = new AbsoluteDate(2003, 9, 16, utc);
69 final Orbit orbit = new CircularOrbit(new PVCoordinates(position, velocity),
70 FramesFactory.getEME2000(), date, mu);
71
72 final boolean[] eventOccurred = new boolean[1];
73 EventHandler<DateDetector> handler = new EventHandler<DateDetector>() {
74 private boolean initCalled;
75 @Override
76 public Action eventOccurred(SpacecraftState s,
77 DateDetector detector,
78 boolean increasing) {
79 if (!initCalled) {
80 throw new RuntimeException("init() not called before eventOccurred()");
81 }
82 eventOccurred[0] = true;
83 return Action.STOP;
84 }
85
86 @Override
87 public void init(SpacecraftState initialState,
88 AbsoluteDate target,
89 DateDetector detector) {
90 initCalled = true;
91 }
92 };
93
94 Propagator propagator = new KeplerianPropagator(orbit);
95 double stepSize = 60.0;
96 propagator.addEventDetector(new DateDetector(date.shiftedBy(5.25 * stepSize)).withHandler(handler));
97 propagator.propagate(date.shiftedBy(10 * stepSize));
98 Assert.assertTrue(eventOccurred[0]);
99
100 }
101
102 @Test
103 public void testBasicScheduling() {
104
105 final TimeScale utc = TimeScalesFactory.getUTC();
106 final Vector3D position = new Vector3D(-6142438.668, 3492467.56, -25767.257);
107 final Vector3D velocity = new Vector3D(505.848, 942.781, 7435.922);
108 final AbsoluteDate date = new AbsoluteDate(2003, 9, 16, utc);
109 final Orbit orbit = new CircularOrbit(new PVCoordinates(position, velocity),
110 FramesFactory.getEME2000(), date, mu);
111
112 Propagator propagator = new KeplerianPropagator(orbit);
113 double stepSize = 60.0;
114 OutOfOrderChecker checker = new OutOfOrderChecker(stepSize);
115 propagator.addEventDetector(new DateDetector(date.shiftedBy(5.25 * stepSize)).withHandler(checker));
116 propagator.setStepHandler(stepSize, checker);
117 propagator.propagate(date.shiftedBy(10 * stepSize));
118 Assert.assertTrue(checker.outOfOrderCallDetected());
119
120 }
121
122 private static class OutOfOrderChecker implements EventHandler<DateDetector>, OrekitFixedStepHandler {
123
124 private AbsoluteDate triggerDate;
125 private boolean outOfOrderCallDetected;
126 private double stepSize;
127
128 public OutOfOrderChecker(final double stepSize) {
129 triggerDate = null;
130 outOfOrderCallDetected = false;
131 this.stepSize = stepSize;
132 }
133
134 public Action eventOccurred(SpacecraftState s, DateDetector detector, boolean increasing) {
135 triggerDate = s.getDate();
136 return Action.CONTINUE;
137 }
138
139 public void handleStep(SpacecraftState currentState) {
140
141
142
143 if (triggerDate != null) {
144 double dt = currentState.getDate().durationFrom(triggerDate);
145 if (dt < 0) {
146 outOfOrderCallDetected = true;
147 Assert.assertTrue(FastMath.abs(dt) < (2 * stepSize));
148 }
149 }
150 }
151
152 public boolean outOfOrderCallDetected() {
153 return outOfOrderCallDetected;
154 }
155
156 @Override
157 public void init(SpacecraftState initialState, AbsoluteDate target, double step) {
158 }
159
160 }
161
162 @Test
163 public void testIssue108Numerical() {
164 final TimeScale utc = TimeScalesFactory.getUTC();
165 final Vector3D position = new Vector3D(-6142438.668, 3492467.56, -25767.257);
166 final Vector3D velocity = new Vector3D(505.848, 942.781, 7435.922);
167 final AbsoluteDate date = new AbsoluteDate(2003, 9, 16, utc);
168 final Orbit orbit = new CircularOrbit(new PVCoordinates(position, velocity),
169 FramesFactory.getEME2000(), date, mu);
170 final double step = 60.0;
171 final int n = 100;
172 NumericalPropagator propagator = new NumericalPropagator(new ClassicalRungeKuttaIntegrator(step));
173 propagator.resetInitialState(new SpacecraftState(orbit));
174 GCallsCounter counter = new GCallsCounter(100000.0, 1.0e-6, 20, new StopOnEvent<GCallsCounter>());
175 propagator.addEventDetector(counter);
176 propagator.propagate(date.shiftedBy(n * step));
177 Assert.assertEquals(n + 1, counter.getCount());
178 }
179
180 @Test
181 public void testIssue108Analytical() {
182 final TimeScale utc = TimeScalesFactory.getUTC();
183 final Vector3D position = new Vector3D(-6142438.668, 3492467.56, -25767.257);
184 final Vector3D velocity = new Vector3D(505.848, 942.781, 7435.922);
185 final AbsoluteDate date = new AbsoluteDate(2003, 9, 16, utc);
186 final Orbit orbit = new CircularOrbit(new PVCoordinates(position, velocity),
187 FramesFactory.getEME2000(), date, mu);
188 final double step = 60.0;
189 final int n = 100;
190 KeplerianPropagator propagator = new KeplerianPropagator(orbit);
191 GCallsCounter counter = new GCallsCounter(100000.0, 1.0e-6, 20, new StopOnEvent<GCallsCounter>());
192 propagator.addEventDetector(counter);
193 propagator.setStepHandler(step, currentState -> {});
194 propagator.propagate(date.shiftedBy(n * step));
195
196 Assert.assertEquals(2, counter.getCount());
197 }
198
199 private static class GCallsCounter extends AbstractDetector<GCallsCounter> {
200
201 private int count;
202
203 public GCallsCounter(final double maxCheck, final double threshold,
204 final int maxIter, final EventHandler<? super GCallsCounter> handler) {
205 super(maxCheck, threshold, maxIter, handler);
206 count = 0;
207 }
208
209 protected GCallsCounter create(final double newMaxCheck, final double newThreshold,
210 final int newMaxIter, final EventHandler<? super GCallsCounter> newHandler) {
211 return new GCallsCounter(newMaxCheck, newThreshold, newMaxIter, newHandler);
212 }
213
214 public int getCount() {
215 return count;
216 }
217
218 public double g(SpacecraftState s) {
219 count++;
220 return 1.0;
221 }
222
223 }
224
225 @Test
226 public void testNoisyGFunction() {
227
228
229 Frame eme2000 = FramesFactory.getEME2000();
230 TimeScale utc = TimeScalesFactory.getUTC();
231 AbsoluteDate initialDate = new AbsoluteDate(2011, 5, 11, utc);
232 AbsoluteDate startDate = new AbsoluteDate(2032, 10, 17, utc);
233 AbsoluteDate interruptDate = new AbsoluteDate(2032, 10, 18, utc);
234 AbsoluteDate targetDate = new AbsoluteDate(2211, 5, 11, utc);
235 KeplerianPropagator k1 =
236 new KeplerianPropagator(new EquinoctialOrbit(new PVCoordinates(new Vector3D(4008462.4706055815, -3155502.5373837613, -5044275.9880020910),
237 new Vector3D(-5012.9298276860990, 1920.3567095973078, -5172.7403501801580)),
238 eme2000, initialDate, Constants.WGS84_EARTH_MU));
239 KeplerianPropagator k2 =
240 new KeplerianPropagator(new EquinoctialOrbit(new PVCoordinates(new Vector3D(4008912.4039522274, -3155453.3125615157, -5044297.6484738905),
241 new Vector3D(-5012.5883854112530, 1920.6332221785074, -5172.2177085540500)),
242 eme2000, initialDate, Constants.WGS84_EARTH_MU));
243 k2.addEventDetector(new CloseApproachDetector(2015.243454166727, 0.0001, 100,
244 new ContinueOnEvent<CloseApproachDetector>(),
245 k1));
246 k2.addEventDetector(new DateDetector(Constants.JULIAN_DAY, 1.0e-6, interruptDate));
247 SpacecraftState s = k2.propagate(startDate, targetDate);
248 Assert.assertEquals(0.0, interruptDate.durationFrom(s.getDate()), 1.1e-6);
249 }
250
251 private static class CloseApproachDetector extends AbstractDetector<CloseApproachDetector> {
252
253 private final PVCoordinatesProvider provider;
254
255 public CloseApproachDetector(double maxCheck, double threshold,
256 final int maxIter, final EventHandler<? super CloseApproachDetector> handler,
257 PVCoordinatesProvider provider) {
258 super(maxCheck, threshold, maxIter, handler);
259 this.provider = provider;
260 }
261
262 public double g(final SpacecraftState s) {
263 PVCoordinates pv1 = provider.getPVCoordinates(s.getDate(), s.getFrame());
264 PVCoordinates pv2 = s.getPVCoordinates();
265 Vector3D deltaP = pv1.getPosition().subtract(pv2.getPosition());
266 Vector3D deltaV = pv1.getVelocity().subtract(pv2.getVelocity());
267 double radialVelocity = Vector3D.dotProduct(deltaP.normalize(), deltaV);
268 return radialVelocity;
269 }
270
271 protected CloseApproachDetector create(final double newMaxCheck, final double newThreshold,
272 final int newMaxIter,
273 final EventHandler<? super CloseApproachDetector> newHandler) {
274 return new CloseApproachDetector(newMaxCheck, newThreshold, newMaxIter, newHandler,
275 provider);
276 }
277
278 }
279
280 @Test
281 public void testWrappedException() {
282 final Throwable dummyCause = new RuntimeException();
283 try {
284
285 Frame eme2000 = FramesFactory.getEME2000();
286 TimeScale utc = TimeScalesFactory.getUTC();
287 final AbsoluteDate initialDate = new AbsoluteDate(2011, 5, 11, utc);
288 final AbsoluteDate exceptionDate = initialDate.shiftedBy(3600.0);
289 KeplerianPropagator k =
290 new KeplerianPropagator(new EquinoctialOrbit(new PVCoordinates(new Vector3D(4008462.4706055815, -3155502.5373837613, -5044275.9880020910),
291 new Vector3D(-5012.9298276860990, 1920.3567095973078, -5172.7403501801580)),
292 eme2000, initialDate, Constants.WGS84_EARTH_MU));
293 k.addEventDetector(new DateDetector(initialDate.shiftedBy(Constants.JULIAN_DAY)) {
294 @Override
295 public double g(final SpacecraftState s) {
296 final double dt = s.getDate().durationFrom(exceptionDate);
297 if (FastMath.abs(dt) < 1.0) {
298 throw new OrekitException(dummyCause, LocalizedCoreFormats.SIMPLE_MESSAGE, "dummy");
299 }
300 return dt;
301 }
302 });
303 k.propagate(initialDate.shiftedBy(Constants.JULIAN_YEAR));
304 Assert.fail("an exception should have been thrown");
305 } catch (OrekitException oe) {
306 Assert.assertSame(OrekitException.class, oe.getClass());
307 Assert.assertSame(dummyCause, oe.getCause().getCause());
308 String expected = "failed to find root between 2011-05-11T00:00:00.000Z " +
309 "(g=-3.6E3) and 2012-05-10T06:00:00.000Z (g=3.1554E7)\n" +
310 "Last iteration at 2011-05-11T00:00:00.000Z (g=-3.6E3)";
311 MatcherAssert.assertThat(oe.getMessage(Locale.US),
312 CoreMatchers.containsString(expected));
313 }
314 }
315
316 @Test
317 public void testDefaultMethods() {
318 EventDetector dummyDetector = new EventDetector() {
319
320 @Override
321 public double getThreshold() {
322 return 1.0e-10;
323 }
324
325 @Override
326 public int getMaxIterationCount() {
327 return 100;
328 }
329
330 @Override
331 public double getMaxCheckInterval() {
332 return 60;
333 }
334
335 @Override
336 public double g(SpacecraftState s) {
337 return s.getDate().durationFrom(AbsoluteDate.J2000_EPOCH);
338 }
339
340 @Override
341 public Action eventOccurred(SpacecraftState s, boolean increasing) {
342 return Action.RESET_STATE;
343 }
344 };
345
346
347 dummyDetector.init(null, null);
348
349
350 SpacecraftState s = new SpacecraftState(new KeplerianOrbit(7e6, 0.01, 0.3, 0, 0, 0,
351 PositionAngle.TRUE, FramesFactory.getEME2000(),
352 AbsoluteDate.J2000_EPOCH, Constants.EIGEN5C_EARTH_MU));
353 Assert.assertSame(s, dummyDetector.resetState(s));
354
355 }
356
357 @Test
358 public void testNumericalNoiseAtIntervalEnd() {
359
360 Frame eme2000 = FramesFactory.getEME2000();
361 TimeScale utc = TimeScalesFactory.getUTC();
362 final AbsoluteDate initialDate = new AbsoluteDate(2011, 5, 11, utc);
363 final Orbit orbit = new EquinoctialOrbit(new PVCoordinates(new Vector3D(4008462.4706055815, -3155502.5373837613, -5044275.9880020910),
364 new Vector3D(-5012.9298276860990, 1920.3567095973078, -5172.7403501801580)),
365 eme2000, initialDate, Constants.WGS84_EARTH_MU);
366 Propagator propagator = new KeplerianPropagator(orbit);
367 double base = 3600.0;
368 double noise = FastMath.scalb(base, -60);
369
370 AbsoluteDate finalTime = initialDate.shiftedBy(base).shiftedBy(2 * noise);
371 AbsoluteDate eventTime = finalTime.shiftedBy(-noise);
372 propagator.addEventDetector(new DateDetector(eventTime).withMaxCheck(base).withThreshold(noise / 2));
373 SpacecraftState finalState = propagator.propagate(finalTime);
374 Assert.assertEquals(0.0, finalState.getDate().durationFrom(eventTime), noise);
375
376 }
377
378 @Test
379 public void testForwardAnalytical() {
380 doTestScheduling(0.0, 1.0, 21, this::buildAnalytical);
381 }
382
383 @Test
384 public void testBackwardAnalytical() {
385 doTestScheduling(1.0, 0.0, 21, this::buildAnalytical);
386 }
387
388 @Test
389 public void testForwardNumerical() {
390 doTestScheduling(0.0, 1.0, 23, this::buildNumerical);
391 }
392
393 @Test
394 public void testBackwardNumerical() {
395 doTestScheduling(1.0, 0.0, 23, this::buildNumerical);
396 }
397
398 private Propagator buildAnalytical(final Orbit orbit) {
399 return new KeplerianPropagator(orbit);
400 }
401
402 private Propagator buildNumerical(final Orbit orbit) {
403 OrbitType type = OrbitType.CARTESIAN;
404 double[][] tol = NumericalPropagator.tolerances(0.0001, orbit, type);
405 ODEIntegrator integrator = new DormandPrince853Integrator(0.0001, 10.0, tol[0], tol[1]);
406 NumericalPropagator propagator = new NumericalPropagator(integrator);
407 propagator.setOrbitType(type);
408 propagator.setInitialState(new SpacecraftState(orbit));
409 return propagator;
410 }
411
412 private void doTestScheduling(final double start, final double stop, final int expectedCalls,
413 final Function<Orbit, Propagator> propagatorBuilder) {
414
415
416 Frame eme2000 = FramesFactory.getEME2000();
417 TimeScale utc = TimeScalesFactory.getUTC();
418 final AbsoluteDate initialDate = new AbsoluteDate(2011, 5, 11, utc);
419 final Orbit orbit = new EquinoctialOrbit(new PVCoordinates(new Vector3D(4008462.4706055815, -3155502.5373837613, -5044275.9880020910),
420 new Vector3D(-5012.9298276860990, 1920.3567095973078, -5172.7403501801580)),
421 eme2000, initialDate, Constants.WGS84_EARTH_MU);
422 Propagator propagator = propagatorBuilder.apply(orbit.shiftedBy(start));
423
424
425
426 final ScheduleChecker checker = new ScheduleChecker(initialDate.shiftedBy(start),
427 initialDate.shiftedBy(stop));
428 propagator.setStepHandler((interpolator) -> {
429 checker.callDate(interpolator.getCurrentState().getDate());
430 });
431
432 for (int i = 0; i < 10; ++i) {
433 propagator.addEventDetector(new DateDetector(initialDate.shiftedBy(0.0625 * (i + 1))).
434 withHandler((state, detector, increasing) -> {
435 checker.callDate(state.getDate());
436 return Action.CONTINUE;
437 }));
438 }
439
440 propagator.propagate(initialDate.shiftedBy(start), initialDate.shiftedBy(stop));
441
442 Assert.assertEquals(expectedCalls, checker.calls);
443
444 }
445
446
447 private static class ScheduleChecker {
448
449 private final AbsoluteDate start;
450 private final AbsoluteDate stop;
451 private AbsoluteDate last;
452 private int calls;
453
454 ScheduleChecker(final AbsoluteDate start, final AbsoluteDate stop) {
455 this.start = start;
456 this.stop = stop;
457 this.last = null;
458 this.calls = 0;
459 }
460
461 void callDate(final AbsoluteDate date) {
462 if (last != null) {
463
464 if (start.isBefore(stop)) {
465
466 Assert.assertTrue(date.isAfterOrEqualTo(start));
467 Assert.assertTrue(date.isBeforeOrEqualTo(stop));
468 Assert.assertTrue(date.isAfterOrEqualTo(last));
469 } else {
470
471 Assert.assertTrue(date.isBeforeOrEqualTo(start));
472 Assert.assertTrue(date.isAfterOrEqualTo(stop));
473 Assert.assertTrue(date.isBeforeOrEqualTo(last));
474 }
475 }
476 last = date;
477 ++calls;
478 }
479
480 }
481
482 @Before
483 public void setUp() {
484 Utils.setDataRoot("regular-data");
485 mu = Constants.EIGEN5C_EARTH_MU;
486 }
487
488 }
489