1   /* Copyright 2002-2012 Space Applications Services
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.files.sp3;
18  
19  import java.util.ArrayList;
20  import java.util.Collection;
21  import java.util.Collections;
22  import java.util.Iterator;
23  import java.util.LinkedHashMap;
24  import java.util.List;
25  import java.util.Map;
26  import java.util.SortedSet;
27  import java.util.TreeSet;
28  
29  import org.hipparchus.geometry.euclidean.threed.Vector3D;
30  import org.hipparchus.util.FastMath;
31  import org.hipparchus.util.Precision;
32  import org.orekit.errors.OrekitException;
33  import org.orekit.errors.OrekitMessages;
34  import org.orekit.files.general.EphemerisFile;
35  import org.orekit.frames.Frame;
36  import org.orekit.frames.Transform;
37  import org.orekit.gnss.IGSUtils;
38  import org.orekit.time.AbsoluteDate;
39  import org.orekit.time.ChronologicalComparator;
40  import org.orekit.utils.PVCoordinates;
41  
42  /**
43   * Represents a parsed SP3 orbit file.
44   * @author Thomas Neidhart
45   * @author Evan Ward
46   */
47  public class SP3 implements EphemerisFile<SP3Coordinate, SP3Segment> {
48  
49      /** Header.
50       * @since 12.0
51       */
52      private final SP3Header header;
53  
54      /** Standard gravitational parameter in m³ / s². */
55      private final double mu;
56  
57      /** Number of samples to use when interpolating. */
58      private final int interpolationSamples;
59  
60      /** Reference frame. */
61      private final Frame frame;
62  
63      /** A map containing satellite information. */
64      private final Map<String, SP3Ephemeris> satellites;
65  
66      /**
67       * Create a new SP3 file object.
68       *
69       * @param mu                   is the standard gravitational parameter in m³ / s².
70       * @param interpolationSamples number of samples to use in interpolation.
71       * @param frame                reference frame
72       */
73      public SP3(final double mu, final int interpolationSamples, final Frame frame) {
74          this(new SP3Header(), mu, interpolationSamples, frame);
75      }
76  
77      /**
78       * Create a new SP3 file object.
79       *
80       * @param header header
81       * @param mu                   is the standard gravitational parameter in m³ / s².
82       * @param interpolationSamples number of samples to use in interpolation.
83       * @param frame                reference frame
84       * @since 12.1
85       */
86      public SP3(final SP3Header header,
87                 final double mu, final int interpolationSamples, final Frame frame) {
88          this.header               = header;
89          this.mu                   = mu;
90          this.interpolationSamples = interpolationSamples;
91          this.frame                = frame;
92          this.satellites           = new LinkedHashMap<>(); // must be linked hash map to preserve order of satellites in the file
93      }
94  
95      /** Check file is valid.
96       * @param parsing if true, we are parsing an existing file, and are more lenient
97       * in order to accept some common errors (like between 86 and 99 satellites
98       * in SP3a, SP3b or SP3c files)
99       * @param fileName file name to generate the error message
100      * @exception OrekitException if file is not valid
101      */
102     public void validate(final boolean parsing, final String fileName) throws OrekitException {
103 
104         // check available data
105         final SortedSet<AbsoluteDate> epochs = new TreeSet<>(new ChronologicalComparator());
106         boolean hasAccuracy = false;
107         for (final Map.Entry<String, SP3Ephemeris> entry : satellites.entrySet()) {
108             SP3Coordinate previous = null;
109             for (final SP3Segment segment : entry.getValue().getSegments()) {
110                 for (final SP3Coordinate coordinate : segment.getCoordinates()) {
111                     final AbsoluteDate previousDate = previous == null ? header.getEpoch() : previous.getDate();
112                     final double       nbSteps      = coordinate.getDate().durationFrom(previousDate) / header.getEpochInterval();
113                     if (FastMath.abs(nbSteps - FastMath.rint(nbSteps)) > 0.001) {
114                         // not an integral number of steps
115                         throw new OrekitException(OrekitMessages.INCONSISTENT_SAMPLING_DATE,
116                                                   previousDate.shiftedBy(FastMath.rint(nbSteps) * header.getEpochInterval()),
117                                                   coordinate.getDate());
118                     }
119                     epochs.add(coordinate.getDate());
120                     previous = coordinate;
121                     hasAccuracy |= !(coordinate.getPositionAccuracy() == null &&
122                                     coordinate.getVelocityAccuracy() == null &&
123                                     Double.isNaN(coordinate.getClockAccuracy()) &&
124                                     Double.isNaN(coordinate.getClockRateAccuracy()));
125                 }
126             }
127         }
128 
129         // check versions limitations
130         if (getSatelliteCount() > getMaxAllowedSatCount(parsing)) {
131             throw new OrekitException(OrekitMessages.SP3_TOO_MANY_SATELLITES_FOR_VERSION,
132                                       header.getVersion(), getMaxAllowedSatCount(parsing), getSatelliteCount(),
133                                       fileName);
134         }
135 
136         header.validate(parsing, hasAccuracy, fileName);
137 
138         // check epochs
139         if (epochs.size() != header.getNumberOfEpochs()) {
140             throw new OrekitException(OrekitMessages.SP3_NUMBER_OF_EPOCH_MISMATCH,
141                                       epochs.size(), fileName, header.getNumberOfEpochs());
142         }
143 
144     }
145 
146     /** Get the header.
147      * @return header
148      * @since 12.0
149      */
150     public SP3Header getHeader() {
151         return header;
152     }
153 
154     /** Get maximum number of satellites allowed for format version.
155      * @param parsing if true, we are parsing an existing file, and are more lenient
156      * in order to accept some common errors (like between 86 and 99 satellites
157      * in SP3a, SP3b or SP3c files)
158      * @return maximum number of satellites allowed for format version
159      * @since 12.0
160      */
161     private int getMaxAllowedSatCount(final boolean parsing) {
162         return header.getVersion() < 'd' ? (parsing ? 99 : 85) : 999;
163     }
164 
165     /** Splice several SP3 files together.
166      * <p>
167      * Splicing SP3 files is intended to be used when continuous computation
168      * covering more than one file is needed. The files should all have the exact same
169      * metadata: {@link SP3Header#getType() type}, {@link SP3Header#getTimeSystem() time system},
170      * {@link SP3Header#getCoordinateSystem() coordinate system}, except for satellite accuracy
171      * which can be different from one file to the next one, and some satellites may
172      * be missing in some files… Once sorted (which is done internally), if the gap between
173      * segments from two file is at most {@link SP3Header#getEpochInterval() epoch interval},
174      * then the segments are merged as one segment, otherwise the segments are kept separated.
175      * </p>
176      * <p>
177      * The spliced file only contains the satellites that were present in all files.
178      * Satellites present in some files and absent from other files are silently
179      * dropped.
180      * </p>
181      * <p>
182      * Depending on producer, successive SP3 files either have a gap between the last
183      * entry of one file and the first entry of the next file (for example files with
184      * a 5 minutes epoch interval may end at 23:55 and the next file start at 00:00),
185      * or both files have one point exactly at the splicing date (i.e. 24:00 one day
186      * and 00:00 next day). In the later case, the last point of the early file is dropped
187      * and the first point of the late file takes precedence, hence only one point remains
188      * in the spliced file ; this design choice is made to enforce continuity and
189      * regular interpolation.
190      * </p>
191      * @param sp3 SP3 files to merge
192      * @return merged SP3
193      * @since 12.0
194      */
195     public static SP3 splice(final Collection<SP3> sp3) {
196 
197         // sort the files
198         final ChronologicalComparator comparator = new ChronologicalComparator();
199         final SortedSet<SP3> sorted = new TreeSet<>((s1, s2) -> comparator.compare(s1.header.getEpoch(), s2.header.getEpoch()));
200         sorted.addAll(sp3);
201 
202         // prepare spliced file
203         final SP3 first   = sorted.first();
204         final SP3 spliced = new SP3(first.mu, first.interpolationSamples, first.frame);
205         spliced.header.setVersion(first.header.getVersion());
206         spliced.header.setFilter(first.header.getFilter());
207         spliced.header.setType(first.header.getType());
208         spliced.header.setTimeSystem(first.header.getTimeSystem());
209         spliced.header.setDataUsed(first.header.getDataUsed());
210         spliced.header.setEpoch(first.header.getEpoch());
211         spliced.header.setGpsWeek(first.header.getGpsWeek());
212         spliced.header.setSecondsOfWeek(first.header.getSecondsOfWeek());
213         spliced.header.setModifiedJulianDay(first.header.getModifiedJulianDay());
214         spliced.header.setDayFraction(first.header.getDayFraction());
215         spliced.header.setEpochInterval(first.header.getEpochInterval());
216         spliced.header.setCoordinateSystem(first.header.getCoordinateSystem());
217         spliced.header.setOrbitTypeKey(first.header.getOrbitTypeKey());
218         spliced.header.setAgency(first.header.getAgency());
219         spliced.header.setPosVelBase(first.header.getPosVelBase());
220         spliced.header.setClockBase(first.header.getClockBase());
221         first.header.getComments().forEach(spliced.header::addComment);
222 
223         // identify the satellites that are present in all files
224         final List<String> commonSats = new ArrayList<>(first.header.getSatIds());
225         for (final SP3 current : sorted) {
226             for (final Iterator<String> iter = commonSats.iterator(); iter.hasNext();) {
227                 final String sat = iter.next();
228                 if (!current.containsSatellite(sat)) {
229                     iter.remove();
230                     break;
231                 }
232             }
233         }
234 
235         // create the spliced list
236         for (final String sat : commonSats) {
237             spliced.addSatellite(sat);
238         }
239 
240         // in order to be conservative, we keep the worst accuracy from all SP3 files for this satellite
241         for (int i = 0; i < commonSats.size(); ++i) {
242             final String sat = commonSats.get(i);
243             double accuracy = 0;
244             for (final SP3 current : sorted) {
245                 accuracy = FastMath.max(accuracy, current.header.getAccuracy(sat));
246             }
247             spliced.header.setAccuracy(i, accuracy);
248         }
249 
250         // splice files
251         SP3 previous = null;
252         int epochCount = 0;
253         for (final SP3 current : sorted) {
254 
255             epochCount += current.header.getNumberOfEpochs();
256             if (previous != null) {
257 
258                 // check metadata and check if we should drop the last entry of previous file
259                 final boolean dropLast = current.checkSplice(previous);
260                 if (dropLast) {
261                     --epochCount;
262                 }
263 
264                 // append the pending data from previous file
265                 for (final Map.Entry<String, SP3Ephemeris> entry : previous.satellites.entrySet()) {
266                     if (commonSats.contains(entry.getKey())) {
267                         final SP3Ephemeris splicedEphemeris = spliced.getEphemeris(entry.getKey());
268                         for (final SP3Segment segment : entry.getValue().getSegments()) {
269                             final List<SP3Coordinate> coordinates = segment.getCoordinates();
270                             for (int i = 0; i < coordinates.size() - (dropLast ? 1 : 0); ++i) {
271                                 splicedEphemeris.addCoordinate(coordinates.get(i), spliced.header.getEpochInterval());
272                             }
273                         }
274                     }
275                 }
276 
277             }
278 
279             previous = current;
280 
281         }
282         spliced.header.setNumberOfEpochs(epochCount);
283 
284         // append the pending data from last file
285         for (final Map.Entry<String, SP3Ephemeris> entry : previous.satellites.entrySet()) {
286             if (commonSats.contains(entry.getKey())) {
287                 final SP3Ephemeris splicedEphemeris = spliced.getEphemeris(entry.getKey());
288                 for (final SP3Segment segment : entry.getValue().getSegments()) {
289                     for (final SP3Coordinate coordinate : segment.getCoordinates()) {
290                         splicedEphemeris.addCoordinate(coordinate, spliced.header.getEpochInterval());
291                     }
292                 }
293             }
294         }
295 
296         return spliced;
297 
298     }
299 
300     /** Change the frame of an SP3 file.
301      * @param original original SP3 file
302      * @param newFrame frame to use for the changed SP3 file
303      * @return changed SP3 file
304      * @since 12.1
305      */
306     public static SP3 changeFrame(final SP3 original, final Frame newFrame) {
307 
308         // copy header
309         final SP3Header header         = new SP3Header();
310         final SP3Header originalHeader = original.header;
311         header.setVersion(originalHeader.getVersion());
312         header.setFilter(originalHeader.getFilter());
313         header.setType(originalHeader.getType());
314         header.setTimeSystem(originalHeader.getTimeSystem());
315         header.setDataUsed(originalHeader.getDataUsed());
316         header.setEpoch(originalHeader.getEpoch());
317         header.setGpsWeek(originalHeader.getGpsWeek());
318         header.setSecondsOfWeek(originalHeader.getSecondsOfWeek());
319         header.setModifiedJulianDay(originalHeader.getModifiedJulianDay());
320         header.setDayFraction(originalHeader.getDayFraction());
321         header.setEpochInterval(originalHeader.getEpochInterval());
322         header.setOrbitTypeKey(originalHeader.getOrbitTypeKey());
323         header.setAgency(originalHeader.getAgency());
324         header.setPosVelBase(originalHeader.getPosVelBase());
325         header.setClockBase(originalHeader.getClockBase());
326         header.setNumberOfEpochs(originalHeader.getNumberOfEpochs());
327 
328         originalHeader.getComments().forEach(header::addComment);
329 
330         // change the frame in the header
331         header.setCoordinateSystem(IGSUtils.frameName(newFrame));
332 
333         // prepare file
334         final SP3 changed = new SP3(header,
335                                     original.mu,
336                                     original.interpolationSamples,
337                                     newFrame);
338 
339         // first loop to add satellite ids only
340         final List<String> ids = originalHeader.getSatIds();
341         ids.forEach(changed::addSatellite);
342 
343         // now that the header knows the number of satellites,
344         // it can lazily allocate the accuracies array,
345         // so we can perform a second loop to set individual accuracies
346         for (int i = 0; i < ids.size(); ++i) {
347             header.setAccuracy(i, originalHeader.getAccuracy(ids.get(i)));
348         }
349 
350         // convert data to new frame
351         for (final Map.Entry<String, SP3Ephemeris> entry : original.satellites.entrySet()) {
352             final SP3Ephemeris originalEphemeris = original.getEphemeris(entry.getKey());
353             final SP3Ephemeris changedEphemeris  = changed.getEphemeris(entry.getKey());
354             for (final SP3Segment segment : originalEphemeris.getSegments()) {
355                 for (SP3Coordinate c : segment.getCoordinates()) {
356                     final Transform t = originalEphemeris.getFrame().getTransformTo(newFrame, c.getDate());
357                     final PVCoordinates newPV = t.transformPVCoordinates(c);
358                     final Vector3D      newP  = newPV.getPosition();
359                     final Vector3D      newPA = c.getPositionAccuracy() == null ?
360                                                 null :
361                                                 t.transformVector(c.getPositionAccuracy());
362                     final Vector3D      newV  = c.getVelocity() == null ? Vector3D.ZERO : newPV.getVelocity();
363                     final Vector3D      newVA = c.getVelocityAccuracy() == null ?
364                                                 null :
365                                                 t.transformVector(c.getVelocityAccuracy());
366                     final SP3Coordinate newC  = new SP3Coordinate(c.getDate(),
367                                                                   newP, newPA, newV, newVA,
368                                                                   c.getClockCorrection(),
369                                                                   c.getClockAccuracy(),
370                                                                   c.getClockRateChange(),
371                                                                   c.getClockRateAccuracy(),
372                                                                   c.hasClockEvent(),
373                                                                   c.hasClockPrediction(),
374                                                                   c.hasOrbitManeuverEvent(),
375                                                                   c.hasOrbitPrediction());
376                     changedEphemeris.addCoordinate(newC, originalHeader.getEpochInterval());
377                 }
378             }
379         }
380 
381         return changed;
382 
383     }
384 
385     /** Check if instance can be spliced after previous one.
386      * @param previous SP3 file (should already be sorted to be before current instance), can be null
387      * @return true if last entry of previous file should be dropped as first entry of current file
388      * is at very close date and will take precedence
389      * @exception OrekitException if metadata are incompatible
390      * @since 12.0
391      */
392     private boolean checkSplice(final SP3 previous) throws OrekitException {
393 
394         if (!(previous.header.getType()             == header.getType()                  &&
395               previous.header.getTimeSystem()       == header.getTimeSystem()            &&
396               previous.header.getOrbitType()        == header.getOrbitType()             &&
397               previous.header.getCoordinateSystem().equals(header.getCoordinateSystem()) &&
398               previous.header.getDataUsed().equals(header.getDataUsed())                 &&
399               previous.header.getAgency().equals(header.getAgency()))) {
400             throw new OrekitException(OrekitMessages.SP3_INCOMPATIBLE_FILE_METADATA);
401         }
402 
403         boolean dropLast = false;
404         for (final Map.Entry<String, SP3Ephemeris> entry : previous.satellites.entrySet()) {
405             final SP3Ephemeris previousEphem = entry.getValue();
406             final SP3Ephemeris currentEphem  = satellites.get(entry.getKey());
407             if (currentEphem != null) {
408                 if (!(previousEphem.getAvailableDerivatives()    == currentEphem.getAvailableDerivatives() &&
409                       previousEphem.getFrame()                   == currentEphem.getFrame()                &&
410                       previousEphem.getInterpolationSamples()    == currentEphem.getInterpolationSamples() &&
411                       Precision.equals(previousEphem.getMu(),       currentEphem.getMu(), 2))) {
412                     throw new OrekitException(OrekitMessages.SP3_INCOMPATIBLE_SATELLITE_MEDATADA,
413                                               entry.getKey());
414                 } else {
415                     final double dt = currentEphem.getStart().durationFrom(previousEphem.getStop());
416                     dropLast = dt < 0.001 * header.getEpochInterval();
417                 }
418             }
419         }
420 
421         return dropLast;
422 
423     }
424 
425     /** Add a new satellite with a given identifier to the list of
426      * stored satellites.
427      * @param satId the satellite identifier
428      */
429     public void addSatellite(final String satId) {
430         header.addSatId(satId);
431         satellites.putIfAbsent(satId, new SP3Ephemeris(satId, mu, frame, interpolationSamples, header.getFilter()));
432     }
433 
434     @Override
435     public Map<String, SP3Ephemeris> getSatellites() {
436         return Collections.unmodifiableMap(satellites);
437     }
438 
439     /** Get an ephemeris.
440      * @param index index of the satellite
441      * @return satellite ephemeris
442      * @since 12.0
443      */
444     public SP3Ephemeris getEphemeris(final int index) {
445         int n = index;
446         for (final Map.Entry<String, SP3Ephemeris> entry : satellites.entrySet()) {
447             if (n == 0) {
448                 return entry.getValue();
449             }
450             n--;
451         }
452 
453         // satellite not found
454         throw new OrekitException(OrekitMessages.INVALID_SATELLITE_ID, index);
455 
456     }
457 
458     /** Get an ephemeris.
459      * @param satId satellite identifier
460      * @return satellite ephemeris, or null if not found
461      * @since 12.0
462      */
463     public SP3Ephemeris getEphemeris(final String satId) {
464         final SP3Ephemeris ephemeris = satellites.get(satId);
465         if (ephemeris == null) {
466             throw new OrekitException(OrekitMessages.INVALID_SATELLITE_ID, satId);
467         } else {
468             return ephemeris;
469         }
470     }
471 
472     /** Get the number of satellites contained in this orbit file.
473      * @return the number of satellites
474      */
475     public int getSatelliteCount() {
476         return satellites.size();
477     }
478 
479     /** Tests whether a satellite with the given id is contained in this orbit
480      * file.
481      * @param satId the satellite id
482      * @return {@code true} if the satellite is contained in the file,
483      *         {@code false} otherwise
484      */
485     public boolean containsSatellite(final String satId) {
486         return header.getSatIds().contains(satId);
487     }
488 
489 }