Examples with CircularOrbit org.orekit.orbits.CircularOrbit used on opensource projects

Example 56 with CircularOrbit

use of org.orekit.orbits.CircularOrbit in project Orekit by CS-SI.

the class EllipsoidTessellatorTest method setUp.

@Before
public void setUp() throws OrekitException {
    Utils.setDataRoot("regular-data");
    // the following orbital parameters have been computed using
    // Orekit tutorial about phasing, using the following configuration:
    // 
    // orbit.date                          = 2012-01-01T00:00:00.000
    // phasing.orbits.number               = 143
    // phasing.days.number                 =  10
    // sun.synchronous.reference.latitude  = 0
    // sun.synchronous.reference.ascending = false
    // sun.synchronous.mean.solar.time     = 10:30:00
    // gravity.field.degree                = 12
    // gravity.field.order                 = 12
    AbsoluteDate date = new AbsoluteDate("2012-01-01T00:00:00.000", TimeScalesFactory.getUTC());
    Frame eme2000 = FramesFactory.getEME2000();
    orbit = new CircularOrbit(7173352.811913891, -4.029194321683225E-4, 0.0013530362644647786, FastMath.toRadians(98.63218182243709), FastMath.toRadians(77.55565567747836), FastMath.PI, PositionAngle.TRUE, eme2000, date, Constants.EIGEN5C_EARTH_MU);
    ellipsoid = new OneAxisEllipsoid(Constants.WGS84_EARTH_EQUATORIAL_RADIUS, Constants.WGS84_EARTH_FLATTENING, FramesFactory.getITRF(IERSConventions.IERS_2010, true));
}

Example 57 with CircularOrbit

use of org.orekit.orbits.CircularOrbit in project Orekit by CS-SI.

the class Phasing method run.

