FieldAuxiliaryElements.java

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

  19. import org.hipparchus.CalculusFieldElement;
  20. import org.hipparchus.geometry.euclidean.threed.FieldVector3D;
  21. import org.hipparchus.util.FastMath;
  22. import org.hipparchus.util.MathUtils;
  23. import org.orekit.frames.Frame;
  24. import org.orekit.orbits.FieldOrbit;
  25. import org.orekit.time.FieldAbsoluteDate;

  27. /** Container class for common parameters used by all DSST forces.
  28.  *  <p>
  29.  *  Most of them are defined in Danielson paper at § 2.1.
  30.  *  </p>
  31.  */
  32. public class FieldAuxiliaryElements<T extends CalculusFieldElement<T>> {

  34.     /** Orbit date. */
  35.     private final FieldAbsoluteDate<T> date;

  37.     /** Orbit frame. */
  38.     private final Frame frame;

  40.     /** Eccentricity. */
  41.     private final T ecc;

  43.     /** Keplerian mean motion. */
  44.     private final T n;

  46.     /** Keplerian period. */
  47.     private final T period;

  49.     /** Semi-major axis. */
  50.     private final T sma;

  52.     /** x component of eccentricity vector. */
  53.     private final T k;

  55.     /** y component of eccentricity vector. */
  56.     private final T h;

  58.     /** x component of inclination vector. */
  59.     private final T q;

  61.     /** y component of inclination vector. */
  62.     private final T p;

  64.     /** Mean longitude. */
  65.     private final T lm;

  67.     /** True longitude. */
  68.     private final T lv;

  70.     /** Eccentric longitude. */
  71.     private final T le;

  73.     /** Retrograde factor I.
  74.      *  <p>
  75.      *  DSST model needs equinoctial orbit as internal representation.
  76.      *  Classical equinoctial elements have discontinuities when inclination
  77.      *  is close to zero. In this representation, I = +1. <br>
  78.      *  To avoid this discontinuity, another representation exists and equinoctial
  79.      *  elements can be expressed in a different way, called "retrograde" orbit.
  80.      *  This implies I = -1. <br>
  81.      *  As Orekit doesn't implement the retrograde orbit, I is always set to +1.
  82.      *  But for the sake of consistency with the theory, the retrograde factor
  83.      *  has been kept in the formulas.
  84.      *  </p>
  85.      */
  86.     private final int    I;

  88.     /** B = sqrt(1 - h² - k²). */
  89.     private final T B;

  91.     /** C = 1 + p² + q². */
  92.     private final T C;

  94.     /** Equinoctial frame f vector. */
  95.     private final FieldVector3D<T> f;

  97.     /** Equinoctial frame g vector. */
  98.     private final FieldVector3D<T> g;

  100.     /** Equinoctial frame w vector. */
  101.     private final FieldVector3D<T> w;

  103.     /** Direction cosine α. */
  104.     private final T alpha;

  106.     /** Direction cosine β. */
  107.     private final T beta;

  109.     /** Direction cosine γ. */
  110.     private final T gamma;

  112.     /** Simple constructor.
  113.      * @param orbit related mean orbit for auxiliary elements
  114.      * @param retrogradeFactor retrograde factor I [Eq. 2.1.2-(2)]
  115.      */
  116.     public FieldAuxiliaryElements(final FieldOrbit<T> orbit, final int retrogradeFactor) {

  118.         final T pi = orbit.getDate().getField().getZero().getPi();

  120.         // Date of the orbit
  121.         date = orbit.getDate();

  123.         // Orbit definition frame
  124.         frame = orbit.getFrame();

  126.         // Eccentricity
  127.         ecc = orbit.getE();

  129.         // Keplerian mean motion
  130.         n = orbit.getKeplerianMeanMotion();

  132.         // Keplerian period
  133.         period = orbit.getKeplerianPeriod();

  135.         // Equinoctial elements [Eq. 2.1.2-(1)]
  136.         sma = orbit.getA();
  137.         k   = orbit.getEquinoctialEx();
  138.         h   = orbit.getEquinoctialEy();
  139.         q   = orbit.getHx();
  140.         p   = orbit.getHy();
  141.         lm  = MathUtils.normalizeAngle(orbit.getLM(), pi);
  142.         lv  = MathUtils.normalizeAngle(orbit.getLv(), pi);
  143.         le  = MathUtils.normalizeAngle(orbit.getLE(), pi);

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

  148.         final T k2 = k.multiply(k);
  149.         final T h2 = h.multiply(h);
  150.         final T q2 = q.multiply(q);
  151.         final T p2 = p.multiply(p);

  153.         // A, B, C parameters [Eq. 2.1.6-(1)]
  154.         B = FastMath.sqrt(k2.add(h2).negate().add(1.));
  155.         C = q2.add(p2).add(1);

  157.         // Equinoctial reference frame [Eq. 2.1.4-(1)]
  158.         final T ooC = C.reciprocal();
  159.         final T px2 = p.multiply(2.);
  160.         final T qx2 = q.multiply(2.);
  161.         final T pq2 = px2.multiply(q);
  162.         f = new FieldVector3D<>(ooC, new FieldVector3D<>(p2.negate().add(1.).add(q2), pq2, px2.multiply(I).negate()));
  163.         g = new FieldVector3D<>(ooC, new FieldVector3D<>(pq2.multiply(I), (p2.add(1.).subtract(q2)).multiply(I), qx2));
  164.         w = new FieldVector3D<>(ooC, new FieldVector3D<>(px2, qx2.negate(), (p2.add(q2).negate().add(1.)).multiply(I)));

  166.         // Direction cosines for central body [Eq. 2.1.9-(1)]
  167.         alpha = (T) f.getZ();
  168.         beta  = (T) g.getZ();
  169.         gamma = (T) w.getZ();
  170.     }

  172.     /** Get the date of the orbit.
  173.      * @return the date
  174.      */
  175.     public FieldAbsoluteDate<T> getDate() {
  176.         return date;
  177.     }

  179.     /** Get the definition frame of the orbit.
  180.      * @return the definition frame
  181.      */
  182.     public Frame getFrame() {
  183.         return frame;
  184.     }

  186.     /** Get the eccentricity.
  187.      * @return ecc
  188.      */
  189.     public T getEcc() {
  190.         return ecc;
  191.     }

  193.     /** Get the Keplerian mean motion.
  194.      * @return n
  195.      */
  196.     public T getMeanMotion() {
  197.         return n;
  198.     }

  200.     /** Get the Keplerian period.
  201.      * @return period
  202.      */
  203.     public T getKeplerianPeriod() {
  204.         return period;
  205.     }

  207.     /** Get the semi-major axis.
  208.      * @return the semi-major axis a
  209.      */
  210.     public T getSma() {
  211.         return sma;
  212.     }

  214.     /** Get the x component of eccentricity vector.
  215.      * <p>
  216.      * This element called k in DSST corresponds to ex
  217.      * for the {@link org.orekit.orbits.EquinoctialOrbit}
  218.      * </p>
  219.      * @return k
  220.      */
  221.     public T getK() {
  222.         return k;
  223.     }

  225.     /** Get the y component of eccentricity vector.
  226.      * <p>
  227.      * This element called h in DSST corresponds to ey
  228.      * for the {@link org.orekit.orbits.EquinoctialOrbit}
  229.      * </p>
  230.      * @return h
  231.      */
  232.     public T getH() {
  233.         return h;
  234.     }

  236.     /** Get the x component of inclination vector.
  237.      * <p>
  238.      * This element called q in DSST corresponds to hx
  239.      * for the {@link org.orekit.orbits.EquinoctialOrbit}
  240.      * </p>
  241.      * @return q
  242.      */
  243.     public T getQ() {
  244.         return q;
  245.     }

  247.     /** Get the y component of inclination vector.
  248.      * <p>
  249.      * This element called p in DSST corresponds to hy
  250.      * for the {@link org.orekit.orbits.EquinoctialOrbit}
  251.      * </p>
  252.      * @return p
  253.      */
  254.     public T getP() {
  255.         return p;
  256.     }

  258.     /** Get the mean longitude.
  259.      * @return lm
  260.      */
  261.     public T getLM() {
  262.         return lm;
  263.     }

  265.     /** Get the true longitude.
  266.      * @return lv
  267.      */
  268.     public T getLv() {
  269.         return lv;
  270.     }

  272.     /** Get the eccentric longitude.
  273.      * @return le
  274.      */
  275.     public T getLe() {
  276.         return le;
  277.     }

  279.     /** Get the retrograde factor.
  280.      * @return the retrograde factor I
  281.      */
  282.     public int getRetrogradeFactor() {
  283.         return I;
  284.     }

  286.     /** Get B = sqrt(1 - e²).
  287.      * @return B
  288.      */
  289.     public T getB() {
  290.         return B;
  291.     }

  293.     /** Get C = 1 + p² + q².
  294.      * @return C
  295.      */
  296.     public T getC() {
  297.         return C;
  298.     }

  300.     /** Get equinoctial frame vector f.
  301.      * @return f vector
  302.      */
  303.     public FieldVector3D<T> getVectorF() {
  304.         return f;
  305.     }

  307.     /** Get equinoctial frame vector g.
  308.      * @return g vector
  309.      */
  310.     public FieldVector3D<T> getVectorG() {
  311.         return g;
  312.     }

  314.     /** Get equinoctial frame vector w.
  315.      * @return w vector
  316.      */
  317.     public FieldVector3D<T> getVectorW() {
  318.         return w;
  319.     }

  321.     /** Get direction cosine α for central body.
  322.      * @return α
  323.      */
  324.     public T getAlpha() {
  325.         return alpha;
  326.     }

  328.     /** Get direction cosine β for central body.
  329.      * @return β
  330.      */
  331.     public T getBeta() {
  332.         return beta;
  333.     }

  335.     /** Get direction cosine γ for central body.
  336.      * @return γ
  337.      */
  338.     public T getGamma() {
  339.         return gamma;
  340.     }

  342. }
Created with JaCoCo 0.8.8.202204050719