use of org.orekit.bodies.OneAxisEllipsoid in project Orekit by CS-SI.
the class OrbitDeterminationTest method createBody.
/**
* Create an orbit from input parameters
* @param parser input file parser
* @param mu central attraction coefficient
* @throws NoSuchElementException if input parameters are missing
* @throws OrekitException if body frame cannot be created
*/
private OneAxisEllipsoid createBody(final KeyValueFileParser<ParameterKey> parser) throws NoSuchElementException, OrekitException {
final Frame bodyFrame;
if (!parser.containsKey(ParameterKey.BODY_FRAME)) {
bodyFrame = FramesFactory.getITRF(IERSConventions.IERS_2010, true);
} else {
bodyFrame = parser.getEarthFrame(ParameterKey.BODY_FRAME);
}
final double equatorialRadius;
if (!parser.containsKey(ParameterKey.BODY_EQUATORIAL_RADIUS)) {
equatorialRadius = Constants.WGS84_EARTH_EQUATORIAL_RADIUS;
} else {
equatorialRadius = parser.getDouble(ParameterKey.BODY_EQUATORIAL_RADIUS);
}
final double flattening;
if (!parser.containsKey(ParameterKey.BODY_INVERSE_FLATTENING)) {
flattening = Constants.WGS84_EARTH_FLATTENING;
} else {
flattening = 1.0 / parser.getDouble(ParameterKey.BODY_INVERSE_FLATTENING);
}
return new OneAxisEllipsoid(equatorialRadius, flattening, bodyFrame);
}
use of org.orekit.bodies.OneAxisEllipsoid in project Orekit by CS-SI.
the class GroundStationTest method doTestCartesianDerivatives.
private void doTestCartesianDerivatives(double latitude, double longitude, double altitude, double stepFactor, double relativeTolerancePositionValue, double relativeTolerancePositionDerivative, double relativeToleranceVelocityValue, double relativeToleranceVelocityDerivative, String... parameterPattern) throws OrekitException {
Utils.setDataRoot("regular-data");
final Frame eme2000 = FramesFactory.getEME2000();
final AbsoluteDate date = AbsoluteDate.J2000_EPOCH;
final AbsoluteDate date0 = date.shiftedBy(50000);
final OneAxisEllipsoid earth = new OneAxisEllipsoid(Constants.WGS84_EARTH_EQUATORIAL_RADIUS, Constants.WGS84_EARTH_FLATTENING, FramesFactory.getITRF(IERSConventions.IERS_2010, true));
final GroundStation station = new GroundStation(new TopocentricFrame(earth, new GeodeticPoint(latitude, longitude, altitude), "dummy"));
final DSFactory factory = new DSFactory(parameterPattern.length, 1);
final FieldAbsoluteDate<DerivativeStructure> dateDS = new FieldAbsoluteDate<>(factory.getDerivativeField(), date);
ParameterDriver[] selectedDrivers = new ParameterDriver[parameterPattern.length];
UnivariateDifferentiableVectorFunction[] dFCartesian = new UnivariateDifferentiableVectorFunction[parameterPattern.length];
final ParameterDriver[] allDrivers = selectAllDrivers(station);
for (ParameterDriver driver : allDrivers) {
driver.setReferenceDate(date0);
}
Map<String, Integer> indices = new HashMap<>();
for (int k = 0; k < dFCartesian.length; ++k) {
for (int i = 0; i < allDrivers.length; ++i) {
if (allDrivers[i].getName().matches(parameterPattern[k])) {
selectedDrivers[k] = allDrivers[i];
dFCartesian[k] = differentiatedStationPV(station, eme2000, date, selectedDrivers[k], stepFactor);
indices.put(selectedDrivers[k].getName(), k);
}
}
}
;
DSFactory factory11 = new DSFactory(1, 1);
RandomGenerator generator = new Well19937a(0x084d58a19c498a54l);
double maxPositionValueRelativeError = 0;
double maxPositionDerivativeRelativeError = 0;
double maxVelocityValueRelativeError = 0;
double maxVelocityDerivativeRelativeError = 0;
for (int i = 0; i < 1000; ++i) {
// randomly change one parameter
ParameterDriver changed = allDrivers[generator.nextInt(allDrivers.length)];
changed.setNormalizedValue(2 * generator.nextDouble() - 1);
// transform to check
FieldTransform<DerivativeStructure> t = station.getOffsetToInertial(eme2000, dateDS, factory, indices);
FieldPVCoordinates<DerivativeStructure> pv = t.transformPVCoordinates(FieldPVCoordinates.getZero(factory.getDerivativeField()));
for (int k = 0; k < dFCartesian.length; ++k) {
// reference values and derivatives computed using finite differences
DerivativeStructure[] refCartesian = dFCartesian[k].value(factory11.variable(0, selectedDrivers[k].getValue()));
// position
final Vector3D refP = new Vector3D(refCartesian[0].getValue(), refCartesian[1].getValue(), refCartesian[2].getValue());
final Vector3D resP = new Vector3D(pv.getPosition().getX().getValue(), pv.getPosition().getY().getValue(), pv.getPosition().getZ().getValue());
maxPositionValueRelativeError = FastMath.max(maxPositionValueRelativeError, Vector3D.distance(refP, resP) / refP.getNorm());
final Vector3D refPD = new Vector3D(refCartesian[0].getPartialDerivative(1), refCartesian[1].getPartialDerivative(1), refCartesian[2].getPartialDerivative(1));
final Vector3D resPD = new Vector3D(pv.getPosition().getX().getAllDerivatives()[k + 1], pv.getPosition().getY().getAllDerivatives()[k + 1], pv.getPosition().getZ().getAllDerivatives()[k + 1]);
maxPositionDerivativeRelativeError = FastMath.max(maxPositionDerivativeRelativeError, Vector3D.distance(refPD, resPD) / refPD.getNorm());
// velocity
final Vector3D refV = new Vector3D(refCartesian[3].getValue(), refCartesian[4].getValue(), refCartesian[5].getValue());
final Vector3D resV = new Vector3D(pv.getVelocity().getX().getValue(), pv.getVelocity().getY().getValue(), pv.getVelocity().getZ().getValue());
maxVelocityValueRelativeError = FastMath.max(maxVelocityValueRelativeError, Vector3D.distance(refV, resV) / refV.getNorm());
final Vector3D refVD = new Vector3D(refCartesian[3].getPartialDerivative(1), refCartesian[4].getPartialDerivative(1), refCartesian[5].getPartialDerivative(1));
final Vector3D resVD = new Vector3D(pv.getVelocity().getX().getAllDerivatives()[k + 1], pv.getVelocity().getY().getAllDerivatives()[k + 1], pv.getVelocity().getZ().getAllDerivatives()[k + 1]);
maxVelocityDerivativeRelativeError = FastMath.max(maxVelocityDerivativeRelativeError, Vector3D.distance(refVD, resVD) / refVD.getNorm());
}
}
Assert.assertEquals(0.0, maxPositionValueRelativeError, relativeTolerancePositionValue);
Assert.assertEquals(0.0, maxPositionDerivativeRelativeError, relativeTolerancePositionDerivative);
Assert.assertEquals(0.0, maxVelocityValueRelativeError, relativeToleranceVelocityValue);
Assert.assertEquals(0.0, maxVelocityDerivativeRelativeError, relativeToleranceVelocityDerivative);
}
use of org.orekit.bodies.OneAxisEllipsoid in project Orekit by CS-SI.
the class GroundStationTest method testNoReferenceDate.
@Test
public void testNoReferenceDate() throws OrekitException {
Utils.setDataRoot("regular-data");
final Frame eme2000 = FramesFactory.getEME2000();
final AbsoluteDate date = AbsoluteDate.J2000_EPOCH;
final OneAxisEllipsoid earth = new OneAxisEllipsoid(Constants.WGS84_EARTH_EQUATORIAL_RADIUS, Constants.WGS84_EARTH_FLATTENING, FramesFactory.getITRF(IERSConventions.IERS_2010, true));
final GroundStation station = new GroundStation(new TopocentricFrame(earth, new GeodeticPoint(0.1, 0.2, 100), "dummy"));
try {
station.getOffsetToInertial(eme2000, date);
Assert.fail("an exception should have been thrown");
} catch (OrekitException oe) {
Assert.assertEquals(OrekitMessages.NO_REFERENCE_DATE_FOR_PARAMETER, oe.getSpecifier());
Assert.assertEquals("prime-meridian-offset", (String) oe.getParts()[0]);
}
try {
DSFactory factory = new DSFactory(9, 1);
Map<String, Integer> indices = new HashMap<>();
for (final ParameterDriver driver : Arrays.asList(station.getPrimeMeridianOffsetDriver(), station.getPrimeMeridianDriftDriver(), station.getPolarOffsetXDriver(), station.getPolarDriftXDriver(), station.getPolarOffsetYDriver(), station.getPolarDriftYDriver(), station.getEastOffsetDriver(), station.getNorthOffsetDriver(), station.getZenithOffsetDriver())) {
indices.put(driver.getName(), indices.size());
}
station.getOffsetToInertial(eme2000, new FieldAbsoluteDate<>(factory.getDerivativeField(), date), factory, indices);
Assert.fail("an exception should have been thrown");
} catch (OrekitException oe) {
Assert.assertEquals(OrekitMessages.NO_REFERENCE_DATE_FOR_PARAMETER, oe.getSpecifier());
Assert.assertEquals("prime-meridian-offset", (String) oe.getParts()[0]);
}
}
use of org.orekit.bodies.OneAxisEllipsoid in project Orekit by CS-SI.
the class GroundStationTest method doTestAngularDerivatives.
private void doTestAngularDerivatives(double latitude, double longitude, double altitude, double stepFactor, double toleranceRotationValue, double toleranceRotationDerivative, double toleranceRotationRateValue, double toleranceRotationRateDerivative, String... parameterPattern) throws OrekitException {
Utils.setDataRoot("regular-data");
final Frame eme2000 = FramesFactory.getEME2000();
final AbsoluteDate date = AbsoluteDate.J2000_EPOCH;
final AbsoluteDate date0 = date.shiftedBy(50000);
final OneAxisEllipsoid earth = new OneAxisEllipsoid(Constants.WGS84_EARTH_EQUATORIAL_RADIUS, Constants.WGS84_EARTH_FLATTENING, FramesFactory.getITRF(IERSConventions.IERS_2010, true));
final GroundStation station = new GroundStation(new TopocentricFrame(earth, new GeodeticPoint(latitude, longitude, altitude), "dummy"));
final DSFactory factory = new DSFactory(parameterPattern.length, 1);
final FieldAbsoluteDate<DerivativeStructure> dateDS = new FieldAbsoluteDate<>(factory.getDerivativeField(), date);
ParameterDriver[] selectedDrivers = new ParameterDriver[parameterPattern.length];
UnivariateDifferentiableVectorFunction[] dFAngular = new UnivariateDifferentiableVectorFunction[parameterPattern.length];
final ParameterDriver[] allDrivers = selectAllDrivers(station);
for (ParameterDriver driver : allDrivers) {
driver.setReferenceDate(date0);
}
Map<String, Integer> indices = new HashMap<>();
for (int k = 0; k < dFAngular.length; ++k) {
for (int i = 0; i < allDrivers.length; ++i) {
if (allDrivers[i].getName().matches(parameterPattern[k])) {
selectedDrivers[k] = allDrivers[i];
dFAngular[k] = differentiatedTransformAngular(station, eme2000, date, selectedDrivers[k], stepFactor);
indices.put(selectedDrivers[k].getName(), k);
}
}
}
;
DSFactory factory11 = new DSFactory(1, 1);
RandomGenerator generator = new Well19937a(0xa01a1d8fe5d80af7l);
double maxRotationValueError = 0;
double maxRotationDerivativeError = 0;
double maxRotationRateValueError = 0;
double maxRotationRateDerivativeError = 0;
for (int i = 0; i < 1000; ++i) {
// randomly change one parameter
ParameterDriver changed = allDrivers[generator.nextInt(allDrivers.length)];
changed.setNormalizedValue(2 * generator.nextDouble() - 1);
// transform to check
FieldTransform<DerivativeStructure> t = station.getOffsetToInertial(eme2000, dateDS, factory, indices);
for (int k = 0; k < dFAngular.length; ++k) {
// reference values and derivatives computed using finite differences
DerivativeStructure[] refAngular = dFAngular[k].value(factory11.variable(0, selectedDrivers[k].getValue()));
// rotation
final Rotation refQ = new Rotation(refAngular[0].getValue(), refAngular[1].getValue(), refAngular[2].getValue(), refAngular[3].getValue(), true);
final Rotation resQ = t.getRotation().toRotation();
maxRotationValueError = FastMath.max(maxRotationValueError, Rotation.distance(refQ, resQ));
double sign = FastMath.copySign(1.0, refAngular[0].getValue() * t.getRotation().getQ0().getValue() + refAngular[1].getValue() * t.getRotation().getQ1().getValue() + refAngular[2].getValue() * t.getRotation().getQ2().getValue() + refAngular[3].getValue() * t.getRotation().getQ3().getValue());
maxRotationDerivativeError = FastMath.max(maxRotationDerivativeError, FastMath.abs(sign * refAngular[0].getPartialDerivative(1) - t.getRotation().getQ0().getAllDerivatives()[k + 1]));
maxRotationDerivativeError = FastMath.max(maxRotationDerivativeError, FastMath.abs(sign * refAngular[1].getPartialDerivative(1) - t.getRotation().getQ1().getAllDerivatives()[k + 1]));
maxRotationDerivativeError = FastMath.max(maxRotationDerivativeError, FastMath.abs(sign * refAngular[2].getPartialDerivative(1) - t.getRotation().getQ2().getAllDerivatives()[k + 1]));
maxRotationDerivativeError = FastMath.max(maxRotationDerivativeError, FastMath.abs(sign * refAngular[3].getPartialDerivative(1) - t.getRotation().getQ3().getAllDerivatives()[k + 1]));
// rotation rate
final Vector3D refRate = new Vector3D(refAngular[4].getValue(), refAngular[5].getValue(), refAngular[6].getValue());
final Vector3D resRate = t.getRotationRate().toVector3D();
final Vector3D refRateD = new Vector3D(refAngular[4].getPartialDerivative(1), refAngular[5].getPartialDerivative(1), refAngular[6].getPartialDerivative(1));
final Vector3D resRateD = new Vector3D(t.getRotationRate().getX().getAllDerivatives()[k + 1], t.getRotationRate().getY().getAllDerivatives()[k + 1], t.getRotationRate().getZ().getAllDerivatives()[k + 1]);
maxRotationRateValueError = FastMath.max(maxRotationRateValueError, Vector3D.distance(refRate, resRate));
maxRotationRateDerivativeError = FastMath.max(maxRotationRateDerivativeError, Vector3D.distance(refRateD, resRateD));
}
}
Assert.assertEquals(0.0, maxRotationValueError, toleranceRotationValue);
Assert.assertEquals(0.0, maxRotationDerivativeError, toleranceRotationDerivative);
Assert.assertEquals(0.0, maxRotationRateValueError, toleranceRotationRateValue);
Assert.assertEquals(0.0, maxRotationRateDerivativeError, toleranceRotationRateDerivative);
}
use of org.orekit.bodies.OneAxisEllipsoid in project Orekit by CS-SI.
the class SolarRadiationPressureTest method RealFieldIsotropicTest.
/**
*Testing if the propagation between the FieldPropagation and the propagation
* is equivalent.
* Also testing if propagating X+dX with the propagation is equivalent to
* propagation X with the FieldPropagation and then applying the taylor
* expansion of dX to the result.
*/
@Test
public void RealFieldIsotropicTest() throws OrekitException {
DSFactory factory = new DSFactory(6, 5);
DerivativeStructure a_0 = factory.variable(0, 7e7);
DerivativeStructure e_0 = factory.variable(1, 0.4);
DerivativeStructure i_0 = factory.variable(2, 85 * FastMath.PI / 180);
DerivativeStructure R_0 = factory.variable(3, 0.7);
DerivativeStructure O_0 = factory.variable(4, 0.5);
DerivativeStructure n_0 = factory.variable(5, 0.1);
Field<DerivativeStructure> field = a_0.getField();
DerivativeStructure zero = field.getZero();
FieldAbsoluteDate<DerivativeStructure> J2000 = FieldAbsoluteDate.getJ2000Epoch(field);
Frame EME = FramesFactory.getEME2000();
FieldKeplerianOrbit<DerivativeStructure> FKO = new FieldKeplerianOrbit<>(a_0, e_0, i_0, R_0, O_0, n_0, PositionAngle.MEAN, EME, J2000, Constants.EIGEN5C_EARTH_MU);
FieldSpacecraftState<DerivativeStructure> initialState = new FieldSpacecraftState<>(FKO);
SpacecraftState iSR = initialState.toSpacecraftState();
final OrbitType type = OrbitType.KEPLERIAN;
double[][] tolerance = NumericalPropagator.tolerances(10.0, FKO.toOrbit(), type);
AdaptiveStepsizeFieldIntegrator<DerivativeStructure> integrator = new DormandPrince853FieldIntegrator<>(field, 0.001, 200, tolerance[0], tolerance[1]);
integrator.setInitialStepSize(zero.add(60));
AdaptiveStepsizeIntegrator RIntegrator = new DormandPrince853Integrator(0.001, 200, tolerance[0], tolerance[1]);
RIntegrator.setInitialStepSize(60);
FieldNumericalPropagator<DerivativeStructure> FNP = new FieldNumericalPropagator<>(field, integrator);
FNP.setOrbitType(type);
FNP.setInitialState(initialState);
NumericalPropagator NP = new NumericalPropagator(RIntegrator);
NP.setOrbitType(type);
NP.setInitialState(iSR);
PVCoordinatesProvider sun = CelestialBodyFactory.getSun();
// creation of the force model
OneAxisEllipsoid earth = new OneAxisEllipsoid(6378136.46, 1.0 / 298.25765, FramesFactory.getITRF(IERSConventions.IERS_2010, true));
SolarRadiationPressure forceModel = new SolarRadiationPressure(sun, earth.getEquatorialRadius(), new IsotropicRadiationCNES95Convention(500.0, 0.7, 0.7));
FNP.addForceModel(forceModel);
NP.addForceModel(forceModel);
FieldAbsoluteDate<DerivativeStructure> target = J2000.shiftedBy(1000.);
FieldSpacecraftState<DerivativeStructure> finalState_DS = FNP.propagate(target);
SpacecraftState finalState_R = NP.propagate(target.toAbsoluteDate());
FieldPVCoordinates<DerivativeStructure> finPVC_DS = finalState_DS.getPVCoordinates();
PVCoordinates finPVC_R = finalState_R.getPVCoordinates();
Assert.assertEquals(0, Vector3D.distance(finPVC_DS.toPVCoordinates().getPosition(), finPVC_R.getPosition()), 4.0e-9);
long number = 23091991;
RandomGenerator RG = new Well19937a(number);
GaussianRandomGenerator NGG = new GaussianRandomGenerator(RG);
UncorrelatedRandomVectorGenerator URVG = new UncorrelatedRandomVectorGenerator(new double[] { 0.0, 0.0, 0.0, 0.0, 0.0, 0.0 }, new double[] { 1e3, 0.01, 0.01, 0.01, 0.01, 0.01 }, NGG);
double a_R = a_0.getReal();
double e_R = e_0.getReal();
double i_R = i_0.getReal();
double R_R = R_0.getReal();
double O_R = O_0.getReal();
double n_R = n_0.getReal();
for (int ii = 0; ii < 1; ii++) {
double[] rand_next = URVG.nextVector();
double a_shift = a_R + rand_next[0];
double e_shift = e_R + rand_next[1];
double i_shift = i_R + rand_next[2];
double R_shift = R_R + rand_next[3];
double O_shift = O_R + rand_next[4];
double n_shift = n_R + rand_next[5];
KeplerianOrbit shiftedOrb = new KeplerianOrbit(a_shift, e_shift, i_shift, R_shift, O_shift, n_shift, PositionAngle.MEAN, EME, J2000.toAbsoluteDate(), Constants.EIGEN5C_EARTH_MU);
SpacecraftState shift_iSR = new SpacecraftState(shiftedOrb);
NumericalPropagator shift_NP = new NumericalPropagator(RIntegrator);
shift_NP.setOrbitType(type);
shift_NP.setInitialState(shift_iSR);
shift_NP.addForceModel(forceModel);
SpacecraftState finalState_shift = shift_NP.propagate(target.toAbsoluteDate());
PVCoordinates finPVC_shift = finalState_shift.getPVCoordinates();
// position check
FieldVector3D<DerivativeStructure> pos_DS = finPVC_DS.getPosition();
double x_DS = pos_DS.getX().taylor(rand_next[0], rand_next[1], rand_next[2], rand_next[3], rand_next[4], rand_next[5]);
double y_DS = pos_DS.getY().taylor(rand_next[0], rand_next[1], rand_next[2], rand_next[3], rand_next[4], rand_next[5]);
double z_DS = pos_DS.getZ().taylor(rand_next[0], rand_next[1], rand_next[2], rand_next[3], rand_next[4], rand_next[5]);
// System.out.println(pos_DS.getX().getPartialDerivative(1));
double x = finPVC_shift.getPosition().getX();
double y = finPVC_shift.getPosition().getY();
double z = finPVC_shift.getPosition().getZ();
Assert.assertEquals(x_DS, x, FastMath.abs(x - pos_DS.getX().getReal()) * 4e-9);
Assert.assertEquals(y_DS, y, FastMath.abs(y - pos_DS.getY().getReal()) * 5e-9);
Assert.assertEquals(z_DS, z, FastMath.abs(z - pos_DS.getZ().getReal()) * 6e-10);
// velocity check
FieldVector3D<DerivativeStructure> vel_DS = finPVC_DS.getVelocity();
double vx_DS = vel_DS.getX().taylor(rand_next[0], rand_next[1], rand_next[2], rand_next[3], rand_next[4], rand_next[5]);
double vy_DS = vel_DS.getY().taylor(rand_next[0], rand_next[1], rand_next[2], rand_next[3], rand_next[4], rand_next[5]);
double vz_DS = vel_DS.getZ().taylor(rand_next[0], rand_next[1], rand_next[2], rand_next[3], rand_next[4], rand_next[5]);
double vx = finPVC_shift.getVelocity().getX();
double vy = finPVC_shift.getVelocity().getY();
double vz = finPVC_shift.getVelocity().getZ();
Assert.assertEquals(vx_DS, vx, FastMath.abs(vx) * 5e-11);
Assert.assertEquals(vy_DS, vy, FastMath.abs(vy) * 3e-10);
Assert.assertEquals(vz_DS, vz, FastMath.abs(vz) * 5e-11);
// acceleration check
FieldVector3D<DerivativeStructure> acc_DS = finPVC_DS.getAcceleration();
double ax_DS = acc_DS.getX().taylor(rand_next[0], rand_next[1], rand_next[2], rand_next[3], rand_next[4], rand_next[5]);
double ay_DS = acc_DS.getY().taylor(rand_next[0], rand_next[1], rand_next[2], rand_next[3], rand_next[4], rand_next[5]);
double az_DS = acc_DS.getZ().taylor(rand_next[0], rand_next[1], rand_next[2], rand_next[3], rand_next[4], rand_next[5]);
double ax = finPVC_shift.getAcceleration().getX();
double ay = finPVC_shift.getAcceleration().getY();
double az = finPVC_shift.getAcceleration().getZ();
Assert.assertEquals(ax_DS, ax, FastMath.abs(ax) * 2e-10);
Assert.assertEquals(ay_DS, ay, FastMath.abs(ay) * 4e-10);
Assert.assertEquals(az_DS, az, FastMath.abs(az) * 7e-10);
}
}
Aggregations