private void run(final File input) throws IOException, IllegalArgumentException, ParseException, OrekitException {
    // read input parameters
    KeyValueFileParser<ParameterKey> parser = new KeyValueFileParser<ParameterKey>(ParameterKey.class);
    try (final FileInputStream fis = new FileInputStream(input)) {
        parser.parseInput(input.getAbsolutePath(), fis);
    }
    TimeScale utc = TimeScalesFactory.getUTC();
    // simulation properties
    AbsoluteDate date = parser.getDate(ParameterKey.ORBIT_DATE, utc);
    int nbOrbits = parser.getInt(ParameterKey.PHASING_ORBITS_NUMBER);
    int nbDays = parser.getInt(ParameterKey.PHASING_DAYS_NUMBER);
    double latitude = parser.getAngle(ParameterKey.SUN_SYNCHRONOUS_REFERENCE_LATITUDE);
    boolean ascending = parser.getBoolean(ParameterKey.SUN_SYNCHRONOUS_REFERENCE_ASCENDING);
    double mst = parser.getTime(ParameterKey.SUN_SYNCHRONOUS_MEAN_SOLAR_TIME).getSecondsInUTCDay() / 3600;
    int degree = parser.getInt(ParameterKey.GRAVITY_FIELD_DEGREE);
    int order = parser.getInt(ParameterKey.GRAVITY_FIELD_ORDER);
    String gridOutput = parser.getString(ParameterKey.GRID_OUTPUT);
    double[] gridLatitudes = new double[] { parser.getAngle(ParameterKey.GRID_LATITUDE_1), parser.getAngle(ParameterKey.GRID_LATITUDE_2), parser.getAngle(ParameterKey.GRID_LATITUDE_3), parser.getAngle(ParameterKey.GRID_LATITUDE_4), parser.getAngle(ParameterKey.GRID_LATITUDE_5) };
    boolean[] gridAscending = new boolean[] { parser.getBoolean(ParameterKey.GRID_ASCENDING_1), parser.getBoolean(ParameterKey.GRID_ASCENDING_2), parser.getBoolean(ParameterKey.GRID_ASCENDING_3), parser.getBoolean(ParameterKey.GRID_ASCENDING_4), parser.getBoolean(ParameterKey.GRID_ASCENDING_5) };
    gravityField = GravityFieldFactory.getNormalizedProvider(degree, order);
    // initial guess for orbit
    CircularOrbit orbit = guessOrbit(date, FramesFactory.getEME2000(), nbOrbits, nbDays, latitude, ascending, mst);
    System.out.println("initial orbit: " + orbit);
    System.out.println("please wait while orbit is adjusted...");
    System.out.println();
    // numerical model for improving orbit
    double[][] tolerances = NumericalPropagator.tolerances(0.1, orbit, OrbitType.CIRCULAR);
    DormandPrince853Integrator integrator = new DormandPrince853Integrator(1.0e-4 * orbit.getKeplerianPeriod(), 1.0e-1 * orbit.getKeplerianPeriod(), tolerances[0], tolerances[1]);
    integrator.setInitialStepSize(1.0e-2 * orbit.getKeplerianPeriod());
    NumericalPropagator propagator = new NumericalPropagator(integrator);
    propagator.addForceModel(new HolmesFeatherstoneAttractionModel(FramesFactory.getGTOD(IERSConventions.IERS_2010, true), gravityField));
    propagator.addForceModel(new ThirdBodyAttraction(CelestialBodyFactory.getSun()));
    propagator.addForceModel(new ThirdBodyAttraction(CelestialBodyFactory.getMoon()));
    double deltaP = Double.POSITIVE_INFINITY;
    double deltaV = Double.POSITIVE_INFINITY;
    int counter = 0;
    DecimalFormat f = new DecimalFormat("0.000E00", new DecimalFormatSymbols(Locale.US));
    while (deltaP > 3.0e-1 || deltaV > 3.0e-4) {
        CircularOrbit previous = orbit;
        CircularOrbit tmp1 = improveEarthPhasing(previous, nbOrbits, nbDays, propagator);
        CircularOrbit tmp2 = improveSunSynchronization(tmp1, nbOrbits * tmp1.getKeplerianPeriod(), latitude, ascending, mst, propagator);
        orbit = improveFrozenEccentricity(tmp2, nbOrbits * tmp2.getKeplerianPeriod(), propagator);
        double da = orbit.getA() - previous.getA();
        double dex = orbit.getCircularEx() - previous.getCircularEx();
        double dey = orbit.getCircularEy() - previous.getCircularEy();
        double di = FastMath.toDegrees(orbit.getI() - previous.getI());
        double dr = FastMath.toDegrees(orbit.getRightAscensionOfAscendingNode() - previous.getRightAscensionOfAscendingNode());
        System.out.println(" iteration " + (++counter) + ": deltaA = " + f.format(da) + " m, deltaEx = " + f.format(dex) + ", deltaEy = " + f.format(dey) + ", deltaI = " + f.format(di) + " deg, deltaRAAN = " + f.format(dr) + " deg");
        PVCoordinates delta = new PVCoordinates(previous.getPVCoordinates(), orbit.getPVCoordinates());
        deltaP = delta.getPosition().getNorm();
        deltaV = delta.getVelocity().getNorm();
    }
    // final orbit
    System.out.println();
    System.out.println("final orbit (osculating): " + orbit);
    // generate the ground track grid file
    try (PrintStream output = new PrintStream(new File(input.getParent(), gridOutput), "UTF-8")) {
        for (int i = 0; i < gridLatitudes.length; ++i) {
            printGridPoints(output, gridLatitudes[i], gridAscending[i], orbit, propagator, nbOrbits);
        }
    }
}

Example 58 with CircularOrbit

use of org.orekit.orbits.CircularOrbit in project Orekit by CS-SI.

the class Phasing method improveFrozenEccentricity.

/**
 * Improve orbit to better match frozen eccentricity property.
 * @param previous previous orbit
 * @param duration sampling duration
 * @param propagator propagator to use
 * @return an improved Earth phased, Sun synchronous orbit with frozen eccentricity
 * @exception OrekitException if orbit cannot be propagated
 */
