Transform.java

  1. /* Copyright 2002-2024 CS GROUP
  2.  * Licensed to CS GROUP (CS) under one or more
  3.  * contributor license agreements.  See the NOTICE file distributed with
  4.  * this work for additional information regarding copyright ownership.
  5.  * CS licenses this file to You under the Apache License, Version 2.0
  6.  * (the "License"); you may not use this file except in compliance with
  7.  * the License.  You may obtain a copy of the License at
  8.  *
  9.  *   http://www.apache.org/licenses/LICENSE-2.0
  10.  *
  11.  * Unless required by applicable law or agreed to in writing, software
  12.  * distributed under the License is distributed on an "AS IS" BASIS,
  13.  * WITHOUT WARRANTIES OR CONDITIONS OF ANY KIND, either express or implied.
  14.  * See the License for the specific language governing permissions and
  15.  * limitations under the License.
  16.  */
  17. package org.orekit.frames;

  18. import java.io.Serializable;
  19. import java.util.ArrayList;
  20. import java.util.Arrays;
  21. import java.util.Collection;
  22. import java.util.List;
  23. import java.util.stream.Collectors;
  24. import java.util.stream.Stream;

  25. import org.hipparchus.CalculusFieldElement;
  26. import org.hipparchus.geometry.euclidean.threed.Line;
  27. import org.hipparchus.geometry.euclidean.threed.Rotation;
  28. import org.hipparchus.geometry.euclidean.threed.Vector3D;
  29. import org.orekit.time.AbsoluteDate;
  30. import org.orekit.time.TimeInterpolator;
  31. import org.orekit.time.TimeShiftable;
  32. import org.orekit.utils.AngularCoordinates;
  33. import org.orekit.utils.AngularDerivativesFilter;
  34. import org.orekit.utils.CartesianDerivativesFilter;
  35. import org.orekit.utils.FieldPVCoordinates;
  36. import org.orekit.utils.PVCoordinates;
  37. import org.orekit.utils.TimeStampedAngularCoordinates;
  38. import org.orekit.utils.TimeStampedAngularCoordinatesHermiteInterpolator;
  39. import org.orekit.utils.TimeStampedFieldPVCoordinates;
  40. import org.orekit.utils.TimeStampedPVCoordinates;
  41. import org.orekit.utils.TimeStampedPVCoordinatesHermiteInterpolator;


  42. /** Transformation class in three dimensional space.
  43.  *
  44.  * <p>This class represents the transformation engine between {@link Frame frames}.
  45.  * It is used both to define the relationship between each frame and its
  46.  * parent frame and to gather all individual transforms into one
  47.  * operation when converting between frames far away from each other.</p>
  48.  * <p> The convention used in OREKIT is vectorial transformation. It means
  49.  * that a transformation is defined as a transform to apply to the
  50.  * coordinates of a vector expressed in the old frame to obtain the
  51.  * same vector expressed in the new frame.
  52.  *
  53.  * <p>Instances of this class are guaranteed to be immutable.</p>
  54.  *
  55.  * <h2> Examples </h2>
  56.  *
  57.  * <h3> Example of translation from R<sub>A</sub> to R<sub>B</sub> </h3>
  58.  *
  59.  * <p> We want to transform the {@link PVCoordinates} PV<sub>A</sub> to
  60.  * PV<sub>B</sub> with :
  61.  * <p> PV<sub>A</sub> = ({1, 0, 0}, {2, 0, 0}, {3, 0, 0}); <br>
  62.  *     PV<sub>B</sub> = ({0, 0, 0}, {0, 0, 0}, {0, 0, 0});
  63.  *
  64.  * <p> The transform to apply then is defined as follows :
  65.  *
  66.  * <pre><code>
  67.  * Vector3D translation  = new Vector3D(-1, 0, 0);
  68.  * Vector3D velocity     = new Vector3D(-2, 0, 0);
  69.  * Vector3D acceleration = new Vector3D(-3, 0, 0);
  70.  *
  71.  * Transform R1toR2 = new Transform(date, translation, velocity, acceleration);
  72.  *
  73.  * PVB = R1toR2.transformPVCoordinates(PVA);
  74.  * </code></pre>
  75.  *
  76.  * <h3> Example of rotation from R<sub>A</sub> to R<sub>B</sub> </h3>
  77.  * <p> We want to transform the {@link PVCoordinates} PV<sub>A</sub> to
  78.  * PV<sub>B</sub> with
  79.  *
  80.  * <p> PV<sub>A</sub> = ({1, 0, 0}, { 1, 0, 0}); <br>
  81.  *     PV<sub>B</sub> = ({0, 1, 0}, {-2, 1, 0});
  82.  *
  83.  * <p> The transform to apply then is defined as follows :
  84.  *
  85.  * <pre><code>
  86.  * Rotation rotation = new Rotation(Vector3D.PLUS_K, FastMath.PI / 2);
  87.  * Vector3D rotationRate = new Vector3D(0, 0, -2);
  88.  *
  89.  * Transform R1toR2 = new Transform(rotation, rotationRate);
  90.  *
  91.  * PVB = R1toR2.transformPVCoordinates(PVA);
  92.  * </code></pre>
  93.  *
  94.  * @author Luc Maisonobe
  95.  * @author Fabien Maussion
  96.  */
  97. public class Transform implements
  98.         TimeShiftable<Transform>,
  99.         Serializable,
  100.         KinematicTransform {

  101.     /** Identity transform. */
  102.     public static final Transform IDENTITY = new IdentityTransform();

  103.     /** Serializable UID. */
  104.     private static final long serialVersionUID = 210140410L;

  105.     /** Date of the transform. */
  106.     private final AbsoluteDate date;

  107.     /** Cartesian coordinates of the target frame with respect to the original frame. */
  108.     private final PVCoordinates cartesian;

  109.     /** Angular coordinates of the target frame with respect to the original frame. */
  110.     private final AngularCoordinates angular;

  111.     /** Build a transform from its primitive operations.
  112.      * @param date date of the transform
  113.      * @param cartesian Cartesian coordinates of the target frame with respect to the original frame
  114.      * @param angular angular coordinates of the target frame with respect to the original frame
  115.      */
  116.     public Transform(final AbsoluteDate date, final PVCoordinates cartesian, final AngularCoordinates angular) {
  117.         this.date      = date;
  118.         this.cartesian = cartesian;
  119.         this.angular   = angular;
  120.     }

  121.     /** Build a translation transform.
  122.      * @param date date of the transform
  123.      * @param translation translation to apply (i.e. coordinates of
  124.      * the transformed origin, or coordinates of the origin of the
  125.      * old frame in the new frame)
  126.      */
  127.     public Transform(final AbsoluteDate date, final Vector3D translation) {
  128.         this(date,
  129.              new PVCoordinates(translation),
  130.              AngularCoordinates.IDENTITY);
  131.     }

  132.     /** Build a rotation transform.
  133.      * @param date date of the transform
  134.      * @param rotation rotation to apply ( i.e. rotation to apply to the
  135.      * coordinates of a vector expressed in the old frame to obtain the
  136.      * same vector expressed in the new frame )
  137.      */
  138.     public Transform(final AbsoluteDate date, final Rotation rotation) {
  139.         this(date,
  140.              PVCoordinates.ZERO,
  141.              new AngularCoordinates(rotation));
  142.     }

  143.     /** Build a combined translation and rotation transform.
  144.      * @param date date of the transform
  145.      * @param translation translation to apply (i.e. coordinates of
  146.      * the transformed origin, or coordinates of the origin of the
  147.      * old frame in the new frame)
  148.      * @param rotation rotation to apply ( i.e. rotation to apply to the
  149.      * coordinates of a vector expressed in the old frame to obtain the
  150.      * same vector expressed in the new frame )
  151.      * @since 12.1
  152.      */
  153.     public Transform(final AbsoluteDate date, final Vector3D translation, final Rotation rotation) {
  154.         this(date, new PVCoordinates(translation), new AngularCoordinates(rotation));
  155.     }

  156.     /** Build a translation transform, with its first time derivative.
  157.      * @param date date of the transform
  158.      * @param translation translation to apply (i.e. coordinates of
  159.      * the transformed origin, or coordinates of the origin of the
  160.      * old frame in the new frame)
  161.      * @param velocity the velocity of the translation (i.e. origin
  162.      * of the old frame velocity in the new frame)
  163.      */
  164.     public Transform(final AbsoluteDate date, final Vector3D translation,
  165.                      final Vector3D velocity) {
  166.         this(date,
  167.              new PVCoordinates(translation, velocity, Vector3D.ZERO),
  168.              AngularCoordinates.IDENTITY);
  169.     }

  170.     /** Build a translation transform, with its first and second time derivatives.
  171.      * @param date date of the transform
  172.      * @param translation translation to apply (i.e. coordinates of
  173.      * the transformed origin, or coordinates of the origin of the
  174.      * old frame in the new frame)
  175.      * @param velocity the velocity of the translation (i.e. origin
  176.      * of the old frame velocity in the new frame)
  177.      * @param acceleration the acceleration of the translation (i.e. origin
  178.      * of the old frame acceleration in the new frame)
  179.      */
  180.     public Transform(final AbsoluteDate date, final Vector3D translation,
  181.                      final Vector3D velocity, final Vector3D acceleration) {
  182.         this(date,
  183.              new PVCoordinates(translation, velocity, acceleration),
  184.              AngularCoordinates.IDENTITY);
  185.     }

  186.     /** Build a translation transform, with its first time derivative.
  187.      * @param date date of the transform
  188.      * @param cartesian Cartesian part of the transformation to apply (i.e. coordinates of
  189.      * the transformed origin, or coordinates of the origin of the
  190.      * old frame in the new frame, with their derivatives)
  191.      */
  192.     public Transform(final AbsoluteDate date, final PVCoordinates cartesian) {
  193.         this(date,
  194.              cartesian,
  195.              AngularCoordinates.IDENTITY);
  196.     }

  197.     /** Build a rotation transform.
  198.      * @param date date of the transform
  199.      * @param rotation rotation to apply ( i.e. rotation to apply to the
  200.      * coordinates of a vector expressed in the old frame to obtain the
  201.      * same vector expressed in the new frame )
  202.      * @param rotationRate the axis of the instant rotation
  203.      * expressed in the new frame. (norm representing angular rate)
  204.      */
  205.     public Transform(final AbsoluteDate date, final Rotation rotation, final Vector3D rotationRate) {
  206.         this(date,
  207.              PVCoordinates.ZERO,
  208.              new AngularCoordinates(rotation, rotationRate, Vector3D.ZERO));
  209.     }

  210.     /** Build a rotation transform.
  211.      * @param date date of the transform
  212.      * @param rotation rotation to apply ( i.e. rotation to apply to the
  213.      * coordinates of a vector expressed in the old frame to obtain the
  214.      * same vector expressed in the new frame )
  215.      * @param rotationRate the axis of the instant rotation
  216.      * @param rotationAcceleration the axis of the instant rotation
  217.      * expressed in the new frame. (norm representing angular rate)
  218.      */
  219.     public Transform(final AbsoluteDate date, final Rotation rotation, final Vector3D rotationRate,
  220.                      final Vector3D rotationAcceleration) {
  221.         this(date,
  222.              PVCoordinates.ZERO,
  223.              new AngularCoordinates(rotation, rotationRate, rotationAcceleration));
  224.     }

  225.     /** Build a rotation transform.
  226.      * @param date date of the transform
  227.      * @param angular angular part of the transformation to apply (i.e. rotation to
  228.      * apply to the coordinates of a vector expressed in the old frame to obtain the
  229.      * same vector expressed in the new frame, with its rotation rate)
  230.      */
  231.     public Transform(final AbsoluteDate date, final AngularCoordinates angular) {
  232.         this(date, PVCoordinates.ZERO, angular);
  233.     }

  234.     /** Build a transform by combining two existing ones.
  235.      * <p>
  236.      * Note that the dates of the two existing transformed are <em>ignored</em>,
  237.      * and the combined transform date is set to the date supplied in this constructor
  238.      * without any attempt to shift the raw transforms. This is a design choice allowing
  239.      * user full control of the combination.
  240.      * </p>
  241.      * @param date date of the transform
  242.      * @param first first transform applied
  243.      * @param second second transform applied
  244.      */
  245.     public Transform(final AbsoluteDate date, final Transform first, final Transform second) {
  246.         this(date,
  247.              new PVCoordinates(StaticTransform.compositeTranslation(first, second),
  248.                                KinematicTransform.compositeVelocity(first, second),
  249.                                compositeAcceleration(first, second)),
  250.              new AngularCoordinates(StaticTransform.compositeRotation(first, second),
  251.                                     KinematicTransform.compositeRotationRate(first, second),
  252.                                     compositeRotationAcceleration(first, second)));
  253.     }

  254.     /** Compute a composite acceleration.
  255.      * @param first first applied transform
  256.      * @param second second applied transform
  257.      * @return acceleration part of the composite transform
  258.      */
  259.     private static Vector3D compositeAcceleration(final Transform first, final Transform second) {

  260.         final Vector3D a1    = first.cartesian.getAcceleration();
  261.         final Rotation r1    = first.angular.getRotation();
  262.         final Vector3D o1    = first.angular.getRotationRate();
  263.         final Vector3D oDot1 = first.angular.getRotationAcceleration();
  264.         final Vector3D p2    = second.cartesian.getPosition();
  265.         final Vector3D v2    = second.cartesian.getVelocity();
  266.         final Vector3D a2    = second.cartesian.getAcceleration();

  267.         final Vector3D crossCrossP = Vector3D.crossProduct(o1,    Vector3D.crossProduct(o1, p2));
  268.         final Vector3D crossV      = Vector3D.crossProduct(o1,    v2);
  269.         final Vector3D crossDotP   = Vector3D.crossProduct(oDot1, p2);

  270.         return a1.add(r1.applyInverseTo(new Vector3D(1, a2, 2, crossV, 1, crossCrossP, 1, crossDotP)));

  271.     }

  272.     /** Compute a composite rotation acceleration.
  273.      * @param first first applied transform
  274.      * @param second second applied transform
  275.      * @return rotation acceleration part of the composite transform
  276.      */
  277.     private static Vector3D compositeRotationAcceleration(final Transform first, final Transform second) {

  278.         final Vector3D o1    = first.angular.getRotationRate();
  279.         final Vector3D oDot1 = first.angular.getRotationAcceleration();
  280.         final Rotation r2    = second.angular.getRotation();
  281.         final Vector3D o2    = second.angular.getRotationRate();
  282.         final Vector3D oDot2 = second.angular.getRotationAcceleration();

  283.         return new Vector3D( 1, oDot2,
  284.                              1, r2.applyTo(oDot1),
  285.                             -1, Vector3D.crossProduct(o2, r2.applyTo(o1)));

  286.     }

  287.     /** {@inheritDoc} */
  288.     public AbsoluteDate getDate() {
  289.         return date;
  290.     }

  291.     /** {@inheritDoc} */
  292.     public Transform shiftedBy(final double dt) {
  293.         return new Transform(date.shiftedBy(dt), cartesian.shiftedBy(dt), angular.shiftedBy(dt));
  294.     }

  295.     /**
  296.      * Shift the transform in time considering all rates, then return only the
  297.      * translation and rotation portion of the transform.
  298.      *
  299.      * @param dt time shift in seconds.
  300.      * @return shifted transform as a static transform. It is static in the
  301.      * sense that it can only be used to transform directions and positions, but
  302.      * not velocities or accelerations.
  303.      * @see #shiftedBy(double)
  304.      */
  305.     public StaticTransform staticShiftedBy(final double dt) {
  306.         return StaticTransform.of(
  307.                 date.shiftedBy(dt),
  308.                 cartesian.positionShiftedBy(dt),
  309.                 angular.rotationShiftedBy(dt));
  310.     }

  311.     /**
  312.      * Create a so-called static transform from the instance.
  313.      *
  314.      * @return static part of the transform. It is static in the
  315.      * sense that it can only be used to transform directions and positions, but
  316.      * not velocities or accelerations.
  317.      * @see StaticTransform
  318.      */
  319.     public StaticTransform toStaticTransform() {
  320.         return StaticTransform.of(date, cartesian.getPosition(), angular.getRotation());
  321.     }

  322.     /** Interpolate a transform from a sample set of existing transforms.
  323.      * <p>
  324.      * Calling this method is equivalent to call {@link #interpolate(AbsoluteDate,
  325.      * CartesianDerivativesFilter, AngularDerivativesFilter, Collection)} with {@code cFilter}
  326.      * set to {@link CartesianDerivativesFilter#USE_PVA} and {@code aFilter} set to
  327.      * {@link AngularDerivativesFilter#USE_RRA}
  328.      * set to true.
  329.      * </p>
  330.      * @param interpolationDate interpolation date
  331.      * @param sample sample points on which interpolation should be done
  332.      * @return a new instance, interpolated at specified date
  333.      */
  334.     public Transform interpolate(final AbsoluteDate interpolationDate, final Stream<Transform> sample) {
  335.         return interpolate(interpolationDate,
  336.                            CartesianDerivativesFilter.USE_PVA, AngularDerivativesFilter.USE_RRA,
  337.                            sample.collect(Collectors.toList()));
  338.     }

  339.     /** Interpolate a transform from a sample set of existing transforms.
  340.      * <p>
  341.      * Note that even if first time derivatives (velocities and rotation rates)
  342.      * from sample can be ignored, the interpolated instance always includes
  343.      * interpolated derivatives. This feature can be used explicitly to
  344.      * compute these derivatives when it would be too complex to compute them
  345.      * from an analytical formula: just compute a few sample points from the
  346.      * explicit formula and set the derivatives to zero in these sample points,
  347.      * then use interpolation to add derivatives consistent with the positions
  348.      * and rotations.
  349.      * </p>
  350.      * <p>
  351.      * As this implementation of interpolation is polynomial, it should be used only
  352.      * with small samples (about 10-20 points) in order to avoid <a
  353.      * href="http://en.wikipedia.org/wiki/Runge%27s_phenomenon">Runge's phenomenon</a>
  354.      * and numerical problems (including NaN appearing).
  355.      * </p>
  356.      * @param date interpolation date
  357.      * @param cFilter filter for derivatives from the sample to use in interpolation
  358.      * @param aFilter filter for derivatives from the sample to use in interpolation
  359.      * @param sample sample points on which interpolation should be done
  360.      * @return a new instance, interpolated at specified date
  361.      * @since 7.0
  362.      */
  363.     public static Transform interpolate(final AbsoluteDate date,
  364.                                         final CartesianDerivativesFilter cFilter,
  365.                                         final AngularDerivativesFilter aFilter,
  366.                                         final Collection<Transform> sample) {

  367.         // Create samples
  368.         final List<TimeStampedPVCoordinates>      datedPV = new ArrayList<>(sample.size());
  369.         final List<TimeStampedAngularCoordinates> datedAC = new ArrayList<>(sample.size());
  370.         for (final Transform t : sample) {
  371.             datedPV.add(new TimeStampedPVCoordinates(t.getDate(), t.getTranslation(), t.getVelocity(), t.getAcceleration()));
  372.             datedAC.add(new TimeStampedAngularCoordinates(t.getDate(), t.getRotation(), t.getRotationRate(), t.getRotationAcceleration()));
  373.         }

  374.         // Create interpolators
  375.         final TimeInterpolator<TimeStampedPVCoordinates> pvInterpolator =
  376.                 new TimeStampedPVCoordinatesHermiteInterpolator(datedPV.size(), cFilter);

  377.         final TimeInterpolator<TimeStampedAngularCoordinates> angularInterpolator =
  378.                 new TimeStampedAngularCoordinatesHermiteInterpolator(datedPV.size(), aFilter);

  379.         // Interpolate
  380.         final TimeStampedPVCoordinates      interpolatedPV = pvInterpolator.interpolate(date, datedPV);
  381.         final TimeStampedAngularCoordinates interpolatedAC = angularInterpolator.interpolate(date, datedAC);
  382.         return new Transform(date, interpolatedPV, interpolatedAC);
  383.     }

  384.     /** Get the inverse transform of the instance.
  385.      * @return inverse transform of the instance
  386.      */
  387.     @Override
  388.     public Transform getInverse() {

  389.         final Rotation r    = angular.getRotation();
  390.         final Vector3D o    = angular.getRotationRate();
  391.         final Vector3D oDot = angular.getRotationAcceleration();
  392.         final Vector3D rp   = r.applyTo(cartesian.getPosition());
  393.         final Vector3D rv   = r.applyTo(cartesian.getVelocity());
  394.         final Vector3D ra   = r.applyTo(cartesian.getAcceleration());

  395.         final Vector3D pInv        = rp.negate();
  396.         final Vector3D crossP      = Vector3D.crossProduct(o, rp);
  397.         final Vector3D vInv        = crossP.subtract(rv);
  398.         final Vector3D crossV      = Vector3D.crossProduct(o, rv);
  399.         final Vector3D crossDotP   = Vector3D.crossProduct(oDot, rp);
  400.         final Vector3D crossCrossP = Vector3D.crossProduct(o, crossP);
  401.         final Vector3D aInv        = new Vector3D(-1, ra,
  402.                                                    2, crossV,
  403.                                                    1, crossDotP,
  404.                                                   -1, crossCrossP);

  405.         return new Transform(getDate(), new PVCoordinates(pInv, vInv, aInv), angular.revert());

  406.     }

  407.     /** Get a frozen transform.
  408.      * <p>
  409.      * This method creates a copy of the instance but frozen in time,
  410.      * i.e. with velocity, acceleration and rotation rate forced to zero.
  411.      * </p>
  412.      * @return a new transform, without any time-dependent parts
  413.      */
  414.     public Transform freeze() {
  415.         return new Transform(date,
  416.                              new PVCoordinates(cartesian.getPosition(), Vector3D.ZERO, Vector3D.ZERO),
  417.                              new AngularCoordinates(angular.getRotation(), Vector3D.ZERO, Vector3D.ZERO));
  418.     }

  419.     /** Transform {@link PVCoordinates} including kinematic effects.
  420.      * @param pva the position-velocity-acceleration triplet to transform.
  421.      * @return transformed position-velocity-acceleration
  422.      */
  423.     public PVCoordinates transformPVCoordinates(final PVCoordinates pva) {
  424.         return angular.applyTo(new PVCoordinates(1, pva, 1, cartesian));
  425.     }

  426.     /** Transform {@link TimeStampedPVCoordinates} including kinematic effects.
  427.      * <p>
  428.      * In order to allow the user more flexibility, this method does <em>not</em> check for
  429.      * consistency between the transform {@link #getDate() date} and the time-stamped
  430.      * position-velocity {@link TimeStampedPVCoordinates#getDate() date}. The returned
  431.      * value will always have the same {@link TimeStampedPVCoordinates#getDate() date} as
  432.      * the input argument, regardless of the instance {@link #getDate() date}.
  433.      * </p>
  434.      * @param pv time-stamped  position-velocity to transform.
  435.      * @return transformed time-stamped position-velocity
  436.      * @since 7.0
  437.      */
  438.     public TimeStampedPVCoordinates transformPVCoordinates(final TimeStampedPVCoordinates pv) {
  439.         return angular.applyTo(new TimeStampedPVCoordinates(pv.getDate(), 1, pv, 1, cartesian));
  440.     }

  441.     /** Transform {@link FieldPVCoordinates} including kinematic effects.
  442.      * @param pv position-velocity to transform.
  443.      * @param <T> type of the field elements
  444.      * @return transformed position-velocity
  445.      */
  446.     public <T extends CalculusFieldElement<T>> FieldPVCoordinates<T> transformPVCoordinates(final FieldPVCoordinates<T> pv) {
  447.         return angular.applyTo(new FieldPVCoordinates<>(pv.getPosition().add(cartesian.getPosition()),
  448.                                                         pv.getVelocity().add(cartesian.getVelocity()),
  449.                                                         pv.getAcceleration().add(cartesian.getAcceleration())));
  450.     }

  451.     /** Transform {@link TimeStampedFieldPVCoordinates} including kinematic effects.
  452.      * <p>
  453.      * In order to allow the user more flexibility, this method does <em>not</em> check for
  454.      * consistency between the transform {@link #getDate() date} and the time-stamped
  455.      * position-velocity {@link TimeStampedFieldPVCoordinates#getDate() date}. The returned
  456.      * value will always have the same {@link TimeStampedFieldPVCoordinates#getDate() date} as
  457.      * the input argument, regardless of the instance {@link #getDate() date}.
  458.      * </p>
  459.      * @param pv time-stamped position-velocity to transform.
  460.      * @param <T> type of the field elements
  461.      * @return transformed time-stamped position-velocity
  462.      * @since 7.0
  463.      */
  464.     public <T extends CalculusFieldElement<T>> TimeStampedFieldPVCoordinates<T> transformPVCoordinates(final TimeStampedFieldPVCoordinates<T> pv) {
  465.         return angular.applyTo(new TimeStampedFieldPVCoordinates<>(pv.getDate(),
  466.                                                                    pv.getPosition().add(cartesian.getPosition()),
  467.                                                                    pv.getVelocity().add(cartesian.getVelocity()),
  468.                                                                    pv.getAcceleration().add(cartesian.getAcceleration())));
  469.     }

  470.     /** Compute the Jacobian of the {@link #transformPVCoordinates(PVCoordinates)}
  471.      * method of the transform.
  472.      * <p>
  473.      * Element {@code jacobian[i][j]} is the derivative of Cartesian coordinate i
  474.      * of the transformed {@link PVCoordinates} with respect to Cartesian coordinate j
  475.      * of the input {@link PVCoordinates} in method {@link #transformPVCoordinates(PVCoordinates)}.
  476.      * </p>
  477.      * <p>
  478.      * This definition implies that if we define position-velocity coordinates
  479.      * <pre>
  480.      * PV₁ = transform.transformPVCoordinates(PV₀), then
  481.      * </pre>
  482.      * <p> their differentials dPV₁ and dPV₀ will obey the following relation
  483.      * where J is the matrix computed by this method:
  484.      * <pre>
  485.      * dPV₁ = J &times; dPV₀
  486.      * </pre>
  487.      *
  488.      * @param selector selector specifying the size of the upper left corner that must be filled
  489.      * (either 3x3 for positions only, 6x6 for positions and velocities, 9x9 for positions,
  490.      * velocities and accelerations)
  491.      * @param jacobian placeholder matrix whose upper-left corner is to be filled with
  492.      * the Jacobian, the rest of the matrix remaining untouched
  493.      */
  494.     public void getJacobian(final CartesianDerivativesFilter selector, final double[][] jacobian) {

  495.         // elementary matrix for rotation
  496.         final double[][] mData = angular.getRotation().getMatrix();

  497.         // dP1/dP0
  498.         System.arraycopy(mData[0], 0, jacobian[0], 0, 3);
  499.         System.arraycopy(mData[1], 0, jacobian[1], 0, 3);
  500.         System.arraycopy(mData[2], 0, jacobian[2], 0, 3);

  501.         if (selector.getMaxOrder() >= 1) {

  502.             // dP1/dV0
  503.             Arrays.fill(jacobian[0], 3, 6, 0.0);
  504.             Arrays.fill(jacobian[1], 3, 6, 0.0);
  505.             Arrays.fill(jacobian[2], 3, 6, 0.0);

  506.             // dV1/dP0
  507.             final Vector3D o = angular.getRotationRate();
  508.             final double ox = o.getX();
  509.             final double oy = o.getY();
  510.             final double oz = o.getZ();
  511.             for (int i = 0; i < 3; ++i) {
  512.                 jacobian[3][i] = -(oy * mData[2][i] - oz * mData[1][i]);
  513.                 jacobian[4][i] = -(oz * mData[0][i] - ox * mData[2][i]);
  514.                 jacobian[5][i] = -(ox * mData[1][i] - oy * mData[0][i]);
  515.             }

  516.             // dV1/dV0
  517.             System.arraycopy(mData[0], 0, jacobian[3], 3, 3);
  518.             System.arraycopy(mData[1], 0, jacobian[4], 3, 3);
  519.             System.arraycopy(mData[2], 0, jacobian[5], 3, 3);

  520.             if (selector.getMaxOrder() >= 2) {

  521.                 // dP1/dA0
  522.                 Arrays.fill(jacobian[0], 6, 9, 0.0);
  523.                 Arrays.fill(jacobian[1], 6, 9, 0.0);
  524.                 Arrays.fill(jacobian[2], 6, 9, 0.0);

  525.                 // dV1/dA0
  526.                 Arrays.fill(jacobian[3], 6, 9, 0.0);
  527.                 Arrays.fill(jacobian[4], 6, 9, 0.0);
  528.                 Arrays.fill(jacobian[5], 6, 9, 0.0);

  529.                 // dA1/dP0
  530.                 final Vector3D oDot = angular.getRotationAcceleration();
  531.                 final double oDotx  = oDot.getX();
  532.                 final double oDoty  = oDot.getY();
  533.                 final double oDotz  = oDot.getZ();
  534.                 for (int i = 0; i < 3; ++i) {
  535.                     jacobian[6][i] = -(oDoty * mData[2][i] - oDotz * mData[1][i]) - (oy * jacobian[5][i] - oz * jacobian[4][i]);
  536.                     jacobian[7][i] = -(oDotz * mData[0][i] - oDotx * mData[2][i]) - (oz * jacobian[3][i] - ox * jacobian[5][i]);
  537.                     jacobian[8][i] = -(oDotx * mData[1][i] - oDoty * mData[0][i]) - (ox * jacobian[4][i] - oy * jacobian[3][i]);
  538.                 }

  539.                 // dA1/dV0
  540.                 for (int i = 0; i < 3; ++i) {
  541.                     jacobian[6][i + 3] = -2 * (oy * mData[2][i] - oz * mData[1][i]);
  542.                     jacobian[7][i + 3] = -2 * (oz * mData[0][i] - ox * mData[2][i]);
  543.                     jacobian[8][i + 3] = -2 * (ox * mData[1][i] - oy * mData[0][i]);
  544.                 }

  545.                 // dA1/dA0
  546.                 System.arraycopy(mData[0], 0, jacobian[6], 6, 3);
  547.                 System.arraycopy(mData[1], 0, jacobian[7], 6, 3);
  548.                 System.arraycopy(mData[2], 0, jacobian[8], 6, 3);

  549.             }

  550.         }

  551.     }

  552.     /** Get the underlying elementary Cartesian part.
  553.      * <p>A transform can be uniquely represented as an elementary
  554.      * translation followed by an elementary rotation. This method
  555.      * returns this unique elementary translation with its derivative.</p>
  556.      * @return underlying elementary Cartesian part
  557.      * @see #getTranslation()
  558.      * @see #getVelocity()
  559.      */
  560.     public PVCoordinates getCartesian() {
  561.         return cartesian;
  562.     }

  563.     /** Get the underlying elementary translation.
  564.      * <p>A transform can be uniquely represented as an elementary
  565.      * translation followed by an elementary rotation. This method
  566.      * returns this unique elementary translation.</p>
  567.      * @return underlying elementary translation
  568.      * @see #getCartesian()
  569.      * @see #getVelocity()
  570.      * @see #getAcceleration()
  571.      */
  572.     public Vector3D getTranslation() {
  573.         return cartesian.getPosition();
  574.     }

  575.     /** Get the first time derivative of the translation.
  576.      * @return first time derivative of the translation
  577.      * @see #getCartesian()
  578.      * @see #getTranslation()
  579.      * @see #getAcceleration()
  580.      */
  581.     public Vector3D getVelocity() {
  582.         return cartesian.getVelocity();
  583.     }

  584.     /** Get the second time derivative of the translation.
  585.      * @return second time derivative of the translation
  586.      * @see #getCartesian()
  587.      * @see #getTranslation()
  588.      * @see #getVelocity()
  589.      */
  590.     public Vector3D getAcceleration() {
  591.         return cartesian.getAcceleration();
  592.     }

  593.     /** Get the underlying elementary angular part.
  594.      * <p>A transform can be uniquely represented as an elementary
  595.      * translation followed by an elementary rotation. This method
  596.      * returns this unique elementary rotation with its derivative.</p>
  597.      * @return underlying elementary angular part
  598.      * @see #getRotation()
  599.      * @see #getRotationRate()
  600.      * @see #getRotationAcceleration()
  601.      */
  602.     public AngularCoordinates getAngular() {
  603.         return angular;
  604.     }

  605.     /** Get the underlying elementary rotation.
  606.      * <p>A transform can be uniquely represented as an elementary
  607.      * translation followed by an elementary rotation. This method
  608.      * returns this unique elementary rotation.</p>
  609.      * @return underlying elementary rotation
  610.      * @see #getAngular()
  611.      * @see #getRotationRate()
  612.      * @see #getRotationAcceleration()
  613.      */
  614.     public Rotation getRotation() {
  615.         return angular.getRotation();
  616.     }

  617.     /** Get the first time derivative of the rotation.
  618.      * <p>The norm represents the angular rate.</p>
  619.      * @return First time derivative of the rotation
  620.      * @see #getAngular()
  621.      * @see #getRotation()
  622.      * @see #getRotationAcceleration()
  623.      */
  624.     public Vector3D getRotationRate() {
  625.         return angular.getRotationRate();
  626.     }

  627.     /** Get the second time derivative of the rotation.
  628.      * @return Second time derivative of the rotation
  629.      * @see #getAngular()
  630.      * @see #getRotation()
  631.      * @see #getRotationRate()
  632.      */
  633.     public Vector3D getRotationAcceleration() {
  634.         return angular.getRotationAcceleration();
  635.     }

  636.     /** Specialized class for identity transform. */
  637.     private static class IdentityTransform extends Transform {

  638.         /** Serializable UID. */
  639.         private static final long serialVersionUID = -9042082036141830517L;

  640.         /** Simple constructor. */
  641.         IdentityTransform() {
  642.             super(AbsoluteDate.ARBITRARY_EPOCH, PVCoordinates.ZERO, AngularCoordinates.IDENTITY);
  643.         }

  644.         /** {@inheritDoc} */
  645.         @Override
  646.         public Transform shiftedBy(final double dt) {
  647.             return this;
  648.         }

  649.         /** {@inheritDoc} */
  650.         @Override
  651.         public Transform getInverse() {
  652.             return this;
  653.         }

  654.         /** {@inheritDoc} */
  655.         @Override
  656.         public Vector3D transformPosition(final Vector3D position) {
  657.             return position;
  658.         }

  659.         /** {@inheritDoc} */
  660.         @Override
  661.         public Vector3D transformVector(final Vector3D vector) {
  662.             return vector;
  663.         }

  664.         /** {@inheritDoc} */
  665.         @Override
  666.         public Line transformLine(final Line line) {
  667.             return line;
  668.         }

  669.         /** {@inheritDoc} */
  670.         @Override
  671.         public PVCoordinates transformPVCoordinates(final PVCoordinates pv) {
  672.             return pv;
  673.         }

  674.         @Override
  675.         public Transform freeze() {
  676.             return this;
  677.         }

  678.         @Override
  679.         public TimeStampedPVCoordinates transformPVCoordinates(
  680.                 final TimeStampedPVCoordinates pv) {
  681.             return pv;
  682.         }

  683.         @Override
  684.         public <T extends CalculusFieldElement<T>> FieldPVCoordinates<T>
  685.             transformPVCoordinates(final FieldPVCoordinates<T> pv) {
  686.             return pv;
  687.         }

  688.         @Override
  689.         public <T extends CalculusFieldElement<T>>
  690.             TimeStampedFieldPVCoordinates<T> transformPVCoordinates(
  691.                     final TimeStampedFieldPVCoordinates<T> pv) {
  692.             return pv;
  693.         }

  694.         /** {@inheritDoc} */
  695.         @Override
  696.         public void getJacobian(final CartesianDerivativesFilter selector, final double[][] jacobian) {
  697.             final int n = 3 * (selector.getMaxOrder() + 1);
  698.             for (int i = 0; i < n; ++i) {
  699.                 Arrays.fill(jacobian[i], 0, n, 0.0);
  700.                 jacobian[i][i] = 1.0;
  701.             }
  702.         }

  703.     }

  704. }