SingleParameterFitter.java

  1. /* Copyright 2022-2025 Luc Maisonobe
  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.frames;

  19. import java.io.Serializable;
  20. import java.util.List;
  21. import java.util.ListIterator;
  22. import java.util.function.ToDoubleFunction;

  24. import org.hipparchus.util.FastMath;
  25. import org.hipparchus.util.MathUtils;
  26. import org.orekit.utils.Constants;
  27. import org.orekit.utils.SecularAndHarmonic;

  29. /** Fitter for one Earth Orientation Parameter.
  30.  * @see PredictedEOPHistory
  31.  * @see EOPFitter
  32.  * @see SecularAndHarmonic
  33.  * @since 12.0
  34.  * @author Luc Maisonobe
  35.  */
  36. public class SingleParameterFitter implements Serializable {

  38.     /** Sun pulsation, one year period. */
  39.     public static final double SUN_PULSATION = MathUtils.TWO_PI / Constants.JULIAN_YEAR;

  41.     /** Moon pulsation (one Moon draconic period). */
  42.     public static final double MOON_DRACONIC_PULSATION = MathUtils.TWO_PI / (27.212221 * Constants.JULIAN_DAY);

  44.     /** Serializable UID. */
  45.     private static final long serialVersionUID = 20230309L;

  47.     /** Time constant of the exponential decay weight. */
  48.     private final double timeConstant;

  50.     /** Convergence on fitted parameter. */
  51.     private final double convergence;

  53.     /** Degree of the polynomial model. */
  54.     private final int degree;

  56.     /** Pulsations of harmonic part (rad/s). */
  57.     private final double[] pulsations;

  59.     /** Simple constructor.
  60.      * @param timeConstant time constant \(\tau\) of the exponential decay weight, point weight is \(e^{\frac{t-t_0}{\tau}}\),
  61.      * i.e. points far in the past before \(t_0\) have smaller weights
  62.      * @param convergence convergence on fitted parameter
  63.      * @param degree degree of the polynomial model
  64.      * @param pulsations pulsations of harmonic part (rad/s)
  65.      * @see #createDefaultDut1FitterShortTermPrediction()
  66.      * @see #createDefaultDut1FitterLongTermPrediction()
  67.      * @see #createDefaultPoleFitterShortTermPrediction()
  68.      * @see #createDefaultPoleFitterLongTermPrediction()
  69.      * @see #createDefaultNutationFitterShortTermPrediction()
  70.      * @see #createDefaultNutationFitterLongTermPrediction()
  71.      * @see SecularAndHarmonic
  72.      * @since 12.0.1
  73.      */
  74.     public SingleParameterFitter(final double timeConstant, final double convergence,
  75.                                  final int degree, final double... pulsations) {
  76.         this.timeConstant    = timeConstant;
  77.         this.convergence     = convergence;
  78.         this.degree          = degree;
  79.         this.pulsations      = pulsations.clone();
  80.     }

  82.     /** Perform secular and harmonic fitting.
  83.      * @param rawHistory EOP history to fit
  84.      * @param extractor extractor for Earth Orientation Parameter
  85.      * @return configured fitter
  86.      */
  87.     public SecularAndHarmonic fit(final EOPHistory rawHistory, final ToDoubleFunction<EOPEntry> extractor) {

  89.         final List<EOPEntry> rawEntries = rawHistory.getEntries();
  90.         final EOPEntry       last       = rawEntries.get(rawEntries.size() - 1);

  92.         // create fitter
  93.         final SecularAndHarmonic sh = new SecularAndHarmonic(degree, pulsations);

  95.         // set up convergence
  96.         sh.setConvergenceRMS(convergence);

  98.         // set up reference date and initial guess to a constant value
  99.         final double[] initialGuess = new double[degree + 1 + 2 * pulsations.length];
  100.         initialGuess[0] = extractor.applyAsDouble(last);
  101.         sh.resetFitting(last.getDate(), initialGuess);

  103.         // sample history
  104.         final ListIterator<EOPEntry> backwardIterator = rawEntries.listIterator(rawEntries.size());
  105.         while (backwardIterator.hasPrevious()) {
  106.             final EOPEntry entry = backwardIterator.previous();
  107.             sh.addWeightedPoint(entry.getDate(), extractor.applyAsDouble(entry),
  108.                                 FastMath.exp(entry.getDate().durationFrom(last.getDate()) / timeConstant));
  109.         }

  111.         // perform fitting
  112.         sh.fit();

  114.         return sh;

  116.     }

  118.     /** Create fitter with default parameters adapted for fitting orientation parameters dUT1 and LOD
  119.      * for short term prediction.
  120.      * <p>
  121.      * The main difference between these settings and {@link #createDefaultDut1FitterLongTermPrediction()
  122.      * the settings for long prediction} is the much smaller \(\tau\). This means more
  123.      * weight is set to the points at the end of the history, hence forcing the fitted prediction
  124.      * model to be closer to these points, hence the prediction error to be smaller just after
  125.      * raw history end. On the other hand, this implies that the model will diverge on long term.
  126.      * These settings are intended when prediction is used for at most 5 days after raw EOP end.
  127.      * </p>
  128.      * <ul>
  129.      *   <li>time constant \(\tau\) of the exponential decay set to 6 {@link Constants#JULIAN_DAY days}</li>
  130.      *   <li>convergence set to 10⁻¹² s</li>
  131.      *   <li>polynomial part set to degree 3</li>
  132.      *   <li>one harmonic term at {@link #SUN_PULSATION}}</li>
  133.      *   <li>one harmonic term at 2 times {@link #SUN_PULSATION}}</li>
  134.      *   <li>one harmonic term at 3 times {@link #SUN_PULSATION}}</li>
  135.      *   <li>one harmonic term at {@link #MOON_DRACONIC_PULSATION}}</li>
  136.      *   <li>one harmonic term at 2 times {@link #MOON_DRACONIC_PULSATION}}</li>
  137.      *   <li>one harmonic term at 3 times {@link #MOON_DRACONIC_PULSATION}}</li>
  138.      * </ul>
  139.      * @return fitter with default configuration for orientation parameters dUT1 and LOD
  140.      * @see #createDefaultDut1FitterShortTermPrediction()
  141.      */
  142.     public static SingleParameterFitter createDefaultDut1FitterShortTermPrediction() {
  143.         return new SingleParameterFitter(6 * Constants.JULIAN_DAY, 1.0e-12, 3,
  144.                                          SUN_PULSATION, 2 * SUN_PULSATION, 3 * SUN_PULSATION,
  145.                                          MOON_DRACONIC_PULSATION, 2 * MOON_DRACONIC_PULSATION, 3 * MOON_DRACONIC_PULSATION);
  146.     }

  148.     /** Create fitter with default parameters adapted for fitting orientation parameters dUT1 and LOD
  149.      * for long term prediction.
  150.      * <p>
  151.      * The main difference between these settings and {@link #createDefaultDut1FitterShortTermPrediction()
  152.      * the settings for short prediction} is the much larger \(\tau\). This means weight
  153.      * is spread throughout history, hence forcing the fitted prediction model to be remain very stable
  154.      * on the long term. On the other hand, this implies that the model will start with already a much
  155.      * larger error just after raw history end.
  156.      * These settings are intended when prediction is used for 5 days after raw EOP end or more.
  157.      * </p>
  158.      * <ul>
  159.      *   <li>time constant \(\tau\) of the exponential decay set to 60 {@link Constants#JULIAN_DAY days}</li>
  160.      *   <li>convergence set to 10⁻¹² s</li>
  161.      *   <li>polynomial part set to degree 3</li>
  162.      *   <li>one harmonic term at {@link #SUN_PULSATION}}</li>
  163.      *   <li>one harmonic term at 2 times {@link #SUN_PULSATION}}</li>
  164.      *   <li>one harmonic term at 3 times {@link #SUN_PULSATION}}</li>
  165.      *   <li>one harmonic term at {@link #MOON_DRACONIC_PULSATION}}</li>
  166.      *   <li>one harmonic term at 2 times {@link #MOON_DRACONIC_PULSATION}}</li>
  167.      *   <li>one harmonic term at 3 times {@link #MOON_DRACONIC_PULSATION}}</li>
  168.      * </ul>
  169.      * @return fitter with default configuration for orientation parameters dUT1 and LOD
  170.      * @see #createDefaultDut1FitterShortTermPrediction()
  171.      */
  172.     public static SingleParameterFitter createDefaultDut1FitterLongTermPrediction() {
  173.         return new SingleParameterFitter(60 * Constants.JULIAN_DAY, 1.0e-12, 3,
  174.                                          SUN_PULSATION, 2 * SUN_PULSATION, 3 * SUN_PULSATION,
  175.                                          MOON_DRACONIC_PULSATION, 2 * MOON_DRACONIC_PULSATION, 3 * MOON_DRACONIC_PULSATION);
  176.     }

  178.     /** Create fitter with default parameters adapted for fitting pole parameters Xp and Yp
  179.      * for long term prediction.
  180.      * <p>
  181.      * The main difference between these settings and {@link #createDefaultPoleFitterLongTermPrediction()
  182.      * the settings for long prediction} is the much smaller \(\tau\). This means more
  183.      * weight is set to the points at the end of the history, hence forcing the fitted prediction
  184.      * model to be closer to these points, hence the prediction error to be smaller just after
  185.      * raw history end. On the other hand, this implies that the model will diverge on long term.
  186.      * These settings are intended when prediction is used for at most 5 days after raw EOP end.
  187.      * </p>
  188.      * <ul>
  189.      *   <li>time constant \(\tau\) of the exponential decay set to 12 {@link Constants#JULIAN_DAY days}</li>
  190.      *   <li>convergence set to 10⁻¹² rad</li>
  191.      *   <li>polynomial part set to degree 3</li>
  192.      *   <li>one harmonic term at {@link #SUN_PULSATION}}</li>
  193.      *   <li>one harmonic term at 2 times {@link #SUN_PULSATION}}</li>
  194.      *   <li>one harmonic term at 3 times {@link #SUN_PULSATION}}</li>
  195.      *   <li>one harmonic term at {@link #MOON_DRACONIC_PULSATION}}</li>
  196.      *   <li>one harmonic term at 2 times {@link #MOON_DRACONIC_PULSATION}}</li>
  197.      *   <li>one harmonic term at 3 times {@link #MOON_DRACONIC_PULSATION}}</li>
  198.      * </ul>
  199.      * @return fitter with default configuration for pole parameters Xp and Yp
  200.      */
  201.     public static SingleParameterFitter createDefaultPoleFitterShortTermPrediction() {
  202.         return new SingleParameterFitter(12 * Constants.JULIAN_DAY, 1.0e-12, 3,
  203.                                          SUN_PULSATION, 2 * SUN_PULSATION, 3 * SUN_PULSATION,
  204.                                          MOON_DRACONIC_PULSATION, 2 * MOON_DRACONIC_PULSATION, 3 * MOON_DRACONIC_PULSATION);
  205.     }

  207.     /** Create fitter with default parameters adapted for fitting pole parameters Xp and Yp
  208.      * for long term prediction.
  209.      * <p>
  210.      * The main difference between these settings and {@link #createDefaultPoleFitterShortTermPrediction()
  211.      * the settings for short prediction} is the much larger \(\tau\). This means weight
  212.      * is spread throughout history, hence forcing the fitted prediction model to be remain very stable
  213.      * on the long term. On the other hand, this implies that the model will start with already a much
  214.      * larger error just after raw history end.
  215.      * These settings are intended when prediction is used for 5 days after raw EOP end or more.
  216.      * </p>
  217.      * <ul>
  218.      *   <li>time constant \(\tau\) of the exponential decay set to 60 {@link Constants#JULIAN_DAY days}</li>
  219.      *   <li>convergence set to 10⁻¹² rad</li>
  220.      *   <li>polynomial part set to degree 3</li>
  221.      *   <li>one harmonic term at {@link #SUN_PULSATION}}</li>
  222.      *   <li>one harmonic term at 2 times {@link #SUN_PULSATION}}</li>
  223.      *   <li>one harmonic term at 3 times {@link #SUN_PULSATION}}</li>
  224.      *   <li>one harmonic term at {@link #MOON_DRACONIC_PULSATION}}</li>
  225.      *   <li>one harmonic term at 2 times {@link #MOON_DRACONIC_PULSATION}}</li>
  226.      *   <li>one harmonic term at 3 times {@link #MOON_DRACONIC_PULSATION}}</li>
  227.      * </ul>
  228.      * @return fitter with default configuration for pole parameters Xp and Yp
  229.      */
  230.     public static SingleParameterFitter createDefaultPoleFitterLongTermPrediction() {
  231.         return new SingleParameterFitter(60 * Constants.JULIAN_DAY, 1.0e-12, 3,
  232.                                          SUN_PULSATION, 2 * SUN_PULSATION, 3 * SUN_PULSATION,
  233.                                          MOON_DRACONIC_PULSATION, 2 * MOON_DRACONIC_PULSATION, 3 * MOON_DRACONIC_PULSATION);
  234.     }

  236.     /** Create fitter with default parameters adapted for fitting nutation parameters dx and dy
  237.      * for long term prediction.
  238.      * <p>
  239.      * The main difference between these settings and {@link #createDefaultNutationFitterLongTermPrediction()
  240.      * the settings for long prediction} is the much smaller \(\tau\). This means more
  241.      * weight is set to the points at the end of the history, hence forcing the fitted prediction
  242.      * model to be closer to these points, hence the prediction error to be smaller just after
  243.      * raw history end. On the other hand, this implies that the model will diverge on long term.
  244.      * These settings are intended when prediction is used for at most 5 days after raw EOP end.
  245.      * </p>
  246.      * <ul>
  247.      *   <li>time constant \(\tau\) of the exponential decay set to 12 {@link Constants#JULIAN_DAY days}</li>
  248.      *   <li>convergence set to 10⁻¹² s</li>
  249.      *   <li>polynomial part set to degree 3</li>
  250.      *   <li>one harmonic term at {@link #SUN_PULSATION}}</li>
  251.      *   <li>one harmonic term at 2 times {@link #SUN_PULSATION}}</li>
  252.      *   <li>one harmonic term at 3 times {@link #SUN_PULSATION}}</li>
  253.      *   <li>one harmonic term at {@link #MOON_DRACONIC_PULSATION}}</li>
  254.      *   <li>one harmonic term at 2 times {@link #MOON_DRACONIC_PULSATION}}</li>
  255.      *   <li>one harmonic term at 3 times {@link #MOON_DRACONIC_PULSATION}}</li>
  256.      * </ul>
  257.      * @return fitter with default configuration for pole nutation parameters dx and dy
  258.      */
  259.     public static SingleParameterFitter createDefaultNutationFitterShortTermPrediction() {
  260.         return new SingleParameterFitter(12 * Constants.JULIAN_DAY, 1.0e-12, 3,
  261.                                          SUN_PULSATION, 2 * SUN_PULSATION, 3 * SUN_PULSATION,
  262.                                          MOON_DRACONIC_PULSATION, 2 * MOON_DRACONIC_PULSATION, 3 * MOON_DRACONIC_PULSATION);
  263.     }

  265.     /** Create fitter with default parameters adapted for fitting nutation parameters dx and dy
  266.      * for long term prediction.
  267.      * <p>
  268.      * The main difference between these settings and {@link #createDefaultNutationFitterShortTermPrediction()
  269.      * the settings for short prediction} is the much larger \(\tau\). This means weight
  270.      * is spread throughout history, hence forcing the fitted prediction model to be remain very stable
  271.      * on the long term. On the other hand, this implies that the model will start with already a much
  272.      * larger error just after raw history end.
  273.      * These settings are intended when prediction is used for 5 days after raw EOP end or more.
  274.      * </p>
  275.      * <ul>
  276.      *   <li>time constant \(\tau\) of the exponential decay set to 60 {@link Constants#JULIAN_DAY days}</li>
  277.      *   <li>convergence set to 10⁻¹² s</li>
  278.      *   <li>polynomial part set to degree 3</li>
  279.      *   <li>one harmonic term at {@link #SUN_PULSATION}}</li>
  280.      *   <li>one harmonic term at 2 times {@link #SUN_PULSATION}}</li>
  281.      *   <li>one harmonic term at 3 times {@link #SUN_PULSATION}}</li>
  282.      *   <li>one harmonic term at {@link #MOON_DRACONIC_PULSATION}}</li>
  283.      *   <li>one harmonic term at 2 times {@link #MOON_DRACONIC_PULSATION}}</li>
  284.      *   <li>one harmonic term at 3 times {@link #MOON_DRACONIC_PULSATION}}</li>
  285.      * </ul>
  286.      * @return fitter with default configuration for pole nutation parameters dx and dy
  287.      */
  288.     public static SingleParameterFitter createDefaultNutationFitterLongTermPrediction() {
  289.         return new SingleParameterFitter(60 * Constants.JULIAN_DAY, 1.0e-12, 3,
  290.                                          SUN_PULSATION, 2 * SUN_PULSATION, 3 * SUN_PULSATION,
  291.                                          MOON_DRACONIC_PULSATION, 2 * MOON_DRACONIC_PULSATION, 3 * MOON_DRACONIC_PULSATION);
  292.     }

  294. }
Created with JaCoCo 0.8.12.202403310830