private CircularOrbit improveFrozenEccentricity(CircularOrbit previous, double duration, Propagator propagator) throws OrekitException {
    propagator.resetInitialState(new SpacecraftState(previous));
    AbsoluteDate start = previous.getDate();
    NormalizedSphericalHarmonicsProvider.NormalizedSphericalHarmonics harmonics = gravityField.onDate(previous.getDate());
    double[][] unnormalization = GravityFieldFactory.getUnnormalizationFactors(2, 0);
    double a = previous.getA();
    double sinI = FastMath.sin(previous.getI());
    double aeOa = gravityField.getAe() / a;
    double mu = gravityField.getMu();
    double n = FastMath.sqrt(mu / a) / a;
    double j2 = -unnormalization[2][0] * harmonics.getNormalizedCnm(2, 0);
    double frozenPulsation = 3 * n * j2 * aeOa * aeOa * (1 - 1.25 * sinI * sinI);
    // fit the eccentricity to an harmonic model with short and medium periods
    // we will only use the medium periods part for the correction
    SecularAndHarmonic exModel = new SecularAndHarmonic(0, frozenPulsation, n, 2 * n);
    SecularAndHarmonic eyModel = new SecularAndHarmonic(0, frozenPulsation, n, 2 * n);
    exModel.resetFitting(start, new double[] { previous.getCircularEx(), -1.0e-10, 1.0e-5, 1.0e-5, 1.0e-5, 1.0e-5, 1.0e-5 });
    eyModel.resetFitting(start, new double[] { previous.getCircularEy(), -1.0e-10, 1.0e-5, 1.0e-5, 1.0e-5, 1.0e-5, 1.0e-5 });
    final double step = previous.getKeplerianPeriod() / 5;
    for (double dt = 0; dt < duration; dt += step) {
        final SpacecraftState state = propagator.propagate(start.shiftedBy(dt));
        final CircularOrbit orbit = (CircularOrbit) OrbitType.CIRCULAR.convertType(state.getOrbit());
        exModel.addPoint(state.getDate(), orbit.getCircularEx());
        eyModel.addPoint(state.getDate(), orbit.getCircularEy());
    }
    // adjust eccentricity
    exModel.fit();
    double dex = -exModel.getFittedParameters()[1];
    eyModel.fit();
    double dey = -eyModel.getFittedParameters()[1];
    // put the eccentricity at center of frozen center
    return new CircularOrbit(previous.getA(), previous.getCircularEx() + dex, previous.getCircularEy() + dey, previous.getI(), previous.getRightAscensionOfAscendingNode(), previous.getAlphaV(), PositionAngle.TRUE, previous.getFrame(), previous.getDate(), previous.getMu());
}

Example 59 with CircularOrbit

use of org.orekit.orbits.CircularOrbit in project Orekit by CS-SI.

the class Phasing method guessOrbit.

/**
 * Guess an initial orbit from theoretical model.
 * @param date orbit date
 * @param frame frame to use for defining orbit
 * @param nbOrbits number of orbits in the phasing cycle
 * @param nbDays number of days in the phasing cycle
 * @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
 * @return an initial guess of Earth phased, Sun synchronous orbit
 * @exception OrekitException if mean solar time cannot be computed
 */
private CircularOrbit guessOrbit(AbsoluteDate date, Frame frame, int nbOrbits, int nbDays, double latitude, boolean ascending, double mst) throws OrekitException {
    double mu = gravityField.getMu();
    NormalizedSphericalHarmonicsProvider.NormalizedSphericalHarmonics harmonics = gravityField.onDate(date);
    // initial semi major axis guess based on Keplerian period
    double period0 = (nbDays * Constants.JULIAN_DAY) / nbOrbits;
    double n0 = 2 * FastMath.PI / period0;
    double a0 = FastMath.cbrt(mu / (n0 * n0));
    // initial inclination guess based on ascending node drift due to J2
    double[][] unnormalization = GravityFieldFactory.getUnnormalizationFactors(3, 0);
    double j2 = -unnormalization[2][0] * harmonics.getNormalizedCnm(2, 0);
    double j3 = -unnormalization[3][0] * harmonics.getNormalizedCnm(3, 0);
    double raanRate = 2 * FastMath.PI / Constants.JULIAN_YEAR;
    double ae = gravityField.getAe();
    double i0 = FastMath.acos(-raanRate * a0 * a0 / (1.5 * ae * ae * j2 * n0));
    // initial eccentricity guess based on J2 and J3
    double ex0 = 0;
    double ey0 = -j3 * ae * FastMath.sin(i0) / (2 * a0 * j2);
    // initial ascending node guess based on mean solar time
    double alpha0 = FastMath.asin(FastMath.sin(latitude) / FastMath.sin(i0));
    if (!ascending) {
        alpha0 = FastMath.PI - alpha0;
    }
    double h = meanSolarTime(new CircularOrbit(a0, ex0, ey0, i0, 0.0, alpha0, PositionAngle.TRUE, frame, date, mu));
    double raan0 = FastMath.PI * (mst - h) / 12.0;
    return new CircularOrbit(a0, ex0, ey0, i0, raan0, alpha0, PositionAngle.TRUE, frame, date, mu);
}

Example 60 with CircularOrbit

use of org.orekit.orbits.CircularOrbit in project Orekit by CS-SI.

the class OrbitDetermination method createOrbit.

/**
 * 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 inertial frame cannot be created
 */
