J2DifferentialEffect.java

  1. /* Copyright 2002-2013 CS Systèmes d'Information
  2.  * Licensed to CS Systèmes d'Information (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;

  19. import org.apache.commons.math3.util.FastMath;
  20. import org.orekit.errors.OrekitException;
  21. import org.orekit.forces.gravity.potential.UnnormalizedSphericalHarmonicsProvider;
  22. import org.orekit.orbits.EquinoctialOrbit;
  23. import org.orekit.orbits.Orbit;
  24. import org.orekit.orbits.OrbitType;
  25. import org.orekit.orbits.PositionAngle;
  26. import org.orekit.propagation.SpacecraftState;
  27. import org.orekit.time.AbsoluteDate;

  29. /** Analytical model for J2 effect.
  30.  * <p>
  31.  * This class computes the differential effect of J2 due to an initial orbit
  32.  * offset. A typical case is when an inclination maneuver changes an orbit
  33.  * inclination at time t<sub>0</sub>. As ascending node drift rate depends on
  34.  * inclination, the change induces a time-dependent change in ascending node
  35.  * for later dates.
  36.  * </p>
  37.  * @see org.orekit.forces.maneuvers.SmallManeuverAnalyticalModel
  38.  * @author Luc Maisonobe
  39.  */
  40. public class J2DifferentialEffect
  41.     implements AdapterPropagator.DifferentialEffect {

  43.     /** Reference date. */
  44.     private final AbsoluteDate referenceDate;

  46.     /** Differential drift on perigee argument. */
  47.     private final double dPaDot;

  49.     /** Differential drift on ascending node. */
  50.     private final double dRaanDot;

  52.     /** Indicator for applying effect before reference date. */
  53.     private final boolean applyBefore;

  55.     /** Simple constructor.
  56.      * <p>
  57.      * The {@code applyBefore} parameter is mainly used when the differential
  58.      * effect is associated with a maneuver. In this case, the parameter must be
  59.      * set to {@code false}.
  60.      * </p>
  61.      * @param original original state at reference date
  62.      * @param directEffect direct effect changing the orbit
  63.      * @param applyBefore if true, effect is applied both before and after
  64.      * reference date, if false it is only applied after reference date
  65.      * @param gravityField gravity field to use
  66.      * @exception OrekitException if gravity field does not contain J2 coefficient
  67.      */
  68.     public J2DifferentialEffect(final SpacecraftState original,
  69.                                 final AdapterPropagator.DifferentialEffect directEffect,
  70.                                 final boolean applyBefore,
  71.                                 final UnnormalizedSphericalHarmonicsProvider gravityField)
  72.         throws OrekitException {
  73.         this(original, directEffect, applyBefore,
  74.              gravityField.getAe(), gravityField.getMu(),
  75.              -gravityField.onDate(original.getDate()).getUnnormalizedCnm(2, 0));
  76.     }

  78.     /** Simple constructor.
  79.          * <p>
  80.          * The {@code applyBefore} parameter is mainly used when the differential
  81.          * effect is associated with a maneuver. In this case, the parameter must be
  82.          * set to {@code false}.
  83.          * </p>
  84.          * @param orbit0 original orbit at reference date
  85.          * @param orbit1 shifted orbit at reference date
  86.          * @param applyBefore if true, effect is applied both before and after
  87.          * reference date, if false it is only applied after reference date
  88.          * @param gravityField gravity field to use
  89.          * @exception OrekitException if gravity field does not contain J2 coefficient
  90.          */
  91.     public J2DifferentialEffect(final Orbit orbit0, final Orbit orbit1, final boolean applyBefore,
  92.                                 final UnnormalizedSphericalHarmonicsProvider gravityField)
  93.         throws OrekitException {
  94.         this(orbit0, orbit1, applyBefore,
  95.              gravityField.getAe(), gravityField.getMu(),
  96.              -gravityField.onDate(orbit0.getDate()).getUnnormalizedCnm(2, 0));
  97.     }

  99.     /** Simple constructor.
  100.      * <p>
  101.      * The {@code applyBefore} parameter is mainly used when the differential
  102.      * effect is associated with a maneuver. In this case, the parameter must be
  103.      * set to {@code false}.
  104.      * </p>
  105.      * @param original original state at reference date
  106.      * @param directEffect direct effect changing the orbit
  107.      * @param applyBefore if true, effect is applied both before and after
  108.      * reference date, if false it is only applied after reference date
  109.      * @param referenceRadius reference radius of the Earth for the potential model (m)
  110.      * @param mu central attraction coefficient (m<sup>3</sup>/s<sup>2</sup>)
  111.      * @param j2 un-normalized zonal coefficient (about +1.08e-3 for Earth)
  112.      * @exception OrekitException if direct effect cannot be applied
  113.      */
  114.     public J2DifferentialEffect(final SpacecraftState original,
  115.                                 final AdapterPropagator.DifferentialEffect directEffect,
  116.                                 final boolean applyBefore,
  117.                                 final double referenceRadius, final double mu, final double j2)
  118.         throws OrekitException {
  119.         this(original.getOrbit(),
  120.              directEffect.apply(original.shiftedBy(0.001)).getOrbit().shiftedBy(-0.001),
  121.              applyBefore, referenceRadius, mu, j2);
  122.     }

  124.     /** Simple constructor.
  125.      * <p>
  126.      * The {@code applyBefore} parameter is mainly used when the differential
  127.      * effect is associated with a maneuver. In this case, the parameter must be
  128.      * set to {@code false}.
  129.      * </p>
  130.      * @param orbit0 original orbit at reference date
  131.      * @param orbit1 shifted orbit at reference date
  132.      * @param applyBefore if true, effect is applied both before and after
  133.      * reference date, if false it is only applied after reference date
  134.      * @param referenceRadius reference radius of the Earth for the potential model (m)
  135.      * @param mu central attraction coefficient (m<sup>3</sup>/s<sup>2</sup>)
  136.      * @param j2 un-normalized zonal coefficient (about +1.08e-3 for Earth)
  137.      */
  138.     public J2DifferentialEffect(final Orbit orbit0, final Orbit orbit1, final boolean applyBefore,
  139.                                 final double referenceRadius, final double mu, final double j2) {

  141.         this.referenceDate = orbit0.getDate();
  142.         this.applyBefore   = applyBefore;

  144.         // extract useful parameters
  145.         final double a0 = orbit0.getA();
  146.         final double e0 = orbit0.getE();
  147.         final double i0 = orbit0.getI();
  148.         final double a1 = orbit1.getA();
  149.         final double e1 = orbit1.getE();
  150.         final double i1 = orbit1.getI();

  152.         // compute reference drifts
  153.         final double oMe2       = 1 - e0 * e0;
  154.         final double ratio      = referenceRadius / (a0 * oMe2);
  155.         final double cosI       = FastMath.cos(i0);
  156.         final double sinI       = FastMath.sin(i0);
  157.         final double n          = FastMath.sqrt(mu / a0) / a0;
  158.         final double c          = ratio * ratio * n * j2;
  159.         final double refPaDot   =  0.75 * c * (4 - 5 * sinI * sinI);
  160.         final double refRaanDot = -1.5  * c * cosI;

  162.         // differential model on perigee argument drift
  163.         final double dPaDotDa = -3.5 * refPaDot / a0;
  164.         final double dPaDotDe = 4 * refPaDot * e0 / oMe2;
  165.         final double dPaDotDi = -7.5 * c * sinI * cosI;
  166.         dPaDot = dPaDotDa * (a1 - a0) + dPaDotDe * (e1 - e0) + dPaDotDi * (i1 - i0);

  168.         // differential model on ascending node drift
  169.         final double dRaanDotDa = -3.5 * refRaanDot / a0;
  170.         final double dRaanDotDe = 4 * refRaanDot * e0 / oMe2;
  171.         final double dRaanDotDi = -refRaanDot * FastMath.tan(i0);
  172.         dRaanDot = dRaanDotDa * (a1 - a0) + dRaanDotDe * (e1 - e0) + dRaanDotDi * (i1 - i0);

  174.     }

  176.     /** Compute the effect of the maneuver on an orbit.
  177.      * @param orbit1 original orbit at t<sub>1</sub>, without maneuver
  178.      * @return orbit at t<sub>1</sub>, taking the maneuver
  179.      * into account if t<sub>1</sub> &gt; t<sub>0</sub>
  180.      * @see #apply(SpacecraftState)
  181.      */
  182.     public Orbit apply(final Orbit orbit1) {

  184.         if (orbit1.getDate().compareTo(referenceDate) <= 0 && !applyBefore) {
  185.             // the orbit change has not occurred yet, don't change anything
  186.             return orbit1;
  187.         }

  189.         return updateOrbit(orbit1);

  191.     }

  193.     /** {@inheritDoc} */
  194.     public SpacecraftState apply(final SpacecraftState state1) {

  196.         if (state1.getDate().compareTo(referenceDate) <= 0 && !applyBefore) {
  197.             // the orbit change has not occurred yet, don't change anything
  198.             return state1;
  199.         }

  201.         return new SpacecraftState(updateOrbit(state1.getOrbit()),
  202.                                    state1.getAttitude(), state1.getMass());

  204.     }

  206.     /** Compute the differential effect of J2 on an orbit.
  207.      * @param orbit1 original orbit at t<sub>1</sub>, without differential J2
  208.      * @return orbit at t<sub>1</sub>, always taking the effect into account
  209.      */
  210.     private Orbit updateOrbit(final Orbit orbit1) {

  212.         // convert current orbital state to equinoctial elements
  213.         final EquinoctialOrbit original =
  214.                 (EquinoctialOrbit) OrbitType.EQUINOCTIAL.convertType(orbit1);

  216.         // compute differential effect
  217.         final AbsoluteDate date = original.getDate();
  218.         final double dt         = date.durationFrom(referenceDate);
  219.         final double dPaRaan    = (dPaDot + dRaanDot) * dt;
  220.         final double cPaRaan    = FastMath.cos(dPaRaan);
  221.         final double sPaRaan    = FastMath.sin(dPaRaan);
  222.         final double dRaan      = dRaanDot * dt;
  223.         final double cRaan      = FastMath.cos(dRaan);
  224.         final double sRaan      = FastMath.sin(dRaan);

  226.         final double ex         = original.getEquinoctialEx() * cPaRaan -
  227.                                   original.getEquinoctialEy() * sPaRaan;
  228.         final double ey         = original.getEquinoctialEx() * sPaRaan +
  229.                                   original.getEquinoctialEy() * cPaRaan;
  230.         final double hx         = original.getHx() * cRaan - original.getHy() * sRaan;
  231.         final double hy         = original.getHx() * sRaan + original.getHy() * cRaan;
  232.         final double lambda     = original.getLv() + dPaRaan;

  234.         // build updated orbit
  235.         final EquinoctialOrbit updated =
  236.                 new EquinoctialOrbit(original.getA(), ex, ey, hx, hy, lambda, PositionAngle.TRUE,
  237.                                      original.getFrame(), date, original.getMu());

  239.         // convert to required type
  240.         return orbit1.getType().convertType(updated);

  242.     }

  244. }
Created with JaCoCo 0.6.3.201306030806