AuxiliaryElements.java

  1. /* Copyright 2002-2024 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.semianalytical.dsst.utilities;

  18. import org.hipparchus.geometry.euclidean.threed.Vector3D;
  19. import org.hipparchus.util.FastMath;
  20. import org.hipparchus.util.MathUtils;
  21. import org.orekit.frames.Frame;
  22. import org.orekit.orbits.Orbit;
  23. import org.orekit.time.AbsoluteDate;


  24. /** Container class for common parameters used by all DSST forces.
  25.  *  <p>
  26.  *  Most of them are defined in Danielson paper at § 2.1.
  27.  *  </p>
  28.  *  @author Pascal Parraud
  29.  */
  30. public class AuxiliaryElements {

  31.     /** Orbit date. */
  32.     private final AbsoluteDate date;

  33.     /** Orbit frame. */
  34.     private final Frame frame;

  35.     /** Eccentricity. */
  36.     private final double ecc;

  37.     /** Keplerian mean motion. */
  38.     private final double n;

  39.     /** Keplerian period. */
  40.     private final double period;

  41.     /** Semi-major axis. */
  42.     private final double sma;

  43.     /** x component of eccentricity vector. */
  44.     private final double k;

  45.     /** y component of eccentricity vector. */
  46.     private final double h;

  47.     /** x component of inclination vector. */
  48.     private final double q;

  49.     /** y component of inclination vector. */
  50.     private final double p;

  51.     /** Mean longitude. */
  52.     private final double lm;

  53.     /** True longitude. */
  54.     private final double lv;

  55.     /** Eccentric longitude. */
  56.     private final double le;

  57.     /** Retrograde factor I.
  58.      *  <p>
  59.      *  DSST model needs equinoctial orbit as internal representation.
  60.      *  Classical equinoctial elements have discontinuities when inclination
  61.      *  is close to zero. In this representation, I = +1. <br>
  62.      *  To avoid this discontinuity, another representation exists and equinoctial
  63.      *  elements can be expressed in a different way, called "retrograde" orbit.
  64.      *  This implies I = -1. <br>
  65.      *  As Orekit doesn't implement the retrograde orbit, I is always set to +1.
  66.      *  But for the sake of consistency with the theory, the retrograde factor
  67.      *  has been kept in the formulas.
  68.      *  </p>
  69.      */
  70.     private final int    I;

  71.     /** Orbit. */
  72.     private Orbit orbit;

  73.     /** B = sqrt(1 - h² - k²). */
  74.     private final double B;

  75.     /** C = 1 + p² + q². */
  76.     private final double C;

  77.     /** Equinoctial frame f vector. */
  78.     private final Vector3D f;

  79.     /** Equinoctial frame g vector. */
  80.     private final Vector3D g;

  81.     /** Equinoctial frame w vector. */
  82.     private final Vector3D w;

  83.     /** Direction cosine α. */
  84.     private final double alpha;

  85.     /** Direction cosine β. */
  86.     private final double beta;

  87.     /** Direction cosine γ. */
  88.     private final double gamma;

  89.     /** Simple constructor.
  90.      * @param orbit related mean orbit for auxiliary elements
  91.      * @param retrogradeFactor retrograde factor I [Eq. 2.1.2-(2)]
  92.      */
  93.     public AuxiliaryElements(final Orbit orbit, final int retrogradeFactor) {

  94.         // Orbit
  95.         this.orbit = orbit;

  96.         // Date of the orbit
  97.         date = orbit.getDate();

  98.         // Orbit definition frame
  99.         frame = orbit.getFrame();

  100.         // Eccentricity
  101.         ecc = orbit.getE();

  102.         // Keplerian mean motion
  103.         n = orbit.getKeplerianMeanMotion();

  104.         // Keplerian period
  105.         period = orbit.getKeplerianPeriod();

  106.         // Equinoctial elements [Eq. 2.1.2-(1)]
  107.         sma = orbit.getA();
  108.         k   = orbit.getEquinoctialEx();
  109.         h   = orbit.getEquinoctialEy();
  110.         q   = orbit.getHx();
  111.         p   = orbit.getHy();
  112.         lm  = MathUtils.normalizeAngle(orbit.getLM(), FastMath.PI);
  113.         lv  = MathUtils.normalizeAngle(orbit.getLv(), FastMath.PI);
  114.         le  = MathUtils.normalizeAngle(orbit.getLE(), FastMath.PI);

  115.         // Retrograde factor [Eq. 2.1.2-(2)]
  116.         I = retrogradeFactor;

  117.         final double k2 = k * k;
  118.         final double h2 = h * h;
  119.         final double q2 = q * q;
  120.         final double p2 = p * p;

  121.         // A, B, C parameters [Eq. 2.1.6-(1)]
  122.         B = FastMath.sqrt(1 - k2 - h2);
  123.         C = 1 + q2 + p2;

  124.         // Equinoctial reference frame [Eq. 2.1.4-(1)]
  125.         final double ooC = 1. / C;
  126.         final double px2 = 2. * p;
  127.         final double qx2 = 2. * q;
  128.         final double pq2 = px2 * q;
  129.         f = new Vector3D(ooC, new Vector3D(1. - p2 + q2, pq2, -px2 * I));
  130.         g = new Vector3D(ooC, new Vector3D(pq2 * I, (1. + p2 - q2) * I, qx2));
  131.         w = new Vector3D(ooC, new Vector3D(px2, -qx2, (1. - p2 - q2) * I));

  132.         // Direction cosines for central body [Eq. 2.1.9-(1)]
  133.         alpha = f.getZ();
  134.         beta  = g.getZ();
  135.         gamma = w.getZ();
  136.     }

  137.     /** Get the orbit.
  138.      * @return the orbit
  139.      */
  140.     public Orbit getOrbit() {
  141.         return orbit;
  142.     }

  143.     /** Get the date of the orbit.
  144.      * @return the date
  145.      */
  146.     public AbsoluteDate getDate() {
  147.         return date;
  148.     }

  149.     /** Get the definition frame of the orbit.
  150.      * @return the definition frame
  151.      */
  152.     public Frame getFrame() {
  153.         return frame;
  154.     }

  155.     /** Get the eccentricity.
  156.      * @return ecc
  157.      */
  158.     public double getEcc() {
  159.         return ecc;
  160.     }

  161.     /** Get the Keplerian mean motion.
  162.      * @return n
  163.      */
  164.     public double getMeanMotion() {
  165.         return n;
  166.     }

  167.     /** Get the Keplerian period.
  168.      * @return period
  169.      */
  170.     public double getKeplerianPeriod() {
  171.         return period;
  172.     }

  173.     /** Get the semi-major axis.
  174.      * @return the semi-major axis a
  175.      */
  176.     public double getSma() {
  177.         return sma;
  178.     }

  179.     /** Get the x component of eccentricity vector.
  180.      * <p>
  181.      * This element called k in DSST corresponds to ex
  182.      * for the {@link org.orekit.orbits.EquinoctialOrbit}
  183.      * </p>
  184.      * @return k
  185.      */
  186.     public double getK() {
  187.         return k;
  188.     }

  189.     /** Get the y component of eccentricity vector.
  190.      * <p>
  191.      * This element called h in DSST corresponds to ey
  192.      * for the {@link org.orekit.orbits.EquinoctialOrbit}
  193.      * </p>
  194.      * @return h
  195.      */
  196.     public double getH() {
  197.         return h;
  198.     }

  199.     /** Get the x component of inclination vector.
  200.      * <p>
  201.      * This element called q in DSST corresponds to hx
  202.      * for the {@link org.orekit.orbits.EquinoctialOrbit}
  203.      * </p>
  204.      * @return q
  205.      */
  206.     public double getQ() {
  207.         return q;
  208.     }

  209.     /** Get the y component of inclination vector.
  210.      * <p>
  211.      * This element called p in DSST corresponds to hy
  212.      * for the {@link org.orekit.orbits.EquinoctialOrbit}
  213.      * </p>
  214.      * @return p
  215.      */
  216.     public double getP() {
  217.         return p;
  218.     }

  219.     /** Get the mean longitude.
  220.      * @return lm
  221.      */
  222.     public double getLM() {
  223.         return lm;
  224.     }

  225.     /** Get the true longitude.
  226.      * @return lv
  227.      */
  228.     public double getLv() {
  229.         return lv;
  230.     }

  231.     /** Get the eccentric longitude.
  232.      * @return lf
  233.      */
  234.     public double getLf() {
  235.         return le;
  236.     }

  237.     /** Get the retrograde factor.
  238.      * @return the retrograde factor I
  239.      */
  240.     public int getRetrogradeFactor() {
  241.         return I;
  242.     }

  243.     /** Get B = sqrt(1 - e²).
  244.      * @return B
  245.      */
  246.     public double getB() {
  247.         return B;
  248.     }

  249.     /** Get C = 1 + p² + q².
  250.      * @return C
  251.      */
  252.     public double getC() {
  253.         return C;
  254.     }

  255.     /** Get equinoctial frame vector f.
  256.      * @return f vector
  257.      */
  258.     public Vector3D getVectorF() {
  259.         return f;
  260.     }

  261.     /** Get equinoctial frame vector g.
  262.      * @return g vector
  263.      */
  264.     public Vector3D getVectorG() {
  265.         return g;
  266.     }

  267.     /** Get equinoctial frame vector w.
  268.      * @return w vector
  269.      */
  270.     public Vector3D getVectorW() {
  271.         return w;
  272.     }

  273.     /** Get direction cosine α for central body.
  274.      * @return α
  275.      */
  276.     public double getAlpha() {
  277.         return alpha;
  278.     }

  279.     /** Get direction cosine β for central body.
  280.      * @return β
  281.      */
  282.     public double getBeta() {
  283.         return beta;
  284.     }

  285.     /** Get direction cosine γ for central body.
  286.      * @return γ
  287.      */
  288.     public double getGamma() {
  289.         return gamma;
  290.     }

  291. }