Impulsive maneuver with TLEPropagator

acapobianchi · November 26, 2021, 10:32pm

Is it possible to propagate a TLE over an impulsive maneuver event?

The following Python code tries to simulate an impulsive maneuver of 1 m/s, in the velocity direction, executed 30 minutes after TLE epoch:

mandv = Vector3D(1.0, 0.0, 0.0)         # delta-v in VNC orbital reference frame
inertialFrame = FramesFactory.getEME2000()
dt = 1.0                                # [s] propagation time step
manT = 30*60.0                          # [s] man time from epoch
simDur = 90*60.0                        # [s] simulatioun duration
# SPOT-5 TLE
tle_line1 = "1 27421U 02021A   02124.48976499 -.00021470  00000-0 -89879-2 0    20"
tle_line2 = "2 27421  98.7490 199.5121 0001333 133.9522 226.1918 14.26113993    62"
tle = TLE(tle_line1,tle_line2)
TLEepoch = tle.getDate()
# attitude law
attitudeOverride = LofOffset(inertialFrame,LOFType.VNC)
# maneuver triggers
triggers = DateDetector(AbsoluteDate(TLEepoch, manT))
# Impulse maneuver definition
impMan = ImpulseManeuver(triggers, attitudeOverride, mandv, isp)
# TLE propagation for simulation duration
propagator = TLEPropagator.selectExtrapolator(tle)
propagator.addEventDetector(impMan)
propagator.propagate(TLEepoch.shiftedBy(simDur))

The impulsive manuever is triggered after TLE epoch and propagation is performed from TLE epoch up to TLEepoch.shiftedBy(simDur), which is after maneuver time. I’d expect to see the maneuver effect in the final state.

But, in last line of the above code, I get the following java error:

    propagator.propagate(TLEepoch.shiftedBy(simDur))
orekit.JavaError: <super: <class 'JavaError'>, <JavaError object>>
    Java stacktrace:
org.orekit.errors.OrekitException: reset state not allowed
        at org.orekit.propagation.analytical.tle.TLEPropagator.resetIntermediateState(TLEPropagator.java:559)
        at org.orekit.propagation.analytical.AbstractAnalyticalPropagator.acceptStep(AbstractAnalyticalPropagator.java:278)
        at org.orekit.propagation.analytical.AbstractAnalyticalPropagator.propagate(AbstractAnalyticalPropagator.java:147)
        at org.orekit.propagation.AbstractPropagator.propagate(AbstractPropagator.java:186)

The message reset state not allowed sounds to me like there is something wrong with this implementation. What is it?

Is there any example of maneuver execution with TLEPropagator I could refer to?

Thank you very much

bcazabonne · November 29, 2021, 1:33pm

Hi @acapobianchi

This behaviour is normal and expected.

The state of the TLEPropagator cannot be reinitialized during an orbit propagation. Below, the code showing the exception:

/** {@inheritDoc} */
protected void resetIntermediateState(final SpacecraftState state, final boolean forward) {
    throw new OrekitException(OrekitMessages.NON_RESETABLE_STATE);
}

TLE propagation is done in Orekit by SGP4 or SDP4 models. The choice of the model depends on the orbital period of the satellite. If the period is greater than 225 minutes, the SDP4 is executed (the D of SDP4 means Deep-space).
SDP4 and SGP4 models are mean elements propagation models. They perform the orbit propagation of the satellite’s state based on the input TLE. If the satellite’s state is reinitialized during the orbit propagation, the input TLE must be reinitialized too because the computed state will not correspond to the input TLE. Otherwise, the propagation results will be wrong. Because it is not recommended to reinitialize the TLE during the propagation, Orekit throws an exception for this method.

What I can recommend you is to build another orbit propagator based on the TLE. It can be an Eckstein-Hechler orbit propagator or a numerical orbit propabator. Both propagator can reset the intermediate state during the orbit propagation. Please find an example for the Eckstein Hechler model (based on your data). Be careful that the state computed by the TLEPropagator will be expressed in TEME frame, not EME2000. It is a frame related to TLE propagation.

# SPOT-5 TLE
tle_line1 = "1 27421U 02021A   02124.48976499 -.00021470  00000-0 -89879-2 0    20"
tle_line2 = "2 27421  98.7490 199.5121 0001333 133.9522 226.1918 14.26113993    62"
tle = TLE(tle_line1,tle_line2)
TLEepoch = tle.getDate()
# TLE propagation at TLE epoch, used to convert the input TLE to a consistent spacecraft state
propagator = TLEPropagator.selectExtrapolator(tle)
# Build a spacecraft state based on the input TLE using TLE propagation model
spacecraftState = propagator.propagate(TLEepoch)
# Intialize the Eckstein-Hechler model
propagatorEck = EcksteinHechlerPropagator(spacecraftState.getOrbit(), GravityFieldFactory.getUnnormalizedProvider(6, 6))
# Create the maneuver, add the maneuver to the Eckstein-Hechler model, and perform the orbit propagation to simDur epoch

The above example can also be done by using a NumericalPropagator instead of the EcksteinHechlerPropagator propagator.

