1   /* Copyright 2002-2026 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.events;
18  
19  import org.hipparchus.ode.events.Action;
20  import org.hipparchus.util.FastMath;
21  import org.hipparchus.util.MathUtils;
22  import org.orekit.frames.Frame;
23  import org.orekit.orbits.KeplerianOrbit;
24  import org.orekit.orbits.Orbit;
25  import org.orekit.orbits.OrbitType;
26  import org.orekit.orbits.PositionAngleType;
27  import org.orekit.propagation.SpacecraftState;
28  import org.orekit.propagation.events.functions.NodeEventFunction;
29  import org.orekit.propagation.events.handlers.EventHandler;
30  import org.orekit.propagation.events.handlers.StopOnIncreasing;
31  import org.orekit.propagation.events.intervals.AdaptableInterval;
32  
33  /** Finder for node crossing events.
34   * <p>This class finds equator crossing events (i.e. ascending
35   * or descending node crossing).</p>
36   * <p>The default implementation behavior is to {@link Action#CONTINUE continue}
37   * propagation at descending node crossing and to {@link Action#STOP stop} propagation
38   * at ascending node crossing. This can be changed by calling
39   * {@link #withHandler(EventHandler)} after construction.</p>
40   * <p>Beware that node detection will fail for almost equatorial orbits. If
41   * for example a node detector is used to trigger an {@link
42   * org.orekit.forces.maneuvers.ImpulseManeuver ImpulseManeuver} and the maneuver
43   * turn the orbit plane to equator, then the detector may completely fail just
44   * after the maneuver has been performed! This is a real case that has been
45   * encountered during validation ...</p>
46   * @see org.orekit.propagation.Propagator#addEventDetector(EventDetector)
47   * @author Luc Maisonobe
48   */
49  public class NodeDetector extends AbstractDetector<NodeDetector> {
50  
51      /** Default max check interval. */
52      private static final double DEFAULT_MAX_CHECK = 1800.0;
53  
54      /** Default convergence threshold. */
55      private static final double DEFAULT_THRESHOLD = 1.0e-3;
56  
57      /** Frame in which the equator is defined. */
58      private final Frame frame;
59  
60      /** Build a new instance.
61       * <p>The default {@link #getMaxCheckInterval() max check interval}
62       * is set to 1800s, it can be changed using {@link #withMaxCheck(double)}
63       * in the fluent API. The default {@link #getThreshold() convergence threshold}
64       * is set to 1.0e-3s, it can be changed using {@link #withThreshold(double)}
65       * in the fluent API.</p>
66       * @param frame frame in which the equator is defined (typical
67       * values are {@link org.orekit.frames.FramesFactory#getEME2000() EME<sub>2000</sub>} or
68       * {@link org.orekit.frames.FramesFactory#getITRF(org.orekit.utils.IERSConventions, boolean) ITRF})
69       * @since 10.3
70       */
71      public NodeDetector(final Frame frame) {
72          this(new EventDetectionSettings(DEFAULT_MAX_CHECK, DEFAULT_THRESHOLD, EventDetectionSettings.DEFAULT_MAX_ITER),
73                  new StopOnIncreasing(), frame);
74      }
75  
76      /** Build a new instance.
77       * <p>The orbit is used only to set an upper bound for the max check interval
78       * to a value related to nodes separation (as computed by a Keplerian model)
79       * and to set the convergence threshold according to orbit size.</p>
80       * @param orbit initial orbit
81       * @param frame frame in which the equator is defined (typical
82       * values are {@link org.orekit.frames.FramesFactory#getEME2000() EME<sub>2000</sub>} or
83       * {@link org.orekit.frames.FramesFactory#getITRF(org.orekit.utils.IERSConventions, boolean) ITRF})
84       */
85      public NodeDetector(final Orbit orbit, final Frame frame) {
86          this(1.0e-13 * orbit.getKeplerianPeriod(), orbit, frame);
87      }
88  
89      /** Build a new instance.
90       * <p>The orbit is used only to set an upper bound for the max check interval
91       * to a value related to nodes separation (as computed by a Keplerian model).</p>
92       * @param threshold convergence threshold (s)
93       * @param orbit initial orbit
94       * @param frame frame in which the equator is defined (typical
95       * values are {@link org.orekit.frames.FramesFactory#getEME2000() EME<sub>2000</sub>} or
96       * {@link org.orekit.frames.FramesFactory#getITRF(org.orekit.utils.IERSConventions, boolean) ITRF})
97       */
98      public NodeDetector(final double threshold, final Orbit orbit, final Frame frame) {
99          this(new EventDetectionSettings(AdaptableInterval.of(2 * estimateNodesTimeSeparation(orbit) / 3), threshold,
100              DEFAULT_MAX_ITER), new StopOnIncreasing(),
101              frame);
102     }
103 
104     /** Constructor with detection settings and handler.
105      * @param detectionSettings detection settings
106      * @param handler event handler to call at event occurrences
107      * @param frame frame in which the equator is defined (typical values are {@link org.orekit.frames.FramesFactory#getEME2000() EME<sub>2000</sub>} or
108      * {@link org.orekit.frames.FramesFactory#getITRF(org.orekit.utils.IERSConventions, boolean) ITRF})
109      */
110     public NodeDetector(final EventDetectionSettings detectionSettings, final EventHandler handler, final Frame frame) {
111         this(new NodeEventFunction(frame), detectionSettings, handler);
112     }
113 
114     /** Protected constructor with full parameters.
115      * <p>
116      * This constructor is not public as users are expected to use the builder
117      * API with the various {@code withXxx()} methods to set up the instance
118      * in a readable manner without using a huge amount of parameters.
119      * </p>
120      * @param nodeEventFunction event function
121      * @param detectionSettings detection settings
122      * @param handler event handler to call at event occurrences
123      * values are {@link org.orekit.frames.FramesFactory#getEME2000() EME<sub>2000</sub>} or
124      * {@link org.orekit.frames.FramesFactory#getITRF(org.orekit.utils.IERSConventions, boolean) ITRF})
125      * @since 14.0
126      */
127     protected NodeDetector(final NodeEventFunction nodeEventFunction,
128                            final EventDetectionSettings detectionSettings, final EventHandler handler) {
129         super(nodeEventFunction, detectionSettings, handler);
130         this.frame = nodeEventFunction.getFrame();
131     }
132 
133     /** {@inheritDoc} */
134     @Override
135     protected NodeDetector create(final EventDetectionSettings detectionSettings, final EventHandler newHandler) {
136         return new NodeDetector((NodeEventFunction) getEventFunction(), detectionSettings, newHandler);
137     }
138 
139     /** Find time separation between nodes.
140      * <p>
141      * The estimation of time separation is based on Keplerian motion, it is only
142      * used as a rough guess for a safe setting of default max check interval for
143      * event detection.
144      * </p>
145      * @param orbit initial orbit
146      * @return minimum time separation between nodes
147      */
148     private static double estimateNodesTimeSeparation(final Orbit orbit) {
149 
150         final KeplerianOrbit keplerian = (KeplerianOrbit) OrbitType.KEPLERIAN.convertType(orbit);
151 
152         // mean anomaly of ascending node
153         final double ascendingM  =  new KeplerianOrbit(keplerian.getA(), keplerian.getE(),
154                                                        keplerian.getI(),
155                                                        keplerian.getPerigeeArgument(),
156                                                        keplerian.getRightAscensionOfAscendingNode(),
157                                                        -keplerian.getPerigeeArgument(), PositionAngleType.TRUE,
158                                                        keplerian.getFrame(), keplerian.getDate(),
159                                                        keplerian.getMu()).getMeanAnomaly();
160 
161         // mean anomaly of descending node
162         final double descendingM =  new KeplerianOrbit(keplerian.getA(), keplerian.getE(),
163                                                        keplerian.getI(),
164                                                        keplerian.getPerigeeArgument(),
165                                                        keplerian.getRightAscensionOfAscendingNode(),
166                                                        FastMath.PI - keplerian.getPerigeeArgument(), PositionAngleType.TRUE,
167                                                        keplerian.getFrame(), keplerian.getDate(),
168                                                        keplerian.getMu()).getMeanAnomaly();
169 
170         // differences between mean anomalies
171         final double delta1 = MathUtils.normalizeAngle(ascendingM, descendingM + FastMath.PI) - descendingM;
172         final double delta2 = 2 * FastMath.PI - delta1;
173 
174         // minimum time separation between the two nodes
175         return FastMath.min(delta1, delta2) / keplerian.getKeplerianMeanMotion();
176 
177     }
178 
179     /** Get the frame in which the equator is defined.
180      * @return the frame in which the equator is defined
181      */
182     public Frame getFrame() {
183         return frame;
184     }
185 
186     /** Compute the value of the switching function.
187      * This function computes the Z position in the defined frame.
188      * @param s the current state information: date, kinematics, attitude
189      * @return value of the switching function
190      */
191     public double g(final SpacecraftState s) {
192         return getEventFunction().value(s);
193     }
194 
195 }