TimeDependentHarmonic.java

  1. /* Copyright 2002-2025 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;

  19. import org.hipparchus.util.SinCos;

  21. /** Time-dependent part of a single component of spherical harmonics.
  22.  * @author Luc Maisonobe
  23.  * @since 11.1
  24.  */
  25. class TimeDependentHarmonic {

  27.     /** Index of the trend reference in the gravity field. */
  28.     private final int trendReferenceIndex;

  30.     /** Base of the cosine coefficient. */
  31.     private final double cBase;

  33.     /** Base of the sine coefficient. */
  34.     private final double sBase;

  36.     /** Secular trend of the cosine coefficient. */
  37.     private double cTrend;

  39.     /** Secular trend of the sine coefficient. */
  40.     private double sTrend;

  42.     /** Indices of the reference dates in the gravity field. */
  43.     private int[] cosReferenceIndices;

  45.     /** Indices of the harmonic pulsations in the gravity field. */
  46.     private int[] cosPulsationIndices;

  48.     /** Cosine component of the cosine coefficient. */
  49.     private double[] cosC;

  51.     /** Cosine component of the sine coefficient. */
  52.     private double[] cosS;

  54.     /** Indices of the reference dates in the gravity field. */
  55.     private int[] sinReferenceIndices;

  57.     /** Indices of the harmonic pulsations in the gravity field. */
  58.     private int[] sinPulsationIndices;

  60.     /** Sine component of the cosine coefficient. */
  61.     private double[] sinC;

  63.     /** Sine component of the sine coefficient. */
  64.     private double[] sinS;

  66.     /** Build a part with only base.
  67.      * @param trendReferenceIndex index of the trend reference in the gravity field
  68.      * @param cBase base of the cosine coefficient
  69.      * @param sBase base of the sine coefficient
  70.      */
  71.     TimeDependentHarmonic(final int trendReferenceIndex, final double cBase, final double sBase) {
  72.         this(trendReferenceIndex, cBase, sBase, 0, 0);
  73.     }

  75.     /** Build a rescaled component.
  76.      * @param scale scaling factor to apply to all coefficients elements
  77.      * @param original original component
  78.      */
  79.     TimeDependentHarmonic(final double scale, final TimeDependentHarmonic original) {

  81.         // rescale base
  82.         this(original.trendReferenceIndex, scale * original.cBase, scale * original.sBase,
  83.              original.cosReferenceIndices.length, original.sinReferenceIndices.length);

  85.         // rescale trend
  86.         cTrend = scale * original.cTrend;
  87.         sTrend = scale * original.sTrend;

  89.         // rescale cosine
  90.         for (int i = 0; i < cosReferenceIndices.length; ++i) {
  91.             cosReferenceIndices[i] = original.cosReferenceIndices[i];
  92.             cosPulsationIndices[i] = original.cosPulsationIndices[i];
  93.             cosC[i]                = scale * original.cosC[i];
  94.             cosS[i]                = scale * original.cosS[i];
  95.         }

  97.         // rescale sine
  98.         for (int i = 0; i < sinReferenceIndices.length; ++i) {
  99.             sinReferenceIndices[i] = original.sinReferenceIndices[i];
  100.             sinPulsationIndices[i] = original.sinPulsationIndices[i];
  101.             sinC[i]                = scale * original.sinC[i];
  102.             sinS[i]                = scale * original.sinS[i];
  103.         }

  105.     }

  107.     /** Build a part with only base.
  108.      * @param trendReferenceIndex index of the trend reference in the gravity field
  109.      * @param cBase base of the cosine coefficient
  110.      * @param sBase base of the sine coefficient
  111.      * @param cSize initial size of the cosine arrays
  112.      * @param sSize initial size of the sine arrays
  113.      */
  114.     private TimeDependentHarmonic(final int trendReferenceIndex, final double cBase, final double sBase,
  115.                                   final int cSize, final int sSize) {

  117.         // linear part
  118.         this.trendReferenceIndex = trendReferenceIndex;
  119.         this.cBase               = cBase;
  120.         this.sBase               = sBase;
  121.         this.cTrend              = 0.0;
  122.         this.sTrend              = 0.0;

  124.         // cosine component
  125.         this.cosReferenceIndices = new int[cSize];
  126.         this.cosPulsationIndices = new int[cSize];
  127.         this.cosC                = new double[cSize];
  128.         this.cosS                = new double[cSize];

  130.         // sine component
  131.         this.sinReferenceIndices = new int[sSize];
  132.         this.sinPulsationIndices = new int[sSize];
  133.         this.sinC                = new double[sSize];
  134.         this.sinS                = new double[sSize];

  136.     }

  138.     /** Set the trend part.
  139.      * @param cDot secular trend of the cosine coefficient (s⁻¹)
  140.      * @param sDot secular trend of the sine coefficient (s⁻¹)
  141.      */
  142.     public void setTrend(final double cDot, final double sDot) {
  143.         this.cTrend = cDot;
  144.         this.sTrend = sDot;
  145.     }

  147.     /** Add a cosine component.
  148.      * @param cosReferenceIndex index of the reference date in the gravity field
  149.      * (if negative, use the trend reference index)
  150.      * @param cosPulsationIndex index of the harmonic pulsation in the gravity field
  151.      * @param cosineC cosine component of the cosine coefficient
  152.      * @param cosineS cosine component of the sine coefficient
  153.      */
  154.     public void addCosine(final int cosReferenceIndex, final int cosPulsationIndex,
  155.                           final double cosineC, final double cosineS) {
  156.         final int refIndex = cosReferenceIndex < 0 ? trendReferenceIndex : cosReferenceIndex;
  157.         this.cosReferenceIndices = addInt(refIndex, this.cosReferenceIndices);
  158.         this.cosPulsationIndices = addInt(cosPulsationIndex, this.cosPulsationIndices);
  159.         this.cosC                = addDouble(cosineC,     this.cosC);
  160.         this.cosS                = addDouble(cosineS,     this.cosS);
  161.     }

  163.     /** Add a sine component.
  164.      * @param sinReferenceIndex index of the reference date in the gravity field
  165.      * (if negative, use the trend reference index)
  166.      * @param sinPulsationIndex index of the harmonic pulsation in the gravity field
  167.      * @param sineC sine component of the cosine coefficient
  168.      * @param sineS sine component of the sine coefficient
  169.      */
  170.     public void addSine(final int sinReferenceIndex, final int sinPulsationIndex,
  171.                         final double sineC, final double sineS) {
  172.         final int refIndex = sinReferenceIndex < 0 ? trendReferenceIndex : sinReferenceIndex;
  173.         this.sinReferenceIndices = addInt(refIndex, this.sinReferenceIndices);
  174.         this.sinPulsationIndices = addInt(sinPulsationIndex, this.sinPulsationIndices);
  175.         this.sinC                = addDouble(sineC,       this.sinC);
  176.         this.sinS                = addDouble(sineS,       this.sinS);
  177.     }

  179.     /** Add an integer to an array.
  180.      * <p>
  181.      * Expanding the array one element at a time may seem a waste of time,
  182.      * but we expect the array to be 0, 1 or 2 elements long only, and this
  183.      * if performed only when reading gravity field, so its is worth doing
  184.      * it this way.
  185.      * </p>
  186.      * @param n integer to add
  187.      * @param array array where to add the integer
  188.      * @return new array
  189.      */
  190.     private static int[] addInt(final int n, final int[] array) {
  191.         final int[] newArray = new int[array.length + 1];
  192.         System.arraycopy(array, 0, newArray, 0, array.length);
  193.         newArray[array.length] = n;
  194.         return newArray;
  195.     }

  197.     /** Add a double number to an array.
  198.      * <p>
  199.      * Expanding the array one element at a time may seem a waste of time,
  200.      * but we expect the array to be 0, 1 or 2 elements long only, and this
  201.      * if performed only when reading gravity field, so its is worth doing
  202.      * it this way.
  203.      * </p>
  204.      * @param d double number to add
  205.      * @param array array where to add the double number
  206.      * @return new array
  207.      */
  208.     private static double[] addDouble(final double d, final double[] array) {
  209.         final double[] newArray = new double[array.length + 1];
  210.         System.arraycopy(array, 0, newArray, 0, array.length);
  211.         newArray[array.length] = d;
  212.         return newArray;
  213.     }

  215.     /** Compute the time-dependent part of a spherical harmonic cosine coefficient.
  216.      * @param offsets offsets to reference dates in the gravity field
  217.      * @param pulsations angular pulsations in the gravity field
  218.      * @return raw coefficient Cnm
  219.      */
  220.     public double computeCnm(final double[] offsets, final SinCos[][] pulsations) {

  222.         // trend effect
  223.         double cnm = cBase + offsets[trendReferenceIndex] * cTrend;

  225.         for (int i = 0; i < cosPulsationIndices.length; ++i) {
  226.             // cosine effect
  227.             cnm += cosC[i] * pulsations[cosReferenceIndices[i]][cosPulsationIndices[i]].cos();
  228.         }

  230.         for (int i = 0; i < sinPulsationIndices.length; ++i) {
  231.             // sine effect
  232.             cnm += sinC[i] * pulsations[sinReferenceIndices[i]][sinPulsationIndices[i]].sin();
  233.         }

  235.         return cnm;

  237.     }

  239.     /** Compute the time-dependent part of a spherical harmonic sine coefficient.
  240.      * @param offsets offsets to reference dates in the gravity field
  241.      * @param pulsations angular pulsations in the gravity field
  242.      * @return raw coefficient Snm
  243.      */
  244.     public double computeSnm(final double[] offsets, final SinCos[][] pulsations) {

  246.         // trend effect
  247.         double snm = sBase + offsets[trendReferenceIndex] * sTrend;

  249.         for (int i = 0; i < cosPulsationIndices.length; ++i) {
  250.             // cosine effect
  251.             snm += cosS[i] * pulsations[cosReferenceIndices[i]][cosPulsationIndices[i]].cos();
  252.         }

  254.         for (int i = 0; i < sinPulsationIndices.length; ++i) {
  255.             // sine effect
  256.             snm += sinS[i] * pulsations[sinReferenceIndices[i]][sinPulsationIndices[i]].sin();
  257.         }

  259.         return snm;

  261.     }

  263. }
Created with JaCoCo 0.8.12.202403310830