1   /* Copyright 2002-2016 CS Systèmes d'Information
2    * Licensed to CS Systèmes d'Information (CS) under one or more
3    * contributor license agreements.  See the NOTICE file distributed with
4    * this work for additional information regarding copyright ownership.
5    * CS licenses this file to You under the Apache License, Version 2.0
6    * (the "License"); you may not use this file except in compliance with
7    * the License.  You may obtain a copy of the License at
8    *
9    *   http://www.apache.org/licenses/LICENSE-2.0
10   *
11   * Unless required by applicable law or agreed to in writing, software
12   * distributed under the License is distributed on an "AS IS" BASIS,
13   * WITHOUT WARRANTIES OR CONDITIONS OF ANY KIND, either express or implied.
14   * See the License for the specific language governing permissions and
15   * limitations under the License.
16   */
17  package org.orekit.frames;
18  
19  import java.io.Serializable;
20  import java.util.ArrayList;
21  import java.util.Collection;
22  import java.util.List;
23  
24  import org.apache.commons.math3.analysis.interpolation.HermiteInterpolator;
25  import org.orekit.errors.OrekitException;
26  import org.orekit.errors.OrekitInternalError;
27  import org.orekit.errors.OrekitMessages;
28  import org.orekit.errors.TimeStampedCacheException;
29  import org.orekit.time.AbsoluteDate;
30  import org.orekit.time.TimeFunction;
31  import org.orekit.time.TimeStamped;
32  import org.orekit.utils.Constants;
33  import org.orekit.utils.GenericTimeStampedCache;
34  import org.orekit.utils.IERSConventions;
35  import org.orekit.utils.ImmutableTimeStampedCache;
36  import org.orekit.utils.OrekitConfiguration;
37  import org.orekit.utils.TimeStampedCache;
38  import org.orekit.utils.TimeStampedGenerator;
39  
40  /** This class loads any kind of Earth Orientation Parameter data throughout a large time range.
41   * @author Pascal Parraud
42   */
43  public class EOPHistory implements Serializable {
44  
45      /** Serializable UID. */
46      private static final long serialVersionUID = 20131010L;
47  
48      /** Number of points to use in interpolation. */
49      private static final int INTERPOLATION_POINTS = 4;
50  
51      /**
52       * If this history has any EOP data.
53       *
54       * @see #hasDataFor(AbsoluteDate)
55       */
56      private final boolean hasData;
57  
58      /** EOP history entries. */
59      private final transient ImmutableTimeStampedCache<EOPEntry> cache;
60  
61      /** IERS conventions to which EOP refers. */
62      private final IERSConventions conventions;
63  
64      /** Correction to apply to EOP (may be null). */
65      private final transient TimeFunction<double[]> tidalCorrection;
66  
67      /** Simple constructor.
68       * @param conventions IERS conventions to which EOP refers
69       * @param data the EOP data to use
70       * @param simpleEOP if true, tidal effects are ignored when interpolating EOP
71       * @exception OrekitException if tidal correction model cannot be loaded
72       */
73      protected EOPHistory(final IERSConventions conventions,
74                           final Collection<EOPEntry> data,
75                           final boolean simpleEOP)
76          throws OrekitException {
77          this(conventions, data, simpleEOP ? null : new CachedCorrection(conventions.getEOPTidalCorrection()));
78      }
79  
80      /** Simple constructor.
81       * @param conventions IERS conventions to which EOP refers
82       * @param data the EOP data to use
83       * @param tidalCorrection correction to apply to EOP
84       * @exception OrekitException if tidal correction model cannot be loaded
85       */
86      private EOPHistory(final IERSConventions conventions,
87                           final Collection<EOPEntry> data,
88                           final TimeFunction<double[]> tidalCorrection)
89          throws OrekitException {
90          this.conventions      = conventions;
91          this.tidalCorrection  = tidalCorrection;
92          if (data.size() >= INTERPOLATION_POINTS) {
93              // enough data to interpolate
94              cache = new ImmutableTimeStampedCache<EOPEntry>(INTERPOLATION_POINTS, data);
95              hasData = true;
96          } else {
97              // not enough data to interpolate -> always use null correction
98              cache = ImmutableTimeStampedCache.emptyCache();
99              hasData = false;
100         }
101     }
102 
103     /** Get non-interpolating version of the instance.
104      * @return non-interpolatig version of the instance
105      * @exception OrekitException if tidal correction model cannot be loaded
106      */
107     public EOPHistory getNonInterpolatingEOPHistory()
108         throws OrekitException {
109         return new EOPHistory(conventions, getEntries(), conventions.getEOPTidalCorrection());
110     }
111 
112     /** Check if the instance uses interpolation on tidal corrections.
113      * @return true if the instance uses interpolation on tidal corrections
114      */
115     public boolean usesInterpolation() {
116         return tidalCorrection != null && tidalCorrection instanceof CachedCorrection;
117     }
118 
119     /** Get the IERS conventions to which these EOP apply.
120      * @return IERS conventions to which these EOP apply
121      */
122     public IERSConventions getConventions() {
123         return conventions;
124     }
125 
126     /** Get the date of the first available Earth Orientation Parameters.
127      * @return the start date of the available data
128      */
129     public AbsoluteDate getStartDate() {
130         return this.cache.getEarliest().getDate();
131     }
132 
133     /** Get the date of the last available Earth Orientation Parameters.
134      * @return the end date of the available data
135      */
136     public AbsoluteDate getEndDate() {
137         return this.cache.getLatest().getDate();
138     }
139 
140     /** Get the UT1-UTC value.
141      * <p>The data provided comes from the IERS files. It is smoothed data.</p>
142      * @param date date at which the value is desired
143      * @return UT1-UTC in seconds (0 if date is outside covered range)
144      */
145     public double getUT1MinusUTC(final AbsoluteDate date) {
146         //check if there is data for date
147         if (!this.hasDataFor(date)) {
148             // no EOP data available for this date, we use a default 0.0 offset
149             return (tidalCorrection == null) ? 0.0 : tidalCorrection.value(date)[2];
150         }
151         //we have EOP data -> interpolate offset
152         try {
153             final List<EOPEntry> neighbors = getNeighbors(date);
154             final HermiteInterpolator interpolator = new HermiteInterpolator();
155             final double firstDUT = neighbors.get(0).getUT1MinusUTC();
156             boolean beforeLeap = true;
157             for (final EOPEntry neighbor : neighbors) {
158                 final double dut;
159                 if (neighbor.getUT1MinusUTC() - firstDUT > 0.9) {
160                     // there was a leap second between the entries
161                     dut = neighbor.getUT1MinusUTC() - 1.0;
162                     if (neighbor.getDate().compareTo(date) <= 0) {
163                         beforeLeap = false;
164                     }
165                 } else {
166                     dut = neighbor.getUT1MinusUTC();
167                 }
168                 interpolator.addSamplePoint(neighbor.getDate().durationFrom(date),
169                                             new double[] {
170                                                 dut
171                                             });
172             }
173             double interpolated = interpolator.value(0)[0];
174             if (tidalCorrection != null) {
175                 interpolated += tidalCorrection.value(date)[2];
176             }
177             return beforeLeap ? interpolated : interpolated + 1.0;
178         } catch (TimeStampedCacheException tce) {
179             //this should not happen because of date check above
180             throw new OrekitInternalError(tce);
181         }
182     }
183 
184     /**
185      * Get the entries surrounding a central date.
186      * <p>
187      * See {@link #hasDataFor(AbsoluteDate)} to determine if the cache has data
188      * for {@code central} without throwing an exception.
189      *
190      * @param central central date
191      * @return array of cached entries surrounding specified date
192      * @exception TimeStampedCacheException if EOP data cannot be retrieved
193      */
194     protected List<EOPEntry> getNeighbors(final AbsoluteDate central) throws TimeStampedCacheException {
195         return cache.getNeighbors(central);
196     }
197 
198     /** Get the LoD (Length of Day) value.
199      * <p>The data provided comes from the IERS files. It is smoothed data.</p>
200      * @param date date at which the value is desired
201      * @return LoD in seconds (0 if date is outside covered range)
202      */
203     public double getLOD(final AbsoluteDate date) {
204         //check if there is data for date
205         if (!this.hasDataFor(date)) {
206             // no EOP data available for this date, we use a default null correction
207             return (tidalCorrection == null) ? 0.0 : tidalCorrection.value(date)[3];
208         }
209         //we have EOP data for date -> interpolate correction
210         try {
211             final HermiteInterpolator interpolator = new HermiteInterpolator();
212             for (final EOPEntry entry : getNeighbors(date)) {
213                 interpolator.addSamplePoint(entry.getDate().durationFrom(date),
214                                             new double[] {
215                                                 entry.getLOD()
216                                             });
217             }
218             double interpolated = interpolator.value(0)[0];
219             if (tidalCorrection != null) {
220                 interpolated += tidalCorrection.value(date)[3];
221             }
222             return interpolated;
223         } catch (TimeStampedCacheException tce) {
224             // this should not happen because of date check above
225             throw new OrekitInternalError(tce);
226         }
227     }
228 
229     /** Get the pole IERS Reference Pole correction.
230      * <p>The data provided comes from the IERS files. It is smoothed data.</p>
231      * @param date date at which the correction is desired
232      * @return pole correction ({@link PoleCorrection#NULL_CORRECTION
233      * PoleCorrection.NULL_CORRECTION} if date is outside covered range)
234      */
235     public PoleCorrection getPoleCorrection(final AbsoluteDate date) {
236         // check if there is data for date
237         if (!this.hasDataFor(date)) {
238             // no EOP data available for this date, we use a default null correction
239             if (tidalCorrection == null) {
240                 return PoleCorrection.NULL_CORRECTION;
241             } else {
242                 final double[] correction = tidalCorrection.value(date);
243                 return new PoleCorrection(correction[0], correction[1]);
244             }
245         }
246         //we have EOP data for date -> interpolate correction
247         try {
248             final HermiteInterpolator interpolator = new HermiteInterpolator();
249             for (final EOPEntry entry : getNeighbors(date)) {
250                 interpolator.addSamplePoint(entry.getDate().durationFrom(date),
251                                             new double[] {
252                                                 entry.getX(), entry.getY()
253                                             });
254             }
255             final double[] interpolated = interpolator.value(0);
256             if (tidalCorrection != null) {
257                 final double[] correction = tidalCorrection.value(date);
258                 interpolated[0] += correction[0];
259                 interpolated[1] += correction[1];
260             }
261             return new PoleCorrection(interpolated[0], interpolated[1]);
262         } catch (TimeStampedCacheException tce) {
263             // this should not happen because of date check above
264             throw new OrekitInternalError(tce);
265         }
266     }
267 
268     /** Get the correction to the nutation parameters for equinox-based paradigm.
269      * <p>The data provided comes from the IERS files. It is smoothed data.</p>
270      * @param date date at which the correction is desired
271      * @return nutation correction in longitude ΔΨ and in obliquity Δε
272      * (zero if date is outside covered range)
273      */
274     public double[] getEquinoxNutationCorrection(final AbsoluteDate date) {
275         // check if there is data for date
276         if (!this.hasDataFor(date)) {
277             // no EOP data available for this date, we use a default null correction
278             return new double[2];
279         }
280         //we have EOP data for date -> interpolate correction
281         try {
282             final HermiteInterpolator interpolator = new HermiteInterpolator();
283             for (final EOPEntry entry : getNeighbors(date)) {
284                 interpolator.addSamplePoint(entry.getDate().durationFrom(date),
285                                             new double[] {
286                                                 entry.getDdPsi(), entry.getDdEps()
287                                             });
288             }
289             return interpolator.value(0);
290         } catch (TimeStampedCacheException tce) {
291             // this should not happen because of date check above
292             throw new OrekitInternalError(tce);
293         }
294     }
295 
296     /** Get the correction to the nutation parameters for Non-Rotating Origin paradigm.
297      * <p>The data provided comes from the IERS files. It is smoothed data.</p>
298      * @param date date at which the correction is desired
299      * @return nutation correction in Celestial Intermediat Pole coordinates
300      * δX and δY (zero if date is outside covered range)
301      */
302     public double[] getNonRotatinOriginNutationCorrection(final AbsoluteDate date) {
303         // check if there is data for date
304         if (!this.hasDataFor(date)) {
305             // no EOP data available for this date, we use a default null correction
306             return new double[2];
307         }
308         //we have EOP data for date -> interpolate correction
309         try {
310             final HermiteInterpolator interpolator = new HermiteInterpolator();
311             for (final EOPEntry entry : getNeighbors(date)) {
312                 interpolator.addSamplePoint(entry.getDate().durationFrom(date),
313                                             new double[] {
314                                                 entry.getDx(), entry.getDy()
315                                             });
316             }
317             return interpolator.value(0);
318         } catch (TimeStampedCacheException tce) {
319             // this should not happen because of date check above
320             throw new OrekitInternalError(tce);
321         }
322     }
323 
324     /** Check Earth orientation parameters continuity.
325      * @param maxGap maximal allowed gap between entries (in seconds)
326      * @exception OrekitException if there are holes in the data sequence
327      */
328     public void checkEOPContinuity(final double maxGap) throws OrekitException {
329         TimeStamped preceding = null;
330         for (final TimeStamped current : this.cache.getAll()) {
331 
332             // compare the dates of preceding and current entries
333             if ((preceding != null) && ((current.getDate().durationFrom(preceding.getDate())) > maxGap)) {
334                 throw new OrekitException(OrekitMessages.MISSING_EARTH_ORIENTATION_PARAMETERS_BETWEEN_DATES,
335                                           preceding.getDate(), current.getDate());
336             }
337 
338             // prepare next iteration
339             preceding = current;
340 
341         }
342     }
343 
344     /**
345      * Check if the cache has data for the given date using
346      * {@link #getStartDate()} and {@link #getEndDate()}.
347      *
348      * @param date the requested date
349      * @return true if the {@link #cache} has data for the requested date, false
350      *         otherwise.
351      */
352     protected boolean hasDataFor(final AbsoluteDate date) {
353         /*
354          * when there is no EOP data, short circuit getStartDate, which will
355          * throw an exception
356          */
357         return this.hasData && this.getStartDate().compareTo(date) <= 0 &&
358                date.compareTo(this.getEndDate()) <= 0;
359     }
360 
361     /** Get a non-modifiable view of the EOP entries.
362      * @return non-modifiable view of the EOP entries
363      */
364     List<EOPEntry> getEntries() {
365         return cache.getAll();
366     }
367 
368     /** Replace the instance with a data transfer object for serialization.
369      * <p>
370      * This intermediate class serializes only the frame key.
371      * </p>
372      * @return data transfer object that will be serialized
373      */
374     private Object writeReplace() {
375         return new DataTransferObject(conventions, getEntries(), tidalCorrection == null);
376     }
377 
378     /** Internal class used only for serialization. */
379     private static class DataTransferObject implements Serializable {
380 
381         /** Serializable UID. */
382         private static final long serialVersionUID = 20131010L;
383 
384         /** IERS conventions. */
385         private final IERSConventions conventions;
386 
387         /** EOP entries. */
388         private final List<EOPEntry> entries;
389 
390         /** Indicator for simple interpolation without tidal effects. */
391         private final boolean simpleEOP;
392 
393         /** Simple constructor.
394          * @param conventions IERS conventions to which EOP refers
395          * @param entries the EOP data to use
396          * @param simpleEOP if true, tidal effects are ignored when interpolating EOP
397          */
398         DataTransferObject(final IERSConventions conventions,
399                                   final List<EOPEntry> entries,
400                                   final boolean simpleEOP) {
401             this.conventions = conventions;
402             this.entries     = entries;
403             this.simpleEOP   = simpleEOP;
404         }
405 
406         /** Replace the deserialized data transfer object with a {@link EOPHistory}.
407          * @return replacement {@link EOPHistory}
408          */
409         private Object readResolve() {
410             try {
411                 // retrieve a managed frame
412                 return new EOPHistory(conventions, entries, simpleEOP);
413             } catch (OrekitException oe) {
414                 throw new OrekitInternalError(oe);
415             }
416         }
417 
418     }
419 
420     /** Internal class for caching tidal correction. */
421     private static class TidalCorrectionEntry implements TimeStamped {
422 
423         /** Entry date. */
424         private final AbsoluteDate date;
425 
426         /** Correction. */
427         private final double[] correction;
428 
429         /** Simple constructor.
430          * @param date entry date
431          * @param correction correction on the EOP parameters (xp, yp, ut1, lod)
432          */
433         TidalCorrectionEntry(final AbsoluteDate date, final double[] correction) {
434             this.date       = date;
435             this.correction = correction;
436         }
437 
438         /** {@inheritDoc} */
439         @Override
440         public AbsoluteDate getDate() {
441             return date;
442         }
443 
444     }
445 
446     /** Local generator for thread-safe cache. */
447     private static class CachedCorrection
448         implements TimeFunction<double[]>, TimeStampedGenerator<TidalCorrectionEntry> {
449 
450         /** Correction to apply to EOP (may be null). */
451         private final TimeFunction<double[]> tidalCorrection;
452 
453         /** Step between generated entries. */
454         private final double step;
455 
456         /** Tidal corrections entries cache. */
457         private final TimeStampedCache<TidalCorrectionEntry> cache;
458 
459         /** Simple constructor.
460          * @param tidalCorrection function computing the tidal correction
461          */
462         CachedCorrection(final TimeFunction<double[]> tidalCorrection) {
463             this.step            = 60 * 60;
464             this.tidalCorrection = tidalCorrection;
465             this.cache           =
466                 new GenericTimeStampedCache<TidalCorrectionEntry>(8,
467                                                                   OrekitConfiguration.getCacheSlotsNumber(),
468                                                                   Constants.JULIAN_DAY * 30,
469                                                                   Constants.JULIAN_DAY,
470                                                                   this,
471                                                                   TidalCorrectionEntry.class);
472         }
473 
474         /** {@inheritDoc} */
475         @Override
476         public double[] value(final AbsoluteDate date) {
477             try {
478                 // set up an interpolator
479                 final HermiteInterpolator interpolator = new HermiteInterpolator();
480                 for (final TidalCorrectionEntry entry : cache.getNeighbors(date)) {
481                     interpolator.addSamplePoint(entry.date.durationFrom(date), entry.correction);
482                 }
483 
484                 // interpolate to specified date
485                 return interpolator.value(0.0);
486             } catch (TimeStampedCacheException tsce) {
487                 // this should never happen
488                 throw new OrekitInternalError(tsce);
489             }
490         }
491 
492         /** {@inheritDoc} */
493         @Override
494         public List<TidalCorrectionEntry> generate(final TidalCorrectionEntry existing, final AbsoluteDate date) {
495 
496             final List<TidalCorrectionEntry> generated = new ArrayList<TidalCorrectionEntry>();
497 
498             if (existing == null) {
499 
500                 // no prior existing entries, just generate a first set
501                 for (int i = -cache.getNeighborsSize() / 2; generated.size() < cache.getNeighborsSize(); ++i) {
502                     final AbsoluteDate t = date.shiftedBy(i * step);
503                     generated.add(new TidalCorrectionEntry(t, tidalCorrection.value(t)));
504                 }
505 
506             } else {
507 
508                 // some entries have already been generated
509                 // add the missing ones up to specified date
510 
511                 AbsoluteDate t = existing.getDate();
512                 if (date.compareTo(t) > 0) {
513                     // forward generation
514                     do {
515                         t = t.shiftedBy(step);
516                         generated.add(new TidalCorrectionEntry(t, tidalCorrection.value(t)));
517                     } while (t.compareTo(date) <= 0);
518                 } else {
519                     // backward generation
520                     do {
521                         t = t.shiftedBy(-step);
522                         generated.add(0, new TidalCorrectionEntry(t, tidalCorrection.value(t)));
523                     } while (t.compareTo(date) >= 0);
524                 }
525             }
526 
527             // return the generated transforms
528             return generated;
529 
530         }
531     }
532 
533 }