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.utils;
18  
19  import java.io.Serializable;
20  import java.util.stream.Stream;
21  
22  import org.hipparchus.analysis.differentiation.Derivative;
23  import org.hipparchus.analysis.interpolation.HermiteInterpolator;
24  import org.hipparchus.geometry.euclidean.threed.FieldVector3D;
25  import org.hipparchus.geometry.euclidean.threed.Vector3D;
26  import org.hipparchus.util.FastMath;
27  import org.orekit.annotation.DefaultDataContext;
28  import org.orekit.data.DataContext;
29  import org.orekit.errors.OrekitException;
30  import org.orekit.errors.OrekitIllegalArgumentException;
31  import org.orekit.errors.OrekitInternalError;
32  import org.orekit.errors.OrekitMessages;
33  import org.orekit.frames.Frame;
34  import org.orekit.frames.Transform;
35  import org.orekit.time.AbsoluteDate;
36  import org.orekit.time.TimeInterpolable;
37  import org.orekit.time.TimeStamped;
38  
39  /** Position - Velocity - Acceleration linked to a date and a frame.
40   */
41  public class AbsolutePVCoordinates extends TimeStampedPVCoordinates
42      implements TimeStamped, TimeInterpolable<AbsolutePVCoordinates>,
43                 Serializable, PVCoordinatesProvider {
44  
45      /** Serializable UID. */
46      private static final long serialVersionUID = 20150824L;
47  
48      /** Frame in which are defined the coordinates. */
49      private final Frame frame;
50  
51      /** Build from position, velocity, acceleration.
52       * @param frame the frame in which the coordinates are defined
53       * @param date coordinates date
54       * @param position the position vector (m)
55       * @param velocity the velocity vector (m/s)
56       * @param acceleration the acceleration vector (m/sÂý)
57       */
58      public AbsolutePVCoordinates(final Frame frame, final AbsoluteDate date,
59                                   final Vector3D position, final Vector3D velocity, final Vector3D acceleration) {
60          super(date, position, velocity, acceleration);
61          this.frame = frame;
62      }
63  
64      /** Build from position and velocity. Acceleration is set to zero.
65       * @param frame the frame in which the coordinates are defined
66       * @param date coordinates date
67       * @param position the position vector (m)
68       * @param velocity the velocity vector (m/s)
69       */
70      public AbsolutePVCoordinates(final Frame frame, final AbsoluteDate date,
71                                   final Vector3D position,
72                                   final Vector3D velocity) {
73          this(frame, date, position, velocity, Vector3D.ZERO);
74      }
75  
76      /** Build from frame, date and PVA coordinates.
77       * @param frame the frame in which the coordinates are defined
78       * @param date date of the coordinates
79       * @param pva TimeStampedPVCoordinates
80       */
81      public AbsolutePVCoordinates(final Frame frame, final AbsoluteDate date, final PVCoordinates pva) {
82          super(date, pva);
83          this.frame = frame;
84      }
85  
86      /** Build from frame and TimeStampedPVCoordinates.
87       * @param frame the frame in which the coordinates are defined
88       * @param pva TimeStampedPVCoordinates
89       */
90      public AbsolutePVCoordinates(final Frame frame, final TimeStampedPVCoordinates pva) {
91          super(pva.getDate(), pva);
92          this.frame = frame;
93      }
94  
95      /** Multiplicative constructor
96       * <p>Build a AbsolutePVCoordinates from another one and a scale factor.</p>
97       * <p>The TimeStampedPVCoordinates built will be a * AbsPva</p>
98       * @param date date of the built coordinates
99       * @param a scale factor
100      * @param AbsPva base (unscaled) AbsolutePVCoordinates
101      */
102     public AbsolutePVCoordinates(final AbsoluteDate date,
103                                  final double a, final AbsolutePVCoordinates AbsPva) {
104         super(date, a, AbsPva);
105         this.frame = AbsPva.frame;
106     }
107 
108     /** Subtractive constructor
109      * <p>Build a relative AbsolutePVCoordinates from a start and an end position.</p>
110      * <p>The AbsolutePVCoordinates built will be end - start.</p>
111      * <p>In case start and end use two different pseudo-inertial frames,
112      * the new AbsolutePVCoordinates arbitrarily be defined in the start frame. </p>
113      * @param date date of the built coordinates
114      * @param start Starting AbsolutePVCoordinates
115      * @param end ending AbsolutePVCoordinates
116      */
117     public AbsolutePVCoordinates(final AbsoluteDate date,
118                                  final AbsolutePVCoordinates start, final AbsolutePVCoordinates end) {
119         super(date, start, end);
120         ensureIdenticalFrames(start, end);
121         this.frame = start.frame;
122     }
123 
124     /** Linear constructor
125      * <p>Build a AbsolutePVCoordinates from two other ones and corresponding scale factors.</p>
126      * <p>The AbsolutePVCoordinates built will be a1 * u1 + a2 * u2</p>
127      * <p>In case the AbsolutePVCoordinates use different pseudo-inertial frames,
128      * the new AbsolutePVCoordinates arbitrarily be defined in the first frame. </p>
129      * @param date date of the built coordinates
130      * @param a1 first scale factor
131      * @param absPv1 first base (unscaled) AbsolutePVCoordinates
132      * @param a2 second scale factor
133      * @param absPv2 second base (unscaled) AbsolutePVCoordinates
134      */
135     public AbsolutePVCoordinates(final AbsoluteDate date,
136                                  final double a1, final AbsolutePVCoordinates absPv1,
137                                  final double a2, final AbsolutePVCoordinates absPv2) {
138         super(date, a1, absPv1.getPVCoordinates(), a2, absPv2.getPVCoordinates());
139         ensureIdenticalFrames(absPv1, absPv2);
140         this.frame = absPv1.getFrame();
141     }
142 
143     /** Linear constructor
144      * <p>Build a AbsolutePVCoordinates from three other ones and corresponding scale factors.</p>
145      * <p>The AbsolutePVCoordinates built will be a1 * u1 + a2 * u2 + a3 * u3</p>
146      * <p>In case the AbsolutePVCoordinates use different pseudo-inertial frames,
147      * the new AbsolutePVCoordinates arbitrarily be defined in the first frame. </p>
148      * @param date date of the built coordinates
149      * @param a1 first scale factor
150      * @param absPv1 first base (unscaled) AbsolutePVCoordinates
151      * @param a2 second scale factor
152      * @param absPv2 second base (unscaled) AbsolutePVCoordinates
153      * @param a3 third scale factor
154      * @param absPv3 third base (unscaled) AbsolutePVCoordinates
155      */
156     public AbsolutePVCoordinates(final AbsoluteDate date,
157                                  final double a1, final AbsolutePVCoordinates absPv1,
158                                  final double a2, final AbsolutePVCoordinates absPv2,
159                                  final double a3, final AbsolutePVCoordinates absPv3) {
160         super(date, a1, absPv1.getPVCoordinates(), a2, absPv2.getPVCoordinates(),
161                 a3, absPv3.getPVCoordinates());
162         ensureIdenticalFrames(absPv1, absPv2);
163         ensureIdenticalFrames(absPv1, absPv3);
164         this.frame = absPv1.getFrame();
165     }
166 
167     /** Linear constructor
168      * <p>Build a AbsolutePVCoordinates from four other ones and corresponding scale factors.</p>
169      * <p>The AbsolutePVCoordinates built will be a1 * u1 + a2 * u2 + a3 * u3 + a4 * u4</p>
170      * <p>In case the AbsolutePVCoordinates use different pseudo-inertial frames,
171      * the new AbsolutePVCoordinates arbitrarily be defined in the first frame. </p>
172      * @param date date of the built coordinates
173      * @param a1 first scale factor
174      * @param absPv1 first base (unscaled) AbsolutePVCoordinates
175      * @param a2 second scale factor
176      * @param absPv2 second base (unscaled) AbsolutePVCoordinates
177      * @param a3 third scale factor
178      * @param absPv3 third base (unscaled) AbsolutePVCoordinates
179      * @param a4 fourth scale factor
180      * @param absPv4 fourth base (unscaled) AbsolutePVCoordinates
181      */
182     public AbsolutePVCoordinates(final AbsoluteDate date,
183                                  final double a1, final AbsolutePVCoordinates absPv1,
184                                  final double a2, final AbsolutePVCoordinates absPv2,
185                                  final double a3, final AbsolutePVCoordinates absPv3,
186                                  final double a4, final AbsolutePVCoordinates absPv4) {
187         super(date, a1, absPv1.getPVCoordinates(), a2, absPv2.getPVCoordinates(),
188                 a3, absPv3.getPVCoordinates(), a4, absPv4.getPVCoordinates());
189         ensureIdenticalFrames(absPv1, absPv2);
190         ensureIdenticalFrames(absPv1, absPv3);
191         ensureIdenticalFrames(absPv1, absPv4);
192         this.frame = absPv1.getFrame();
193     }
194 
195     /** Builds a AbsolutePVCoordinates triplet from  a {@link FieldVector3D}&lt;{@link Derivative}&gt;.
196      * <p>
197      * The vector components must have time as their only derivation parameter and
198      * have consistent derivation orders.
199      * </p>
200      * @param frame the frame in which the parameters are defined
201      * @param date date of the built coordinates
202      * @param p vector with time-derivatives embedded within the coordinates
203      * @param <U> type of the derivative
204      */
205     public <U extends Derivative<U>> AbsolutePVCoordinates(final Frame frame, final AbsoluteDate date,
206                                                            final FieldVector3D<U> p) {
207         super(date, p);
208         this.frame = frame;
209     }
210 
211     /** Ensure that the frames from two AbsolutePVCoordinates are identical.
212      * @param absPv1 first AbsolutePVCoordinates
213      * @param absPv2 first AbsolutePVCoordinates
214      * @throws OrekitIllegalArgumentException if frames are different
215      */
216     private static void ensureIdenticalFrames(final AbsolutePVCoordinates absPv1, final AbsolutePVCoordinates absPv2)
217         throws OrekitIllegalArgumentException {
218         if (!absPv1.frame.equals(absPv2.frame)) {
219             throw new OrekitIllegalArgumentException(OrekitMessages.INCOMPATIBLE_FRAMES,
220                                                      absPv1.frame.getName(), absPv2.frame.getName());
221         }
222     }
223 
224     /** Get a time-shifted state.
225      * <p>
226      * The state can be slightly shifted to close dates. This shift is based on
227      * a simple Taylor expansion. It is <em>not</em> intended as a replacement for
228      * proper orbit propagation (it is not even Keplerian!) but should be sufficient
229      * for either small time shifts or coarse accuracy.
230      * </p>
231      * @param dt time shift in seconds
232      * @return a new state, shifted with respect to the instance (which is immutable)
233      */
234     public AbsolutePVCoordinates shiftedBy(final double dt) {
235         final TimeStampedPVCoordinates spv = super.shiftedBy(dt);
236         return new AbsolutePVCoordinates(frame, spv);
237     }
238 
239     /** Create a local provider using simply Taylor expansion through {@link #shiftedBy(double)}.
240      * <p>
241      * The time evolution is based on a simple Taylor expansion. It is <em>not</em> intended as a
242      * replacement for proper orbit propagation (it is not even Keplerian!) but should be sufficient
243      * for either small time shifts or coarse accuracy.
244      * </p>
245      * @return provider based on Taylor expansion, for small time shifts around instance date
246      */
247     public PVCoordinatesProvider toTaylorProvider() {
248         return new PVCoordinatesProvider() {
249             /** {@inheritDoc} */
250             public TimeStampedPVCoordinates getPVCoordinates(final AbsoluteDate d,  final Frame f) {
251                 final TimeStampedPVCoordinates shifted   = shiftedBy(d.durationFrom(getDate()));
252                 final Transform                transform = frame.getTransformTo(f, d);
253                 return transform.transformPVCoordinates(shifted);
254             }
255         };
256     }
257 
258     /** Get the frame in which the coordinates are defined.
259      * @return frame in which the coordinates are defined
260      */
261     public Frame getFrame() {
262         return frame;
263     }
264 
265     /** Get the TimeStampedPVCoordinates.
266      * @return TimeStampedPVCoordinates
267      */
268     public TimeStampedPVCoordinates getPVCoordinates() {
269         return this;
270     }
271 
272     /** Get the TimeStampedPVCoordinates in a specified frame.
273      * @param outputFrame frame in which the position/velocity coordinates shall be computed
274      * @return TimeStampedPVCoordinates
275      * @exception OrekitException if transformation between frames cannot be computed
276      * @see #getPVCoordinates()
277      */
278     public TimeStampedPVCoordinates getPVCoordinates(final Frame outputFrame) {
279         // If output frame requested is the same as definition frame,
280         // PV coordinates are returned directly
281         if (outputFrame == frame) {
282             return getPVCoordinates();
283         }
284 
285         // Else, PV coordinates are transformed to output frame
286         final Transform t = frame.getTransformTo(outputFrame, getDate());
287         return t.transformPVCoordinates(getPVCoordinates());
288     }
289 
290     @Override
291     public TimeStampedPVCoordinates getPVCoordinates(final AbsoluteDate otherDate, final Frame outputFrame) {
292         return shiftedBy(otherDate.durationFrom(getDate())).getPVCoordinates(outputFrame);
293     }
294 
295     @Override
296     public AbsolutePVCoordinates interpolate(final AbsoluteDate date, final Stream<AbsolutePVCoordinates> sample) {
297         return interpolate(getFrame(), date, CartesianDerivativesFilter.USE_PVA, sample);
298     }
299 
300     /** Interpolate position-velocity.
301      * <p>
302      * The interpolated instance is created by polynomial Hermite interpolation
303      * ensuring velocity remains the exact derivative of position.
304      * </p>
305      * <p>
306      * Note that even if first time derivatives (velocities)
307      * from sample can be ignored, the interpolated instance always includes
308      * interpolated derivatives. This feature can be used explicitly to
309      * compute these derivatives when it would be too complex to compute them
310      * from an analytical formula: just compute a few sample points from the
311      * explicit formula and set the derivatives to zero in these sample points,
312      * then use interpolation to add derivatives consistent with the positions.
313      * </p>
314      * @param frame frame for the interpolted instance
315      * @param date interpolation date
316      * @param filter filter for derivatives from the sample to use in interpolation
317      * @param sample sample points on which interpolation should be done
318      * @return a new position-velocity, interpolated at specified date
319      * @exception OrekitIllegalArgumentException if some elements in the sample do not
320      * have the same defining frame as other
321      */
322     public static AbsolutePVCoordinates interpolate(final Frame frame, final AbsoluteDate date,
323                                                     final CartesianDerivativesFilter filter,
324                                                     final Stream<AbsolutePVCoordinates> sample) {
325 
326 
327         // set up an interpolator taking derivatives into account
328         final HermiteInterpolator interpolator = new HermiteInterpolator();
329 
330         // add sample points
331         switch (filter) {
332             case USE_P :
333                 // populate sample with position data, ignoring velocity
334                 sample.forEach(pv -> {
335                     final Vector3D position = pv.getPosition();
336                     interpolator.addSamplePoint(pv.getDate().durationFrom(date),
337                                                 position.toArray());
338                 });
339                 break;
340             case USE_PV :
341                 // populate sample with position and velocity data
342                 sample.forEach(pv -> {
343                     final Vector3D position = pv.getPosition();
344                     final Vector3D velocity = pv.getVelocity();
345                     interpolator.addSamplePoint(pv.getDate().durationFrom(date),
346                                                 position.toArray(), velocity.toArray());
347                 });
348                 break;
349             case USE_PVA :
350                 // populate sample with position, velocity and acceleration data
351                 sample.forEach(pv -> {
352                     final Vector3D position     = pv.getPosition();
353                     final Vector3D velocity     = pv.getVelocity();
354                     final Vector3D acceleration = pv.getAcceleration();
355                     interpolator.addSamplePoint(pv.getDate().durationFrom(date),
356                                                 position.toArray(), velocity.toArray(), acceleration.toArray());
357                 });
358                 break;
359             default :
360                 // this should never happen
361                 throw new OrekitInternalError(null);
362         }
363 
364         // interpolate
365         final double[][] p = interpolator.derivatives(0.0, 2);
366 
367         // build a new interpolated instance
368         return new AbsolutePVCoordinates(frame, date, new Vector3D(p[0]), new Vector3D(p[1]), new Vector3D(p[2]));
369 
370     }
371 
372     /** Replace the instance with a data transfer object for serialization.
373      * @return data transfer object that will be serialized
374      */
375     @DefaultDataContext
376     private Object writeReplace() {
377         return new DTO(this);
378     }
379 
380     /** Internal class used only for serialization. */
381     @DefaultDataContext
382     private static class DTO implements Serializable {
383 
384         /** Serializable UID. */
385         private static final long serialVersionUID = 20150916L;
386 
387         /** Double values. */
388         private double[] d;
389 
390         /** Frame in which acoordinates are defined. */
391         private final Frame frame;
392 
393         /** Simple constructor.
394          * @param absPva instance to serialize
395          */
396         private DTO(final AbsolutePVCoordinates absPva) {
397 
398             // decompose date
399             final AbsoluteDate j2000Epoch =
400                     DataContext.getDefault().getTimeScales().getJ2000Epoch();
401             final double epoch  = FastMath.floor(absPva.getDate().durationFrom(j2000Epoch));
402             final double offset = absPva.getDate().durationFrom(j2000Epoch.shiftedBy(epoch));
403 
404             this.d = new double[] {
405                 epoch, offset,
406                 absPva.getPosition().getX(),     absPva.getPosition().getY(),     absPva.getPosition().getZ(),
407                 absPva.getVelocity().getX(),     absPva.getVelocity().getY(),     absPva.getVelocity().getZ(),
408                 absPva.getAcceleration().getX(), absPva.getAcceleration().getY(), absPva.getAcceleration().getZ()
409             };
410             this.frame = absPva.frame;
411 
412         }
413 
414         /** Replace the deserialized data transfer object with a {@link AbsolutePVCoordinates}.
415          * @return replacement {@link AbsolutePVCoordinates}
416          */
417         private Object readResolve() {
418             final AbsoluteDate j2000Epoch =
419                     DataContext.getDefault().getTimeScales().getJ2000Epoch();
420             return new AbsolutePVCoordinates(frame,
421                                              j2000Epoch.shiftedBy(d[0]).shiftedBy(d[1]),
422                                              new Vector3D(d[2], d[3], d[ 4]),
423                                              new Vector3D(d[5], d[6], d[ 7]),
424                                              new Vector3D(d[8], d[9], d[10]));
425         }
426 
427     }
428 
429 }
430 
431 
432