CartesianOrbit.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.orbits;

  18. import java.io.Serializable;

  19. import org.hipparchus.analysis.differentiation.UnivariateDerivative2;
  20. import org.hipparchus.geometry.euclidean.threed.FieldVector3D;
  21. import org.hipparchus.geometry.euclidean.threed.Vector3D;
  22. import org.hipparchus.linear.MatrixUtils;
  23. import org.hipparchus.util.FastMath;
  24. import org.orekit.annotation.DefaultDataContext;
  25. import org.orekit.data.DataContext;
  26. import org.orekit.frames.Frame;
  27. import org.orekit.time.AbsoluteDate;
  28. import org.orekit.utils.FieldPVCoordinates;
  29. import org.orekit.utils.PVCoordinates;
  30. import org.orekit.utils.TimeStampedPVCoordinates;


  31. /** This class holds Cartesian orbital parameters.

  32.  * <p>
  33.  * The parameters used internally are the Cartesian coordinates:
  34.  *   <ul>
  35.  *     <li>x</li>
  36.  *     <li>y</li>
  37.  *     <li>z</li>
  38.  *     <li>xDot</li>
  39.  *     <li>yDot</li>
  40.  *     <li>zDot</li>
  41.  *   </ul>
  42.  * contained in {@link PVCoordinates}.
  43.  *

  44.  * <p>
  45.  * Note that the implementation of this class delegates all non-Cartesian related
  46.  * computations ({@link #getA()}, {@link #getEquinoctialEx()}, ...) to an underlying
  47.  * instance of the {@link EquinoctialOrbit} class. This implies that using this class
  48.  * only for analytical computations which are always based on non-Cartesian
  49.  * parameters is perfectly possible but somewhat sub-optimal.
  50.  * </p>
  51.  * <p>
  52.  * The instance <code>CartesianOrbit</code> is guaranteed to be immutable.
  53.  * </p>
  54.  * @see    Orbit
  55.  * @see    KeplerianOrbit
  56.  * @see    CircularOrbit
  57.  * @see    EquinoctialOrbit
  58.  * @author Luc Maisonobe
  59.  * @author Guylaine Prat
  60.  * @author Fabien Maussion
  61.  * @author V&eacute;ronique Pommier-Maurussane
  62.  * @author Andrew Goetz
  63.  */
  64. public class CartesianOrbit extends Orbit {

  65.     /** Serializable UID. */
  66.     private static final long serialVersionUID = 20170414L;

  67.     /** 6x6 identity matrix. */
  68.     private static final double[][] SIX_BY_SIX_IDENTITY = MatrixUtils.createRealIdentityMatrix(6).getData();

  69.     /** Indicator for non-Keplerian derivatives. */
  70.     private final transient boolean hasNonKeplerianAcceleration;

  71.     /** Underlying equinoctial orbit to which high-level methods are delegated. */
  72.     private transient EquinoctialOrbit equinoctial;

  73.     /** Constructor from Cartesian parameters.
  74.      *
  75.      * <p> The acceleration provided in {@code pvCoordinates} is accessible using
  76.      * {@link #getPVCoordinates()} and {@link #getPVCoordinates(Frame)}. All other methods
  77.      * use {@code mu} and the position to compute the acceleration, including
  78.      * {@link #shiftedBy(double)} and {@link #getPVCoordinates(AbsoluteDate, Frame)}.
  79.      *
  80.      * @param pvaCoordinates the position, velocity and acceleration of the satellite.
  81.      * @param frame the frame in which the {@link PVCoordinates} are defined
  82.      * (<em>must</em> be a {@link Frame#isPseudoInertial pseudo-inertial frame})
  83.      * @param mu central attraction coefficient (m³/s²)
  84.      * @exception IllegalArgumentException if frame is not a {@link
  85.      * Frame#isPseudoInertial pseudo-inertial frame}
  86.      */
  87.     public CartesianOrbit(final TimeStampedPVCoordinates pvaCoordinates,
  88.                           final Frame frame, final double mu)
  89.         throws IllegalArgumentException {
  90.         super(pvaCoordinates, frame, mu);
  91.         hasNonKeplerianAcceleration = hasNonKeplerianAcceleration(pvaCoordinates, mu);
  92.         equinoctial = null;
  93.     }

  94.     /** Constructor from Cartesian parameters.
  95.      *
  96.      * <p> The acceleration provided in {@code pvCoordinates} is accessible using
  97.      * {@link #getPVCoordinates()} and {@link #getPVCoordinates(Frame)}. All other methods
  98.      * use {@code mu} and the position to compute the acceleration, including
  99.      * {@link #shiftedBy(double)} and {@link #getPVCoordinates(AbsoluteDate, Frame)}.
  100.      *
  101.      * @param pvaCoordinates the position and velocity of the satellite.
  102.      * @param frame the frame in which the {@link PVCoordinates} are defined
  103.      * (<em>must</em> be a {@link Frame#isPseudoInertial pseudo-inertial frame})
  104.      * @param date date of the orbital parameters
  105.      * @param mu central attraction coefficient (m³/s²)
  106.      * @exception IllegalArgumentException if frame is not a {@link
  107.      * Frame#isPseudoInertial pseudo-inertial frame}
  108.      */
  109.     public CartesianOrbit(final PVCoordinates pvaCoordinates, final Frame frame,
  110.                           final AbsoluteDate date, final double mu)
  111.         throws IllegalArgumentException {
  112.         this(new TimeStampedPVCoordinates(date, pvaCoordinates), frame, mu);
  113.     }

  114.     /** Constructor from any kind of orbital parameters.
  115.      * @param op orbital parameters to copy
  116.      */
  117.     public CartesianOrbit(final Orbit op) {
  118.         super(op.getPVCoordinates(), op.getFrame(), op.getMu());
  119.         hasNonKeplerianAcceleration = op.hasDerivatives();
  120.         if (op instanceof EquinoctialOrbit) {
  121.             equinoctial = (EquinoctialOrbit) op;
  122.         } else if (op instanceof CartesianOrbit) {
  123.             equinoctial = ((CartesianOrbit) op).equinoctial;
  124.         } else {
  125.             equinoctial = null;
  126.         }
  127.     }

  128.     /** {@inheritDoc} */
  129.     public OrbitType getType() {
  130.         return OrbitType.CARTESIAN;
  131.     }

  132.     /** Lazy evaluation of equinoctial parameters. */
  133.     private void initEquinoctial() {
  134.         if (equinoctial == null) {
  135.             if (hasDerivatives()) {
  136.                 // getPVCoordinates includes accelerations that will be interpreted as derivatives
  137.                 equinoctial = new EquinoctialOrbit(getPVCoordinates(), getFrame(), getDate(), getMu());
  138.             } else {
  139.                 // get rid of Keplerian acceleration so we don't assume
  140.                 // we have derivatives when in fact we don't have them
  141.                 equinoctial = new EquinoctialOrbit(new PVCoordinates(getPosition(),
  142.                                                                      getPVCoordinates().getVelocity()),
  143.                                                    getFrame(), getDate(), getMu());
  144.             }
  145.         }
  146.     }

  147.     /** Get the position/velocity with derivatives.
  148.      * @return position/velocity with derivatives
  149.      * @since 10.2
  150.      */
  151.     private FieldPVCoordinates<UnivariateDerivative2> getPVDerivatives() {
  152.         // PVA coordinates
  153.         final PVCoordinates pva = getPVCoordinates();
  154.         final Vector3D      p   = pva.getPosition();
  155.         final Vector3D      v   = pva.getVelocity();
  156.         final Vector3D      a   = pva.getAcceleration();
  157.         // Field coordinates
  158.         final FieldVector3D<UnivariateDerivative2> pG = new FieldVector3D<>(new UnivariateDerivative2(p.getX(), v.getX(), a.getX()),
  159.                                                                new UnivariateDerivative2(p.getY(), v.getY(), a.getY()),
  160.                                                                new UnivariateDerivative2(p.getZ(), v.getZ(), a.getZ()));
  161.         final FieldVector3D<UnivariateDerivative2> vG = new FieldVector3D<>(new UnivariateDerivative2(v.getX(), a.getX(), 0.0),
  162.                                                                new UnivariateDerivative2(v.getY(), a.getY(), 0.0),
  163.                                                                new UnivariateDerivative2(v.getZ(), a.getZ(), 0.0));
  164.         return new FieldPVCoordinates<>(pG, vG);
  165.     }

  166.     /** {@inheritDoc} */
  167.     public double getA() {
  168.         final double r  = getPosition().getNorm();
  169.         final double V2 = getPVCoordinates().getVelocity().getNormSq();
  170.         return r / (2 - r * V2 / getMu());
  171.     }

  172.     /** {@inheritDoc} */
  173.     public double getADot() {
  174.         if (hasDerivatives()) {
  175.             final FieldPVCoordinates<UnivariateDerivative2> pv = getPVDerivatives();
  176.             final UnivariateDerivative2 r  = pv.getPosition().getNorm();
  177.             final UnivariateDerivative2 V2 = pv.getVelocity().getNormSq();
  178.             final UnivariateDerivative2 a  = r.divide(r.multiply(V2).divide(getMu()).subtract(2).negate());
  179.             return a.getDerivative(1);
  180.         } else {
  181.             return Double.NaN;
  182.         }
  183.     }

  184.     /** {@inheritDoc} */
  185.     public double getE() {
  186.         final double muA = getMu() * getA();
  187.         if (isElliptical()) {
  188.             // elliptic or circular orbit
  189.             final Vector3D pvP   = getPosition();
  190.             final Vector3D pvV   = getPVCoordinates().getVelocity();
  191.             final double rV2OnMu = pvP.getNorm() * pvV.getNormSq() / getMu();
  192.             final double eSE     = Vector3D.dotProduct(pvP, pvV) / FastMath.sqrt(muA);
  193.             final double eCE     = rV2OnMu - 1;
  194.             return FastMath.sqrt(eCE * eCE + eSE * eSE);
  195.         } else {
  196.             // hyperbolic orbit
  197.             final Vector3D pvM = getPVCoordinates().getMomentum();
  198.             return FastMath.sqrt(1 - pvM.getNormSq() / muA);
  199.         }
  200.     }

  201.     /** {@inheritDoc} */
  202.     public double getEDot() {
  203.         if (hasDerivatives()) {
  204.             final FieldPVCoordinates<UnivariateDerivative2> pv = getPVDerivatives();
  205.             final FieldVector3D<UnivariateDerivative2> pvP   = pv.getPosition();
  206.             final FieldVector3D<UnivariateDerivative2> pvV   = pv.getVelocity();
  207.             final UnivariateDerivative2 r       = pvP.getNorm();
  208.             final UnivariateDerivative2 V2      = pvV.getNormSq();
  209.             final UnivariateDerivative2 rV2OnMu = r.multiply(V2).divide(getMu());
  210.             final UnivariateDerivative2 a       = r.divide(rV2OnMu.negate().add(2));
  211.             final UnivariateDerivative2 eSE     = FieldVector3D.dotProduct(pvP, pvV).divide(a.multiply(getMu()).sqrt());
  212.             final UnivariateDerivative2 eCE     = rV2OnMu.subtract(1);
  213.             final UnivariateDerivative2 e       = eCE.multiply(eCE).add(eSE.multiply(eSE)).sqrt();
  214.             return e.getDerivative(1);
  215.         } else {
  216.             return Double.NaN;
  217.         }
  218.     }

  219.     /** {@inheritDoc} */
  220.     public double getI() {
  221.         return Vector3D.angle(Vector3D.PLUS_K, getPVCoordinates().getMomentum());
  222.     }

  223.     /** {@inheritDoc} */
  224.     public double getIDot() {
  225.         if (hasDerivatives()) {
  226.             final FieldPVCoordinates<UnivariateDerivative2> pv = getPVDerivatives();
  227.             final FieldVector3D<UnivariateDerivative2> momentum =
  228.                             FieldVector3D.crossProduct(pv.getPosition(), pv.getVelocity());
  229.             final UnivariateDerivative2 i = FieldVector3D.angle(Vector3D.PLUS_K, momentum);
  230.             return i.getDerivative(1);
  231.         } else {
  232.             return Double.NaN;
  233.         }
  234.     }

  235.     /** {@inheritDoc} */
  236.     public double getEquinoctialEx() {
  237.         initEquinoctial();
  238.         return equinoctial.getEquinoctialEx();
  239.     }

  240.     /** {@inheritDoc} */
  241.     public double getEquinoctialExDot() {
  242.         initEquinoctial();
  243.         return equinoctial.getEquinoctialExDot();
  244.     }

  245.     /** {@inheritDoc} */
  246.     public double getEquinoctialEy() {
  247.         initEquinoctial();
  248.         return equinoctial.getEquinoctialEy();
  249.     }

  250.     /** {@inheritDoc} */
  251.     public double getEquinoctialEyDot() {
  252.         initEquinoctial();
  253.         return equinoctial.getEquinoctialEyDot();
  254.     }

  255.     /** {@inheritDoc} */
  256.     public double getHx() {
  257.         final Vector3D w = getPVCoordinates().getMomentum().normalize();
  258.         // Check for equatorial retrograde orbit
  259.         if ((w.getX() * w.getX() + w.getY() * w.getY()) == 0 && w.getZ() < 0) {
  260.             return Double.NaN;
  261.         }
  262.         return -w.getY() / (1 + w.getZ());
  263.     }

  264.     /** {@inheritDoc} */
  265.     public double getHxDot() {
  266.         if (hasDerivatives()) {
  267.             final FieldPVCoordinates<UnivariateDerivative2> pv = getPVDerivatives();
  268.             final FieldVector3D<UnivariateDerivative2> w =
  269.                             FieldVector3D.crossProduct(pv.getPosition(), pv.getVelocity()).normalize();
  270.             // Check for equatorial retrograde orbit
  271.             final double x = w.getX().getValue();
  272.             final double y = w.getY().getValue();
  273.             final double z = w.getZ().getValue();
  274.             if ((x * x + y * y) == 0 && z < 0) {
  275.                 return Double.NaN;
  276.             }
  277.             final UnivariateDerivative2 hx = w.getY().negate().divide(w.getZ().add(1));
  278.             return hx.getDerivative(1);
  279.         } else {
  280.             return Double.NaN;
  281.         }
  282.     }

  283.     /** {@inheritDoc} */
  284.     public double getHy() {
  285.         final Vector3D w = getPVCoordinates().getMomentum().normalize();
  286.         // Check for equatorial retrograde orbit
  287.         if ((w.getX() * w.getX() + w.getY() * w.getY()) == 0 && w.getZ() < 0) {
  288.             return Double.NaN;
  289.         }
  290.         return  w.getX() / (1 + w.getZ());
  291.     }

  292.     /** {@inheritDoc} */
  293.     public double getHyDot() {
  294.         if (hasDerivatives()) {
  295.             final FieldPVCoordinates<UnivariateDerivative2> pv = getPVDerivatives();
  296.             final FieldVector3D<UnivariateDerivative2> w =
  297.                             FieldVector3D.crossProduct(pv.getPosition(), pv.getVelocity()).normalize();
  298.             // Check for equatorial retrograde orbit
  299.             final double x = w.getX().getValue();
  300.             final double y = w.getY().getValue();
  301.             final double z = w.getZ().getValue();
  302.             if ((x * x + y * y) == 0 && z < 0) {
  303.                 return Double.NaN;
  304.             }
  305.             final UnivariateDerivative2 hy = w.getX().divide(w.getZ().add(1));
  306.             return hy.getDerivative(1);
  307.         } else {
  308.             return Double.NaN;
  309.         }
  310.     }

  311.     /** {@inheritDoc} */
  312.     public double getLv() {
  313.         initEquinoctial();
  314.         return equinoctial.getLv();
  315.     }

  316.     /** {@inheritDoc} */
  317.     public double getLvDot() {
  318.         initEquinoctial();
  319.         return equinoctial.getLvDot();
  320.     }

  321.     /** {@inheritDoc} */
  322.     public double getLE() {
  323.         initEquinoctial();
  324.         return equinoctial.getLE();
  325.     }

  326.     /** {@inheritDoc} */
  327.     public double getLEDot() {
  328.         initEquinoctial();
  329.         return equinoctial.getLEDot();
  330.     }

  331.     /** {@inheritDoc} */
  332.     public double getLM() {
  333.         initEquinoctial();
  334.         return equinoctial.getLM();
  335.     }

  336.     /** {@inheritDoc} */
  337.     public double getLMDot() {
  338.         initEquinoctial();
  339.         return equinoctial.getLMDot();
  340.     }

  341.     /** {@inheritDoc} */
  342.     @Override
  343.     public boolean hasDerivatives() {
  344.         return hasNonKeplerianAcceleration;
  345.     }

  346.     /** {@inheritDoc} */
  347.     protected Vector3D initPosition() {
  348.         // nothing to do here, as the canonical elements are already the Cartesian ones
  349.         return getPVCoordinates().getPosition();
  350.     }

  351.     /** {@inheritDoc} */
  352.     protected TimeStampedPVCoordinates initPVCoordinates() {
  353.         // nothing to do here, as the canonical elements are already the Cartesian ones
  354.         return getPVCoordinates();
  355.     }

  356.     /** {@inheritDoc} */
  357.     public CartesianOrbit shiftedBy(final double dt) {
  358.         final PVCoordinates shiftedPV = shiftPV(dt);
  359.         return new CartesianOrbit(shiftedPV, getFrame(), getDate().shiftedBy(dt), getMu());
  360.     }

  361.     /** Compute shifted position and velocity.
  362.      * @param dt time shift
  363.      * @return shifted position and velocity
  364.      */
  365.     private PVCoordinates shiftPV(final double dt) {

  366.         final Vector3D pvP = getPosition();
  367.         final PVCoordinates shiftedPV = KeplerianMotionCartesianUtility.predictPositionVelocity(dt, pvP,
  368.             getPVCoordinates().getVelocity(), getMu());

  369.         if (hasNonKeplerianAcceleration) {

  370.             // extract non-Keplerian part of the initial acceleration
  371.             final double r2 = pvP.getNormSq();
  372.             final double r = FastMath.sqrt(r2);
  373.             final Vector3D nonKeplerianAcceleration = new Vector3D(1, getPVCoordinates().getAcceleration(),
  374.                                                                    getMu() / (r2 * r), pvP);

  375.             // add the quadratic motion due to the non-Keplerian acceleration to the Keplerian motion
  376.             final Vector3D shiftedP = shiftedPV.getPosition();
  377.             final Vector3D shiftedV = shiftedPV.getVelocity();
  378.             final Vector3D fixedP   = new Vector3D(1, shiftedP,
  379.                                                    0.5 * dt * dt, nonKeplerianAcceleration);
  380.             final double   fixedR2 = fixedP.getNormSq();
  381.             final double   fixedR  = FastMath.sqrt(fixedR2);
  382.             final Vector3D fixedV  = new Vector3D(1, shiftedV,
  383.                                                   dt, nonKeplerianAcceleration);
  384.             final Vector3D fixedA  = new Vector3D(-getMu() / (fixedR2 * fixedR), shiftedP,
  385.                                                   1, nonKeplerianAcceleration);

  386.             return new PVCoordinates(fixedP, fixedV, fixedA);

  387.         } else {
  388.             // don't include acceleration,
  389.             // so the shifted orbit is not considered to have derivatives
  390.             return shiftedPV;
  391.         }

  392.     }

  393.     @Override
  394.     protected double[][] computeJacobianMeanWrtCartesian() {
  395.         return SIX_BY_SIX_IDENTITY;
  396.     }

  397.     @Override
  398.     protected double[][] computeJacobianEccentricWrtCartesian() {
  399.         return SIX_BY_SIX_IDENTITY;
  400.     }

  401.     @Override
  402.     protected double[][] computeJacobianTrueWrtCartesian() {
  403.         return SIX_BY_SIX_IDENTITY;
  404.     }

  405.     /** {@inheritDoc} */
  406.     public void addKeplerContribution(final PositionAngleType type, final double gm,
  407.                                       final double[] pDot) {

  408.         final PVCoordinates pv = getPVCoordinates();

  409.         // position derivative is velocity
  410.         final Vector3D velocity = pv.getVelocity();
  411.         pDot[0] += velocity.getX();
  412.         pDot[1] += velocity.getY();
  413.         pDot[2] += velocity.getZ();

  414.         // velocity derivative is Newtonian acceleration
  415.         final Vector3D position = pv.getPosition();
  416.         final double r2         = position.getNormSq();
  417.         final double coeff      = -gm / (r2 * FastMath.sqrt(r2));
  418.         pDot[3] += coeff * position.getX();
  419.         pDot[4] += coeff * position.getY();
  420.         pDot[5] += coeff * position.getZ();

  421.     }

  422.     /**  Returns a string representation of this Orbit object.
  423.      * @return a string representation of this object
  424.      */
  425.     public String toString() {
  426.         // use only the six defining elements, like the other Orbit.toString() methods
  427.         final String comma = ", ";
  428.         final PVCoordinates pv = getPVCoordinates();
  429.         final Vector3D p = pv.getPosition();
  430.         final Vector3D v = pv.getVelocity();
  431.         return "Cartesian parameters: {P(" +
  432.                 p.getX() + comma +
  433.                 p.getY() + comma +
  434.                 p.getZ() + "), V(" +
  435.                 v.getX() + comma +
  436.                 v.getY() + comma +
  437.                 v.getZ() + ")}";
  438.     }

  439.     /** Replace the instance with a data transfer object for serialization.
  440.      * <p>
  441.      * This intermediate class serializes all needed parameters,
  442.      * including position-velocity which are <em>not</em> serialized by parent
  443.      * {@link Orbit} class.
  444.      * </p>
  445.      * @return data transfer object that will be serialized
  446.      */
  447.     @DefaultDataContext
  448.     private Object writeReplace() {
  449.         return new DTO(this);
  450.     }

  451.     /** Internal class used only for serialization. */
  452.     @DefaultDataContext
  453.     private static class DTO implements Serializable {

  454.         /** Serializable UID. */
  455.         private static final long serialVersionUID = 20170414L;

  456.         /** Double values. */
  457.         private double[] d;

  458.         /** Frame in which are defined the orbital parameters. */
  459.         private final Frame frame;

  460.         /** Simple constructor.
  461.          * @param orbit instance to serialize
  462.          */
  463.         private DTO(final CartesianOrbit orbit) {

  464.             final TimeStampedPVCoordinates pv = orbit.getPVCoordinates();

  465.             // decompose date
  466.             final AbsoluteDate j2000Epoch =
  467.                     DataContext.getDefault().getTimeScales().getJ2000Epoch();
  468.             final double epoch  = FastMath.floor(pv.getDate().durationFrom(j2000Epoch));
  469.             final double offset = pv.getDate().durationFrom(j2000Epoch.shiftedBy(epoch));

  470.             if (orbit.hasDerivatives()) {
  471.                 this.d = new double[] {
  472.                     epoch, offset, orbit.getMu(),
  473.                     pv.getPosition().getX(),     pv.getPosition().getY(),     pv.getPosition().getZ(),
  474.                     pv.getVelocity().getX(),     pv.getVelocity().getY(),     pv.getVelocity().getZ(),
  475.                     pv.getAcceleration().getX(), pv.getAcceleration().getY(), pv.getAcceleration().getZ()
  476.                 };
  477.             } else {
  478.                 this.d = new double[] {
  479.                     epoch, offset, orbit.getMu(),
  480.                     pv.getPosition().getX(),     pv.getPosition().getY(),     pv.getPosition().getZ(),
  481.                     pv.getVelocity().getX(),     pv.getVelocity().getY(),     pv.getVelocity().getZ()
  482.                 };
  483.             }

  484.             this.frame = orbit.getFrame();

  485.         }

  486.         /** Replace the deserialized data transfer object with a {@link CartesianOrbit}.
  487.          * @return replacement {@link CartesianOrbit}
  488.          */
  489.         private Object readResolve() {
  490.             final AbsoluteDate j2000Epoch =
  491.                     DataContext.getDefault().getTimeScales().getJ2000Epoch();
  492.             if (d.length >= 12) {
  493.                 // we have derivatives
  494.                 return new CartesianOrbit(new TimeStampedPVCoordinates(j2000Epoch.shiftedBy(d[0]).shiftedBy(d[1]),
  495.                                                                        new Vector3D(d[3], d[ 4], d[ 5]),
  496.                                                                        new Vector3D(d[6], d[ 7], d[ 8]),
  497.                                                                        new Vector3D(d[9], d[10], d[11])),
  498.                                           frame, d[2]);
  499.             } else {
  500.                 // we don't have derivatives
  501.                 return new CartesianOrbit(new TimeStampedPVCoordinates(j2000Epoch.shiftedBy(d[0]).shiftedBy(d[1]),
  502.                                                                        new Vector3D(d[3], d[ 4], d[ 5]),
  503.                                                                        new Vector3D(d[6], d[ 7], d[ 8])),
  504.                                           frame, d[2]);
  505.             }
  506.         }

  507.     }

  508. }