use of org.orekit.propagation.numerical.NumericalPropagator in project Orekit by CS-SI.
the class AttitudesSequenceTest method testResetDuringTransitionForward.
@Test
public void testResetDuringTransitionForward() throws OrekitException {
// Initial state definition : date, orbit
final AbsoluteDate initialDate = new AbsoluteDate(2004, 01, 01, 23, 30, 00.000, TimeScalesFactory.getUTC());
final Vector3D position = new Vector3D(-6142438.668, 3492467.560, -25767.25680);
final Vector3D velocity = new Vector3D(505.8479685, 942.7809215, 7435.922231);
final Orbit initialOrbit = new KeplerianOrbit(new PVCoordinates(position, velocity), FramesFactory.getEME2000(), initialDate, Constants.EIGEN5C_EARTH_MU);
final OneAxisEllipsoid earth = new OneAxisEllipsoid(Constants.WGS84_EARTH_EQUATORIAL_RADIUS, Constants.WGS84_EARTH_FLATTENING, FramesFactory.getITRF(IERSConventions.IERS_2010, true));
final TopocentricFrame volgograd = new TopocentricFrame(earth, new GeodeticPoint(FastMath.toRadians(48.7), FastMath.toRadians(44.5), 24.0), "Волгоград");
final AttitudesSequence attitudesSequence = new AttitudesSequence();
final double transitionTime = 250.0;
final AttitudeProvider nadirPointing = new NadirPointing(initialOrbit.getFrame(), earth);
final AttitudeProvider targetPointing = new TargetPointing(initialOrbit.getFrame(), volgograd.getPoint(), earth);
final ElevationDetector eventDetector = new ElevationDetector(volgograd).withConstantElevation(FastMath.toRadians(5.0)).withHandler(new ContinueOnEvent<>());
final List<AbsoluteDate> nadirToTarget = new ArrayList<>();
attitudesSequence.addSwitchingCondition(nadirPointing, targetPointing, eventDetector, true, false, transitionTime, AngularDerivativesFilter.USE_RR, (previous, next, state) -> nadirToTarget.add(state.getDate()));
final double[][] tolerance = NumericalPropagator.tolerances(10.0, initialOrbit, initialOrbit.getType());
final AdaptiveStepsizeIntegrator integrator = new DormandPrince853Integrator(0.001, 300.0, tolerance[0], tolerance[1]);
final NumericalPropagator propagator = new NumericalPropagator(integrator);
GravityFieldFactory.addPotentialCoefficientsReader(new ICGEMFormatReader("g007_eigen_05c_coef", false));
propagator.addForceModel(new HolmesFeatherstoneAttractionModel(earth.getBodyFrame(), GravityFieldFactory.getNormalizedProvider(8, 8)));
propagator.setInitialState(new SpacecraftState(initialOrbit, nadirPointing.getAttitude(initialOrbit, initialOrbit.getDate(), initialOrbit.getFrame())));
propagator.setAttitudeProvider(attitudesSequence);
attitudesSequence.registerSwitchEvents(propagator);
propagator.propagate(initialDate.shiftedBy(6000));
// check that if we restart a forward propagation from an intermediate state
// we properly get an interpolated attitude despite we missed the event trigger
final AbsoluteDate midTransition = nadirToTarget.get(0).shiftedBy(0.5 * transitionTime);
SpacecraftState state = propagator.propagate(midTransition.shiftedBy(-60), midTransition);
Rotation nadirR = nadirPointing.getAttitude(state.getOrbit(), state.getDate(), state.getFrame()).getRotation();
Rotation targetR = targetPointing.getAttitude(state.getOrbit(), state.getDate(), state.getFrame()).getRotation();
final double reorientationAngle = Rotation.distance(nadirR, targetR);
Assert.assertEquals(0.5 * reorientationAngle, Rotation.distance(state.getAttitude().getRotation(), nadirR), 0.03 * reorientationAngle);
}
use of org.orekit.propagation.numerical.NumericalPropagator in project Orekit by CS-SI.
the class AttitudesSequenceTest method testOutOfSyncCalls.
@Test
public void testOutOfSyncCalls() throws OrekitException {
// Initial state definition : date, orbit
final AbsoluteDate initialDate = new AbsoluteDate(2004, 01, 01, 23, 30, 00.000, TimeScalesFactory.getUTC());
final Vector3D position = new Vector3D(-6142438.668, 3492467.560, -25767.25680);
final Vector3D velocity = new Vector3D(505.8479685, 942.7809215, 7435.922231);
final Orbit initialOrbit = new KeplerianOrbit(new PVCoordinates(position, velocity), FramesFactory.getEME2000(), initialDate, Constants.EIGEN5C_EARTH_MU);
final OneAxisEllipsoid earth = new OneAxisEllipsoid(Constants.WGS84_EARTH_EQUATORIAL_RADIUS, Constants.WGS84_EARTH_FLATTENING, FramesFactory.getITRF(IERSConventions.IERS_2010, true));
final TopocentricFrame volgograd = new TopocentricFrame(earth, new GeodeticPoint(FastMath.toRadians(48.7), FastMath.toRadians(44.5), 24.0), "Волгоград");
final AttitudesSequence attitudesSequence = new AttitudesSequence();
final double transitionTime = 250.0;
final AttitudeProvider nadirPointing = new NadirPointing(initialOrbit.getFrame(), earth);
final AttitudeProvider targetPointing = new TargetPointing(initialOrbit.getFrame(), volgograd.getPoint(), earth);
final ElevationDetector eventDetector = new ElevationDetector(volgograd).withConstantElevation(FastMath.toRadians(5.0)).withHandler(new ContinueOnEvent<>());
final Handler nadirToTarget = new Handler(nadirPointing, targetPointing);
attitudesSequence.addSwitchingCondition(nadirPointing, targetPointing, eventDetector, true, false, transitionTime, AngularDerivativesFilter.USE_RR, nadirToTarget);
final Handler targetToNadir = new Handler(targetPointing, nadirPointing);
attitudesSequence.addSwitchingCondition(targetPointing, nadirPointing, eventDetector, false, true, transitionTime, AngularDerivativesFilter.USE_RR, targetToNadir);
final double[][] tolerance = NumericalPropagator.tolerances(10.0, initialOrbit, initialOrbit.getType());
final AdaptiveStepsizeIntegrator integrator = new DormandPrince853Integrator(0.001, 300.0, tolerance[0], tolerance[1]);
final NumericalPropagator propagator = new NumericalPropagator(integrator);
GravityFieldFactory.addPotentialCoefficientsReader(new ICGEMFormatReader("g007_eigen_05c_coef", false));
propagator.addForceModel(new HolmesFeatherstoneAttractionModel(earth.getBodyFrame(), GravityFieldFactory.getNormalizedProvider(8, 8)));
propagator.setInitialState(new SpacecraftState(initialOrbit, nadirPointing.getAttitude(initialOrbit, initialOrbit.getDate(), initialOrbit.getFrame())));
propagator.setAttitudeProvider(attitudesSequence);
attitudesSequence.registerSwitchEvents(propagator);
propagator.setMasterMode(10, (state, isLast) -> {
Attitude nadirAttitude = nadirPointing.getAttitude(state.getOrbit(), state.getDate(), state.getFrame());
Attitude targetAttitude = targetPointing.getAttitude(state.getOrbit(), state.getDate(), state.getFrame());
Attitude stateAttitude = state.getAttitude();
if (nadirToTarget.dates.isEmpty() || state.getDate().durationFrom(nadirToTarget.dates.get(0)) < 0) {
// we are stabilized in nadir pointing, before first switch
checkEqualAttitudes(nadirAttitude, stateAttitude);
} else if (state.getDate().durationFrom(nadirToTarget.dates.get(0)) <= transitionTime) {
// we are in transition from nadir to target
checkBetweenAttitudes(nadirAttitude, targetAttitude, stateAttitude);
} else if (targetToNadir.dates.isEmpty() || state.getDate().durationFrom(targetToNadir.dates.get(0)) < 0) {
// we are stabilized in target pointing between the two switches
checkEqualAttitudes(targetAttitude, stateAttitude);
} else if (state.getDate().durationFrom(targetToNadir.dates.get(0)) <= transitionTime) {
// we are in transition from target to nadir
checkBetweenAttitudes(targetAttitude, nadirAttitude, stateAttitude);
} else {
// we are stabilized back in nadir pointing, after second switch
checkEqualAttitudes(nadirAttitude, stateAttitude);
}
});
propagator.propagate(initialDate.shiftedBy(6000));
}
use of org.orekit.propagation.numerical.NumericalPropagator in project SpriteOrbits by ProjectPersephone.
the class SpriteProp method createPropagator.
/**
* Create a numerical propagator for a state.
* @param state state to propagate
* @param attitudeProvider provider for the attitude
* @param crossSection cross section of the object
* @param dragCoeff drag coefficient
*/
private Propagator createPropagator(final SpacecraftState state, final AttitudeProvider attitudeProvider, final double crossSection, final double dragCoeff) throws OrekitException {
// see https://www.orekit.org/static/architecture/propagation.html
// steps limits
final double minStep = 0.001;
final double maxStep = 1000;
final double initStep = 60;
// error control parameters (absolute and relative)
final double positionError = 10.0;
// we will propagate in Cartesian coordinates
final OrbitType orbitType = OrbitType.CARTESIAN;
final double[][] tolerances = NumericalPropagator.tolerances(positionError, state.getOrbit(), orbitType);
// set up mathematical integrator
AdaptiveStepsizeIntegrator integrator = new DormandPrince853Integrator(minStep, maxStep, tolerances[0], tolerances[1]);
integrator.setInitialStepSize(initStep);
// set up space dynamics propagator
NumericalPropagator propagator = new NumericalPropagator(integrator);
propagator.setOrbitType(orbitType);
// add gravity field force model
final NormalizedSphericalHarmonicsProvider gravityProvider = GravityFieldFactory.getNormalizedProvider(8, 8);
propagator.addForceModel(new HolmesFeatherstoneAttractionModel(earth.getBodyFrame(), gravityProvider));
// add atmospheric drag force model
propagator.addForceModel(new DragForce(new HarrisPriester(sun, earth), new SphericalSpacecraft(crossSection, dragCoeff, 0.0, 0.0)));
// set attitude mode
propagator.setAttitudeProvider(attitudeProvider);
propagator.setInitialState(state);
return propagator;
}
use of org.orekit.propagation.numerical.NumericalPropagator in project Orekit by CS-SI.
the class KalmanEstimatorTest method testCartesianRangeRate.
/**
* Perfect range rate measurements with a perfect start
* Cartesian formalism
* @throws OrekitException
*/
@Test
public void testCartesianRangeRate() throws OrekitException {
// Create context
Context context = EstimationTestUtils.eccentricContext("regular-data:potential:tides");
// Create initial orbit and propagator builder
final OrbitType orbitType = OrbitType.CARTESIAN;
final PositionAngle positionAngle = PositionAngle.TRUE;
final boolean perfectStart = true;
final double minStep = 1.e-6;
final double maxStep = 60.;
final double dP = 1.;
final NumericalPropagatorBuilder propagatorBuilder = context.createBuilder(orbitType, positionAngle, perfectStart, minStep, maxStep, dP);
// Create perfect range measurements
final Propagator propagator = EstimationTestUtils.createPropagator(context.initialOrbit, propagatorBuilder);
final List<ObservedMeasurement<?>> measurements = EstimationTestUtils.createMeasurements(propagator, new RangeRateMeasurementCreator(context, false), 1.0, 3.0, 300.0);
// Reference propagator for estimation performances
final NumericalPropagator referencePropagator = propagatorBuilder.buildPropagator(propagatorBuilder.getSelectedNormalizedParameters());
// Reference position/velocity at last measurement date
final Orbit refOrbit = referencePropagator.propagate(measurements.get(measurements.size() - 1).getDate()).getOrbit();
// Cartesian covariance matrix initialization
// 100m on position / 1e-2m/s on velocity
final RealMatrix cartesianP = MatrixUtils.createRealDiagonalMatrix(new double[] { 1e-4, 1e-4, 1e-4, 1e-10, 1e-10, 1e-10 });
// Jacobian of the orbital parameters w/r to Cartesian
final Orbit initialOrbit = orbitType.convertType(context.initialOrbit);
final double[][] dYdC = new double[6][6];
initialOrbit.getJacobianWrtCartesian(PositionAngle.TRUE, dYdC);
final RealMatrix Jac = MatrixUtils.createRealMatrix(dYdC);
// Initial covariance matrix
final RealMatrix initialP = Jac.multiply(cartesianP.multiply(Jac.transpose()));
// Process noise matrix
final RealMatrix cartesianQ = MatrixUtils.createRealDiagonalMatrix(new double[] { 1.e-6, 1.e-6, 1.e-6, 1.e-12, 1.e-12, 1.e-12 });
final RealMatrix Q = Jac.multiply(cartesianQ.multiply(Jac.transpose()));
// Build the Kalman filter
final KalmanEstimatorBuilder kalmanBuilder = new KalmanEstimatorBuilder();
kalmanBuilder.builder(propagatorBuilder);
kalmanBuilder.estimatedMeasurementsParameters(new ParameterDriversList());
kalmanBuilder.initialCovarianceMatrix(initialP);
kalmanBuilder.processNoiseMatrixProvider(new ConstantProcessNoise(Q));
final KalmanEstimator kalman = kalmanBuilder.build();
// Filter the measurements and check the results
final double expectedDeltaPos = 0.;
final double posEps = 9.50e-4;
final double expectedDeltaVel = 0.;
final double velEps = 3.49e-7;
final double[] expectedSigmasPos = { 0.324398, 1.347031, 1.743310 };
final double sigmaPosEps = 1e-6;
final double[] expectedSigmasVel = { 2.856883e-4, 5.765844e-4, 5.056186e-4 };
final double sigmaVelEps = 1e-10;
EstimationTestUtils.checkKalmanFit(context, kalman, measurements, refOrbit, positionAngle, expectedDeltaPos, posEps, expectedDeltaVel, velEps, expectedSigmasPos, sigmaPosEps, expectedSigmasVel, sigmaVelEps);
}
use of org.orekit.propagation.numerical.NumericalPropagator in project Orekit by CS-SI.
the class KalmanEstimatorTest method testCircularAzimuthElevation.
/**
* Perfect azimuth/elevation measurements with a perfect start
* Circular formalism
* @throws OrekitException
*/
@Test
public void testCircularAzimuthElevation() throws OrekitException {
// Create context
Context context = EstimationTestUtils.eccentricContext("regular-data:potential:tides");
// Create initial orbit and propagator builder
final OrbitType orbitType = OrbitType.CIRCULAR;
final PositionAngle positionAngle = PositionAngle.TRUE;
final boolean perfectStart = true;
final double minStep = 1.e-6;
final double maxStep = 60.;
final double dP = 1.;
final NumericalPropagatorBuilder propagatorBuilder = context.createBuilder(orbitType, positionAngle, perfectStart, minStep, maxStep, dP);
// Create perfect range measurements
final Propagator propagator = EstimationTestUtils.createPropagator(context.initialOrbit, propagatorBuilder);
final List<ObservedMeasurement<?>> measurements = EstimationTestUtils.createMeasurements(propagator, new AngularAzElMeasurementCreator(context), 1.0, 4.0, 60.0);
// Reference propagator for estimation performances
final NumericalPropagator referencePropagator = propagatorBuilder.buildPropagator(propagatorBuilder.getSelectedNormalizedParameters());
// Reference position/velocity at last measurement date
final Orbit refOrbit = referencePropagator.propagate(measurements.get(measurements.size() - 1).getDate()).getOrbit();
// Cartesian covariance matrix initialization
final RealMatrix cartesianP = MatrixUtils.createRealDiagonalMatrix(new double[] { 1e-4, 1e-4, 1e-4, 1e-10, 1e-10, 1e-10 });
// Jacobian of the orbital parameters w/r to Cartesian
final Orbit initialOrbit = orbitType.convertType(context.initialOrbit);
final double[][] dYdC = new double[6][6];
initialOrbit.getJacobianWrtCartesian(PositionAngle.TRUE, dYdC);
final RealMatrix Jac = MatrixUtils.createRealMatrix(dYdC);
// Initial covariance matrix
final RealMatrix initialP = Jac.multiply(cartesianP.multiply(Jac.transpose()));
// Process noise matrix
final RealMatrix cartesianQ = MatrixUtils.createRealDiagonalMatrix(new double[] { 1.e-6, 1.e-6, 1.e-6, 1.e-12, 1.e-12, 1.e-12 });
final RealMatrix Q = Jac.multiply(cartesianQ.multiply(Jac.transpose()));
// Build the Kalman filter
final KalmanEstimatorBuilder kalmanBuilder = new KalmanEstimatorBuilder();
kalmanBuilder.builder(propagatorBuilder);
kalmanBuilder.estimatedMeasurementsParameters(new ParameterDriversList());
kalmanBuilder.initialCovarianceMatrix(initialP);
kalmanBuilder.processNoiseMatrixProvider(new ConstantProcessNoise(Q));
final KalmanEstimator kalman = kalmanBuilder.build();
// Filter the measurements and check the results
final double expectedDeltaPos = 0.;
final double posEps = 4.78e-7;
final double expectedDeltaVel = 0.;
final double velEps = 1.54e-10;
final double[] expectedSigmasPos = { 0.356902, 1.297507, 1.798551 };
final double sigmaPosEps = 1e-6;
final double[] expectedSigmasVel = { 2.468745e-4, 5.810027e-4, 3.887394e-4 };
final double sigmaVelEps = 1e-10;
EstimationTestUtils.checkKalmanFit(context, kalman, measurements, refOrbit, positionAngle, expectedDeltaPos, posEps, expectedDeltaVel, velEps, expectedSigmasPos, sigmaPosEps, expectedSigmasVel, sigmaVelEps);
}
Aggregations