1   /* Copyright 2002-2024 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  import java.util.ArrayList;
20  import java.util.Collection;
21  import java.util.Collections;
22  import java.util.Comparator;
23  import java.util.List;
24  import java.util.PriorityQueue;
25  import java.util.Queue;
26  
27  import org.hipparchus.exception.MathRuntimeException;
28  import org.hipparchus.ode.events.Action;
29  import org.orekit.attitudes.Attitude;
30  import org.orekit.attitudes.AttitudeProvider;
31  import org.orekit.errors.OrekitException;
32  import org.orekit.errors.OrekitIllegalArgumentException;
33  import org.orekit.errors.OrekitInternalError;
34  import org.orekit.errors.OrekitMessages;
35  import org.orekit.frames.Frame;
36  import org.orekit.orbits.Orbit;
37  import org.orekit.propagation.AbstractPropagator;
38  import org.orekit.propagation.AdditionalStateProvider;
39  import org.orekit.propagation.BoundedPropagator;
40  import org.orekit.propagation.EphemerisGenerator;
41  import org.orekit.propagation.MatricesHarvester;
42  import org.orekit.propagation.SpacecraftState;
43  import org.orekit.propagation.events.EventDetector;
44  import org.orekit.propagation.events.EventState;
45  import org.orekit.propagation.events.EventState.EventOccurrence;
46  import org.orekit.propagation.sampling.OrekitStepInterpolator;
47  import org.orekit.time.AbsoluteDate;
48  import org.orekit.utils.PVCoordinatesProvider;
49  import org.orekit.utils.TimeStampedPVCoordinates;
50  
51  /** Common handling of {@link org.orekit.propagation.Propagator} methods for analytical propagators.
52   * <p>
53   * This abstract class allows to provide easily the full set of {@link
54   * org.orekit.propagation.Propagator Propagator} methods, including all propagation
55   * modes support and discrete events support for any simple propagation method. Only
56   * two methods must be implemented by derived classes: {@link #propagateOrbit(AbsoluteDate)}
57   * and {@link #getMass(AbsoluteDate)}. The first method should perform straightforward
58   * propagation starting from some internally stored initial state up to the specified target date.
59   * </p>
60   * @author Luc Maisonobe
61   */
62  public abstract class AbstractAnalyticalPropagator extends AbstractPropagator {
63  
64      /** Provider for attitude computation. */
65      private PVCoordinatesProvider pvProvider;
66  
67      /** Start date of last propagation. */
68      private AbsoluteDate lastPropagationStart;
69  
70      /** End date of last propagation. */
71      private AbsoluteDate lastPropagationEnd;
72  
73      /** Initialization indicator of events states. */
74      private boolean statesInitialized;
75  
76      /** Indicator for last step. */
77      private boolean isLastStep;
78  
79      /** Event steps. */
80      private final Collection<EventState<?>> eventsStates;
81  
82      /** Build a new instance.
83       * @param attitudeProvider provider for attitude computation
84       */
85      protected AbstractAnalyticalPropagator(final AttitudeProvider attitudeProvider) {
86          setAttitudeProvider(attitudeProvider);
87          pvProvider           = new LocalPVProvider();
88          lastPropagationStart = AbsoluteDate.PAST_INFINITY;
89          lastPropagationEnd   = AbsoluteDate.FUTURE_INFINITY;
90          statesInitialized    = false;
91          eventsStates         = new ArrayList<>();
92      }
93  
94      /** {@inheritDoc} */
95      @Override
96      public EphemerisGenerator getEphemerisGenerator() {
97          return () -> new BoundedPropagatorView(lastPropagationStart, lastPropagationEnd);
98      }
99  
100     /** {@inheritDoc} */
101     public <T extends EventDetector> void addEventDetector(final T detector) {
102         eventsStates.add(new EventState<>(detector));
103     }
104 
105     /** {@inheritDoc} */
106     public Collection<EventDetector> getEventsDetectors() {
107         final List<EventDetector> list = new ArrayList<>();
108         for (final EventState<?> state : eventsStates) {
109             list.add(state.getEventDetector());
110         }
111         return Collections.unmodifiableCollection(list);
112     }
113 
114     /** {@inheritDoc} */
115     public void clearEventsDetectors() {
116         eventsStates.clear();
117     }
118 
119     /** {@inheritDoc} */
120     public SpacecraftState propagate(final AbsoluteDate start, final AbsoluteDate target) {
121         checkStartDateIsNotInfinity(start);
122         try {
123             initializePropagation();
124 
125             lastPropagationStart = start;
126 
127             // Initialize additional states
128             initializeAdditionalStates(target);
129 
130             final boolean isForward = target.compareTo(start) >= 0;
131             SpacecraftState state   = updateAdditionalStates(basicPropagate(start));
132 
133             // initialize event detectors
134             for (final EventState<?> es : eventsStates) {
135                 es.init(state, target);
136             }
137 
138             // initialize step handlers
139             getMultiplexer().init(state, target);
140 
141             // iterate over the propagation range, need loop due to reset events
142             statesInitialized = false;
143             isLastStep = false;
144             do {
145 
146                 // attempt to advance to the target date
147                 final SpacecraftState previous = state;
148                 final SpacecraftState current = updateAdditionalStates(basicPropagate(target));
149                 final OrekitStepInterpolator interpolator =
150                         new BasicStepInterpolator(isForward, previous, current);
151 
152                 // accept the step, trigger events and step handlers
153                 state = acceptStep(interpolator, target);
154 
155                 // Update the potential changes in the spacecraft state due to the events
156                 // especially the potential attitude transition
157                 state = updateAdditionalStates(basicPropagate(state.getDate()));
158 
159             } while (!isLastStep);
160 
161             // Finalize event detectors
162             for (final EventState<?> es : eventsStates) {
163                 es.finish(state);
164             }
165 
166             // return the last computed state
167             lastPropagationEnd = state.getDate();
168             setStartDate(state.getDate());
169             return state;
170 
171         } catch (MathRuntimeException mrte) {
172             throw OrekitException.unwrap(mrte);
173         }
174     }
175 
176     /**
177      * Check the starting date is not {@code AbsoluteDate.PAST_INFINITY} or {@code AbsoluteDate.FUTURE_INFINITY}.
178      * @param start propagation starting date
179      */
180     private void checkStartDateIsNotInfinity(final AbsoluteDate start) {
181         if (start.isEqualTo(AbsoluteDate.PAST_INFINITY) || start.isEqualTo(AbsoluteDate.FUTURE_INFINITY)) {
182             throw new OrekitIllegalArgumentException(OrekitMessages.CANNOT_START_PROPAGATION_FROM_INFINITY);
183         }
184     }
185 
186     /** Accept a step, triggering events and step handlers.
187      * @param interpolator interpolator for the current step
188      * @param target final propagation time
189      * @return state at the end of the step
190      * @exception MathRuntimeException if an event cannot be located
191      */
192     protected SpacecraftState acceptStep(final OrekitStepInterpolator interpolator,
193                                          final AbsoluteDate target)
194         throws MathRuntimeException {
195 
196         SpacecraftState        previous   = interpolator.getPreviousState();
197         final SpacecraftState  current    = interpolator.getCurrentState();
198         OrekitStepInterpolator restricted = interpolator;
199 
200 
201         // initialize the events states if needed
202         if (!statesInitialized) {
203 
204             if (!eventsStates.isEmpty()) {
205                 // initialize the events states
206                 for (final EventState<?> state : eventsStates) {
207                     state.reinitializeBegin(interpolator);
208                 }
209             }
210 
211             statesInitialized = true;
212 
213         }
214 
215         // search for next events that may occur during the step
216         final int orderingSign = interpolator.isForward() ? +1 : -1;
217         final Queue<EventState<?>> occurringEvents = new PriorityQueue<>(new Comparator<EventState<?>>() {
218             /** {@inheritDoc} */
219             @Override
220             public int compare(final EventState<?> es0, final EventState<?> es1) {
221                 return orderingSign * es0.getEventDate().compareTo(es1.getEventDate());
222             }
223         });
224 
225         boolean doneWithStep = false;
226         resetEvents:
227         do {
228 
229             // Evaluate all event detectors for events
230             occurringEvents.clear();
231             for (final EventState<?> state : eventsStates) {
232                 if (state.evaluateStep(interpolator)) {
233                     // the event occurs during the current step
234                     occurringEvents.add(state);
235                 }
236             }
237 
238             do {
239 
240                 eventLoop:
241                 while (!occurringEvents.isEmpty()) {
242 
243                     // handle the chronologically first event
244                     final EventState<?> currentEvent = occurringEvents.poll();
245 
246                     // get state at event time
247                     SpacecraftState eventState = restricted.getInterpolatedState(currentEvent.getEventDate());
248 
249                     // restrict the interpolator to the first part of the step, up to the event
250                     restricted = restricted.restrictStep(previous, eventState);
251 
252                     // try to advance all event states to current time
253                     for (final EventState<?> state : eventsStates) {
254                         if (state != currentEvent && state.tryAdvance(eventState, interpolator)) {
255                             // we need to handle another event first
256                             // remove event we just updated to prevent heap corruption
257                             occurringEvents.remove(state);
258                             // add it back to update its position in the heap
259                             occurringEvents.add(state);
260                             // re-queue the event we were processing
261                             occurringEvents.add(currentEvent);
262                             continue eventLoop;
263                         }
264                     }
265                     // all event detectors agree we can advance to the current event time
266 
267                     // handle the first part of the step, up to the event
268                     getMultiplexer().handleStep(restricted);
269 
270                     // acknowledge event occurrence
271                     final EventOccurrence occurrence = currentEvent.doEvent(eventState);
272                     final Action action = occurrence.getAction();
273                     isLastStep = action == Action.STOP;
274 
275                     if (isLastStep) {
276 
277                         // ensure the event is after the root if it is returned STOP
278                         // this lets the user integrate to a STOP event and then restart
279                         // integration from the same time.
280                         final SpacecraftState savedState = eventState;
281                         eventState = interpolator.getInterpolatedState(occurrence.getStopDate());
282                         restricted = restricted.restrictStep(savedState, eventState);
283 
284                         // handle the almost zero size last part of the final step, at event time
285                         getMultiplexer().handleStep(restricted);
286                         getMultiplexer().finish(restricted.getCurrentState());
287 
288                     }
289 
290                     if (isLastStep) {
291                         // the event asked to stop integration
292                         return eventState;
293                     }
294 
295                     if (action == Action.RESET_DERIVATIVES || action == Action.RESET_STATE) {
296                         // some event handler has triggered changes that
297                         // invalidate the derivatives, we need to recompute them
298                         final SpacecraftState resetState = occurrence.getNewState();
299                         resetIntermediateState(resetState, interpolator.isForward());
300                         return resetState;
301                     }
302                     // at this point action == Action.CONTINUE or Action.RESET_EVENTS
303 
304                     // prepare handling of the remaining part of the step
305                     previous = eventState;
306                     restricted = new BasicStepInterpolator(restricted.isForward(), eventState, current);
307 
308                     if (action == Action.RESET_EVENTS) {
309                         continue resetEvents;
310                     }
311 
312                     // at this point action == Action.CONTINUE
313                     // check if the same event occurs again in the remaining part of the step
314                     if (currentEvent.evaluateStep(restricted)) {
315                         // the event occurs during the current step
316                         occurringEvents.add(currentEvent);
317                     }
318 
319                 }
320 
321                 // last part of the step, after the last event. Advance all detectors to
322                 // the end of the step. Should only detect a new event here if an event
323                 // modified the g function of another detector. Detecting such events here
324                 // is unreliable and RESET_EVENTS should be used instead. Might as well
325                 // re-check here because we have to loop through all the detectors anyway
326                 // and the alternative is to throw an exception.
327                 for (final EventState<?> state : eventsStates) {
328                     if (state.tryAdvance(current, interpolator)) {
329                         occurringEvents.add(state);
330                     }
331                 }
332 
333             } while (!occurringEvents.isEmpty());
334 
335             doneWithStep = true;
336         } while (!doneWithStep);
337 
338         isLastStep = target.equals(current.getDate());
339 
340         // handle the remaining part of the step, after all events if any
341         getMultiplexer().handleStep(restricted);
342         if (isLastStep) {
343             getMultiplexer().finish(restricted.getCurrentState());
344         }
345 
346         return current;
347 
348     }
349 
350     /** Get the mass.
351      * @param date target date for the orbit
352      * @return mass mass
353      */
354     protected abstract double getMass(AbsoluteDate date);
355 
356     /** Get PV coordinates provider.
357      * @return PV coordinates provider
358      */
359     public PVCoordinatesProvider getPvProvider() {
360         return pvProvider;
361     }
362 
363     /** Reset an intermediate state.
364      * @param state new intermediate state to consider
365      * @param forward if true, the intermediate state is valid for
366      * propagations after itself
367      */
368     protected abstract void resetIntermediateState(SpacecraftState state, boolean forward);
369 
370     /** Extrapolate an orbit up to a specific target date.
371      * @param date target date for the orbit
372      * @return extrapolated parameters
373      */
374     protected abstract Orbit propagateOrbit(AbsoluteDate date);
375 
376     /** Propagate an orbit without any fancy features.
377      * <p>This method is similar in spirit to the {@link #propagate} method,
378      * except that it does <strong>not</strong> call any handler during
379      * propagation, nor any discrete events, not additional states. It always
380      * stop exactly at the specified date.</p>
381      * @param date target date for propagation
382      * @return state at specified date
383      */
384     protected SpacecraftState basicPropagate(final AbsoluteDate date) {
385         try {
386 
387             // evaluate orbit
388             final Orbit orbit = propagateOrbit(date);
389 
390             // evaluate attitude
391             final Attitude attitude =
392                 getAttitudeProvider().getAttitude(pvProvider, date, orbit.getFrame());
393 
394             // build raw state
395             return new SpacecraftState(orbit, attitude, getMass(date));
396 
397         } catch (OrekitException oe) {
398             throw new OrekitException(oe);
399         }
400     }
401 
402     /**
403      * Get the names of the parameters in the matrix returned by {@link MatricesHarvester#getParametersJacobian}.
404      * @return names of the parameters (i.e. columns) of the Jacobian matrix
405      * @since 11.1
406      */
407     protected List<String> getJacobiansColumnsNames() {
408         return Collections.emptyList();
409     }
410 
411     /** Internal PVCoordinatesProvider for attitude computation. */
412     private class LocalPVProvider implements PVCoordinatesProvider {
413 
414         /** {@inheritDoc} */
415         public TimeStampedPVCoordinates getPVCoordinates(final AbsoluteDate date, final Frame frame) {
416             return propagateOrbit(date).getPVCoordinates(frame);
417         }
418 
419     }
420 
421     /** {@link BoundedPropagator} view of the instance. */
422     private class BoundedPropagatorView extends AbstractAnalyticalPropagator implements BoundedPropagator {
423 
424         /** Min date. */
425         private final AbsoluteDate minDate;
426 
427         /** Max date. */
428         private final AbsoluteDate maxDate;
429 
430         /** Simple constructor.
431          * @param startDate start date of the propagation
432          * @param endDate end date of the propagation
433          */
434         BoundedPropagatorView(final AbsoluteDate startDate, final AbsoluteDate endDate) {
435             super(AbstractAnalyticalPropagator.this.getAttitudeProvider());
436             super.resetInitialState(AbstractAnalyticalPropagator.this.getInitialState());
437             if (startDate.compareTo(endDate) <= 0) {
438                 minDate = startDate;
439                 maxDate = endDate;
440             } else {
441                 minDate = endDate;
442                 maxDate = startDate;
443             }
444 
445             try {
446                 // copy the same additional state providers as the original propagator
447                 for (AdditionalStateProvider provider : AbstractAnalyticalPropagator.this.getAdditionalStateProviders()) {
448                     addAdditionalStateProvider(provider);
449                 }
450             } catch (OrekitException oe) {
451                 // as the generators are already compatible with each other,
452                 // this should never happen
453                 throw new OrekitInternalError(null);
454             }
455 
456         }
457 
458         /** {@inheritDoc} */
459         public AbsoluteDate getMinDate() {
460             return minDate;
461         }
462 
463         /** {@inheritDoc} */
464         public AbsoluteDate getMaxDate() {
465             return maxDate;
466         }
467 
468         /** {@inheritDoc} */
469         protected Orbit propagateOrbit(final AbsoluteDate target) {
470             return AbstractAnalyticalPropagator.this.propagateOrbit(target);
471         }
472 
473         /** {@inheritDoc} */
474         public double getMass(final AbsoluteDate date) {
475             return AbstractAnalyticalPropagator.this.getMass(date);
476         }
477 
478         /** {@inheritDoc} */
479         public void resetInitialState(final SpacecraftState state) {
480             super.resetInitialState(state);
481             AbstractAnalyticalPropagator.this.resetInitialState(state);
482         }
483 
484         /** {@inheritDoc} */
485         protected void resetIntermediateState(final SpacecraftState state, final boolean forward) {
486             AbstractAnalyticalPropagator.this.resetIntermediateState(state, forward);
487         }
488 
489         /** {@inheritDoc} */
490         public SpacecraftState getInitialState() {
491             return AbstractAnalyticalPropagator.this.getInitialState();
492         }
493 
494         /** {@inheritDoc} */
495         public Frame getFrame() {
496             return AbstractAnalyticalPropagator.this.getFrame();
497         }
498 
499     }
500 
501     /** Internal class for local propagation. */
502     private class BasicStepInterpolator implements OrekitStepInterpolator {
503 
504         /** Previous state. */
505         private final SpacecraftState previousState;
506 
507         /** Current state. */
508         private final SpacecraftState currentState;
509 
510         /** Forward propagation indicator. */
511         private final boolean forward;
512 
513         /** Simple constructor.
514          * @param isForward integration direction indicator
515          * @param previousState start of the step
516          * @param currentState end of the step
517          */
518         BasicStepInterpolator(final boolean isForward,
519                               final SpacecraftState previousState,
520                               final SpacecraftState currentState) {
521             this.forward         = isForward;
522             this.previousState   = previousState;
523             this.currentState    = currentState;
524         }
525 
526         /** {@inheritDoc} */
527         @Override
528         public SpacecraftState getPreviousState() {
529             return previousState;
530         }
531 
532         /** {@inheritDoc} */
533         @Override
534         public boolean isPreviousStateInterpolated() {
535             // no difference in analytical propagators
536             return false;
537         }
538 
539         /** {@inheritDoc} */
540         @Override
541         public SpacecraftState getCurrentState() {
542             return currentState;
543         }
544 
545         /** {@inheritDoc} */
546         @Override
547         public boolean isCurrentStateInterpolated() {
548             // no difference in analytical propagators
549             return false;
550         }
551 
552         /** {@inheritDoc} */
553         @Override
554         public SpacecraftState getInterpolatedState(final AbsoluteDate date) {
555 
556             // compute the basic spacecraft state
557             final SpacecraftState basicState = basicPropagate(date);
558 
559             // add the additional states
560             return updateAdditionalStates(basicState);
561 
562         }
563 
564         /** {@inheritDoc} */
565         @Override
566         public boolean isForward() {
567             return forward;
568         }
569 
570         /** {@inheritDoc} */
571         @Override
572         public BasicStepInterpolator restrictStep(final SpacecraftState newPreviousState,
573                                                   final SpacecraftState newCurrentState) {
574             return new BasicStepInterpolator(forward, newPreviousState, newCurrentState);
575         }
576 
577     }
578 
579 }