FieldIntelsatElevenElementsPropagator.java

  1. /* Copyright 2002-2025 Airbus Defence and Space
  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.  * Airbus Defence and Space 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.intelsat;

  18. import java.util.Collections;
  19. import java.util.List;
  20. import org.hipparchus.CalculusFieldElement;
  21. import org.hipparchus.Field;
  22. import org.hipparchus.analysis.differentiation.FieldUnivariateDerivative2;
  23. import org.hipparchus.geometry.euclidean.threed.FieldVector3D;
  24. import org.hipparchus.util.FastMath;
  25. import org.hipparchus.util.FieldSinCos;
  26. import org.orekit.annotation.DefaultDataContext;
  27. import org.orekit.attitudes.AttitudeProvider;
  28. import org.orekit.attitudes.FieldAttitude;
  29. import org.orekit.attitudes.FrameAlignedProvider;
  30. import org.orekit.data.DataContext;
  31. import org.orekit.errors.OrekitException;
  32. import org.orekit.errors.OrekitMessages;
  33. import org.orekit.frames.Frame;
  34. import org.orekit.frames.FramesFactory;
  35. import org.orekit.orbits.FieldCartesianOrbit;
  36. import org.orekit.orbits.FieldOrbit;
  37. import org.orekit.propagation.FieldSpacecraftState;
  38. import org.orekit.propagation.Propagator;
  39. import org.orekit.propagation.analytical.FieldAbstractAnalyticalPropagator;
  40. import org.orekit.time.FieldAbsoluteDate;
  41. import org.orekit.utils.Constants;
  42. import org.orekit.utils.FieldPVCoordinates;
  43. import org.orekit.utils.IERSConventions;
  44. import org.orekit.utils.ParameterDriver;
  45. import org.orekit.utils.units.Unit;

  46. /**
  47.  * This class provides elements to propagate Intelsat's 11 elements.
  48.  * <p>
  49.  * Intelsat's 11 elements propagation is defined in ITU-R S.1525 standard.
  50.  * </p>
  51.  *
  52.  * @author Bryan Cazabonne
  53.  * @since 12.1
  54.  */
  55. public class FieldIntelsatElevenElementsPropagator<T extends CalculusFieldElement<T>> extends FieldAbstractAnalyticalPropagator<T> {

  56.     /**
  57.      * Intelsat's 11 elements.
  58.      */
  59.     private final FieldIntelsatElevenElements<T> elements;

  60.     /**
  61.      * Inertial frame for the output orbit.
  62.      */
  63.     private final Frame inertialFrame;

  64.     /**
  65.      * ECEF frame related to the Intelsat's 11 elements.
  66.      */
  67.     private final Frame ecefFrame;

  68.     /**
  69.      * Spacecraft mass in kilograms.
  70.      */
  71.     private final T mass;

  72.     /**
  73.      * Compute spacecraft's east longitude.
  74.      */
  75.     private FieldUnivariateDerivative2<T> eastLongitudeDegrees;

  76.     /**
  77.      * Compute spacecraft's geocentric latitude.
  78.      */
  79.     private FieldUnivariateDerivative2<T> geocentricLatitudeDegrees;

  80.     /**
  81.      * Compute spacecraft's orbit radius.
  82.      */
  83.     private FieldUnivariateDerivative2<T> orbitRadius;

  84.     /**
  85.      * Default constructor.
  86.      * <p>
  87.      * This constructor uses the {@link DataContext#getDefault() default data context}.
  88.      * </p>
  89.      * <p> The attitude provider is set by default to be aligned with the inertial frame.<br>
  90.      * The mass is set by default to the
  91.      * {@link org.orekit.propagation.Propagator#DEFAULT_MASS DEFAULT_MASS}.<br>
  92.      * The inertial frame is set by default to the
  93.      * {@link org.orekit.frames.Predefined#TOD_CONVENTIONS_2010_SIMPLE_EOP TOD frame} in the default data
  94.      * context.<br>
  95.      * The ECEF frame is set by default to the
  96.      * {@link org.orekit.frames.Predefined#ITRF_CIO_CONV_2010_SIMPLE_EOP
  97.      * CIO/2010-based ITRF simple EOP} in the default data context.
  98.      * </p>
  99.      *
  100.      * @param elements Intelsat's 11 elements
  101.      */
  102.     @DefaultDataContext
  103.     public FieldIntelsatElevenElementsPropagator(final FieldIntelsatElevenElements<T> elements) {
  104.         this(elements, FramesFactory.getTOD(IERSConventions.IERS_2010, true), FramesFactory.getITRF(IERSConventions.IERS_2010, true));
  105.     }

  106.     /**
  107.      * Constructor.
  108.      *
  109.      * <p> The attitude provider is set by default to be aligned with the inertial frame.<br>
  110.      * The mass is set by default to the
  111.      * {@link org.orekit.propagation.Propagator#DEFAULT_MASS DEFAULT_MASS}.<br>
  112.      * </p>
  113.      *
  114.      * @param elements      Intelsat's 11 elements
  115.      * @param inertialFrame inertial frame for the output orbit
  116.      * @param ecefFrame     ECEF frame related to the Intelsat's 11 elements
  117.      */
  118.     public FieldIntelsatElevenElementsPropagator(final FieldIntelsatElevenElements<T> elements, final Frame inertialFrame, final Frame ecefFrame) {
  119.         this(elements, inertialFrame, ecefFrame, FrameAlignedProvider.of(inertialFrame), elements.getEpoch().getField().getZero().add(Propagator.DEFAULT_MASS));
  120.     }

  121.     /**
  122.      * Constructor.
  123.      *
  124.      * @param elements         Intelsat's 11 elements
  125.      * @param inertialFrame    inertial frame for the output orbit
  126.      * @param ecefFrame        ECEF frame related to the Intelsat's 11 elements
  127.      * @param attitudeProvider attitude provider
  128.      * @param mass             spacecraft mass
  129.      */
  130.     public FieldIntelsatElevenElementsPropagator(final FieldIntelsatElevenElements<T> elements, final Frame inertialFrame, final Frame ecefFrame,
  131.                                                  final AttitudeProvider attitudeProvider, final T mass) {
  132.         super(elements.getEpoch().getField(), attitudeProvider);
  133.         this.elements = elements;
  134.         this.inertialFrame = inertialFrame;
  135.         this.ecefFrame = ecefFrame;
  136.         this.mass = mass;
  137.         setStartDate(elements.getEpoch());
  138.         final FieldOrbit<T> orbitAtElementsDate = propagateOrbit(elements.getEpoch(), getParameters(elements.getEpoch().getField()));
  139.         final FieldAttitude<T> attitude = attitudeProvider.getAttitude(orbitAtElementsDate, elements.getEpoch(), inertialFrame);
  140.         super.resetInitialState(new FieldSpacecraftState<>(orbitAtElementsDate, attitude).withMass(mass));
  141.     }

  142.     /**
  143.      * Converts the Intelsat's 11 elements into Position/Velocity coordinates in ECEF.
  144.      *
  145.      * @param date computation epoch
  146.      * @return Position/Velocity coordinates in ECEF
  147.      */
  148.     public FieldPVCoordinates<T> propagateInEcef(final FieldAbsoluteDate<T> date) {
  149.         final Field<T> field = date.getField();
  150.         final FieldUnivariateDerivative2<T> tDays = new FieldUnivariateDerivative2<>(date.durationFrom(elements.getEpoch()), field.getOne(), field.getZero()).divide(
  151.                 Constants.JULIAN_DAY);
  152.         final T wDegreesPerDay = elements.getLm1().add(IntelsatElevenElements.DRIFT_RATE_SHIFT_DEG_PER_DAY);
  153.         final FieldUnivariateDerivative2<T> wt = FastMath.toRadians(tDays.multiply(wDegreesPerDay));
  154.         final FieldSinCos<FieldUnivariateDerivative2<T>> scWt = FastMath.sinCos(wt);
  155.         final FieldSinCos<FieldUnivariateDerivative2<T>> sc2Wt = FastMath.sinCos(wt.multiply(2.0));
  156.         final FieldUnivariateDerivative2<T> satelliteEastLongitudeDegrees = computeSatelliteEastLongitudeDegrees(tDays, scWt, sc2Wt);
  157.         final FieldUnivariateDerivative2<T> satelliteGeocentricLatitudeDegrees = computeSatelliteGeocentricLatitudeDegrees(tDays, scWt);
  158.         final FieldUnivariateDerivative2<T> satelliteRadius = computeSatelliteRadiusKilometers(wDegreesPerDay, scWt).multiply(Unit.KILOMETRE.getScale());
  159.         this.eastLongitudeDegrees = satelliteEastLongitudeDegrees;
  160.         this.geocentricLatitudeDegrees = satelliteGeocentricLatitudeDegrees;
  161.         this.orbitRadius = satelliteRadius;
  162.         final FieldSinCos<FieldUnivariateDerivative2<T>> scLongitude = FastMath.sinCos(FastMath.toRadians(satelliteEastLongitudeDegrees));
  163.         final FieldSinCos<FieldUnivariateDerivative2<T>> scLatitude = FastMath.sinCos(FastMath.toRadians(satelliteGeocentricLatitudeDegrees));
  164.         final FieldVector3D<FieldUnivariateDerivative2<T>> positionWithDerivatives = new FieldVector3D<>(satelliteRadius.multiply(scLatitude.cos()).multiply(scLongitude.cos()),
  165.                                                                                                          satelliteRadius.multiply(scLatitude.cos()).multiply(scLongitude.sin()),
  166.                                                                                                          satelliteRadius.multiply(scLatitude.sin()));
  167.         return new FieldPVCoordinates<>(new FieldVector3D<>(positionWithDerivatives.getX().getValue(), //
  168.                                                             positionWithDerivatives.getY().getValue(), //
  169.                                                             positionWithDerivatives.getZ().getValue()), //
  170.                                         new FieldVector3D<>(positionWithDerivatives.getX().getFirstDerivative(), //
  171.                                                             positionWithDerivatives.getY().getFirstDerivative(), //
  172.                                                             positionWithDerivatives.getZ().getFirstDerivative()), //
  173.                                         new FieldVector3D<>(positionWithDerivatives.getX().getSecondDerivative(), //
  174.                                                             positionWithDerivatives.getY().getSecondDerivative(), //
  175.                                                             positionWithDerivatives.getZ().getSecondDerivative()));
  176.     }

  177.     /**
  178.      * {@inheritDoc}.
  179.      */
  180.     @Override
  181.     public void resetInitialState(final FieldSpacecraftState<T> state) {
  182.         throw new OrekitException(OrekitMessages.NON_RESETABLE_STATE);
  183.     }

  184.     /**
  185.      * {@inheritDoc}.
  186.      */
  187.     @Override
  188.     protected T getMass(final FieldAbsoluteDate<T> date) {
  189.         return mass;
  190.     }

  191.     /**
  192.      * {@inheritDoc}.
  193.      */
  194.     @Override
  195.     protected void resetIntermediateState(final FieldSpacecraftState<T> state, final boolean forward) {
  196.         throw new OrekitException(OrekitMessages.NON_RESETABLE_STATE);
  197.     }

  198.     /**
  199.      * {@inheritDoc}.
  200.      */
  201.     @Override
  202.     public FieldOrbit<T> propagateOrbit(final FieldAbsoluteDate<T> date,
  203.                                         final T[] parameters) {
  204.         return new FieldCartesianOrbit<>(ecefFrame.getTransformTo(inertialFrame, date).transformPVCoordinates(propagateInEcef(date)), inertialFrame, date,
  205.                                          date.getField().getZero().add(Constants.WGS84_EARTH_MU));
  206.     }

  207.     /**
  208.      * Computes the satellite's east longitude.
  209.      *
  210.      * @param tDays delta time in days
  211.      * @param scW   sin/cos of the W angle
  212.      * @param sc2W  sin/cos of the 2xW angle
  213.      * @return the satellite's east longitude in degrees
  214.      */
  215.     private FieldUnivariateDerivative2<T> computeSatelliteEastLongitudeDegrees(final FieldUnivariateDerivative2<T> tDays, final FieldSinCos<FieldUnivariateDerivative2<T>> scW,
  216.                                                                                final FieldSinCos<FieldUnivariateDerivative2<T>> sc2W) {
  217.         final FieldUnivariateDerivative2<T> longitude = tDays.multiply(tDays).multiply(elements.getLm2()) //
  218.                                                              .add(tDays.multiply(elements.getLm1())) //
  219.                                                              .add(elements.getLm0());
  220.         final FieldUnivariateDerivative2<T> cosineLongitudeTerm = scW.cos().multiply(tDays.multiply(elements.getLonC1()).add(elements.getLonC()));
  221.         final FieldUnivariateDerivative2<T> sineLongitudeTerm = scW.sin().multiply(tDays.multiply(elements.getLonS1()).add(elements.getLonS()));
  222.         final FieldUnivariateDerivative2<T> latitudeTerm = sc2W.sin()
  223.                                                                .multiply(elements.getLatC()
  224.                                                                                  .multiply(elements.getLatC())
  225.                                                                                  .subtract(elements.getLatS().multiply(elements.getLatS()))
  226.                                                                                  .multiply(0.5)) //
  227.                                                                .subtract(sc2W.cos().multiply(elements.getLatC().multiply(elements.getLatS()))) //
  228.                                                                .multiply(IntelsatElevenElements.K);
  229.         return longitude.add(cosineLongitudeTerm).add(sineLongitudeTerm).add(latitudeTerm);
  230.     }

  231.     /**
  232.      * Computes the satellite's geocentric latitude.
  233.      *
  234.      * @param tDays delta time in days
  235.      * @param scW   sin/cos of the W angle
  236.      * @return he satellite geocentric latitude in degrees
  237.      */
  238.     private FieldUnivariateDerivative2<T> computeSatelliteGeocentricLatitudeDegrees(final FieldUnivariateDerivative2<T> tDays,
  239.                                                                                     final FieldSinCos<FieldUnivariateDerivative2<T>> scW) {
  240.         final FieldUnivariateDerivative2<T> cosineTerm = scW.cos().multiply(tDays.multiply(elements.getLatC1()).add(elements.getLatC()));
  241.         final FieldUnivariateDerivative2<T> sineTerm = scW.sin().multiply(tDays.multiply(elements.getLatS1()).add(elements.getLatS()));
  242.         return cosineTerm.add(sineTerm);
  243.     }

  244.     /**
  245.      * Computes the satellite's orbit radius.
  246.      *
  247.      * @param wDegreesPerDay W angle in degrees/day
  248.      * @param scW            sin/cos of the W angle
  249.      * @return the satellite's orbit radius in kilometers
  250.      */
  251.     private FieldUnivariateDerivative2<T> computeSatelliteRadiusKilometers(final T wDegreesPerDay, final FieldSinCos<FieldUnivariateDerivative2<T>> scW) {
  252.         final T coefficient = elements.getLm1()
  253.                                       .multiply(2.0)
  254.                                       .divide(wDegreesPerDay.subtract(elements.getLm1()).multiply(3.0))
  255.                                       .negate()
  256.                                       .add(1.0)
  257.                                       .multiply(IntelsatElevenElements.SYNCHRONOUS_RADIUS_KM);
  258.         return scW.sin()
  259.                   .multiply(elements.getLonC().multiply(IntelsatElevenElements.K))
  260.                   .add(1.0)
  261.                   .subtract(scW.cos().multiply(elements.getLonS().multiply(IntelsatElevenElements.K)))
  262.                   .multiply(coefficient);
  263.     }

  264.     /**
  265.      * Get the computed satellite's east longitude.
  266.      *
  267.      * @return the satellite's east longitude in degrees
  268.      */
  269.     public FieldUnivariateDerivative2<T> getEastLongitudeDegrees() {
  270.         return eastLongitudeDegrees;
  271.     }

  272.     /**
  273.      * Get the computed satellite's geocentric latitude.
  274.      *
  275.      * @return the satellite's geocentric latitude in degrees
  276.      */
  277.     public FieldUnivariateDerivative2<T> getGeocentricLatitudeDegrees() {
  278.         return geocentricLatitudeDegrees;
  279.     }

  280.     /**
  281.      * Get the computed satellite's orbit.
  282.      *
  283.      * @return satellite's orbit radius in meters
  284.      */
  285.     public FieldUnivariateDerivative2<T> getOrbitRadius() {
  286.         return orbitRadius;
  287.     }

  288.     /**
  289.      * {@inheritDoc}.
  290.      */
  291.     @Override
  292.     public Frame getFrame() {
  293.         return inertialFrame;
  294.     }

  295.     /**
  296.      * {@inheritDoc}.
  297.      */
  298.     @Override
  299.     public List<ParameterDriver> getParametersDrivers() {
  300.         return Collections.emptyList();
  301.     }

  302.     /**
  303.      * Get the Intelsat's 11 elements used by the propagator.
  304.      *
  305.      * @return the Intelsat's 11 elements used by the propagator
  306.      */
  307.     public FieldIntelsatElevenElements<T> getIntelsatElevenElements() {
  308.         return elements;
  309.     }
  310. }