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