diff --git a/.gitignore b/.gitignore
index 85a319caa3bfbf6087430f4ee25da709a9626041..08e36c564480fa8fc5498c4102102663e977ae00 100644
--- a/.gitignore
+++ b/.gitignore
@@ -8,4 +8,4 @@ target
/lib
*.class
.idea/
-orekit.iml
\ No newline at end of file
+orekit.iml
diff --git a/src/changes/changes.xml b/src/changes/changes.xml
index 5c973b7bdfe43be0e5f7619e3586f1f562955e29..62ed35696a1a5ac239eae98b8c8653e8e473f2f7 100644
--- a/src/changes/changes.xml
+++ b/src/changes/changes.xml
@@ -33,6 +33,9 @@
+ * This class finds extremum approach events (i.e. closest or farthest approach).
+ *
+ * The default implementation behavior is to {@link Action#CONTINUE continue} propagation at farthest approach and to
+ * {@link Action#STOP stop} propagation at closest approach. This can be changed by calling
+ * {@link #withHandler(EventHandler)} after construction (go to the end of the documentation to see an example).
+ *
+ * As this detector needs two objects (moving relative to each other), it embeds one
+ * {@link PVCoordinatesProvider coordinates provider} for the secondary object and is registered as an event detector in
+ * the propagator of the primary object. The secondary object {@link PVCoordinatesProvider coordinates provider} will
+ * therefore be driven by this detector (and hence by the propagator in which this detector is registered).
+ *
+ * In order to avoid infinite recursion, care must be taken to have the secondary object provider being completely
+ * independent from anything else. In particular, if the provider is a propagator, it should not be run
+ * together in a {@link PropagatorsParallelizer propagators parallelizer} with the propagator this detector is
+ * registered in. It is fine however to configure two separate propagators PsA and PsB with similar settings for the
+ * secondary object and one propagator Pm for the primary object and then use Psa in this detector registered within Pm
+ * while Pm and Psb are run in the context of a {@link PropagatorsParallelizer propagators parallelizer}.
+ *
+ * For efficiency reason during the event search loop, it is recommended to have the secondary provider be an analytical
+ * propagator or an ephemeris. A numerical propagator as a secondary propagator works but is expected to be
+ * computationally costly.
+ *
+ * Also, it is possible to detect solely one type of event using an {@link EventSlopeFilter event slope filter}. For
+ * example in order to only detect closest approach, one should type the following :
+ * {@code
+ * ExtremumApproachDetector extremumApproachDetector = new ExtremumApproachDetector(secondaryPVProvider);
+ * EventDetector closeApproachDetector = new EventSlopeFilter
+ *
+ * By default, the implemented behavior is to {@link Action#CONTINUE continue} propagation at farthest approach and + * to {@link Action#STOP stop} propagation at closest approach. + *
+ * + * @param secondaryPVProvider PVCoordinates provider of the other object with which we want to find out the extremum + * approach. + */ + public ExtremumApproachDetector(final PVCoordinatesProvider secondaryPVProvider) { + this(DEFAULT_MAXCHECK, DEFAULT_THRESHOLD, DEFAULT_MAX_ITER, new StopOnIncreasing<>(), secondaryPVProvider); + } + + /** + * Constructor. + *+ * This constructor is to be used if the user wants to change the default behavior of the detector. + *
+ * + * @param maxCheck Maximum checking interval (s). + * @param threshold Convergence threshold (s). + * @param maxIter Maximum number of iterations in the event time search. + * @param handler Event handler to call at event occurrences. + * @param secondaryPVProvider PVCoordinates provider of the other object with which we want to find out the extremum + * * approach. + * @see EventHandler + */ + public ExtremumApproachDetector( + final double maxCheck, final double threshold, final int maxIter, + final EventHandler super ExtremumApproachDetector> handler, + final PVCoordinatesProvider secondaryPVProvider) { + super(maxCheck, threshold, maxIter, handler); + this.secondaryPVProvider = secondaryPVProvider; + } + + /** + * The {@code g} is positive when the primary object is getting further away from the secondary object and is + * negative when it is getting closer to it. + * + * @param s the current state information: date, kinematics, attitude + * @return value of the switching function + */ + public double g(final SpacecraftState s) { + final PVCoordinates deltaPV = computeDeltaPV(s); + return Vector3D.dotProduct(deltaPV.getPosition(), deltaPV.getVelocity()); + } + + /** + * Compute the relative PV between primary and secondary objects. + * + * @param s Spacecraft state. + * @return Relative position between primary (=s) and secondaryPVProvider. + */ + protected PVCoordinates computeDeltaPV(final SpacecraftState s) { + return new PVCoordinates(s.getPVCoordinates(), + secondaryPVProvider.getPVCoordinates(s.getDate(), s.getFrame())); + } + + /** {@inheritDoc} */ + @Override + protected ExtremumApproachDetector create(final double newMaxCheck, final double newThreshold, final int newMaxIter, + final EventHandler super ExtremumApproachDetector> newHandler) { + return new ExtremumApproachDetector(newMaxCheck, newThreshold, newMaxIter, newHandler, secondaryPVProvider); + } +} diff --git a/src/test/java/org/orekit/propagation/events/ExtremumApproachDetectorTest.java b/src/test/java/org/orekit/propagation/events/ExtremumApproachDetectorTest.java new file mode 100644 index 0000000000000000000000000000000000000000..faa45e496c5e3275b2bb27ae5395cd2ba9fdb79d --- /dev/null +++ b/src/test/java/org/orekit/propagation/events/ExtremumApproachDetectorTest.java @@ -0,0 +1,131 @@ +/* Copyright 2002-2022 CS GROUP + * Licensed to CS GROUP (CS) under one or more + * contributor license agreements. See the NOTICE file distributed with + * this work for additional information regarding copyright ownership. + * CS licenses this file to You under the Apache License, Version 2.0 + * (the "License"); you may not use this file except in compliance with + * the License. You may obtain a copy of the License at + * + * http://www.apache.org/licenses/LICENSE-2.0 + * + * Unless required by applicable law or agreed to in writing, software + * distributed under the License is distributed on an "AS IS" BASIS, + * WITHOUT WARRANTIES OR CONDITIONS OF ANY KIND, either express or implied. + * See the License for the specific language governing permissions and + * limitations under the License. + */ +package org.orekit.propagation.events; + +import org.hipparchus.util.FastMath; +import org.junit.jupiter.api.Assertions; +import org.junit.jupiter.api.Test; +import org.orekit.Utils; +import org.orekit.bodies.CelestialBodyFactory; +import org.orekit.frames.Frame; +import org.orekit.frames.FramesFactory; +import org.orekit.orbits.KeplerianOrbit; +import org.orekit.orbits.Orbit; +import org.orekit.orbits.PositionAngle; +import org.orekit.propagation.Propagator; +import org.orekit.propagation.SpacecraftState; +import org.orekit.propagation.analytical.KeplerianPropagator; +import org.orekit.propagation.events.handlers.StopOnEvent; +import org.orekit.time.AbsoluteDate; +import org.orekit.utils.PVCoordinatesProvider; + +public class ExtremumApproachDetectorTest { + + /** + * Test the detector on a keplerian orbit and detect extremum approach with Earth. + */ + @Test + public void testStopPropagationClosestApproachByDefault() { + // Given + // Loading Orekit data + Utils.setDataRoot("regular-data"); + + // Generating orbit + final AbsoluteDate initialDate = new AbsoluteDate(); + final Frame frame = FramesFactory.getEME2000(); + final double mu = 398600e9; //m**3/s**2 + + final double rp = (6378 + 400) * 1000; //m + final double ra = (6378 + 800) * 1000; //m + + final double a = (ra + rp) / 2; //m + final double e = (ra - rp) / (ra + rp); //m + final double i = 0; //rad + final double pa = 0; //rad + final double raan = 0; //rad + final double anomaly = FastMath.toRadians(0); //rad + final Orbit orbit = + new KeplerianOrbit(a, e, i, pa, raan, anomaly, PositionAngle.TRUE, frame, initialDate, mu); + System.out.println("Keplerian period is : "); + System.out.println(orbit.getKeplerianPeriod()); + + // Will detect extremum approaches with Earth + final PVCoordinatesProvider earthPVProvider = CelestialBodyFactory.getEarth(); + + // Initializing detector + final ExtremumApproachDetector detector = new ExtremumApproachDetector(earthPVProvider); + + // Initializing propagator + final Propagator propagator = new KeplerianPropagator(orbit); + propagator.addEventDetector(detector); + + // When + final SpacecraftState stateAtEvent = + propagator.propagate(initialDate.shiftedBy(orbit.getKeplerianPeriod() * 2)); + + // Then + Assertions.assertEquals(stateAtEvent.getDate().durationFrom(initialDate),orbit.getKeplerianPeriod(),1e-9); + + } + + /** + * Test the detector on a keplerian orbit and detect extremum approach with Earth. + */ + @Test + public void testStopPropagationFarthestApproachWithHandler() { + + // Given + // Loading Orekit data + Utils.setDataRoot("regular-data"); + + // Generating orbit + final AbsoluteDate initialDate = new AbsoluteDate(); + final Frame frame = FramesFactory.getEME2000(); + final double mu = 398600e9; //m**3/s**2 + + final double rp = (6378 + 400) * 1000; //m + final double ra = (6378 + 800) * 1000; //m + + final double a = (ra + rp) / 2; //m + final double e = (ra - rp) / (ra + rp); //m + final double i = 0; //rad + final double pa = 0; //rad + final double raan = 0; //rad + final double anomaly = FastMath.toRadians(0); //rad + final Orbit orbit = + new KeplerianOrbit(a, e, i, pa, raan, anomaly, PositionAngle.TRUE, frame, initialDate, mu); + + // Will detect extremum approaches with Earth + final PVCoordinatesProvider earthPVProvider = CelestialBodyFactory.getEarth(); + + // Initializing detector with custom handler + final ExtremumApproachDetector detector = + new ExtremumApproachDetector(earthPVProvider).withHandler(new StopOnEvent<>()); + + // Initializing propagator + final Propagator propagator = new KeplerianPropagator(orbit); + propagator.addEventDetector(detector); + + // When + final SpacecraftState stateAtEvent = + propagator.propagate(initialDate.shiftedBy(orbit.getKeplerianPeriod() * 2)); + + // Then + Assertions.assertEquals(stateAtEvent.getDate().durationFrom(initialDate),orbit.getKeplerianPeriod() / 2,1e-7); + + } +}