TLEPropagator.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.propagation.analytical.tle;

  19. import java.util.ArrayList;
  20. import java.util.Collections;
  21. import java.util.List;

  23. import org.hipparchus.geometry.euclidean.threed.Vector3D;
  24. import org.hipparchus.linear.RealMatrix;
  25. import org.hipparchus.util.FastMath;
  26. import org.hipparchus.util.MathUtils;
  27. import org.hipparchus.util.Pair;
  28. import org.hipparchus.util.SinCos;
  29. import org.orekit.annotation.DefaultDataContext;
  30. import org.orekit.attitudes.Attitude;
  31. import org.orekit.attitudes.AttitudeProvider;
  32. import org.orekit.attitudes.FrameAlignedProvider;
  33. import org.orekit.data.DataContext;
  34. import org.orekit.errors.OrekitException;
  35. import org.orekit.errors.OrekitMessages;
  36. import org.orekit.frames.Frame;
  37. import org.orekit.orbits.CartesianOrbit;
  38. import org.orekit.orbits.Orbit;
  39. import org.orekit.propagation.AbstractMatricesHarvester;
  40. import org.orekit.propagation.MatricesHarvester;
  41. import org.orekit.propagation.SpacecraftState;
  42. import org.orekit.propagation.analytical.AbstractAnalyticalPropagator;
  43. import org.orekit.propagation.analytical.tle.generation.FixedPointTleGenerationAlgorithm;
  44. import org.orekit.propagation.analytical.tle.generation.TleGenerationAlgorithm;
  45. import org.orekit.time.AbsoluteDate;
  46. import org.orekit.time.TimeScale;
  47. import org.orekit.utils.DoubleArrayDictionary;
  48. import org.orekit.utils.PVCoordinates;
  49. import org.orekit.utils.ParameterDriver;
  50. import org.orekit.utils.TimeSpanMap;
  51. import org.orekit.utils.TimeSpanMap.Span;

  53. /** This class provides elements to propagate TLE's.
  54.  * <p>
  55.  * The models used are SGP4 and SDP4, initially proposed by NORAD as the unique convenient
  56.  * propagator for TLE's. Inputs and outputs of this propagator are only suited for
  57.  * NORAD two lines elements sets, since it uses estimations and mean values appropriate
  58.  * for TLE's only.
  59.  * </p>
  60.  * <p>
  61.  * Deep- or near- space propagator is selected internally according to NORAD recommendations
  62.  * so that the user has not to worry about the used computation methods. One instance is created
  63.  * for each TLE (this instance can only be get using {@link #selectExtrapolator(TLE)} method,
  64.  * and can compute {@link PVCoordinates position and velocity coordinates} at any
  65.  * time. Maximum accuracy is guaranteed in a 24h range period before and after the provided
  66.  * TLE epoch (of course this accuracy is not really measurable nor predictable: according to
  67.  * <a href="https://www.celestrak.com/">CelesTrak</a>, the precision is close to one kilometer
  68.  * and error won't probably rise above 2 km).
  69.  * </p>
  70.  * <p>This implementation is largely inspired from the paper and source code <a
  71.  * href="https://www.celestrak.com/publications/AIAA/2006-6753/">Revisiting Spacetrack
  72.  * Report #3</a> and is fully compliant with its results and tests cases.</p>
  73.  * @author Felix R. Hoots, Ronald L. Roehrich, December 1980 (original fortran)
  74.  * @author David A. Vallado, Paul Crawford, Richard Hujsak, T.S. Kelso (C++ translation and improvements)
  75.  * @author Fabien Maussion (java translation)
  76.  * @see TLE
  77.  */
  78. public abstract class TLEPropagator extends AbstractAnalyticalPropagator {

  80.     // CHECKSTYLE: stop VisibilityModifier check

  82.     /** Initial state. */
  83.     protected TLE tle;

  85.     /** UTC time scale. */
  86.     protected final TimeScale utc;

  88.     /** final RAAN. */
  89.     protected double xnode;

  91.     /** final semi major axis. */
  92.     protected double a;

  94.     /** final eccentricity. */
  95.     protected double e;

  97.     /** final inclination. */
  98.     protected double i;

  100.     /** final perigee argument. */
  101.     protected double omega;

  103.     /** L from SPTRCK #3. */
  104.     protected double xl;

  106.     /** original recovered semi major axis. */
  107.     protected double a0dp;

  109.     /** original recovered mean motion. */
  110.     protected double xn0dp;

  112.     /** cosinus original inclination. */
  113.     protected double cosi0;

  115.     /** cos io squared. */
  116.     protected double theta2;

  118.     /** sinus original inclination. */
  119.     protected double sini0;

  121.     /** common parameter for mean anomaly (M) computation. */
  122.     protected double xmdot;

  124.     /** common parameter for perigee argument (omega) computation. */
  125.     protected double omgdot;

  127.     /** common parameter for raan (OMEGA) computation. */
  128.     protected double xnodot;

  130.     /** original eccentricity squared. */
  131.     protected double e0sq;
  132.     /** 1 - e2. */
  133.     protected double beta02;

  135.     /** sqrt (1 - e2). */
  136.     protected double beta0;

  138.     /** perigee, expressed in KM and ALTITUDE. */
  139.     protected double perige;

  141.     /** eta squared. */
  142.     protected double etasq;

  144.     /** original eccentricity * eta. */
  145.     protected double eeta;

  147.     /** s* new value for the contant s. */
  148.     protected double s4;

  150.     /** tsi from SPTRCK #3. */
  151.     protected double tsi;

  153.     /** eta from SPTRCK #3. */
  154.     protected double eta;

  156.     /** coef for SGP C3 computation. */
  157.     protected double coef;

  159.     /** coef for SGP C5 computation. */
  160.     protected double coef1;

  162.     /** C1 from SPTRCK #3. */
  163.     protected double c1;

  165.     /** C2 from SPTRCK #3. */
  166.     protected double c2;

  168.     /** C4 from SPTRCK #3. */
  169.     protected double c4;

  171.     /** common parameter for raan (OMEGA) computation. */
  172.     protected double xnodcf;

  174.     /** 3/2 * C1. */
  175.     protected double t2cof;

  177.     // CHECKSTYLE: resume VisibilityModifier check

  179.     /** TLE frame. */
  180.     private final Frame teme;

  182.     /** All TLEs and masses. */
  183.     private TimeSpanMap<Pair<TLE, Double>> tlesAndMasses;

  185.     /** Protected constructor for derived classes.
  186.      *
  187.      * <p>This constructor uses the {@link DataContext#getDefault() default data context}.
  188.      *
  189.      * @param initialTLE the unique TLE to propagate
  190.      * @param attitudeProvider provider for attitude computation
  191.      * @param mass spacecraft mass (kg)
  192.      * @see #TLEPropagator(TLE, AttitudeProvider, double, Frame)
  193.      */
  194.     @DefaultDataContext
  195.     protected TLEPropagator(final TLE initialTLE, final AttitudeProvider attitudeProvider,
  196.                             final double mass) {
  197.         this(initialTLE, attitudeProvider, mass,
  198.                 DataContext.getDefault().getFrames().getTEME());
  199.     }

  201.     /** Protected constructor for derived classes.
  202.      * @param initialTLE the unique TLE to propagate
  203.      * @param attitudeProvider provider for attitude computation
  204.      * @param mass spacecraft mass (kg)
  205.      * @param teme the TEME frame to use for propagation.
  206.      * @since 10.1
  207.      */
  208.     protected TLEPropagator(final TLE initialTLE,
  209.                             final AttitudeProvider attitudeProvider,
  210.                             final double mass,
  211.                             final Frame teme) {
  212.         super(attitudeProvider);
  213.         setStartDate(initialTLE.getDate());
  214.         this.utc       = initialTLE.getUtc();
  215.         initializeTle(initialTLE);
  216.         this.teme      = teme;
  217.         this.tlesAndMasses = new TimeSpanMap<>(new Pair<>(tle, mass));

  219.         // set the initial state
  220.         final Orbit orbit = propagateOrbit(initialTLE.getDate());
  221.         final Attitude attitude = attitudeProvider.getAttitude(orbit, orbit.getDate(), orbit.getFrame());
  222.         super.resetInitialState(new SpacecraftState(orbit, attitude).withMass(mass));
  223.     }

  225.     /** Selects the extrapolator to use with the selected TLE.
  226.      *
  227.      * <p>This method uses the {@link DataContext#getDefault() default data context}.
  228.      *
  229.      * @param tle the TLE to propagate.
  230.      * @return the correct propagator.
  231.      * @see #selectExtrapolator(TLE, Frame)
  232.      */
  233.     @DefaultDataContext
  234.     public static TLEPropagator selectExtrapolator(final TLE tle) {
  235.         return selectExtrapolator(tle, DataContext.getDefault().getFrames().getTEME());
  236.     }

  238.     /** Selects the extrapolator to use with the selected TLE.
  239.      * @param tle the TLE to propagate.
  240.      * @param teme TEME frame.
  241.      * @return the correct propagator.
  242.      * @since 10.1
  243.      * @see #selectExtrapolator(TLE, Frame, AttitudeProvider)
  244.      */
  245.     public static TLEPropagator selectExtrapolator(final TLE tle, final Frame teme) {
  246.         return selectExtrapolator(tle, teme, FrameAlignedProvider.of(teme));
  247.     }

  249.     /** Selects the extrapolator to use with the selected TLE.
  250.      * @param tle the TLE to propagate.
  251.      * @param teme TEME frame.
  252.      * @param attitudeProvider provider for attitude computation
  253.      * @return the correct propagator.
  254.      * @since 12.2
  255.      */
  256.     public static TLEPropagator selectExtrapolator(final TLE tle, final Frame teme, final AttitudeProvider attitudeProvider) {
  257.         return selectExtrapolator(tle, attitudeProvider, DEFAULT_MASS, teme);
  258.     }

  260.     /** Selects the extrapolator to use with the selected TLE.
  261.      *
  262.      * <p>This method uses the {@link DataContext#getDefault() default data context}.
  263.      *
  264.      * @param tle the TLE to propagate.
  265.      * @param attitudeProvider provider for attitude computation
  266.      * @param mass spacecraft mass (kg)
  267.      * @return the correct propagator.
  268.      * @see #selectExtrapolator(TLE, AttitudeProvider, double, Frame)
  269.      */
  270.     @DefaultDataContext
  271.     public static TLEPropagator selectExtrapolator(final TLE tle, final AttitudeProvider attitudeProvider,
  272.                                                    final double mass) {
  273.         return selectExtrapolator(tle, attitudeProvider, mass,
  274.                                   DataContext.getDefault().getFrames().getTEME());
  275.     }

  277.     /** Selects the extrapolator to use with the selected TLE.
  278.      * @param tle the TLE to propagate.
  279.      * @param attitudeProvider provider for attitude computation
  280.      * @param mass spacecraft mass (kg)
  281.      * @param teme the TEME frame to use for propagation.
  282.      * @return the correct propagator.
  283.      * @since 10.1
  284.      */
  285.     public static TLEPropagator selectExtrapolator(final TLE tle,
  286.                                                    final AttitudeProvider attitudeProvider,
  287.                                                    final double mass,
  288.                                                    final Frame teme) {

  290.         final double a1 = FastMath.pow( TLEConstants.XKE / (tle.getMeanMotion() * 60.0), TLEConstants.TWO_THIRD);
  291.         final double cosi0 = FastMath.cos(tle.getI());
  292.         final double temp = TLEConstants.CK2 * 1.5 * (3 * cosi0 * cosi0 - 1.0) *
  293.                             FastMath.pow(1.0 - tle.getE() * tle.getE(), -1.5);
  294.         final double delta1 = temp / (a1 * a1);
  295.         final double a0 = a1 * (1.0 - delta1 * (TLEConstants.ONE_THIRD + delta1 * (delta1 * 134.0 / 81.0 + 1.0)));
  296.         final double delta0 = temp / (a0 * a0);

  298.         // recover original mean motion :
  299.         final double xn0dp = tle.getMeanMotion() * 60.0 / (delta0 + 1.0);

  301.         // Period >= 225 minutes is deep space
  302.         if (MathUtils.TWO_PI / (xn0dp * TLEConstants.MINUTES_PER_DAY) >= (1.0 / 6.4)) {
  303.             return new DeepSDP4(tle, attitudeProvider, mass, teme);
  304.         } else {
  305.             return new SGP4(tle, attitudeProvider, mass, teme);
  306.         }
  307.     }

  309.     /** Get the Earth gravity coefficient used for TLE propagation.
  310.      * @return the Earth gravity coefficient.
  311.      */
  312.     public static double getMU() {
  313.         return TLEConstants.MU;
  314.     }

  316.     /** Get the extrapolated position and velocity from an initial TLE.
  317.      * @param date the final date
  318.      * @return the final PVCoordinates
  319.      */
  320.     public PVCoordinates getPVCoordinates(final AbsoluteDate date) {

  322.         sxpPropagate(date.durationFrom(tle.getDate()) / 60.0);

  324.         // Compute PV with previous calculated parameters
  325.         return computePVCoordinates();
  326.     }

  328.     /** Computation of the first commons parameters.
  329.      */
  330.     private void initializeCommons() {

  332.         // Sine and cosine of inclination
  333.         final SinCos scI0 = FastMath.sinCos(tle.getI());

  335.         final double a1 = FastMath.pow(TLEConstants.XKE / (tle.getMeanMotion() * 60.0), TLEConstants.TWO_THIRD);
  336.         cosi0 = scI0.cos();
  337.         theta2 = cosi0 * cosi0;
  338.         final double x3thm1 = 3.0 * theta2 - 1.0;
  339.         e0sq = tle.getE() * tle.getE();
  340.         beta02 = 1.0 - e0sq;
  341.         beta0 = FastMath.sqrt(beta02);
  342.         final double tval = TLEConstants.CK2 * 1.5 * x3thm1 / (beta0 * beta02);
  343.         final double delta1 = tval / (a1 * a1);
  344.         final double a0 = a1 * (1.0 - delta1 * (TLEConstants.ONE_THIRD + delta1 * (1.0 + 134.0 / 81.0 * delta1)));
  345.         final double delta0 = tval / (a0 * a0);

  347.         // recover original mean motion and semi-major axis :
  348.         xn0dp = tle.getMeanMotion() * 60.0 / (delta0 + 1.0);
  349.         a0dp = a0 / (1.0 - delta0);

  351.         // Values of s and qms2t :
  352.         s4 = TLEConstants.S;  // unmodified value for s
  353.         double q0ms24 = TLEConstants.QOMS2T; // unmodified value for q0ms2T

  355.         perige = (a0dp * (1 - tle.getE()) - TLEConstants.NORMALIZED_EQUATORIAL_RADIUS) * TLEConstants.EARTH_RADIUS; // perige

  357.         //  For perigee below 156 km, the values of s and qoms2t are changed :
  358.         if (perige < 156.0) {
  359.             if (perige <= 98.0) {
  360.                 s4 = 20.0;
  361.             } else {
  362.                 s4 = perige - 78.0;
  363.             }
  364.             final double temp_val = (120.0 - s4) * TLEConstants.NORMALIZED_EQUATORIAL_RADIUS / TLEConstants.EARTH_RADIUS;
  365.             final double temp_val_squared = temp_val * temp_val;
  366.             q0ms24 = temp_val_squared * temp_val_squared;
  367.             s4 = s4 / TLEConstants.EARTH_RADIUS + TLEConstants.NORMALIZED_EQUATORIAL_RADIUS; // new value for q0ms2T and s
  368.         }

  370.         final double pinv = 1.0 / (a0dp * beta02);
  371.         final double pinvsq = pinv * pinv;
  372.         tsi = 1.0 / (a0dp - s4);
  373.         eta = a0dp * tle.getE() * tsi;
  374.         etasq = eta * eta;
  375.         eeta = tle.getE() * eta;

  377.         final double psisq = FastMath.abs(1.0 - etasq); // abs because pow 3.5 needs positive value
  378.         final double tsi_squared = tsi * tsi;
  379.         coef = q0ms24 * tsi_squared * tsi_squared;
  380.         coef1 = coef / FastMath.pow(psisq, 3.5);

  382.         // C2 and C1 coefficients computation :
  383.         c2 = coef1 * xn0dp * (a0dp * (1.0 + 1.5 * etasq + eeta * (4.0 + etasq)) +
  384.              0.75 * TLEConstants.CK2 * tsi / psisq * x3thm1 * (8.0 + 3.0 * etasq * (8.0 + etasq)));
  385.         c1 = tle.getBStar() * c2;
  386.         sini0 = scI0.sin();

  388.         final double x1mth2 = 1.0 - theta2;

  390.         // C4 coefficient computation :
  391.         c4 = 2.0 * xn0dp * coef1 * a0dp * beta02 * (eta * (2.0 + 0.5 * etasq) +
  392.              tle.getE() * (0.5 + 2.0 * etasq) -
  393.              2 * TLEConstants.CK2 * tsi / (a0dp * psisq) *
  394.              (-3.0 * x3thm1 * (1.0 - 2.0 * eeta + etasq * (1.5 - 0.5 * eeta)) +
  395.               0.75 * x1mth2 * (2.0 * etasq - eeta * (1.0 + etasq)) * FastMath.cos(2.0 * tle.getPerigeeArgument())));

  397.         final double theta4 = theta2 * theta2;
  398.         final double temp1 = 3 * TLEConstants.CK2 * pinvsq * xn0dp;
  399.         final double temp2 = temp1 * TLEConstants.CK2 * pinvsq;
  400.         final double temp3 = 1.25 * TLEConstants.CK4 * pinvsq * pinvsq * xn0dp;

  402.         // atmospheric and gravitation coefs :(Mdf and OMEGAdf)
  403.         xmdot = xn0dp +
  404.                 0.5 * temp1 * beta0 * x3thm1 +
  405.                 0.0625 * temp2 * beta0 * (13.0 - 78.0 * theta2 + 137.0 * theta4);

  407.         final double x1m5th = 1.0 - 5.0 * theta2;

  409.         omgdot = -0.5 * temp1 * x1m5th +
  410.                  0.0625 * temp2 * (7.0 - 114.0 * theta2 + 395.0 * theta4) +
  411.                  temp3 * (3.0 - 36.0 * theta2 + 49.0 * theta4);

  413.         final double xhdot1 = -temp1 * cosi0;

  415.         xnodot = xhdot1 + (0.5 * temp2 * (4.0 - 19.0 * theta2) + 2.0 * temp3 * (3.0 - 7.0 * theta2)) * cosi0;
  416.         xnodcf = 3.5 * beta02 * xhdot1 * c1;
  417.         t2cof = 1.5 * c1;

  419.     }

  421.     /** Retrieves the position and velocity.
  422.      * @return the computed PVCoordinates.
  423.      */
  424.     private PVCoordinates computePVCoordinates() {

  426.         // Sine and cosine of final perigee argument
  427.         final SinCos scOmega = FastMath.sinCos(omega);

  429.         // Long period periodics
  430.         final double axn = e * scOmega.cos();
  431.         double temp = 1.0 / (a * (1.0 - e * e));
  432.         final double xlcof = 0.125 * TLEConstants.A3OVK2 * sini0 * (3.0 + 5.0 * cosi0) / (1.0 + cosi0);
  433.         final double aycof = 0.25 * TLEConstants.A3OVK2 * sini0;
  434.         final double xll = temp * xlcof * axn;
  435.         final double aynl = temp * aycof;
  436.         final double xlt = xl + xll;
  437.         final double ayn = e * scOmega.sin() + aynl;
  438.         final double elsq = axn * axn + ayn * ayn;
  439.         final double capu = MathUtils.normalizeAngle(xlt - xnode, FastMath.PI);
  440.         double epw = capu;
  441.         double ecosE = 0;
  442.         double esinE = 0;
  443.         double sinEPW = 0;
  444.         double cosEPW = 0;

  446.         // Dundee changes:  items dependent on cosio get recomputed:
  447.         final double cosi0Sq = cosi0 * cosi0;
  448.         final double x3thm1 = 3.0 * cosi0Sq - 1.0;
  449.         final double x1mth2 = 1.0 - cosi0Sq;
  450.         final double x7thm1 = 7.0 * cosi0Sq - 1.0;

  452.         if (e > (1 - 1e-6)) {
  453.             throw new OrekitException(OrekitMessages.TOO_LARGE_ECCENTRICITY_FOR_PROPAGATION_MODEL, e);
  454.         }

  456.         // Solve Kepler's' Equation.
  457.         final double newtonRaphsonEpsilon = 1e-12;
  458.         for (int j = 0; j < 10; j++) {

  460.             boolean doSecondOrderNewtonRaphson = true;

  462.             final SinCos scEPW = FastMath.sinCos(epw);
  463.             sinEPW = scEPW.sin();
  464.             cosEPW = scEPW.cos();
  465.             ecosE = axn * cosEPW + ayn * sinEPW;
  466.             esinE = axn * sinEPW - ayn * cosEPW;
  467.             final double f = capu - epw + esinE;
  468.             if (FastMath.abs(f) < newtonRaphsonEpsilon) {
  469.                 break;
  470.             }
  471.             final double fdot = 1.0 - ecosE;
  472.             double delta_epw = f / fdot;
  473.             if (j == 0) {
  474.                 final double maxNewtonRaphson = 1.25 * FastMath.abs(e);
  475.                 doSecondOrderNewtonRaphson = false;
  476.                 if (delta_epw > maxNewtonRaphson) {
  477.                     delta_epw = maxNewtonRaphson;
  478.                 } else if (delta_epw < -maxNewtonRaphson) {
  479.                     delta_epw = -maxNewtonRaphson;
  480.                 } else {
  481.                     doSecondOrderNewtonRaphson = true;
  482.                 }
  483.             }
  484.             if (doSecondOrderNewtonRaphson) {
  485.                 delta_epw = f / (fdot + 0.5 * esinE * delta_epw);
  486.             }
  487.             epw += delta_epw;
  488.         }

  490.         // Short period preliminary quantities
  491.         temp = 1.0 - elsq;
  492.         final double pl = a * temp;
  493.         final double r = a * (1.0 - ecosE);
  494.         double temp2 = a / r;
  495.         final double betal = FastMath.sqrt(temp);
  496.         temp = esinE / (1.0 + betal);
  497.         final double cosu = temp2 * (cosEPW - axn + ayn * temp);
  498.         final double sinu = temp2 * (sinEPW - ayn - axn * temp);
  499.         final double u = FastMath.atan2(sinu, cosu);
  500.         final double sin2u = 2.0 * sinu * cosu;
  501.         final double cos2u = 2.0 * cosu * cosu - 1.0;
  502.         final double temp1 = TLEConstants.CK2 / pl;
  503.         temp2 = temp1 / pl;

  505.         // Update for short periodics
  506.         final double rk = r * (1.0 - 1.5 * temp2 * betal * x3thm1) + 0.5 * temp1 * x1mth2 * cos2u;
  507.         final double uk = u - 0.25 * temp2 * x7thm1 * sin2u;
  508.         final double xnodek = xnode + 1.5 * temp2 * cosi0 * sin2u;
  509.         final double xinck = i + 1.5 * temp2 * cosi0 * sini0 * cos2u;

  511.         // Orientation vectors
  512.         final SinCos scuk   = FastMath.sinCos(uk);
  513.         final SinCos scik   = FastMath.sinCos(xinck);
  514.         final SinCos scnok  = FastMath.sinCos(xnodek);
  515.         final double sinuk  = scuk.sin();
  516.         final double cosuk  = scuk.cos();
  517.         final double sinik  = scik.sin();
  518.         final double cosik  = scik.cos();
  519.         final double sinnok = scnok.sin();
  520.         final double cosnok = scnok.cos();
  521.         final double xmx = -sinnok * cosik;
  522.         final double xmy = cosnok * cosik;
  523.         final double ux  = xmx * sinuk + cosnok * cosuk;
  524.         final double uy  = xmy * sinuk + sinnok * cosuk;
  525.         final double uz  = sinik * sinuk;

  527.         // Position and velocity
  528.         final double cr = 1000 * rk * TLEConstants.EARTH_RADIUS;
  529.         final Vector3D pos = new Vector3D(cr * ux, cr * uy, cr * uz);

  531.         final double rdot   = TLEConstants.XKE * FastMath.sqrt(a) * esinE / r;
  532.         final double rfdot  = TLEConstants.XKE * FastMath.sqrt(pl) / r;
  533.         final double xn     = TLEConstants.XKE / (a * FastMath.sqrt(a));
  534.         final double rdotk  = rdot - xn * temp1 * x1mth2 * sin2u;
  535.         final double rfdotk = rfdot + xn * temp1 * (x1mth2 * cos2u + 1.5 * x3thm1);
  536.         final double vx     = xmx * cosuk - cosnok * sinuk;
  537.         final double vy     = xmy * cosuk - sinnok * sinuk;
  538.         final double vz     = sinik * cosuk;

  540.         final double cv = 1000.0 * TLEConstants.EARTH_RADIUS / 60.0;
  541.         final Vector3D vel = new Vector3D(cv * (rdotk * ux + rfdotk * vx),
  542.                                           cv * (rdotk * uy + rfdotk * vy),
  543.                                           cv * (rdotk * uz + rfdotk * vz));

  545.         return new PVCoordinates(pos, vel);

  547.     }

  549.     /** Initialization proper to each propagator (SGP or SDP).
  550.      */
  551.     protected abstract void sxpInitialize();

  553.     /** Propagation proper to each propagator (SGP or SDP).
  554.      * @param t the offset from initial epoch (min)
  555.      */
  556.     protected abstract void sxpPropagate(double t);

  558.     /** {@inheritDoc}
  559.      * <p>
  560.      * For TLE propagator, calling this method is only recommended
  561.      * for covariance propagation when the new <code>state</code>
  562.      * differs from the previous one by only adding the additional
  563.      * state containing the derivatives.
  564.      * </p>
  565.      */
  566.     public void resetInitialState(final SpacecraftState state) {
  567.         super.resetInitialState(state);
  568.         resetTle(state);
  569.         tlesAndMasses = new TimeSpanMap<>(new Pair<>(tle, state.getMass()));
  570.     }

  572.     /** {@inheritDoc} */
  573.     protected void resetIntermediateState(final SpacecraftState state, final boolean forward) {
  574.         resetTle(state);
  575.         final Pair<TLE, Double> tleAndMass = new Pair<>(tle, state.getMass());
  576.         if (forward) {
  577.             tlesAndMasses.addValidAfter(tleAndMass, state.getDate(), false);
  578.         } else {
  579.             tlesAndMasses.addValidBefore(tleAndMass, state.getDate(), false);
  580.         }
  581.         stateChanged(state);
  582.     }

  584.     /** Reset internal TLE from a SpacecraftState.
  585.      * @param state spacecraft state on which to base new TLE
  586.      */
  587.     private void resetTle(final SpacecraftState state) {
  588.         final TleGenerationAlgorithm algorithm = getDefaultTleGenerationAlgorithm(utc, teme);
  589.         final TLE newTle = algorithm.generate(state, tle);
  590.         initializeTle(newTle);
  591.     }

  593.     /** Initialize internal TLE.
  594.      * @param newTle tle to replace current one
  595.      */
  596.     private void initializeTle(final TLE newTle) {
  597.         tle = newTle;
  598.         initializeCommons();
  599.         sxpInitialize();
  600.     }

  602.     /** {@inheritDoc} */
  603.     protected double getMass(final AbsoluteDate date) {
  604.         return tlesAndMasses.get(date).getValue();
  605.     }

  607.     /** {@inheritDoc} */
  608.     public Orbit propagateOrbit(final AbsoluteDate date) {
  609.         final TLE closestTle = tlesAndMasses.get(date).getKey();
  610.         if (!tle.equals(closestTle)) {
  611.             initializeTle(closestTle);
  612.         }
  613.         return new CartesianOrbit(getPVCoordinates(date), teme, date, TLEConstants.MU);
  614.     }

  616.     /** Get the underlying TLE.
  617.      * If there has been calls to #resetInitialState or #resetIntermediateState,
  618.      * it will not be the same as given to the constructor.
  619.      * @return underlying TLE
  620.      */
  621.     public TLE getTLE() {
  622.         return tle;
  623.     }

  625.     /** {@inheritDoc} */
  626.     public Frame getFrame() {
  627.         return teme;
  628.     }

  630.     /** {@inheritDoc} */
  631.     @Override
  632.     protected AbstractMatricesHarvester createHarvester(final String stmName, final RealMatrix initialStm,
  633.                                                         final DoubleArrayDictionary initialJacobianColumns) {
  634.         // Create the harvester
  635.         final TLEHarvester harvester = new TLEHarvester(this, stmName, initialStm, initialJacobianColumns);
  636.         // Update the list of additional state provider
  637.         addAdditionalDataProvider(harvester);
  638.         // Return the configured harvester
  639.         return harvester;
  640.     }

  642.     /**
  643.      * Get the names of the parameters in the matrix returned by {@link MatricesHarvester#getParametersJacobian}.
  644.      * @return names of the parameters (i.e. columns) of the Jacobian matrix
  645.      */
  646.     protected List<String> getJacobiansColumnsNames() {
  647.         final List<String> columnsNames = new ArrayList<>();
  648.         for (final ParameterDriver driver : tle.getParametersDrivers()) {

  650.             if (driver.isSelected() && !columnsNames.contains(driver.getNamesSpanMap().getFirstSpan().getData())) {
  651.                 // As driver with same name should have same NamesSpanMap we only check the if condition on the
  652.                 // first span map and then if the condition is OK all the span names are added to the jacobian column names
  653.                 for (Span<String> span = driver.getNamesSpanMap().getFirstSpan(); span != null; span = span.next()) {
  654.                     columnsNames.add(span.getData());
  655.                 }
  656.             }
  657.         }
  658.         Collections.sort(columnsNames);
  659.         return columnsNames;
  660.     }

  662.     /**
  663.      * Get the default TLE generation algorithm.
  664.      * @param utc UTC time scale
  665.      * @param teme TEME frame
  666.      * @return a TLE generation algorithm
  667.      * @since 12.0
  668.      */
  669.     public static TleGenerationAlgorithm getDefaultTleGenerationAlgorithm(final TimeScale utc, final Frame teme) {
  670.         return new FixedPointTleGenerationAlgorithm(FixedPointTleGenerationAlgorithm.EPSILON_DEFAULT,
  671.                                                     FixedPointTleGenerationAlgorithm.MAX_ITERATIONS_DEFAULT,
  672.                                                     FixedPointTleGenerationAlgorithm.SCALE_DEFAULT, utc, teme);
  673.     }

  675. }
Created with JaCoCo 0.8.12.202403310830