List of usage examples for org.apache.commons.math3.ode.nonstiff AdaptiveStepsizeIntegrator setInitialStepSize
public void setInitialStepSize(final double initialStepSize)
From source file:Ephemeris.mySlavePropagator.java
public static void propagate(double a, double e, double i, double omega, double raan, double lM, double stepS, double duration, boolean solarBL, boolean aeroBL, JTextArea TF) { try {//from w ww . ja va 2s.co m // configure Orekit Autoconfiguration.configureOrekit(); // Inertial frame Frame inertialFrame = FramesFactory.getEME2000(); // Initial date in UTC time scale TimeScale utc = TimeScalesFactory.getUTC(); AbsoluteDate initialDate = new AbsoluteDate(2016, 01, 01, 12, 00, 00.000, utc); // gravitation coefficient double mu = 3.986004415e+14; // Orbit construction as Keplerian Orbit initialOrbit = new KeplerianOrbit(a, e, i, omega, raan, lM, PositionAngle.MEAN, inertialFrame, initialDate, mu); // Initial state definition SpacecraftState initialState = new SpacecraftState(initialOrbit); // Adaptive step integrator with a minimum step of 0.001 and a maximum step of 1000 final double minStep = 0.001; final double maxstep = 1000.0; final double positionTolerance = 0.1; final double initStep = stepS; final OrbitType propagationType = OrbitType.KEPLERIAN; final double[][] tolerances = NumericalPropagator.tolerances(positionTolerance, initialOrbit, propagationType); AdaptiveStepsizeIntegrator integrator = new DormandPrince853Integrator(minStep, maxstep, tolerances[0], tolerances[1]); integrator.setInitialStepSize(initStep); NumericalPropagator propagator = new NumericalPropagator(integrator); propagator.setOrbitType(propagationType); // Force Model (reduced to perturbing gravity field) final NormalizedSphericalHarmonicsProvider provider = GravityFieldFactory.getNormalizedProvider(10, 10); ForceModel holmesFeatherstone = new HolmesFeatherstoneAttractionModel( FramesFactory.getITRF(IERSConventions.IERS_2010, true), provider); // Aero drag final Frame earthFrame = CelestialBodyFactory.getEarth().getBodyOrientedFrame(); final double ae = provider.getAe(); final SphericalSpacecraft ssc = new SphericalSpacecraft(1, 0.47, 0., 1.2); final OneAxisEllipsoid earth = new OneAxisEllipsoid(ae, Constants.WGS84_EARTH_FLATTENING, earthFrame); Atmosphere atmosphere = new HarrisPriester(CelestialBodyFactory.getSun(), earth); ForceModel aerodrag = new DragForce(atmosphere, ssc); //Solar radiation pressure ForceModel solarPressure = new SolarRadiationPressure(CelestialBodyFactory.getSun(), ae, ssc); // Add force model to the propagator propagator.addForceModel(holmesFeatherstone); if (aeroBL) propagator.addForceModel(aerodrag); if (solarBL) propagator.addForceModel(solarPressure); // Set up initial state in the propagator propagator.setInitialState(initialState); // Set the propagator to slave mode (could be omitted as it is the default mode) propagator.setSlaveMode(); // Overall duration in seconds for extrapolation //double duration = 630.; // Stop date final AbsoluteDate finalDate = initialDate.shiftedBy(duration); // Step duration in seconds double stepT = stepS; // Extrapolation loop //int cpt = 1; for (AbsoluteDate extrapDate = initialDate; extrapDate .compareTo(finalDate) <= 0; extrapDate = extrapDate.shiftedBy(60 * stepT)) { SpacecraftState currentState = propagator.propagate(extrapDate); //System.out.println(/*"step " + */cpt++); //System.out.println(/*" time : " + */currentState.getDate()); //System.out.println(" " + currentState.getOrbit()); KeplerianOrbit o = (KeplerianOrbit) currentState.getOrbit(); String tempSTR; tempSTR = String.format(Locale.US, "%s %12.3f %10.8f %10.6f %10.6f %10.6f %10.6f%n", currentState.getDate(), o.getA(), o.getE(), FastMath.toDegrees(o.getI()), FastMath.toDegrees(o.getPerigeeArgument()), FastMath.toDegrees(o.getRightAscensionOfAscendingNode()), FastMath.toDegrees(o.getTrueAnomaly())); TF.append(tempSTR); } } catch (OrekitException oe) { System.err.println(oe.getMessage()); } }
From source file:gbt.ubt.tool.Orthorectifier.java
public static BoundedPropagator generateEphemeris(InputParameters params) throws IOException, OrekitException { /*// www .j a v a 2 s . c o m This function generates the satellite ephemeris for the duration of the acquisition using the orbital state vector that is embedded in the product header (MPH). This state vector is normally the restituted orbit from the Flight Operations Segment. The open source Orekit library (https://www.orekit.org/) is used to perform the numerical propagation. The generated ephemeris compares well against the available DORIS Precise Orbit files (maximum 10 metre displacement vector magnitude during propagation duration). */ if (params.orthorectify) { System.out.println("Orthorectification is selected, performing orbit propagation"); /* Get the acquisition duration and state vector from the L1b product */ DataProvidersManager.getInstance().addProvider(new ZipJarCrawler(new File("orekit-data.zip"))); Product readProduct = ProductIO.readProduct(params.inputFileLocation); Date startTime = readProduct.getStartTime().getAsDate(); Date stopTime = readProduct.getEndTime().getAsDate(); MetadataElement metadataRoot = readProduct.getMetadataRoot(); Date vectorTime = metadataRoot.getElementAt(0).getAttributeUTC("STATE_VECTOR_TIME").getAsDate(); Double xPos = metadataRoot.getElementAt(0).getAttributeDouble("X_POSITION"); Double yPos = metadataRoot.getElementAt(0).getAttributeDouble("Y_POSITION"); Double zPos = metadataRoot.getElementAt(0).getAttributeDouble("Z_POSITION"); Double xVel = metadataRoot.getElementAt(0).getAttributeDouble("X_VELOCITY"); Double yVel = metadataRoot.getElementAt(0).getAttributeDouble("Y_VELOCITY"); Double zVel = metadataRoot.getElementAt(0).getAttributeDouble("Z_VELOCITY"); /* Product state vector in fixed frame and UTC Note state vector and product/acquisition start time are not coincident */ Frame fixedFrame = FramesFactory.getITRF(IERSConventions.IERS_2010, false); TimeScale utc = TimeScalesFactory.getUTC(); GregorianCalendar cal = new GregorianCalendar(TimeZone.getTimeZone("UTC")); cal.setTime(vectorTime); AbsoluteDate initialDate = new AbsoluteDate(cal.get(Calendar.YEAR), cal.get(Calendar.MONTH) + 1, cal.get(Calendar.DATE), cal.get(Calendar.HOUR_OF_DAY), cal.get(Calendar.MINUTE), cal.get(Calendar.SECOND), utc); Vector3D position = new Vector3D(xPos, yPos, zPos); Vector3D velocity = new Vector3D(xVel, yVel, zVel); /* Orekit can only propate orbits defined in inertial reference frame Convert state vector to J2000 (EME2000) frame. */ Frame inertialFrame = FramesFactory.getEME2000(); Transform frameTransform = fixedFrame.getTransformTo(inertialFrame, initialDate); PVCoordinates transformPVCoordinates = frameTransform .transformPVCoordinates(new PVCoordinates(position, velocity)); /* Set initial spacecraft state */ double mu = 3.986004415e+14; // Central attraction coefficient Orbit initialOrbit = new CartesianOrbit(transformPVCoordinates, inertialFrame, initialDate, mu); SpacecraftState initialState = new SpacecraftState(initialOrbit, 7892.0); // Orbital parameters and Envisat dry mass /* Set up numerical integrator for propagation */ double minStep = 1; double maxStep = 1000; double initialStep = 60; OrbitType propagationType = OrbitType.CARTESIAN; double[][] tolerance = NumericalPropagator.tolerances(0.001, initialOrbit, propagationType); AdaptiveStepsizeIntegrator integrator = new DormandPrince853Integrator(minStep, maxStep, tolerance[0], tolerance[1]); integrator.setInitialStepSize(initialStep); NumericalPropagator propagator = new NumericalPropagator(integrator); propagator.setOrbitType(propagationType); /* Add force models to the propagator: Consider central body (Earth) gravity pertubation Consider third body (Sun & Moon) gravity pertubation Other pertubations are not considered (computation speed/accuracy trade-off) */ NormalizedSphericalHarmonicsProvider provider = GravityFieldFactory.getNormalizedProvider(21, 21); ForceModel holmesFeatherstone = new HolmesFeatherstoneAttractionModel( FramesFactory.getITRF(IERSConventions.IERS_2010, true), provider); propagator.addForceModel(holmesFeatherstone); propagator.addForceModel(new ThirdBodyAttraction(CelestialBodyFactory.getSun())); propagator.addForceModel(new ThirdBodyAttraction(CelestialBodyFactory.getMoon())); /* Run the propagator to generate the ephemeris */ propagator.setInitialState(initialState); propagator.setEphemerisMode(); cal.setTime(startTime); AbsoluteDate startDate = new AbsoluteDate(cal.get(Calendar.YEAR), cal.get(Calendar.MONTH) + 1, cal.get(Calendar.DATE), cal.get(Calendar.HOUR_OF_DAY), cal.get(Calendar.MINUTE), cal.get(Calendar.SECOND), utc); cal.setTime(stopTime); AbsoluteDate stopDate = new AbsoluteDate(cal.get(Calendar.YEAR), cal.get(Calendar.MONTH) + 1, cal.get(Calendar.DATE), cal.get(Calendar.HOUR_OF_DAY), cal.get(Calendar.MINUTE), cal.get(Calendar.SECOND), utc); SpacecraftState finalState = propagator.propagate(startDate.shiftedBy(-600.0), stopDate.shiftedBy(600.0)); BoundedPropagator ephemeris = propagator.getGeneratedEphemeris(); return ephemeris; } else { return null; } }
From source file:fr.cs.examples.propagation.DSSTPropagation.java
/** Create the numerical propagator * * @param orbit initial orbit/*from w w w.jav a2 s . c om*/ * @param mass S/C mass (kg) * @throws OrekitException */ private NumericalPropagator createNumProp(final Orbit orbit, final double mass) throws OrekitException { final double[][] tol = NumericalPropagator.tolerances(1.0, orbit, orbit.getType()); final double minStep = 1.e-3; final double maxStep = 1.e+3; AdaptiveStepsizeIntegrator integrator = new DormandPrince853Integrator(minStep, maxStep, tol[0], tol[1]); integrator.setInitialStepSize(100.); NumericalPropagator numProp = new NumericalPropagator(integrator); numProp.setInitialState(new SpacecraftState(orbit, mass)); return numProp; }
From source file:org.orekit.forces.gravity.CunninghamAttractionModelTest.java
@Before public void setUp() { itrf2008 = null;//from w w w . j a v a2 s .c om propagator = null; Utils.setDataRoot("regular-data"); try { // Eigen c1 model truncated to degree 6 mu = 3.986004415e+14; ae = 6378136.460; c20 = -1.08262631303e-3; c30 = 2.53248017972e-6; c40 = 1.61994537014e-6; c50 = 2.27888264414e-7; c60 = -5.40618601332e-7; itrf2008 = FramesFactory.getITRF(IERSConventions.IERS_2010, true); double[] absTolerance = { 0.001, 1.0e-9, 1.0e-9, 1.0e-6, 1.0e-6, 1.0e-6, 0.001 }; double[] relTolerance = { 1.0e-7, 1.0e-4, 1.0e-4, 1.0e-7, 1.0e-7, 1.0e-7, 1.0e-7 }; AdaptiveStepsizeIntegrator integrator = new DormandPrince853Integrator(0.001, 1000, absTolerance, relTolerance); integrator.setInitialStepSize(60); propagator = new NumericalPropagator(integrator); } catch (OrekitException oe) { Assert.fail(oe.getMessage()); } }
From source file:org.orekit.forces.gravity.DrozinerAttractionModelTest.java
@Before public void setUp() { itrf2008 = null;//from www. ja v a2 s.c om propagator = null; Utils.setDataRoot("regular-data"); try { mu = 3.986004415e+14; ae = 6378136.460; c20 = -1.08262631303e-3; c30 = 2.53248017972e-6; c40 = 1.61994537014e-6; c50 = 2.27888264414e-7; c60 = -5.40618601332e-7; itrf2008 = FramesFactory.getITRF(IERSConventions.IERS_2010, true); double[] absTolerance = { 0.001, 1.0e-9, 1.0e-9, 1.0e-6, 1.0e-6, 1.0e-6, 0.001 }; double[] relTolerance = { 1.0e-7, 1.0e-4, 1.0e-4, 1.0e-7, 1.0e-7, 1.0e-7, 1.0e-7 }; AdaptiveStepsizeIntegrator integrator = new DormandPrince853Integrator(0.001, 1000, absTolerance, relTolerance); integrator.setInitialStepSize(60); propagator = new NumericalPropagator(integrator); } catch (OrekitException oe) { Assert.fail(oe.getMessage()); } }
From source file:org.orekit.forces.gravity.HolmesFeatherstoneAttractionModelTest.java
@Before public void setUp() { itrf = null;/*from w ww . j a v a2 s . c o m*/ propagator = null; Utils.setDataRoot("regular-data"); try { // Eigen 6s model truncated to degree 6 mu = 3.986004415e+14; ae = 6378136.460; normalizedC20 = -4.84165299820e-04; normalizedC30 = 9.57211326674e-07; normalizedC40 = 5.39990167207e-07; normalizedC50 = 6.86846073356e-08; normalizedC60 = -1.49953256913e-07; unnormalizedC20 = FastMath.sqrt(5) * normalizedC20; unnormalizedC30 = FastMath.sqrt(7) * normalizedC30; unnormalizedC40 = FastMath.sqrt(9) * normalizedC40; unnormalizedC50 = FastMath.sqrt(11) * normalizedC50; unnormalizedC60 = FastMath.sqrt(13) * normalizedC60; itrf = FramesFactory.getITRF(IERSConventions.IERS_2010, true); double[] absTolerance = { 0.001, 1.0e-9, 1.0e-9, 1.0e-6, 1.0e-6, 1.0e-6, 0.001 }; double[] relTolerance = { 1.0e-7, 1.0e-4, 1.0e-4, 1.0e-7, 1.0e-7, 1.0e-7, 1.0e-7 }; AdaptiveStepsizeIntegrator integrator = new DormandPrince853Integrator(0.001, 1000, absTolerance, relTolerance); integrator.setInitialStepSize(60); propagator = new NumericalPropagator(integrator); } catch (OrekitException oe) { Assert.fail(oe.getMessage()); } }
From source file:org.orekit.forces.maneuvers.ConstantThrustManeuverTest.java
@Test public void testRoughBehaviour() throws OrekitException { final double isp = 318; final double mass = 2500; final double a = 24396159; final double e = 0.72831215; final double i = FastMath.toRadians(7); final double omega = FastMath.toRadians(180); final double OMEGA = FastMath.toRadians(261); final double lv = 0; final double duration = 3653.99; final double f = 420; final double delta = FastMath.toRadians(-7.4978); final double alpha = FastMath.toRadians(351); final AttitudeProvider law = new InertialProvider( new Rotation(new Vector3D(alpha, delta), Vector3D.PLUS_I)); final AbsoluteDate initDate = new AbsoluteDate(new DateComponents(2004, 01, 01), new TimeComponents(23, 30, 00.000), TimeScalesFactory.getUTC()); final Orbit orbit = new KeplerianOrbit(a, e, i, omega, OMEGA, lv, PositionAngle.TRUE, FramesFactory.getEME2000(), initDate, mu); final SpacecraftState initialState = new SpacecraftState(orbit, law.getAttitude(orbit, orbit.getDate(), orbit.getFrame()), mass); final AbsoluteDate fireDate = new AbsoluteDate(new DateComponents(2004, 01, 02), new TimeComponents(04, 15, 34.080), TimeScalesFactory.getUTC()); final ConstantThrustManeuver maneuver = new ConstantThrustManeuver(fireDate, duration, f, isp, Vector3D.PLUS_I);/* w ww. j a v a 2 s . c o m*/ Assert.assertEquals(f, maneuver.getThrust(), 1.0e-10); Assert.assertEquals(isp, maneuver.getISP(), 1.0e-10); double[] absTolerance = { 0.001, 1.0e-9, 1.0e-9, 1.0e-6, 1.0e-6, 1.0e-6, 0.001 }; double[] relTolerance = { 1.0e-7, 1.0e-4, 1.0e-4, 1.0e-7, 1.0e-7, 1.0e-7, 1.0e-7 }; AdaptiveStepsizeIntegrator integrator = new DormandPrince853Integrator(0.001, 1000, absTolerance, relTolerance); integrator.setInitialStepSize(60); final NumericalPropagator propagator = new NumericalPropagator(integrator); propagator.setInitialState(initialState); propagator.setAttitudeProvider(law); propagator.addForceModel(maneuver); final SpacecraftState finalorb = propagator.propagate(fireDate.shiftedBy(3800)); final double massTolerance = FastMath.abs(maneuver.getFlowRate()) * maneuver.getEventsDetectors()[0].getThreshold(); Assert.assertEquals(2007.8824544261233, finalorb.getMass(), massTolerance); Assert.assertEquals(2.6872, FastMath.toDegrees(MathUtils.normalizeAngle(finalorb.getI(), FastMath.PI)), 1e-4); Assert.assertEquals(28970, finalorb.getA() / 1000, 1); }
From source file:org.orekit.forces.maneuvers.SmallManeuverAnalyticalModelTest.java
private BoundedPropagator getEphemeris(final Orbit orbit, final double mass, final AbsoluteDate t0, final Vector3D dV, final double f, final double isp) throws OrekitException { AttitudeProvider law = new LofOffset(orbit.getFrame(), LOFType.LVLH); final SpacecraftState initialState = new SpacecraftState(orbit, law.getAttitude(orbit, orbit.getDate(), orbit.getFrame()), mass); // set up numerical propagator final double dP = 1.0; double[][] tolerances = NumericalPropagator.tolerances(dP, orbit, orbit.getType()); AdaptiveStepsizeIntegrator integrator = new DormandPrince853Integrator(0.001, 1000, tolerances[0], tolerances[1]);/* w ww . j a va 2s . c o m*/ integrator.setInitialStepSize(orbit.getKeplerianPeriod() / 100.0); final NumericalPropagator propagator = new NumericalPropagator(integrator); propagator.setOrbitType(orbit.getType()); propagator.setInitialState(initialState); propagator.setAttitudeProvider(law); if (dV.getNorm() > 1.0e-6) { // set up maneuver final double vExhaust = Constants.G0_STANDARD_GRAVITY * isp; final double dt = -(mass * vExhaust / f) * FastMath.expm1(-dV.getNorm() / vExhaust); final ConstantThrustManeuver maneuver = new ConstantThrustManeuver(t0, dt, f, isp, dV.normalize()); propagator.addForceModel(maneuver); } propagator.setEphemerisMode(); propagator.propagate(t0.shiftedBy(5 * orbit.getKeplerianPeriod())); return propagator.getGeneratedEphemeris(); }
From source file:org.orekit.forces.radiation.SolarRadiationPressureTest.java
@Test public void testRoughOrbitalModifs() throws ParseException, OrekitException, FileNotFoundException { // initialization AbsoluteDate date = new AbsoluteDate(new DateComponents(1970, 7, 1), new TimeComponents(13, 59, 27.816), TimeScalesFactory.getUTC()); Orbit orbit = new EquinoctialOrbit(42164000, 10e-3, 10e-3, FastMath.tan(0.001745329) * FastMath.cos(2 * FastMath.PI / 3), FastMath.tan(0.001745329) * FastMath.sin(2 * FastMath.PI / 3), 0.1, PositionAngle.TRUE, FramesFactory.getEME2000(), date, mu); final double period = orbit.getKeplerianPeriod(); Assert.assertEquals(86164, period, 1); 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 SRP = new SolarRadiationPressure(sun, earth.getEquatorialRadius(), (RadiationSensitive) new SphericalSpacecraft(500.0, 0.7, 0.7, 0.7)); // creation of the propagator double[] absTolerance = { 0.1, 1.0e-9, 1.0e-9, 1.0e-5, 1.0e-5, 1.0e-5, 0.001 }; double[] relTolerance = { 1.0e-4, 1.0e-4, 1.0e-4, 1.0e-6, 1.0e-6, 1.0e-6, 1.0e-7 }; AdaptiveStepsizeIntegrator integrator = new DormandPrince853Integrator(900.0, 60000, absTolerance, relTolerance);// ww w. jav a 2 s.c o m integrator.setInitialStepSize(3600); final NumericalPropagator calc = new NumericalPropagator(integrator); calc.addForceModel(SRP); // Step Handler calc.setMasterMode(FastMath.floor(period), new SolarStepHandler()); AbsoluteDate finalDate = date.shiftedBy(10 * period); calc.setInitialState(new SpacecraftState(orbit, 1500.0)); calc.propagate(finalDate); Assert.assertTrue(calc.getCalls() < 7100); }
From source file:org.orekit.propagation.analytical.AdapterPropagatorTest.java
private BoundedPropagator getEphemeris(final Orbit orbit, final double mass, final int nbOrbits, final AttitudeProvider law, final AbsoluteDate t0, final Vector3D dV, final double f, final double isp, final boolean sunAttraction, final boolean moonAttraction, final NormalizedSphericalHarmonicsProvider gravityField) throws OrekitException, ParseException, IOException { SpacecraftState initialState = new SpacecraftState(orbit, law.getAttitude(orbit, orbit.getDate(), orbit.getFrame()), mass); // add some dummy additional states initialState = initialState.addAdditionalState("dummy 1", 1.25, 2.5); initialState = initialState.addAdditionalState("dummy 2", 5.0); // set up numerical propagator final double dP = 1.0; double[][] tolerances = NumericalPropagator.tolerances(dP, orbit, orbit.getType()); AdaptiveStepsizeIntegrator integrator = new DormandPrince853Integrator(0.001, 1000, tolerances[0], tolerances[1]);// ww w. ja va2s .c o m integrator.setInitialStepSize(orbit.getKeplerianPeriod() / 100.0); final NumericalPropagator propagator = new NumericalPropagator(integrator); propagator.addAdditionalStateProvider(new AdditionalStateProvider() { public String getName() { return "dummy 2"; } public double[] getAdditionalState(SpacecraftState state) { return new double[] { 5.0 }; } }); propagator.setInitialState(initialState); propagator.setAttitudeProvider(law); if (dV.getNorm() > 1.0e-6) { // set up maneuver final double vExhaust = Constants.G0_STANDARD_GRAVITY * isp; final double dt = -(mass * vExhaust / f) * FastMath.expm1(-dV.getNorm() / vExhaust); final ConstantThrustManeuver maneuver = new ConstantThrustManeuver(t0.shiftedBy(-0.5 * dt), dt, f, isp, dV.normalize()); propagator.addForceModel(maneuver); } if (sunAttraction) { propagator.addForceModel(new ThirdBodyAttraction(CelestialBodyFactory.getSun())); } if (moonAttraction) { propagator.addForceModel(new ThirdBodyAttraction(CelestialBodyFactory.getMoon())); } if (gravityField != null) { propagator.addForceModel( new HolmesFeatherstoneAttractionModel(FramesFactory.getGTOD(false), gravityField)); } propagator.setEphemerisMode(); propagator.propagate(t0.shiftedBy(nbOrbits * orbit.getKeplerianPeriod())); final BoundedPropagator ephemeris = propagator.getGeneratedEphemeris(); // both the initial propagator and generated ephemeris manage one of the two // additional states, but they also contain unmanaged copies of the other one Assert.assertFalse(propagator.isAdditionalStateManaged("dummy 1")); Assert.assertTrue(propagator.isAdditionalStateManaged("dummy 2")); Assert.assertFalse(ephemeris.isAdditionalStateManaged("dummy 1")); Assert.assertTrue(ephemeris.isAdditionalStateManaged("dummy 2")); Assert.assertEquals(2, ephemeris.getInitialState().getAdditionalState("dummy 1").length); Assert.assertEquals(1, ephemeris.getInitialState().getAdditionalState("dummy 2").length); return ephemeris; }