Examples with SpacecraftState org.orekit.propagation.SpacecraftState used on opensource projects

Example 86 with SpacecraftState

use of org.orekit.propagation.SpacecraftState in project Orekit by CS-SI.

the class DSSTPropagatorTest method testPropagationWithSolarRadiationPressure.

@Test
public void testPropagationWithSolarRadiationPressure() throws OrekitException {
    // Central Body geopotential 2x0
    final UnnormalizedSphericalHarmonicsProvider provider = GravityFieldFactory.getUnnormalizedProvider(2, 0);
    DSSTForceModel zonal = new DSSTZonal(provider, 2, 1, 5);
    DSSTForceModel tesseral = new DSSTTesseral(CelestialBodyFactory.getEarth().getBodyOrientedFrame(), Constants.WGS84_EARTH_ANGULAR_VELOCITY, provider, 2, 0, 0, 2, 2, 0, 0);
    // SRP Force Model
    DSSTForceModel srp = new DSSTSolarRadiationPressure(1.2, 100., CelestialBodyFactory.getSun(), Constants.WGS84_EARTH_EQUATORIAL_RADIUS);
    // GEO Orbit
    final AbsoluteDate initDate = new AbsoluteDate(2003, 9, 16, 0, 0, 00.000, TimeScalesFactory.getUTC());
    final Orbit orbit = new KeplerianOrbit(42166258., 0.0001, FastMath.toRadians(0.001), FastMath.toRadians(315.4985), FastMath.toRadians(130.7562), FastMath.toRadians(44.2377), PositionAngle.MEAN, FramesFactory.getGCRF(), initDate, provider.getMu());
    // Set propagator with state and force model
    dsstProp = new DSSTPropagator(new ClassicalRungeKuttaIntegrator(86400.));
    dsstProp.setInitialState(new SpacecraftState(orbit), false);
    dsstProp.addForceModel(zonal);
    dsstProp.addForceModel(tesseral);
    dsstProp.addForceModel(srp);
    // 10 days propagation
    final SpacecraftState state = dsstProp.propagate(initDate.shiftedBy(10. * 86400.));
    // Ref Standalone_DSST:
    // a    = 42166257.99807995 m
    // h/ey = -0.1191876027555493D-03
    // k/ex = -0.1781865038201885D-05
    // p/hy =  0.6618387121369373D-05
    // q/hx = -0.5624363171289686D-05
    // lM   = 140°3496229467104
    Assert.assertEquals(42166257.99807995, state.getA(), 0.8);
    Assert.assertEquals(-0.1781865038201885e-05, state.getEquinoctialEx(), 3.e-7);
    Assert.assertEquals(-0.1191876027555493e-03, state.getEquinoctialEy(), 4.e-6);
    Assert.assertEquals(-0.5624363171289686e-05, state.getHx(), 4.e-9);
    Assert.assertEquals(0.6618387121369373e-05, state.getHy(), 3.e-10);
    Assert.assertEquals(140.3496229467104, FastMath.toDegrees(MathUtils.normalizeAngle(state.getLM(), FastMath.PI)), 2.e-4);
}

Example 87 with SpacecraftState

use of org.orekit.propagation.SpacecraftState 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);
}

Example 88 with SpacecraftState

use of org.orekit.propagation.SpacecraftState in project Orekit by CS-SI.

the class Phasing method printGridPoints.

/**
 * Print ground track grid point
 * @param out output stream
 * @param latitude point latitude
 * @param ascending indicator for latitude crossing direction
 * @param orbit phased orbit
 * @param propagator propagator for orbit
 * @param nbOrbits number of orbits in the cycle
 * @exception OrekitException if orbit cannot be propagated
 */
