Looking for Umbra and Penumbra Event Detection Example on Existing Ephemeris

Hello,

First time posting here so I hope this is the appropriate place to request some support. I am looking to define Umbra and Penumbra events on various ephemerides. I am using the python wrapper for orekit and it’s been working great.

It’s been very helpful to find example usages for other analysis I’ve completed with orekit but I can’t find any for implementing the event detection capability for an existing ephemeris. It appears I need to define an Ephemeris object which involves a list of spacecraft states, which come from the Orbit class, which is generated from TimeStampedPVCoordinates. Then I can potentially apply an EclipseDetector via the addEventDetector with the withPenumbra() and withUmbra() sub-methods.

Can anyone with more experience using the event detection functions on an existing ephemeris provide any feedback on if this is the correct way to go? Or is there a more direct way for this? An example usage would be much appreciated.

Thanks in advance!
Regards,
Justin

If you have an EphemerisSegment from loading an ephemeris file use an EphemerisSegmentPropagator. Or you could turn your PV into SpacecraftState by adding a Frame. E.g. new SpacecraftState(new CartesianOrbit(pv, frame)) Or if you want to implement your own for a special case just copy EphemerisSegmentPropagator and adapt to your needs. It is fairly short.

2 Likes

Per what Evan said, here is an example in Java of doing that sort of conversion in a streamlined way. It really just behaves like any other propagator in terms of adding event detectors etc.:

List<TimeStampedPVCoordinates> ephemeris = buildEphemeris();
Frame frame = FramesFactory.getGCRF();
double gm = CelestialBodyFactory.getEarth().getGM();
List<SpacecraftState> states = ephemeris.stream().map(e -> new SpacecraftState(new CartesianOrbit(e, frame ,gm))).collect(Collectors.toList());
int interpolationPoints = 5;
Ephemeris ephemerisPropagator = new Ephemeris(states, interpolationPoints);
ApsideDetector apsideDetector2 = new ApsideDetector(states.get(0).getOrbit());
ephemerisPropagator.addEventDetector(apsideDetector2);
SpacecraftState periapsisCrossingState2 = ephemerisPropagator.propagate(states.get(0).getDate().shiftedBy(43200.0));
System.out.println(periapsisCrossingState2.getDate() + " TA = " + new KeplerianOrbit(periapsisCrossingState2.getOrbit()).getAnomaly(PositionAngle.TRUE));
1 Like

Thanks to both for the ideas. This is what I have so far in psuedo-code (using Python wrapper):

GCRF = FramesFactory.getGCRF()
UTC = TimeScalesFactory.getUTC()
orekit_states = ArrayList()

for state in ephemeris_data:
    epoch = AbsoluteDate(*time data from state[0] input here*, UTC)
    rx, ry, rz = state[1:4] 
    vx, vy, vz = state [4:7]
    pos = Vector3D(rx, ry, rz)
    vel = Vector3D(vx, vy, vz)
    coords_GCRF = PVCoordinates(pos, vel)
    abs_coords_GCRF = AbsolutePVCoordinates(GCRF, epoch, coords_GCRF)
    orekit_states.add(SpacecraftState(abs_coords_GCRF))

orekit_ephemeris = Ephemeris(orekit_states, 10)

That is how I saw I could form the Ephemeris object. Now it appears I must use an ephemerisPropagator or maybe an EphemerisSegmentPropagator with an attached event detector to get the umbra/penumbra events I am looking for.

Also thinking my eclipse detector will look something like this:

eclipse_detector = EclipseDetector(sun, Constants.WGS84_EARTH_EQUATORIAL_RADIUS, earth)

I will attempt this and let you know my results. Thanks!

Another question: I assume the Ephemeris Propagator is not actually propagating states? I would hope it is simply using a given Ephemeris to interpolate between points. I am not interested in an actual force model propagation, just want it to find the umbra/penumbra events given an Ephemeris generated outside of orekit.

I now have the Ephemeris object defined and an eclipse detector added. I then “propagate” to the end of the ephemeris. The result is a SpacecraftState object with AbsolutePVACoordinates I can access with .getAbsPVA().

orekit_ephemeris = Ephemeris(orekit_states, 10)
eclipse_detector = EclipseDetector(sun, Constants.WGS84_EARTH_EQUATORIAL_RADIUS, earth)
orekit_ephemeris.addEventDetector(eclipse_detector)
result = orekit_ephemeris.propagate(orekit_ephemeris.getMaxDate()) 

The method above using getOrbit() to find the detected events does not work since the resulting object does not have an orbit defined. Does this mean I can’t used AbsolutePVCoordinates and instead must work with Orbit to allow me to access the detected events?

Edit: I think what is happening is it’s only propagating to the first event and stopping. I need to determine how to allow it to not stop at the events and continue logging them for the entire ephemeris length. Then I need to determine how to access the event log - any tips for this? Not having any luck finding this in docs…

