FieldEquinoctialLongitudeArgumentUtility.java

  1. /* Copyright 2022-2024 Romain Serra
  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.orbits;

  18. import org.hipparchus.CalculusFieldElement;
  19. import org.hipparchus.util.FastMath;
  20. import org.hipparchus.util.FieldSinCos;
  21. import org.orekit.errors.OrekitException;
  22. import org.orekit.errors.OrekitMessages;

  23. /**
  24.  * Utility methods for converting between different longitude arguments used by {@link FieldEquinoctialOrbit}.
  25.  * @author Romain Serra
  26.  * @see FieldEquinoctialOrbit
  27.  * @since 12.1
  28.  */
  29. public class FieldEquinoctialLongitudeArgumentUtility {

  30.     /** Tolerance for stopping criterion in iterative conversion from mean to eccentric angle. */
  31.     private static final double TOLERANCE_CONVERGENCE = 1.0e-12;

  32.     /** Maximum number of iterations in iterative conversion from mean to eccentric angle. */
  33.     private static final int MAXIMUM_ITERATION = 50;

  34.     /** Private constructor for utility class. */
  35.     private FieldEquinoctialLongitudeArgumentUtility() {
  36.         // nothing here (utils class)
  37.     }

  38.     /**
  39.      * Computes the true longitude argument from the eccentric longitude argument.
  40.      *
  41.      * @param <T> Type of the field elements
  42.      * @param ex  e cos(ω), first component of eccentricity vector
  43.      * @param ey  e sin(ω), second component of eccentricity vector
  44.      * @param lE  = E + ω + Ω eccentric longitude argument (rad)
  45.      * @return the true longitude argument.
  46.      */
  47.     public static <T extends CalculusFieldElement<T>> T eccentricToTrue(final T ex, final T ey, final T lE) {
  48.         final T epsilon           = eccentricAndTrueEpsilon(ex, ey);
  49.         final FieldSinCos<T> scLE = FastMath.sinCos(lE);
  50.         final T cosLE             = scLE.cos();
  51.         final T sinLE             = scLE.sin();
  52.         final T num               = ex.multiply(sinLE).subtract(ey.multiply(cosLE));
  53.         final T den               = epsilon.add(1).subtract(ex.multiply(cosLE)).subtract(ey.multiply(sinLE));
  54.         return lE.add(eccentricAndTrueAtan(num, den));
  55.     }

  56.     /**
  57.      * Computes the eccentric longitude argument from the true longitude argument.
  58.      *
  59.      * @param <T> Type of the field elements
  60.      * @param ex  e cos(ω), first component of eccentricity vector
  61.      * @param ey  e sin(ω), second component of eccentricity vector
  62.      * @param lV  = V + ω + Ω true longitude argument (rad)
  63.      * @return the eccentric longitude argument.
  64.      */
  65.     public static <T extends CalculusFieldElement<T>> T trueToEccentric(final T ex, final T ey, final T lV) {
  66.         final T epsilon           = eccentricAndTrueEpsilon(ex, ey);
  67.         final FieldSinCos<T> scLv = FastMath.sinCos(lV);
  68.         final T cosLv             = scLv.cos();
  69.         final T sinLv             = scLv.sin();
  70.         final T num               = ey.multiply(cosLv).subtract(ex.multiply(sinLv));
  71.         final T den               = epsilon.add(1).add(ex.multiply(cosLv)).add(ey.multiply(sinLv));
  72.         return lV.add(eccentricAndTrueAtan(num, den));
  73.     }

  74.     /**
  75.      * Computes an intermediate quantity for conversions between true and eccentric.
  76.      *
  77.      * @param <T>    Type of the field elements
  78.      * @param ex e cos(ω), first component of eccentricity vector
  79.      * @param ey e sin(ω), second component of eccentricity vector
  80.      * @return intermediate variable referred to as epsilon.
  81.      */
  82.     private static <T extends CalculusFieldElement<T>> T eccentricAndTrueEpsilon(final T ex, final T ey) {
  83.         return (ex.square().negate().subtract(ey.square()).add(1.)).sqrt();
  84.     }

  85.     /**
  86.      * Computes another intermediate quantity for conversions between true and eccentric.
  87.      *
  88.      * @param <T>    Type of the field elements
  89.      * @param num numerator for angular conversion
  90.      * @param den denominator for angular conversion
  91.      * @return arc-tangent of ratio of inputs times two.
  92.      */
  93.     private static <T extends CalculusFieldElement<T>> T eccentricAndTrueAtan(final T num, final T den) {
  94.         return (num.divide(den)).atan().multiply(2);
  95.     }

  96.     /**
  97.      * Computes the eccentric longitude argument from the mean longitude argument.
  98.      *
  99.      * @param <T> Type of the field elements
  100.      * @param ex  e cos(ω), first component of eccentricity vector
  101.      * @param ey  e sin(ω), second component of eccentricity vector
  102.      * @param lM  = M + ω + Ω  mean longitude argument (rad)
  103.      * @return the eccentric longitude argument.
  104.      */
  105.     public static <T extends CalculusFieldElement<T>> T meanToEccentric(final T ex, final T ey, final T lM) {
  106.         // Generalization of Kepler equation to equinoctial parameters
  107.         // with lE = PA + RAAN + E and
  108.         //      lM = PA + RAAN + M = lE - ex.sin(lE) + ey.cos(lE)
  109.         T lE = lM;
  110.         T shift;
  111.         T lEmlM = lM.getField().getZero();
  112.         boolean hasConverged;
  113.         int iter = 0;
  114.         do {
  115.             final FieldSinCos<T> scLE = FastMath.sinCos(lE);
  116.             final T f2 = ex.multiply(scLE.sin()).subtract(ey.multiply(scLE.cos()));
  117.             final T f1 = ex.multiply(scLE.cos()).add(ey.multiply(scLE.sin())).negate().add(1);
  118.             final T f0 = lEmlM.subtract(f2);

  119.             final T f12 = f1.multiply(2.0);
  120.             shift = f0.multiply(f12).divide(f1.multiply(f12).subtract(f0.multiply(f2)));

  121.             lEmlM = lEmlM.subtract(shift);
  122.             lE     = lM.add(lEmlM);

  123.             hasConverged = FastMath.abs(shift.getReal()) <= TOLERANCE_CONVERGENCE;
  124.         } while (++iter < MAXIMUM_ITERATION && !hasConverged);

  125.         if (!hasConverged) {
  126.             throw new OrekitException(OrekitMessages.UNABLE_TO_COMPUTE_ECCENTRIC_LONGITUDE_ARGUMENT, iter);
  127.         }
  128.         return lE;

  129.     }

  130.     /**
  131.      * Computes the mean longitude argument from the eccentric longitude argument.
  132.      *
  133.      * @param <T> Type of the field elements
  134.      * @param ex  e cos(ω), first component of eccentricity vector
  135.      * @param ey  e sin(ω), second component of eccentricity vector
  136.      * @param lE  = E + ω + Ω  mean longitude argument (rad)
  137.      * @return the mean longitude argument.
  138.      */
  139.     public static <T extends CalculusFieldElement<T>> T eccentricToMean(final T ex, final T ey, final T lE) {
  140.         final FieldSinCos<T> scLE = FastMath.sinCos(lE);
  141.         return lE.subtract(ex.multiply(scLE.sin())).add(ey.multiply(scLE.cos()));
  142.     }

  143.     /**
  144.      * Computes the mean longitude argument from the eccentric longitude argument.
  145.      *
  146.      * @param <T> Type of the field elements
  147.      * @param ex  e cos(ω), first component of eccentricity vector
  148.      * @param ey  e sin(ω), second component of eccentricity vector
  149.      * @param lV  = V + ω + Ω  true longitude argument (rad)
  150.      * @return the mean longitude argument.
  151.      */
  152.     public static <T extends CalculusFieldElement<T>> T trueToMean(final T ex, final T ey, final T lV) {
  153.         final T alphaE = trueToEccentric(ex, ey, lV);
  154.         return eccentricToMean(ex, ey, alphaE);
  155.     }

  156.     /**
  157.      * Computes the true longitude argument from the eccentric longitude argument.
  158.      *
  159.      * @param <T> Type of the field elements
  160.      * @param ex  e cos(ω), first component of eccentricity vector
  161.      * @param ey  e sin(ω), second component of eccentricity vector
  162.      * @param lM  = M + ω + Ω  mean longitude argument (rad)
  163.      * @return the true longitude argument.
  164.      */
  165.     public static <T extends CalculusFieldElement<T>> T meanToTrue(final T ex, final T ey, final T lM) {
  166.         final T alphaE = meanToEccentric(ex, ey, lM);
  167.         return eccentricToTrue(ex, ey, alphaE);
  168.     }

  169. }