diff --git a/src/main/java/org/orekit/errors/OrekitMessages.java b/src/main/java/org/orekit/errors/OrekitMessages.java index f8cf78f06ff107c4d0ad46b04e5c11105dd6a4c7..adeaa7fd875286e759c8b58fd01fb9d9e97a21a8 100644 --- a/src/main/java/org/orekit/errors/OrekitMessages.java +++ b/src/main/java/org/orekit/errors/OrekitMessages.java @@ -16,6 +16,8 @@ */ package org.orekit.errors; +import org.hipparchus.exception.Localizable; + import java.io.IOException; import java.io.InputStream; import java.io.InputStreamReader; @@ -27,8 +29,6 @@ import java.util.MissingResourceException; import java.util.PropertyResourceBundle; import java.util.ResourceBundle; -import org.hipparchus.exception.Localizable; - /** * Enumeration for localized messages formats. *

@@ -330,12 +330,15 @@ public enum OrekitMessages implements Localizable { INCOMPATIBLE_UNITS("units {0} and {1} are not compatible"), MISSING_VELOCITY("missing velocity data"), ATTEMPT_TO_GENERATE_MALFORMED_FILE("attempt to generate file {0} with a formatting error"), - FIND_ROOT("{0} failed to find root between {1} (g={2,number,0.0##############E0}) and {3} (g={4,number,0.0##############E0})\nLast iteration at {5} (g={6,number,0.0##############E0})"), + FIND_ROOT( + "{0} failed to find root between {1} (g={2,number,0.0##############E0}) and {3} (g={4,number,0.0##############E0})\nLast iteration at {5} (g={6,number,0.0##############E0})"), BACKWARD_PROPAGATION_NOT_ALLOWED("backward propagation not allowed here"), - NO_STATION_ECCENTRICITY_FOR_EPOCH("no station eccentricity values for the given epoch {0}, validity interval is between {1} and {2}"), + NO_STATION_ECCENTRICITY_FOR_EPOCH( + "no station eccentricity values for the given epoch {0}, validity interval is between {1} and {2}"), INCONSISTENT_SELECTION("inconsistent parameters selection between pairs {0}/{1} and {2}/{3}"), NO_UNSCENTED_TRANSFORM_CONFIGURED("no unscented transform configured"), - NOT_STRICTLY_POSITIVE("value is not strictly positive: {0}"); + NOT_STRICTLY_POSITIVE("value is not strictly positive: {0}"), + WRONG_ORBIT_PARAMETERS_TYPE("orbital parameters type: {0} is different from expected orbital type : {1}"); // CHECKSTYLE: resume JavadocVariable check diff --git a/src/main/java/org/orekit/frames/LOFType.java b/src/main/java/org/orekit/frames/LOFType.java index 54e4f8c6ce6befec7c05691d81f60e69c59f3125..988611338063167d32caf81a42dd97858b98d208 100644 --- a/src/main/java/org/orekit/frames/LOFType.java +++ b/src/main/java/org/orekit/frames/LOFType.java @@ -16,18 +16,21 @@ */ package org.orekit.frames; -import org.hipparchus.Field; import org.hipparchus.CalculusFieldElement; +import org.hipparchus.Field; import org.hipparchus.geometry.euclidean.threed.FieldRotation; import org.hipparchus.geometry.euclidean.threed.FieldVector3D; import org.hipparchus.geometry.euclidean.threed.Rotation; +import org.hipparchus.geometry.euclidean.threed.RotationConvention; import org.hipparchus.geometry.euclidean.threed.Vector3D; import org.orekit.time.AbsoluteDate; import org.orekit.time.FieldAbsoluteDate; import org.orekit.utils.FieldPVCoordinates; import org.orekit.utils.PVCoordinates; -/** Enumerate for different types of Local Orbital Frames. +/** + * Enumerate for different types of Local Orbital Frames. + * * @author Luc Maisonobe */ public enum LOFType { @@ -47,16 +50,23 @@ public enum LOFType { * @see #NTW */ TNW { - /** {@inheritDoc} */ + @Override public Rotation rotationFromInertial(final PVCoordinates pv) { return new Rotation(pv.getVelocity(), pv.getMomentum(), Vector3D.PLUS_I, Vector3D.PLUS_K); } /** {@inheritDoc} */ + @Override + public Rotation rotationFromLOFType(final LOFType fromLOF, final PVCoordinates pv) { + return rotationFromLOFInToLOFOut(fromLOF, TNW, pv); + } + + /** {@inheritDoc} */ + @Override public > FieldRotation rotationFromInertial(final Field field, - final FieldPVCoordinates pv) { + final FieldPVCoordinates pv) { return new FieldRotation<>(pv.getVelocity(), pv.getMomentum(), new FieldVector3D<>(field, Vector3D.PLUS_I), new FieldVector3D<>(field, Vector3D.PLUS_K)); @@ -81,16 +91,81 @@ public enum LOFType { * @see #VVLH */ QSW { - /** {@inheritDoc} */ + @Override public Rotation rotationFromInertial(final PVCoordinates pv) { return new Rotation(pv.getPosition(), pv.getMomentum(), Vector3D.PLUS_I, Vector3D.PLUS_K); } /** {@inheritDoc} */ + @Override + public Rotation rotationFromLOFType(final LOFType fromLOF, final PVCoordinates pv) { + + // Compute the common rotation to all cases + final Rotation rotationFromECIToQSW = LOFType.QSW.rotationFromInertial(pv); + + switch (fromLOF) { + case EQW: + final Rotation rotationFromEQWToECI = LOFType.EQW.rotationFromInertial(pv); + + final Rotation rotationFromEQWToQSW = + rotationFromEQWToECI.compose(rotationFromECIToQSW, RotationConvention.FRAME_TRANSFORM); + + return rotationFromEQWToQSW; + case NTW: + final Rotation rotationFromNTWToECI = LOFType.NTW.rotationFromInertial(pv).revert(); + + final Rotation rotationFromNTWToQSW = + rotationFromNTWToECI.compose(rotationFromECIToQSW, RotationConvention.FRAME_TRANSFORM); + + return rotationFromNTWToQSW; + + case VNC: + final Rotation rotationFromVNCToECI = LOFType.VNC.rotationFromInertial(pv); + + final Rotation rotationFromVNCToQSW = + rotationFromVNCToECI.compose(rotationFromECIToQSW, RotationConvention.FRAME_TRANSFORM); + + return rotationFromVNCToQSW; + + case LVLH: + return Rotation.IDENTITY; + + case VVLH: + final Rotation rotationFromVVLHToECI = LOFType.VVLH.rotationFromInertial(pv); + + final Rotation rotationFromVVLHToQSW = + rotationFromVVLHToECI.compose(rotationFromECIToQSW, RotationConvention.FRAME_TRANSFORM); + + return rotationFromVVLHToQSW; + + case LVLH_CCSDS: + final Rotation rotationFromLVLH_CCSDSToECI = LOFType.LVLH_CCSDS.rotationFromInertial(pv); + + final Rotation rotationFromLVLH_CCSDSToQSW = + rotationFromLVLH_CCSDSToECI.compose(rotationFromECIToQSW, + RotationConvention.FRAME_TRANSFORM); + + return rotationFromLVLH_CCSDSToQSW; + + case TNW: + final Rotation rotationFromTNWToECI = LOFType.TNW.rotationFromInertial(pv); + + final Rotation rotationFromTNWToQSW = + rotationFromTNWToECI.compose(rotationFromECIToQSW, RotationConvention.FRAME_TRANSFORM); + + return rotationFromTNWToQSW; + + default: + return Rotation.IDENTITY; + } + } + + /** {@inheritDoc} */ + @Override public > FieldRotation rotationFromInertial(final Field field, - final FieldPVCoordinates pv) { + final FieldPVCoordinates pv) { return new FieldRotation<>(pv.getPosition(), pv.getMomentum(), new FieldVector3D<>(field, Vector3D.PLUS_I), new FieldVector3D<>(field, Vector3D.PLUS_K)); @@ -121,16 +196,23 @@ public enum LOFType { * @see #VVLH */ LVLH { - /** {@inheritDoc} */ + @Override public Rotation rotationFromInertial(final PVCoordinates pv) { return new Rotation(pv.getPosition(), pv.getMomentum(), Vector3D.PLUS_I, Vector3D.PLUS_K); } /** {@inheritDoc} */ + @Override + public Rotation rotationFromLOFType(final LOFType fromLOF, final PVCoordinates pv) { + return rotationFromLOFInToLOFOut(fromLOF, LVLH, pv); + } + + /** {@inheritDoc} */ + @Override public > FieldRotation rotationFromInertial(final Field field, - final FieldPVCoordinates pv) { + final FieldPVCoordinates pv) { return new FieldRotation<>(pv.getPosition(), pv.getMomentum(), new FieldVector3D<>(field, Vector3D.PLUS_I), new FieldVector3D<>(field, Vector3D.PLUS_K)); @@ -159,16 +241,23 @@ public enum LOFType { * @since 11.0 */ LVLH_CCSDS { - /** {@inheritDoc} */ + @Override public Rotation rotationFromInertial(final PVCoordinates pv) { return new Rotation(pv.getPosition(), pv.getMomentum(), Vector3D.MINUS_K, Vector3D.MINUS_J); } /** {@inheritDoc} */ + @Override + public Rotation rotationFromLOFType(final LOFType fromLOF, final PVCoordinates pv) { + return rotationFromLOFInToLOFOut(fromLOF, LVLH_CCSDS, pv); + } + + /** {@inheritDoc} */ + @Override public > FieldRotation rotationFromInertial(final Field field, - final FieldPVCoordinates pv) { + final FieldPVCoordinates pv) { return new FieldRotation<>(pv.getPosition(), pv.getMomentum(), new FieldVector3D<>(field, Vector3D.MINUS_K), new FieldVector3D<>(field, Vector3D.MINUS_J)); @@ -196,15 +285,22 @@ public enum LOFType { * @see #LVLH_CCSDS */ VVLH { - /** {@inheritDoc} */ + @Override public Rotation rotationFromInertial(final PVCoordinates pv) { return LVLH_CCSDS.rotationFromInertial(pv); } /** {@inheritDoc} */ + @Override + public Rotation rotationFromLOFType(final LOFType fromLOF, final PVCoordinates pv) { + return rotationFromLOFInToLOFOut(fromLOF, VVLH, pv); + } + + /** {@inheritDoc} */ + @Override public > FieldRotation rotationFromInertial(final Field field, - final FieldPVCoordinates pv) { + final FieldPVCoordinates pv) { return LVLH_CCSDS.rotationFromInertial(field, pv); } @@ -225,16 +321,23 @@ public enum LOFType { * @see #NTW */ VNC { - /** {@inheritDoc} */ + @Override public Rotation rotationFromInertial(final PVCoordinates pv) { return new Rotation(pv.getVelocity(), pv.getMomentum(), Vector3D.PLUS_I, Vector3D.PLUS_J); } + /** {@inheritDoc} */ + @Override + public Rotation rotationFromLOFType(final LOFType fromLOF, final PVCoordinates pv) { + return rotationFromLOFInToLOFOut(fromLOF, VNC, pv); + } + + /** {@inheritDoc} */ @Override public > FieldRotation rotationFromInertial(final Field field, - final FieldPVCoordinates pv) { + final FieldPVCoordinates pv) { return new FieldRotation<>(pv.getVelocity(), pv.getMomentum(), new FieldVector3D<>(field, Vector3D.PLUS_I), new FieldVector3D<>(field, Vector3D.PLUS_J)); @@ -247,17 +350,24 @@ public enum LOFType { * @since 11.0 */ EQW { - /** {@inheritDoc} */ + @Override public Rotation rotationFromInertial(final PVCoordinates pv) { final Vector3D m = pv.getMomentum(); return new Rotation(new Vector3D(-m.getY(), m.getX(), 0), m, Vector3D.PLUS_I, Vector3D.PLUS_J); } + /** {@inheritDoc} */ + @Override + public Rotation rotationFromLOFType(final LOFType fromLOF, final PVCoordinates pv) { + return rotationFromLOFInToLOFOut(fromLOF, EQW, pv); + } + + /** {@inheritDoc} */ @Override public > FieldRotation rotationFromInertial(final Field field, - final FieldPVCoordinates pv) { + final FieldPVCoordinates pv) { final FieldVector3D m = pv.getMomentum(); return new FieldRotation<>(new FieldVector3D<>(m.getY().negate(), m.getX(), field.getZero()), m, @@ -282,16 +392,23 @@ public enum LOFType { * @since 11.0 */ NTW { - /** {@inheritDoc} */ + @Override public Rotation rotationFromInertial(final PVCoordinates pv) { return new Rotation(pv.getVelocity(), pv.getMomentum(), Vector3D.PLUS_J, Vector3D.PLUS_K); } + /** {@inheritDoc} */ + @Override + public Rotation rotationFromLOFType(final LOFType fromLOF, final PVCoordinates pv) { + return rotationFromLOFInToLOFOut(fromLOF, NTW, pv); + } + + /** {@inheritDoc} */ @Override public > FieldRotation rotationFromInertial(final Field field, - final FieldPVCoordinates pv) { + final FieldPVCoordinates pv) { return new FieldRotation<>(pv.getVelocity(), pv.getMomentum(), new FieldVector3D<>(field, Vector3D.PLUS_J), new FieldVector3D<>(field, Vector3D.PLUS_K)); @@ -299,9 +416,31 @@ public enum LOFType { }; - /** Get the transform from an inertial frame defining position-velocity and the local orbital frame. - * @param date current date + /** + * Get the rotation from input to output {@link LOFType commonly used local orbital frame}. + * + * @param in input commonly used local orbital frame + * @param out output commonly used local orbital frame * @param pv position-velocity of the spacecraft in some inertial frame + * @return rotation from input to output {@link LOFType commonly used local orbital frame}. + */ + public static Rotation rotationFromLOFInToLOFOut(final LOFType in, final LOFType out, final PVCoordinates pv) { + + // First compute the rotation from the input LOF to the pivot LOF QSW + final Rotation inToQSW = LOFType.QSW.rotationFromLOFType(in, pv); + + // Then compute the rotation from the pivot LOF QSW to the output LOF + final Rotation QSWToOut = LOFType.QSW.rotationFromLOFType(out, pv).revert(); + + // Output composed rotation + return inToQSW.compose(QSWToOut, RotationConvention.FRAME_TRANSFORM); + } + + /** + * Get the transform from an inertial frame defining position-velocity and the local orbital frame. + * + * @param date current date + * @param pv position-velocity of the spacecraft in some inertial frame * @return transform from the frame where position-velocity are defined to local orbital frame */ public Transform transformFromInertial(final AbsoluteDate date, final PVCoordinates pv) { @@ -310,8 +449,8 @@ public enum LOFType { final Transform translation = new Transform(date, pv.negate()); // compute the rotation part of the transform - final Rotation r = rotationFromInertial(pv); - final Vector3D p = pv.getPosition(); + final Rotation r = rotationFromInertial(pv); + final Vector3D p = pv.getPosition(); final Vector3D momentum = pv.getMomentum(); final Transform rotation = new Transform(date, r, new Vector3D(1.0 / p.getNormSq(), r.applyTo(momentum))); @@ -320,22 +459,37 @@ public enum LOFType { } - /** Get the transform from an inertial frame defining position-velocity and the local orbital frame. - * @param date current date + /** + * Get the rotation from inertial frame to local orbital frame. + *

+ * This rotation does not include any time derivatives. If first time derivatives (i.e. rotation rate) is needed as + * well, the full {@link #transformFromInertial(AbsoluteDate, PVCoordinates) transformFromInertial} method must be + * called and the complete rotation transform must be extracted from it. + *

+ * * @param pv position-velocity of the spacecraft in some inertial frame - * @param type of the fiels elements + * @return rotation from inertial frame to local orbital frame + */ + public abstract Rotation rotationFromInertial(PVCoordinates pv); + + /** + * Get the transform from an inertial frame defining position-velocity and the local orbital frame. + * + * @param date current date + * @param pv position-velocity of the spacecraft in some inertial frame + * @param type of the fields elements * @return transform from the frame where position-velocity are defined to local orbital frame * @since 9.0 */ public > FieldTransform transformFromInertial(final FieldAbsoluteDate date, - final FieldPVCoordinates pv) { + final FieldPVCoordinates pv) { // compute the translation part of the transform final FieldTransform translation = new FieldTransform<>(date, pv.negate()); // compute the rotation part of the transform - final FieldRotation r = rotationFromInertial(date.getField(), pv); - final FieldVector3D p = pv.getPosition(); + final FieldRotation r = rotationFromInertial(date.getField(), pv); + final FieldVector3D p = pv.getPosition(); final FieldVector3D momentum = pv.getMomentum(); final FieldTransform rotation = new FieldTransform(date, r, new FieldVector3D<>(p.getNormSq().reciprocal(), r.applyTo(momentum))); @@ -344,18 +498,7 @@ public enum LOFType { } - /** Get the rotation from inertial frame to local orbital frame. - *

- * This rotation does not include any time derivatives. If first - * time derivatives (i.e. rotation rate) is needed as well, the full - * {@link #transformFromInertial(AbsoluteDate, PVCoordinates) transformFromInertial} - * method must be called and the complete rotation transform must be extracted - * from it. - *

- * @param pv position-velocity of the spacecraft in some inertial frame - * @return rotation from inertial frame to local orbital frame - */ - public abstract Rotation rotationFromInertial(PVCoordinates pv); + public abstract Rotation rotationFromLOFType(LOFType fromLOF, PVCoordinates pv); /** Get the rotation from inertial frame to local orbital frame. *

diff --git a/src/main/java/org/orekit/propagation/StateCovarianceMatrixProvider.java b/src/main/java/org/orekit/propagation/StateCovarianceMatrixProvider.java index c342d76301686e577a68ea6fcd35214a0f3fbac1..4de49bcd263507d66ec045763b8c6643cf0e19a6 100644 --- a/src/main/java/org/orekit/propagation/StateCovarianceMatrixProvider.java +++ b/src/main/java/org/orekit/propagation/StateCovarianceMatrixProvider.java @@ -16,12 +16,15 @@ */ package org.orekit.propagation; +import org.hipparchus.geometry.euclidean.threed.Rotation; import org.hipparchus.linear.Array2DRowRealMatrix; import org.hipparchus.linear.MatrixUtils; import org.hipparchus.linear.RealMatrix; +import org.orekit.errors.OrekitException; +import org.orekit.errors.OrekitMessages; import org.orekit.frames.Frame; +import org.orekit.frames.LOFType; import org.orekit.frames.Transform; -import org.orekit.orbits.CartesianOrbit; import org.orekit.orbits.Orbit; import org.orekit.orbits.OrbitType; import org.orekit.orbits.PositionAngle; @@ -31,40 +34,32 @@ import org.orekit.utils.CartesianDerivativesFilter; /** * Additional state provider for state covariance matrix. *

- * This additional state provider allows computing a propagated - * covariance matrix based on a user defined input state covariance - * matrix. The computation of the propagated covariance matrix uses - * the State Transition Matrix between the propagated spacecraft state - * and the initial state. As a result, the user must define the name + * This additional state provider allows computing a propagated covariance matrix based on a user defined input state + * covariance matrix. The computation of the propagated covariance matrix uses the State Transition Matrix between the + * propagated spacecraft state and the initial state. As a result, the user must define the name * {@link #stmName of the provider for the State Transition Matrix}. *

- * As the State Transition Matrix and the input state covariance - * matrix can be expressed in different orbit types, the user must - * specify both orbit types when building the covariance provider. - * In addition, the position angle used in both matrices must also - * be specified. + * As the State Transition Matrix and the input state covariance matrix can be expressed in different orbit types, the + * user must specify both orbit types when building the covariance provider. In addition, the position angle used in + * both matrices must also be specified. *

- * In order to add this additional state provider to an orbit - * propagator, user must use the - * {@link Propagator#addAdditionalStateProvider(AdditionalStateProvider)} - * method. + * In order to add this additional state provider to an orbit propagator, user must use the + * {@link Propagator#addAdditionalStateProvider(AdditionalStateProvider)} method. *

- * For a given propagated spacecraft {@code state}, the propagated state - * covariance matrix is accessible through the method - * {@link #getStateCovariance(SpacecraftState)} + * For a given propagated spacecraft {@code state}, the propagated state covariance matrix is accessible through the + * method {@link #getStateCovariance(SpacecraftState)} *

* It is possible to change the covariance frame by using the - * {@link #changeCovarianceFrame(Orbit, Frame, Frame, RealMatrix, OrbitType, PositionAngle)} - * method. This method is based on Equation (18) of Covariance - * Transformations for Satellite Flight Dynamics Operations - * by David A. Vallado. It is important to highlight that the frames - * must be inertial frames. + * {@link #changeCovarianceFrame(Orbit, Frame, Frame, RealMatrix, OrbitType, PositionAngle)} method. This method is + * based on Equation (18) of Covariance Transformations for Satellite Flight Dynamics Operations by David A. + * Vallado. *

* Finally, covariance orbit type can be changed using the - * {@link #changeCovarianceType(Orbit, OrbitType, PositionAngle, OrbitType, PositionAngle, RealMatrix)} - * method. + * {@link #changeCovarianceType(Orbit, OrbitType, PositionAngle, OrbitType, PositionAngle, RealMatrix)} method. *

+ * * @author Bryan Cazabonne + * @author Vincent Cucchietti * @since 11.3 */ public class StateCovarianceMatrixProvider implements AdditionalStateProvider { @@ -98,16 +93,18 @@ public class StateCovarianceMatrixProvider implements AdditionalStateProvider { /** * Constructor. + * * @param additionalName name of the additional state * @param stmName name of the state for State Transition Matrix - * @param harvester matrix harvester as returned by {@code propagator.setupMatricesComputation(stmName, null, null)} + * @param harvester matrix harvester as returned by + * {@code propagator.setupMatricesComputation(stmName, null, null)} * @param stmOrbitType orbit type used for the State Transition Matrix computation - * @param stmAngleType position angle used for State Transition Matrix computation - * (not used if stmOrbitType equals {@code CARTESIAN}) + * @param stmAngleType position angle used for State Transition Matrix computation (not used if stmOrbitType equals + * {@code CARTESIAN}) * @param covInit initial state covariance matrix (6x6 dimension) * @param covOrbitType orbit type for the covariance matrix - * @param covAngleType position angle used for the covariance matrix - * (not used if covOrbitType equals {@code CARTESIAN}) + * @param covAngleType position angle used for the covariance matrix (not used if covOrbitType equals + * {@code CARTESIAN}) */ public StateCovarianceMatrixProvider(final String additionalName, final String stmName, final MatricesHarvester harvester, @@ -116,20 +113,51 @@ public class StateCovarianceMatrixProvider implements AdditionalStateProvider { final OrbitType covOrbitType, final PositionAngle covAngleType) { // Initialize fields this.additionalName = additionalName; - this.stmName = stmName; - this.harvester = harvester; - this.covInit = covInit; - this.stmOrbitType = stmOrbitType; - this.stmAngleType = stmAngleType; - this.covOrbitType = covOrbitType; - this.covAngleType = covAngleType; + this.stmName = stmName; + this.harvester = harvester; + this.covInit = covInit; + this.stmOrbitType = stmOrbitType; + this.stmAngleType = stmAngleType; + this.covOrbitType = covOrbitType; + this.covAngleType = covAngleType; } - /** Get the orbit type in which the covariance matrix is expressed. - * @return the orbit type + /** + * Convert the covariance matrix from a {@link LOFType commonly used local orbital frame} to another + * {@link LOFType commonly used local orbital frame}. + *

+ * The transformation is based on Equation (20) of "Covariance Transformations for Satellite Flight Dynamics + * Operations" by David A. Vallado". + *

+ * As this method transforms from and to a {@link LOFType commonly used local orbital frame}, it necessarily takes + * in a covariance matrix expressed in cartesian elements and output a covariance matrix also expressed in + * the + * cartesian elements. + * + * @param orbit orbit to which the covariance matrix should correspond + * @param lofIn the local orbital frame in which the input covariance matrix is expressed + * @param lofOut the target local orbital frame + * @param inputCartesianCov input covariance {@code CARTESIAN}) + * @return the covariance matrix expressed in the target commonly used local orbital frame in cartesian elements */ - public OrbitType getCovarianceOrbitType() { - return covOrbitType; + public static RealMatrix changeCovarianceFrame(final Orbit orbit, + final LOFType lofIn, final LOFType lofOut, + final RealMatrix inputCartesianCov) { + + // Compute rotation matrix from lofIn to lofOut + final Rotation rotationFromLofInToLofOut = LOFType.rotationFromLOFInToLOFOut(lofIn, lofOut, + orbit.getPVCoordinates()); + + // Builds the matrix to perform covariance transformation + final RealMatrix transformationMatrix = buildTransformationMatrixFromRotation(rotationFromLofInToLofOut); + + // Get the Cartesian covariance matrix converted to frameOut + final RealMatrix cartesianCovarianceOut = + transformationMatrix.multiply(inputCartesianCov.multiplyTransposed(transformationMatrix)); + + // Output converted covariance + return cartesianCovarianceOut; + } /** {@inheritDoc} */ @@ -148,14 +176,18 @@ public class StateCovarianceMatrixProvider implements AdditionalStateProvider { covInit); } - /** {@inheritDoc} - *

- * The covariance matrix can be computed only if the State Transition Matrix state is available. - *

- */ - @Override - public boolean yield(final SpacecraftState state) { - return !state.hasAdditionalState(stmName); + private static RealMatrix buildTransformationMatrixFromRotation(final Rotation rotation) { + + final double[][] rotationMatrixData = rotation.getMatrix(); + + final RealMatrix transformationMatrix = MatrixUtils.createRealMatrix(6, 6); + + // Fills in the upper left and lower right blocks with the rotation + transformationMatrix.setSubMatrix(rotationMatrixData, 0, 0); + transformationMatrix.setSubMatrix(rotationMatrixData, 3, 3); + + return transformationMatrix; + } /** {@inheritDoc} */ @@ -178,105 +210,101 @@ public class StateCovarianceMatrixProvider implements AdditionalStateProvider { } - /** Get the state covariance matrix (6x6 dimension). + /** + * Convert the covariance matrix from a {@link LOFType commonly used local orbital frame} to a {@link Frame frame}. *

- * The output covariance matrix is expressed in the - * same orbit type as {@link #getCovarianceOrbitType()}. + * The transformation is based on Equation (20) of "Covariance Transformations for Satellite Flight Dynamics + * Operations" by David A. Vallado". *

- * It is possible to change the covariance frame by using the - * {@link #changeCovarianceFrame(Orbit, Frame, Frame, RealMatrix, OrbitType, PositionAngle)} - * method. + * The input covariance matrix is necessarily expressed in cartesian elements. *

- * It is also possible to change the covariance orbit type by - * using the {@link #changeCovarianceType(Orbit, OrbitType, PositionAngle, OrbitType, PositionAngle, RealMatrix)} - * method. - * @param state spacecraft state to which the covariance matrix should correspond - * @return the state covariance matrix - * @see #getStateCovariance(SpacecraftState, Frame) - * @see #getStateCovariance(SpacecraftState, OrbitType, PositionAngle) + * The output covariance matrix will be expressed in cartesian elements. + * + * @param orbit orbit to which the covariance matrix should correspond + * @param lofIn the local orbital frame in which the input covariance matrix is expressed + * @param frameOut the target frame + * @param inputCartesianCov input covariance ({@code CARTESIAN}) + * @return the covariance matrix expressed in the target frame in cartesian elements */ - public RealMatrix getStateCovariance(final SpacecraftState state) { + public static RealMatrix changeCovarianceFrame(final Orbit orbit, + final LOFType lofIn, final Frame frameOut, + final RealMatrix inputCartesianCov) { - // Get the current propagated covariance - final RealMatrix covariance = toRealMatrix(state.getAdditionalState(additionalName)); + // Initialize output variable + final RealMatrix cartesianCovarianceOut; - // Return the converted covariance - return covariance; + // Get the orbit inertial frame + final Frame orbitInertialFrame = orbit.getFrame(); - } + if (frameOut.isPseudoInertial()) { + // Compute rotation matrix from lofIn to frameOut + final Rotation rotationFromLofInToFrameOut = + lofIn.rotationFromInertial(orbit.getPVCoordinates(frameOut)).revert(); - /** Get the state covariance matrix (6x6 dimension) expressed in a given frame. - *

- * The output covariance matrix is expressed in the - * same orbit type as {@link #getCovarianceOrbitType()}. - *

- * It is also possible to change the covariance orbit type by - * using the {@link #changeCovarianceType(Orbit, OrbitType, PositionAngle, OrbitType, PositionAngle, RealMatrix)} - * method. - * @param state spacecraft state to which the covariance matrix should correspond - * @param frame output frame for which the output covariance matrix must be expressed - * (must be inertial) - * @return the state covariance matrix expressed in frame - * @see #getStateCovariance(SpacecraftState) - * @see #getStateCovariance(SpacecraftState, OrbitType, PositionAngle) - */ - public RealMatrix getStateCovariance(final SpacecraftState state, final Frame frame) { + // Builds the matrix to perform covariance transformation + final RealMatrix transformationMatrix = + buildTransformationMatrixFromRotation(rotationFromLofInToFrameOut); - // Get the current propagated covariance - final RealMatrix covariance = toRealMatrix(state.getAdditionalState(additionalName)); - - // Return the converted covariance - return changeCovarianceFrame(state.getOrbit(), state.getFrame(), frame, covariance, covOrbitType, covAngleType); - - } - - /** Get the state covariance matrix (6x6 dimension) expressed in a given orbit type. - * @param state spacecraft state to which the covariance matrix should correspond - * @param orbitType output orbit type - * @param angleType output position angle - * (not used if orbitType equals {@code CARTESIAN}) - * @return the state covariance matrix in orbitType and angleType - * @see #getStateCovariance(SpacecraftState) - * @see #getStateCovariance(SpacecraftState, Frame) - */ - public RealMatrix getStateCovariance(final SpacecraftState state, final OrbitType orbitType, final PositionAngle angleType) { + // Get the Cartesian covariance matrix converted to frameOut + cartesianCovarianceOut = + transformationMatrix.multiply(inputCartesianCov.multiplyTransposed(transformationMatrix)); - // Get the current propagated covariance - final RealMatrix covariance = toRealMatrix(state.getAdditionalState(additionalName)); + } + else { + // Compute rotation matrix from lofIn to orbit inertial frame + final Rotation rotationFromLofInToOrbitFrame = + lofIn.rotationFromInertial(orbit.getPVCoordinates()).revert(); + + // Builds the matrix to perform covariance transformation + final RealMatrix transformationMatrixFromLofInToOrbitFrame = + buildTransformationMatrixFromRotation(rotationFromLofInToOrbitFrame); + + // Get the Cartesian covariance matrix converted to orbit inertial frame + final RealMatrix cartesianCovarianceInOrbitFrame = transformationMatrixFromLofInToOrbitFrame.multiply( + inputCartesianCov.multiplyTransposed(transformationMatrixFromLofInToOrbitFrame)); + + // Get the Cartesian covariance matrix converted to frameOut + cartesianCovarianceOut = + changeCovarianceFrame(orbit, orbitInertialFrame, frameOut, cartesianCovarianceInOrbitFrame, + OrbitType.CARTESIAN, PositionAngle.MEAN); + } - // Return the converted covariance - return changeCovarianceType(state.getOrbit(), covOrbitType, covAngleType, orbitType, angleType, covariance); + // Output converted covariance + return cartesianCovarianceOut; } - /** Get the covariance matrix in another frame. + /** + * Get the covariance matrix in another frame. *

- * The transformation is based on Equation (18) of - * "Covariance Transformations for Satellite Flight Dynamics Operations" - * by David A. Vallado". + * The transformation is based on Equation (18) of "Covariance Transformations for Satellite Flight Dynamics + * Operations" by David A. Vallado". *

- * As the frame transformation must be performed with the covariance - * expressed in Cartesian elements, both the orbit and position angle - * types of the input covariance must be provided. + * As the frame transformation must be performed with the covariance expressed in Cartesian elements, both the orbit + * and position angle types of the input covariance must be provided. + *

+ * In case the input frame is not pseudo-inertial, the input covariance matrix is necessarily + * expressed in cartesian elements. + *

+ * In case the output frame is not pseudo-inertial, the output covariance matrix will be + * expressed in cartesian elements. + * * @param orbit orbit to which the covariance matrix should correspond - * @param frameIn the frame in which the input covariance matrix is expressed (must be inertial) - * @param frameOut the target frame (must be inertial) + * @param frameIn the frame in which the input covariance matrix is expressed + * @param frameOut the target frame * @param inputCov input covariance - * @param covOrbitType orbit type of the covariance matrix - * @param covAngleType position angle type of the covariance matrix - * (not used if covOrbitType equals {@code CARTESIAN}) + * @param covOrbitType orbit type of the covariance matrix (used if frameIn is pseudo-inertial) + * @param covAngleType position angle type of the covariance matrix (used if frameIn is pseudo-inertial) (not + * used if covOrbitType equals {@code CARTESIAN}) * @return the covariance matrix expressed in the target frame + * @throws OrekitException if input frame is not pseudo-inertial and the input covariance is + * not expressed in cartesian elements. */ public static RealMatrix changeCovarianceFrame(final Orbit orbit, final Frame frameIn, final Frame frameOut, final RealMatrix inputCov, final OrbitType covOrbitType, final PositionAngle covAngleType) { - // Convert input matrix to Cartesian parameters in input frame - final RealMatrix cartesianCovarianceIn = changeCovarianceType(orbit, covOrbitType, covAngleType, - OrbitType.CARTESIAN, PositionAngle.MEAN, - inputCov); - // Get the transform from the covariance frame to the output frame final Transform inToOut = frameIn.getTransformTo(frameOut, orbit.getDate()); @@ -287,27 +315,62 @@ public class StateCovarianceMatrixProvider implements AdditionalStateProvider { // Matrix to perform the covariance transformation final RealMatrix j = new Array2DRowRealMatrix(jacobian, false); - // Get the Cartesian covariance matrix converted to frameOut - final RealMatrix cartesianCovarianceOut = j.multiply(cartesianCovarianceIn.multiplyTransposed(j)); + // Input frame pseudo-inertial + if (frameIn.isPseudoInertial()) { + + // Convert input matrix to Cartesian parameters in input frame + final RealMatrix cartesianCovarianceIn = changeCovarianceType(orbit, covOrbitType, covAngleType, + OrbitType.CARTESIAN, PositionAngle.MEAN, + inputCov); + + // Get the Cartesian covariance matrix converted to frameOut + final RealMatrix cartesianCovarianceOut = j.multiply(cartesianCovarianceIn.multiplyTransposed(j)); + + // Output frame is pseudo-inertial + if (frameOut.isPseudoInertial()) { + + // Convert output Cartesian matrix to initial orbit type and position angle + return changeCovarianceType(orbit, OrbitType.CARTESIAN, PositionAngle.MEAN, + covOrbitType, covAngleType, cartesianCovarianceOut); + } + // Output frame is not pseudo-inertial + else { + + // Output cartesian matrix + return cartesianCovarianceOut; + } + + } + + // Input frame is not pseudo-inertial + else { + // Covariance is expressed in cartesian elements + if (covOrbitType.equals(OrbitType.CARTESIAN)) { + + // Output the Cartesian covariance matrix converted to frameOut + return j.multiply(inputCov.multiplyTransposed(j)); - // Convert orbit frame to output frame - final Orbit outOrbit = new CartesianOrbit(orbit.getPVCoordinates(frameOut), frameOut, orbit.getMu()); + } + // Covariance is not expressed in cartesian elements + else { + throw new OrekitException(OrekitMessages.WRONG_ORBIT_PARAMETERS_TYPE, covOrbitType, + OrbitType.CARTESIAN); + } - // Convert output Cartesian matrix to initial orbit type and position angle - return changeCovarianceType(outOrbit, OrbitType.CARTESIAN, PositionAngle.MEAN, - covOrbitType, covAngleType, cartesianCovarianceOut); + } } - /** Get the covariance matrix in another orbit type. + /** + * Get the covariance matrix in another orbit type. + * * @param orbit orbit to which the covariance matrix should correspond * @param inOrbitType initial orbit type of the state covariance matrix * @param inAngleType initial position angle type of the state covariance matrix * @param outOrbitType target orbit type of the state covariance matrix * @param outAngleType target position angle type of the state covariance matrix * @param inputCov input covariance - * @return the covariance expressed in the target orbit type with the - * target position angle + * @return the covariance expressed in the target orbit type with the target position angle */ public static RealMatrix changeCovarianceType(final Orbit orbit, final OrbitType inOrbitType, final PositionAngle inAngleType, @@ -320,14 +383,14 @@ public class StateCovarianceMatrixProvider implements AdditionalStateProvider { // C: Cartesian parameters // Compute dOutputdCartesian - final double[][] aOC = new double[STATE_DIMENSION][STATE_DIMENSION]; - final Orbit orbitInOutputType = outOrbitType.convertType(orbit); + final double[][] aOC = new double[STATE_DIMENSION][STATE_DIMENSION]; + final Orbit orbitInOutputType = outOrbitType.convertType(orbit); orbitInOutputType.getJacobianWrtCartesian(outAngleType, aOC); final RealMatrix dOdC = new Array2DRowRealMatrix(aOC, false); // Compute dCartesiandInput - final double[][] aCI = new double[STATE_DIMENSION][STATE_DIMENSION]; - final Orbit orbitInInputType = inOrbitType.convertType(orbit); + final double[][] aCI = new double[STATE_DIMENSION][STATE_DIMENSION]; + final Orbit orbitInInputType = inOrbitType.convertType(orbit); orbitInInputType.getJacobianWrtParameters(inAngleType, aCI); final RealMatrix dCdI = new Array2DRowRealMatrix(aCI, false); @@ -342,28 +405,204 @@ public class StateCovarianceMatrixProvider implements AdditionalStateProvider { } - /** Convert an array to a matrix (6x6 dimension). + /** + * Convert the covariance matrix from a {@link Frame frame} to a {@link LOFType commonly used local orbital frame}. + *

+ * The transformation is based on Equation (20) of "Covariance Transformations for Satellite Flight Dynamics + * Operations" by David A. Vallado". + *

+ * As the frame transformation must be performed with the covariance expressed in Cartesian elements, both the orbit + * and position angle types of the input covariance must be provided. + *

+ * The output covariance matrix will necessarily be expressed in cartesian elements and not converted back to + * its original expression (if input different from cartesian elements). + *

+ * + * @param orbit orbit to which the covariance matrix should correspond + * @param frameIn the frame in which the input covariance matrix is expressed. In case the frame is not + * pseudo-inertial, the input covariance matrix is expected to be expressed in cartesian elements. + * @param lofOut the target local orbital frame + * @param inputCov input covariance + * @param covOrbitType orbit type of the covariance matrix (used if frameIn is pseudo-inertial) + * @param covAngleType position angle type of the covariance matrix (used if frameIn is pseudo-inertial) (not used + * if covOrbitType equals {@code CARTESIAN}) + * @return the covariance matrix expressed in the target local orbital frame in cartesian elements + * @throws OrekitException if input frame is not pseudo-inertial and the input covariance is + * not expressed in cartesian elements. + */ + public static RealMatrix changeCovarianceFrame(final Orbit orbit, + final Frame frameIn, final LOFType lofOut, + final RealMatrix inputCov, + final OrbitType covOrbitType, final PositionAngle covAngleType) { + + // Input frame is inertial + if (frameIn.isPseudoInertial()) { + + // Convert input matrix to Cartesian parameters in input frame + final RealMatrix cartesianCovarianceIn = changeCovarianceType(orbit, covOrbitType, covAngleType, + OrbitType.CARTESIAN, PositionAngle.MEAN, + inputCov); + + // Compute rotation matrix from frameIn to lofOut + final Rotation rotationFromFrameInToLofOut = lofOut.rotationFromInertial(orbit.getPVCoordinates(frameIn)); + + // Builds the matrix to perform covariance transformation + final RealMatrix transformationMatrix = buildTransformationMatrixFromRotation(rotationFromFrameInToLofOut); + + // Get the Cartesian covariance matrix converted to frameOut + final RealMatrix cartesianCovarianceOut = + transformationMatrix.multiply(cartesianCovarianceIn.multiplyTransposed(transformationMatrix)); + + // Return converted covariance matrix expressed in cartesian elements + return cartesianCovarianceOut; + + } + + // Input frame is not inertial so the covariance matrix is expected to be in cartesian elements + else { + if (covOrbitType.equals(OrbitType.CARTESIAN)) { + final Frame orbitInertialFrame = orbit.getFrame(); + + // Compute rotation matrix from frameIn to orbit inertial frame + final RealMatrix cartesianCovarianceInOrbitFrame = + changeCovarianceFrame(orbit, frameIn, orbitInertialFrame, inputCov, + OrbitType.CARTESIAN, + PositionAngle.MEAN); + + // Convert from orbit inertial frame to lofOut + return changeCovarianceFrame(orbit, orbitInertialFrame, lofOut, cartesianCovarianceInOrbitFrame, + OrbitType.CARTESIAN, + PositionAngle.MEAN); + + } + else { + throw new OrekitException(OrekitMessages.WRONG_ORBIT_PARAMETERS_TYPE, covOrbitType, + OrbitType.CARTESIAN); + + } + + } + + } + + /** + * Get the orbit type in which the covariance matrix is expressed. + * + * @return the orbit type + */ + public OrbitType getCovarianceOrbitType() { + return covOrbitType; + } + + /** + * {@inheritDoc} + *

+ * The covariance matrix can be computed only if the State Transition Matrix state is available. + *

+ */ + @Override + public boolean yield(final SpacecraftState state) { + return !state.hasAdditionalState(stmName); + } + + /** + * Get the state covariance matrix (6x6 dimension). + *

+ * The output covariance matrix is expressed in the same orbit type as {@link #getCovarianceOrbitType()}. + *

+ * It is possible to change the covariance frame by using the + * {@link #changeCovarianceFrame(Orbit, Frame, Frame, RealMatrix, OrbitType, PositionAngle)} method. + *

+ * It is also possible to change the covariance orbit type by using the + * {@link #changeCovarianceType(Orbit, OrbitType, PositionAngle, OrbitType, PositionAngle, RealMatrix)} method. + * + * @param state spacecraft state to which the covariance matrix should correspond + * @return the state covariance matrix + * @see #getStateCovariance(SpacecraftState, Frame) + * @see #getStateCovariance(SpacecraftState, OrbitType, PositionAngle) + */ + public RealMatrix getStateCovariance(final SpacecraftState state) { + + // Get the current propagated covariance + final RealMatrix covariance = toRealMatrix(state.getAdditionalState(additionalName)); + + // Return the converted covariance + return covariance; + + } + + /** + * Convert an array to a matrix (6x6 dimension). + * * @param array input array * @return the corresponding matrix */ private RealMatrix toRealMatrix(final double[] array) { final RealMatrix matrix = MatrixUtils.createRealMatrix(STATE_DIMENSION, STATE_DIMENSION); - int index = 0; + int index = 0; for (int i = 0; i < STATE_DIMENSION; ++i) { for (int j = 0; j < STATE_DIMENSION; ++j) { matrix.setEntry(i, j, array[index++]); } } return matrix; + + } + + /** + * Get the state covariance matrix (6x6 dimension) expressed in a given frame. + *

+ * The output covariance matrix is expressed in the same orbit type as {@link #getCovarianceOrbitType()}. + *

+ * It is also possible to change the covariance orbit type by using the + * {@link #changeCovarianceType(Orbit, OrbitType, PositionAngle, OrbitType, PositionAngle, RealMatrix)} method. + * + * @param state spacecraft state to which the covariance matrix should correspond + * @param frame output frame for which the output covariance matrix must be expressed (must be inertial) + * @return the state covariance matrix expressed in frame + * @see #getStateCovariance(SpacecraftState) + * @see #getStateCovariance(SpacecraftState, OrbitType, PositionAngle) + */ + public RealMatrix getStateCovariance(final SpacecraftState state, final Frame frame) { + + // Get the current propagated covariance + final RealMatrix covariance = toRealMatrix(state.getAdditionalState(additionalName)); + + // Return the converted covariance + return changeCovarianceFrame(state.getOrbit(), state.getFrame(), frame, covariance, covOrbitType, covAngleType); + + } + + /** + * Get the state covariance matrix (6x6 dimension) expressed in a given orbit type. + * + * @param state spacecraft state to which the covariance matrix should correspond + * @param orbitType output orbit type + * @param angleType output position angle (not used if orbitType equals {@code CARTESIAN}) + * @return the state covariance matrix in orbitType and angleType + * @see #getStateCovariance(SpacecraftState) + * @see #getStateCovariance(SpacecraftState, Frame) + */ + public RealMatrix getStateCovariance(final SpacecraftState state, final OrbitType orbitType, + final PositionAngle angleType) { + + // Get the current propagated covariance + final RealMatrix covariance = toRealMatrix(state.getAdditionalState(additionalName)); + + // Return the converted covariance + return changeCovarianceType(state.getOrbit(), covOrbitType, covAngleType, orbitType, angleType, covariance); + } - /** Set the covariance data into an array. + /** + * Set the covariance data into an array. + * * @param covariance covariance matrix * @return an array containing the covariance data */ private double[] toArray(final RealMatrix covariance) { final double[] array = new double[STATE_DIMENSION * STATE_DIMENSION]; - int index = 0; + int index = 0; for (int i = 0; i < STATE_DIMENSION; ++i) { for (int j = 0; j < STATE_DIMENSION; ++j) { array[index++] = covariance.getEntry(i, j); diff --git a/src/main/resources/assets/org/orekit/localization/OrekitMessages_da.utf8 b/src/main/resources/assets/org/orekit/localization/OrekitMessages_da.utf8 index 896542aa4f954c8930fbf5f1109c8fe1c3f0d4ab..73e273c2a1c972573324578fba8eebddb4a9055c 100644 --- a/src/main/resources/assets/org/orekit/localization/OrekitMessages_da.utf8 +++ b/src/main/resources/assets/org/orekit/localization/OrekitMessages_da.utf8 @@ -761,3 +761,6 @@ NO_UNSCENTED_TRANSFORM_CONFIGURED = # value is not strictly positive: {0} NOT_STRICTLY_POSITIVE = + +# orbital parameters type: {0} is different from expected orbital type : {1} +WRONG_ORBIT_PARAMETERS_TYPE = diff --git a/src/main/resources/assets/org/orekit/localization/OrekitMessages_de.utf8 b/src/main/resources/assets/org/orekit/localization/OrekitMessages_de.utf8 index b33ca77dec4c347a3e56c4db18e2627eea577617..7ef66367ad2232f46ef3d113209f66657c6f9d2e 100644 --- a/src/main/resources/assets/org/orekit/localization/OrekitMessages_de.utf8 +++ b/src/main/resources/assets/org/orekit/localization/OrekitMessages_de.utf8 @@ -761,3 +761,6 @@ NO_UNSCENTED_TRANSFORM_CONFIGURED = # value is not strictly positive: {0} NOT_STRICTLY_POSITIVE = + +# orbital parameters type: {0} is different from expected orbital type : {1} +WRONG_ORBIT_PARAMETERS_TYPE = diff --git a/src/main/resources/assets/org/orekit/localization/OrekitMessages_el.utf8 b/src/main/resources/assets/org/orekit/localization/OrekitMessages_el.utf8 index eb913f12effa1ea8d0989a78ef6429578cdab68f..2cc462a6d9896a28d6a6ad49888ec330b934abaa 100644 --- a/src/main/resources/assets/org/orekit/localization/OrekitMessages_el.utf8 +++ b/src/main/resources/assets/org/orekit/localization/OrekitMessages_el.utf8 @@ -761,3 +761,6 @@ NO_UNSCENTED_TRANSFORM_CONFIGURED = # value is not strictly positive: {0} NOT_STRICTLY_POSITIVE = + +# orbital parameters type: {0} is different from expected orbital type : {1} +WRONG_ORBIT_PARAMETERS_TYPE = diff --git a/src/main/resources/assets/org/orekit/localization/OrekitMessages_en.utf8 b/src/main/resources/assets/org/orekit/localization/OrekitMessages_en.utf8 index f535b31df09752312fc3ea4c0afdfc4d78c6e867..1f120fe32fb7267941ddaafdf7495fb3b3abeaed 100644 --- a/src/main/resources/assets/org/orekit/localization/OrekitMessages_en.utf8 +++ b/src/main/resources/assets/org/orekit/localization/OrekitMessages_en.utf8 @@ -761,3 +761,6 @@ NO_UNSCENTED_TRANSFORM_CONFIGURED = no unscented transform configured # value is not strictly positive: {0} NOT_STRICTLY_POSITIVE = value is not strictly positive: {0} + +# orbital parameters type: {0} is different from expected orbital type : {1} +WRONG_ORBIT_PARAMETERS_TYPE = orbital parameters type: {0} is different from expected orbital type : {1} diff --git a/src/main/resources/assets/org/orekit/localization/OrekitMessages_es.utf8 b/src/main/resources/assets/org/orekit/localization/OrekitMessages_es.utf8 index f90f57ffe3814773c4e87aa73cb6a1382f64c194..7546e0818e8a8e22006daa59285176869d541f49 100644 --- a/src/main/resources/assets/org/orekit/localization/OrekitMessages_es.utf8 +++ b/src/main/resources/assets/org/orekit/localization/OrekitMessages_es.utf8 @@ -761,3 +761,6 @@ NO_UNSCENTED_TRANSFORM_CONFIGURED = # value is not strictly positive: {0} NOT_STRICTLY_POSITIVE = + +# orbital parameters type: {0} is different from expected orbital type : {1} +WRONG_ORBIT_PARAMETERS_TYPE = diff --git a/src/main/resources/assets/org/orekit/localization/OrekitMessages_fr.utf8 b/src/main/resources/assets/org/orekit/localization/OrekitMessages_fr.utf8 index 34104a890222f760aa6bb41e8998aaa3ac25932f..36e36e8babfbd7f76d676877581ecd17b9273ce1 100644 --- a/src/main/resources/assets/org/orekit/localization/OrekitMessages_fr.utf8 +++ b/src/main/resources/assets/org/orekit/localization/OrekitMessages_fr.utf8 @@ -761,3 +761,6 @@ NO_UNSCENTED_TRANSFORM_CONFIGURED = aucune transformation unscented configurée # value is not strictly positive: {0} NOT_STRICTLY_POSITIVE = la valeur {0} n''est pas strictement positive + +# orbital parameters type: {0} is different from expected orbital type : {1} +WRONG_ORBIT_PARAMETERS_TYPE = Le type de paramètres orbitaux: {0} est différent du type de paramètres orbitaux attendus: {1} diff --git a/src/main/resources/assets/org/orekit/localization/OrekitMessages_gl.utf8 b/src/main/resources/assets/org/orekit/localization/OrekitMessages_gl.utf8 index ad989ef00b867ba1d88ad64a413372a677168fb1..0d616bfad85905e28e475c8183e01f35dbf3dc9c 100644 --- a/src/main/resources/assets/org/orekit/localization/OrekitMessages_gl.utf8 +++ b/src/main/resources/assets/org/orekit/localization/OrekitMessages_gl.utf8 @@ -761,3 +761,6 @@ NO_UNSCENTED_TRANSFORM_CONFIGURED = # value is not strictly positive: {0} NOT_STRICTLY_POSITIVE = + +# orbital parameters type: {0} is different from expected orbital type : {1} +WRONG_ORBIT_PARAMETERS_TYPE = diff --git a/src/main/resources/assets/org/orekit/localization/OrekitMessages_it.utf8 b/src/main/resources/assets/org/orekit/localization/OrekitMessages_it.utf8 index 0ea505e69a670f10a0505a54b1ada70b80a2a75d..29c0d3f29b8c5ab1285f8576e70dee92025e5ff2 100644 --- a/src/main/resources/assets/org/orekit/localization/OrekitMessages_it.utf8 +++ b/src/main/resources/assets/org/orekit/localization/OrekitMessages_it.utf8 @@ -761,3 +761,6 @@ NO_UNSCENTED_TRANSFORM_CONFIGURED = # value is not strictly positive: {0} NOT_STRICTLY_POSITIVE = + +# orbital parameters type: {0} is different from expected orbital type : {1} +WRONG_ORBIT_PARAMETERS_TYPE = diff --git a/src/main/resources/assets/org/orekit/localization/OrekitMessages_no.utf8 b/src/main/resources/assets/org/orekit/localization/OrekitMessages_no.utf8 index dc717e547f14c128e2eb991148dca87f4c03fd39..ad25f3850bface0ed6938efa62cbe36a3a4ad896 100644 --- a/src/main/resources/assets/org/orekit/localization/OrekitMessages_no.utf8 +++ b/src/main/resources/assets/org/orekit/localization/OrekitMessages_no.utf8 @@ -761,3 +761,6 @@ NO_UNSCENTED_TRANSFORM_CONFIGURED = # value is not strictly positive: {0} NOT_STRICTLY_POSITIVE = + +# orbital parameters type: {0} is different from expected orbital type : {1} +WRONG_ORBIT_PARAMETERS_TYPE = diff --git a/src/main/resources/assets/org/orekit/localization/OrekitMessages_ro.utf8 b/src/main/resources/assets/org/orekit/localization/OrekitMessages_ro.utf8 index 148eeb975fdbb3189aa911acab3a57c4cb5532e9..5dfe61e05fb5308eb3e2079a859dcb0a2a8d4fe7 100644 --- a/src/main/resources/assets/org/orekit/localization/OrekitMessages_ro.utf8 +++ b/src/main/resources/assets/org/orekit/localization/OrekitMessages_ro.utf8 @@ -761,3 +761,6 @@ NO_UNSCENTED_TRANSFORM_CONFIGURED = # value is not strictly positive: {0} NOT_STRICTLY_POSITIVE = + +# orbital parameters type: {0} is different from expected orbital type : {1} +WRONG_ORBIT_PARAMETERS_TYPE = diff --git a/src/test/java/org/orekit/errors/OrekitMessagesTest.java b/src/test/java/org/orekit/errors/OrekitMessagesTest.java index 240958f9eb2bda5966f1417f5cb391827deb2501..1bae7b79f6aa596f204eb5067bda7d321d059389 100644 --- a/src/test/java/org/orekit/errors/OrekitMessagesTest.java +++ b/src/test/java/org/orekit/errors/OrekitMessagesTest.java @@ -30,7 +30,7 @@ public class OrekitMessagesTest { @Test public void testMessageNumber() { - Assertions.assertEquals(254, OrekitMessages.values().length); + Assertions.assertEquals(255, OrekitMessages.values().length); } @Test diff --git a/src/test/java/org/orekit/propagation/StateCovarianceMatrixProviderTest.java b/src/test/java/org/orekit/propagation/StateCovarianceMatrixProviderTest.java index 73a37b091be86c9d84dcc903957b1449fa075b38..b0720bf943e681d9daf970f9e641dc2de161fd3d 100644 --- a/src/test/java/org/orekit/propagation/StateCovarianceMatrixProviderTest.java +++ b/src/test/java/org/orekit/propagation/StateCovarianceMatrixProviderTest.java @@ -17,15 +17,17 @@ package org.orekit.propagation; import org.hipparchus.geometry.euclidean.threed.Vector3D; +import org.hipparchus.linear.BlockRealMatrix; import org.hipparchus.linear.MatrixUtils; import org.hipparchus.linear.RealMatrix; import org.hipparchus.ode.ODEIntegrator; import org.hipparchus.ode.nonstiff.ClassicalRungeKuttaIntegrator; import org.hipparchus.util.FastMath; import org.junit.jupiter.api.Assertions; -import org.junit.jupiter.api.BeforeEach; +import org.junit.jupiter.api.DisplayName; import org.junit.jupiter.api.Test; import org.orekit.Utils; +import org.orekit.errors.OrekitException; import org.orekit.forces.ForceModel; import org.orekit.forces.gravity.HolmesFeatherstoneAttractionModel; import org.orekit.forces.gravity.potential.GravityFieldFactory; @@ -33,6 +35,7 @@ import org.orekit.forces.gravity.potential.ICGEMFormatReader; import org.orekit.forces.gravity.potential.NormalizedSphericalHarmonicsProvider; import org.orekit.frames.Frame; import org.orekit.frames.FramesFactory; +import org.orekit.frames.LOFType; import org.orekit.orbits.CartesianOrbit; import org.orekit.orbits.Orbit; import org.orekit.orbits.OrbitType; @@ -49,26 +52,7 @@ public class StateCovarianceMatrixProviderTest { private SpacecraftState initialState; private double[][] initCov; - - @BeforeEach - public void setUp() { - Utils.setDataRoot("orbit-determination/february-2016:potential/icgem-format"); - GravityFieldFactory.addPotentialCoefficientsReader(new ICGEMFormatReader("eigen-6s-truncated", true)); - Orbit initialOrbit = new CartesianOrbit(new PVCoordinates(new Vector3D(7526993.581890527, -9646310.10026971, 1464110.4928112086), - new Vector3D(3033.79456099698, 1715.265069098717, -4447.658745923895)), - FramesFactory.getEME2000(), - new AbsoluteDate("2016-02-13T16:00:00.000", TimeScalesFactory.getUTC()), - Constants.WGS84_EARTH_MU); - initialState = new SpacecraftState(initialOrbit); - initCov = new double[][] { - {8.651816029e+01, 5.689987127e+01, -2.763870764e+01, -2.435617201e-02, 2.058274137e-02, -5.872883051e-03}, - {5.689987127e+01, 7.070624321e+01, 1.367120909e+01, -6.112622013e-03, 7.623626008e-03, -1.239413190e-02}, - {-2.763870764e+01, 1.367120909e+01, 1.811858898e+02, 3.143798992e-02, -4.963106559e-02, -7.420114385e-04}, - {-2.435617201e-02, -6.112622013e-03, 3.143798992e-02, 4.657077389e-05, 1.469943634e-05, 3.328475593e-05}, - {2.058274137e-02, 7.623626008e-03, -4.963106559e-02, 1.469943634e-05, 3.950715934e-05, 2.516044258e-05}, - {-5.872883051e-03, -1.239413190e-02, -7.420114385e-04, 3.328475593e-05, 2.516044258e-05, 3.547466120e-05} - }; - } + final private double DEFAULT_VALLADO_THRESHOLD = 1e-6; /** * Unit test for the covariance frame transformation. @@ -76,6 +60,9 @@ public class StateCovarianceMatrixProviderTest { @Test public void testFrameConversion() { + // Initialization + setUp(); + // Reference final RealMatrix referenceCov = MatrixUtils.createRealMatrix(initCov); @@ -84,22 +71,128 @@ public class StateCovarianceMatrixProviderTest { final Frame frameB = FramesFactory.getTEME(); // First transformation - RealMatrix transformedCov = StateCovarianceMatrixProvider.changeCovarianceFrame(initialState.getOrbit(), frameA, frameB, referenceCov, OrbitType.CARTESIAN, PositionAngle.MEAN); + RealMatrix transformedCov = + StateCovarianceMatrixProvider.changeCovarianceFrame(initialState.getOrbit(), frameA, frameB, + referenceCov, OrbitType.CARTESIAN, + PositionAngle.MEAN); // Second transformation - transformedCov = StateCovarianceMatrixProvider.changeCovarianceFrame(initialState.getOrbit(), frameB, frameA, transformedCov, OrbitType.CARTESIAN, PositionAngle.MEAN); + transformedCov = StateCovarianceMatrixProvider.changeCovarianceFrame(initialState.getOrbit(), frameB, frameA, + transformedCov, OrbitType.CARTESIAN, + PositionAngle.MEAN); // Verify compareCovariance(referenceCov, transformedCov, 5.2e-15); } + public void setUp() { + Utils.setDataRoot("orbit-determination/february-2016:potential/icgem-format"); + GravityFieldFactory.addPotentialCoefficientsReader(new ICGEMFormatReader("eigen-6s-truncated", true)); + Orbit initialOrbit = new CartesianOrbit( + new PVCoordinates(new Vector3D(7526993.581890527, -9646310.10026971, 1464110.4928112086), + new Vector3D(3033.79456099698, 1715.265069098717, -4447.658745923895)), + FramesFactory.getEME2000(), + new AbsoluteDate("2016-02-13T16:00:00.000", TimeScalesFactory.getUTC()), + Constants.WGS84_EARTH_MU); + initialState = new SpacecraftState(initialOrbit); + initCov = new double[][] { + { 8.651816029e+01, 5.689987127e+01, -2.763870764e+01, -2.435617201e-02, 2.058274137e-02, + -5.872883051e-03 }, + { 5.689987127e+01, 7.070624321e+01, 1.367120909e+01, -6.112622013e-03, 7.623626008e-03, + -1.239413190e-02 }, + { -2.763870764e+01, 1.367120909e+01, 1.811858898e+02, 3.143798992e-02, -4.963106559e-02, + -7.420114385e-04 }, + { -2.435617201e-02, -6.112622013e-03, 3.143798992e-02, 4.657077389e-05, 1.469943634e-05, + 3.328475593e-05 }, + { 2.058274137e-02, 7.623626008e-03, -4.963106559e-02, 1.469943634e-05, 3.950715934e-05, + 2.516044258e-05 }, + { -5.872883051e-03, -1.239413190e-02, -7.420114385e-04, 3.328475593e-05, 2.516044258e-05, + 3.547466120e-05 } + }; + } + + /** + * Compare two covariance matrices + * + * @param reference reference covariance + * @param computed computed covariance + * @param threshold threshold for comparison + */ + private void compareCovariance(final RealMatrix reference, final RealMatrix computed, final double threshold) { + for (int row = 0; row < reference.getRowDimension(); row++) { + for (int column = 0; column < reference.getColumnDimension(); column++) { + if (reference.getEntry(row, column) == 0) { + Assertions.assertEquals(reference.getEntry(row, column), computed.getEntry(row, column), + threshold); + } + else { + Assertions.assertEquals(reference.getEntry(row, column), computed.getEntry(row, column), + FastMath.abs(threshold * reference.getEntry(row, column))); + } + } + } + } + + /** + * In that particular case, the RTN local orbital frame is perfectly aligned with the inertial frame. Hence, the + * frame conversion should return the same covariance matrix as input. + */ + @Test + @DisplayName("Test conversion from inertial frame to RTN local orbital frame") + void should_return_same_covariance_matrix() { + + // Given + final AbsoluteDate initialDate = new AbsoluteDate(); + final Frame initialInertialFrame = FramesFactory.getGCRF(); + final double mu = 398600e9; + + final PVCoordinates initialPV = new PVCoordinates( + new Vector3D(6778000, 0, 0), + new Vector3D(0, 7668.63, 0)); + + final Orbit initialOrbit = new CartesianOrbit(initialPV, initialInertialFrame, initialDate, mu); + + final RealMatrix initialCovarianceInInertialFrame = new BlockRealMatrix(new double[][] { + { 1, 0, 0, 0, 1e-5, 1e-4 }, + { 0, 1, 0, 0, 0, 1e-5 }, + { 0, 0, 1, 0, 0, 0 }, + { 0, 0, 0, 1e-3, 0, 0 }, + { 1e-5, 0, 0, 0, 1e-3, 0 }, + { 1e-4, 1e-5, 0, 0, 0, 1e-3 } + }); + + // When + final RealMatrix covarianceMatrixInRTN = StateCovarianceMatrixProvider.changeCovarianceFrame(initialOrbit, + initialInertialFrame, + LOFType.QSW, + initialCovarianceInInertialFrame, + OrbitType.CARTESIAN, + PositionAngle.MEAN); + + // Then + final RealMatrix expectedMatrixInRTN = new BlockRealMatrix(new double[][] { + { 1, 0, 0, 0, 1e-5, 1e-4 }, + { 0, 1, 0, 0, 0, 1e-5 }, + { 0, 0, 1, 0, 0, 0 }, + { 0, 0, 0, 1e-3, 0, 0 }, + { 1e-5, 0, 0, 0, 1e-3, 0 }, + { 1e-4, 1e-5, 0, 0, 0, 1e-3 } + }); + + compareCovariance(covarianceMatrixInRTN, expectedMatrixInRTN, 1e-20); + + } + /** * Unit test for the covariance type transformation. */ @Test public void testTypeConversion() { + // Initialization + setUp(); + // Reference final RealMatrix referenceCov = MatrixUtils.createRealMatrix(initCov); @@ -108,12 +201,14 @@ public class StateCovarianceMatrixProviderTest { final OrbitType kep = OrbitType.KEPLERIAN; // First transformation - RealMatrix transformedCov = StateCovarianceMatrixProvider.changeCovarianceType(initialState.getOrbit(), cart, PositionAngle.MEAN, - kep, PositionAngle.MEAN, referenceCov); + RealMatrix transformedCov = + StateCovarianceMatrixProvider.changeCovarianceType(initialState.getOrbit(), cart, PositionAngle.MEAN, + kep, PositionAngle.MEAN, referenceCov); // Second transformation - transformedCov = StateCovarianceMatrixProvider.changeCovarianceType(initialState.getOrbit(), kep, PositionAngle.MEAN, - cart, PositionAngle.MEAN, transformedCov); + transformedCov = + StateCovarianceMatrixProvider.changeCovarianceType(initialState.getOrbit(), kep, PositionAngle.MEAN, + cart, PositionAngle.MEAN, transformedCov); // Verify compareCovariance(referenceCov, transformedCov, 3.5e-12); @@ -126,18 +221,22 @@ public class StateCovarianceMatrixProviderTest { @Test public void testWithNumericalPropagatorCartesian() { + // Initialization + setUp(); + // Integrator final double step = 60.0; final ODEIntegrator integrator = new ClassicalRungeKuttaIntegrator(step); // Numerical propagator - final String stmName = "STM"; - final OrbitType propType = OrbitType.CARTESIAN; - final PositionAngle angleType = PositionAngle.MEAN; + final String stmName = "STM"; + final OrbitType propType = OrbitType.CARTESIAN; + final PositionAngle angleType = PositionAngle.MEAN; final NumericalPropagator propagator = new NumericalPropagator(integrator); // Add a force model final NormalizedSphericalHarmonicsProvider gravity = GravityFieldFactory.getNormalizedProvider(2, 0); - final ForceModel holmesFeatherstone = new HolmesFeatherstoneAttractionModel(FramesFactory.getITRF(IERSConventions.IERS_2010, true), gravity); + final ForceModel holmesFeatherstone = + new HolmesFeatherstoneAttractionModel(FramesFactory.getITRF(IERSConventions.IERS_2010, true), gravity); propagator.addForceModel(holmesFeatherstone); // Finalize setting final MatricesHarvester harvester = propagator.setupMatricesComputation(stmName, null, null); @@ -145,12 +244,13 @@ public class StateCovarianceMatrixProviderTest { propagator.setPositionAngleType(angleType); // Create additional state - final String additionalName = "cartCov"; - final RealMatrix initialCov = MatrixUtils.createRealMatrix(initCov); - final StateCovarianceMatrixProvider provider = new StateCovarianceMatrixProvider(additionalName, stmName, harvester, - propagator.getOrbitType(), propagator.getPositionAngleType(), - initialCov, - OrbitType.CARTESIAN, PositionAngle.MEAN); + final String additionalName = "cartCov"; + final RealMatrix initialCov = MatrixUtils.createRealMatrix(initCov); + final StateCovarianceMatrixProvider provider = + new StateCovarianceMatrixProvider(additionalName, stmName, harvester, + propagator.getOrbitType(), propagator.getPositionAngleType(), + initialCov, + OrbitType.CARTESIAN, PositionAngle.MEAN); propagator.setInitialState(initialState); propagator.addAdditionalStateProvider(provider); @@ -162,12 +262,18 @@ public class StateCovarianceMatrixProviderTest { // Reference (computed using a different solution) final double[][] ref = new double[][] { - { 5.770543135e+02, 2.316979550e+02, -5.172369105e+02, -2.585893247e-01, 2.113809017e-01, -1.759509343e-01}, - { 2.316979550e+02, 1.182942930e+02, -1.788422178e+02, -9.570305681e-02, 7.792155309e-02, -7.435822327e-02}, - {-5.172369105e+02, -1.788422178e+02, 6.996248500e+02, 2.633605389e-01, -2.480144888e-01, 1.908427233e-01}, - {-2.585893247e-01, -9.570305681e-02, 2.633605389e-01, 1.419148897e-04, -8.715858320e-05, 1.024944399e-04}, - { 2.113809017e-01, 7.792155309e-02, -2.480144888e-01, -8.715858320e-05, 1.069566588e-04, -5.667563856e-05}, - {-1.759509343e-01, -7.435822327e-02, 1.908427233e-01, 1.024944399e-04, -5.667563856e-05, 8.178356868e-05} + { 5.770543135e+02, 2.316979550e+02, -5.172369105e+02, -2.585893247e-01, 2.113809017e-01, + -1.759509343e-01 }, + { 2.316979550e+02, 1.182942930e+02, -1.788422178e+02, -9.570305681e-02, 7.792155309e-02, + -7.435822327e-02 }, + { -5.172369105e+02, -1.788422178e+02, 6.996248500e+02, 2.633605389e-01, -2.480144888e-01, + 1.908427233e-01 }, + { -2.585893247e-01, -9.570305681e-02, 2.633605389e-01, 1.419148897e-04, -8.715858320e-05, + 1.024944399e-04 }, + { 2.113809017e-01, 7.792155309e-02, -2.480144888e-01, -8.715858320e-05, 1.069566588e-04, + -5.667563856e-05 }, + { -1.759509343e-01, -7.435822327e-02, 1.908427233e-01, 1.024944399e-04, -5.667563856e-05, + 8.178356868e-05 } }; final RealMatrix referenceCov = MatrixUtils.createRealMatrix(ref); @@ -186,7 +292,10 @@ public class StateCovarianceMatrixProviderTest { RealMatrix transformedCovA = provider.getStateCovariance(propagated, frameB); // Second transformation - RealMatrix transformedCovB = StateCovarianceMatrixProvider.changeCovarianceFrame(propagated.getOrbit(), propagated.getFrame(), frameB, propagatedCov, OrbitType.CARTESIAN, PositionAngle.MEAN); + RealMatrix transformedCovB = + StateCovarianceMatrixProvider.changeCovarianceFrame(propagated.getOrbit(), propagated.getFrame(), + frameB, propagatedCov, OrbitType.CARTESIAN, + PositionAngle.MEAN); // Verify compareCovariance(transformedCovA, transformedCovB, 1.0e-15); @@ -196,14 +305,17 @@ public class StateCovarianceMatrixProviderTest { /////////// // Define a new output frame - final OrbitType outOrbitType = OrbitType.KEPLERIAN; + final OrbitType outOrbitType = OrbitType.KEPLERIAN; final PositionAngle outAngleType = PositionAngle.MEAN; // Transformation using getStateJacobian() method RealMatrix transformedCovC = provider.getStateCovariance(propagated, outOrbitType, outAngleType); // Second transformation - RealMatrix transformedCovD = StateCovarianceMatrixProvider.changeCovarianceType(propagated.getOrbit(), OrbitType.CARTESIAN, PositionAngle.MEAN, outOrbitType, outAngleType, propagatedCov); + RealMatrix transformedCovD = + StateCovarianceMatrixProvider.changeCovarianceType(propagated.getOrbit(), OrbitType.CARTESIAN, + PositionAngle.MEAN, outOrbitType, outAngleType, + propagatedCov); // Verify compareCovariance(transformedCovC, transformedCovD, 1.0e-15); @@ -211,33 +323,38 @@ public class StateCovarianceMatrixProviderTest { } /** - * Unit test for covariance propagation in Cartesian elements. - * The difference here is that the propagator uses its default orbit type: EQUINOCTIAL + * Unit test for covariance propagation in Cartesian elements. The difference here is that the propagator uses its + * default orbit type: EQUINOCTIAL */ @Test public void testWithNumericalPropagatorDefault() { + // Initialization + setUp(); + // Integrator final double step = 60.0; final ODEIntegrator integrator = new ClassicalRungeKuttaIntegrator(step); // Numerical propagator - final String stmName = "STM"; + final String stmName = "STM"; final NumericalPropagator propagator = new NumericalPropagator(integrator); // Add a force model final NormalizedSphericalHarmonicsProvider gravity = GravityFieldFactory.getNormalizedProvider(2, 0); - final ForceModel holmesFeatherstone = new HolmesFeatherstoneAttractionModel(FramesFactory.getITRF(IERSConventions.IERS_2010, true), gravity); + final ForceModel holmesFeatherstone = + new HolmesFeatherstoneAttractionModel(FramesFactory.getITRF(IERSConventions.IERS_2010, true), gravity); propagator.addForceModel(holmesFeatherstone); // Finalize setting final MatricesHarvester harvester = propagator.setupMatricesComputation(stmName, null, null); // Create additional state - final String additionalName = "cartCov"; - final RealMatrix initialCov = MatrixUtils.createRealMatrix(initCov); - final StateCovarianceMatrixProvider provider = new StateCovarianceMatrixProvider(additionalName, stmName, harvester, - propagator.getOrbitType(), propagator.getPositionAngleType(), - initialCov, - OrbitType.CARTESIAN, PositionAngle.MEAN); + final String additionalName = "cartCov"; + final RealMatrix initialCov = MatrixUtils.createRealMatrix(initCov); + final StateCovarianceMatrixProvider provider = + new StateCovarianceMatrixProvider(additionalName, stmName, harvester, + propagator.getOrbitType(), propagator.getPositionAngleType(), + initialCov, + OrbitType.CARTESIAN, PositionAngle.MEAN); propagator.setInitialState(initialState); propagator.addAdditionalStateProvider(provider); @@ -249,12 +366,18 @@ public class StateCovarianceMatrixProviderTest { // Reference (computed using a different solution) final double[][] ref = new double[][] { - { 5.770543135e+02, 2.316979550e+02, -5.172369105e+02, -2.585893247e-01, 2.113809017e-01, -1.759509343e-01}, - { 2.316979550e+02, 1.182942930e+02, -1.788422178e+02, -9.570305681e-02, 7.792155309e-02, -7.435822327e-02}, - {-5.172369105e+02, -1.788422178e+02, 6.996248500e+02, 2.633605389e-01, -2.480144888e-01, 1.908427233e-01}, - {-2.585893247e-01, -9.570305681e-02, 2.633605389e-01, 1.419148897e-04, -8.715858320e-05, 1.024944399e-04}, - { 2.113809017e-01, 7.792155309e-02, -2.480144888e-01, -8.715858320e-05, 1.069566588e-04, -5.667563856e-05}, - {-1.759509343e-01, -7.435822327e-02, 1.908427233e-01, 1.024944399e-04, -5.667563856e-05, 8.178356868e-05} + { 5.770543135e+02, 2.316979550e+02, -5.172369105e+02, -2.585893247e-01, 2.113809017e-01, + -1.759509343e-01 }, + { 2.316979550e+02, 1.182942930e+02, -1.788422178e+02, -9.570305681e-02, 7.792155309e-02, + -7.435822327e-02 }, + { -5.172369105e+02, -1.788422178e+02, 6.996248500e+02, 2.633605389e-01, -2.480144888e-01, + 1.908427233e-01 }, + { -2.585893247e-01, -9.570305681e-02, 2.633605389e-01, 1.419148897e-04, -8.715858320e-05, + 1.024944399e-04 }, + { 2.113809017e-01, 7.792155309e-02, -2.480144888e-01, -8.715858320e-05, 1.069566588e-04, + -5.667563856e-05 }, + { -1.759509343e-01, -7.435822327e-02, 1.908427233e-01, 1.024944399e-04, -5.667563856e-05, + 8.178356868e-05 } }; final RealMatrix referenceCov = MatrixUtils.createRealMatrix(ref); @@ -273,20 +396,26 @@ public class StateCovarianceMatrixProviderTest { RealMatrix transformedCovA = provider.getStateCovariance(propagated, frameB); // Second transformation - RealMatrix transformedCovB = StateCovarianceMatrixProvider.changeCovarianceFrame(propagated.getOrbit(), propagated.getFrame(), frameB, propagatedCov, OrbitType.CARTESIAN, PositionAngle.MEAN); + RealMatrix transformedCovB = + StateCovarianceMatrixProvider.changeCovarianceFrame(propagated.getOrbit(), propagated.getFrame(), + frameB, propagatedCov, OrbitType.CARTESIAN, + PositionAngle.MEAN); // Verify compareCovariance(transformedCovA, transformedCovB, 1.0e-15); // Define a new output frame - final OrbitType outOrbitType = OrbitType.KEPLERIAN; + final OrbitType outOrbitType = OrbitType.KEPLERIAN; final PositionAngle outAngleType = PositionAngle.MEAN; // Transformation using getStateJacobian() method RealMatrix transformedCovC = provider.getStateCovariance(propagated, outOrbitType, outAngleType); // Second transformation - RealMatrix transformedCovD = StateCovarianceMatrixProvider.changeCovarianceType(propagated.getOrbit(), OrbitType.CARTESIAN, PositionAngle.MEAN, outOrbitType, outAngleType, propagatedCov); + RealMatrix transformedCovD = + StateCovarianceMatrixProvider.changeCovarianceType(propagated.getOrbit(), OrbitType.CARTESIAN, + PositionAngle.MEAN, outOrbitType, outAngleType, + propagatedCov); // Verify compareCovariance(transformedCovC, transformedCovD, 1.0e-15); @@ -299,22 +428,28 @@ public class StateCovarianceMatrixProviderTest { @Test public void testWithAnalyticalPropagator() { + // Initialization + setUp(); + // Numerical propagator - final String stmName = "STM"; - final OrbitType propType = OrbitType.CARTESIAN; + final String stmName = "STM"; + final OrbitType propType = OrbitType.CARTESIAN; final PositionAngle angleType = PositionAngle.MEAN; final EcksteinHechlerPropagator propagator = new EcksteinHechlerPropagator(initialState.getOrbit(), - GravityFieldFactory.getUnnormalizedProvider(6, 0)); + GravityFieldFactory.getUnnormalizedProvider( + 6, 0)); + // Finalize setting final MatricesHarvester harvester = propagator.setupMatricesComputation(stmName, null, null); // Create additional state - final String additionalName = "cartCov"; - final RealMatrix initialCov = MatrixUtils.createRealMatrix(initCov); - final StateCovarianceMatrixProvider provider = new StateCovarianceMatrixProvider(additionalName, stmName, harvester, - propType, angleType, - initialCov, - OrbitType.CARTESIAN, PositionAngle.MEAN); + final String additionalName = "cartCov"; + final RealMatrix initialCov = MatrixUtils.createRealMatrix(initCov); + final StateCovarianceMatrixProvider provider = + new StateCovarianceMatrixProvider(additionalName, stmName, harvester, + propType, angleType, + initialCov, + OrbitType.CARTESIAN, PositionAngle.MEAN); propagator.addAdditionalStateProvider(provider); // Propagate @@ -325,12 +460,18 @@ public class StateCovarianceMatrixProviderTest { // Reference (computed using a numerical solution) final double[][] ref = new double[][] { - { 5.770543135e+02, 2.316979550e+02, -5.172369105e+02, -2.585893247e-01, 2.113809017e-01, -1.759509343e-01}, - { 2.316979550e+02, 1.182942930e+02, -1.788422178e+02, -9.570305681e-02, 7.792155309e-02, -7.435822327e-02}, - {-5.172369105e+02, -1.788422178e+02, 6.996248500e+02, 2.633605389e-01, -2.480144888e-01, 1.908427233e-01}, - {-2.585893247e-01, -9.570305681e-02, 2.633605389e-01, 1.419148897e-04, -8.715858320e-05, 1.024944399e-04}, - { 2.113809017e-01, 7.792155309e-02, -2.480144888e-01, -8.715858320e-05, 1.069566588e-04, -5.667563856e-05}, - {-1.759509343e-01, -7.435822327e-02, 1.908427233e-01, 1.024944399e-04, -5.667563856e-05, 8.178356868e-05} + { 5.770543135e+02, 2.316979550e+02, -5.172369105e+02, -2.585893247e-01, 2.113809017e-01, + -1.759509343e-01 }, + { 2.316979550e+02, 1.182942930e+02, -1.788422178e+02, -9.570305681e-02, 7.792155309e-02, + -7.435822327e-02 }, + { -5.172369105e+02, -1.788422178e+02, 6.996248500e+02, 2.633605389e-01, -2.480144888e-01, + 1.908427233e-01 }, + { -2.585893247e-01, -9.570305681e-02, 2.633605389e-01, 1.419148897e-04, -8.715858320e-05, + 1.024944399e-04 }, + { 2.113809017e-01, 7.792155309e-02, -2.480144888e-01, -8.715858320e-05, 1.069566588e-04, + -5.667563856e-05 }, + { -1.759509343e-01, -7.435822327e-02, 1.908427233e-01, 1.024944399e-04, -5.667563856e-05, + 8.178356868e-05 } }; final RealMatrix referenceCov = MatrixUtils.createRealMatrix(ref); @@ -349,38 +490,589 @@ public class StateCovarianceMatrixProviderTest { RealMatrix transformedCovA = provider.getStateCovariance(propagated, frameB); // Second transformation - RealMatrix transformedCovB = StateCovarianceMatrixProvider.changeCovarianceFrame(propagated.getOrbit(), propagated.getFrame(), frameB, propagatedCov, OrbitType.CARTESIAN, PositionAngle.MEAN); + RealMatrix transformedCovB = + StateCovarianceMatrixProvider.changeCovarianceFrame(propagated.getOrbit(), propagated.getFrame(), + frameB, propagatedCov, OrbitType.CARTESIAN, + PositionAngle.MEAN); // Verify compareCovariance(transformedCovA, transformedCovB, 1.0e-15); // Define a new output frame - final OrbitType outOrbitType = OrbitType.KEPLERIAN; + final OrbitType outOrbitType = OrbitType.KEPLERIAN; final PositionAngle outAngleType = PositionAngle.MEAN; // Transformation using getStateJacobian() method RealMatrix transformedCovC = provider.getStateCovariance(propagated, outOrbitType, outAngleType); // Second transformation - RealMatrix transformedCovD = StateCovarianceMatrixProvider.changeCovarianceType(propagated.getOrbit(), OrbitType.CARTESIAN, PositionAngle.MEAN, outOrbitType, outAngleType, propagatedCov); + RealMatrix transformedCovD = + StateCovarianceMatrixProvider.changeCovarianceType(propagated.getOrbit(), OrbitType.CARTESIAN, + PositionAngle.MEAN, outOrbitType, outAngleType, + propagatedCov); // Verify compareCovariance(transformedCovC, transformedCovD, 1.0e-15); } + @Test + @DisplayName("Test covariance conversion from inertial frame to RTN local orbital frame") + void should_rotate_covariance_matrix_by_ninety_degrees() { + + // Given + final AbsoluteDate initialDate = new AbsoluteDate(); + final Frame initialInertialFrame = FramesFactory.getGCRF(); + final double mu = 398600e9; + + final PVCoordinates initialPV = new PVCoordinates( + new Vector3D(0, 6778000, 0), + new Vector3D(-7668.63, 0, 0)); + + final Orbit initialOrbit = new CartesianOrbit(initialPV, initialInertialFrame, initialDate, mu); + + final RealMatrix initialCovarianceInInertialFrame = new BlockRealMatrix(new double[][] { + { 1, 0, 0, 0, 0, 1e-5 }, + { 0, 1, 0, 0, 0, 0 }, + { 0, 0, 1, 0, 0, 0 }, + { 0, 0, 0, 1e-3, 0, 0 }, + { 0, 0, 0, 0, 1e-3, 0 }, + { 1e-5, 0, 0, 0, 0, 1e-3 } + }); + + // When + final RealMatrix convertedCovarianceMatrixInRTN = + StateCovarianceMatrixProvider.changeCovarianceFrame(initialOrbit, + initialInertialFrame, LOFType.QSW, + initialCovarianceInInertialFrame, + OrbitType.CARTESIAN, + PositionAngle.MEAN); + + // Then + // Expected covariance matrix obtained by rotation initial covariance matrix by 90 degrees + final RealMatrix expectedMatrixInRTN = new BlockRealMatrix(new double[][] { + { 1, 0, 0, 0, 0, 0 }, + { 0, 1, 0, 0, 0, -1e-5 }, + { 0, 0, 1, 0, 0, 0 }, + { 0, 0, 0, 1e-3, 0, 0 }, + { 0, 0, 0, 0, 1e-3, 0 }, + { 0, -1e-5, 0, 0, 0, 1e-3 } + }); + + compareCovariance(expectedMatrixInRTN, convertedCovarianceMatrixInRTN, 1e-20); + } + + @Test + @DisplayName("Test covariance conversion from RTN local orbital frame to inertial frame") + void should_rotate_covariance_matrix_by_minus_ninety_degrees() { + + // Given + final AbsoluteDate initialDate = new AbsoluteDate(); + final Frame initialInertialFrame = FramesFactory.getGCRF(); + final double mu = 398600e9; + + final PVCoordinates initialPV = new PVCoordinates( + new Vector3D(0, 6778000, 0), + new Vector3D(-7668.63, 0, 0)); + + final Orbit initialOrbit = new CartesianOrbit(initialPV, initialInertialFrame, initialDate, mu); + + final RealMatrix initialCovarianceInRTN = new BlockRealMatrix(new double[][] { + { 1, 0, 0, 0, 0, 0 }, + { 0, 1, 0, 0, 0, -1e-5 }, + { 0, 0, 1, 0, 0, 0 }, + { 0, 0, 0, 1e-3, 0, 0 }, + { 0, 0, 0, 0, 1e-3, 0 }, + { 0, -1e-5, 0, 0, 0, 1e-3 } + }); + + // When + final RealMatrix convertedCovarianceMatrixInInertialFrame = + StateCovarianceMatrixProvider.changeCovarianceFrame(initialOrbit, + LOFType.QSW, initialInertialFrame, + initialCovarianceInRTN); + + // Then + + // Expected covariance matrix obtained by rotation initial covariance matrix by -90 degrees + final RealMatrix expectedMatrixInInertialFrame = new BlockRealMatrix(new double[][] { + { 1, 0, 0, 0, 0, 1e-5 }, + { 0, 1, 0, 0, 0, 0 }, + { 0, 0, 1, 0, 0, 0 }, + { 0, 0, 0, 1e-3, 0, 0 }, + { 0, 0, 0, 0, 1e-3, 0 }, + { 1e-5, 0, 0, 0, 0, 1e-3 } + }); + + compareCovariance(expectedMatrixInInertialFrame, convertedCovarianceMatrixInInertialFrame, 1e-20); + + } + /** - * Compare two covariance matrices - * @param reference reference covariance - * @param computed computed covariance - * @param threshold threshold for comparison + * This test is based on the following paper : Covariance Transformations for Satellite Flight Dynamics Operations + * from David A. Vallado. + *

+ * More specifically, we're using the initial covariance matrix from p.14 and compare the computed result with the + * cartesian covariance in RSW from p.19. + *

+ * Note that the followings local orbital frame are equivalent RSW=RTN=QSW. */ - private void compareCovariance(final RealMatrix reference, final RealMatrix computed, final double threshold) { - for (int row = 0; row < reference.getRowDimension(); row++) { - for (int column = 0; column < reference.getColumnDimension(); column++) { - Assertions.assertEquals(reference.getEntry(row, column), computed.getEntry(row, column), FastMath.abs(threshold * reference.getEntry(row, column))); - } - } + @Test + @DisplayName("Test covariance conversion Vallado test case : ECI cartesian to RTN") + void should_return_Vallado_RSW_covariance_matrix_from_ECI() { + + // Initialize Orekit + Utils.setDataRoot("regular-data"); + + // Given + final AbsoluteDate initialDate = getValladoInitialDate(); + final PVCoordinates initialPV = getValladoInitialPV(); + final Frame initialInertialFrame = FramesFactory.getGCRF(); + final Orbit initialOrbit = + new CartesianOrbit(initialPV, initialInertialFrame, initialDate, getValladoMu()); + + final RealMatrix initialCovarianceMatrix = getValladoInitialCovarianceMatrix(); + + // When + final RealMatrix convertedCovarianceMatrixInRTN = + StateCovarianceMatrixProvider.changeCovarianceFrame(initialOrbit, + initialInertialFrame, LOFType.QSW, + initialCovarianceMatrix, + OrbitType.CARTESIAN, PositionAngle.MEAN); + + // Then + final RealMatrix expectedCovarianceMatrixInRTN = new BlockRealMatrix(new double[][] { + { 9.918921e-001, 6.700644e-003, -2.878187e-003, 1.892086e-005, 6.700644e-005, -2.878187e-005 }, + { 6.700644e-003, 1.013730e+000, -1.019283e-002, 6.700644e-005, 2.372970e-004, -1.019283e-004 }, + { -2.878187e-003, -1.019283e-002, 9.943782e-001, -2.878187e-005, -1.019283e-004, 4.378217e-005 }, + { 1.892086e-005, 6.700644e-005, -2.878187e-005, 1.892086e-007, 6.700644e-007, -2.878187e-007 }, + { 6.700644e-005, 2.372970e-004, -1.019283e-004, 6.700644e-007, 2.372970e-006, -1.019283e-006 }, + { -2.878187e-005, -1.019283e-004, 4.378217e-005, -2.878187e-007, -1.019283e-006, 4.378217e-007 } + }); + + compareCovariance(expectedCovarianceMatrixInRTN, convertedCovarianceMatrixInRTN, DEFAULT_VALLADO_THRESHOLD); + + } + + private AbsoluteDate getValladoInitialDate() { + return new AbsoluteDate(2000, 12, 15, 16, 58, 50.208, TimeScalesFactory.getUTC()); + } + + private PVCoordinates getValladoInitialPV() { + return new PVCoordinates( + new Vector3D(-605792.21660, -5870229.51108, 3493053.19896), + new Vector3D(-1568.25429, -3702.34891, -6479.48395)); + } + + private double getValladoMu() { + return Constants.IERS2010_EARTH_MU; + } + + private RealMatrix getValladoInitialCovarianceMatrix() { + return new BlockRealMatrix(new double[][] { + { 1, 1e-2, 1e-2, 1e-4, 1e-4, 1e-4 }, + { 1e-2, 1, 1e-2, 1e-4, 1e-4, 1e-4 }, + { 1e-2, 1e-2, 1, 1e-4, 1e-4, 1e-4 }, + { 1e-4, 1e-4, 1e-4, 1e-6, 1e-6, 1e-6 }, + { 1e-4, 1e-4, 1e-4, 1e-6, 1e-6, 1e-6 }, + { 1e-4, 1e-4, 1e-4, 1e-6, 1e-6, 1e-6 } + }); + } + + /** + * This test is based on the following paper : Covariance Transformations for Satellite Flight Dynamics Operations + * from David A. Vallado. + *

+ * More specifically, we're using the initial covariance matrix from p.14 and compare the computed result with the + * cartesian covariance in NTW from p.19. + */ + @Test + @DisplayName("Test covariance conversion Vallado test case : ECI cartesian to NTW") + void should_return_Vallado_NTW_covariance_matrix_from_ECI() { + + // Initialize orekit + Utils.setDataRoot("regular-data"); + + // Given + final AbsoluteDate initialDate = getValladoInitialDate(); + final PVCoordinates initialPV = getValladoInitialPV(); + final Frame initialInertialFrame = FramesFactory.getGCRF(); + final Orbit initialOrbit = + new CartesianOrbit(initialPV, initialInertialFrame, initialDate, getValladoMu()); + + final RealMatrix initialCovarianceMatrix = getValladoInitialCovarianceMatrix(); + + // When + final RealMatrix convertedCovarianceMatrixInNTW = + StateCovarianceMatrixProvider.changeCovarianceFrame(initialOrbit, + initialInertialFrame, LOFType.NTW, + initialCovarianceMatrix, + OrbitType.CARTESIAN, PositionAngle.MEAN); + + // Then + final RealMatrix expectedCovarianceMatrixInNTW = new BlockRealMatrix(new double[][] { + { 9.918792e-001, 6.679546e-003, -2.868345e-003, 1.879167e-005, 6.679546e-005, -2.868345e-005 }, + { 6.679546e-003, 1.013743e+000, -1.019560e-002, 6.679546e-005, 2.374262e-004, -1.019560e-004 }, + { -2.868345e-003, -1.019560e-002, 9.943782e-001, -2.868345e-005, -1.019560e-004, 4.378217e-005 }, + { 1.879167e-005, 6.679546e-005, -2.868345e-005, 1.879167e-007, 6.679546e-007, -2.868345e-007 }, + { 6.679546e-005, 2.374262e-004, -1.019560e-004, 6.679546e-007, 2.374262e-006, -1.019560e-006 }, + { -2.868345e-005, -1.019560e-004, 4.378217e-005, -2.868345e-007, -1.019560e-006, 4.378217e-007 } + }); + + compareCovariance(expectedCovarianceMatrixInNTW, convertedCovarianceMatrixInNTW, DEFAULT_VALLADO_THRESHOLD); + + } + + /** + * This test is based on the following paper : Covariance Transformations for Satellite Flight Dynamics Operations + * from David A. Vallado. + *

+ * More specifically, we're using the initial covariance matrix from p.14 as a reference to test multiple + * conversions. + *

+ * This test aims to verify the numerical precision after various conversions and serves as a non regression test + * for future updates. + *

+ * Also, note that the conversion from the RTN to TEME tests the fact that the orbit is initially expressed in GCRF + * while we want the covariance expressed in TEME. Hence, it tests that the rotation from RTN to TEME needs to be + * obtained by expressing the orbit PVCoordinates in the TEME frame (hence the use of orbit.gtPVCoordinates(frameOut) + * ,see relevant changeCovarianceFrame method). + */ + @Test + @DisplayName("Test custom covariance conversion Vallado test case : GCRF -> TEME -> IRTF -> NTW -> RTN -> ITRF -> GCRF") + void should_return_initial_covariance_after_multiple_conversion() { + + // Initialize orekit + Utils.setDataRoot("regular-data"); + + // Given + final AbsoluteDate initialDate = getValladoInitialDate(); + final PVCoordinates initialPV = getValladoInitialPV(); + final Frame initialInertialFrame = FramesFactory.getGCRF(); + final Orbit initialOrbit = + new CartesianOrbit(initialPV, initialInertialFrame, initialDate, getValladoMu()); + + final RealMatrix initialCovarianceMatrixInGCRF = getValladoInitialCovarianceMatrix(); + + final Frame teme = FramesFactory.getTEME(); + + final Frame itrf = FramesFactory.getITRF(IERSConventions.IERS_2010, false); + + // When + // GCRF -> TEME + final RealMatrix convertedCovarianceMatrixInTEME = + StateCovarianceMatrixProvider.changeCovarianceFrame(initialOrbit, + initialInertialFrame, teme, + initialCovarianceMatrixInGCRF, + OrbitType.CARTESIAN, PositionAngle.MEAN); + + // TEME -> ITRF + final RealMatrix convertedCovarianceMatrixInITRF = + StateCovarianceMatrixProvider.changeCovarianceFrame(initialOrbit, + teme, itrf, + convertedCovarianceMatrixInTEME, + OrbitType.CARTESIAN, PositionAngle.MEAN); + // ITRF -> NTW + final RealMatrix convertedCovarianceMatrixInNTW = + StateCovarianceMatrixProvider.changeCovarianceFrame(initialOrbit, + itrf, LOFType.NTW, + convertedCovarianceMatrixInITRF, + OrbitType.CARTESIAN, PositionAngle.MEAN); + // NTW -> RTN + final RealMatrix convertedCovarianceMatrixInRTN = + StateCovarianceMatrixProvider.changeCovarianceFrame(initialOrbit, + LOFType.NTW, LOFType.QSW, + convertedCovarianceMatrixInNTW); + // RTN -> ITRF + final RealMatrix convertedCovarianceMatrixBackInITRF = + StateCovarianceMatrixProvider.changeCovarianceFrame(initialOrbit, + LOFType.QSW, itrf, + convertedCovarianceMatrixInRTN); + + // ITRF -> TEME + final RealMatrix convertedCovarianceMatrixBackInTEME = + StateCovarianceMatrixProvider.changeCovarianceFrame(initialOrbit, + itrf, teme, + convertedCovarianceMatrixBackInITRF, + OrbitType.CARTESIAN, PositionAngle.MEAN); + // TEME -> GCRF + final RealMatrix convertedCovarianceMatrixInGCRF = + StateCovarianceMatrixProvider.changeCovarianceFrame(initialOrbit, + teme, initialInertialFrame, + convertedCovarianceMatrixBackInTEME, + OrbitType.CARTESIAN, PositionAngle.MEAN); + + // Then + final RealMatrix expectedCovarianceMatrixInGCRF = getValladoInitialCovarianceMatrix(); + + compareCovariance(expectedCovarianceMatrixInGCRF, convertedCovarianceMatrixInGCRF, 1e-12); + + } + + /** + * This test is based on the following paper : Covariance Transformations for Satellite Flight Dynamics Operations + * from David A. Vallado. + *

+ * More specifically, we're using the initial covariance matrix from p.14 and compare the computed result with the + * cartesian covariance in ECEF from p.18. + *

+ *

+ * BEWARE: It has been found that the earth rotation in this Vallado's case is given 1 million times slower than + * the expected value. This has been corrected and the expected covariance matrix is now the covariance matrix + * computed by Orekit given the similarities with Vallado's results. In addition, the small differences potentially + * come from the custom EOP that Vallado uses. Hence, this test can be considered as a non regression + * test. + *

+ */ + @Test + @DisplayName("Test covariance conversion Vallado test case : ECI cartesian to PEF") + void should_return_Vallado_PEF_covariance_matrix_from_ECI() { + + // Initialize orekit + Utils.setDataRoot("regular-data"); + + // Given + // Initialize orbit + final AbsoluteDate initialDate = getValladoInitialDate(); + final PVCoordinates initialPV = getValladoInitialPV(); + final Frame initialInertialFrame = FramesFactory.getGCRF(); + + final Orbit initialOrbit = new CartesianOrbit(initialPV, initialInertialFrame, initialDate, getValladoMu()); + + // Initialize input covariance matrix + final RealMatrix initialCovarianceMatrix = getValladoInitialCovarianceMatrix(); + + // Initialize output frame + final Frame outputFrame = FramesFactory.getITRF(IERSConventions.IERS_2010, true); + + // When + final RealMatrix convertedCovarianceMatrixInITRF = + StateCovarianceMatrixProvider.changeCovarianceFrame(initialOrbit, + initialInertialFrame, outputFrame, + initialCovarianceMatrix, + OrbitType.CARTESIAN, PositionAngle.MEAN); + + // Then + final RealMatrix expectedCovarianceMatrixInITRF = new BlockRealMatrix(new double[][] { + { 9.9340005761276870e-01, 7.5124999798868530e-03, 5.8312675007359050e-03, 3.4548396261054936e-05, + 2.6851237046859200e-06, 5.8312677693153940e-05 }, + { 7.5124999798868025e-03, 1.0065990293034541e+00, 1.2884310200351924e-02, 1.4852736004690684e-04, + 1.6544247282904867e-04, 1.2884310644320954e-04 }, + { 5.8312675007359040e-03, 1.2884310200351924e-02, 1.0000009130837746e+00, 5.9252211072590390e-05, + 1.2841787487219444e-04, 1.0000913090989617e-04 }, + { 3.4548396261054936e-05, 1.4852736004690686e-04, 5.9252211072590403e-05, 3.5631474857130520e-07, + 7.6083489184819870e-07, 5.9252213790760030e-07 }, + { 2.6851237046859150e-06, 1.6544247282904864e-04, 1.2841787487219447e-04, 7.6083489184819880e-07, + 1.6542289254142709e-06, 1.2841787929229964e-06 }, + { 5.8312677693153934e-05, 1.2884310644320950e-04, 1.0000913090989616e-04, 5.9252213790760020e-07, + 1.2841787929229960e-06, 1.0000913098203875e-06 } + }); + + compareCovariance(expectedCovarianceMatrixInITRF, convertedCovarianceMatrixInITRF, 1e-20); + + } + + /** + * This test is based on the following paper : Covariance Transformations for Satellite Flight Dynamics Operations + * from David A. Vallado. + *

+ * More specifically, we're using the initial covariance matrix from p.14 and compare the computed result with the + * cartesian covariance in PEF from p.18. + */ + @Test + @DisplayName("Test covariance conversion Vallado test case : PEF cartesian to ECI") + void should_return_Vallado_ECI_covariance_matrix_from_PEF() { + + // Initialize orekit + Utils.setDataRoot("regular-data"); + + // Given + final AbsoluteDate initialDate = getValladoInitialDate(); + final PVCoordinates initialPV = getValladoInitialPV(); + final Frame initialInertialFrame = FramesFactory.getGCRF(); + final Orbit initialOrbit = + new CartesianOrbit(initialPV, initialInertialFrame, initialDate, getValladoMu()); + + final RealMatrix initialCovarianceMatrixInPEF = new BlockRealMatrix(new double[][] { + { 9.9340005761276870e-01, 7.5124999798868530e-03, 5.8312675007359050e-03, 3.4548396261054936e-05, + 2.6851237046859200e-06, 5.8312677693153940e-05 }, + { 7.5124999798868025e-03, 1.0065990293034541e+00, 1.2884310200351924e-02, 1.4852736004690684e-04, + 1.6544247282904867e-04, 1.2884310644320954e-04 }, + { 5.8312675007359040e-03, 1.2884310200351924e-02, 1.0000009130837746e+00, 5.9252211072590390e-05, + 1.2841787487219444e-04, 1.0000913090989617e-04 }, + { 3.4548396261054936e-05, 1.4852736004690686e-04, 5.9252211072590403e-05, 3.5631474857130520e-07, + 7.6083489184819870e-07, 5.9252213790760030e-07 }, + { 2.6851237046859150e-06, 1.6544247282904864e-04, 1.2841787487219447e-04, 7.6083489184819880e-07, + 1.6542289254142709e-06, 1.2841787929229964e-06 }, + { 5.8312677693153934e-05, 1.2884310644320950e-04, 1.0000913090989616e-04, 5.9252213790760020e-07, + 1.2841787929229960e-06, 1.0000913098203875e-06 } + }); + + final Frame inputFrame = FramesFactory.getITRF(IERSConventions.IERS_2010, false); + + // When + final RealMatrix convertedCovarianceMatrixInECI = + StateCovarianceMatrixProvider.changeCovarianceFrame(initialOrbit, + inputFrame, initialInertialFrame, + initialCovarianceMatrixInPEF, + OrbitType.CARTESIAN, PositionAngle.MEAN); + + // Then + final RealMatrix expectedCovarianceMatrixInECI = getValladoInitialCovarianceMatrix(); + + compareCovariance(expectedCovarianceMatrixInECI, convertedCovarianceMatrixInECI, 1e-7); + + } + + /** + * This test is based on the following paper : Covariance Transformations for Satellite Flight Dynamics Operations + * from David A. Vallado. + *

+ * More specifically, we're using the initial covariance matrix from p.14 and compare the computed result with the + * cartesian covariance in MOD from p.17. + *

+ */ + @Test + @DisplayName("Test covariance conversion Vallado test case : ECI cartesian to MOD") + void should_return_Vallado_MOD_covariance_matrix_from_ECI() { + + // Initialize orekit + Utils.setDataRoot("regular-data"); + + // Given + final AbsoluteDate initialDate = getValladoInitialDate(); + final PVCoordinates initialPV = getValladoInitialPV(); + final Frame initialInertialFrame = FramesFactory.getGCRF(); + final Orbit initialOrbit = + new CartesianOrbit(initialPV, initialInertialFrame, initialDate, getValladoMu()); + + final RealMatrix initialCovarianceMatrix = getValladoInitialCovarianceMatrix(); + + final Frame outputFrame = FramesFactory.getMOD(IERSConventions.IERS_2010); + + // When + final RealMatrix convertedCovarianceMatrixInMOD = + StateCovarianceMatrixProvider.changeCovarianceFrame(initialOrbit, + initialInertialFrame, outputFrame, + initialCovarianceMatrix, + OrbitType.CARTESIAN, PositionAngle.MEAN); + + // Then + final RealMatrix expectedCovarianceMatrixInMOD = new BlockRealMatrix(new double[][] { + { 9.999939e-001, 9.999070e-003, 9.997861e-003, 9.993866e-005, 9.999070e-005, 9.997861e-005 }, + { 9.999070e-003, 1.000004e+000, 1.000307e-002, 9.999070e-005, 1.000428e-004, 1.000307e-004 }, + { 9.997861e-003, 1.000307e-002, 1.000002e+000, 9.997861e-005, 1.000307e-004, 1.000186e-004 }, + { 9.993866e-005, 9.999070e-005, 9.997861e-005, 9.993866e-007, 9.999070e-007, 9.997861e-007 }, + { 9.999070e-005, 1.000428e-004, 1.000307e-004, 9.999070e-007, 1.000428e-006, 1.000307e-006 }, + { 9.997861e-005, 1.000307e-004, 1.000186e-004, 9.997861e-007, 1.000307e-006, 1.000186e-006 }, + }); + + compareCovariance(expectedCovarianceMatrixInMOD, convertedCovarianceMatrixInMOD, DEFAULT_VALLADO_THRESHOLD); + + } + + /** + * This test is based on the following paper : Covariance Transformations for Satellite Flight Dynamics Operations + * from David A. Vallado. + *

+ * More specifically, we're using the initial NTW covariance matrix from p.19 and compare the computed result with + * the cartesian covariance in RSW from the same page. + *

+ */ + @Test + @DisplayName("Test covariance conversion from Vallado test case NTW to RSW") + void should_convert_Vallado_NTW_to_RSW() { + + // Initialize orekit + Utils.setDataRoot("regular-data"); + + // Given + final AbsoluteDate initialDate = getValladoInitialDate(); + final PVCoordinates initialPV = getValladoInitialPV(); + final Frame initialInertialFrame = FramesFactory.getGCRF(); + final Orbit initialOrbit = + new CartesianOrbit(initialPV, initialInertialFrame, initialDate, getValladoMu()); + + final RealMatrix initialCovarianceMatrixInNTW = new BlockRealMatrix(new double[][] { + { 9.918792e-001, 6.679546e-003, -2.868345e-003, 1.879167e-005, 6.679546e-005, -2.868345e-005 }, + { 6.679546e-003, 1.013743e+000, -1.019560e-002, 6.679546e-005, 2.374262e-004, -1.019560e-004 }, + { -2.868345e-003, -1.019560e-002, 9.943782e-001, -2.868345e-005, -1.019560e-004, 4.378217e-005 }, + { 1.879167e-005, 6.679546e-005, -2.868345e-005, 1.879167e-007, 6.679546e-007, -2.868345e-007 }, + { 6.679546e-005, 2.374262e-004, -1.019560e-004, 6.679546e-007, 2.374262e-006, -1.019560e-006 }, + { -2.868345e-005, -1.019560e-004, 4.378217e-005, -2.868345e-007, -1.019560e-006, 4.378217e-007 } + }); + + // When + final RealMatrix convertedCovarianceMatrixInRTN = + StateCovarianceMatrixProvider.changeCovarianceFrame(initialOrbit, LOFType.NTW, LOFType.QSW, + initialCovarianceMatrixInNTW); + + // Then + final RealMatrix expectedCovarianceMatrixInRTN = new BlockRealMatrix(new double[][] { + { 9.918921e-001, 6.700644e-003, -2.878187e-003, 1.892086e-005, 6.700644e-005, -2.878187e-005 }, + { 6.700644e-003, 1.013730e+000, -1.019283e-002, 6.700644e-005, 2.372970e-004, -1.019283e-004 }, + { -2.878187e-003, -1.019283e-002, 9.943782e-001, -2.878187e-005, -1.019283e-004, 4.378217e-005 }, + { 1.892086e-005, 6.700644e-005, -2.878187e-005, 1.892086e-007, 6.700644e-007, -2.878187e-007 }, + { 6.700644e-005, 2.372970e-004, -1.019283e-004, 6.700644e-007, 2.372970e-006, -1.019283e-006 }, + { -2.878187e-005, -1.019283e-004, 4.378217e-005, -2.878187e-007, -1.019283e-006, 4.378217e-007 } + }); + + compareCovariance(expectedCovarianceMatrixInRTN, convertedCovarianceMatrixInRTN, DEFAULT_VALLADO_THRESHOLD); + + } + + /** + * This test is based on the following paper : Covariance Transformations for Satellite Flight Dynamics Operations + * from David A. Vallado. + *

+ * More specifically, we're using the initial NTW covariance matrix from p.19 and compare the computed result with + * the cartesian covariance in RSW from the same page. + *

+ */ + @Test + @DisplayName("Test thrown error if input frame is not pseudo-inertial and the covariance matrix is not expressed in cartesian elements") + void should_return_orekit_exception() { + + // Initialize orekit + Utils.setDataRoot("regular-data"); + + // Given + final AbsoluteDate initialDate = getValladoInitialDate(); + final PVCoordinates initialPV = getValladoInitialPV(); + final Frame initialInertialFrame = FramesFactory.getGCRF(); + final Orbit initialOrbit = + new CartesianOrbit(initialPV, initialInertialFrame, initialDate, getValladoMu()); + + final RealMatrix randomCovarianceMatrix = new BlockRealMatrix(new double[][] { + { 9.918792e-001, 6.679546e-003, -2.868345e-003, 1.879167e-005, 6.679546e-005, -2.868345e-005 }, + { 6.679546e-003, 1.013743e+000, -1.019560e-002, 6.679546e-005, 2.374262e-004, -1.019560e-004 }, + { -2.868345e-003, -1.019560e-002, 9.943782e-001, -2.868345e-005, -1.019560e-004, 4.378217e-005 }, + { 1.879167e-005, 6.679546e-005, -2.868345e-005, 1.879167e-007, 6.679546e-007, -2.868345e-007 }, + { 6.679546e-005, 2.374262e-004, -1.019560e-004, 6.679546e-007, 2.374262e-006, -1.019560e-006 }, + { -2.868345e-005, -1.019560e-004, 4.378217e-005, -2.868345e-007, -1.019560e-006, 4.378217e-007 } + }); + + final Frame nonInertialFrame = FramesFactory.getITRF(IERSConventions.IERS_2010, true); + + final Frame inertialFrame = FramesFactory.getGCRF(); + + // When & Then + Assertions.assertThrows(OrekitException.class, + () -> { + StateCovarianceMatrixProvider.changeCovarianceFrame( + initialOrbit, nonInertialFrame, inertialFrame, + randomCovarianceMatrix, OrbitType.CIRCULAR, + PositionAngle.MEAN); + }); + + Assertions.assertThrows(OrekitException.class, + () -> { + StateCovarianceMatrixProvider.changeCovarianceFrame( + initialOrbit, nonInertialFrame, LOFType.QSW, + randomCovarianceMatrix, OrbitType.EQUINOCTIAL, + PositionAngle.MEAN); + }); + } }