private void printGridPoints(final PrintStream out, final double latitude, final boolean ascending, final Orbit orbit, final Propagator propagator, int nbOrbits) throws OrekitException {
    propagator.resetInitialState(new SpacecraftState(orbit));
    AbsoluteDate start = orbit.getDate();
    // find the first latitude crossing
    double period = orbit.getKeplerianPeriod();
    double stepSize = period / 100;
    SpacecraftState crossing = findFirstCrossing(latitude, ascending, start, start.shiftedBy(2 * period), stepSize, propagator);
    // find all other latitude crossings from regular schedule
    DecimalFormat fTime = new DecimalFormat("0000000.000", new DecimalFormatSymbols(Locale.US));
    DecimalFormat fAngle = new DecimalFormat("000.00000", new DecimalFormatSymbols(Locale.US));
    while (nbOrbits-- > 0) {
        GeodeticPoint gp = earth.transform(crossing.getPVCoordinates().getPosition(), crossing.getFrame(), crossing.getDate());
        out.println(fTime.format(crossing.getDate().durationFrom(start)) + " " + fAngle.format(FastMath.toDegrees(gp.getLatitude())) + " " + fAngle.format(FastMath.toDegrees(gp.getLongitude())) + " " + ascending);
        final AbsoluteDate previousDate = crossing.getDate();
        crossing = findLatitudeCrossing(latitude, previousDate.shiftedBy(period), previousDate.shiftedBy(2 * period), stepSize, period / 8, propagator);
        period = crossing.getDate().durationFrom(previousDate);
    }
}

Example 89 with SpacecraftState

use of org.orekit.propagation.SpacecraftState in project Orekit by CS-SI.

the class Phasing method improveEarthPhasing.

/**
 * Improve orbit to better match Earth phasing parameters.
 * @param previous previous orbit
 * @param nbOrbits number of orbits in the phasing cycle
 * @param nbDays number of days in the phasing cycle
 * @param propagator propagator to use
 * @return an improved Earth phased orbit
 * @exception OrekitException if orbit cannot be propagated
 */
private CircularOrbit improveEarthPhasing(CircularOrbit previous, int nbOrbits, int nbDays, Propagator propagator) throws OrekitException {
    propagator.resetInitialState(new SpacecraftState(previous));
    // find first ascending node
    double period = previous.getKeplerianPeriod();
    SpacecraftState firstState = findFirstCrossing(0.0, true, previous.getDate(), previous.getDate().shiftedBy(2 * period), 0.01 * period, propagator);
    // go to next cycle, one orbit at a time
    SpacecraftState state = firstState;
    for (int i = 0; i < nbOrbits; ++i) {
        final AbsoluteDate previousDate = state.getDate();
        state = findLatitudeCrossing(0.0, previousDate.shiftedBy(period), previousDate.shiftedBy(2 * period), 0.01 * period, period, propagator);
        period = state.getDate().durationFrom(previousDate);
    }
    double cycleDuration = state.getDate().durationFrom(firstState.getDate());
    double deltaT;
    if (((int) FastMath.rint(cycleDuration / Constants.JULIAN_DAY)) != nbDays) {
        // we are very far from expected duration
        deltaT = nbDays * Constants.JULIAN_DAY - cycleDuration;
    } else {
        // we are close to expected duration
        GeodeticPoint startPoint = earth.transform(firstState.getPVCoordinates().getPosition(), firstState.getFrame(), firstState.getDate());
        GeodeticPoint endPoint = earth.transform(state.getPVCoordinates().getPosition(), state.getFrame(), state.getDate());
        double deltaL = MathUtils.normalizeAngle(endPoint.getLongitude() - startPoint.getLongitude(), 0.0);
        deltaT = deltaL * Constants.JULIAN_DAY / (2 * FastMath.PI);
    }
    double deltaA = 2 * previous.getA() * deltaT / (3 * nbOrbits * previous.getKeplerianPeriod());
    return new CircularOrbit(previous.getA() + deltaA, previous.getCircularEx(), previous.getCircularEy(), previous.getI(), previous.getRightAscensionOfAscendingNode(), previous.getAlphaV(), PositionAngle.TRUE, previous.getFrame(), previous.getDate(), previous.getMu());
}