Hi,

In Orekit ‘Propagator’ means a trajectory that supports event detection and step handling. Orekit supports integrated propagators and analytical propagators. Ephemeris is in the latter category. It interpolates some quantities (orbit) and can compute others (attitude) depending on how you configure it. You can take a look at the source code to see what it is actually doing.

If your orbit has a central body it is more common to use an Orbit instead of AbsolutePVCoordinates. In Orekit an Orbit is an osculating orbit, most of the time. That means it is just position and velocity (or another equivalent set of six numbers) w.r.t. a central body at a point in time.

Easiest way is to use an EventHandler, either a custom one you write or perhaps RecordAndContinue. The default event handler stops propagation when exiting eclipse.[1] See [2] for an example of printing events. You could also use EventsLogger.

[1] EclipseDetector (ORbit Extrapolation KIT 10.1 API)
[1] Orekit Tutorials – Propagation

Thanks again for the helpful info, the support is much appreciated. Looking at the Orbit objects… it appears I cannot give it a list of states (ephem) as it only takes a single state. This means it must be propagating which is what I do not want as I am looking for interpolation of an existing ephemeris to get umbra events and thus I must stick with an ArrayList() of SpacecraftStates to define the Ephemeris object.

I will take your advice and try to figure out how to implement the RecordAndContinue function (or custom as needed) to allow for proper events logging. Thanks for pointing me in the right direction!

I believe I have finially figured it out. My code is now:

umbra_logger = EventsLogger()
penumbra_logger = EventsLogger()
orekit_ephemeris = Ephemeris(orekit_sc_states, 10)  # orekit_sc_states is an ArrayList() of SpacecraftStates using AbsolutePVCoordinates() from pre-existing ephemeris
umbra_detector = EclipseDetector(sun, Constants.WGS84_EARTH_EQUATORIAL_RADIUS, earth).withUmbra().withHandler(RecordAndContinue())
penumbra_detector = EclipseDetector(sun, Constants.WGS84_EARTH_EQUATORIAL_RADIUS, earth).withPenumbra().withHandler(RecordAndContinue())
orekit_ephemeris.addEventDetector(umbra_detector)
orekit_ephemeris.addEventDetector(umbra_logger.monitorDetector(umbra_detector))
orekit_ephemeris.addEventDetector(penumbra_detector)
orekit_ephemeris.addEventDetector(penumbra_logger.monitorDetector(penumbra_detector))
result = orekit_ephemeris.propagate(orekit_ephemeris.getMaxDate())
umbra_events = umbra_logger.getLoggedEvents()
penumbra_events = penumbra_logger.getLoggedEvents()

And to access the events:

# get first event position
umbra_events.get(0).getState().getPVCoordinates().getPosition()
# get first event epoch
umbra_events.get(0).getDate()
# get inc/dec status
umbra_events.get(0).isIncreasing()  # True or False

I am a bit confused on how to interpret whether each event is either an eclipse entry or exit event. I would’ve assumed that if isIncreasing() = True it would mean an entry and vice versa for False…but compared to another umbra source it suggested that False = entry and True = exit. Is this correct? I see that in the tutorial you provided the usage is such that isIncreasing() = True means an event is beginning (for a visibility event). Can you help me understand how isIncreasing() correlates to entry/exit?

The increasing/decreasing is about the shape of the function not the beginning/entry point. So you have to look at the g function to see if increasing is the beginning or ending of the event depending on whether you are trying to trigger on an up or a downswing of the curve. Per the documentation, “This function becomes negative when entering the region of shadow and positive when exiting.”

1 Like

Yes I did see the g function as well as getTotalEclipse() previously, thanks for redirecting me to those. It appears I need to couple both the g function and isIncreasing() to directly determine an exit/entry event.

For example:

In [18]: umbra_detector.g(umbra_events.get(0).getState())                                                                                                                                                          
Out[18]: -4.5382885546351587e-13

which suggests that the event is an umbra entrance event since it is very slightly negative (really is the zero point). However, use of the direction:

In [20]: umbra_events.get(0).isIncreasing()                                                                                                                                                                        
Out[20]: True

demonstrates the value is slightly negative/zero and about to become positive since the above command yielded true, thus it is actually an eclipe exit event (this agrees with other eclipse event sources) per the documentation in the g function. Would you agree with my assessment?

Thanks again!

Yes. isIncreasing() indicates if it is an entry or exit event. Sign of g function at a detected event is not significant. Event may be offset from the true root in either direction by up to the tolerance you set in the event detector.

1 Like

Got it, isIncreasing() = True indicates the g function is increasing at the event epoch, a positive value will occur next/soon, and thus it’s an eclipse exit event per the docs (and vice versa for the False case.) Thanks again for all the help!