1   /* Copyright 2013-2019 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.api;
18  
19  import java.io.IOException;
20  import java.io.InputStream;
21  import java.io.ObjectInputStream;
22  import java.io.ObjectOutputStream;
23  import java.io.OutputStream;
24  import java.util.ArrayList;
25  import java.util.Collections;
26  import java.util.List;
27  
28  import org.hipparchus.exception.LocalizedCoreFormats;
29  import org.hipparchus.geometry.euclidean.threed.Rotation;
30  import org.hipparchus.geometry.euclidean.threed.Vector3D;
31  import org.orekit.bodies.OneAxisEllipsoid;
32  import org.orekit.frames.Frame;
33  import org.orekit.frames.FramesFactory;
34  import org.orekit.propagation.Propagator;
35  import org.orekit.propagation.SpacecraftState;
36  import org.orekit.propagation.sampling.OrekitFixedStepHandler;
37  import org.orekit.rugged.errors.RuggedException;
38  import org.orekit.rugged.errors.RuggedInternalError;
39  import org.orekit.rugged.errors.RuggedMessages;
40  import org.orekit.rugged.intersection.BasicScanAlgorithm;
41  import org.orekit.rugged.intersection.ConstantElevationAlgorithm;
42  import org.orekit.rugged.intersection.IgnoreDEMAlgorithm;
43  import org.orekit.rugged.intersection.IntersectionAlgorithm;
44  import org.orekit.rugged.intersection.duvenhage.DuvenhageAlgorithm;
45  import org.orekit.rugged.linesensor.LineSensor;
46  import org.orekit.rugged.raster.TileUpdater;
47  import org.orekit.rugged.refraction.AtmosphericRefraction;
48  import org.orekit.rugged.utils.ExtendedEllipsoid;
49  import org.orekit.rugged.utils.SpacecraftToObservedBody;
50  import org.orekit.time.AbsoluteDate;
51  import org.orekit.utils.AngularDerivativesFilter;
52  import org.orekit.utils.CartesianDerivativesFilter;
53  import org.orekit.utils.Constants;
54  import org.orekit.utils.IERSConventions;
55  import org.orekit.utils.PVCoordinates;
56  import org.orekit.utils.TimeStampedAngularCoordinates;
57  import org.orekit.utils.TimeStampedPVCoordinates;
58  
59  /** Builder for {@link Rugged} instances.
60   * <p>
61   * This class implements the <em>builder pattern</em> to create {@link Rugged} instances.
62   * It does so by using a <em>fluent API</em> in order to clarify reading and allow
63   * later extensions with new configuration parameters.
64   * </p>
65   * <p>
66   * A typical use would be:
67   * </p>
68   * <pre>
69   *   Rugged rugged = new RuggedBuilder().
70   *                   setEllipsoid(EllipsoidId.WGS84, BodyRotatingFrameId.ITRF).
71   *                   setAlgorithmID(AlgorithmId.Duvenhage).
72   *                   setDigitalElevationModel(tileUpdater, maxCachedTiles).
73   *                   setTimeSpan(minDate, maxDate, tStep, overshootTolerance).
74   *                   setTrajectory(positionsVelocities, pvInterpolationNumber, pvFilter,
75   *                                 quaternions, aInterpolationNumber, aFilter).
76   *                   addLineSensor(sensor1).
77   *                   addLineSensor(sensor2).
78   *                   addLineSensor(sensor3).
79   *                   build();
80   * </pre>
81   * <p>
82   * If a configuration parameter has not been set prior to the call to {]link #build()}, then
83   * an exception will be triggered with an explicit error message.
84   * </p>
85   * @see <a href="https://en.wikipedia.org/wiki/Builder_pattern">Builder pattern (wikipedia)</a>
86   * @see <a href="https://en.wikipedia.org/wiki/Fluent_interface">Fluent interface (wikipedia)</a>
87   * @author Luc Maisonobe
88   * @author Guylaine Prat
89   */
90  public class RuggedBuilder {
91  
92      /** Reference ellipsoid. */
93      private ExtendedEllipsoid ellipsoid;
94  
95      /** DEM intersection algorithm. */
96      private AlgorithmId algorithmID;
97  
98      /** Updater used to load Digital Elevation Model tiles. */
99      private TileUpdater tileUpdater;
100 
101     /** Constant elevation over ellipsoid (m).
102      * used only with {@link AlgorithmId#CONSTANT_ELEVATION_OVER_ELLIPSOID. */
103     private double constantElevation;
104 
105     /** Maximum number of tiles stored in the cache. */
106     private int maxCachedTiles;
107 
108     /** Start of search time span. */
109     private AbsoluteDate minDate;
110 
111     /** End of search time span. */
112     private AbsoluteDate maxDate;
113 
114     /** Step to use for inertial frame to body frame transforms cache computations (s). */
115     private double tStep;
116 
117     /** OvershootTolerance tolerance in seconds allowed for {@code minDate} and {@code maxDate} overshooting (s). */
118     private double overshootTolerance;
119 
120     /** Inertial frame for position/velocity/attitude. */
121     private Frame inertial;
122 
123     /** Satellite position and velocity (m and m/s in inertial frame). */
124     private List<TimeStampedPVCoordinates> pvSample;
125 
126     /** Number of points to use for position/velocity interpolation. */
127     private int pvNeighborsSize;
128 
129     /** Filter for derivatives from the sample to use in position/velocity interpolation. */
130     private CartesianDerivativesFilter pvDerivatives;
131 
132     /** Satellite quaternions with respect to inertial frame. */
133     private List<TimeStampedAngularCoordinates> aSample;
134 
135     /** Number of points to use for attitude interpolation. */
136     private int aNeighborsSize;
137 
138     /** Filter for derivatives from the sample to use in attitude interpolation. */
139     private AngularDerivativesFilter aDerivatives;
140 
141     /** Propagator for position/velocity/attitude. */
142     private Propagator pvaPropagator;
143 
144     /** Step to use for inertial/Earth/spacecraft transforms interpolations (s). */
145     private double iStep;
146 
147     /** Number of points to use for inertial/Earth/spacecraft transforms interpolations. */
148     private int iN;
149 
150     /** Converter between spacecraft and body. */
151     private SpacecraftToObservedBody scToBody;
152 
153     /** Flag for light time correction. */
154     private boolean lightTimeCorrection;
155 
156     /** Flag for aberration of light correction. */
157     private boolean aberrationOfLightCorrection;
158 
159     /** Atmospheric refraction to use for line of sight correction. */
160     private AtmosphericRefraction atmosphericRefraction;
161 
162     /** Sensors. */
163     private final List<LineSensor> sensors;
164 
165     /** Rugged name. */
166     private String name;
167 
168     /** Create a non-configured builder.
169      * <p>
170      * The builder <em>must</em> be configured before calling the
171      * {@link #build} method, otherwise an exception will be triggered
172      * at build time.
173      * </p>
174      */
175     public RuggedBuilder() {
176         sensors                     = new ArrayList<LineSensor>();
177         constantElevation           = Double.NaN;
178         lightTimeCorrection         = true;
179         aberrationOfLightCorrection = true;
180         name                        = "Rugged";
181     }
182 
183     /** Set the reference ellipsoid.
184      * @param ellipsoidID reference ellipsoid
185      * @param bodyRotatingFrameID body rotating frame identifier
186      * from an earlier run and frames mismatch
187      * @return the builder instance
188      * @see #setEllipsoid(OneAxisEllipsoid)
189      * @see #getEllipsoid()
190      */
191     public RuggedBuilder setEllipsoid(final EllipsoidId ellipsoidID, final BodyRotatingFrameId bodyRotatingFrameID) {
192         return setEllipsoid(selectEllipsoid(ellipsoidID, selectBodyRotatingFrame(bodyRotatingFrameID)));
193     }
194 
195     /** Set the reference ellipsoid.
196      * @param newEllipsoid reference ellipsoid
197      * @return the builder instance
198      * @see #setEllipsoid(EllipsoidId, BodyRotatingFrameId)
199      * @see #getEllipsoid()
200      */
201     public RuggedBuilder setEllipsoid(final OneAxisEllipsoid newEllipsoid) {
202         this.ellipsoid = new ExtendedEllipsoid(newEllipsoid.getEquatorialRadius(),
203                                                newEllipsoid.getFlattening(),
204                                                newEllipsoid.getBodyFrame());
205         checkFramesConsistency();
206         return this;
207     }
208 
209     /** Get the ellipsoid.
210      * @return the ellipsoid
211      * @see #setEllipsoid(EllipsoidId, BodyRotatingFrameId)
212      * @see #setEllipsoid(OneAxisEllipsoid)
213      */
214     public ExtendedEllipsoid getEllipsoid() {
215         return ellipsoid;
216     }
217 
218     /** Get the Rugged name.
219      * @return the Rugged name
220      * @since 2.0
221      */
222     public String getName() {
223         return name;
224     }
225 
226     /** Set the Rugged name.
227      * @param name the Rugged name
228      * @since 2.0
229      */
230     public void setName(final String name) {
231         this.name = name;
232     }
233 
234     /** Set the algorithm to use for Digital Elevation Model intersection.
235      * <p>
236      * Note that some algorithms require specific other methods to be called too:
237      * <ul>
238      *   <li>{@link AlgorithmId#DUVENHAGE DUVENHAGE},
239      *   {@link AlgorithmId#DUVENHAGE_FLAT_BODY DUVENHAGE_FLAT_BODY}
240      *   and {@link AlgorithmId#BASIC_SLOW_EXHAUSTIVE_SCAN_FOR_TESTS_ONLY
241      *   BASIC_SLOW_EXHAUSTIVE_SCAN_FOR_TESTS_ONLY} all
242      *   require {@link #setDigitalElevationModel(TileUpdater, int) setDigitalElevationModel}
243      *   to be called,</li>
244      *   <li>{@link AlgorithmId#CONSTANT_ELEVATION_OVER_ELLIPSOID
245      *   CONSTANT_ELEVATION_OVER_ELLIPSOID} requires
246      *   {@link #setConstantElevation(double) setConstantElevation} to be called,</li>
247      *   <li>{@link AlgorithmId#IGNORE_DEM_USE_ELLIPSOID
248      *   IGNORE_DEM_USE_ELLIPSOID} does not require
249      *   any methods tobe called.</li>
250      * </ul>
251      *
252      * @param newAlgorithmId identifier of algorithm to use for Digital Elevation Model intersection
253      * @return the builder instance
254      * @see #setDigitalElevationModel(TileUpdater, int)
255      * @see #getAlgorithm()
256      */
257     public RuggedBuilder setAlgorithm(final AlgorithmId newAlgorithmId) {
258         this.algorithmID = newAlgorithmId;
259         return this;
260     }
261 
262     /** Get the algorithm to use for Digital Elevation Model intersection.
263      * @return algorithm to use for Digital Elevation Model intersection
264      * @see #setAlgorithm(AlgorithmId)
265      */
266     public AlgorithmId getAlgorithm() {
267         return algorithmID;
268     }
269 
270     /** Set the user-provided {@link TileUpdater tile updater}.
271      * <p>
272      * Note that when the algorithm specified in {@link #setAlgorithm(AlgorithmId)}
273      * is either {@link AlgorithmId#CONSTANT_ELEVATION_OVER_ELLIPSOID
274      * CONSTANT_ELEVATION_OVER_ELLIPSOID} or {@link
275      * AlgorithmId#IGNORE_DEM_USE_ELLIPSOID IGNORE_DEM_USE_ELLIPSOID},
276      * then this method becomes irrelevant and can either be not called at all,
277      * or it can be called with an updater set to {@code null}. For all other
278      * algorithms, the updater must be properly configured.
279      * </p>
280      * @param newTileUpdater updater used to load Digital Elevation Model tiles
281      * @param newMaxCachedTiles maximum number of tiles stored in the cache
282      * @return the builder instance
283      * @see #setAlgorithm(AlgorithmId)
284      * @see #getTileUpdater()
285      * @see #getMaxCachedTiles()
286      */
287     public RuggedBuilder setDigitalElevationModel(final TileUpdater newTileUpdater, final int newMaxCachedTiles) {
288         this.tileUpdater    = newTileUpdater;
289         this.maxCachedTiles = newMaxCachedTiles;
290         return this;
291     }
292 
293     /** Get the updater used to load Digital Elevation Model tiles.
294      * @return updater used to load Digital Elevation Model tiles
295      * @see #setDigitalElevationModel(TileUpdater, int)
296      * @see #getMaxCachedTiles()
297      */
298     public TileUpdater getTileUpdater() {
299         return tileUpdater;
300     }
301 
302     /** Set the user-provided constant elevation model.
303      * <p>
304      * Note that this method is relevant <em>only</em> if the algorithm specified
305      * in {@link #setAlgorithm(AlgorithmId)} is {@link
306      * AlgorithmId#CONSTANT_ELEVATION_OVER_ELLIPSOID CONSTANT_ELEVATION_OVER_ELLIPSOID}.
307      * If it is called for another algorithm, the elevation set here will be ignored.
308      * </p>
309      * @param newConstantElevation constant elevation to use (m)
310      * @return the builder instance
311      * @see #setAlgorithm(AlgorithmId)
312      * @see #getConstantElevation()
313      */
314     public RuggedBuilder setConstantElevation(final double newConstantElevation) {
315         this.constantElevation = newConstantElevation;
316         return this;
317     }
318 
319     /** Get the constant elevation over ellipsoid to use with {@link AlgorithmId#CONSTANT_ELEVATION_OVER_ELLIPSOID}.
320      * @return updater used to load Digital Elevation Model tiles
321      * @see #setConstantElevation(double)
322      */
323     public double getConstantElevation() {
324         return constantElevation;
325     }
326 
327     /** Get the maximum number of tiles stored in the cache.
328      * @return maximum number of tiles stored in the cache
329      * @see #setDigitalElevationModel(TileUpdater, int)
330      * @see #getTileUpdater()
331      */
332     public int getMaxCachedTiles() {
333         return maxCachedTiles;
334     }
335 
336     /** Set the time span to be covered for direct and inverse location calls.
337      * <p>
338      * This method set only the time span and not the trajectory, therefore it
339      * <em>must</em> be used together with either
340      * {@link #setTrajectory(InertialFrameId, List, int, CartesianDerivativesFilter, List, int, AngularDerivativesFilter)},
341      * {@link #setTrajectory(Frame, List, int, CartesianDerivativesFilter, List, int, AngularDerivativesFilter)},
342      * or {@link #setTrajectory(double, int, CartesianDerivativesFilter, AngularDerivativesFilter, Propagator)}
343      * but should <em>not</em> be mixed with {@link #setTrajectoryAndTimeSpan(InputStream)}.
344      * </p>
345      * @param newMinDate start of search time span
346      * @param newMaxDate end of search time span
347      * @param newTstep step to use for inertial frame to body frame transforms cache computations (s)
348      * @param newOvershootTolerance tolerance in seconds allowed for {@code minDate} and {@code maxDate} overshooting (s)
349      * @return the builder instance
350      * @see #setTrajectoryAndTimeSpan(InputStream)
351      * @see #getMinDate()
352      * @see #getMaxDate()
353      * @see #getTStep()
354      * @see #getOvershootTolerance()
355      */
356     public RuggedBuilder setTimeSpan(final AbsoluteDate newMinDate, final AbsoluteDate newMaxDate,
357                                      final double newTstep, final double newOvershootTolerance) {
358         this.minDate            = newMinDate;
359         this.maxDate            = newMaxDate;
360         this.tStep              = newTstep;
361         this.overshootTolerance = newOvershootTolerance;
362         this.scToBody           = null;
363         return this;
364     }
365 
366     /** Get the start of search time span.
367      * @return start of search time span
368      * @see #setTimeSpan(AbsoluteDate, AbsoluteDate, double, double)
369      */
370     public AbsoluteDate getMinDate() {
371         return minDate;
372     }
373 
374     /** Get the end of search time span.
375      * @return end of search time span
376      * @see #setTimeSpan(AbsoluteDate, AbsoluteDate, double, double)
377      */
378     public AbsoluteDate getMaxDate() {
379         return maxDate;
380     }
381 
382     /** Get the step to use for inertial frame to body frame transforms cache computations.
383      * @return step to use for inertial frame to body frame transforms cache computations
384      * @see #setTimeSpan(AbsoluteDate, AbsoluteDate, double, double)
385      */
386     public double getTStep() {
387         return tStep;
388     }
389 
390     /** Get the tolerance in seconds allowed for {@link #getMinDate()} and {@link #getMaxDate()} overshooting.
391      * @return tolerance in seconds allowed for {@link #getMinDate()} and {@link #getMaxDate()} overshooting
392      * @see #setTimeSpan(AbsoluteDate, AbsoluteDate, double, double)
393      */
394     public double getOvershootTolerance() {
395         return overshootTolerance;
396     }
397 
398     /** Set the spacecraft trajectory.
399      * <p>
400      * This method set only the trajectory and not the time span, therefore it
401      * <em>must</em> be used together with the {@link #setTimeSpan(AbsoluteDate, AbsoluteDate, double, double)}
402      * but should <em>not</em> be mixed with {@link #setTrajectoryAndTimeSpan(InputStream)}.
403      * </p>
404      * @param inertialFrameId inertial frame Id used for spacecraft positions/velocities/quaternions
405      * @param positionsVelocities satellite position and velocity (m and m/s in inertial frame)
406      * @param pvInterpolationNumber number of points to use for position/velocity interpolation
407      * @param pvFilter filter for derivatives from the sample to use in position/velocity interpolation
408      * @param quaternions satellite quaternions with respect to inertial frame
409      * @param aInterpolationNumber number of points to use for attitude interpolation
410      * @param aFilter filter for derivatives from the sample to use in attitude interpolation
411      * @return the builder instance
412      * @see #setTrajectory(Frame, List, int, CartesianDerivativesFilter, List, int, AngularDerivativesFilter)
413      * @see #setTrajectory(double, int, CartesianDerivativesFilter, AngularDerivativesFilter, Propagator)
414      * @see #setTrajectoryAndTimeSpan(InputStream)
415      * @see #getInertialFrame()
416      * @see #getPositionsVelocities()
417      * @see #getPVInterpolationNumber()
418      * @see #getPVInterpolationNumber()
419      * @see #getQuaternions()
420      * @see #getAInterpolationNumber()
421      * @see #getAFilter()
422      */
423     public RuggedBuilder setTrajectory(final InertialFrameId inertialFrameId,
424                                        final List<TimeStampedPVCoordinates> positionsVelocities, final int pvInterpolationNumber,
425                                        final CartesianDerivativesFilter pvFilter,
426                                        final List<TimeStampedAngularCoordinates> quaternions, final int aInterpolationNumber,
427                                        final AngularDerivativesFilter aFilter) {
428 
429         return setTrajectory(selectInertialFrame(inertialFrameId),
430                              positionsVelocities, pvInterpolationNumber, pvFilter,
431                              quaternions, aInterpolationNumber, aFilter);
432     }
433 
434     /** Set the spacecraft trajectory.
435      * <p>
436      * This method set only the trajectory and not the time span, therefore it
437      * <em>must</em> be used together with the {@link #setTimeSpan(AbsoluteDate, AbsoluteDate, double, double)}
438      * but should <em>not</em> be mixed with {@link #setTrajectoryAndTimeSpan(InputStream)}.
439      * </p>
440      * @param inertialFrame inertial frame used for spacecraft positions/velocities/quaternions
441      * @param positionsVelocities satellite position and velocity (m and m/s in inertial frame)
442      * @param pvInterpolationNumber number of points to use for position/velocity interpolation
443      * @param pvFilter filter for derivatives from the sample to use in position/velocity interpolation
444      * @param quaternions satellite quaternions with respect to inertial frame
445      * @param aInterpolationNumber number of points to use for attitude interpolation
446      * @param aFilter filter for derivatives from the sample to use in attitude interpolation
447      * @return the builder instance
448      * @see #setTrajectory(InertialFrameId, List, int, CartesianDerivativesFilter, List, int, AngularDerivativesFilter)
449      * @see #setTrajectory(double, int, CartesianDerivativesFilter, AngularDerivativesFilter, Propagator)
450      * @see #setTrajectoryAndTimeSpan(InputStream)
451      * @see #getPositionsVelocities()
452      * @see #getPVInterpolationNumber()
453      * @see #getPVInterpolationNumber()
454      * @see #getQuaternions()
455      * @see #getAInterpolationNumber()
456      * @see #getAFilter()
457      */
458     public RuggedBuilder setTrajectory(final Frame inertialFrame,
459                                        final List<TimeStampedPVCoordinates> positionsVelocities, final int pvInterpolationNumber,
460                                        final CartesianDerivativesFilter pvFilter,
461                                        final List<TimeStampedAngularCoordinates> quaternions, final int aInterpolationNumber,
462                                        final AngularDerivativesFilter aFilter) {
463 
464         this.inertial        = inertialFrame;
465         this.pvSample        = positionsVelocities;
466         this.pvNeighborsSize = pvInterpolationNumber;
467         this.pvDerivatives   = pvFilter;
468         this.aSample         = quaternions;
469         this.aNeighborsSize  = aInterpolationNumber;
470         this.aDerivatives    = aFilter;
471         this.pvaPropagator   = null;
472         this.iStep           = Double.NaN;
473         this.iN              = -1;
474         this.scToBody        = null;
475         return this;
476     }
477 
478     /** Set the spacecraft trajectory.
479      * <p>
480      * This method set only the trajectory and not the time span, therefore it
481      * <em>must</em> be used together with the {@link #setTimeSpan(AbsoluteDate, AbsoluteDate, double, double)}
482      * but should <em>not</em> be mixed with {@link #setTrajectoryAndTimeSpan(InputStream)}.
483      * </p>
484      * @param interpolationStep step to use for inertial/Earth/spacecraft transforms interpolations (s)
485      * @param interpolationNumber number of points to use for inertial/Earth/spacecraft transforms interpolations
486      * @param pvFilter filter for derivatives from the sample to use in position/velocity interpolation
487      * @param aFilter filter for derivatives from the sample to use in attitude interpolation
488      * @param propagator global propagator
489      * @return the builder instance
490      * @see #setTrajectory(InertialFrameId, List, int, CartesianDerivativesFilter, List, int, AngularDerivativesFilter)
491      * @see #setTrajectory(Frame, List, int, CartesianDerivativesFilter, List, int, AngularDerivativesFilter)
492      * @see #setTrajectoryAndTimeSpan(InputStream)
493      */
494     public RuggedBuilder setTrajectory(final double interpolationStep, final int interpolationNumber,
495                                        final CartesianDerivativesFilter pvFilter, final AngularDerivativesFilter aFilter,
496                                        final Propagator propagator) {
497         this.inertial        = propagator.getFrame();
498         this.pvSample        = null;
499         this.pvNeighborsSize = -1;
500         this.pvDerivatives   = pvFilter;
501         this.aSample         = null;
502         this.aNeighborsSize  = -1;
503         this.aDerivatives    = aFilter;
504         this.pvaPropagator   = propagator;
505         this.iStep           = interpolationStep;
506         this.iN              = interpolationNumber;
507         this.scToBody        = null;
508         return this;
509     }
510 
511     /** Get the inertial frame.
512      * @return inertial frame
513      */
514     public Frame getInertialFrame() {
515         return inertial;
516     }
517 
518     /** Get the satellite position and velocity (m and m/s in inertial frame).
519      * @return satellite position and velocity (m and m/s in inertial frame)
520      * @see #setTrajectory(double, int, CartesianDerivativesFilter, AngularDerivativesFilter, Propagator)
521      */
522     public List<TimeStampedPVCoordinates> getPositionsVelocities() {
523         return pvSample;
524     }
525 
526     /** Get the number of points to use for position/velocity interpolation.
527      * @return number of points to use for position/velocity interpolation
528      * @see #setTrajectory(double, int, CartesianDerivativesFilter, AngularDerivativesFilter, Propagator)
529      */
530     public int getPVInterpolationNumber() {
531         return pvNeighborsSize;
532     }
533 
534     /** Get the filter for derivatives from the sample to use in position/velocity interpolation.
535      * @return filter for derivatives from the sample to use in position/velocity interpolation
536      * @see #setTrajectory(double, int, CartesianDerivativesFilter, AngularDerivativesFilter, Propagator)
537      */
538     public CartesianDerivativesFilter getPVFilter() {
539         return pvDerivatives;
540     }
541 
542     /** Get the satellite quaternions with respect to inertial frame.
543      * @return satellite quaternions with respect to inertial frame
544      * @see #setTrajectory(double, int, CartesianDerivativesFilter, AngularDerivativesFilter, Propagator)
545      */
546     public List<TimeStampedAngularCoordinates> getQuaternions() {
547         return aSample;
548     }
549 
550     /** Get the number of points to use for attitude interpolation.
551      * @return number of points to use for attitude interpolation
552      * @see #setTrajectory(double, int, CartesianDerivativesFilter, AngularDerivativesFilter, Propagator)
553      */
554     public int getAInterpolationNumber() {
555         return aNeighborsSize;
556     }
557 
558     /** Get the filter for derivatives from the sample to use in attitude interpolation.
559      * @return filter for derivatives from the sample to use in attitude interpolation
560      * @see #setTrajectory(double, int, CartesianDerivativesFilter, AngularDerivativesFilter, Propagator)
561      */
562     public AngularDerivativesFilter getAFilter() {
563         return aDerivatives;
564     }
565 
566     /** Set both the spacecraft trajectory and the time span.
567      * <p>
568      * This method set both the trajectory and the time span in a tightly coupled
569      * way, therefore it should <em>not</em> be mixed with the individual methods
570      * {@link #setTrajectory(double, int, CartesianDerivativesFilter, AngularDerivativesFilter, Propagator)},
571      * {@link #setTrajectory(Frame, List, int, CartesianDerivativesFilter, List, int, AngularDerivativesFilter)},
572      * {@link #setTrajectory(InertialFrameId, List, int, CartesianDerivativesFilter, List, int, AngularDerivativesFilter)},
573      * or {@link #setTimeSpan(AbsoluteDate, AbsoluteDate, double, double)}.
574      * </p>
575      * @param storageStream stream from where to read previous instance {@link #storeInterpolator(OutputStream)
576      * stored interpolator} (caller opened it and remains responsible for closing it)
577      * @return the builder instance
578      * @see #setTrajectory(InertialFrameId, List, int, CartesianDerivativesFilter, List, int, AngularDerivativesFilter)
579      * @see #setTrajectory(Frame, List, int, CartesianDerivativesFilter, List, int, AngularDerivativesFilter)
580      * @see #setTrajectory(double, int, CartesianDerivativesFilter, AngularDerivativesFilter, Propagator)
581      * @see #storeInterpolator(OutputStream)
582      */
583     public RuggedBuilder setTrajectoryAndTimeSpan(final InputStream storageStream) {
584 
585         try {
586             this.inertial           = null;
587             this.pvSample           = null;
588             this.pvNeighborsSize    = -1;
589             this.pvDerivatives      = null;
590             this.aSample            = null;
591             this.aNeighborsSize     = -1;
592             this.aDerivatives       = null;
593             this.pvaPropagator      = null;
594             this.iStep              = Double.NaN;
595             this.iN                 = -1;
596             this.scToBody           = (SpacecraftToObservedBody) new ObjectInputStream(storageStream).readObject();
597             this.minDate            = scToBody.getMinDate();
598             this.maxDate            = scToBody.getMaxDate();
599             this.tStep              = scToBody.getTStep();
600             this.overshootTolerance = scToBody.getOvershootTolerance();
601             checkFramesConsistency();
602             return this;
603 
604         } catch (ClassNotFoundException cnfe) {
605             throw new RuggedException(cnfe, RuggedMessages.NOT_INTERPOLATOR_DUMP_DATA);
606         } catch (ClassCastException cce) {
607             throw new RuggedException(cce, RuggedMessages.NOT_INTERPOLATOR_DUMP_DATA);
608         } catch (IOException ioe) {
609             throw new RuggedException(ioe, RuggedMessages.NOT_INTERPOLATOR_DUMP_DATA);
610         }
611     }
612 
613     /** Store frames transform interpolator.
614      * <p>
615      * This method allows to reuse the interpolator built in one instance, to build
616      * another instance by calling {@link #setTrajectoryAndTimeSpan(InputStream)}.
617      * This reduces the builder initialization time as setting up the interpolator can be long, it is
618      * mainly intended to be used when several runs are done (for example in an image processing chain)
619      * with the same configuration.
620      * </p>
621      * <p>
622      * This method must be called <em>after</em> both the ellipsoid and trajectory have been set.
623      * </p>
624      * @param storageStream stream where to store the interpolator
625      * (caller opened it and remains responsible for closing it)
626      * @see #setEllipsoid(EllipsoidId, BodyRotatingFrameId)
627      * @see #setEllipsoid(OneAxisEllipsoid)
628      * @see #setTrajectory(InertialFrameId, List, int, CartesianDerivativesFilter, List, int, AngularDerivativesFilter)
629      * @see #setTrajectory(Frame, List, int, CartesianDerivativesFilter, List, int, AngularDerivativesFilter)
630      * @see #setTrajectory(double, int, CartesianDerivativesFilter, AngularDerivativesFilter, Propagator)
631      * @see #setTrajectoryAndTimeSpan(InputStream)
632      */
633     public void storeInterpolator(final OutputStream storageStream) {
634         try {
635             createInterpolatorIfNeeded();
636             new ObjectOutputStream(storageStream).writeObject(scToBody);
637         } catch (IOException ioe) {
638             throw new RuggedException(ioe, LocalizedCoreFormats.SIMPLE_MESSAGE, ioe.getMessage());
639         }
640     }
641 
642     /** Check frames consistency.
643      */
644     private void checkFramesConsistency() {
645         if (ellipsoid != null && scToBody != null &&
646             !ellipsoid.getBodyFrame().getName().equals(scToBody.getBodyFrame().getName())) {
647             // if frames have been set both by direct calls and by deserializing an interpolator dump and a mismatch occurs
648             throw new RuggedException(RuggedMessages.FRAMES_MISMATCH_WITH_INTERPOLATOR_DUMP,
649                                       ellipsoid.getBodyFrame().getName(), scToBody.getBodyFrame().getName());
650         }
651     }
652 
653     /** Create a transform interpolator if needed.
654      */
655     private void createInterpolatorIfNeeded() {
656 
657         if (ellipsoid == null) {
658             throw new RuggedException(RuggedMessages.UNINITIALIZED_CONTEXT, "RuggedBuilder.setEllipsoid()");
659         }
660 
661         if (scToBody == null) {
662             if (pvSample != null) {
663                 scToBody = createInterpolator(inertial, ellipsoid.getBodyFrame(),
664                                               minDate, maxDate, tStep, overshootTolerance,
665                                               pvSample, pvNeighborsSize, pvDerivatives,
666                                               aSample, aNeighborsSize, aDerivatives);
667             } else if (pvaPropagator != null) {
668                 scToBody = createInterpolator(inertial, ellipsoid.getBodyFrame(),
669                                               minDate, maxDate, tStep, overshootTolerance,
670                                               iStep, iN, pvDerivatives, aDerivatives, pvaPropagator);
671             } else {
672                 throw new RuggedException(RuggedMessages.UNINITIALIZED_CONTEXT, "RuggedBuilder.setTrajectory()");
673             }
674         }
675     }
676 
677     /** Create a transform interpolator from positions and quaternions lists.
678      * @param inertialFrame inertial frame
679      * @param bodyFrame observed body frame
680      * @param minDate start of search time span
681      * @param maxDate end of search time span
682      * @param tStep step to use for inertial frame to body frame transforms cache computations
683      * @param overshootTolerance tolerance in seconds allowed for {@code minDate} and {@code maxDate} overshooting
684      * @param positionsVelocities satellite position and velocity
685      * @param pvInterpolationNumber number of points to use for position/velocity interpolation
686      * @param pvFilter filter for derivatives from the sample to use in position/velocity interpolation
687      * @param quaternions satellite quaternions
688      * @param aInterpolationNumber number of points to use for attitude interpolation
689      * @param aFilter filter for derivatives from the sample to use in attitude interpolation
690      * @return transforms interpolator
691      */
692     private static SpacecraftToObservedBody createInterpolator(final Frame inertialFrame, final Frame bodyFrame,
693                                                                final AbsoluteDate minDate, final AbsoluteDate maxDate,
694                                                                final double tStep, final double overshootTolerance,
695                                                                final List<TimeStampedPVCoordinates> positionsVelocities,
696                                                                final int pvInterpolationNumber,
697                                                                final CartesianDerivativesFilter pvFilter,
698                                                                final List<TimeStampedAngularCoordinates> quaternions,
699                                                                final int aInterpolationNumber,
700                                                                final AngularDerivativesFilter aFilter) {
701 
702         return new SpacecraftToObservedBody(inertialFrame, bodyFrame,
703                                             minDate, maxDate, tStep, overshootTolerance,
704                                             positionsVelocities, pvInterpolationNumber,
705                                             pvFilter, quaternions, aInterpolationNumber,
706                                             aFilter);
707     }
708 
709     /** Create a transform interpolator from a propagator.
710      * @param inertialFrame inertial frame
711      * @param bodyFrame observed body frame
712      * @param minDate start of search time span
713      * @param maxDate end of search time span
714      * @param tStep step to use for inertial frame to body frame transforms cache computations
715      * @param overshootTolerance tolerance in seconds allowed for {@code minDate} and {@code maxDate} overshooting
716      * @param interpolationStep step to use for inertial/Earth/spacecraft transforms interpolations
717      * @param interpolationNumber number of points of to use for inertial/Earth/spacecraft transforms interpolations
718      * @param pvFilter filter for derivatives from the sample to use in position/velocity interpolation
719      * @param aFilter filter for derivatives from the sample to use in attitude interpolation
720      * @param propagator global propagator
721      * @return transforms interpolator
722      */
723     private static SpacecraftToObservedBody createInterpolator(final Frame inertialFrame, final Frame bodyFrame,
724                                                                final AbsoluteDate minDate, final AbsoluteDate maxDate,
725                                                                final double tStep, final double overshootTolerance,
726                                                                final double interpolationStep, final int interpolationNumber,
727                                                                final CartesianDerivativesFilter pvFilter,
728                                                                final AngularDerivativesFilter aFilter,
729                                                                final Propagator propagator) {
730 
731         // extract position/attitude samples from propagator
732         final List<TimeStampedPVCoordinates> positionsVelocities =
733                 new ArrayList<TimeStampedPVCoordinates>();
734         final List<TimeStampedAngularCoordinates> quaternions =
735                 new ArrayList<TimeStampedAngularCoordinates>();
736         propagator.setMasterMode(interpolationStep, new OrekitFixedStepHandler() {
737 
738             /** {@inheritDoc} */
739             @Override
740             public void handleStep(final SpacecraftState currentState, final boolean isLast) {
741                 final AbsoluteDate  date = currentState.getDate();
742                 final PVCoordinates pv   = currentState.getPVCoordinates(inertialFrame);
743                 final Rotation      q    = currentState.getAttitude().getRotation();
744                 positionsVelocities.add(new TimeStampedPVCoordinates(date, pv.getPosition(), pv.getVelocity(), Vector3D.ZERO));
745                 quaternions.add(new TimeStampedAngularCoordinates(date, q, Vector3D.ZERO, Vector3D.ZERO));
746             }
747 
748         });
749         propagator.propagate(minDate.shiftedBy(-interpolationStep), maxDate.shiftedBy(interpolationStep));
750 
751         // orbit/attitude to body converter
752         return createInterpolator(inertialFrame, bodyFrame,
753                 minDate, maxDate, tStep, overshootTolerance,
754                 positionsVelocities, interpolationNumber,
755                 pvFilter, quaternions, interpolationNumber,
756                 aFilter);
757     }
758 
759     /** Set flag for light time correction.
760      * <p>
761      * This methods set the flag for compensating or not light time between
762      * ground and spacecraft. Compensating this delay improves location
763      * accuracy and is <em>enabled</em> by default (i.e. not calling this
764      * method before building is therefore equivalent to calling it with
765      * a parameter set to {@code true}). Not compensating it is mainly useful
766      * for validation purposes against system that do not compensate it.
767      * </p>
768      * @param newLightTimeCorrection if true, the light travel time between ground
769      * and spacecraft is compensated for more accurate location
770      * @return the builder instance
771      * @see #setAberrationOfLightCorrection(boolean)
772      * @see #getLightTimeCorrection()
773      */
774     public RuggedBuilder setLightTimeCorrection(final boolean newLightTimeCorrection) {
775         this.lightTimeCorrection = newLightTimeCorrection;
776         return this;
777     }
778 
779     /** Get the light time correction flag.
780      * @return light time correction flag
781      * @see #setLightTimeCorrection(boolean)
782      */
783     public boolean getLightTimeCorrection() {
784         return lightTimeCorrection;
785     }
786 
787     /** Set flag for aberration of light correction.
788      * <p>
789      * This methods set the flag for compensating or not aberration of light,
790      * which is velocity composition between light and spacecraft when the
791      * light from ground points reaches the sensor.
792      * Compensating this velocity composition improves location
793      * accuracy and is <em>enabled</em> by default (i.e. not calling this
794      * method before building is therefore equivalent to calling it with
795      * a parameter set to {@code true}). Not compensating it is useful
796      * in two cases: for validation purposes against system that do not
797      * compensate it or when the pixels line of sight already include the
798      * correction.
799      * </p>
800      * @param newAberrationOfLightCorrection if true, the aberration of light
801      * is corrected for more accurate location
802      * @return the builder instance
803      * @see #setLightTimeCorrection(boolean)
804      * @see #getAberrationOfLightCorrection()
805      */
806     public RuggedBuilder setAberrationOfLightCorrection(final boolean newAberrationOfLightCorrection) {
807         this.aberrationOfLightCorrection = newAberrationOfLightCorrection;
808         return this;
809     }
810 
811     /** Get the aberration of light correction flag.
812      * @return aberration of light correction flag
813      * @see #setAberrationOfLightCorrection(boolean)
814      */
815     public boolean getAberrationOfLightCorrection() {
816         return aberrationOfLightCorrection;
817     }
818 
819     /** Set atmospheric refraction for line of sight correction.
820      * <p>
821      * This method sets an atmospheric refraction model to be used between
822      * spacecraft and ground for the correction of intersected points on ground.
823      * Compensating for the effect of atmospheric refraction improves location
824      * accuracy.
825      * </p>
826      * @param newAtmosphericRefraction the atmospheric refraction model to be used for more accurate location
827      * @return the builder instance
828      * @see #getRefractionCorrection()
829      */
830     public RuggedBuilder setRefractionCorrection(final AtmosphericRefraction newAtmosphericRefraction) {
831         this.atmosphericRefraction = newAtmosphericRefraction;
832         return this;
833     }
834 
835     /** Get the atmospheric refraction model.
836      * @return atmospheric refraction model
837      * @see #setRefractionCorrection(AtmosphericRefraction)
838      */
839     public AtmosphericRefraction getRefractionCorrection() {
840         return atmosphericRefraction;
841     }
842 
843     /** Set up line sensor model.
844      * @param lineSensor line sensor model
845      * @return the builder instance
846      */
847     public RuggedBuilder addLineSensor(final LineSensor lineSensor) {
848         sensors.add(lineSensor);
849         return this;
850     }
851 
852     /** Remove all line sensors.
853      * @return the builder instance
854      */
855     public RuggedBuilder clearLineSensors() {
856         sensors.clear();
857         return this;
858     }
859 
860     /** Get all line sensors.
861      * @return all line sensors (in an unmodifiable list)
862      */
863     public List<LineSensor> getLineSensors() {
864         return Collections.unmodifiableList(sensors);
865     }
866 
867     /** Select inertial frame.
868      * @param inertialFrameId inertial frame identifier
869      * @return inertial frame
870      */
871     private static Frame selectInertialFrame(final InertialFrameId inertialFrameId) {
872 
873         // set up the inertial frame
874         switch (inertialFrameId) {
875             case GCRF :
876                 return FramesFactory.getGCRF();
877             case EME2000 :
878                 return FramesFactory.getEME2000();
879             case MOD :
880                 return FramesFactory.getMOD(IERSConventions.IERS_1996);
881             case TOD :
882                 return FramesFactory.getTOD(IERSConventions.IERS_1996, true);
883             case VEIS1950 :
884                 return FramesFactory.getVeis1950();
885             default :
886                 // this should never happen
887                 throw new RuggedInternalError(null);
888         }
889     }
890 
891     /** Select body rotating frame.
892      * @param bodyRotatingFrame body rotating frame identifier
893      * @return body rotating frame
894      */
895     private static Frame selectBodyRotatingFrame(final BodyRotatingFrameId bodyRotatingFrame) {
896 
897         // set up the rotating frame
898         switch (bodyRotatingFrame) {
899             case ITRF :
900                 return FramesFactory.getITRF(IERSConventions.IERS_2010, true);
901             case ITRF_EQUINOX :
902                 return FramesFactory.getITRFEquinox(IERSConventions.IERS_1996, true);
903             case GTOD :
904                 return FramesFactory.getGTOD(IERSConventions.IERS_1996, true);
905             default :
906                 // this should never happen
907                 throw new RuggedInternalError(null);
908         }
909     }
910 
911     /** Select ellipsoid.
912      * @param ellipsoidID reference ellipsoid identifier
913      * @param bodyFrame body rotating frame
914      * @return selected ellipsoid
915      */
916     private static OneAxisEllipsoid selectEllipsoid(final EllipsoidId ellipsoidID, final Frame bodyFrame) {
917 
918         // set up the ellipsoid
919         switch (ellipsoidID) {
920             case GRS80 :
921                 return new OneAxisEllipsoid(6378137.0, 1.0 / 298.257222101, bodyFrame);
922             case WGS84 :
923                 return new OneAxisEllipsoid(Constants.WGS84_EARTH_EQUATORIAL_RADIUS,
924                                             Constants.WGS84_EARTH_FLATTENING,
925                                             bodyFrame);
926             case IERS96 :
927                 return new OneAxisEllipsoid(6378136.49, 1.0 / 298.25645, bodyFrame);
928             case IERS2003 :
929                 return new OneAxisEllipsoid(6378136.6, 1.0 / 298.25642, bodyFrame);
930             default :
931                 // this should never happen
932                 throw new RuggedInternalError(null);
933         }
934 
935     }
936 
937     /** Create DEM intersection algorithm.
938      * @param algorithmID intersection algorithm identifier
939      * @param updater updater used to load Digital Elevation Model tiles
940      * @param maxCachedTiles maximum number of tiles stored in the cache
941      * @param constantElevation constant elevation over ellipsoid
942      * @return selected algorithm
943      */
944     private static IntersectionAlgorithm createAlgorithm(final AlgorithmId algorithmID,
945                                                          final TileUpdater updater, final int maxCachedTiles,
946                                                          final double constantElevation) {
947 
948         // set up the algorithm
949         switch (algorithmID) {
950             case DUVENHAGE :
951                 return new DuvenhageAlgorithm(updater, maxCachedTiles, false);
952             case DUVENHAGE_FLAT_BODY :
953                 return new DuvenhageAlgorithm(updater, maxCachedTiles, true);
954             case BASIC_SLOW_EXHAUSTIVE_SCAN_FOR_TESTS_ONLY :
955                 return new BasicScanAlgorithm(updater, maxCachedTiles);
956             case CONSTANT_ELEVATION_OVER_ELLIPSOID :
957                 return new ConstantElevationAlgorithm(constantElevation);
958             case IGNORE_DEM_USE_ELLIPSOID :
959                 return new IgnoreDEMAlgorithm();
960             default :
961                 // this should never happen
962                 throw new RuggedInternalError(null);
963         }
964     }
965 
966     /** Build a {@link Rugged} instance.
967      * @return built instance
968      */
969     public Rugged build() {
970 
971         if (algorithmID == null) {
972             throw new RuggedException(RuggedMessages.UNINITIALIZED_CONTEXT, "RuggedBuilder.setAlgorithmID()");
973         }
974         if (algorithmID == AlgorithmId.CONSTANT_ELEVATION_OVER_ELLIPSOID) {
975             if (Double.isNaN(constantElevation)) {
976                 throw new RuggedException(RuggedMessages.UNINITIALIZED_CONTEXT, "RuggedBuilder.setConstantElevation()");
977             }
978         } else if (algorithmID != AlgorithmId.IGNORE_DEM_USE_ELLIPSOID) {
979             if (tileUpdater == null) {
980                 throw new RuggedException(RuggedMessages.UNINITIALIZED_CONTEXT, "RuggedBuilder.setDigitalElevationModel()");
981             }
982         }
983         createInterpolatorIfNeeded();
984         return new Rugged(createAlgorithm(algorithmID, tileUpdater, maxCachedTiles, constantElevation), ellipsoid,
985                           lightTimeCorrection, aberrationOfLightCorrection, atmosphericRefraction, scToBody, sensors, name);
986     }
987 }