private Orbit createOrbit(final KeyValueFileParser<ParameterKey> parser, final double mu) throws NoSuchElementException, OrekitException {
    final Frame frame;
    if (!parser.containsKey(ParameterKey.INERTIAL_FRAME)) {
        frame = FramesFactory.getEME2000();
    } else {
        frame = parser.getInertialFrame(ParameterKey.INERTIAL_FRAME);
    }
    // Orbit definition
    PositionAngle angleType = PositionAngle.MEAN;
    if (parser.containsKey(ParameterKey.ORBIT_ANGLE_TYPE)) {
        angleType = PositionAngle.valueOf(parser.getString(ParameterKey.ORBIT_ANGLE_TYPE).toUpperCase());
    }
    if (parser.containsKey(ParameterKey.ORBIT_KEPLERIAN_A)) {
        return new KeplerianOrbit(parser.getDouble(ParameterKey.ORBIT_KEPLERIAN_A), parser.getDouble(ParameterKey.ORBIT_KEPLERIAN_E), parser.getAngle(ParameterKey.ORBIT_KEPLERIAN_I), parser.getAngle(ParameterKey.ORBIT_KEPLERIAN_PA), parser.getAngle(ParameterKey.ORBIT_KEPLERIAN_RAAN), parser.getAngle(ParameterKey.ORBIT_KEPLERIAN_ANOMALY), angleType, frame, parser.getDate(ParameterKey.ORBIT_DATE, TimeScalesFactory.getUTC()), mu);
    } else if (parser.containsKey(ParameterKey.ORBIT_EQUINOCTIAL_A)) {
        return new EquinoctialOrbit(parser.getDouble(ParameterKey.ORBIT_EQUINOCTIAL_A), parser.getDouble(ParameterKey.ORBIT_EQUINOCTIAL_EX), parser.getDouble(ParameterKey.ORBIT_EQUINOCTIAL_EY), parser.getDouble(ParameterKey.ORBIT_EQUINOCTIAL_HX), parser.getDouble(ParameterKey.ORBIT_EQUINOCTIAL_HY), parser.getAngle(ParameterKey.ORBIT_EQUINOCTIAL_LAMBDA), angleType, frame, parser.getDate(ParameterKey.ORBIT_DATE, TimeScalesFactory.getUTC()), mu);
    } else if (parser.containsKey(ParameterKey.ORBIT_CIRCULAR_A)) {
        return new CircularOrbit(parser.getDouble(ParameterKey.ORBIT_CIRCULAR_A), parser.getDouble(ParameterKey.ORBIT_CIRCULAR_EX), parser.getDouble(ParameterKey.ORBIT_CIRCULAR_EY), parser.getAngle(ParameterKey.ORBIT_CIRCULAR_I), parser.getAngle(ParameterKey.ORBIT_CIRCULAR_RAAN), parser.getAngle(ParameterKey.ORBIT_CIRCULAR_ALPHA), angleType, frame, parser.getDate(ParameterKey.ORBIT_DATE, TimeScalesFactory.getUTC()), mu);
    } else if (parser.containsKey(ParameterKey.ORBIT_TLE_LINE_1)) {
        final String line1 = parser.getString(ParameterKey.ORBIT_TLE_LINE_1);
        final String line2 = parser.getString(ParameterKey.ORBIT_TLE_LINE_2);
        final TLE tle = new TLE(line1, line2);
        TLEPropagator propagator = TLEPropagator.selectExtrapolator(tle);
        // propagator.setEphemerisMode();
        AbsoluteDate initDate = tle.getDate();
        SpacecraftState initialState = propagator.getInitialState();
        // Transformation from TEME to frame.
        return new CartesianOrbit(initialState.getPVCoordinates(FramesFactory.getEME2000()), frame, initDate, mu);
    } else {
        final double[] pos = { parser.getDouble(ParameterKey.ORBIT_CARTESIAN_PX), parser.getDouble(ParameterKey.ORBIT_CARTESIAN_PY), parser.getDouble(ParameterKey.ORBIT_CARTESIAN_PZ) };
        final double[] vel = { parser.getDouble(ParameterKey.ORBIT_CARTESIAN_VX), parser.getDouble(ParameterKey.ORBIT_CARTESIAN_VY), parser.getDouble(ParameterKey.ORBIT_CARTESIAN_VZ) };
        return new CartesianOrbit(new PVCoordinates(new Vector3D(pos), new Vector3D(vel)), frame, parser.getDate(ParameterKey.ORBIT_DATE, TimeScalesFactory.getUTC()), mu);
    }
}