1   /* Copyright 2022-2024 Romain Serra
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.intervals;
18  
19  import org.hipparchus.util.FastMath;
20  import org.hipparchus.util.MathUtils;
21  import org.orekit.orbits.KeplerianOrbit;
22  import org.orekit.orbits.Orbit;
23  import org.orekit.orbits.OrbitType;
24  import org.orekit.propagation.events.AdaptableInterval;
25  
26  /**
27   * Factory class for {@link AdaptableInterval} suitable for apside detection on eccentric orbits.
28   * It requires {@link org.orekit.propagation.SpacecraftState} to be based on {@link Orbit} in order to work.
29   * @see org.orekit.propagation.events.AdaptableInterval
30   * @see org.orekit.propagation.events.ApsideDetector
31   * @see org.orekit.propagation.events.EventSlopeFilter
32   * @author Romain Serra
33   * @since 12.1
34   */
35  public class ApsideDetectionAdaptableIntervalFactory {
36  
37      /**
38       * Private constructor.
39       */
40      private ApsideDetectionAdaptableIntervalFactory() {
41          // factory class
42      }
43  
44      /**
45       * Method providing a candidate {@link AdaptableInterval} for arbitrary apside detection with forward propagation.
46       * It uses a Keplerian, eccentric approximation.
47       * @return adaptable interval for forward apside detection
48       */
49      public static AdaptableInterval getForwardApsideDetectionAdaptableInterval() {
50          return state -> {
51              final Orbit orbit = state.getOrbit();
52              final KeplerianOrbit keplerianOrbit = convertOrbitIntoKeplerianOne(orbit);
53              final double meanMotion = keplerianOrbit.getKeplerianMeanMotion();
54              final double meanAnomaly = keplerianOrbit.getMeanAnomaly();
55              final double durationToNextPeriapsis = computeKeplerianDurationToNextPeriapsis(meanAnomaly, meanMotion);
56              final double durationToNextApoapsis = computeKeplerianDurationToNextApoapsis(meanAnomaly, meanMotion);
57              return FastMath.min(durationToNextPeriapsis, durationToNextApoapsis);
58          };
59      }
60  
61      /**
62       * Method providing a candidate {@link AdaptableInterval} for arbitrary apside detection with backward propagation.
63       * It uses a Keplerian, eccentric approximation.
64       * @return adaptable interval for backward apside detection
65       */
66      public static AdaptableInterval getBackwardApsideDetectionAdaptableInterval() {
67          return state -> {
68              final Orbit orbit = state.getOrbit();
69              final KeplerianOrbit keplerianOrbit = convertOrbitIntoKeplerianOne(orbit);
70              final double meanMotion = keplerianOrbit.getKeplerianMeanMotion();
71              final double meanAnomaly = keplerianOrbit.getMeanAnomaly();
72              final double durationFromPreviousPeriapsis = computeKeplerianDurationFromPreviousPeriapsis(meanAnomaly,
73                      meanMotion);
74              final double durationFromPreviousApoapsis = computeKeplerianDurationFromPreviousApoapsis(meanAnomaly,
75                      meanMotion);
76              return FastMath.min(durationFromPreviousApoapsis, durationFromPreviousPeriapsis);
77          };
78      }
79  
80      /**
81       * Method providing a candidate {@link AdaptableInterval} for periapsis detection with forward propagation.
82       * It uses a Keplerian, eccentric approximation.
83       * @return adaptable interval for forward periaspsis detection
84       */
85      public static AdaptableInterval getForwardPeriapsisDetectionAdaptableInterval() {
86          return state -> {
87              final Orbit orbit = state.getOrbit();
88              final KeplerianOrbit keplerianOrbit = convertOrbitIntoKeplerianOne(orbit);
89              final double meanMotion = keplerianOrbit.getKeplerianMeanMotion();
90              final double meanAnomaly = keplerianOrbit.getMeanAnomaly();
91              return computeKeplerianDurationToNextPeriapsis(meanAnomaly, meanMotion);
92          };
93      }
94  
95      /**
96       * Method providing a candidate {@link AdaptableInterval} for periapsis detection with backward propagation.
97       * It uses a Keplerian, eccentric approximation.
98       * @return adaptable interval for backward periaspsis detection
99       */
100     public static AdaptableInterval getBackwardPeriapsisDetectionAdaptableInterval() {
101         return state -> {
102             final Orbit orbit = state.getOrbit();
103             final KeplerianOrbit keplerianOrbit = convertOrbitIntoKeplerianOne(orbit);
104             final double meanMotion = keplerianOrbit.getKeplerianMeanMotion();
105             final double meanAnomaly = keplerianOrbit.getMeanAnomaly();
106             return computeKeplerianDurationFromPreviousPeriapsis(meanAnomaly, meanMotion);
107         };
108     }
109 
110     /**
111      * Method providing a candidate {@link AdaptableInterval} for apoapsis detection with forward propagation.
112      * It uses a Keplerian, eccentric approximation.
113      * @return adaptable interval for forward apoapsis detection
114      */
115     public static AdaptableInterval getForwardApoapsisDetectionAdaptableInterval() {
116         return state -> {
117             final Orbit orbit = state.getOrbit();
118             final KeplerianOrbit keplerianOrbit = convertOrbitIntoKeplerianOne(orbit);
119             final double meanMotion = keplerianOrbit.getKeplerianMeanMotion();
120             final double meanAnomaly = keplerianOrbit.getMeanAnomaly();
121             return computeKeplerianDurationToNextApoapsis(meanAnomaly, meanMotion);
122         };
123     }
124 
125     /**
126      * Method providing a candidate {@link AdaptableInterval} for apoapsis detection with backward propagation.
127      * It uses a Keplerian, eccentric approximation.
128      * @return adaptable interval for backward apoapsis detection
129      */
130     public static AdaptableInterval getBackwardApoapsisDetectionAdaptableInterval() {
131         return state -> {
132             final Orbit orbit = state.getOrbit();
133             final KeplerianOrbit keplerianOrbit = convertOrbitIntoKeplerianOne(orbit);
134             final double meanMotion = keplerianOrbit.getKeplerianMeanMotion();
135             final double meanAnomaly = keplerianOrbit.getMeanAnomaly();
136             return computeKeplerianDurationFromPreviousApoapsis(meanAnomaly, meanMotion);
137         };
138     }
139 
140     /**
141      * Convert a generic {@link Orbit} into a {@link KeplerianOrbit}.
142      * @param orbit orbit to convert
143      * @return Keplerian orbit
144      */
145     private static KeplerianOrbit convertOrbitIntoKeplerianOne(final Orbit orbit) {
146         return (KeplerianOrbit) OrbitType.KEPLERIAN.convertType(orbit);
147     }
148 
149     /**
150      * Method computing time to go until next periapsis, assuming Keplerian motion.
151      * @param meanAnomaly mean anomaly
152      * @param meanMotion Keplerian mean motion
153      * @return duration to next periapsis
154      */
155     private static double computeKeplerianDurationToNextPeriapsis(final double meanAnomaly,
156                                                                   final double meanMotion) {
157         final double normalizedMeanAnomaly = MathUtils.normalizeAngle(meanAnomaly, FastMath.PI);
158         return (MathUtils.TWO_PI - normalizedMeanAnomaly) / meanMotion;
159     }
160 
161     /**
162      * Method computing time elapsed since last periapsis, assuming Keplerian motion.
163      * @param meanAnomaly mean anomaly
164      * @param meanMotion Keplerian mean motion
165      * @return duration elapsed since last periapsis
166      */
167     public static double computeKeplerianDurationFromPreviousPeriapsis(final double meanAnomaly,
168                                                                        final double meanMotion) {
169         final double normalizedMeanAnomaly = MathUtils.normalizeAngle(meanAnomaly, FastMath.PI);
170         return normalizedMeanAnomaly / meanMotion;
171     }
172 
173     /**
174      * Method computing time to go until next apoapsis, assuming Keplerian motion.
175      * @param meanAnomaly mean anomaly
176      * @param meanMotion Keplerian mean motion
177      * @return duration to next apoapsis
178      */
179     private static double computeKeplerianDurationToNextApoapsis(final double meanAnomaly,
180                                                                  final double meanMotion) {
181         final double normalizedMeanAnomaly = MathUtils.normalizeAngle(meanAnomaly, MathUtils.TWO_PI);
182         return (MathUtils.TWO_PI + FastMath.PI - normalizedMeanAnomaly) / meanMotion;
183     }
184 
185     /**
186      * Method computing time elapsed since last apoapsis, assuming Keplerian motion.
187      * @param meanAnomaly mean anomaly
188      * @param meanMotion Keplerian mean motion
189      * @return duration elapsed since last apoapsis
190      */
191     public static double computeKeplerianDurationFromPreviousApoapsis(final double meanAnomaly,
192                                                                       final double meanMotion) {
193         final double normalizedMeanAnomaly = MathUtils.normalizeAngle(meanAnomaly, MathUtils.TWO_PI);
194         return (normalizedMeanAnomaly - FastMath.PI) / meanMotion;
195     }
196 }