FieldAbsolutePVCoordinatesHermiteInterpolator.java
/* Copyright 2002-2024 CS GROUP
* Licensed to CS GROUP (CS) under one or more
* contributor license agreements. See the NOTICE file distributed with
* this work for additional information regarding copyright ownership.
* CS licenses this file to You under the Apache License, Version 2.0
* (the "License"); you may not use this file except in compliance with
* the License. You may obtain a copy of the License at
*
* http://www.apache.org/licenses/LICENSE-2.0
*
* Unless required by applicable law or agreed to in writing, software
* distributed under the License is distributed on an "AS IS" BASIS,
* WITHOUT WARRANTIES OR CONDITIONS OF ANY KIND, either express or implied.
* See the License for the specific language governing permissions and
* limitations under the License.
*/
package org.orekit.utils;
import org.hipparchus.CalculusFieldElement;
import org.hipparchus.analysis.interpolation.FieldHermiteInterpolator;
import org.hipparchus.geometry.euclidean.threed.FieldVector3D;
import org.orekit.errors.OrekitInternalError;
import org.orekit.frames.Frame;
import org.orekit.time.AbstractFieldTimeInterpolator;
import org.orekit.time.FieldAbsoluteDate;
import java.util.List;
/**
* Class using a Hermite interpolator to interpolate absolute position-velocity-acceleration coordinates.
* <p>
* As this implementation of interpolation is polynomial, it should be used only with small number of interpolation points
* (about 10-20 points) in order to avoid <a href="http://en.wikipedia.org/wiki/Runge%27s_phenomenon">Runge's phenomenon</a>
* and numerical problems (including NaN appearing).
*
* @author Luc Maisonobe
* @author Vincent Cucchietti
* @see FieldHermiteInterpolator
* @see FieldAbsolutePVCoordinates
* @param <KK> type of the field elements
*/
public class FieldAbsolutePVCoordinatesHermiteInterpolator<KK extends CalculusFieldElement<KK>>
extends AbstractFieldTimeInterpolator<FieldAbsolutePVCoordinates<KK>, KK> {
/** Filter for derivatives from the sample to use in interpolation. */
private final CartesianDerivativesFilter filter;
/** Output frame for the interpolated instance. */
private final Frame outputFrame;
/**
* Constructor with :
* <ul>
* <li>Default number of interpolation points of {@code DEFAULT_INTERPOLATION_POINTS}</li>
* <li>Default extrapolation threshold value ({@code DEFAULT_EXTRAPOLATION_THRESHOLD_SEC} s)</li>
* <li>Use of position and two time derivatives during interpolation</li>
* </ul>
* As this implementation of interpolation is polynomial, it should be used only with small number of interpolation
* points (about 10-20 points) in order to avoid <a href="http://en.wikipedia.org/wiki/Runge%27s_phenomenon">Runge's
* phenomenon</a> and numerical problems (including NaN appearing).
*
* @param outputFrame frame for the interpolated instance
*/
public FieldAbsolutePVCoordinatesHermiteInterpolator(final Frame outputFrame) {
this(DEFAULT_INTERPOLATION_POINTS, outputFrame);
}
/**
* Constructor with :
* <ul>
* <li>Default extrapolation threshold value ({@code DEFAULT_EXTRAPOLATION_THRESHOLD_SEC} s)</li>
* <li>Use of position and two time derivatives during interpolation</li>
* </ul>
* As this implementation of interpolation is polynomial, it should be used only with small number of interpolation
* points (about 10-20 points) in order to avoid <a href="http://en.wikipedia.org/wiki/Runge%27s_phenomenon">Runge's
* phenomenon</a> and numerical problems (including NaN appearing).
*
* @param interpolationPoints number of interpolation points
* @param outputFrame frame for the interpolated instance
*/
public FieldAbsolutePVCoordinatesHermiteInterpolator(final int interpolationPoints, final Frame outputFrame) {
this(interpolationPoints, outputFrame, CartesianDerivativesFilter.USE_PVA);
}
/**
* Constructor with default extrapolation threshold value ({@code DEFAULT_EXTRAPOLATION_THRESHOLD_SEC} s).
* <p>
* As this implementation of interpolation is polynomial, it should be used only with small number of interpolation
* points (about 10-20 points) in order to avoid <a href="http://en.wikipedia.org/wiki/Runge%27s_phenomenon">Runge's
* phenomenon</a> and numerical problems (including NaN appearing).
*
* @param interpolationPoints number of interpolation points
* @param outputFrame frame for the interpolated instance
* @param filter filter for derivatives from the sample to use in interpolation
*/
public FieldAbsolutePVCoordinatesHermiteInterpolator(final int interpolationPoints, final Frame outputFrame,
final CartesianDerivativesFilter filter) {
super(interpolationPoints, DEFAULT_EXTRAPOLATION_THRESHOLD_SEC);
this.outputFrame = outputFrame;
this.filter = filter;
}
/**
* Constructor.
* <p>
* As this implementation of interpolation is polynomial, it should be used only with small number of interpolation
* points (about 10-20 points) in order to avoid <a href="http://en.wikipedia.org/wiki/Runge%27s_phenomenon">Runge's
* phenomenon</a> and numerical problems (including NaN appearing).
*
* @param interpolationPoints number of interpolation points
* @param extrapolationThreshold extrapolation threshold beyond which the propagation will fail
* @param outputFrame frame for the interpolated instance
* @param filter filter for derivatives from the sample to use in interpolation
*/
public FieldAbsolutePVCoordinatesHermiteInterpolator(final int interpolationPoints, final double extrapolationThreshold,
final Frame outputFrame, final CartesianDerivativesFilter filter) {
super(interpolationPoints, extrapolationThreshold);
this.outputFrame = outputFrame;
this.filter = filter;
}
/** Get filter for derivatives from the sample to use in interpolation.
* @return filter for derivatives from the sample to use in interpolation
*/
public CartesianDerivativesFilter getFilter() {
return filter;
}
/** Get output frame for the interpolated instance.
* @return output frame for the interpolated instance
*/
public Frame getOutputFrame() {
return outputFrame;
}
/**
* {@inheritDoc}
* <p>
* The interpolated instance is created by polynomial Hermite interpolation ensuring velocity remains the exact
* derivative of position.
* <p>
* Note that even if first time derivatives (velocities) from sample can be ignored, the interpolated instance always
* includes interpolated derivatives. This feature can be used explicitly to compute these derivatives when it would be
* too complex to compute them from an analytical formula: just compute a few sample points from the explicit formula and
* set the derivatives to zero in these sample points, then use interpolation to add derivatives consistent with the
* positions.
*/
@Override
protected FieldAbsolutePVCoordinates<KK> interpolate(final InterpolationData interpolationData) {
// Get interpolation date
final FieldAbsoluteDate<KK> date = interpolationData.getInterpolationDate();
// Get sample
final List<FieldAbsolutePVCoordinates<KK>> sample = interpolationData.getNeighborList();
// Set up an interpolator taking derivatives into account
final FieldHermiteInterpolator<KK> interpolator = new FieldHermiteInterpolator<>();
// Add sample points
switch (filter) {
case USE_P:
// Populate sample with position data, ignoring velocity
sample.forEach(pv -> {
final FieldVector3D<KK> position = pv.getPosition();
interpolator.addSamplePoint(pv.getDate().durationFrom(date),
position.toArray());
});
break;
case USE_PV:
// Populate sample with position and velocity data
sample.forEach(pv -> {
final FieldVector3D<KK> position = pv.getPosition();
final FieldVector3D<KK> velocity = pv.getVelocity();
interpolator.addSamplePoint(pv.getDate().durationFrom(date),
position.toArray(), velocity.toArray());
});
break;
case USE_PVA:
// Populate sample with position, velocity and acceleration data
sample.forEach(pv -> {
final FieldVector3D<KK> position = pv.getPosition();
final FieldVector3D<KK> velocity = pv.getVelocity();
final FieldVector3D<KK> acceleration = pv.getAcceleration();
interpolator.addSamplePoint(pv.getDate().durationFrom(date),
position.toArray(), velocity.toArray(), acceleration.toArray());
});
break;
default:
// this should never happen
throw new OrekitInternalError(null);
}
// interpolate
final KK zero = interpolationData.getZero();
final KK[][] p = interpolator.derivatives(zero, 2);
// build a new interpolated instance
return new FieldAbsolutePVCoordinates<>(outputFrame, date, new FieldVector3D<>(p[0]), new FieldVector3D<>(p[1]),
new FieldVector3D<>(p[2]));
}
}