1   /* Copyright 2002-2022 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.forces.empirical;
18  
19  import java.util.ArrayList;
20  import java.util.Collections;
21  import java.util.List;
22  import java.util.NavigableSet;
23  import java.util.stream.Stream;
24  
25  import org.hipparchus.Field;
26  import org.hipparchus.CalculusFieldElement;
27  import org.hipparchus.geometry.euclidean.threed.FieldVector3D;
28  import org.hipparchus.geometry.euclidean.threed.Vector3D;
29  import org.hipparchus.util.MathArrays;
30  import org.orekit.attitudes.Attitude;
31  import org.orekit.attitudes.AttitudeProvider;
32  import org.orekit.attitudes.FieldAttitude;
33  import org.orekit.forces.AbstractForceModel;
34  import org.orekit.propagation.FieldSpacecraftState;
35  import org.orekit.propagation.SpacecraftState;
36  import org.orekit.propagation.events.EventDetector;
37  import org.orekit.propagation.events.FieldEventDetector;
38  import org.orekit.time.AbsoluteDate;
39  import org.orekit.time.FieldAbsoluteDate;
40  import org.orekit.utils.ParameterDriver;
41  import org.orekit.utils.TimeSpanMap;
42  import org.orekit.utils.TimeSpanMap.Span;
43  import org.orekit.utils.TimeSpanMap.Transition;
44  
45  /** Time span parametric acceleration model.
46   *  <p>
47   *  This class is closely related to {@link org.orekit.forces.empirical.ParametricAcceleration ParametricAcceleration} class.<br>
48   *  The difference is that it has a {@link TimeSpanMap} of {@link AccelerationModel} objects as attribute
49   *  instead of a single {@link AccelerationModel} object. <br>
50   *  The idea behind this model is to allow the user to design a parametric acceleration model that can see its physical parameters
51   *  change with time, at dates chosen by the user. <br>
52   *  </p>
53   *  <p>
54   *  This is a behavior that can be sought in precise orbit determination.<br>
55   *  Indeed for this type of application, the empirical parameters must be revalued at
56   *  each new orbit.
57   *  </p>
58   *  <b>Usage</b>:<ul>
59   *  <li><u>Construction</u>: constructor takes an acceleration direction, an attitude mode (or an inertial flag) and
60   *  an AccelerationModel model.<br>
61   *  This last model will be your initial AccelerationModel model and it will be initially valid for the whole time line.<br>
62   *  The real validity of this first entry will be truncated as other AccelerationModel models are added.
63   *  <li><u>Time spans</u>: AccelerationModel models are added using methods {@link #addAccelerationModelValidAfter(AccelerationModel, AbsoluteDate)}
64   *   or {@link #addAccelerationModelValidBefore(AccelerationModel, AbsoluteDate)}.<br>
65   *   Recommendations are the same than the ones in {@link TimeSpanMap}, meaning: <ul>
66   *   <li>As an entry is added, it truncates the validity of the neighboring entries already present in the map;
67   *   <li><b>The transition dates should be entered only once</b>. Repeating a transition date will lead to unexpected result and is not supported;
68   *   <li>It is advised to order your AccelerationModel models chronologically when adding them to avoid any confusion.
69   *   </ul>
70   *   <li><u>Naming the parameter drivers</u>: It is strongly advised to give a custom name to the {@link ParameterDriver}(s)
71   *   of each AccelerationModel model that is added to the object. This will allow you keeping track of the evolution of your models.<br>
72   *   Different names are mandatory to differentiate the different drivers.<br>
73   *   Since there is no default name for acceleration model parameters, you must handle the driver names to consider
74   *   different names when adding a new acceleration model.
75   *   </ul>
76   * @author Bryan Cazabonne
77   * @since 10.3
78   */
79  public class TimeSpanParametricAcceleration extends AbstractForceModel {
80  
81      /** Prefix for dates before in the parameter drivers' name. */
82      public static final String DATE_BEFORE = " - Before ";
83  
84      /** Prefix for dates after in the parameter drivers' name. */
85      public static final String DATE_AFTER = " - After ";
86  
87      /** Direction of the acceleration in defining frame. */
88      private final Vector3D direction;
89  
90      /** Flag for inertial acceleration direction. */
91      private final boolean isInertial;
92  
93      /** The attitude to override, if set. */
94      private final AttitudeProvider attitudeOverride;
95  
96      /** TimeSpanMap of AccelerationModel objects. */
97      private final TimeSpanMap<AccelerationModel> accelerationModelTimeSpanMap;
98  
99      /** Simple constructor.
100      * @param direction acceleration direction in overridden spacecraft frame
101      * @param isInertial if true, direction is defined in the same inertial
102      * frame used for propagation (i.e. {@link SpacecraftState#getFrame()}),
103      * otherwise direction is defined in spacecraft frame (i.e. using the
104      * propagation {@link
105      * org.orekit.propagation.Propagator#setAttitudeProvider(AttitudeProvider)
106      * attitude law})
107      * @param accelerationModel acceleration model used to compute the contribution of the empirical acceleration
108      */
109     public TimeSpanParametricAcceleration(final Vector3D direction,
110                                           final boolean isInertial,
111                                           final AccelerationModel accelerationModel) {
112         this(direction, isInertial, null, accelerationModel);
113     }
114 
115     /** Simple constructor.
116      * @param direction acceleration direction in overridden spacecraft frame
117      * frame used for propagation (i.e. {@link SpacecraftState#getFrame()}),
118      * otherwise direction is defined in spacecraft frame (i.e. using the
119      * propagation {@link
120      * org.orekit.propagation.Propagator#setAttitudeProvider(AttitudeProvider)
121      * attitude law})
122      * @param attitudeOverride provider for attitude used to compute acceleration
123      * @param accelerationModel acceleration model used to compute the contribution of the empirical acceleration
124      */
125     public TimeSpanParametricAcceleration(final Vector3D direction,
126                                           final AttitudeProvider attitudeOverride,
127                                           final AccelerationModel accelerationModel) {
128         this(direction, false, attitudeOverride, accelerationModel);
129     }
130 
131     /** Simple constructor.
132      * @param direction acceleration direction in overridden spacecraft frame
133      * @param isInertial if true, direction is defined in the same inertial
134      * frame used for propagation (i.e. {@link SpacecraftState#getFrame()}),
135      * otherwise direction is defined in spacecraft frame (i.e. using the
136      * propagation {@link
137      * org.orekit.propagation.Propagator#setAttitudeProvider(AttitudeProvider)
138      * attitude law})
139      * @param attitudeOverride provider for attitude used to compute acceleration
140      * @param accelerationModel acceleration model used to compute the contribution of the empirical acceleration
141      */
142     private TimeSpanParametricAcceleration(final Vector3D direction,
143                                            final boolean isInertial,
144                                            final AttitudeProvider attitudeOverride,
145                                            final AccelerationModel accelerationModel) {
146         this.direction                    = direction;
147         this.isInertial                   = isInertial;
148         this.attitudeOverride             = attitudeOverride;
149         this.accelerationModelTimeSpanMap = new TimeSpanMap<>(accelerationModel);
150     }
151 
152     /** {@inheritDoc} */
153     @Override
154     public void init(final SpacecraftState initialState, final AbsoluteDate target) {
155         accelerationModelTimeSpanMap.forEach(accelerationModel -> accelerationModel.init(initialState, target));
156     }
157 
158     /** Add an AccelerationModel entry valid before a limit date.<br>
159      * <p>
160      * Using <code>addAccelerationModelValidBefore(entry, t)</code> will make <code>entry</code>
161      * valid in ]-∞, t[ (note the open bracket).
162      * <p>
163      * <b>WARNING</b>: Since there is no default name for acceleration model parameters,
164      * the user must handle itself the driver names to consider different names
165      * (i.e. different parameters) when adding a new acceleration model.
166      * @param accelerationModel AccelerationModel entry
167      * @param latestValidityDate date before which the entry is valid
168      * (must be different from <b>all</b> dates already used for transitions)
169      */
170     public void addAccelerationModelValidBefore(final AccelerationModel accelerationModel, final AbsoluteDate latestValidityDate) {
171         accelerationModelTimeSpanMap.addValidBefore(accelerationModel, latestValidityDate, false);
172     }
173 
174     /** Add a AccelerationModel entry valid after a limit date.<br>
175      * <p>
176      * Using <code>addAccelerationModelValidAfter(entry, t)</code> will make <code>entry</code>
177      * valid in [t, +∞[ (note the closed bracket).
178      * <p>
179      * <b>WARNING</b>: Since there is no default name for acceleration model parameters,
180      * the user must handle itself the driver names to consider different names
181      * (i.e. different parameters) when adding a new acceleration model.
182      * @param accelerationModel AccelerationModel entry
183      * @param earliestValidityDate date after which the entry is valid
184      * (must be different from <b>all</b> dates already used for transitions)
185      */
186     public void addAccelerationModelValidAfter(final AccelerationModel accelerationModel, final AbsoluteDate earliestValidityDate) {
187         accelerationModelTimeSpanMap.addValidAfter(accelerationModel, earliestValidityDate, false);
188     }
189 
190     /** Get the {@link AccelerationModel} model valid at a date.
191      * @param date the date of validity
192      * @return the AccelerationModel model valid at date
193      */
194     public AccelerationModel getAccelerationModel(final AbsoluteDate date) {
195         return accelerationModelTimeSpanMap.get(date);
196     }
197 
198     /** Get the {@link AccelerationModel} {@link Span} containing a specified date.
199      * @param date date belonging to the desired time span
200      * @return the AccelerationModel time span containing the specified date
201      */
202     public Span<AccelerationModel> getAccelerationModelSpan(final AbsoluteDate date) {
203         return accelerationModelTimeSpanMap.getSpan(date);
204     }
205 
206     /** Extract a range of the {@link AccelerationModel} map.
207      * <p>
208      * The object returned will be a new independent instance that will contain
209      * only the transitions that lie in the specified range.
210      * </p>
211      * See the {@link TimeSpanMap#extractRange TimeSpanMap.extractRange method} for more.
212      * @param start earliest date at which a transition is included in the range
213      * (may be set to {@link AbsoluteDate#PAST_INFINITY} to keep all early transitions)
214      * @param end latest date at which a transition is included in the r
215      * (may be set to {@link AbsoluteDate#FUTURE_INFINITY} to keep all late transitions)
216      * @return a new TimeSpanMap instance of AccelerationModel with all transitions restricted to the specified range
217      */
218     public TimeSpanMap<AccelerationModel> extractAccelerationModelRange(final AbsoluteDate start, final AbsoluteDate end) {
219         return accelerationModelTimeSpanMap.extractRange(start, end);
220     }
221 
222     /** Get the {@link Transition}s of the acceleration model time span map.
223      * @return the {@link Transition}s for the acceleration model time span map
224      * @deprecated as of 11.1, replace by {@link #getFirstSpan()}
225      */
226     @Deprecated
227     public NavigableSet<Transition<AccelerationModel>> getTransitions() {
228         return accelerationModelTimeSpanMap.getTransitions();
229     }
230 
231     /** Get the first {@link Span time span} of the acceleration model time span map.
232      * @return the first {@link Span time span} of the acceleration model time span map
233      * @since 11.1
234      */
235     public Span<AccelerationModel> getFirstSpan() {
236         return accelerationModelTimeSpanMap.getFirstSpan();
237     }
238 
239     /** {@inheritDoc} */
240     @Override
241     public boolean dependsOnPositionOnly() {
242         return isInertial;
243     }
244 
245     /** {@inheritDoc} */
246     @Override
247     public Vector3D acceleration(final SpacecraftState state,
248                                  final double[] parameters) {
249 
250         // Date
251         final AbsoluteDate date = state.getDate();
252 
253         // Compute inertial direction
254         final Vector3D inertialDirection;
255         if (isInertial) {
256             // the acceleration direction is already defined in the inertial frame
257             inertialDirection = direction;
258         } else {
259             final Attitude attitude;
260             if (attitudeOverride == null) {
261                 // the acceleration direction is defined in spacecraft frame as set by the propagator
262                 attitude = state.getAttitude();
263             } else {
264                 // the acceleration direction is defined in a dedicated frame
265                 attitude = attitudeOverride.getAttitude(state.getOrbit(), date, state.getFrame());
266             }
267             inertialDirection = attitude.getRotation().applyInverseTo(direction);
268         }
269 
270         // Extract the proper parameters valid at date from the input array
271         final double[] extractedParameters = extractParameters(parameters, date);
272 
273         // Compute and return the parametric acceleration
274         return new Vector3D(getAccelerationModel(date).signedAmplitude(state, extractedParameters), inertialDirection);
275 
276     }
277 
278     /** {@inheritDoc} */
279     @Override
280     public <T extends CalculusFieldElement<T>> FieldVector3D<T> acceleration(final FieldSpacecraftState<T> state,
281                                                                          final T[] parameters) {
282 
283         // Date
284         final FieldAbsoluteDate<T> date = state.getDate();
285 
286         // Compute inertial direction
287         final FieldVector3D<T> inertialDirection;
288         if (isInertial) {
289             // the acceleration direction is already defined in the inertial frame
290             inertialDirection = new FieldVector3D<>(state.getDate().getField(), direction);
291         } else {
292             final FieldAttitude<T> attitude;
293             if (attitudeOverride == null) {
294                 // the acceleration direction is defined in spacecraft frame as set by the propagator
295                 attitude = state.getAttitude();
296             } else {
297                 // the acceleration direction is defined in a dedicated frame
298                 attitude = attitudeOverride.getAttitude(state.getOrbit(), date, state.getFrame());
299             }
300             inertialDirection = attitude.getRotation().applyInverseTo(direction);
301         }
302 
303         // Extract the proper parameters valid at date from the input array
304         final T[] extractedParameters = extractParameters(parameters, date);
305 
306         // Compute and return the parametric acceleration
307         return new FieldVector3D<>(getAccelerationModel(date.toAbsoluteDate()).signedAmplitude(state, extractedParameters), inertialDirection);
308 
309     }
310 
311     /** {@inheritDoc} */
312     @Override
313     public Stream<EventDetector> getEventsDetectors() {
314         return Stream.empty();
315     }
316 
317     /** {@inheritDoc} */
318     @Override
319     public <T extends CalculusFieldElement<T>> Stream<FieldEventDetector<T>> getFieldEventsDetectors(final Field<T> field) {
320         return Stream.empty();
321     }
322 
323     /** {@inheritDoc}
324      * <p>
325      * All the parameter drivers of all AccelerationModel models are returned in an array.
326      * Models are ordered chronologically.
327      * </p>
328      */
329     @Override
330     public List<ParameterDriver> getParametersDrivers() {
331 
332         // Get all transitions from the TimeSpanMap
333         final List<ParameterDriver> listParameterDrivers = new ArrayList<>();
334 
335         // Loop on the spans
336         for (Span<AccelerationModel> span = getFirstSpan(); span != null; span = span.next()) {
337             // Add all the parameter drivers of the time span
338             for (ParameterDriver driver : span.getData().getParametersDrivers()) {
339                 // Add the driver only if the name does not exist already
340                 if (!findByName(listParameterDrivers, driver.getName())) {
341                     listParameterDrivers.add(driver);
342                 }
343             }
344         }
345 
346         // Return an array of parameter drivers with no duplicated name
347         return Collections.unmodifiableList(listParameterDrivers);
348 
349     }
350 
351     /** Extract the proper parameter drivers' values from the array in input of the
352      * {@link #acceleration(SpacecraftState, double[]) acceleration} method.
353      *  Parameters are filtered given an input date.
354      * @param parameters the input parameters array
355      * @param date the date
356      * @return the parameters given the date
357      */
358     public double[] extractParameters(final double[] parameters, final AbsoluteDate date) {
359 
360         // Get the acceleration model parameter drivers of the date
361         final List<ParameterDriver> empiricalParameterDriver = getAccelerationModel(date).getParametersDrivers();
362 
363         // Find out the indexes of the parameters in the whole array of parameters
364         final List<ParameterDriver> allParameters = getParametersDrivers();
365         final double[] outParameters = new double[empiricalParameterDriver.size()];
366         int index = 0;
367         for (int i = 0; i < allParameters.size(); i++) {
368             final String driverName = allParameters.get(i).getName();
369             for (ParameterDriver accDriver : empiricalParameterDriver) {
370                 if (accDriver.getName().equals(driverName)) {
371                     outParameters[index++] = parameters[i];
372                 }
373             }
374         }
375         return outParameters;
376     }
377 
378     /** Extract the proper parameter drivers' values from the array in input of the
379      * {@link #acceleration(FieldSpacecraftState, CalculusFieldElement[]) acceleration} method.
380      *  Parameters are filtered given an input date.
381      * @param parameters the input parameters array
382      * @param date the date
383      * @param <T> extends CalculusFieldElement
384      * @return the parameters given the date
385      */
386     public <T extends CalculusFieldElement<T>> T[] extractParameters(final T[] parameters,
387                                                                  final FieldAbsoluteDate<T> date) {
388 
389         // Get the acceleration parameter drivers of the date
390         final List<ParameterDriver> empiricalParameterDriver = getAccelerationModel(date.toAbsoluteDate()).getParametersDrivers();
391 
392         // Find out the indexes of the parameters in the whole array of parameters
393         final List<ParameterDriver> allParameters = getParametersDrivers();
394         final T[] outParameters = MathArrays.buildArray(date.getField(), empiricalParameterDriver.size());
395         int index = 0;
396         for (int i = 0; i < allParameters.size(); i++) {
397             final String driverName = allParameters.get(i).getName();
398             for (ParameterDriver accDriver : empiricalParameterDriver) {
399                 if (accDriver.getName().equals(driverName)) {
400                     outParameters[index++] = parameters[i];
401                 }
402             }
403         }
404         return outParameters;
405     }
406 
407     /** Find if a parameter driver with a given name already exists in a list of parameter drivers.
408      * @param driversList the list of parameter drivers
409      * @param name the parameter driver's name to filter with
410      * @return true if the name was found, false otherwise
411      */
412     private boolean findByName(final List<ParameterDriver> driversList, final String name) {
413         for (final ParameterDriver d : driversList) {
414             if (d.getName().equals(name)) {
415                 return true;
416             }
417         }
418         return false;
419     }
420 
421 }