1   /* Copyright 2013-2016 CS Systèmes d'Information
2    * Licensed to CS Systèmes d'Information (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.rugged.utils;
18  
19  import java.io.Serializable;
20  import java.util.ArrayList;
21  import java.util.List;
22  
23  import org.apache.commons.math3.util.FastMath;
24  import org.orekit.errors.OrekitException;
25  import org.orekit.frames.Frame;
26  import org.orekit.frames.Transform;
27  import org.orekit.rugged.errors.DumpManager;
28  import org.orekit.rugged.errors.RuggedException;
29  import org.orekit.rugged.errors.RuggedMessages;
30  import org.orekit.time.AbsoluteDate;
31  import org.orekit.utils.AngularDerivativesFilter;
32  import org.orekit.utils.CartesianDerivativesFilter;
33  import org.orekit.utils.ImmutableTimeStampedCache;
34  import org.orekit.utils.TimeStampedAngularCoordinates;
35  import org.orekit.utils.TimeStampedCache;
36  import org.orekit.utils.TimeStampedPVCoordinates;
37  
38  /** Provider for observation transforms.
39   * @author Luc Maisonobe
40   */
41  public class SpacecraftToObservedBody implements Serializable {
42  
43      /** Serializable UID. */
44      private static final long serialVersionUID = 20140909L;
45  
46      /** Inertial frame. */
47      private final Frame inertialFrame;
48  
49      /** Body frame. */
50      private final Frame bodyFrame;
51  
52      /** Start of search time span. */
53      private final AbsoluteDate minDate;
54  
55      /** End of search time span. */
56      private final AbsoluteDate maxDate;
57  
58      /** Step to use for inertial frame to body frame transforms cache computations. */
59      private final double tStep;
60  
61      /** Tolerance in seconds allowed for {@code minDate} and {@code maxDate} overshooting. */
62      private final double overshootTolerance;
63  
64      /** Transforms sample from observed body frame to inertial frame. */
65      private final List<Transform> bodyToInertial;
66  
67      /** Transforms sample from inertial frame to observed body frame. */
68      private final List<Transform> inertialToBody;
69  
70      /** Transforms sample from spacecraft frame to inertial frame. */
71      private final List<Transform> scToInertial;
72  
73      /** Simple constructor.
74       * @param inertialFrame inertial frame
75       * @param bodyFrame observed body frame
76       * @param minDate start of search time span
77       * @param maxDate end of search time span
78       * @param tStep step to use for inertial frame to body frame transforms cache computations
79       * @param overshootTolerance tolerance in seconds allowed for {@code minDate} and {@code maxDate} overshooting
80       * slightly the position, velocity and quaternions ephemerides
81       * @param positionsVelocities satellite position and velocity
82       * @param pvInterpolationNumber number of points to use for position/velocity interpolation
83       * @param pvFilter filter for derivatives from the sample to use in position/velocity interpolation
84       * @param quaternions satellite quaternions
85       * @param aInterpolationNumber number of points to use for attitude interpolation
86       * @param aFilter filter for derivatives from the sample to use in attitude interpolation
87       * @exception RuggedException if [{@code minDate}, {@code maxDate}] search time span overshoots
88       * position or attitude samples by more than {@code overshootTolerance}
89       * ,
90       */
91      public SpacecraftToObservedBody(final Frame inertialFrame, final Frame bodyFrame,
92                                      final AbsoluteDate minDate, final AbsoluteDate maxDate, final double tStep,
93                                      final double overshootTolerance,
94                                      final List<TimeStampedPVCoordinates> positionsVelocities, final int pvInterpolationNumber,
95                                      final CartesianDerivativesFilter pvFilter,
96                                      final List<TimeStampedAngularCoordinates> quaternions, final int aInterpolationNumber,
97                                      final AngularDerivativesFilter aFilter)
98          throws RuggedException {
99          try {
100 
101             this.inertialFrame      = inertialFrame;
102             this.bodyFrame          = bodyFrame;
103             this.minDate            = minDate;
104             this.maxDate            = maxDate;
105             this.overshootTolerance = overshootTolerance;
106 
107             // safety checks
108             final AbsoluteDate minPVDate = positionsVelocities.get(0).getDate();
109             final AbsoluteDate maxPVDate = positionsVelocities.get(positionsVelocities.size() - 1).getDate();
110             if (minPVDate.durationFrom(minDate) > overshootTolerance) {
111                 throw new RuggedException(RuggedMessages.OUT_OF_TIME_RANGE, minDate, minPVDate, maxPVDate);
112             }
113             if (maxDate.durationFrom(maxDate) > overshootTolerance) {
114                 throw new RuggedException(RuggedMessages.OUT_OF_TIME_RANGE, maxDate, minPVDate, maxPVDate);
115             }
116 
117             final AbsoluteDate minQDate  = quaternions.get(0).getDate();
118             final AbsoluteDate maxQDate  = quaternions.get(quaternions.size() - 1).getDate();
119             if (minQDate.durationFrom(minDate) > overshootTolerance) {
120                 throw new RuggedException(RuggedMessages.OUT_OF_TIME_RANGE, minDate, minQDate, maxQDate);
121             }
122             if (maxDate.durationFrom(maxQDate) > overshootTolerance) {
123                 throw new RuggedException(RuggedMessages.OUT_OF_TIME_RANGE, maxDate, minQDate, maxQDate);
124             }
125 
126             // set up the cache for position-velocities
127             final TimeStampedCache<TimeStampedPVCoordinates> pvCache =
128                     new ImmutableTimeStampedCache<TimeStampedPVCoordinates>(pvInterpolationNumber, positionsVelocities);
129 
130             // set up the cache for attitudes
131             final TimeStampedCache<TimeStampedAngularCoordinates> aCache =
132                     new ImmutableTimeStampedCache<TimeStampedAngularCoordinates>(aInterpolationNumber, quaternions);
133 
134             final int n = (int) FastMath.ceil(maxDate.durationFrom(minDate) / tStep);
135             this.tStep          = tStep;
136             this.bodyToInertial = new ArrayList<Transform>(n);
137             this.inertialToBody = new ArrayList<Transform>(n);
138             this.scToInertial   = new ArrayList<Transform>(n);
139             for (AbsoluteDate date = minDate; bodyToInertial.size() < n; date = date.shiftedBy(tStep)) {
140 
141                 // interpolate position-velocity, allowing slight extrapolation near the boundaries
142                 final AbsoluteDate pvInterpolationDate;
143                 if (date.compareTo(pvCache.getEarliest().getDate()) < 0) {
144                     pvInterpolationDate = pvCache.getEarliest().getDate();
145                 } else if (date.compareTo(pvCache.getLatest().getDate()) > 0) {
146                     pvInterpolationDate = pvCache.getLatest().getDate();
147                 } else {
148                     pvInterpolationDate = date;
149                 }
150                 final TimeStampedPVCoordinates interpolatedPV =
151                         TimeStampedPVCoordinates.interpolate(pvInterpolationDate, pvFilter,
152                                                              pvCache.getNeighbors(pvInterpolationDate));
153                 final TimeStampedPVCoordinates pv = interpolatedPV.shiftedBy(date.durationFrom(pvInterpolationDate));
154 
155                 // interpolate attitude, allowing slight extrapolation near the boundaries
156                 final AbsoluteDate aInterpolationDate;
157                 if (date.compareTo(aCache.getEarliest().getDate()) < 0) {
158                     aInterpolationDate = aCache.getEarliest().getDate();
159                 } else if (date.compareTo(aCache.getLatest().getDate()) > 0) {
160                     aInterpolationDate = aCache.getLatest().getDate();
161                 } else {
162                     aInterpolationDate = date;
163                 }
164                 final TimeStampedAngularCoordinates interpolatedQuaternion =
165                         TimeStampedAngularCoordinates.interpolate(aInterpolationDate, aFilter,
166                                                                   aCache.getNeighbors(aInterpolationDate));
167                 final TimeStampedAngularCoordinates quaternion = interpolatedQuaternion.shiftedBy(date.durationFrom(aInterpolationDate));
168 
169                 // store transform from spacecraft frame to inertial frame
170                 scToInertial.add(new Transform(date,
171                                                new Transform(date, quaternion.revert()),
172                                                new Transform(date, pv)));
173 
174                 // store transform from body frame to inertial frame
175                 final Transform b2i = bodyFrame.getTransformTo(inertialFrame, date);
176                 bodyToInertial.add(b2i);
177                 inertialToBody.add(b2i.getInverse());
178 
179             }
180 
181         } catch (OrekitException oe) {
182             throw new RuggedException(oe, oe.getSpecifier(), oe.getParts());
183         }
184     }
185 
186     /** Simple constructor.
187      * @param inertialFrame inertial frame
188      * @param bodyFrame observed body frame
189      * @param minDate start of search time span
190      * @param maxDate end of search time span
191      * @param tStep step to use for inertial frame to body frame transforms cache computations
192      * @param overshootTolerance tolerance in seconds allowed for {@code minDate} and {@code maxDate} overshooting
193      * slightly the position, velocity and quaternions ephemerides
194      * @param bodyToInertial transforms sample from observed body frame to inertial frame
195      * @param scToInertial transforms sample from spacecraft frame to inertial frame
196      */
197     public SpacecraftToObservedBody(final Frame inertialFrame, final Frame bodyFrame,
198                                     final AbsoluteDate minDate, final AbsoluteDate maxDate, final double tStep,
199                                     final double overshootTolerance,
200                                     final List<Transform> bodyToInertial, final List<Transform> scToInertial) {
201 
202         this.inertialFrame      = inertialFrame;
203         this.bodyFrame          = bodyFrame;
204         this.minDate            = minDate;
205         this.maxDate            = maxDate;
206         this.tStep              = tStep;
207         this.overshootTolerance = overshootTolerance;
208         this.bodyToInertial     = bodyToInertial;
209         this.scToInertial       = scToInertial;
210 
211         this.inertialToBody = new ArrayList<Transform>(bodyToInertial.size());
212         for (final Transform b2i : bodyToInertial) {
213             inertialToBody.add(b2i.getInverse());
214         }
215 
216     }
217 
218     /** Get the inertial frame.
219      * @return inertial frame
220      */
221     public Frame getInertialFrame() {
222         return inertialFrame;
223     }
224 
225     /** Get the body frame.
226      * @return body frame
227      */
228     public Frame getBodyFrame() {
229         return bodyFrame;
230     }
231 
232     /** Get the start of search time span.
233      * @return start of search time span
234      */
235     public AbsoluteDate getMinDate() {
236         return minDate;
237     }
238 
239     /** Get the end of search time span.
240      * @return end of search time span
241      */
242     public AbsoluteDate getMaxDate() {
243         return maxDate;
244     }
245 
246     /** Get the step to use for inertial frame to body frame transforms cache computations.
247      * @return step to use for inertial frame to body frame transforms cache computations
248      */
249     public double getTStep() {
250         return tStep;
251     }
252 
253     /** Get the tolerance in seconds allowed for {@link #getMinDate()} and {@link #getMaxDate()} overshooting.
254      * @return tolerance in seconds allowed for {@link #getMinDate()} and {@link #getMaxDate()} overshooting
255      */
256     public double getOvershootTolerance() {
257         return overshootTolerance;
258     }
259 
260     /** Get transform from spacecraft to inertial frame.
261      * @param date date of the transform
262      * @return transform from spacecraft to inertial frame
263      * @exception RuggedException if spacecraft position or attitude cannot be computed at date
264      */
265     public Transform getScToInertial(final AbsoluteDate date)
266         throws RuggedException {
267         return interpolate(date, scToInertial);
268     }
269 
270     /** Get transform from inertial frame to observed body frame.
271      * @param date date of the transform
272      * @return transform from inertial frame to observed body frame
273      * @exception RuggedException if frames cannot be computed at date
274      */
275     public Transform getInertialToBody(final AbsoluteDate date)
276         throws RuggedException {
277         return interpolate(date, inertialToBody);
278     }
279 
280     /** Get transform from observed body frame to inertial frame.
281      * @param date date of the transform
282      * @return transform from observed body frame to inertial frame
283      * @exception RuggedException if frames cannot be computed at date
284      */
285     public Transform getBodyToInertial(final AbsoluteDate date)
286         throws RuggedException {
287         return interpolate(date, bodyToInertial);
288     }
289 
290     /** Interpolate transform.
291      * @param date date of the transform
292      * @param list transforms list to interpolate from
293      * @return interpolated transform
294      * @exception RuggedException if frames cannot be computed at date
295      */
296     private Transform interpolate(final AbsoluteDate date, final List<Transform> list)
297         throws RuggedException {
298 
299         // check date range
300         if (!isInRange(date)) {
301             throw new RuggedException(RuggedMessages.OUT_OF_TIME_RANGE, date, minDate, maxDate);
302         }
303 
304         final double    s     = date.durationFrom(list.get(0).getDate()) / tStep;
305         final int       index = FastMath.max(0, FastMath.min(list.size() - 1, (int) FastMath.rint(s)));
306 
307         DumpManager.dumpTransform(this, index, bodyToInertial.get(index), scToInertial.get(index));
308 
309         final Transform close = list.get(index);
310         return close.shiftedBy(date.durationFrom(close.getDate()));
311 
312     }
313 
314     /** Check if a date is in the supported range.
315      * @param date date to check
316      * @return true if date is in the supported range
317      */
318     public boolean isInRange(final AbsoluteDate date) {
319         return (minDate.durationFrom(date) <= overshootTolerance) &&
320                (date.durationFrom(maxDate) <= overshootTolerance);
321     }
322 
323 }