LibrationOrbit.java
/* Copyright 2002-2021 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.orbits;
import org.hipparchus.geometry.euclidean.threed.Vector3D;
import org.hipparchus.linear.EigenDecomposition;
import org.hipparchus.linear.RealMatrix;
import org.hipparchus.linear.RealVector;
import org.orekit.annotation.DefaultDataContext;
import org.orekit.attitudes.AttitudeProvider;
import org.orekit.attitudes.InertialProvider;
import org.orekit.bodies.CR3BPSystem;
import org.orekit.data.DataContext;
import org.orekit.propagation.SpacecraftState;
import org.orekit.propagation.numerical.cr3bp.STMEquations;
import org.orekit.time.TimeScale;
import org.orekit.utils.PVCoordinates;
/**
* Base class for libration orbits.
* @see HaloOrbit
* @see LyapunovOrbit
* @author Vincent Mouraux
* @author Bryan Cazabonne
* @since 10.2
*/
public abstract class LibrationOrbit {
/** CR3BP System of the libration Orbit. */
private final CR3BPSystem syst;
/** Position-Velocity initial position on a libration Orbit. */
private PVCoordinates initialPV;
/** Orbital Period of the libration Orbit. */
private double orbitalPeriod;
/**
* Constructor.
* @param system CR3BP System considered
* @param initialPV initial position on a libration Orbit
* @param orbitalPeriod initial orbital period of the libration Orbit
*/
protected LibrationOrbit(final CR3BPSystem system,
final PVCoordinates initialPV,
final double orbitalPeriod) {
this.syst = system;
this.initialPV = initialPV;
this.orbitalPeriod = orbitalPeriod;
}
/** Return the orbital period of the libration orbit.
* @return orbitalPeriod orbital period of the libration orbit
*/
public double getOrbitalPeriod() {
return orbitalPeriod;
}
/** Return the initialPV on the libration orbit.
* <p>
* This will return the exact initialPV only if you applied a prior
* differential correction. If you did not, you can use the method
* {@link #applyCorrectionOnPV(CR3BPDifferentialCorrection)}
* </p>
* @return initialPV initialPV on the libration orbit
*/
public PVCoordinates getInitialPV() {
return initialPV;
}
/** Apply differential correction.
* <p>
* This will update {@link #initialPV} and
* {@link #orbitalPeriod} parameters.
* </p>
*/
@DefaultDataContext
public void applyDifferentialCorrection() {
applyDifferentialCorrection(InertialProvider.of(syst.getRotatingFrame()),
DataContext.getDefault().getTimeScales().getUTC());
}
/** Apply differential correction.
* <p>
* This will update {@link #initialPV} and
* {@link #orbitalPeriod} parameters.
* </p>
* @param attitudeProvider the attitude law for the numerocal propagator
* @param utc UTC time scale
*/
public void applyDifferentialCorrection(final AttitudeProvider attitudeProvider,
final TimeScale utc) {
final CR3BPDifferentialCorrection diff = new CR3BPDifferentialCorrection(initialPV, syst, orbitalPeriod, attitudeProvider, utc);
initialPV = applyCorrectionOnPV(diff);
orbitalPeriod = diff.getOrbitalPeriod();
}
/** Return a manifold direction from one position on a libration Orbit.
* @param s SpacecraftState with additional equations
* @param isStable true if the manifold is stable
* @return manifold first guess Position-Velocity of a point on the libration Orbit
*/
public PVCoordinates getManifolds(final SpacecraftState s,
final boolean isStable) {
return isStable ? getStableManifolds(s) : getUnstableManifolds(s);
}
/** Return the stable manifold direction for several positions on a libration Orbit.
* @param s SpacecraftStates (with STM equations) to compute from
* @return Stable manifold first direction from a point on the libration Orbit
*/
private PVCoordinates getStableManifolds(final SpacecraftState s) {
// Small delta, linked to the characteristic velocity of the CR3BP system
final double epsilon = syst.getVdim() * 1E2 / syst.getDdim();
// Get Normalize eigen vector linked to the stability of the manifold
final RealMatrix phi = new STMEquations(syst).getStateTransitionMatrix(s);
final RealVector eigenVector = new EigenDecomposition(phi).getEigenvector(1).unitVector();
// New PVCoordinates following the manifold
return new PVCoordinates(s.getPVCoordinates().getPosition()
.add(new Vector3D(eigenVector.getEntry(0), eigenVector
.getEntry(1), eigenVector.getEntry(2))
.scalarMultiply(epsilon)), s.getPVCoordinates()
.getVelocity()
.add(new Vector3D(eigenVector.getEntry(3),
eigenVector.getEntry(4),
eigenVector.getEntry(5))
.scalarMultiply(epsilon)));
}
/** Get the Unstable manifold direction for several positions on a libration Orbit.
* @param s spacecraft state (with STM equations) to compute from
* @return pv coordinates representing the unstable manifold first direction
* from a point on the libration Orbit
*/
private PVCoordinates getUnstableManifolds(final SpacecraftState s) {
// Small delta, linked to the characteristic velocity of the CR3BP system
final double epsilon =
syst.getVdim() * 1E2 / syst.getDdim();
// Get Normalize eigen vector linked to the stability of the manifold
final RealMatrix phi = new STMEquations(syst).getStateTransitionMatrix(s);
final RealVector eigenVector = new EigenDecomposition(phi).getEigenvector(0).unitVector();
// New PVCoordinates following the manifold
return new PVCoordinates(s.getPVCoordinates().getPosition()
.add(new Vector3D(eigenVector.getEntry(0), eigenVector
.getEntry(1), eigenVector.getEntry(2))
.scalarMultiply(epsilon)), s.getPVCoordinates()
.getVelocity()
.add(new Vector3D(eigenVector.getEntry(3),
eigenVector.getEntry(4),
eigenVector.getEntry(5))
.scalarMultiply(epsilon)));
}
/**
* Apply the differential correction to compute more accurate initial PV.
* @param diff cr3bp differential correction
* @return corrected PV coordinates
*/
protected abstract PVCoordinates applyCorrectionOnPV(CR3BPDifferentialCorrection diff);
}