List of usage examples for org.apache.commons.math3.ode.nonstiff DormandPrince853Integrator DormandPrince853Integrator
public DormandPrince853Integrator(final double minStep, final double maxStep, final double[] vecAbsoluteTolerance, final double[] vecRelativeTolerance)
From source file:jat.examples.TwoBodyExample.TwoBodyExample.java
public static void main(String[] args) { TwoBodyExample x = new TwoBodyExample(); // create a TwoBody orbit using orbit elements TwoBodyAPL sat = new TwoBodyAPL(7000.0, 0.3, 0.0, 0.0, 0.0, 0.0); double[] y = sat.randv(); ArrayRealVector v = new ArrayRealVector(y); DecimalFormat df2 = new DecimalFormat("#,###,###,##0.00"); RealVectorFormat format = new RealVectorFormat(df2); System.out.println(format.format(v)); // find out the period of the orbit double period = sat.period(); // set the final time = one orbit period double tf = period; // set the initial time to zero double t0 = 0.0; // propagate the orbit FirstOrderIntegrator dp853 = new DormandPrince853Integrator(1.0e-8, 100.0, 1.0e-10, 1.0e-10); dp853.addStepHandler(sat.stepHandler); // double[] y = new double[] { 7000.0, 0, 0, .0, 8, 0 }; // initial // state//from www. j av a 2s. com dp853.integrate(sat, 0.0, y, 8000, y); // now y contains final state at // tf Double[] objArray = sat.time.toArray(new Double[sat.time.size()]); double[] timeArray = ArrayUtils.toPrimitive(objArray); double[] xsolArray = ArrayUtils.toPrimitive(sat.xsol.toArray(new Double[sat.time.size()])); double[] ysolArray = ArrayUtils.toPrimitive(sat.ysol.toArray(new Double[sat.time.size()])); double[][] XY = new double[timeArray.length][2]; // int a=0; // System.arraycopy(timeArray,0,XY[a],0,timeArray.length); // System.arraycopy(ysolArray,0,XY[1],0,ysolArray.length); for (int i = 0; i < timeArray.length; i++) { XY[i][0] = xsolArray[i]; XY[i][1] = ysolArray[i]; } Plot2DPanel p = new Plot2DPanel(); // Plot2DPanel p = new Plot2DPanel(min, max, axesScales, axesLabels); ScatterPlot s = new ScatterPlot("orbit", Color.RED, XY); // LinePlot l = new LinePlot("sin", Color.RED, XY); // l.closed_curve = false; // l.draw_dot = true; p.addPlot(s); p.setLegendOrientation(PlotPanel.SOUTH); double plotSize = 10000.; double[] min = { -plotSize, -plotSize }; double[] max = { plotSize, plotSize }; p.setFixedBounds(min, max); new FrameView(p).setDefaultCloseOperation(JFrame.EXIT_ON_CLOSE); System.out.println("end"); }
From source file:fr.cs.examples.propagation.MasterMode.java
/** Program entry point. * @param args program arguments (unused here) *//* w ww. j a v a 2s. c om*/ public static void main(String[] args) { try { // configure Orekit Autoconfiguration.configureOrekit(); // gravitation coefficient double mu = 3.986004415e+14; // inertial frame Frame inertialFrame = FramesFactory.getEME2000(); // Initial date AbsoluteDate initialDate = new AbsoluteDate(2004, 01, 01, 23, 30, 00.000, TimeScalesFactory.getUTC()); // Initial orbit double a = 24396159; // semi major axis in meters double e = 0.72831215; // eccentricity double i = FastMath.toRadians(7); // inclination double omega = FastMath.toRadians(180); // perigee argument double raan = FastMath.toRadians(261); // right ascention of ascending node double lM = 0; // mean anomaly 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 = 10.0; final OrbitType propagationType = OrbitType.KEPLERIAN; final double[][] tolerances = NumericalPropagator.tolerances(positionTolerance, initialOrbit, propagationType); AdaptiveStepsizeIntegrator integrator = new DormandPrince853Integrator(minStep, maxstep, tolerances[0], tolerances[1]); // Propagator 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); // Add force model to the propagator propagator.addForceModel(holmesFeatherstone); // Set up initial state in the propagator propagator.setInitialState(initialState); // Set up operating mode for the propagator as master mode // with fixed step and specialized step handler propagator.setMasterMode(60., new TutorialStepHandler()); // Extrapolate from the initial to the final date SpacecraftState finalState = propagator.propagate(initialDate.shiftedBy(630.)); KeplerianOrbit o = (KeplerianOrbit) OrbitType.KEPLERIAN.convertType(finalState.getOrbit()); System.out.format(Locale.US, "Final state:%n%s %12.3f %10.8f %10.6f %10.6f %10.6f %10.6f%n", finalState.getDate(), o.getA(), o.getE(), FastMath.toDegrees(o.getI()), FastMath.toDegrees(o.getPerigeeArgument()), FastMath.toDegrees(o.getRightAscensionOfAscendingNode()), FastMath.toDegrees(o.getTrueAnomaly())); } catch (OrekitException oe) { System.err.println(oe.getMessage()); } }
From source file:fr.cs.examples.conversion.PropagatorConversion.java
/** Program entry point. * @param args program arguments (unused here) *//*from w ww . j ava2s .com*/ public static void main(String[] args) { try { // configure Orekit Autoconfiguration.configureOrekit(); // gravity field NormalizedSphericalHarmonicsProvider provider = GravityFieldFactory.getNormalizedProvider(2, 0); double mu = provider.getMu(); // inertial frame Frame inertialFrame = FramesFactory.getEME2000(); // Initial date AbsoluteDate initialDate = new AbsoluteDate(2004, 01, 01, 23, 30, 00.000, TimeScalesFactory.getUTC()); // Initial orbit (GTO) final double a = 24396159; // semi major axis in meters final double e = 0.72831215; // eccentricity final double i = FastMath.toRadians(7); // inclination final double omega = FastMath.toRadians(180); // perigee argument final double raan = FastMath.toRadians(261); // right ascention of ascending node final double lM = 0; // mean anomaly Orbit initialOrbit = new KeplerianOrbit(a, e, i, omega, raan, lM, PositionAngle.MEAN, inertialFrame, initialDate, mu); final double period = initialOrbit.getKeplerianPeriod(); // Initial state definition final 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.; final double dP = 1.e-2; final OrbitType orbType = OrbitType.CARTESIAN; final double[][] tol = NumericalPropagator.tolerances(dP, initialOrbit, orbType); final AbstractIntegrator integrator = new DormandPrince853Integrator(minStep, maxStep, tol[0], tol[1]); // Propagator NumericalPropagator numProp = new NumericalPropagator(integrator); numProp.setInitialState(initialState); numProp.setOrbitType(orbType); // Force Models: // 1 - Perturbing gravity field (only J2 is considered here) ForceModel gravity = new HolmesFeatherstoneAttractionModel( FramesFactory.getITRF(IERSConventions.IERS_2010, true), provider); // Add force models to the propagator numProp.addForceModel(gravity); // Propagator factory PropagatorBuilder builder = new KeplerianPropagatorBuilder(mu, inertialFrame, OrbitType.KEPLERIAN, PositionAngle.TRUE); // Propagator converter PropagatorConverter fitter = new FiniteDifferencePropagatorConverter(builder, 1.e-6, 5000); // Resulting propagator KeplerianPropagator kepProp = (KeplerianPropagator) fitter.convert(numProp, 2 * period, 251); // Step handlers StatesHandler numStepHandler = new StatesHandler(); StatesHandler kepStepHandler = new StatesHandler(); // Set up operating mode for the propagator as master mode // with fixed step and specialized step handler numProp.setMasterMode(60., numStepHandler); kepProp.setMasterMode(60., kepStepHandler); // Extrapolate from the initial to the final date numProp.propagate(initialDate.shiftedBy(10. * period)); kepProp.propagate(initialDate.shiftedBy(10. * period)); // retrieve the states List<SpacecraftState> numStates = numStepHandler.getStates(); List<SpacecraftState> kepStates = kepStepHandler.getStates(); // Print the results on the output file File output = new File(new File(System.getProperty("user.home")), "elements.dat"); PrintStream stream = new PrintStream(output); stream.println("# date Anum Akep Enum Ekep Inum Ikep LMnum LMkep"); for (SpacecraftState numState : numStates) { for (SpacecraftState kepState : kepStates) { if (numState.getDate().compareTo(kepState.getDate()) == 0) { stream.println(numState.getDate() + " " + numState.getA() + " " + kepState.getA() + " " + numState.getE() + " " + kepState.getE() + " " + FastMath.toDegrees(numState.getI()) + " " + FastMath.toDegrees(kepState.getI()) + " " + FastMath.toDegrees(MathUtils.normalizeAngle(numState.getLM(), FastMath.PI)) + " " + FastMath.toDegrees(MathUtils.normalizeAngle(kepState.getLM(), FastMath.PI))); break; } } } stream.close(); System.out.println("Results saved as file " + output); File output1 = new File(new File(System.getProperty("user.home")), "elts_pv.dat"); PrintStream stream1 = new PrintStream(output1); stream.println("# date pxn pyn pzn vxn vyn vzn pxk pyk pzk vxk vyk vzk"); for (SpacecraftState numState : numStates) { for (SpacecraftState kepState : kepStates) { if (numState.getDate().compareTo(kepState.getDate()) == 0) { final double pxn = numState.getPVCoordinates().getPosition().getX(); final double pyn = numState.getPVCoordinates().getPosition().getY(); final double pzn = numState.getPVCoordinates().getPosition().getZ(); final double vxn = numState.getPVCoordinates().getVelocity().getX(); final double vyn = numState.getPVCoordinates().getVelocity().getY(); final double vzn = numState.getPVCoordinates().getVelocity().getZ(); final double pxk = kepState.getPVCoordinates().getPosition().getX(); final double pyk = kepState.getPVCoordinates().getPosition().getY(); final double pzk = kepState.getPVCoordinates().getPosition().getZ(); final double vxk = kepState.getPVCoordinates().getVelocity().getX(); final double vyk = kepState.getPVCoordinates().getVelocity().getY(); final double vzk = kepState.getPVCoordinates().getVelocity().getZ(); stream1.println(numState.getDate() + " " + pxn + " " + pyn + " " + pzn + " " + vxn + " " + vyn + " " + vzn + " " + pxk + " " + pyk + " " + pzk + " " + vxk + " " + vyk + " " + vzk); break; } } } stream1.close(); System.out.println("Results saved as file " + output1); } catch (OrekitException oe) { System.err.println(oe.getLocalizedMessage()); System.exit(1); } catch (FileNotFoundException fnfe) { System.err.println(fnfe.getLocalizedMessage()); System.exit(1); } }
From source file:gentracklets.Propagate.java
public static void main(String[] args) throws OrekitException { // load the data files File data = new File("/home/zittersteijn/Documents/java/libraries/orekit-data.zip"); DataProvidersManager DM = DataProvidersManager.getInstance(); ZipJarCrawler crawler = new ZipJarCrawler(data); DM.clearProviders();// w w w. j a v a2s. c om DM.addProvider(crawler); // Read in TLE elements File tleFile = new File("/home/zittersteijn/Documents/TLEs/ASTRA20151207.tle"); FileReader TLEfr; Vector<TLE> tles = new Vector<>(); tles.setSize(30); try { // read and save TLEs to a vector TLEfr = new FileReader("/home/zittersteijn/Documents/TLEs/ASTRA20151207.tle"); BufferedReader readTLE = new BufferedReader(TLEfr); Scanner s = new Scanner(tleFile); String line1, line2; TLE2 tle = new TLE2(); int nrOfObj = 4; for (int ii = 1; ii < nrOfObj + 1; ii++) { System.out.println(ii); line1 = s.nextLine(); line2 = s.nextLine(); if (TLE.isFormatOK(line1, line2)) { tles.setElementAt(new TLE(line1, line2), ii); System.out.println(tles.get(ii).toString()); } else { System.out.println("format problem"); } } readTLE.close(); // define a groundstation Frame inertialFrame = FramesFactory.getEME2000(); TimeScale utc = TimeScalesFactory.getUTC(); double longitude = FastMath.toRadians(7.465); double latitude = FastMath.toRadians(46.87); double altitude = 950.; GeodeticPoint station = new GeodeticPoint(latitude, longitude, altitude); Frame earthFrame = FramesFactory.getITRF(IERSConventions.IERS_2010, true); BodyShape earth = new OneAxisEllipsoid(Constants.WGS84_EARTH_EQUATORIAL_RADIUS, Constants.WGS84_EARTH_FLATTENING, earthFrame); TopocentricFrame staF = new TopocentricFrame(earth, station, "station"); Vector<Orbit> eles = new Vector<>(); eles.setSize(tles.size()); for (int ii = 1; ii < nrOfObj + 1; ii++) { double a = FastMath.pow(Constants.WGS84_EARTH_MU / FastMath.pow(tles.get(ii).getMeanMotion(), 2), (1.0 / 3)); // convert them to orbits Orbit kep = new KeplerianOrbit(a, tles.get(ii).getE(), tles.get(ii).getI(), tles.get(ii).getPerigeeArgument(), tles.get(ii).getRaan(), tles.get(ii).getMeanAnomaly(), PositionAngle.MEAN, inertialFrame, tles.get(ii).getDate(), Constants.WGS84_EARTH_MU); eles.setElementAt(kep, ii); // set up propagators KeplerianPropagator kepler = new KeplerianPropagator(eles.get(ii)); System.out.println("a: " + a); // Initial state definition double mass = 1000.0; SpacecraftState initialState = new SpacecraftState(kep, mass); // Adaptive step integrator // with a minimum step of 0.001 and a maximum step of 1000 double minStep = 0.001; double maxstep = 1000.0; double positionTolerance = 10.0; OrbitType propagationType = OrbitType.KEPLERIAN; double[][] tolerances = NumericalPropagator.tolerances(positionTolerance, kep, propagationType); AdaptiveStepsizeIntegrator integrator = new DormandPrince853Integrator(minStep, maxstep, tolerances[0], tolerances[1]); NumericalPropagator propagator = new NumericalPropagator(integrator); propagator.setOrbitType(propagationType); // set up and add force models double AMR = 4.0; double crossSection = mass * AMR; double Cd = 0.01; double Cr = 0.5; double Co = 0.8; NormalizedSphericalHarmonicsProvider provider = GravityFieldFactory.getNormalizedProvider(4, 4); ForceModel holmesFeatherstone = new HolmesFeatherstoneAttractionModel( FramesFactory.getITRF(IERSConventions.IERS_2010, true), provider); SphericalSpacecraft ssc = new SphericalSpacecraft(crossSection, Cd, Cr, Co); PVCoordinatesProvider sun = CelestialBodyFactory.getSun(); SolarRadiationPressure srp = new SolarRadiationPressure(sun, Constants.WGS84_EARTH_EQUATORIAL_RADIUS, ssc); // propagator.addForceModel(srp); // propagator.addForceModel(holmesFeatherstone); propagator.setInitialState(initialState); // propagate the orbits with steps size and tracklet lenght at several epochs (tracklets) Vector<AbsoluteDate> startDates = new Vector<>(); startDates.setSize(1); startDates.setElementAt(new AbsoluteDate(2016, 1, 26, 20, 00, 00, utc), 0); // set the step size [s] and total length double tstep = 100; double ld = 3; double ls = FastMath.floor(ld * (24 * 60 * 60) / tstep); System.out.println(ls); SpacecraftState currentStateKep = kepler.propagate(startDates.get(0)); SpacecraftState currentStatePer = propagator.propagate(startDates.get(0)); for (int tt = 0; tt < startDates.size(); tt++) { // set up output file String app = tles.get(ii).getSatelliteNumber() + "_" + startDates.get(tt) + ".txt"; // with formatted output File file1 = new File("/home/zittersteijn/Documents/propagate/keplerian/MEO/" + app); File file2 = new File("/home/zittersteijn/Documents/propagate/perturbed/MEO/" + app); file1.createNewFile(); file2.createNewFile(); Formatter fmt1 = new Formatter(file1); Formatter fmt2 = new Formatter(file2); for (int kk = 0; kk < (int) ls; kk++) { AbsoluteDate propDate = startDates.get(tt).shiftedBy(tstep * kk); currentStateKep = kepler.propagate(propDate); currentStatePer = propagator.propagate(propDate); System.out.println(currentStateKep.getPVCoordinates().getPosition() + "\t" + currentStateKep.getDate()); // convert to RADEC coordinates double[] radecKep = conversions.geo2radec(currentStateKep.getPVCoordinates(), staF, inertialFrame, propDate); double[] radecPer = conversions.geo2radec(currentStatePer.getPVCoordinates(), staF, inertialFrame, propDate); // write the orbit to seperate files with the RA, DEC, epoch and fence given AbsoluteDate year = new AbsoluteDate(YEAR, utc); fmt1.format("%.12f %.12f %.12f %d%n", radecKep[0], radecKep[2], (currentStateKep.getDate().durationFrom(year) / (24 * 3600)), (tt + 1)); fmt2.format("%.12f %.12f %.12f %d%n", radecPer[0], radecPer[2], (currentStateKep.getDate().durationFrom(year) / (24 * 3600)), (tt + 1)); } fmt1.flush(); fmt1.close(); fmt2.flush(); fmt2.close(); } double[] radecKep = conversions.geo2radec(currentStateKep.getPVCoordinates(), staF, inertialFrame, new AbsoluteDate(startDates.get(0), ls * tstep)); double[] radecPer = conversions.geo2radec(currentStatePer.getPVCoordinates(), staF, inertialFrame, new AbsoluteDate(startDates.get(0), ls * tstep)); double sig0 = 1.0 / 3600.0 / 180.0 * FastMath.PI; double dRA = radecKep[0] - radecPer[0] / (sig0 * sig0); double dDEC = radecKep[2] - radecPer[2] / (sig0 * sig0); System.out.println(dRA + "\t" + dDEC); } } catch (FileNotFoundException ex) { Logger.getLogger(GenTracklets.class.getName()).log(Level.SEVERE, null, ex); } catch (IOException iox) { Logger.getLogger(GenTracklets.class.getName()).log(Level.SEVERE, null, iox); } }
From source file:gentracklets.GenTracklets.java
public static void main(String[] args) throws OrekitException { // load the data files File data = new File("/home/zittersteijn/Documents/java/libraries/orekit-data.zip"); DataProvidersManager DM = DataProvidersManager.getInstance(); ZipJarCrawler crawler = new ZipJarCrawler(data); DM.clearProviders();/*from www . j a v a 2 s . c o m*/ DM.addProvider(crawler); // Read in TLE elements File tleFile = new File("/home/zittersteijn/Documents/TLEs/ASTRA20151207.tle"); FileReader TLEfr; Vector<TLE> tles = new Vector<>(); tles.setSize(30); try { // read and save TLEs to a vector TLEfr = new FileReader("/home/zittersteijn/Documents/TLEs/ASTRA20151207.tle"); BufferedReader readTLE = new BufferedReader(TLEfr); Scanner s = new Scanner(tleFile); String line1, line2; TLE2 tle = new TLE2(); int nrOfObj = 4; for (int ii = 1; ii < nrOfObj + 1; ii++) { System.out.println(ii); line1 = s.nextLine(); line2 = s.nextLine(); if (TLE.isFormatOK(line1, line2)) { tles.setElementAt(new TLE(line1, line2), ii); System.out.println(tles.get(ii).toString()); } else { System.out.println("format problem"); } } readTLE.close(); // define a groundstation Frame inertialFrame = FramesFactory.getEME2000(); TimeScale utc = TimeScalesFactory.getUTC(); double longitude = FastMath.toRadians(7.465); double latitude = FastMath.toRadians(46.87); double altitude = 950.; GeodeticPoint station = new GeodeticPoint(latitude, longitude, altitude); Frame earthFrame = FramesFactory.getITRF(IERSConventions.IERS_2010, true); BodyShape earth = new OneAxisEllipsoid(Constants.WGS84_EARTH_EQUATORIAL_RADIUS, Constants.WGS84_EARTH_FLATTENING, earthFrame); TopocentricFrame staF = new TopocentricFrame(earth, station, "station"); Vector<Orbit> eles = new Vector<>(); eles.setSize(tles.size()); for (int ii = 1; ii < nrOfObj + 1; ii++) { double a = FastMath.pow(Constants.WGS84_EARTH_MU / FastMath.pow(tles.get(ii).getMeanMotion(), 2), (1.0 / 3)); // convert them to orbits Orbit kep = new KeplerianOrbit(a, tles.get(ii).getE(), tles.get(ii).getI(), tles.get(ii).getPerigeeArgument(), tles.get(ii).getRaan(), tles.get(ii).getMeanAnomaly(), PositionAngle.MEAN, inertialFrame, tles.get(ii).getDate(), Constants.WGS84_EARTH_MU); eles.setElementAt(kep, ii); // set up propagators KeplerianPropagator kepler = new KeplerianPropagator(eles.get(ii)); System.out.println("a: " + a); // Initial state definition double mass = 1000.0; SpacecraftState initialState = new SpacecraftState(kep, mass); // Adaptive step integrator // with a minimum step of 0.001 and a maximum step of 1000 double minStep = 0.001; double maxstep = 1000.0; double positionTolerance = 10.0; OrbitType propagationType = OrbitType.KEPLERIAN; double[][] tolerances = NumericalPropagator.tolerances(positionTolerance, kep, propagationType); AdaptiveStepsizeIntegrator integrator = new DormandPrince853Integrator(minStep, maxstep, tolerances[0], tolerances[1]); NumericalPropagator propagator = new NumericalPropagator(integrator); propagator.setOrbitType(propagationType); // set up and add force models double AMR = 0.4; double crossSection = mass * AMR; double Cd = 0.01; double Cr = 0.5; double Co = 0.8; NormalizedSphericalHarmonicsProvider provider = GravityFieldFactory.getNormalizedProvider(4, 4); ForceModel holmesFeatherstone = new HolmesFeatherstoneAttractionModel( FramesFactory.getITRF(IERSConventions.IERS_2010, true), provider); SphericalSpacecraft ssc = new SphericalSpacecraft(crossSection, Cd, Cr, Co); PVCoordinatesProvider sun = CelestialBodyFactory.getSun(); SolarRadiationPressure srp = new SolarRadiationPressure(sun, Constants.WGS84_EARTH_EQUATORIAL_RADIUS, ssc); propagator.addForceModel(srp); propagator.addForceModel(holmesFeatherstone); propagator.setInitialState(initialState); // propagate the orbits with steps size and tracklet lenght at several epochs (tracklets) Vector<AbsoluteDate> startDates = new Vector<>(); startDates.setSize(3); startDates.setElementAt(new AbsoluteDate(2015, 12, 8, 20, 00, 00, utc), 0); startDates.setElementAt(new AbsoluteDate(2015, 12, 9, 21, 00, 00, utc), 1); startDates.setElementAt(new AbsoluteDate(2015, 12, 10, 22, 00, 00, utc), 2); double tstep = 30; int l = 7; for (int tt = 0; tt < startDates.size(); tt++) { // set up output file String app = "S_" + tles.get(ii).getSatelliteNumber() + "_" + startDates.get(tt) + ".txt"; // FileWriter trackletsOutKep = new FileWriter("/home/zittersteijn/Documents/tracklets/simulated/keplerian/ASTRA/dt1h/AMR040/" + app); // FileWriter trackletsOutPer = new FileWriter("/home/zittersteijn/Documents/tracklets/simulated/perturbed/ASTRA/dt1h/AMR040/" + app); // BufferedWriter trackletsKepBW = new BufferedWriter(trackletsOutKep); // BufferedWriter trackletsPerBW = new BufferedWriter(trackletsOutPer); // with formatted output File file1 = new File( "/home/zittersteijn/Documents/tracklets/simulated/keplerian/ASTRA/dt1d/AMR040/" + app); File file2 = new File( "/home/zittersteijn/Documents/tracklets/simulated/perturbed/ASTRA/dt1d/AMR040/" + app); file1.createNewFile(); file2.createNewFile(); Formatter fmt1 = new Formatter(file1); Formatter fmt2 = new Formatter(file2); for (int kk = 0; kk < l; kk++) { AbsoluteDate propDate = startDates.get(tt).shiftedBy(tstep * kk); SpacecraftState currentStateKep = kepler.propagate(propDate); SpacecraftState currentStatePer = propagator.propagate(propDate); System.out.println(currentStateKep.getPVCoordinates().getPosition() + "\t" + currentStateKep.getDate()); // convert to RADEC coordinates double[] radecKep = conversions.geo2radec(currentStateKep.getPVCoordinates(), staF, inertialFrame, propDate); double[] radecPer = conversions.geo2radec(currentStatePer.getPVCoordinates(), staF, inertialFrame, propDate); // write the tracklets to seperate files with the RA, DEC, epoch and fence given // System.out.println(tles.get(kk).getSatelliteNumber() + "\t" + radec[0] / (2 * FastMath.PI) * 180 + "\t" + currentState.getDate()); AbsoluteDate year = new AbsoluteDate(YEAR, utc); fmt1.format("%.12f %.12f %.12f %d%n", radecKep[0], radecKep[2], (currentStateKep.getDate().durationFrom(year) / (24 * 3600)), (tt + 1)); fmt2.format("%.12f %.12f %.12f %d%n", radecPer[0], radecPer[2], (currentStateKep.getDate().durationFrom(year) / (24 * 3600)), (tt + 1)); } fmt1.flush(); fmt1.close(); fmt2.flush(); fmt2.close(); } } } catch (FileNotFoundException ex) { Logger.getLogger(GenTracklets.class.getName()).log(Level.SEVERE, null, ex); } catch (IOException iox) { Logger.getLogger(GenTracklets.class.getName()).log(Level.SEVERE, null, iox); } }
From source file:jat.examples.ephemeris.DE405PropagateText.java
void doExample() { //double tf = 10000000.; double tf = 3600 * 24 * 500; double[] y = { 2e8, 0, 0, 0, 24.2, 0 }; // initial state PathUtil path = new PathUtil(); DE405Plus Eph = new DE405Plus(path); Eph.setFrame(DE405Frame.frame.HEE); Eph.printSteps = true;// www .j av a 2 s . c o m TimeAPL myTime = new TimeAPL(2003, 3, 1, 12, 0, 0); Eph.setIntegrationStartTime(myTime); Eph.bodyGravOnOff[body.SUN.ordinal()] = true; FirstOrderIntegrator dp853 = new DormandPrince853Integrator(1.0e-8, tf / 10.0, 1.0e-10, 1.0e-10); dp853.addStepHandler(Eph.stepHandler); FirstOrderDifferentialEquations ode = Eph; dp853.integrate(ode, 0.0, y, tf, y); // now y contains final state // at // time tf if (print) { String nf = "%10.3f "; String format = nf + nf + nf + nf + nf; System.out.printf(format, tf, y[0], y[1], y[2], Eph.energy(tf, y)); System.out.println(); } }
From source file:edu.gcsc.vrl.commons.math.ode.ODESolver.java
public Trajectory solve(@ParamInfo(name = "Label", options = "value=\"Label 1\"") String label, @ParamInfo(name = "x0", options = "value=0.0D") double t0, @ParamInfo(name = "xn", options = "value=3.0D") double tn, @ParamInfo(name = "y0", options = "value=0.0") double y0, @ParamInfo(name = "Min Step", options = "value=1e-6D") double minStep, @ParamInfo(name = "Max Step", options = "value=1e-2D") double maxStep, @ParamInfo(name = "Abs.Tol.", options = "value=1e-10") double absTol, @ParamInfo(name = "Rel.Tol.", options = "value=1e-10") double relTol, @ParamInfo(name = "RHS") FirstOrderDifferentialEquations rhs) { FirstOrderIntegrator integrator = new DormandPrince853Integrator(minStep, maxStep, absTol, relTol); final Trajectory result = new Trajectory(label); StepHandler stepHandler = new StepHandler() { @Override/*from w w w . java 2s.c o m*/ public void init(double t0, double[] y0, double t) { result.add(t, y0[0]); } @Override public void handleStep(StepInterpolator interpolator, boolean isLast) { double t = interpolator.getCurrentTime(); double[] y = interpolator.getInterpolatedState(); result.add(t, y[0]); } }; integrator.addStepHandler(stepHandler); double[] y = new double[] { y0 }; // initial state integrator.integrate(rhs, t0, y, tn, y); return result; }
From source file:jat.examples.CRTBP.CRTBPPlot.java
void doExample() { double mu = 0.15; double[] y0 = { .11, 0, 0, 1.35, 1.33, 0 }; // initial state // double mu = 0.1; // double mu = 3.035909999e-6; // double mu = 0.012277471; // double[] y0 = { .1, 0, 0, 2.69, 2.69, 0 }; // initial state // double mu = 0.2; CRTBP myCRTBP = new CRTBP(mu); FirstOrderIntegrator dp853 = new DormandPrince853Integrator(1.0e-8, 100.0, 1.0e-10, 1.0e-10); dp853.addStepHandler(myCRTBP.stepHandler); FirstOrderDifferentialEquations ode = myCRTBP; double tf;//from w ww . ja v a 2 s . c o m double[] y = new double[6]; // initial state // for (int i = 1; i < 2; i++) { // tf = i * 20.; tf = 40.; System.arraycopy(y0, 0, y, 0, 6); dp853.integrate(ode, 0.0, y, tf, y); // now y contains final state // at // time tf if (print) { System.out.printf("%9.6f %9.6f %9.6f %9.6f %9.6f", tf, y[0], y[1], y[2], myCRTBP.JacobiIntegral(y)); System.out.println(); } int arraySize = myCRTBP.time.size(); double[] timeArray = ArrayUtils.toPrimitive(myCRTBP.time.toArray(new Double[arraySize])); double[] xsolArray = ArrayUtils.toPrimitive(myCRTBP.xsol.toArray(new Double[arraySize])); double[] ysolArray = ArrayUtils.toPrimitive(myCRTBP.ysol.toArray(new Double[arraySize])); double[][] XY = new double[timeArray.length][2]; for (int i = 0; i < timeArray.length; i++) { XY[i][0] = xsolArray[i]; XY[i][1] = ysolArray[i]; } Plot2DPanel p = new Plot2DPanel(); LinePlot l = new LinePlot("spacecraft", Color.RED, XY); l.closed_curve = false; l.draw_dot = true; p.addPlot(l); double plotSize = 1.2; myCRTBP.findLibrationPoints(); Color darkGreen = new java.awt.Color(0, 190, 0); addPoint(p, "Earth", java.awt.Color.BLUE, -mu, 0); addPoint(p, "Moon", java.awt.Color.gray, 1 - mu, 0); addPoint(p, "L1", darkGreen, myCRTBP.LibPoints[0].getX(), 0); addPoint(p, "L2", darkGreen, myCRTBP.LibPoints[1].getX(), 0); addPoint(p, "L3", darkGreen, myCRTBP.LibPoints[2].getX(), 0); String Labelmu = "mu = " + myCRTBP.mu; p.addLabel(Labelmu, java.awt.Color.black, 1, .9 * plotSize); String initial = "initial x,v = (" + y0[0] + "," + y0[1] + "),(" + y0[3] + "," + y0[4] + ")"; p.addLabel(initial, java.awt.Color.black, 1, .8 * plotSize); String Jacobi = "spacecraft C = " + myCRTBP.C; p.addLabel(Jacobi, java.awt.Color.black, 1, .7 * plotSize); String L1C = "L1 C = " + myCRTBP.C1; p.addLabel(L1C, java.awt.Color.black, 1, .6 * plotSize); myCRTBP.findZeroVelocity(); int size = myCRTBP.xzv.size(); double[] xzvArray = ArrayUtils.toPrimitive(myCRTBP.xzv.toArray(new Double[size])); double[] yzvArray = ArrayUtils.toPrimitive(myCRTBP.yzv.toArray(new Double[size])); double[][] XYzv = new double[size][2]; for (int i = 0; i < size; i++) { XYzv[i][0] = xzvArray[i]; XYzv[i][1] = yzvArray[i]; } LinePlot lzv = new LinePlot("zero vel", Color.blue, XYzv); lzv.closed_curve = false; lzv.draw_dot = true; p.addPlot(lzv); p.setLegendOrientation(PlotPanel.SOUTH); p.setFixedBounds(0, -plotSize, plotSize); p.setFixedBounds(1, -plotSize, plotSize); new FrameView(p).setDefaultCloseOperation(JFrame.EXIT_ON_CLOSE); }
From source file:jat.examples.ephemeris.DE405PropagatePlot.java
void doExample() { double tf = 3600 * 24 * 300; double[] y0 = { 2e8, 0, 0, 0, 24.2, 0 }; // initial state double[] y = new double[6]; PathUtil path = new PathUtil(); DE405Plus Eph = new DE405Plus(path); Eph.setFrame(DE405Frame.frame.HEE); Eph.printSteps = true;/*from ww w .ja va 2s . co m*/ TimeAPL myTime = new TimeAPL(2003, 3, 1, 12, 0, 0); Eph.setIntegrationStartTime(myTime); Eph.bodyGravOnOff[body.SUN.ordinal()] = true; // Eph.planetOnOff[body.JUPITER.ordinal()] = true; FirstOrderIntegrator dp853 = new DormandPrince853Integrator(1.0e-8, tf / 10.0, 1.0e-10, 1.0e-10); dp853.addStepHandler(Eph.stepHandler); FirstOrderDifferentialEquations ode = Eph; dp853.integrate(ode, 0.0, y0, tf, y); // now y contains final state at // time tf if (print) { String nf = "%10.3f "; String format = nf + nf + nf + nf + nf; System.out.printf(format, tf, y[0], y[1], y[2], Eph.energy(tf, y)); System.out.println(); } Plot2DPanel p = new Plot2DPanel(); LinePlot l1 = new LinePlot("Jup. off", Color.RED, getXYforPlot(Eph.xsol, Eph.ysol)); l1.closed_curve = false; p.addPlot(l1); Eph.reset(); Eph.bodyGravOnOff[body.JUPITER.ordinal()] = true; dp853.integrate(ode, 0.0, y0, tf, y); // now y contains final state at LinePlot l2 = new LinePlot("Jup. on", Color.BLUE, getXYforPlot(Eph.xsol, Eph.ysol)); l2.closed_curve = false; p.addPlot(l2); VectorN EarthPos = null; try { EarthPos = Eph.get_planet_pos(body.EARTH, myTime); } catch (IOException e) { // TODO Auto-generated catch block e.printStackTrace(); } addPoint(p, "Earth", java.awt.Color.BLUE, EarthPos.x[0], EarthPos.x[1]); p.setLegendOrientation(PlotPanel.SOUTH); double plotSize = 2e8; p.setFixedBounds(0, -plotSize, plotSize); p.setFixedBounds(1, -plotSize, plotSize); new FrameView(p).setDefaultCloseOperation(JFrame.EXIT_ON_CLOSE); }
From source file:jat.application.DE405Propagator.DE405PropagatorPlot.java
public void add_scene() { // Update Ephemeris to current user parameters for (body b : body.values()) { Eph.bodyGravOnOff[b.ordinal()] = dpMain.dpParam.bodyGravOnOff[b.ordinal()]; }//from ww w . ja v a2s . c o m Eph.setIntegrationStartTime(dpMain.dpParam.simulationDate); try { Eph.setEarthMoonPlaneNormal(dpMain.dpParam.simulationDate); } catch (IOException e1) { e1.printStackTrace(); } Eph.setFrame(dpMain.dpParam.Frame); Eph.reset(); // Spacecraft Trajectory FirstOrderIntegrator dp853 = new DormandPrince853Integrator(1.0e-8, dpMain.dpParam.tf, 1.0e-10, 1.0e-10); dp853.addStepHandler(Eph.stepHandler); FirstOrderDifferentialEquations ode = Eph; double[] y = new double[6]; dp853.integrate(ode, 0.0, dpMain.dpParam.y0, dpMain.dpParam.tf, y); if (print) { String nf = "%10.3f "; String format = nf + nf + nf + nf + nf; System.out.printf(format, dpMain.dpParam.tf, y[0], y[1], y[2], Eph.energy(dpMain.dpParam.tf, y)); System.out.println(); } LinePlot l1 = new LinePlot("spacecraft", Color.RED, getXYZforPlot(Eph.xsol, Eph.ysol, Eph.zsol)); l1.closed_curve = false; plot.addPlot(l1); addBodies(); // Vector3D y0v=new Vector3D(dpParam.y0[0],dpParam.y0[1],dpParam.y0[2]); // double plotBounds = 2*y0v.getNorm(); plot.setFixedBounds(0, -plotBounds, plotBounds); plot.setFixedBounds(1, -plotBounds, plotBounds); plot.setFixedBounds(2, -plotBounds, plotBounds); plot.setLegendOrientation(PlotPanel.SOUTH); }