use of org.orekit.propagation.semianalytical.dsst.forces.DSSTZonal in project Orekit by CS-SI.
the class DSSTPropagatorTest method testOsculatingToMeanState.
@Test
public void testOsculatingToMeanState() throws IllegalArgumentException, OrekitException {
final SpacecraftState meanState = getGEOState();
final UnnormalizedSphericalHarmonicsProvider provider = GravityFieldFactory.getUnnormalizedProvider(2, 0);
final Frame earthFrame = CelestialBodyFactory.getEarth().getBodyOrientedFrame();
DSSTForceModel zonal = new DSSTZonal(provider, 2, 1, 5);
DSSTForceModel tesseral = new DSSTTesseral(earthFrame, Constants.WGS84_EARTH_ANGULAR_VELOCITY, provider, 2, 0, 0, 2, 2, 0, 0);
final Collection<DSSTForceModel> forces = new ArrayList<DSSTForceModel>();
forces.add(zonal);
forces.add(tesseral);
final SpacecraftState osculatingState = DSSTPropagator.computeOsculatingState(meanState, null, forces);
// there are no Gaussian force models, we don't need an attitude provider
final SpacecraftState computedMeanState = DSSTPropagator.computeMeanState(osculatingState, null, forces);
Assert.assertEquals(meanState.getA(), computedMeanState.getA(), 2.0e-8);
Assert.assertEquals(0.0, Vector3D.distance(meanState.getPVCoordinates().getPosition(), computedMeanState.getPVCoordinates().getPosition()), 2.0e-8);
}
use of org.orekit.propagation.semianalytical.dsst.forces.DSSTZonal in project Orekit by CS-SI.
the class DSSTPropagation method setForceModel.
/**
* Set DSST propagator force models
*
* @param parser input file parser
* @param unnormalized spherical harmonics provider
* @param earthFrame Earth rotating frame
* @param rotationRate central body rotation rate (rad/s)
* @param dsstProp DSST propagator
* @throws IOException
* @throws OrekitException
*/
private void setForceModel(final KeyValueFileParser<ParameterKey> parser, final UnnormalizedSphericalHarmonicsProvider unnormalized, final Frame earthFrame, final double rotationRate, final DSSTPropagator dsstProp) throws IOException, OrekitException {
final double ae = unnormalized.getAe();
final int degree = parser.getInt(ParameterKey.CENTRAL_BODY_DEGREE);
final int order = parser.getInt(ParameterKey.CENTRAL_BODY_ORDER);
if (order > degree) {
throw new IOException("Potential order cannot be higher than potential degree");
}
// Central Body Force Model with un-normalized coefficients
dsstProp.addForceModel(new DSSTZonal(unnormalized, parser.getInt(ParameterKey.MAX_DEGREE_ZONAL_SHORT_PERIODS), parser.getInt(ParameterKey.MAX_ECCENTRICITY_POWER_ZONAL_SHORT_PERIODS), parser.getInt(ParameterKey.MAX_FREQUENCY_TRUE_LONGITUDE_ZONAL_SHORT_PERIODS)));
dsstProp.addForceModel(new DSSTTesseral(earthFrame, rotationRate, unnormalized, parser.getInt(ParameterKey.MAX_DEGREE_TESSERAL_SHORT_PERIODS), parser.getInt(ParameterKey.MAX_ORDER_TESSERAL_SHORT_PERIODS), parser.getInt(ParameterKey.MAX_ECCENTRICITY_POWER_TESSERAL_SHORT_PERIODS), parser.getInt(ParameterKey.MAX_FREQUENCY_MEAN_LONGITUDE_TESSERAL_SHORT_PERIODS), parser.getInt(ParameterKey.MAX_DEGREE_TESSERAL_M_DAILIES_SHORT_PERIODS), parser.getInt(ParameterKey.MAX_ORDER_TESSERAL_M_DAILIES_SHORT_PERIODS), parser.getInt(ParameterKey.MAX_ECCENTRICITY_POWER_TESSERAL_M_DAILIES_SHORT_PERIODS)));
// 3rd body (SUN)
if (parser.containsKey(ParameterKey.THIRD_BODY_SUN) && parser.getBoolean(ParameterKey.THIRD_BODY_SUN)) {
dsstProp.addForceModel(new DSSTThirdBody(CelestialBodyFactory.getSun()));
}
// 3rd body (MOON)
if (parser.containsKey(ParameterKey.THIRD_BODY_MOON) && parser.getBoolean(ParameterKey.THIRD_BODY_MOON)) {
dsstProp.addForceModel(new DSSTThirdBody(CelestialBodyFactory.getMoon()));
}
// Drag
if (parser.containsKey(ParameterKey.DRAG) && parser.getBoolean(ParameterKey.DRAG)) {
final OneAxisEllipsoid earth = new OneAxisEllipsoid(ae, Constants.WGS84_EARTH_FLATTENING, earthFrame);
final Atmosphere atm = new HarrisPriester(CelestialBodyFactory.getSun(), earth, 6);
dsstProp.addForceModel(new DSSTAtmosphericDrag(atm, parser.getDouble(ParameterKey.DRAG_CD), parser.getDouble(ParameterKey.DRAG_SF)));
}
// Solar Radiation Pressure
if (parser.containsKey(ParameterKey.SOLAR_RADIATION_PRESSURE) && parser.getBoolean(ParameterKey.SOLAR_RADIATION_PRESSURE)) {
dsstProp.addForceModel(new DSSTSolarRadiationPressure(parser.getDouble(ParameterKey.SOLAR_RADIATION_PRESSURE_CR), parser.getDouble(ParameterKey.SOLAR_RADIATION_PRESSURE_SF), CelestialBodyFactory.getSun(), ae));
}
}
use of org.orekit.propagation.semianalytical.dsst.forces.DSSTZonal in project Orekit by CS-SI.
the class DSSTPropagatorTest method testIssue157.
@Test
public void testIssue157() throws OrekitException {
Utils.setDataRoot("regular-data:potential/icgem-format");
GravityFieldFactory.addPotentialCoefficientsReader(new ICGEMFormatReader("^eigen-6s-truncated$", false));
UnnormalizedSphericalHarmonicsProvider nshp = GravityFieldFactory.getUnnormalizedProvider(8, 8);
Orbit orbit = new KeplerianOrbit(13378000, 0.05, 0, 0, FastMath.PI, 0, PositionAngle.MEAN, FramesFactory.getTOD(false), new AbsoluteDate(2003, 5, 6, TimeScalesFactory.getUTC()), nshp.getMu());
double period = orbit.getKeplerianPeriod();
double[][] tolerance = DSSTPropagator.tolerances(1.0, orbit);
AdaptiveStepsizeIntegrator integrator = new DormandPrince853Integrator(period / 100, period * 100, tolerance[0], tolerance[1]);
integrator.setInitialStepSize(10 * period);
DSSTPropagator propagator = new DSSTPropagator(integrator, true);
OneAxisEllipsoid earth = new OneAxisEllipsoid(Constants.WGS84_EARTH_EQUATORIAL_RADIUS, Constants.WGS84_EARTH_FLATTENING, FramesFactory.getGTOD(false));
CelestialBody sun = CelestialBodyFactory.getSun();
CelestialBody moon = CelestialBodyFactory.getMoon();
propagator.addForceModel(new DSSTZonal(nshp, 8, 7, 17));
propagator.addForceModel(new DSSTTesseral(earth.getBodyFrame(), Constants.WGS84_EARTH_ANGULAR_VELOCITY, nshp, 8, 8, 4, 12, 8, 8, 4));
propagator.addForceModel(new DSSTThirdBody(sun));
propagator.addForceModel(new DSSTThirdBody(moon));
propagator.addForceModel(new DSSTAtmosphericDrag(new HarrisPriester(sun, earth), 2.1, 180));
propagator.addForceModel(new DSSTSolarRadiationPressure(1.2, 180, sun, earth.getEquatorialRadius()));
propagator.setInitialState(new SpacecraftState(orbit, 45.0), true);
SpacecraftState finalState = propagator.propagate(orbit.getDate().shiftedBy(30 * Constants.JULIAN_DAY));
// the following comparison is in fact meaningless
// the initial orbit is osculating the final orbit is a mean orbit
// and they are not considered at the same epoch
// we keep it only as is was an historical test
Assert.assertEquals(2189.4, orbit.getA() - finalState.getA(), 1.0);
propagator.setInitialState(new SpacecraftState(orbit, 45.0), false);
finalState = propagator.propagate(orbit.getDate().shiftedBy(30 * Constants.JULIAN_DAY));
// the following comparison is realistic
// both the initial orbit and final orbit are mean orbits
Assert.assertEquals(1478.05, orbit.getA() - finalState.getA(), 1.0);
}
use of org.orekit.propagation.semianalytical.dsst.forces.DSSTZonal in project Orekit by CS-SI.
the class DSSTPropagatorTest method testShortPeriodCoefficients.
@Test
public void testShortPeriodCoefficients() throws OrekitException {
Utils.setDataRoot("regular-data:potential/icgem-format");
GravityFieldFactory.addPotentialCoefficientsReader(new ICGEMFormatReader("^eigen-6s-truncated$", false));
UnnormalizedSphericalHarmonicsProvider nshp = GravityFieldFactory.getUnnormalizedProvider(4, 4);
Orbit orbit = new KeplerianOrbit(13378000, 0.05, 0, 0, FastMath.PI, 0, PositionAngle.MEAN, FramesFactory.getTOD(false), new AbsoluteDate(2003, 5, 6, TimeScalesFactory.getUTC()), nshp.getMu());
double period = orbit.getKeplerianPeriod();
double[][] tolerance = DSSTPropagator.tolerances(1.0, orbit);
AdaptiveStepsizeIntegrator integrator = new DormandPrince853Integrator(period / 100, period * 100, tolerance[0], tolerance[1]);
integrator.setInitialStepSize(10 * period);
DSSTPropagator propagator = new DSSTPropagator(integrator, false);
OneAxisEllipsoid earth = new OneAxisEllipsoid(Constants.WGS84_EARTH_EQUATORIAL_RADIUS, Constants.WGS84_EARTH_FLATTENING, FramesFactory.getGTOD(false));
CelestialBody sun = CelestialBodyFactory.getSun();
CelestialBody moon = CelestialBodyFactory.getMoon();
propagator.addForceModel(new DSSTZonal(nshp, 4, 3, 9));
propagator.addForceModel(new DSSTTesseral(earth.getBodyFrame(), Constants.WGS84_EARTH_ANGULAR_VELOCITY, nshp, 4, 4, 4, 8, 4, 4, 2));
propagator.addForceModel(new DSSTThirdBody(sun));
propagator.addForceModel(new DSSTThirdBody(moon));
propagator.addForceModel(new DSSTAtmosphericDrag(new HarrisPriester(sun, earth), 2.1, 180));
propagator.addForceModel(new DSSTSolarRadiationPressure(1.2, 180, sun, earth.getEquatorialRadius()));
final AbsoluteDate finalDate = orbit.getDate().shiftedBy(30 * Constants.JULIAN_DAY);
propagator.resetInitialState(new SpacecraftState(orbit, 45.0));
final SpacecraftState stateNoConfig = propagator.propagate(finalDate);
Assert.assertEquals(0, stateNoConfig.getAdditionalStates().size());
propagator.setSelectedCoefficients(new HashSet<String>());
propagator.resetInitialState(new SpacecraftState(orbit, 45.0));
final SpacecraftState stateConfigEmpty = propagator.propagate(finalDate);
Assert.assertEquals(234, stateConfigEmpty.getAdditionalStates().size());
final Set<String> selected = new HashSet<String>();
selected.add("DSST-3rd-body-Moon-s[7]");
selected.add("DSST-central-body-tesseral-c[-2][3]");
propagator.setSelectedCoefficients(selected);
propagator.resetInitialState(new SpacecraftState(orbit, 45.0));
final SpacecraftState stateConfigeSelected = propagator.propagate(finalDate);
Assert.assertEquals(selected.size(), stateConfigeSelected.getAdditionalStates().size());
propagator.setSelectedCoefficients(null);
propagator.resetInitialState(new SpacecraftState(orbit, 45.0));
final SpacecraftState stateConfigNull = propagator.propagate(finalDate);
Assert.assertEquals(0, stateConfigNull.getAdditionalStates().size());
}
use of org.orekit.propagation.semianalytical.dsst.forces.DSSTZonal in project Orekit by CS-SI.
the class DSSTPropagatorTest method testPropagationWithDrag.
@Test
public void testPropagationWithDrag() throws OrekitException {
// Central Body geopotential 2x0
final UnnormalizedSphericalHarmonicsProvider provider = GravityFieldFactory.getUnnormalizedProvider(2, 0);
final Frame earthFrame = CelestialBodyFactory.getEarth().getBodyOrientedFrame();
DSSTForceModel zonal = new DSSTZonal(provider, 2, 0, 5);
DSSTForceModel tesseral = new DSSTTesseral(earthFrame, Constants.WGS84_EARTH_ANGULAR_VELOCITY, provider, 2, 0, 0, 2, 2, 0, 0);
// Drag Force Model
final OneAxisEllipsoid earth = new OneAxisEllipsoid(provider.getAe(), Constants.WGS84_EARTH_FLATTENING, earthFrame);
final Atmosphere atm = new HarrisPriester(CelestialBodyFactory.getSun(), earth, 6);
final double cd = 2.0;
final double area = 25.0;
DSSTForceModel drag = new DSSTAtmosphericDrag(atm, cd, area);
// LEO Orbit
final AbsoluteDate initDate = new AbsoluteDate(2003, 7, 1, 0, 0, 00.000, TimeScalesFactory.getUTC());
final Orbit orbit = new KeplerianOrbit(7204535.848109440, 0.0012402238462686, FastMath.toRadians(98.74341600466740), FastMath.toRadians(111.1990175076630), FastMath.toRadians(43.32990110790340), FastMath.toRadians(68.66852509725620), PositionAngle.MEAN, FramesFactory.getEME2000(), initDate, provider.getMu());
// Set propagator with state and force model
setDSSTProp(new SpacecraftState(orbit));
dsstProp.addForceModel(zonal);
dsstProp.addForceModel(tesseral);
dsstProp.addForceModel(drag);
// 5 days propagation
final SpacecraftState state = dsstProp.propagate(initDate.shiftedBy(5. * 86400.));
// Ref Standalone_DSST:
// a = 7204521.657141485 m
// h/ey = 0.0007093755541595772
// k/ex = -0.001016800430994036
// p/hy = 0.8698955648709271
// q/hx = 0.7757573478894775
// lM = 193°0939742953394
Assert.assertEquals(7204521.657141485, state.getA(), 6.e-1);
Assert.assertEquals(-0.001016800430994036, state.getEquinoctialEx(), 5.e-8);
Assert.assertEquals(0.0007093755541595772, state.getEquinoctialEy(), 2.e-8);
Assert.assertEquals(0.7757573478894775, state.getHx(), 5.e-8);
Assert.assertEquals(0.8698955648709271, state.getHy(), 5.e-8);
Assert.assertEquals(193.0939742953394, FastMath.toDegrees(MathUtils.normalizeAngle(state.getLM(), FastMath.PI)), 2.e-3);
// Assert.assertEquals(((DSSTAtmosphericDrag)drag).getCd(), cd, 1e-9);
// Assert.assertEquals(((DSSTAtmosphericDrag)drag).getArea(), area, 1e-9);
Assert.assertEquals(((DSSTAtmosphericDrag) drag).getAtmosphere(), atm);
// DSSTAtmosphericDrag.ATMOSPHERE_ALTITUDE_MAX
final double atmosphericMaxConstant = 1000000.0;
Assert.assertEquals(((DSSTAtmosphericDrag) drag).getRbar(), atmosphericMaxConstant + Constants.WGS84_EARTH_EQUATORIAL_RADIUS, 1e-9);
}
Aggregations