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;
18
19 import java.io.Serializable;
20 import java.util.List;
21 import java.util.ListIterator;
22 import java.util.function.ToDoubleFunction;
23
24 import org.hipparchus.util.FastMath;
25 import org.hipparchus.util.MathUtils;
26 import org.orekit.utils.Constants;
27 import org.orekit.utils.SecularAndHarmonic;
28
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 {
37
38 /** Sun pulsation, one year period. */
39 public static final double SUN_PULSATION = MathUtils.TWO_PI / Constants.JULIAN_YEAR;
40
41 /** Moon pulsation (one Moon draconic period). */
42 public static final double MOON_DRACONIC_PULSATION = MathUtils.TWO_PI / (27.212221 * Constants.JULIAN_DAY);
43
44 /** Serializable UID. */
45 private static final long serialVersionUID = 20230309L;
46
47 /** Time constant of the exponential decay weight. */
48 private final double timeConstant;
49
50 /** Convergence on fitted parameter. */
51 private final double convergence;
52
53 /** Degree of the polynomial model. */
54 private final int degree;
55
56 /** Pulsations of harmonic part (rad/s). */
57 private final double[] pulsations;
58
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 }
81
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) {
88
89 final List<EOPEntry> rawEntries = rawHistory.getEntries();
90 final EOPEntry last = rawEntries.get(rawEntries.size() - 1);
91
92 // create fitter
93 final SecularAndHarmonic sh = new SecularAndHarmonic(degree, pulsations);
94
95 // set up convergence
96 sh.setConvergenceRMS(convergence);
97
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);
102
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 }
110
111 // perform fitting
112 sh.fit();
113
114 return sh;
115
116 }
117
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 }
147
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 }
177
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 }
206
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 }
235
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 }
264
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 }
293
294 }