I hop this will help you.

Best regards,
Bryan

acapobianchi · November 29, 2021, 4:11pm

Thanks Bryan for your reply, it helped me!

I have modified my code, switching to the Eckstein Hechler propagator, and made this test:

isp = 280.0
mandv = Vector3D(0.0, 2000.0, 0.0)    # delta-v in VNC orbital reference frame
inertialFrame = FramesFactory.getEME2000()
dt = 60.0                               # [s] propagation time step
startDate = AbsoluteDate(2002, 5, 7, 12, 0, 0.0, TimeScalesFactory.getUTC())
manT = 60*60.0                          # [s] man time from simulation start (it's at node)
simDur = 3.5*60*60.0                        # [s] simulatioun duration
satelliteName = 'SPOT-5'
tle_line1 = "1 27421U 02021A   02124.48976499 -.00021470  00000-0 -89879-2 0    20"
tle_line2 = "2 27421  98.7490 199.5121 0001333 133.9522 226.1918 14.26113993    62"
finalDate = startDate.shiftedBy(simDur)
tle = TLE(tle_line1,tle_line2)
TLEepoch = tle.getDate()
print(tle.getLine1())
print(tle.getLine2())
# set up an attitude law dedicated to the maneuver
attitudeOverride = LofOffset(inertialFrame,LOFType.VNC)
# maneuver will start at a known date and stop after a known duration
firingDate = AbsoluteDate(startDate, manT)
print('TLE Epoch: '+TLEepoch.toString())
print('Sim Start: '+startDate.toString())
print(' Man time: '+firingDate.toString())
print(' Sim Stop: '+finalDate.toString())
triggers = DateDetector(firingDate)
# Impulse maneuver definition
impMan = ImpulseManeuver(triggers, attitudeOverride, mandv, isp)
# wrap-up everything in a propagator
propagatorTLE = TLEPropagator.selectExtrapolator(tle)
# Build a spacecraft state based on the input TLE using TLE propagation model
spacecraftState = propagatorTLE.propagate(TLEepoch)
# Intialize the Eckstein-Hechler model
propagatorEck = EcksteinHechlerPropagator(spacecraftState.getOrbit(), GravityFieldFactory.getUnnormalizedProvider(6, 6))
# Create the maneuver, add the maneuver to the Eckstein-Hechler model, and perform the orbit propagation to simDur epoch
propagatorEck.addEventDetector(impMan)
# run simulation
propagatorEck.propagate(startDate)

I can see the inclination maneuver effect that I expected! So, the maneuver is being taken into account during propagation.

Now, if I replace the maneuver direction from Normal to Velocity:

mandv = Vector3D(2000.0, 0.0, 0.0) # delta-v in VNC orbital reference frame

I get the following error in the last line (propagation):

orekit.JavaError: <super: <class 'JavaError'>, <JavaError object>>
    Java stacktrace:
org.orekit.errors.OrekitException: too large eccentricity for propagation model: e = 0.61
        at org.orekit.propagation.analytical.EcksteinHechlerPropagator$EHModel.<init>(EcksteinHechlerPropagator.java:668)
        at org.orekit.propagation.analytical.EcksteinHechlerPropagator.computeMeanParameters(EcksteinHechlerPropagator.java:491)
        at org.orekit.propagation.analytical.EcksteinHechlerPropagator.resetIntermediateState(EcksteinHechlerPropagator.java:467)
        at org.orekit.propagation.analytical.AbstractAnalyticalPropagator.acceptStep(AbstractAnalyticalPropagator.java:278)
        at org.orekit.propagation.analytical.AbstractAnalyticalPropagator.propagate(AbstractAnalyticalPropagator.java:147)
        at org.orekit.propagation.AbstractPropagator.propagate(AbstractPropagator.java:186)

I understand that this is happening because the (big) Along-Track maneuver modifies the orbit eccentricity too much, pushing it out of the limits that the Eckstein Hechler propagator can manage. In fact, from its class description I can read:

The Eckstein-Hechler model is suited for near circular orbits (e < 0.1, with poor accuracy between 0.005 and 0.1) and inclination neither equatorial (direct or retrograde) nor critical (direct or retrograde).

I’d like to use an analytical propagator that can handle maneuvers of any size and any post-maneuver orbital elements. Since I’m performing a preliminary mission analysis, I don’t have precision requirements, but, instead, I have computational speed requirements. That’s why my first guess was using TLEs and the TLEPropagator (for its speed), but, as you explained before,

the state of the TLEPropagator cannot be reinitialized during an orbit propagation

Although it is not correct, is it possible to consider the following workaround for TLE propagation and impulsive maneuvers?

propagate a TLE up to the maneuver time
convert it to SpacecraftState
applying the maneuver
get a post-maneuver SpacecraftState
generate a new TLE from post-maneuver SpacecraftState (how?)
propagate the new TLE up to the next maneuver

If the above workaround is not an option, which Orekit analytical propagator would you recommend to handle any maneuver and any orbit?

Should I consider switching to a semi-analytical propagator, like the DSSTPropagator? Which limitations does it have?

