1 /* Copyright 2002-2022 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.forces.gravity.potential;
18
19 import java.io.IOException;
20 import java.io.InputStream;
21 import java.text.ParseException;
22 import java.util.ArrayList;
23 import java.util.Arrays;
24 import java.util.List;
25 import java.util.Locale;
26
27 import org.hipparchus.util.FastMath;
28 import org.hipparchus.util.Precision;
29 import org.orekit.annotation.DefaultDataContext;
30 import org.orekit.data.DataContext;
31 import org.orekit.data.DataLoader;
32 import org.orekit.errors.OrekitException;
33 import org.orekit.errors.OrekitMessages;
34 import org.orekit.time.AbsoluteDate;
35 import org.orekit.time.DateComponents;
36 import org.orekit.time.TimeComponents;
37 import org.orekit.time.TimeScale;
38
39 /**This abstract class represents a Gravitational Potential Coefficients file reader.
40 *
41 * <p> As it exits many different coefficients models and containers this
42 * interface represents all the methods that should be implemented by a reader.
43 * The proper way to use this interface is to call the {@link GravityFieldFactory}
44 * which will determine which reader to use with the selected potential
45 * coefficients file.<p>
46 *
47 * @see GravityFields
48 * @author Fabien Maussion
49 */
50 public abstract class PotentialCoefficientsReader implements DataLoader {
51
52 /** Maximal degree to parse. */
53 private int maxParseDegree;
54
55 /** Maximal order to parse. */
56 private int maxParseOrder;
57
58 /** Regular expression for supported files names. */
59 private final String supportedNames;
60
61 /** Allow missing coefficients in the input data. */
62 private final boolean missingCoefficientsAllowed;
63
64 /** Time scale for parsing dates. */
65 private final TimeScale timeScale;
66
67 /** Indicator for complete read. */
68 private boolean readComplete;
69
70 /** Central body reference radius. */
71 private double ae;
72
73 /** Central body attraction coefficient. */
74 private double mu;
75
76 /** Converter from triangular to flat form. */
77 private Flattener flattener;
78
79 /** Raw tesseral-sectorial coefficients matrix. */
80 private double[] rawC;
81
82 /** Raw tesseral-sectorial coefficients matrix. */
83 private double[] rawS;
84
85 /** Indicator for normalized raw coefficients. */
86 private boolean normalized;
87
88 /** Tide system. */
89 private TideSystem tideSystem;
90
91 /** Simple constructor.
92 * <p>Build an uninitialized reader.</p>
93 *
94 * <p>This constructor uses the {@link DataContext#getDefault() default data context}.
95 *
96 * @param supportedNames regular expression for supported files names
97 * @param missingCoefficientsAllowed allow missing coefficients in the input data
98 * @see #PotentialCoefficientsReader(String, boolean, TimeScale)
99 */
100 @DefaultDataContext
101 protected PotentialCoefficientsReader(final String supportedNames,
102 final boolean missingCoefficientsAllowed) {
103 this(supportedNames, missingCoefficientsAllowed,
104 DataContext.getDefault().getTimeScales().getTT());
105 }
106
107 /** Simple constructor.
108 * <p>Build an uninitialized reader.</p>
109 * @param supportedNames regular expression for supported files names
110 * @param missingCoefficientsAllowed allow missing coefficients in the input data
111 * @param timeScale to use when parsing dates.
112 * @since 10.1
113 */
114 protected PotentialCoefficientsReader(final String supportedNames,
115 final boolean missingCoefficientsAllowed,
116 final TimeScale timeScale) {
117 this.supportedNames = supportedNames;
118 this.missingCoefficientsAllowed = missingCoefficientsAllowed;
119 this.maxParseDegree = Integer.MAX_VALUE;
120 this.maxParseOrder = Integer.MAX_VALUE;
121 this.readComplete = false;
122 this.ae = Double.NaN;
123 this.mu = Double.NaN;
124 this.flattener = null;
125 this.rawC = null;
126 this.rawS = null;
127 this.normalized = false;
128 this.tideSystem = TideSystem.UNKNOWN;
129 this.timeScale = timeScale;
130 }
131
132 /** Get the regular expression for supported files names.
133 * @return regular expression for supported files names
134 */
135 public String getSupportedNames() {
136 return supportedNames;
137 }
138
139 /** Check if missing coefficients are allowed in the input data.
140 * @return true if missing coefficients are allowed in the input data
141 */
142 public boolean missingCoefficientsAllowed() {
143 return missingCoefficientsAllowed;
144 }
145
146 /** Set the degree limit for the next file parsing.
147 * @param maxParseDegree maximal degree to parse (may be safely
148 * set to {@link Integer#MAX_VALUE} to parse all available coefficients)
149 * @since 6.0
150 */
151 public void setMaxParseDegree(final int maxParseDegree) {
152 this.maxParseDegree = maxParseDegree;
153 }
154
155 /** Get the degree limit for the next file parsing.
156 * @return degree limit for the next file parsing
157 * @since 6.0
158 */
159 public int getMaxParseDegree() {
160 return maxParseDegree;
161 }
162
163 /** Set the order limit for the next file parsing.
164 * @param maxParseOrder maximal order to parse (may be safely
165 * set to {@link Integer#MAX_VALUE} to parse all available coefficients)
166 * @since 6.0
167 */
168 public void setMaxParseOrder(final int maxParseOrder) {
169 this.maxParseOrder = maxParseOrder;
170 }
171
172 /** Get the order limit for the next file parsing.
173 * @return order limit for the next file parsing
174 * @since 6.0
175 */
176 public int getMaxParseOrder() {
177 return maxParseOrder;
178 }
179
180 /** {@inheritDoc} */
181 public boolean stillAcceptsData() {
182 return !(readComplete &&
183 getMaxAvailableDegree() >= getMaxParseDegree() &&
184 getMaxAvailableOrder() >= getMaxParseOrder());
185 }
186
187 /** Set the indicator for completed read.
188 * @param readComplete if true, a gravity field has been completely read
189 */
190 protected void setReadComplete(final boolean readComplete) {
191 this.readComplete = readComplete;
192 }
193
194 /** Set the central body reference radius.
195 * @param ae central body reference radius
196 */
197 protected void setAe(final double ae) {
198 this.ae = ae;
199 }
200
201 /** Get the central body reference radius.
202 * @return central body reference radius
203 */
204 protected double getAe() {
205 return ae;
206 }
207
208 /** Set the central body attraction coefficient.
209 * @param mu central body attraction coefficient
210 */
211 protected void setMu(final double mu) {
212 this.mu = mu;
213 }
214
215 /** Get the central body attraction coefficient.
216 * @return central body attraction coefficient
217 */
218 protected double getMu() {
219 return mu;
220 }
221
222 /** Set the {@link TideSystem} used in the gravity field.
223 * @param tideSystem tide system used in the gravity field
224 */
225 protected void setTideSystem(final TideSystem tideSystem) {
226 this.tideSystem = tideSystem;
227 }
228
229 /** Get the {@link TideSystem} used in the gravity field.
230 * @return tide system used in the gravity field
231 */
232 protected TideSystem getTideSystem() {
233 return tideSystem;
234 }
235
236 /** Set the tesseral-sectorial coefficients matrix.
237 * @param rawNormalized if true, raw coefficients are normalized
238 * @param c raw tesseral-sectorial coefficients matrix
239 * @param s raw tesseral-sectorial coefficients matrix
240 * @param name name of the file (or zip entry)
241 * @deprecated as of 11.1, replaced by {@link #setRawCoefficients(boolean,
242 * Flattener, double[], double[], String)}
243 */
244 @Deprecated
245 protected void setRawCoefficients(final boolean rawNormalized,
246 final double[][] c, final double[][] s,
247 final String name) {
248 setRawCoefficients(rawNormalized, buildFlattener(c),
249 buildFlattener(c).flatten(c), buildFlattener(s).flatten(s),
250 name);
251 }
252
253 /** Set the tesseral-sectorial coefficients matrix.
254 * @param rawNormalized if true, raw coefficients are normalized
255 * @param f converter from triangular to flat form
256 * @param c raw tesseral-sectorial coefficients matrix
257 * @param s raw tesseral-sectorial coefficients matrix
258 * @param name name of the file (or zip entry)
259 * @since 11.1
260 */
261 protected void setRawCoefficients(final boolean rawNormalized, final Flattener f,
262 final double[] c, final double[] s,
263 final String name) {
264
265 this.flattener = f;
266
267 // normalization indicator
268 normalized = rawNormalized;
269
270 // set known constant values, if they were not defined in the file.
271 // See Hofmann-Wellenhof and Moritz, "Physical Geodesy",
272 // section 2.6 Harmonics of Lower Degree.
273 // All S_i,0 are irrelevant because they are multiplied by zero.
274 // C0,0 is 1, the central part, since all coefficients are normalized by GM.
275 setIfUnset(c, flattener.index(0, 0), 1);
276 setIfUnset(s, flattener.index(0, 0), 0);
277
278 if (flattener.getDegree() >= 1) {
279 // C1,0, C1,1, and S1,1 are the x,y,z coordinates of the center of mass,
280 // which are 0 since all coefficients are given in an Earth centered frame
281 setIfUnset(c, flattener.index(1, 0), 0);
282 setIfUnset(s, flattener.index(1, 0), 0);
283 if (flattener.getOrder() >= 1) {
284 setIfUnset(c, flattener.index(1, 1), 0);
285 setIfUnset(s, flattener.index(1, 1), 0);
286 }
287 }
288
289 // cosine part
290 for (int i = 0; i <= flattener.getDegree(); ++i) {
291 for (int j = 0; j <= FastMath.min(i, flattener.getOrder()); ++j) {
292 if (Double.isNaN(c[flattener.index(i, j)])) {
293 throw new OrekitException(OrekitMessages.MISSING_GRAVITY_FIELD_COEFFICIENT_IN_FILE,
294 'C', i, j, name);
295 }
296 }
297 }
298 rawC = c.clone();
299
300 // sine part
301 for (int i = 0; i <= flattener.getDegree(); ++i) {
302 for (int j = 0; j <= FastMath.min(i, flattener.getOrder()); ++j) {
303 if (Double.isNaN(s[flattener.index(i, j)])) {
304 throw new OrekitException(OrekitMessages.MISSING_GRAVITY_FIELD_COEFFICIENT_IN_FILE,
305 'S', i, j, name);
306 }
307 }
308 }
309 rawS = s.clone();
310
311 }
312
313 /**
314 * Set a coefficient if it has not been set already.
315 * <p>
316 * If {@code array[i]} is 0 or NaN this method sets it to {@code value} and returns
317 * {@code true}. Otherwise the original value of {@code array[i]} is preserved and
318 * this method return {@code false}.
319 * <p>
320 * If {@code array[i]} does not exist then this method returns {@code false}.
321 *
322 * @param array the coefficient array.
323 * @param i index in array.
324 * @param value the new value to set.
325 * @return {@code true} if the coefficient was set to {@code value}, {@code false} if
326 * the coefficient was not set to {@code value}. A {@code false} return indicates the
327 * coefficient has previously been set to a non-NaN, non-zero value.
328 */
329 private boolean setIfUnset(final double[] array, final int i, final double value) {
330 if (array.length > i && (Double.isNaN(array[i]) || Precision.equals(array[i], 0.0, 0))) {
331 // the coefficient was not already initialized
332 array[i] = value;
333 return true;
334 } else {
335 return false;
336 }
337 }
338
339 /** Get the maximal degree available in the last file parsed.
340 * @return maximal degree available in the last file parsed
341 * @since 6.0
342 */
343 public int getMaxAvailableDegree() {
344 return flattener.getDegree();
345 }
346
347 /** Get the maximal order available in the last file parsed.
348 * @return maximal order available in the last file parsed
349 * @since 6.0
350 */
351 public int getMaxAvailableOrder() {
352 return flattener.getOrder();
353 }
354
355 /** {@inheritDoc} */
356 public abstract void loadData(InputStream input, String name)
357 throws IOException, ParseException, OrekitException;
358
359 /** Get a provider for read spherical harmonics coefficients.
360 * @param wantNormalized if true, the provider will provide normalized coefficients,
361 * otherwise it will provide un-normalized coefficients
362 * @param degree maximal degree
363 * @param order maximal order
364 * @return a new provider
365 * @see #getConstantProvider(boolean, int, int)
366 * @since 6.0
367 */
368 public abstract RawSphericalHarmonicsProvider getProvider(boolean wantNormalized, int degree, int order);
369
370 /** Get a time-independent provider containing base harmonics coefficients.
371 * <p>
372 * Beware that some coeefficients may be missing here, if they are managed as time-dependent
373 * piecewise models (as in ICGEM V2.0).
374 * </p>
375 * @param wantNormalized if true, the raw provider must provide normalized coefficients,
376 * otherwise it will provide un-normalized coefficients
377 * @param degree maximal degree
378 * @param order maximal order
379 * @return a new provider, with no time-dependent parts
380 * @see #getProvider(boolean, int, int)
381 * @since 11.1
382 */
383 protected ConstantSphericalHarmonics getBaseProvider(final boolean wantNormalized,
384 final int degree, final int order) {
385
386 if (!readComplete) {
387 throw new OrekitException(OrekitMessages.NO_GRAVITY_FIELD_DATA_LOADED);
388 }
389
390 final Flattener truncatedFlattener = new Flattener(degree, order);
391 return new ConstantSphericalHarmonics(ae, mu, tideSystem, truncatedFlattener,
392 rescale(1.0, wantNormalized, truncatedFlattener, flattener, rawC),
393 rescale(1.0, wantNormalized, truncatedFlattener, flattener, rawS));
394
395 }
396
397 /** Get a time-independent provider for read spherical harmonics coefficients.
398 * @param wantNormalized if true, the raw provider must provide normalized coefficients,
399 * otherwise it will provide un-normalized coefficients
400 * @param degree maximal degree
401 * @param order maximal order
402 * @return a new provider, with no time-dependent parts
403 * @see #getProvider(boolean, int, int)
404 * @since 6.0
405 * @deprecated as of 11.1, not used anymore
406 */
407 @Deprecated
408 protected ConstantSphericalHarmonics getConstantProvider(final boolean wantNormalized,
409 final int degree, final int order) {
410 return getBaseProvider(wantNormalized, degree, order);
411 }
412
413 /** Get a flattener for a triangular array.
414 * @param triangular triangular array to flatten
415 * @return flattener suited for triangular array dimensions
416 * @since 11.1
417 */
418 private static Flattener buildFlattener(final double[][] triangular) {
419 return new Flattener(triangular.length - 1, triangular[triangular.length - 1].length - 1);
420 }
421
422 /** Build a coefficients triangular array.
423 * @param degree array degree
424 * @param order array order
425 * @param value initial value to put in array elements
426 * @return built array
427 * @deprecated as of 11.1, replaced by {@link #buildFlatArray(Flattener, double)}
428 */
429 @Deprecated
430 protected static double[][] buildTriangularArray(final int degree, final int order, final double value) {
431 final int rows = degree + 1;
432 final double[][] array = new double[rows][];
433 for (int k = 0; k < array.length; ++k) {
434 array[k] = buildRow(k, order, value);
435 }
436 return array;
437 }
438
439 /** Build a coefficients array in flat form.
440 * @param flattener converter from triangular to flat form
441 * @param value initial value to put in array elements
442 * @return built array
443 * @since 11.1
444 */
445 protected static double[] buildFlatArray(final Flattener flattener, final double value) {
446 final double[] array = new double[flattener.arraySize()];
447 Arrays.fill(array, value);
448 return array;
449 }
450
451 /**
452 * Parse a double from a string. Accept the Fortran convention of using a 'D' or
453 * 'd' instead of an 'E' or 'e'.
454 *
455 * @param string to be parsed.
456 * @return the double value of {@code string}.
457 */
458 protected static double parseDouble(final String string) {
459 return Double.parseDouble(string.toUpperCase(Locale.ENGLISH).replace('D', 'E'));
460 }
461
462 /** Build a coefficients row.
463 * @param degree row degree
464 * @param order row order
465 * @param value initial value to put in array elements
466 * @return built row
467 */
468 protected static double[] buildRow(final int degree, final int order, final double value) {
469 final double[] row = new double[FastMath.min(order, degree) + 1];
470 Arrays.fill(row, value);
471 return row;
472 }
473
474 /** Extend a list of lists of coefficients if needed.
475 * @param list list of lists of coefficients
476 * @param degree degree required to be present
477 * @param order order required to be present
478 * @param value initial value to put in list elements
479 * @deprecated as of 11.1, not used anymore
480 */
481 @Deprecated
482 protected void extendListOfLists(final List<List<Double>> list, final int degree, final int order,
483 final double value) {
484 for (int i = list.size(); i <= degree; ++i) {
485 list.add(new ArrayList<>());
486 }
487 final List<Double> listN = list.get(degree);
488 final Double v = value;
489 for (int j = listN.size(); j <= order; ++j) {
490 listN.add(v);
491 }
492 }
493
494 /** Convert a list of list into an array.
495 * @param list list of lists of coefficients
496 * @return a new array
497 * @deprecated as of 11.1, not used anymore
498 */
499 @Deprecated
500 protected double[][] toArray(final List<List<Double>> list) {
501 final double[][] array = new double[list.size()][];
502 for (int i = 0; i < array.length; ++i) {
503 array[i] = new double[list.get(i).size()];
504 for (int j = 0; j < array[i].length; ++j) {
505 array[i][j] = list.get(i).get(j);
506 }
507 }
508 return array;
509 }
510
511 /** Parse a coefficient.
512 * @param field text field to parse
513 * @param list list where to put the coefficient
514 * @param i first index in the list
515 * @param j second index in the list
516 * @param cName name of the coefficient
517 * @param name name of the file
518 * @deprecated as of 11.1, replaced by {@link #parseCoefficient(String,
519 * Flattener, double[], int, int, String, String)}
520 */
521 @Deprecated
522 protected void parseCoefficient(final String field, final List<List<Double>> list,
523 final int i, final int j,
524 final String cName, final String name) {
525 final double value = parseDouble(field);
526 final double oldValue = list.get(i).get(j);
527 if (Double.isNaN(oldValue) || Precision.equals(oldValue, 0.0, 0)) {
528 // the coefficient was not already initialized
529 list.get(i).set(j, value);
530 } else {
531 throw new OrekitException(OrekitMessages.DUPLICATED_GRAVITY_FIELD_COEFFICIENT_IN_FILE,
532 name, i, j, name);
533 }
534 }
535
536 /** Parse a coefficient.
537 * @param field text field to parse
538 * @param array array where to put the coefficient
539 * @param i first index in the list
540 * @param j second index in the list
541 * @param cName name of the coefficient
542 * @param name name of the file
543 * @deprecated as of 11.1, replaced by {@link #parseCoefficient(String,
544 * Flattener, double[], int, int, String, String)}
545 */
546 @Deprecated
547 protected void parseCoefficient(final String field, final double[][] array,
548 final int i, final int j,
549 final String cName, final String name) {
550 final double value = parseDouble(field);
551 final double oldValue = array[i][j];
552 if (Double.isNaN(oldValue) || Precision.equals(oldValue, 0.0, 0)) {
553 // the coefficient was not already initialized
554 array[i][j] = value;
555 } else {
556 throw new OrekitException(OrekitMessages.DUPLICATED_GRAVITY_FIELD_COEFFICIENT_IN_FILE,
557 name, i, j, name);
558 }
559 }
560
561 /** Parse a coefficient.
562 * @param field text field to parse
563 * @param f converter from triangular to flat form
564 * @param array array where to put the coefficient
565 * @param i first index in the list
566 * @param j second index in the list
567 * @param cName name of the coefficient
568 * @param name name of the file
569 * @since 11.1
570 */
571 protected void parseCoefficient(final String field, final Flattener f,
572 final double[] array, final int i, final int j,
573 final String cName, final String name) {
574 final int index = f.index(i, j);
575 final double value = parseDouble(field);
576 final double oldValue = array[index];
577 if (Double.isNaN(oldValue) || Precision.equals(oldValue, 0.0, 0)) {
578 // the coefficient was not already initialized
579 array[index] = value;
580 } else {
581 throw new OrekitException(OrekitMessages.DUPLICATED_GRAVITY_FIELD_COEFFICIENT_IN_FILE,
582 name, i, j, name);
583 }
584 }
585
586 /** Rescale coefficients arrays.
587 * @param scale general scaling factor to apply to all elements
588 * @param normalizedOrigin if true, the origin coefficients are normalized
589 * @param originC cosine part of the original coefficients
590 * @param originS sine part of the origin coefficients
591 * @param wantNormalized if true, the rescaled coefficients must be normalized
592 * @param rescaledC cosine part of the rescaled coefficients to fill in (may be the originC array)
593 * @param rescaledS sine part of the rescaled coefficients to fill in (may be the originS array)
594 * @deprecated as of 11.1, replaced by {@link #rescale(double, boolean, Flattener, Flattener, double[])}
595 */
596 @Deprecated
597 protected static void rescale(final double scale,
598 final boolean normalizedOrigin, final double[][] originC,
599 final double[][] originS, final boolean wantNormalized,
600 final double[][] rescaledC, final double[][] rescaledS) {
601
602 if (wantNormalized == normalizedOrigin) {
603 // apply only the general scaling factor
604 for (int i = 0; i < rescaledC.length; ++i) {
605 final double[] rCi = rescaledC[i];
606 final double[] rSi = rescaledS[i];
607 final double[] oCi = originC[i];
608 final double[] oSi = originS[i];
609 for (int j = 0; j < rCi.length; ++j) {
610 rCi[j] = oCi[j] * scale;
611 rSi[j] = oSi[j] * scale;
612 }
613 }
614 } else {
615
616 // we have to re-scale the coefficients
617 // (we use rescaledC.length - 1 for the order instead of rescaledC[rescaledC.length - 1].length - 1
618 // because typically trend or pulsation arrays are irregular, some test cases have
619 // order 2 elements at degree 2, but only order 1 elements for higher degrees for example)
620 final double[][] factors = GravityFieldFactory.getUnnormalizationFactors(rescaledC.length - 1,
621 rescaledC.length - 1);
622
623 if (wantNormalized) {
624 // normalize the coefficients
625 for (int i = 0; i < rescaledC.length; ++i) {
626 final double[] rCi = rescaledC[i];
627 final double[] rSi = rescaledS[i];
628 final double[] oCi = originC[i];
629 final double[] oSi = originS[i];
630 final double[] fi = factors[i];
631 for (int j = 0; j < rCi.length; ++j) {
632 final double factor = scale / fi[j];
633 rCi[j] = oCi[j] * factor;
634 rSi[j] = oSi[j] * factor;
635 }
636 }
637 } else {
638 // un-normalize the coefficients
639 for (int i = 0; i < rescaledC.length; ++i) {
640 final double[] rCi = rescaledC[i];
641 final double[] rSi = rescaledS[i];
642 final double[] oCi = originC[i];
643 final double[] oSi = originS[i];
644 final double[] fi = factors[i];
645 for (int j = 0; j < rCi.length; ++j) {
646 final double factor = scale * fi[j];
647 rCi[j] = oCi[j] * factor;
648 rSi[j] = oSi[j] * factor;
649 }
650 }
651 }
652
653 }
654 }
655
656 /** Rescale coefficients arrays.
657 * <p>
658 * The normalized/unnormalized nature of original coefficients is inherited from previous parsing.
659 * </p>
660 * @param scale general scaling factor to apply to all elements
661 * @param wantNormalized if true, the rescaled coefficients must be normalized,
662 * otherwise they must be un-normalized
663 * @param rescaledFlattener converter from triangular to flat form
664 * @param originalFlattener converter from triangular to flat form
665 * @param original original coefficients
666 * @return rescaled coefficients
667 * @since 11.1
668 */
669 protected double[] rescale(final double scale, final boolean wantNormalized, final Flattener rescaledFlattener,
670 final Flattener originalFlattener, final double[] original) {
671
672 if (rescaledFlattener.getDegree() > originalFlattener.getDegree()) {
673 throw new OrekitException(OrekitMessages.TOO_LARGE_DEGREE_FOR_GRAVITY_FIELD,
674 rescaledFlattener.getDegree(), flattener.getDegree());
675 }
676
677 if (rescaledFlattener.getOrder() > originalFlattener.getOrder()) {
678 throw new OrekitException(OrekitMessages.TOO_LARGE_ORDER_FOR_GRAVITY_FIELD,
679 rescaledFlattener.getOrder(), flattener.getOrder());
680 }
681
682 // scaling and normalization factors
683 final FactorsGenerator generator;
684 if (wantNormalized == normalized) {
685 // the parsed coefficients already match the specified normalization
686 generator = (n, m) -> scale;
687 } else {
688 // we need to normalize/unnormalize parsed coefficients
689 final double[][] unnormalizationFactors =
690 GravityFieldFactory.getUnnormalizationFactors(rescaledFlattener.getDegree(),
691 rescaledFlattener.getOrder());
692 generator = wantNormalized ?
693 (n, m) -> scale / unnormalizationFactors[n][m] :
694 (n, m) -> scale * unnormalizationFactors[n][m];
695 }
696
697 // perform rescaling
698 final double[] rescaled = buildFlatArray(rescaledFlattener, 0.0);
699 for (int n = 0; n <= rescaledFlattener.getDegree(); ++n) {
700 for (int m = 0; m <= FastMath.min(n, rescaledFlattener.getOrder()); ++m) {
701 final int rescaledndex = rescaledFlattener.index(n, m);
702 final int originalndex = originalFlattener.index(n, m);
703 rescaled[rescaledndex] = original[originalndex] * generator.factor(n, m);
704 }
705 }
706
707 return rescaled;
708
709 }
710
711 /** Rescale coefficients arrays.
712 * <p>
713 * The normalized/unnormalized nature of original coefficients is inherited from previous parsing.
714 * </p>
715 * @param wantNormalized if true, the rescaled coefficients must be normalized,
716 * otherwise they must be un-normalized
717 * @param rescaledFlattener converter from triangular to flat form
718 * @param originalFlattener converter from triangular to flat form
719 * @param original original coefficients
720 * @return rescaled coefficients
721 * @since 11.1
722 */
723 protected TimeDependentHarmonic[] rescale(final boolean wantNormalized, final Flattener rescaledFlattener,
724 final Flattener originalFlattener, final TimeDependentHarmonic[] original) {
725
726 if (rescaledFlattener.getDegree() > originalFlattener.getDegree()) {
727 throw new OrekitException(OrekitMessages.TOO_LARGE_DEGREE_FOR_GRAVITY_FIELD,
728 rescaledFlattener.getDegree(), flattener.getDegree());
729 }
730
731 if (rescaledFlattener.getOrder() > originalFlattener.getOrder()) {
732 throw new OrekitException(OrekitMessages.TOO_LARGE_ORDER_FOR_GRAVITY_FIELD,
733 rescaledFlattener.getOrder(), flattener.getOrder());
734 }
735
736 // scaling and normalization factors
737 final FactorsGenerator generator;
738 if (wantNormalized == normalized) {
739 // the parsed coefficients already match the specified normalization
740 generator = (n, m) -> 1.0;
741 } else {
742 // we need to normalize/unnormalize parsed coefficients
743 final double[][] unnormalizationFactors =
744 GravityFieldFactory.getUnnormalizationFactors(rescaledFlattener.getDegree(),
745 rescaledFlattener.getOrder());
746 generator = wantNormalized ?
747 (n, m) -> 1.0 / unnormalizationFactors[n][m] :
748 (n, m) -> unnormalizationFactors[n][m];
749 }
750
751 // perform rescaling
752 final TimeDependentHarmonic[] rescaled = new TimeDependentHarmonic[rescaledFlattener.arraySize()];
753 for (int n = 0; n <= rescaledFlattener.getDegree(); ++n) {
754 for (int m = 0; m <= FastMath.min(n, rescaledFlattener.getOrder()); ++m) {
755 final int originalndex = originalFlattener.index(n, m);
756 if (original[originalndex] != null) {
757 final int rescaledndex = rescaledFlattener.index(n, m);
758 final double factor = generator.factor(n, m);
759 rescaled[rescaledndex] = new TimeDependentHarmonic(factor, original[originalndex]);
760 }
761 }
762 }
763
764 return rescaled;
765
766 }
767
768 /**
769 * Create a date from components. Assumes the time part is noon.
770 *
771 * @param components year, month, day.
772 * @return date.
773 */
774 protected AbsoluteDate toDate(final DateComponents components) {
775 return toDate(components, TimeComponents.H12);
776 }
777
778 /**
779 * Create a date from components.
780 *
781 * @param dc dates components.
782 * @param tc time components
783 * @return date.
784 * @since 11.1
785 */
786 protected AbsoluteDate toDate(final DateComponents dc, final TimeComponents tc) {
787 return new AbsoluteDate(dc, tc, timeScale);
788 }
789
790 /** Generator for normalization/unnormalization factors.
791 * @since 11.1
792 */
793 private interface FactorsGenerator {
794 /** Generator the normalization/unnormalization factors.
795 * @param n degree of the gravity field component
796 * @param m order of the gravity field component
797 * @return factor to apply to term
798 */
799 double factor(int n, int m);
800 }
801
802 }