Example 90 with SpacecraftState

use of org.orekit.propagation.SpacecraftState in project Orekit by CS-SI.

the class Phasing method improveSunSynchronization.

/**
 * Improve orbit to better match phasing parameters.
 * @param previous previous orbit
 * @param duration sampling duration
 * @param latitude reference latitude for Sun synchronous orbit
 * @param ascending if true, crossing latitude is from South to North
 * @param mst desired mean solar time at reference latitude crossing
 * @param propagator propagator to use
 * @return an improved Earth phased, Sun synchronous orbit
 * @exception OrekitException if orbit cannot be propagated
 */
private CircularOrbit improveSunSynchronization(CircularOrbit previous, double duration, double latitude, boolean ascending, double mst, Propagator propagator) throws OrekitException {
    propagator.resetInitialState(new SpacecraftState(previous));
    AbsoluteDate start = previous.getDate();
    // find the first latitude crossing
    double period = previous.getKeplerianPeriod();
    double stepSize = period / 100;
    SpacecraftState crossing = findFirstCrossing(latitude, ascending, start, start.shiftedBy(2 * period), stepSize, propagator);
    // find all other latitude crossings from regular schedule
    SecularAndHarmonic mstModel = new SecularAndHarmonic(2, 2.0 * FastMath.PI / Constants.JULIAN_YEAR, 4.0 * FastMath.PI / Constants.JULIAN_YEAR, 2.0 * FastMath.PI / Constants.JULIAN_DAY, 4.0 * FastMath.PI / Constants.JULIAN_DAY);
    mstModel.resetFitting(start, new double[] { mst, -1.0e-10, -1.0e-17, 1.0e-3, 1.0e-3, 1.0e-5, 1.0e-5, 1.0e-5, 1.0e-5, 1.0e-5, 1.0e-5 });
    while (crossing != null && crossing.getDate().durationFrom(start) < duration) {
        final AbsoluteDate previousDate = crossing.getDate();
        crossing = findLatitudeCrossing(latitude, previousDate.shiftedBy(period), previousDate.shiftedBy(2 * period), stepSize, period / 8, propagator);
        if (crossing != null) {
            // store current point
            mstModel.addPoint(crossing.getDate(), meanSolarTime(crossing.getOrbit()));
            // use the same time separation to pinpoint next crossing
            period = crossing.getDate().durationFrom(previousDate);
        }
    }
    // fit the mean solar time to a parabolic plus medium periods model
    // we will only use the linear part for the correction
    mstModel.fit();
    final double[] fittedH = mstModel.approximateAsPolynomialOnly(1, start, 2, 2, start, start.shiftedBy(duration), stepSize);
    // solar time bias must be compensated by shifting ascending node
    double deltaRaan = FastMath.PI * (mst - fittedH[0]) / 12;
    // solar time slope must be compensated by changing inclination
    // linearized relationship between hDot and inclination:
    // hDot = alphaDot - raanDot where alphaDot is the angular rate of Sun right ascension
    // and raanDot is the angular rate of ascending node right ascension. So hDot evolution
    // is the opposite of raan evolution, which itself is proportional to cos(i) due to J2
    // effect. So hDot = alphaDot - k cos(i) and hence Delta hDot = -k sin(i) Delta i
    // so Delta hDot / Delta i = (alphaDot - hDot) tan(i)
    double dhDotDi = (24.0 / Constants.JULIAN_YEAR - fittedH[1]) * FastMath.tan(previous.getI());
    // compute inclination offset needed to achieve station-keeping target
    final double deltaI = fittedH[1] / dhDotDi;
    return new CircularOrbit(previous.getA(), previous.getCircularEx(), previous.getCircularEy(), previous.getI() + deltaI, previous.getRightAscensionOfAscendingNode() + deltaRaan, previous.getAlphaV(), PositionAngle.TRUE, previous.getFrame(), previous.getDate(), previous.getMu());
}