Thank you very much,
Alfredo

bcazabonne · November 29, 2021, 4:39pm

Yes … I forgot that Eckstein-Hechler model doesn’t like eccentric orbits …
In 11.1 we will add a new analytical orbit propagator which is not limited in eccentricity or inclination: the Brouwer-Lyddane model. It will be very useful for your test. Are you using the development version of Orekit or a stable version (e.g., 11.0 or 11.0.1 or 11.0.2)?

DSSTPropagator has one limitation: it doesn’t like maneuvers.
DSST is a mean elements satellite theory. It computes the osculating orbital elements by summing the short periodic terms to the integrated mean orbital elements (orb_osc = orb_mean + short_period). The mean orbital elements are computed numerically while the short periodic terms are computed analytically. That’s why the DSST is a semi-analytical orbit propagator. The fact that only the mean elements equation of motions are numerically integrated allows to use big integration steps for mean elements integration: about 12 hours.

An integration step of the DSST propagator costs a lot of computation time compare to an integration step of the numerical propagator. However, the DSST orbit propagator performs much less integration steps than a numerical one during an orbit propagation. This is why it is really faster.

As a result resetting the DSST propagator is possible, but it involves to pay a computation time cost.
What could be interesting is to test. You can test using both a DSST propagator and a numerical propagator and compare the computation time.

Best regards,
Bryan

acapobianchi · November 29, 2021, 5:23pm

Thanks Bryan,
I’m using python wrapper for Orekit 11.0.1.

The Brower-Lyddane propagator looks promising! When do you expect to release it?

For what you say, the DSST propagator doesn’t seem to fit my needs.

I don’t need a very sophisticated gravitational model. Is there the possibility to use a simple Keplerian propagator plus something to take into account the J2 effect? How should I configure it? (or should I just set up a simple numerical propagator, with gravity field of degree 2 and order 0?)

Thank you,
Alfredo

bcazabonne · November 30, 2021, 8:34am

From an optimistic point of view, by the end of the year. From a pessimistic point of view, by the end of January 2022. I’m sure that the Brouwer-Lyddane method will fully answer your problem

I think it is interesting to try.

I think the best is indeed to initialize a numerical orbit propagator with a simple propagation model. The following example shows how to initialize it with a gravity field of degree 2 and order 0 (I try to do it in Python).

integrationStep = 60
integrator      = ClassicalRungeKuttaIntegrator(integrationStep)
propagator      = NumericalPropagator(integrator)
initialState    = SpacecraftState(initialOrbit)
propagator.setInitialState(spacecraftState)
# Add a force model
gravityField = HolmesFeatherstoneAttractionModel(ecef, GravityFieldFactory.getNormalizedProvider(2, 0))
propagator.addForceModel(gravityField)
# Create the maneuver, add the maneuver to the numerical propagator using .addEventDetector() method

Do not hesitate to increase to degree and order of the gravity field. It can be increase it until the compromise between computation time and accuracy is acceptable. You can also try to add other force model. There are all available in the org.orekit.forces package.

At the end, I think the numerical propagator solution is interesting while waiting for the addition of Brouwer-Lyddane.

Best regards,
Bryan

acapobianchi · November 30, 2021, 6:26pm

Thanks for your answer, Bryan, it really helps me.

Here’s my (working) Python code:

mandv = Vector3D(1.0, 0.0, 0.0)    # delta-v in VNC orbital reference frame
orbit = KeplerianOrbit(pv, inertialFrame, startDate, Constants.EIGEN5C_EARTH_MU)
initialState = SpacecraftState(orbit, mass)
# set up propagator
orbitType = OrbitType.EQUINOCTIAL
integrator = ClassicalRungeKuttaIntegrator(dt)
NumProp = NumericalPropagator(integrator)
NumProp.setOrbitType(orbitType)
NumProp.setInitialState(initialState)
NumProp.setAttitudeProvider(LofOffset(inertialFrame, LOFType.VNC))
# Add a force model (gravity field of degree 2 and order 0)
gravityField = HolmesFeatherstoneAttractionModel(FramesFactory.getGCRF(), GravityFieldFactory.getNormalizedProvider(2, 0))
NumProp.addForceModel(gravityField)
# set up an attitude law dedicated to the maneuver
attitudeOverride = LofOffset(inertialFrame,LOFType.VNC)
# maneuver will start at a known date and stop after a known duration
# Maneuver (impulsive)
# note that ImpulseManeuver is a event detector, not a force model!
triggers = DateDetector(firingDate)
impMan = ImpulseManeuver(triggers, attitudeOverride, mandv, isp)
NumProp.addEventDetector(impMan)

Do you see any inconsistency in it? Or does it look good to you?

Thank you,
Alfredo

bcazabonne · November 30, 2021, 7:29pm

There is one error. When you initialize the Holmes-Featherstone model, you use GCRF as frame, which is an inertial frame. However, this model needs an Earth Centered, Earth Fixed (ECEF) frame. For instance you can use an ITRF frame that can be initialized like that: FramesFactory.getITRF(IERSConvention.IERS_2010, true)

Except this error, everything look good

Bryan