List of usage examples for org.apache.commons.math3.geometry.euclidean.threed Vector3D getNorm
public double getNorm()
From source file:org.orekit.propagation.analytical.AdapterPropagatorTest.java
@Test public void testLowEarthOrbit() throws OrekitException, ParseException, IOException { Orbit leo = new CircularOrbit( 7200000.0, -1.0e-5, 2.0e-4, FastMath.toRadians(98.0), FastMath.toRadians(123.456), 0.0, PositionAngle.MEAN, FramesFactory.getEME2000(), new AbsoluteDate(new DateComponents(2004, 01, 01), new TimeComponents(23, 30, 00.000), TimeScalesFactory.getUTC()), Constants.EIGEN5C_EARTH_MU); double mass = 5600.0; AbsoluteDate t0 = leo.getDate().shiftedBy(1000.0); Vector3D dV = new Vector3D(-0.1, 0.2, 0.3); double f = 20.0; double isp = 315.0; double vExhaust = Constants.G0_STANDARD_GRAVITY * isp; double dt = -(mass * vExhaust / f) * FastMath.expm1(-dV.getNorm() / vExhaust); BoundedPropagator withoutManeuver = getEphemeris(leo, mass, 5, new LofOffset(leo.getFrame(), LOFType.LVLH), t0, Vector3D.ZERO, f, isp, false, false, null); BoundedPropagator withManeuver = getEphemeris(leo, mass, 5, new LofOffset(leo.getFrame(), LOFType.LVLH), t0, dV, f, isp, false, false, null); // we set up a model that reverts the maneuvers AdapterPropagator adapterPropagator = new AdapterPropagator(withManeuver); AdapterPropagator.DifferentialEffect effect = new SmallManeuverAnalyticalModel( adapterPropagator.propagate(t0), dV.negate(), isp); adapterPropagator.addEffect(effect); adapterPropagator.addAdditionalStateProvider(new AdditionalStateProvider() { public String getName() { return "dummy 3"; }/*from w w w. jav a2s . c o m*/ public double[] getAdditionalState(SpacecraftState state) { return new double[3]; } }); // the adapted propagators do not manage the additional states from the reference, // they simply forward them Assert.assertFalse(adapterPropagator.isAdditionalStateManaged("dummy 1")); Assert.assertFalse(adapterPropagator.isAdditionalStateManaged("dummy 2")); Assert.assertTrue(adapterPropagator.isAdditionalStateManaged("dummy 3")); for (AbsoluteDate t = t0.shiftedBy(0.5 * dt); t.compareTo(withoutManeuver.getMaxDate()) < 0; t = t .shiftedBy(60.0)) { PVCoordinates pvWithout = withoutManeuver.getPVCoordinates(t, leo.getFrame()); PVCoordinates pvReverted = adapterPropagator.getPVCoordinates(t, leo.getFrame()); double revertError = new PVCoordinates(pvWithout, pvReverted).getPosition().getNorm(); Assert.assertEquals(0, revertError, 0.45); Assert.assertEquals(2, adapterPropagator.propagate(t).getAdditionalState("dummy 1").length); Assert.assertEquals(1, adapterPropagator.propagate(t).getAdditionalState("dummy 2").length); Assert.assertEquals(3, adapterPropagator.propagate(t).getAdditionalState("dummy 3").length); } }
From source file:org.orekit.propagation.analytical.AdapterPropagatorTest.java
@Test public void testEccentricOrbit() throws OrekitException, ParseException, IOException { Orbit heo = new KeplerianOrbit(90000000.0, 0.92, FastMath.toRadians(98.0), FastMath.toRadians(12.3456), FastMath.toRadians(123.456), FastMath.toRadians(1.23456), PositionAngle.MEAN, FramesFactory.getEME2000(), new AbsoluteDate(new DateComponents(2004, 01, 01), new TimeComponents(23, 30, 00.000), TimeScalesFactory.getUTC()), Constants.EIGEN5C_EARTH_MU); double mass = 5600.0; AbsoluteDate t0 = heo.getDate().shiftedBy(1000.0); Vector3D dV = new Vector3D(-0.01, 0.02, 0.03); double f = 20.0; double isp = 315.0; double vExhaust = Constants.G0_STANDARD_GRAVITY * isp; double dt = -(mass * vExhaust / f) * FastMath.expm1(-dV.getNorm() / vExhaust); BoundedPropagator withoutManeuver = getEphemeris(heo, mass, 5, new LofOffset(heo.getFrame(), LOFType.LVLH), t0, Vector3D.ZERO, f, isp, false, false, null); BoundedPropagator withManeuver = getEphemeris(heo, mass, 5, new LofOffset(heo.getFrame(), LOFType.LVLH), t0, dV, f, isp, false, false, null); // we set up a model that reverts the maneuvers AdapterPropagator adapterPropagator = new AdapterPropagator(withManeuver); AdapterPropagator.DifferentialEffect effect = new SmallManeuverAnalyticalModel( adapterPropagator.propagate(t0), dV.negate(), isp); adapterPropagator.addEffect(effect); adapterPropagator.addAdditionalStateProvider(new AdditionalStateProvider() { public String getName() { return "dummy 3"; }//from w ww . ja v a 2s . co m public double[] getAdditionalState(SpacecraftState state) { return new double[3]; } }); // the adapted propagators do not manage the additional states from the reference, // they simply forward them Assert.assertFalse(adapterPropagator.isAdditionalStateManaged("dummy 1")); Assert.assertFalse(adapterPropagator.isAdditionalStateManaged("dummy 2")); Assert.assertTrue(adapterPropagator.isAdditionalStateManaged("dummy 3")); for (AbsoluteDate t = t0.shiftedBy(0.5 * dt); t.compareTo(withoutManeuver.getMaxDate()) < 0; t = t .shiftedBy(300.0)) { PVCoordinates pvWithout = withoutManeuver.getPVCoordinates(t, heo.getFrame()); PVCoordinates pvReverted = adapterPropagator.getPVCoordinates(t, heo.getFrame()); double revertError = new PVCoordinates(pvWithout, pvReverted).getPosition().getNorm(); Assert.assertEquals(0, revertError, 180.0); Assert.assertEquals(2, adapterPropagator.propagate(t).getAdditionalState("dummy 1").length); Assert.assertEquals(1, adapterPropagator.propagate(t).getAdditionalState("dummy 2").length); Assert.assertEquals(3, adapterPropagator.propagate(t).getAdditionalState("dummy 3").length); } }
From source file:org.orekit.propagation.analytical.AdapterPropagatorTest.java
@Test public void testNonKeplerian() throws OrekitException, ParseException, IOException { Orbit leo = new CircularOrbit(7204319.233600575, 4.434564637450575E-4, 0.0011736728299091088, 1.7211611441767323, 5.5552084166959474, 24950.321259193086, PositionAngle.TRUE, FramesFactory.getEME2000(), new AbsoluteDate(new DateComponents(2003, 9, 16), new TimeComponents(23, 11, 20.264), TimeScalesFactory.getUTC()), Constants.EIGEN5C_EARTH_MU); double mass = 4093.0; AbsoluteDate t0 = new AbsoluteDate(new DateComponents(2003, 9, 16), new TimeComponents(23, 14, 40.264), TimeScalesFactory.getUTC()); Vector3D dV = new Vector3D(0.0, 3.0, 0.0); double f = 40.0; double isp = 300.0; double vExhaust = Constants.G0_STANDARD_GRAVITY * isp; double dt = -(mass * vExhaust / f) * FastMath.expm1(-dV.getNorm() / vExhaust); // setup a specific coefficient file for gravity potential as it will also // try to read a corrupted one otherwise GravityFieldFactory.addPotentialCoefficientsReader(new ICGEMFormatReader("g007_eigen_05c_coef", false)); NormalizedSphericalHarmonicsProvider gravityField = GravityFieldFactory.getNormalizedProvider(8, 8); BoundedPropagator withoutManeuver = getEphemeris(leo, mass, 10, new LofOffset(leo.getFrame(), LOFType.VNC), t0, Vector3D.ZERO, f, isp, true, true, gravityField); BoundedPropagator withManeuver = getEphemeris(leo, mass, 10, new LofOffset(leo.getFrame(), LOFType.VNC), t0, dV, f, isp, true, true, gravityField); // we set up a model that reverts the maneuvers AdapterPropagator adapterPropagator = new AdapterPropagator(withManeuver); SpacecraftState state0 = adapterPropagator.propagate(t0); AdapterPropagator.DifferentialEffect directEffect = new SmallManeuverAnalyticalModel(state0, dV.negate(), isp);//from w w w.j a va 2 s. c o m AdapterPropagator.DifferentialEffect derivedEffect = new J2DifferentialEffect(state0, directEffect, false, GravityFieldFactory.getUnnormalizedProvider(gravityField)); adapterPropagator.addEffect(directEffect); adapterPropagator.addEffect(derivedEffect); adapterPropagator.addAdditionalStateProvider(new AdditionalStateProvider() { public String getName() { return "dummy 3"; } public double[] getAdditionalState(SpacecraftState state) { return new double[3]; } }); // the adapted propagators do not manage the additional states from the reference, // they simply forward them Assert.assertFalse(adapterPropagator.isAdditionalStateManaged("dummy 1")); Assert.assertFalse(adapterPropagator.isAdditionalStateManaged("dummy 2")); Assert.assertTrue(adapterPropagator.isAdditionalStateManaged("dummy 3")); double maxDelta = 0; double maxNominal = 0; for (AbsoluteDate t = t0.shiftedBy(0.5 * dt); t.compareTo(withoutManeuver.getMaxDate()) < 0; t = t .shiftedBy(600.0)) { PVCoordinates pvWithout = withoutManeuver.getPVCoordinates(t, leo.getFrame()); PVCoordinates pvWith = withManeuver.getPVCoordinates(t, leo.getFrame()); PVCoordinates pvReverted = adapterPropagator.getPVCoordinates(t, leo.getFrame()); double nominal = new PVCoordinates(pvWithout, pvWith).getPosition().getNorm(); double revertError = new PVCoordinates(pvWithout, pvReverted).getPosition().getNorm(); maxDelta = FastMath.max(maxDelta, revertError); maxNominal = FastMath.max(maxNominal, nominal); Assert.assertEquals(2, adapterPropagator.propagate(t).getAdditionalState("dummy 1").length); Assert.assertEquals(1, adapterPropagator.propagate(t).getAdditionalState("dummy 2").length); Assert.assertEquals(3, adapterPropagator.propagate(t).getAdditionalState("dummy 3").length); } Assert.assertTrue(maxDelta < 120); Assert.assertTrue(maxNominal > 2800); }
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]);/* w w w.ja v a 2 s. 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; }
From source file:org.orekit.propagation.analytical.EcksteinHechlerPropagatorTest.java
@Test public void propagatedCartesian() throws OrekitException { // Definition of initial conditions with position and velocity // ------------------------------------------------------------ // with e around e = 1.4e-4 and i = 1.7 rad Vector3D position = new Vector3D(3220103., 69623., 6449822.); Vector3D velocity = new Vector3D(6414.7, -2006., -3180.); AbsoluteDate initDate = AbsoluteDate.J2000_EPOCH.shiftedBy(584.); Orbit initialOrbit = new EquinoctialOrbit(new PVCoordinates(position, velocity), FramesFactory.getEME2000(), initDate, provider.getMu()); // Extrapolator definition // ----------------------- EcksteinHechlerPropagator extrapolator = new EcksteinHechlerPropagator(initialOrbit, new LofOffset(initialOrbit.getFrame(), LOFType.VNC, RotationOrder.XYZ, 0, 0, 0), provider); // Extrapolation at a final date different from initial date // --------------------------------------------------------- double delta_t = 100000.0; // extrapolation duration in seconds AbsoluteDate extrapDate = initDate.shiftedBy(delta_t); SpacecraftState finalOrbit = extrapolator.propagate(extrapDate); Assert.assertEquals(0.0, finalOrbit.getDate().durationFrom(extrapDate), 1.0e-9); // computation of M final orbit double LM = finalOrbit.getLE() - finalOrbit.getEquinoctialEx() * FastMath.sin(finalOrbit.getLE()) + finalOrbit.getEquinoctialEy() * FastMath.cos(finalOrbit.getLE()); Assert.assertEquals(LM, finalOrbit.getLM(), Utils.epsilonAngle * FastMath.abs(finalOrbit.getLM())); // test of tan ((LE - Lv)/2) : Assert.assertEquals(FastMath.tan((finalOrbit.getLE() - finalOrbit.getLv()) / 2.), tangLEmLv(finalOrbit.getLv(), finalOrbit.getEquinoctialEx(), finalOrbit.getEquinoctialEy()), Utils.epsilonAngle);//from w w w .jav a 2 s . c o m // test of evolution of M vs E: LM = LE - ex*sin(LE) + ey*cos(LE) double deltaM = finalOrbit.getLM() - initialOrbit.getLM(); double deltaE = finalOrbit.getLE() - initialOrbit.getLE(); double delta = finalOrbit.getEquinoctialEx() * FastMath.sin(finalOrbit.getLE()) - initialOrbit.getEquinoctialEx() * FastMath.sin(initialOrbit.getLE()) - finalOrbit.getEquinoctialEy() * FastMath.cos(finalOrbit.getLE()) + initialOrbit.getEquinoctialEy() * FastMath.cos(initialOrbit.getLE()); Assert.assertEquals(deltaM, deltaE - delta, Utils.epsilonAngle * FastMath.abs(deltaE - delta)); // for final orbit double ex = finalOrbit.getEquinoctialEx(); double ey = finalOrbit.getEquinoctialEy(); double hx = finalOrbit.getHx(); double hy = finalOrbit.getHy(); double LE = finalOrbit.getLE(); double ex2 = ex * ex; double ey2 = ey * ey; double hx2 = hx * hx; double hy2 = hy * hy; double h2p1 = 1. + hx2 + hy2; double beta = 1. / (1. + FastMath.sqrt(1. - ex2 - ey2)); double x3 = -ex + (1. - beta * ey2) * FastMath.cos(LE) + beta * ex * ey * FastMath.sin(LE); double y3 = -ey + (1. - beta * ex2) * FastMath.sin(LE) + beta * ex * ey * FastMath.cos(LE); Vector3D U = new Vector3D((1. + hx2 - hy2) / h2p1, (2. * hx * hy) / h2p1, (-2. * hy) / h2p1); Vector3D V = new Vector3D((2. * hx * hy) / h2p1, (1. - hx2 + hy2) / h2p1, (2. * hx) / h2p1); Vector3D r = new Vector3D(finalOrbit.getA(), (new Vector3D(x3, U, y3, V))); Assert.assertEquals(finalOrbit.getPVCoordinates().getPosition().getNorm(), r.getNorm(), Utils.epsilonTest * r.getNorm()); }
From source file:org.orekit.propagation.analytical.EcksteinHechlerPropagatorTest.java
@Test public void propagatedKeplerian() throws OrekitException { // Definition of initial conditions with keplerian parameters // ----------------------------------------------------------- AbsoluteDate initDate = AbsoluteDate.J2000_EPOCH.shiftedBy(584.); Orbit initialOrbit = new KeplerianOrbit(7209668.0, 0.5e-4, 1.7, 2.1, 2.9, 6.2, PositionAngle.TRUE, FramesFactory.getEME2000(), initDate, provider.getMu()); // Extrapolator definition // ----------------------- EcksteinHechlerPropagator extrapolator = new EcksteinHechlerPropagator(initialOrbit, new LofOffset(initialOrbit.getFrame(), LOFType.VNC, RotationOrder.XYZ, 0, 0, 0), 2000.0, provider); // Extrapolation at a final date different from initial date // --------------------------------------------------------- double delta_t = 100000.0; // extrapolation duration in seconds AbsoluteDate extrapDate = initDate.shiftedBy(delta_t); SpacecraftState finalOrbit = extrapolator.propagate(extrapDate); Assert.assertEquals(0.0, finalOrbit.getDate().durationFrom(extrapDate), 1.0e-9); // computation of M final orbit double LM = finalOrbit.getLE() - finalOrbit.getEquinoctialEx() * FastMath.sin(finalOrbit.getLE()) + finalOrbit.getEquinoctialEy() * FastMath.cos(finalOrbit.getLE()); Assert.assertEquals(LM, finalOrbit.getLM(), Utils.epsilonAngle); // test of tan((LE - Lv)/2) : Assert.assertEquals(FastMath.tan((finalOrbit.getLE() - finalOrbit.getLv()) / 2.), tangLEmLv(finalOrbit.getLv(), finalOrbit.getEquinoctialEx(), finalOrbit.getEquinoctialEy()), Utils.epsilonAngle);//from w w w .j av a 2 s. c om // test of evolution of M vs E: LM = LE - ex*sin(LE) + ey*cos(LE) // with ex and ey the same for initial and final orbit double deltaM = finalOrbit.getLM() - initialOrbit.getLM(); double deltaE = finalOrbit.getLE() - initialOrbit.getLE(); double delta = finalOrbit.getEquinoctialEx() * FastMath.sin(finalOrbit.getLE()) - initialOrbit.getEquinoctialEx() * FastMath.sin(initialOrbit.getLE()) - finalOrbit.getEquinoctialEy() * FastMath.cos(finalOrbit.getLE()) + initialOrbit.getEquinoctialEy() * FastMath.cos(initialOrbit.getLE()); Assert.assertEquals(deltaM, deltaE - delta, Utils.epsilonAngle * FastMath.abs(deltaE - delta)); // for final orbit double ex = finalOrbit.getEquinoctialEx(); double ey = finalOrbit.getEquinoctialEy(); double hx = finalOrbit.getHx(); double hy = finalOrbit.getHy(); double LE = finalOrbit.getLE(); double ex2 = ex * ex; double ey2 = ey * ey; double hx2 = hx * hx; double hy2 = hy * hy; double h2p1 = 1. + hx2 + hy2; double beta = 1. / (1. + FastMath.sqrt(1. - ex2 - ey2)); double x3 = -ex + (1. - beta * ey2) * FastMath.cos(LE) + beta * ex * ey * FastMath.sin(LE); double y3 = -ey + (1. - beta * ex2) * FastMath.sin(LE) + beta * ex * ey * FastMath.cos(LE); Vector3D U = new Vector3D((1. + hx2 - hy2) / h2p1, (2. * hx * hy) / h2p1, (-2. * hy) / h2p1); Vector3D V = new Vector3D((2. * hx * hy) / h2p1, (1. - hx2 + hy2) / h2p1, (2. * hx) / h2p1); Vector3D r = new Vector3D(finalOrbit.getA(), (new Vector3D(x3, U, y3, V))); Assert.assertEquals(finalOrbit.getPVCoordinates().getPosition().getNorm(), r.getNorm(), Utils.epsilonTest * r.getNorm()); }
From source file:org.orekit.propagation.analytical.KeplerianPropagatorTest.java
@Test public void propagatedCartesian() throws OrekitException { // Definition of initial conditions with position and velocity //------------------------------------------------------------ Vector3D position = new Vector3D(7.0e6, 1.0e6, 4.0e6); Vector3D velocity = new Vector3D(-500.0, 8000.0, 1000.0); double mu = 3.9860047e14; AbsoluteDate initDate = AbsoluteDate.J2000_EPOCH.shiftedBy(584.); Orbit initialOrbit = new EquinoctialOrbit(new PVCoordinates(position, velocity), FramesFactory.getEME2000(), initDate, mu);//from w ww . j av a 2 s .c o m // Extrapolator definition // ----------------------- KeplerianPropagator extrapolator = new KeplerianPropagator(initialOrbit); // Extrapolation at a final date different from initial date // --------------------------------------------------------- double delta_t = 100000.0; // extrapolation duration in seconds AbsoluteDate extrapDate = initDate.shiftedBy(delta_t); SpacecraftState finalOrbit = extrapolator.propagate(extrapDate); // computation of (M final - M initial) with another method double a = finalOrbit.getA(); // another way to compute n double n = FastMath.sqrt(finalOrbit.getMu() / FastMath.pow(a, 3)); Assert.assertEquals(n * delta_t, finalOrbit.getLM() - initialOrbit.getLM(), Utils.epsilonAngle); // computation of M final orbit double LM = finalOrbit.getLE() - finalOrbit.getEquinoctialEx() * FastMath.sin(finalOrbit.getLE()) + finalOrbit.getEquinoctialEy() * FastMath.cos(finalOrbit.getLE()); Assert.assertEquals(LM, finalOrbit.getLM(), Utils.epsilonAngle); // test of tan ((LE - Lv)/2) : Assert.assertEquals(FastMath.tan((finalOrbit.getLE() - finalOrbit.getLv()) / 2.), tangLEmLv(finalOrbit.getLv(), finalOrbit.getEquinoctialEx(), finalOrbit.getEquinoctialEy()), Utils.epsilonAngle); // test of evolution of M vs E: LM = LE - ex*sin(LE) + ey*cos(LE) // with ex and ey the same for initial and final orbit double deltaM = finalOrbit.getLM() - initialOrbit.getLM(); double deltaE = finalOrbit.getLE() - initialOrbit.getLE(); double delta = finalOrbit.getEquinoctialEx() * (FastMath.sin(finalOrbit.getLE()) - FastMath.sin(initialOrbit.getLE())) - finalOrbit.getEquinoctialEy() * (FastMath.cos(finalOrbit.getLE()) - FastMath.cos(initialOrbit.getLE())); Assert.assertEquals(deltaM, deltaE - delta, Utils.epsilonAngle); // the orbital elements except for Mean/True/Eccentric latitude arguments are the same Assert.assertEquals(finalOrbit.getA(), initialOrbit.getA(), Utils.epsilonTest * initialOrbit.getA()); Assert.assertEquals(finalOrbit.getEquinoctialEx(), initialOrbit.getEquinoctialEx(), Utils.epsilonE); Assert.assertEquals(finalOrbit.getEquinoctialEy(), initialOrbit.getEquinoctialEy(), Utils.epsilonE); Assert.assertEquals(finalOrbit.getHx(), initialOrbit.getHx(), Utils.epsilonAngle); Assert.assertEquals(finalOrbit.getHy(), initialOrbit.getHy(), Utils.epsilonAngle); // for final orbit double ex = finalOrbit.getEquinoctialEx(); double ey = finalOrbit.getEquinoctialEy(); double hx = finalOrbit.getHx(); double hy = finalOrbit.getHy(); double LE = finalOrbit.getLE(); double ex2 = ex * ex; double ey2 = ey * ey; double hx2 = hx * hx; double hy2 = hy * hy; double h2p1 = 1. + hx2 + hy2; double beta = 1. / (1. + FastMath.sqrt(1. - ex2 - ey2)); double x3 = -ex + (1. - beta * ey2) * FastMath.cos(LE) + beta * ex * ey * FastMath.sin(LE); double y3 = -ey + (1. - beta * ex2) * FastMath.sin(LE) + beta * ex * ey * FastMath.cos(LE); Vector3D U = new Vector3D((1. + hx2 - hy2) / h2p1, (2. * hx * hy) / h2p1, (-2. * hy) / h2p1); Vector3D V = new Vector3D((2. * hx * hy) / h2p1, (1. - hx2 + hy2) / h2p1, (2. * hx) / h2p1); Vector3D r = new Vector3D(finalOrbit.getA(), (new Vector3D(x3, U, y3, V))); Assert.assertEquals(finalOrbit.getPVCoordinates().getPosition().getNorm(), r.getNorm(), Utils.epsilonTest * r.getNorm()); }
From source file:org.orekit.propagation.analytical.KeplerianPropagatorTest.java
@Test public void propagatedKeplerian() throws OrekitException { // Definition of initial conditions with keplerian parameters //----------------------------------------------------------- AbsoluteDate initDate = AbsoluteDate.J2000_EPOCH.shiftedBy(584.); Orbit initialOrbit = new KeplerianOrbit(7209668.0, 0.5e-4, 1.7, 2.1, 2.9, 6.2, PositionAngle.TRUE, FramesFactory.getEME2000(), initDate, mu); // Extrapolator definition // ----------------------- KeplerianPropagator extrapolator = new KeplerianPropagator(initialOrbit); // Extrapolation at a final date different from initial date // --------------------------------------------------------- double delta_t = 100000.0; // extrapolation duration in seconds AbsoluteDate extrapDate = initDate.shiftedBy(delta_t); SpacecraftState finalOrbit = extrapolator.propagate(extrapDate); Assert.assertEquals(6092.3362422560844633, finalOrbit.getKeplerianPeriod(), 1.0e-12); Assert.assertEquals(0.001031326088602888358, finalOrbit.getKeplerianMeanMotion(), 1.0e-16); // computation of (M final - M initial) with another method double a = finalOrbit.getA(); // another way to compute n double n = FastMath.sqrt(finalOrbit.getMu() / FastMath.pow(a, 3)); Assert.assertEquals(n * delta_t, finalOrbit.getLM() - initialOrbit.getLM(), Utils.epsilonAngle); // computation of M final orbit double LM = finalOrbit.getLE() - finalOrbit.getEquinoctialEx() * FastMath.sin(finalOrbit.getLE()) + finalOrbit.getEquinoctialEy() * FastMath.cos(finalOrbit.getLE()); Assert.assertEquals(LM, finalOrbit.getLM(), Utils.epsilonAngle); // test of tan ((LE - Lv)/2) : Assert.assertEquals(FastMath.tan((finalOrbit.getLE() - finalOrbit.getLv()) / 2.), tangLEmLv(finalOrbit.getLv(), finalOrbit.getEquinoctialEx(), finalOrbit.getEquinoctialEy()), Utils.epsilonAngle);/*from w w w . j a v a 2 s .co m*/ // test of evolution of M vs E: LM = LE - ex*sin(LE) + ey*cos(LE) // with ex and ey the same for initial and final orbit double deltaM = finalOrbit.getLM() - initialOrbit.getLM(); double deltaE = finalOrbit.getLE() - initialOrbit.getLE(); double delta = finalOrbit.getEquinoctialEx() * (FastMath.sin(finalOrbit.getLE()) - FastMath.sin(initialOrbit.getLE())) - finalOrbit.getEquinoctialEy() * (FastMath.cos(finalOrbit.getLE()) - FastMath.cos(initialOrbit.getLE())); Assert.assertEquals(deltaM, deltaE - delta, Utils.epsilonAngle); // the orbital elements except for Mean/True/Eccentric latitude arguments are the same Assert.assertEquals(finalOrbit.getA(), initialOrbit.getA(), Utils.epsilonTest * initialOrbit.getA()); Assert.assertEquals(finalOrbit.getEquinoctialEx(), initialOrbit.getEquinoctialEx(), Utils.epsilonE); Assert.assertEquals(finalOrbit.getEquinoctialEy(), initialOrbit.getEquinoctialEy(), Utils.epsilonE); Assert.assertEquals(finalOrbit.getHx(), initialOrbit.getHx(), Utils.epsilonAngle); Assert.assertEquals(finalOrbit.getHy(), initialOrbit.getHy(), Utils.epsilonAngle); // for final orbit double ex = finalOrbit.getEquinoctialEx(); double ey = finalOrbit.getEquinoctialEy(); double hx = finalOrbit.getHx(); double hy = finalOrbit.getHy(); double LE = finalOrbit.getLE(); double ex2 = ex * ex; double ey2 = ey * ey; double hx2 = hx * hx; double hy2 = hy * hy; double h2p1 = 1. + hx2 + hy2; double beta = 1. / (1. + FastMath.sqrt(1. - ex2 - ey2)); double x3 = -ex + (1. - beta * ey2) * FastMath.cos(LE) + beta * ex * ey * FastMath.sin(LE); double y3 = -ey + (1. - beta * ex2) * FastMath.sin(LE) + beta * ex * ey * FastMath.cos(LE); Vector3D U = new Vector3D((1. + hx2 - hy2) / h2p1, (2. * hx * hy) / h2p1, (-2. * hy) / h2p1); Vector3D V = new Vector3D((2. * hx * hy) / h2p1, (1. - hx2 + hy2) / h2p1, (2. * hx) / h2p1); Vector3D r = new Vector3D(finalOrbit.getA(), (new Vector3D(x3, U, y3, V))); Assert.assertEquals(finalOrbit.getPVCoordinates().getPosition().getNorm(), r.getNorm(), Utils.epsilonTest * r.getNorm()); }
From source file:org.orekit.propagation.events.EclipseDetector.java
/** Compute the value of the switching function. * This function becomes negative when entering the region of shadow * and positive when exiting./* w w w . j a va 2s. com*/ * @param s the current state information: date, kinematics, attitude * @return value of the switching function * @exception OrekitException if some specific error occurs */ public double g(final SpacecraftState s) throws OrekitException { final Vector3D pted = occulted.getPVCoordinates(s.getDate(), s.getFrame()).getPosition(); final Vector3D ping = occulting.getPVCoordinates(s.getDate(), s.getFrame()).getPosition(); final Vector3D psat = s.getPVCoordinates().getPosition(); final Vector3D ps = pted.subtract(psat); final Vector3D po = ping.subtract(psat); final double angle = Vector3D.angle(ps, po); final double rs = FastMath.asin(occultedRadius / ps.getNorm()); if (Double.isNaN(rs)) { return FastMath.PI; } final double ro = FastMath.asin(occultingRadius / po.getNorm()); if (Double.isNaN(ro)) { return -FastMath.PI; } return totalEclipse ? (angle - ro + rs) : (angle - ro - rs); }
From source file:org.orekit.propagation.numerical.Jacobianizer.java
/** Compute acceleration and derivatives with respect to state. * @param date current date/*from w ww . j ava2 s .c om*/ * @param frame inertial reference frame for state (both orbit and attitude) * @param position position of spacecraft in reference frame * @param velocity velocity of spacecraft in reference frame * @param rotation orientation (attitude) of the spacecraft with respect to reference frame * @param mass spacecraft mass * @return acceleration with derivatives * @exception OrekitException if the underlying force models cannot compute the acceleration */ public FieldVector3D<DerivativeStructure> accelerationDerivatives(final AbsoluteDate date, final Frame frame, final FieldVector3D<DerivativeStructure> position, final FieldVector3D<DerivativeStructure> velocity, final FieldRotation<DerivativeStructure> rotation, final DerivativeStructure mass) throws OrekitException { final int parameters = mass.getFreeParameters(); final int order = mass.getOrder(); // estimate the scalar velocity step, assuming energy conservation // and hence differentiating equation V = sqrt(mu (2/r - 1/a)) final Vector3D p0 = position.toVector3D(); final Vector3D v0 = velocity.toVector3D(); final double r2 = p0.getNormSq(); final double hVel = mu * hPos / (v0.getNorm() * r2); // estimate mass step, applying the same relative value as position final double hMass = mass.getValue() * hPos / FastMath.sqrt(r2); // compute nominal acceleration final AccelerationRetriever nominal = new AccelerationRetriever(); computeShiftedAcceleration(nominal, date, frame, p0, v0, rotation.toRotation(), mass.getValue()); final double[] a0 = nominal.getAcceleration().toArray(); // compute accelerations with shifted states final AccelerationRetriever shifted = new AccelerationRetriever(); // shift position by hPos alon x, y and z computeShiftedAcceleration(shifted, date, frame, shift(position, 0, hPos), v0, shift(rotation, 0, hPos), shift(mass, 0, hPos)); final double[] derPx = new Vector3D(1 / hPos, shifted.getAcceleration(), -1 / hPos, nominal.getAcceleration()).toArray(); computeShiftedAcceleration(shifted, date, frame, shift(position, 1, hPos), v0, shift(rotation, 1, hPos), shift(mass, 1, hPos)); final double[] derPy = new Vector3D(1 / hPos, shifted.getAcceleration(), -1 / hPos, nominal.getAcceleration()).toArray(); computeShiftedAcceleration(shifted, date, frame, shift(position, 2, hPos), v0, shift(rotation, 2, hPos), shift(mass, 2, hPos)); final double[] derPz = new Vector3D(1 / hPos, shifted.getAcceleration(), -1 / hPos, nominal.getAcceleration()).toArray(); // shift velocity by hVel alon x, y and z computeShiftedAcceleration(shifted, date, frame, p0, shift(velocity, 3, hVel), shift(rotation, 3, hVel), shift(mass, 3, hVel)); final double[] derVx = new Vector3D(1 / hVel, shifted.getAcceleration(), -1 / hVel, nominal.getAcceleration()).toArray(); computeShiftedAcceleration(shifted, date, frame, p0, shift(velocity, 4, hVel), shift(rotation, 4, hVel), shift(mass, 4, hVel)); final double[] derVy = new Vector3D(1 / hVel, shifted.getAcceleration(), -1 / hVel, nominal.getAcceleration()).toArray(); computeShiftedAcceleration(shifted, date, frame, p0, shift(velocity, 5, hVel), shift(rotation, 5, hVel), shift(mass, 5, hVel)); final double[] derVz = new Vector3D(1 / hVel, shifted.getAcceleration(), -1 / hVel, nominal.getAcceleration()).toArray(); final double[] derM; if (parameters < 7) { derM = null; } else { // shift mass by hMass computeShiftedAcceleration(shifted, date, frame, p0, v0, shift(rotation, 6, hMass), shift(mass, 6, hMass)); derM = new Vector3D(1 / hMass, shifted.getAcceleration(), -1 / hMass, nominal.getAcceleration()) .toArray(); } final double[] derivatives = new double[1 + parameters]; final DerivativeStructure[] accDer = new DerivativeStructure[3]; for (int i = 0; i < 3; ++i) { // first element is value of acceleration derivatives[0] = a0[i]; // next three elements are derivatives with respect to position derivatives[1] = derPx[i]; derivatives[2] = derPy[i]; derivatives[3] = derPz[i]; // next three elements are derivatives with respect to velocity derivatives[4] = derVx[i]; derivatives[5] = derVy[i]; derivatives[6] = derVz[i]; if (derM != null) { derivatives[7] = derM[i]; } accDer[i] = new DerivativeStructure(parameters, order, derivatives); } return new FieldVector3D<DerivativeStructure>(accDer); }