StateCovariance.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;

  18. import org.hipparchus.linear.Array2DRowRealMatrix;
  19. import org.hipparchus.linear.MatrixUtils;
  20. import org.hipparchus.linear.RealMatrix;
  21. import org.hipparchus.util.FastMath;
  22. import org.orekit.errors.OrekitException;
  23. import org.orekit.errors.OrekitMessages;
  24. import org.orekit.frames.Frame;
  25. import org.orekit.frames.LOFType;
  26. import org.orekit.frames.Transform;
  27. import org.orekit.orbits.Orbit;
  28. import org.orekit.orbits.OrbitType;
  29. import org.orekit.orbits.PositionAngle;
  30. import org.orekit.time.AbsoluteDate;
  31. import org.orekit.time.TimeStamped;
  32. import org.orekit.utils.CartesianDerivativesFilter;

  33. /** This class is the representation of a covariance matrix at a given date.
  34.  * <p>
  35.  * Currently, the covariance only represents the orbital elements.
  36.  * <p>
  37.  * It is possible to change the covariance frame by using the
  38.  * {@link #changeCovarianceFrame(Orbit, Frame)} or {@link #changeCovarianceFrame(Orbit, LOFType)} method.
  39.  * These methods are based on Equations (18) and (20) of <i>Covariance Transformations for Satellite
  40.  * Flight Dynamics Operations</i> by David A. SVallado.
  41.  * <p>
  42.  * Finally, covariance orbit type can be changed using the
  43.  * {@link #changeCovarianceType(Orbit, OrbitType, PositionAngle) changeCovarianceType(Orbit, OrbitType, PositionAngle)}
  44.  * method.
  45.  * </p>
  46.  * @author Bryan Cazabonne
  47.  * @author Vincent Cucchietti
  48.  * @since 11.3
  49.  */
  50. public class StateCovariance implements TimeStamped {

  51.     /** State dimension. */
  52.     private static final int STATE_DIMENSION = 6;

  53.     /** Default position angle for covariance expressed in cartesian elements. */
  54.     private static final PositionAngle DEFAULT_POSITION_ANGLE = PositionAngle.TRUE;

  55.     /** Orbital covariance [6x6]. */
  56.     private final RealMatrix orbitalCovariance;

  57.     /** Covariance frame (can be null if lofType is defined). */
  58.     private final Frame frame;

  59.     /** Covariance LOF type (can be null if frame is defined). */
  60.     private final LOFType lofType;

  61.     /** Covariance epoch. */
  62.     private final AbsoluteDate epoch;

  63.     /** Covariance orbit type. */
  64.     private final OrbitType orbitType;

  65.     /** Covariance position angle type (not used if orbitType is CARTESIAN). */
  66.     private final PositionAngle angleType;

  67.     /**
  68.      * Constructor.
  69.      * @param orbitalCovariance 6x6 orbital parameters covariance
  70.      * @param epoch epoch of the covariance
  71.      * @param lofType covariance LOF type
  72.      */
  73.     public StateCovariance(final RealMatrix orbitalCovariance, final AbsoluteDate epoch, final LOFType lofType) {
  74.         this(orbitalCovariance, epoch, null, lofType, OrbitType.CARTESIAN, DEFAULT_POSITION_ANGLE);
  75.     }

  76.     /**
  77.      * Constructor.
  78.      * @param orbitalCovariance 6x6 orbital parameters covariance
  79.      * @param epoch epoch of the covariance
  80.      * @param covarianceFrame covariance frame (inertial or Earth fixed)
  81.      * @param orbitType orbit type of the covariance (CARTESIAN if covarianceFrame is not pseudo-inertial)
  82.      * @param angleType position angle type of the covariance (not used if orbitType is CARTESIAN)
  83.      */
  84.     public StateCovariance(final RealMatrix orbitalCovariance, final AbsoluteDate epoch,
  85.                            final Frame covarianceFrame,
  86.                            final OrbitType orbitType, final PositionAngle angleType) {
  87.         this(orbitalCovariance, epoch, covarianceFrame, null, orbitType, angleType);
  88.     }

  89.     /**
  90.      * Private constructor.
  91.      * @param orbitalCovariance 6x6 orbital parameters covariance
  92.      * @param epoch epoch of the covariance
  93.      * @param covarianceFrame covariance frame (inertial or Earth fixed)
  94.      * @param lofType covariance LOF type
  95.      * @param orbitType orbit type of the covariance
  96.      * @param angleType position angle type of the covariance (not used if orbitType is CARTESIAN)
  97.      */
  98.     private StateCovariance(final RealMatrix orbitalCovariance, final AbsoluteDate epoch,
  99.                             final Frame covarianceFrame, final LOFType lofType,
  100.                             final OrbitType orbitType, final PositionAngle angleType) {

  101.         checkInputConsistency(covarianceFrame, orbitType);

  102.         this.orbitalCovariance = orbitalCovariance;
  103.         this.epoch = epoch;
  104.         this.frame = covarianceFrame;
  105.         this.lofType = lofType;
  106.         this.orbitType = orbitType;
  107.         this.angleType = angleType;

  108.     }

  109.     /** Check constructor's inputs consistency.
  110.      *
  111.      * @param covarianceFrame covariance frame (inertial or Earth fixed)
  112.      * @param inputType orbit type of the covariance
  113.      */
  114.     private void checkInputConsistency(final Frame covarianceFrame, final OrbitType inputType) {

  115.         // State covariance expressed in a celestial body frame
  116.         if (covarianceFrame != null) {

  117.             // Input frame is pseudo-inertial
  118.             if (!covarianceFrame.isPseudoInertial() && inputType != OrbitType.CARTESIAN) {
  119.                 throw new OrekitException(OrekitMessages.WRONG_ORBIT_PARAMETERS_TYPE,
  120.                                           inputType.name(),
  121.                                           OrbitType.CARTESIAN.name());
  122.             }

  123.         }
  124.     }

  125.     /** {@inheritDoc}. */
  126.     @Override
  127.     public AbsoluteDate getDate() {
  128.         return epoch;
  129.     }

  130.     /**
  131.      * Get the covariance matrix.
  132.      * @return the covariance matrix
  133.      */
  134.     public RealMatrix getMatrix() {
  135.         return orbitalCovariance;
  136.     }

  137.     /**
  138.      * Get the covariance orbit type.
  139.      * @return the covariance orbit type
  140.      */
  141.     public OrbitType getOrbitType() {
  142.         return orbitType;
  143.     }

  144.     /**
  145.      * Get the covariance angle type.
  146.      * @return the covariance angle type
  147.      */
  148.     public PositionAngle getPositionAngle() {
  149.         return angleType;
  150.     }

  151.     /**
  152.      * Get the covariance frame.
  153.      * @return the covariance frame (can be null)
  154.      * @see #getLOFType()
  155.      */
  156.     public Frame getFrame() {
  157.         return frame;
  158.     }

  159.     /**
  160.      * Get the covariance LOF type.
  161.      * @return the covariance LOF type (can be null)
  162.      * @see #getFrame()
  163.      */
  164.     public LOFType getLOFType() {
  165.         return lofType;
  166.     }

  167.     /**
  168.      * Get the covariance matrix in another orbit type.
  169.      * <p>
  170.      * The covariance orbit type <b>cannot</b> be changed if the covariance
  171.      * matrix is expressed in a {@link LOFType local orbital frame} or a
  172.      * non-pseudo inertial frame.
  173.      * </p>
  174.      *
  175.      * @param orbit orbit to which the covariance matrix should correspond
  176.      * @param outOrbitType target orbit type of the state covariance matrix
  177.      * @param outAngleType target position angle type of the state covariance matrix
  178.      * @return a new covariance state, expressed in the target orbit type with the target position angle
  179.      * @see #changeCovarianceFrame(Orbit, Frame)
  180.      */
  181.     public StateCovariance changeCovarianceType(final Orbit orbit, final OrbitType outOrbitType,
  182.                                                 final PositionAngle outAngleType) {

  183.         // Check if the covariance expressed in a celestial body frame
  184.         if (frame != null) {

  185.             // Check if the covarianc is defined in inertial frame
  186.             if (frame.isPseudoInertial()) {
  187.                 return changeTypeAndCreate(orbit, epoch, frame, orbitType, angleType, outOrbitType, outAngleType,
  188.                                            orbitalCovariance);
  189.             }

  190.             // The covariance is not defined in inertial frame. The orbit type cannot be changes
  191.             throw new OrekitException(OrekitMessages.CANNOT_CHANGE_COVARIANCE_TYPE_IF_DEFINED_IN_NON_INERTIAL_FRAME);

  192.         }

  193.         // The covariance is not expressed in a celestial body frame. The orbit type cannot be changed
  194.         throw new OrekitException(OrekitMessages.CANNOT_CHANGE_COVARIANCE_TYPE_IF_DEFINED_IN_LOF);

  195.     }

  196.     /**
  197.      * Get the covariance in a given local orbital frame.
  198.      *
  199.      * @param orbit orbit to which the covariance matrix should correspond
  200.      * @param lofOut output local orbital frame
  201.      * @return a new covariance state, expressed in the output local orbital frame
  202.      */
  203.     public StateCovariance changeCovarianceFrame(final Orbit orbit, final LOFType lofOut) {

  204.         // Verify current covariance frame
  205.         if (lofType != null) {

  206.             // Change the covariance local orbital frame
  207.             return changeFrameAndCreate(orbit, epoch, lofType, lofOut, orbitalCovariance);

  208.         } else {

  209.             // Covariance is expressed in celestial body frame
  210.             return changeFrameAndCreate(orbit, epoch, frame, lofOut, orbitalCovariance, orbitType, angleType);

  211.         }

  212.     }

  213.     /**
  214.      * Get the covariance in the output frame.
  215.      *
  216.      * @param orbit orbit to which the covariance matrix should correspond
  217.      * @param frameOut output frame
  218.      * @return a new covariance state, expressed in the output frame
  219.      */
  220.     public StateCovariance changeCovarianceFrame(final Orbit orbit, final Frame frameOut) {

  221.         // Verify current covariance frame
  222.         if (lofType != null) {

  223.             // Covariance is expressed in local orbital frame
  224.             return changeFrameAndCreate(orbit, epoch, lofType, frameOut, orbitalCovariance);

  225.         } else {

  226.             // Change covariance frame
  227.             return changeFrameAndCreate(orbit, epoch, frame, frameOut, orbitalCovariance, orbitType, angleType);

  228.         }

  229.     }

  230.     /**
  231.      * Get a time-shifted covariance matrix.
  232.      * <p>
  233.      * The shifting model is a Keplerian one. In other words, the state transition matrix is computed supposing
  234.      * Keplerian motion.
  235.      * <p>
  236.      * Shifting is <em>not</em> intended as a replacement for proper covariance propagation, but should be sufficient
  237.      * for small time shifts or coarse accuracy.
  238.      *
  239.      * @param orbit orbit to which the covariance matrix should correspond
  240.      * @param dt time shift in seconds
  241.      * @return a new covariance state, shifted with respect to the instance
  242.      */
  243.     public StateCovariance shiftedBy(final Orbit orbit, final double dt) {

  244.         // Shifted orbit
  245.         final Orbit shifted = orbit.shiftedBy(dt);

  246.         // State covariance expressed in celestial body frame
  247.         if (frame != null) {

  248.             // State covariance expressed in a pseudo-inertial frame
  249.             if (frame.isPseudoInertial()) {

  250.                 // Compute STM
  251.                 final RealMatrix stm = getStm(orbit, dt);

  252.                 // Convert covariance in STM type (i.e., Equinoctial elements)
  253.                 final StateCovariance inStmType = changeTypeAndCreate(orbit, epoch, frame, orbitType, angleType,
  254.                                                                       OrbitType.EQUINOCTIAL, PositionAngle.MEAN,
  255.                                                                       orbitalCovariance);

  256.                 // Shift covariance by applying the STM
  257.                 final RealMatrix shiftedCov = stm.multiply(inStmType.getMatrix().multiplyTransposed(stm));

  258.                 // Restore the initial covariance type
  259.                 return changeTypeAndCreate(shifted, shifted.getDate(), frame,
  260.                                            OrbitType.EQUINOCTIAL, PositionAngle.MEAN,
  261.                                            orbitType, angleType, shiftedCov);
  262.             }

  263.             // State covariance expressed in a non pseudo-inertial frame
  264.             else {

  265.                 // Convert state covariance to orbit pseudo-inertial frame
  266.                 final StateCovariance inOrbitFrame = changeCovarianceFrame(orbit, orbit.getFrame());

  267.                 // Shift the state covariance
  268.                 final StateCovariance shiftedCovariance = inOrbitFrame.shiftedBy(orbit, dt);

  269.                 // Restore the initial covariance frame
  270.                 return shiftedCovariance.changeCovarianceFrame(shifted, frame);
  271.             }
  272.         }

  273.         // State covariance expressed in a commonly used local orbital frame (LOFType)
  274.         else {

  275.             // Convert state covariance to orbit pseudo-inertial frame
  276.             final StateCovariance inOrbitFrame = changeCovarianceFrame(orbit, orbit.getFrame());

  277.             // Shift the state covariance
  278.             final StateCovariance shiftedCovariance = inOrbitFrame.shiftedBy(orbit, dt);

  279.             // Restore the initial covariance frame
  280.             return shiftedCovariance.changeCovarianceFrame(shifted, lofType);
  281.         }

  282.     }

  283.     /**
  284.      * Create a covariance matrix in another orbit type.
  285.      *
  286.      * @param orbit orbit to which the covariance matrix should correspond
  287.      * @param date covariance epoch
  288.      * @param covFrame covariance frame
  289.      * @param inOrbitType initial orbit type of the state covariance matrix
  290.      * @param inAngleType initial position angle type of the state covariance matrix
  291.      * @param outOrbitType target orbit type of the state covariance matrix
  292.      * @param outAngleType target position angle type of the state covariance matrix
  293.      * @param inputCov input covariance
  294.      * @return the covariance expressed in the target orbit type with the target position angle
  295.      */
  296.     private static StateCovariance changeTypeAndCreate(final Orbit orbit, final AbsoluteDate date,
  297.                                                        final Frame covFrame,
  298.                                                        final OrbitType inOrbitType, final PositionAngle inAngleType,
  299.                                                        final OrbitType outOrbitType, final PositionAngle outAngleType,
  300.                                                        final RealMatrix inputCov) {

  301.         // Notations:
  302.         // I: Input orbit type
  303.         // O: Output orbit type
  304.         // C: Cartesian parameters

  305.         // Compute dOutputdCartesian
  306.         final double[][] aOC               = new double[STATE_DIMENSION][STATE_DIMENSION];
  307.         final Orbit      orbitInOutputType = outOrbitType.convertType(orbit);
  308.         orbitInOutputType.getJacobianWrtCartesian(outAngleType, aOC);
  309.         final RealMatrix dOdC = new Array2DRowRealMatrix(aOC, false);

  310.         // Compute dCartesiandInput
  311.         final double[][] aCI              = new double[STATE_DIMENSION][STATE_DIMENSION];
  312.         final Orbit      orbitInInputType = inOrbitType.convertType(orbit);
  313.         orbitInInputType.getJacobianWrtParameters(inAngleType, aCI);
  314.         final RealMatrix dCdI = new Array2DRowRealMatrix(aCI, false);

  315.         // Compute dOutputdInput
  316.         final RealMatrix dOdI = dOdC.multiply(dCdI);

  317.         // Conversion of the state covariance matrix in target orbit type
  318.         final RealMatrix outputCov = dOdI.multiply(inputCov.multiplyTransposed(dOdI));

  319.         // Return the converted covariance
  320.         return new StateCovariance(outputCov, date, covFrame, outOrbitType, outAngleType);

  321.     }

  322.     /**
  323.      * Create a covariance matrix from a {@link LOFType local orbital frame} to another
  324.      * {@link LOFType local orbital frame}.
  325.      * <p>
  326.      * The transformation is based on Equation (20) of "Covariance Transformations for Satellite Flight Dynamics
  327.      * Operations" by David A. Vallado.
  328.      *
  329.      * @param orbit orbit to which the covariance matrix should correspond
  330.      * @param date covariance epoch
  331.      * @param lofIn the local orbital frame in which the input covariance matrix is expressed
  332.      * @param lofOut the target local orbital frame
  333.      * @param inputCartesianCov input covariance {@code CARTESIAN})
  334.      * @return the covariance matrix expressed in the target commonly used local orbital frame in cartesian elements
  335.      */
  336.     private static StateCovariance changeFrameAndCreate(final Orbit orbit, final AbsoluteDate date,
  337.                                                         final LOFType lofIn, final LOFType lofOut,
  338.                                                         final RealMatrix inputCartesianCov) {

  339.         // Builds the matrix to perform covariance transformation
  340.         final RealMatrix jacobianFromLofInToLofOut =
  341.                         getJacobian(LOFType.transformFromLOFInToLOFOut(lofIn, lofOut, date, orbit.getPVCoordinates()));

  342.         // Get the Cartesian covariance matrix converted to frameOut
  343.         final RealMatrix cartesianCovarianceOut =
  344.                 jacobianFromLofInToLofOut.multiply(inputCartesianCov.multiplyTransposed(jacobianFromLofInToLofOut));

  345.         // Output converted covariance
  346.         return new StateCovariance(cartesianCovarianceOut, date, lofOut);

  347.     }

  348.     /**
  349.      * Convert the covariance matrix from a {@link Frame frame} to a {@link LOFType local orbital frame}.
  350.      * <p>
  351.      * The transformation is based on Equation (20) of "Covariance Transformations for Satellite Flight Dynamics
  352.      * Operations" by David A. Vallado.
  353.      * <p>
  354.      * As the frame transformation must be performed with the covariance expressed in Cartesian elements, both the orbit
  355.      * and position angle types of the input covariance must be provided.
  356.      * <p>
  357.      * <b>The output covariance matrix will provided in cartesian orbit type and not converted back to
  358.      * its original expression (if input different from cartesian elements).</b>
  359.      * <p>
  360.      *
  361.      * @param orbit orbit to which the covariance matrix should correspond
  362.      * @param date covariance epoch
  363.      * @param frameIn the frame in which the input covariance matrix is expressed. In case the frame is <b>not</b>
  364.      * pseudo-inertial, the input covariance matrix is expected to be expressed in <b>cartesian elements</b>.
  365.      * @param lofOut the target local orbital frame
  366.      * @param inputCov input covariance
  367.      * @param covOrbitType orbit type of the covariance matrix (used if frameIn is pseudo-inertial)
  368.      * @param covAngleType position angle type of the covariance matrix (used if frameIn is pseudo-inertial) (not used
  369.      * if covOrbitType equals {@code CARTESIAN})
  370.      * @return the covariance matrix expressed in the target local orbital frame in cartesian elements
  371.      * @throws OrekitException if input frame is <b>not</b> pseudo-inertial <b>and</b> the input covariance is
  372.      * <b>not</b> expressed in cartesian elements.
  373.      */
  374.     private static StateCovariance changeFrameAndCreate(final Orbit orbit, final AbsoluteDate date,
  375.                                                         final Frame frameIn, final LOFType lofOut,
  376.                                                         final RealMatrix inputCov,
  377.                                                         final OrbitType covOrbitType,
  378.                                                         final PositionAngle covAngleType) {

  379.         // Input frame is inertial
  380.         if (frameIn.isPseudoInertial()) {

  381.             // Convert input matrix to Cartesian parameters in input frame
  382.             final RealMatrix cartesianCovarianceIn =
  383.                     changeTypeAndCreate(orbit, date, frameIn, covOrbitType, covAngleType,
  384.                                         OrbitType.CARTESIAN, PositionAngle.MEAN,
  385.                                         inputCov).getMatrix();

  386.             // Builds the matrix to perform covariance transformation
  387.             final RealMatrix jacobianFromFrameInToLofOut =
  388.                             getJacobian(lofOut.transformFromInertial(date, orbit.getPVCoordinates(frameIn)));

  389.             // Get the Cartesian covariance matrix converted to frameOut
  390.             final RealMatrix cartesianCovarianceOut =
  391.                     jacobianFromFrameInToLofOut.multiply(cartesianCovarianceIn.multiplyTransposed(jacobianFromFrameInToLofOut));

  392.             // Return converted covariance matrix expressed in cartesian elements
  393.             return new StateCovariance(cartesianCovarianceOut, date, lofOut);

  394.         }

  395.         // Input frame is not inertial so the covariance matrix is expected to be in cartesian elements
  396.         else {

  397.             // Method checkInputConsistency already verify that input covariance is defined in CARTESIAN type
  398.             final Frame orbitInertialFrame = orbit.getFrame();

  399.             // Compute rotation matrix from frameIn to orbit inertial frame
  400.             final RealMatrix cartesianCovarianceInOrbitFrame =
  401.                    changeFrameAndCreate(orbit, date, frameIn, orbitInertialFrame, inputCov,
  402.                                          OrbitType.CARTESIAN, PositionAngle.MEAN).getMatrix();

  403.             // Convert from orbit inertial frame to lofOut
  404.             return changeFrameAndCreate(orbit, date, orbitInertialFrame, lofOut, cartesianCovarianceInOrbitFrame,
  405.                                         OrbitType.CARTESIAN, PositionAngle.MEAN);

  406.         }

  407.     }

  408.     /**
  409.      * Convert the covariance matrix from a {@link LOFType  local orbital frame} to a {@link Frame frame}.
  410.      * <p>
  411.      * The transformation is based on Equation (20) of "Covariance Transformations for Satellite Flight Dynamics
  412.      * Operations" by David A. Vallado.
  413.      * <p>
  414.      * The <u>input</u> covariance matrix is necessarily provided in <b>cartesian orbit type</b>.
  415.      * <p>
  416.      * The <u>output</u> covariance matrix will be expressed in <b>cartesian elements</b>.
  417.      *
  418.      * @param orbit orbit to which the covariance matrix should correspond
  419.      * @param date covariance epoch
  420.      * @param lofIn the local orbital frame in which the input covariance matrix is expressed
  421.      * @param frameOut the target frame
  422.      * @param inputCartesianCov input covariance ({@code CARTESIAN})
  423.      * @return the covariance matrix expressed in the target frame in cartesian elements
  424.      */
  425.     private static StateCovariance changeFrameAndCreate(final Orbit orbit, final AbsoluteDate date,
  426.                                                         final LOFType lofIn, final Frame frameOut,
  427.                                                         final RealMatrix inputCartesianCov) {

  428.         // Output frame is pseudo-inertial
  429.         if (frameOut.isPseudoInertial()) {

  430.             // Builds the matrix to perform covariance transformation
  431.             final RealMatrix jacobianFromLofInToFrameOut =
  432.                     getJacobian(lofIn.transformFromInertial(date, orbit.getPVCoordinates(frameOut)).getInverse());

  433.             // Transform covariance
  434.             final RealMatrix transformedCovariance =
  435.                     jacobianFromLofInToFrameOut.multiply(inputCartesianCov.multiplyTransposed(jacobianFromLofInToFrameOut));

  436.             // Get the Cartesian covariance matrix converted to frameOut
  437.             return new StateCovariance(transformedCovariance, date, frameOut, OrbitType.CARTESIAN,
  438.                                        DEFAULT_POSITION_ANGLE);

  439.         }

  440.         // Output frame is not pseudo-inertial
  441.         else {

  442.             // Builds the matrix to perform covariance transformation
  443.             final RealMatrix jacobianFromLofInToOrbitFrame =
  444.                     getJacobian(lofIn.transformFromInertial(date, orbit.getPVCoordinates()).getInverse());

  445.             // Get the Cartesian covariance matrix converted to orbit inertial frame
  446.             final RealMatrix cartesianCovarianceInOrbitFrame = jacobianFromLofInToOrbitFrame.multiply(
  447.                     inputCartesianCov.multiplyTransposed(jacobianFromLofInToOrbitFrame));

  448.             // Get the Cartesian covariance matrix converted to frameOut
  449.             return changeFrameAndCreate(orbit, date, orbit.getFrame(), frameOut, cartesianCovarianceInOrbitFrame,
  450.                                         OrbitType.CARTESIAN, PositionAngle.MEAN);
  451.         }

  452.     }

  453.     /**
  454.      * Get the covariance matrix in another frame.
  455.      * <p>
  456.      * The transformation is based on Equation (18) of "Covariance Transformations for Satellite Flight Dynamics
  457.      * Operations" by David A. Vallado.
  458.      * <p>
  459.      * As the frame transformation must be performed with the covariance expressed in Cartesian elements, both the orbit
  460.      * and position angle types of the input covariance must be provided.
  461.      * <p>
  462.      * In case the <u>input</u> frame is <b>not</b> pseudo-inertial, the <u>input</u> covariance matrix is necessarily
  463.      * expressed in <b>cartesian elements</b>.
  464.      * <p>
  465.      * In case the <u>output</u> frame is <b>not</b> pseudo-inertial, the <u>output</u> covariance matrix will be
  466.      * expressed in <b>cartesian elements</b>.
  467.      *
  468.      * @param orbit orbit to which the covariance matrix should correspond
  469.      * @param date covariance epoch
  470.      * @param frameIn the frame in which the input covariance matrix is expressed
  471.      * @param frameOut the target frame
  472.      * @param inputCov input covariance
  473.      * @param covOrbitType orbit type of the covariance matrix (used if frameIn is pseudo-inertial)
  474.      * @param covAngleType position angle type of the covariance matrix (used if frameIn is pseudo-inertial) (<b>not</b>
  475.      * used if covOrbitType equals {@code CARTESIAN})
  476.      * @return the covariance matrix expressed in the target frame
  477.      * @throws OrekitException if input frame is <b>not</b> pseudo-inertial <b>and</b> the input covariance is
  478.      * <b>not</b> expressed in cartesian elements.
  479.      */
  480.     private static StateCovariance changeFrameAndCreate(final Orbit orbit, final AbsoluteDate date,
  481.                                                         final Frame frameIn, final Frame frameOut,
  482.                                                         final RealMatrix inputCov,
  483.                                                         final OrbitType covOrbitType,
  484.                                                         final PositionAngle covAngleType) {

  485.         // Get the transform from the covariance frame to the output frame
  486.         final Transform inToOut = frameIn.getTransformTo(frameOut, orbit.getDate());

  487.         // Matrix to perform the covariance transformation
  488.         final RealMatrix j = getJacobian(inToOut);

  489.         // Input frame pseudo-inertial
  490.         if (frameIn.isPseudoInertial()) {

  491.             // Convert input matrix to Cartesian parameters in input frame
  492.             final RealMatrix cartesianCovarianceIn =
  493.                     changeTypeAndCreate(orbit, date, frameIn, covOrbitType, covAngleType,
  494.                                         OrbitType.CARTESIAN, PositionAngle.MEAN,
  495.                                         inputCov).getMatrix();

  496.             // Get the Cartesian covariance matrix converted to frameOut
  497.             final RealMatrix cartesianCovarianceOut = j.multiply(cartesianCovarianceIn.multiplyTransposed(j));

  498.             // Output frame is pseudo-inertial
  499.             if (frameOut.isPseudoInertial()) {

  500.                 // Convert output Cartesian matrix to initial orbit type and position angle
  501.                 return changeTypeAndCreate(orbit, date, frameOut, OrbitType.CARTESIAN, PositionAngle.MEAN,
  502.                                            covOrbitType, covAngleType, cartesianCovarianceOut);

  503.             }

  504.             // Output frame is not pseudo-inertial
  505.             else {

  506.                 // Output cartesian matrix
  507.                 return new StateCovariance(cartesianCovarianceOut, date, frameOut, OrbitType.CARTESIAN,
  508.                                            DEFAULT_POSITION_ANGLE);

  509.             }
  510.         }

  511.         // Input frame is not pseudo-inertial
  512.         else {

  513.             // Method checkInputConsistency already verify that input covariance is defined in CARTESIAN type

  514.             // Convert covariance matrix to frameOut
  515.             final RealMatrix covInFrameOut = j.multiply(inputCov.multiplyTransposed(j));

  516.             // Output the Cartesian covariance matrix converted to frameOut
  517.             return new StateCovariance(covInFrameOut, date, frameOut, OrbitType.CARTESIAN, DEFAULT_POSITION_ANGLE);

  518.         }

  519.     }

  520.     /**
  521.      * Builds the matrix to perform covariance frame transformation.
  522.      * @param transform input transformation
  523.      * @return the matrix to perform the covariance frame transformation
  524.      */
  525.     private static RealMatrix getJacobian(final Transform transform) {
  526.         // Get the Jacobian of the transform
  527.         final double[][] jacobian = new double[STATE_DIMENSION][STATE_DIMENSION];
  528.         transform.getJacobian(CartesianDerivativesFilter.USE_PV, jacobian);
  529.         // Return
  530.         return new Array2DRowRealMatrix(jacobian, false);

  531.     }

  532.     /**
  533.      * Get the state transition matrix considering Keplerian contribution only.
  534.      *
  535.      * @param initialOrbit orbit to which the initial covariance matrix should correspond
  536.      * @param dt time difference between the two orbits
  537.      * @return the state transition matrix used to shift the covariance matrix
  538.      */
  539.     private RealMatrix getStm(final Orbit initialOrbit, final double dt) {

  540.         // initialize the STM
  541.         final RealMatrix stm = MatrixUtils.createRealIdentityMatrix(STATE_DIMENSION);

  542.         // State transition matrix using Keplerian contribution only
  543.         final double mu           = initialOrbit.getMu();
  544.         final double sma          = initialOrbit.getA();
  545.         final double contribution = -1.5 * dt * FastMath.sqrt(mu / FastMath.pow(sma, 5));
  546.         stm.setEntry(5, 0, contribution);

  547.         // Return
  548.         return stm;

  549.     }

  550. }