@MaximeJ
Thank you for this tip Maxime, I’ve implemented it that way and it worked. I tried to replicate the procedure with the DSST propagator using the JB2008 atmospheric model but i get the following error.
“Exception in thread “main” org.orekit.errors.OrekitException: altitude (-1,031,992.212 m) is below the 90,000 m allowed threshold”
I have an altitude detector implemented that’s meant to stop the propagation when the altitude gets below 120km but it doesn’t seem to work with the DSST propagator. My code is shown below.
// Set path to orekit data folder
File orekitData = new File("Path to Orekit Data folder");
DataProvidersManager manager = DataContext.getDefault().getDataProvidersManager();
manager.addProvider(new DirectoryCrawler(orekitData));
// TimeScale
TimeScale timescale = TimeScalesFactory.getUTC();
// Enter TLE
String line1 = "1 23757U 95074A 12012.49666365 .00007180 00000-0 20301-3 0 9995";
String line2 = "2 23757 022.9885 165.8605 0007544 261.0176 098.9451 15.32981482889705";
TLE inputTLE = new TLE(line1, line2, timescale);
System.out.println(inputTLE.getDate());
// Retrieve Propagator used for TLE
Propagator tlePropagator = TLEPropagator.selectExtrapolator(inputTLE);
Orbit initialOrbit = new CartesianOrbit(tlePropagator.getInitialState().getOrbit());
SpacecraftState initialState = new SpacecraftState(initialOrbit, mass);
// Initial reference Orbit for the builder
Orbit referenceOrbit = tlePropagator.getInitialState().getOrbit();
System.out.println("Initial orbit");
System.out.println(initialState.getDate());
System.out.println(initialState.getOrbit());
// Set up Integrator Parameters
double minStep = 360; // Minimum TimeStep
double maxStep = 86400; // Maximum TimeStep
double dP = 10; // Position Tolerance
// Builder
DormandPrince853IntegratorBuilder builder = new DormandPrince853IntegratorBuilder(minStep, maxStep, dP); // Propagator Factory
DSSTPropagatorBuilder dsstPropbuilder = new DSSTPropagatorBuilder(referenceOrbit,builder,dP, PropagationType.MEAN,PropagationType.MEAN);
dsstPropbuilder.setMass(mass);
// Add Force models to Fitted Propagator Builder
// Initialize Solar Activity Provider Data Provider
final JB2008SpaceEnvironmentData jbdata = new JB2008SpaceEnvironmentData(JB2008SpaceEnvironmentData.DEFAULT_SUPPORTED_NAMES_SOLFSMY,JB2008SpaceEnvironmentData.DEFAULT_SUPPORTED_NAMES_DTC,manager,timescale);
// Initialize JB2008 Atmospheric Model
final JB2008 jb2008 = new JB2008(jbdata,sun, earth);
dsstPropagatorBuilder.addForceModel(new DSSTAtmosphericDrag(jb2008, Cd, area, provider.getMu()));
// SGP4 Propagator Converter
FiniteDifferencePropagatorConverter dsstFitter = new FiniteDifferencePropagatorConverter(dsstPropBuilder, 1.0e-3, 1000);
// Fitting Parameters
double duration = 86400.0; // Fitting Timespan
double stepSize = 300.0; // Fitting Stepsize
Double points = duration / stepSize; // Number of fitting points
// Convert SGP4 Propagator to DSST Propagator
System.out.println("START DSST PROPAGATOR FITTING");
dsstFitter.convert(tlePropagator, duration, points.intValue()); // Conversion from SGP4 to DSST Propagator
System.out.println("END DSST PROPAGATOR FITTING");
DSSTPropagator dsstPropagator = (DSSTPropagator) dsstFitter.getAdaptedPropagator(); // Retrieving the converted numerical propagator
System.out.println(dsstPropagator.getInitialState().getOrbit());
// Converting Initial State to EME2000 Frame
TimeStampedPVCoordinates tspvc1 = dsstPropagator.getInitialState().getPVCoordinates(FramesFactory.getEME2000());
CartesianOrbit corbit = new CartesianOrbit(tspvc1,FramesFactory.getEME2000(),dsstPropagator.getMu());
SpacecraftState newState = new SpacecraftState(corbit,mass);
dsstPropagator.resetInitialState(newState);
// Create Earth Model
OneAxisEllipsoid earth = new OneAxisEllipsoid(Constants.WGS84_EARTH_EQUATORIAL_RADIUS, Constants.WGS84_EARTH_FLATTENING,
FramesFactory.getITRF(IERSConventions.IERS_2010, false));
// Re-entry Detector
EventDetector altdet = new AltitudeDetector(maxCheck, threshold, alt, earth).withHandler(((s, detector, decreasing) -> {
System.out.println("Altitude is " + FastMath.min(((s.getA() * (1 + s.getE())) - Constants.WGS84_EARTH_EQUATORIAL_RADIUS) / 1000, ((s.getA() * (1 - s.getE())) - Constants.WGS84_EARTH_EQUATORIAL_RADIUS) / 1000) + " km and the spacecraft will re-enter soon");
System.out.println("Re-entry event occurs on " + s.getDate().toString());
return Action.STOP;
}));
// Add Event detector
dsstPropagator.addEventDetector(altdet);
// Setup Step Handlers
final StatesHandler FittedStepHandler = new StatesHandler();
// Add fixed step handlers to the propagator
dsstPropagator.getMultiplexer().add(3600., FittedStepHandler);
// Get Final States
final SpacecraftState finalState = dsstPropagator.propagate(initialState.getDate().shiftedBy(10. * Constants.JULIAN_YEAR));
// Retrieve the states
final List<SpacecraftState> FittedStates = FittedStepHandler.getStates();
final File output = new File("Output path" + name + ".txt");
try (PrintStream stream = new PrintStream(output, StandardCharsets.UTF_8)) {
stream.println(" Date SMA Ecc Incl RAAN AoP MAN ");
stream.println(" [YYYY-MM-DD hh:mm:ss] (km) (deg) (deg) (deg) (deg) ");
for (SpacecraftState numstate: FittedStates){
KeplerianOrbit Orb = (KeplerianOrbit) OrbitType.KEPLERIAN.convertType(numstate.getOrbit());
String text = " %s %.4f %.7f %.5f %.5f %.5f %.5f";
text = String.format(text,Orb.getDate(),Orb.getA()/1000,Orb.getE(),FastMath.toDegrees(Orb.getI()),FastMath.toDegrees(Orb.getRightAscensionOfAscendingNode()),FastMath.toDegrees(Orb.getPerigeeArgument()),FastMath.toDegrees(Orb.getMeanAnomaly()));
stream.println(text);
}
}
final String resultSavedAs = "Output file saved to file " + name + ".txt in specified path";
System.out.println(resultSavedAs);
KeplerianOrbit Neworbit = (KeplerianOrbit) OrbitType.KEPLERIAN.convertType(finalState.getOrbit());
System.out.println("Final Fitted orbit keplerian elements");
System.out.println(Neworbit.toString());
double years = (finalState.getDate().durationFrom((initialState.getDate()))) / (Constants.JULIAN_YEAR);
System.out.println("The remaining Lifetime is " + years + "years");
System.out.println("The remaining Lifetime with a 5% margin is " + years * 1.05 + "years");
}
private static class StatesHandler implements OrekitFixedStepHandler {
// Points
private final List<SpacecraftState> states;
// Simple Constructor//
StatesHandler() {
// prepare an empty list of SpacecraftStates
states = new ArrayList<SpacecraftState>();
}
public void handleStep(final SpacecraftState currentState) {
// add the current state
states.add(currentState);
}
/**
* Get the list of handled orbital elements.
*
* @return orbital elements list
*/
public List<SpacecraftState> getStates() {
return states;
}
}