Custom Event Condition

In this notebook, we will:

1. Create two custom event conditions

   1. Beta angle condition (real value based condition)

   2. Eclipse condition (boolean value based condition)

import numpy as np

from ostk.physics import Environment
from ostk.physics.environment.object.celestial import Earth
from ostk.physics.environment.object.celestial import Sun
from ostk.physics.environment.gravitational import Earth as EarthGravitationalModel
from ostk.physics.environment.atmospheric import Earth as EarthAtmosphericModel
from ostk.physics.environment.magnetic import Earth as EarthMagneticModel
from ostk.physics.time import Instant
from ostk.physics.time import Duration
from ostk.physics.coordinate import Frame
from ostk.physics.coordinate import Position
from ostk.physics.coordinate import Velocity

from ostk.astrodynamics.trajectory import State
from ostk.astrodynamics.trajectory import Propagator
from ostk.astrodynamics.trajectory.state import NumericalSolver
from ostk.astrodynamics.event_condition import BooleanCondition
from ostk.astrodynamics.event_condition import RealCondition

Setup the initial state and propagation interval

start_instant = Instant.J2000()
initial_state = State(
    instant=start_instant,
    position=Position.meters([7000000.0, 0.0, 0.0], Frame.GCRF()),
    velocity=Velocity.meters_per_second(
        [0.0, 5335.865450622126, 5335.865450622126], Frame.GCRF()
    ),
)

earth = Earth.from_models(
    EarthGravitationalModel(EarthGravitationalModel.Type.Spherical),
    EarthMagneticModel(EarthMagneticModel.Type.Undefined),
    EarthAtmosphericModel(EarthAtmosphericModel.Type.Undefined),
)
environment = Environment(Instant.J2000(), [earth, Sun.default()])

Eclipse condition

A boolean condition is defined by:

• A name

• A criteria that defines what kind of root crossing is desired

• The evaluation function which accepts 2 arguments

  • a state_vector

  • a time

• If the condition is to be inverted

def in_eclipse(state):
    # global environment
    environment.set_instant(state.get_instant())
    return environment.is_position_in_eclipse(state.get_position())


eclipse_condition = BooleanCondition(
    "Eclipse Condition", BooleanCondition.Criterion.AnyCrossing, in_eclipse, False
)

Beta angle condition

A real condition is similarly defined by:

• A name

• A criteria that defines what kind of root crossing is desired

• The evaluation function

• The target value

sun = Sun.default()


def get_beta_angle(state):
    instant = state.get_instant()
    sun_vector = sun.get_position_in(Frame.GCRF(), instant).get_coordinates()
    orbit_plane_normal = np.cross(
        state.get_position().get_coordinates(), state.get_velocity().get_coordinates()
    )
    return 90 - np.rad2deg(
        np.arccos(
            np.clip(
                np.dot(
                    sun_vector / np.linalg.norm(sun_vector),
                    orbit_plane_normal / np.linalg.norm(orbit_plane_normal),
                ),
                -1.0,
                1.0,
            )
        )
    )


beta_angle_condition = RealCondition(
    "Beta Angle Condition", RealCondition.Criterion.AnyCrossing, get_beta_angle, 21.17
)

propagator = Propagator.from_environment(
    NumericalSolver.default_conditional(), environment
)

Let's find the first time the satellite is in eclipse

solution = propagator.calculate_state_to_condition(
    initial_state, initial_state.get_instant() + Duration.hours(24.0), eclipse_condition
)
print(
    f"Satellite eclipse state [{in_eclipse(solution.state)}] at {solution.state.get_instant().to_string()}"
)

Satellite eclipse state [True] at 2000-01-01 12:09:00.876.134.713 [UTC]

Let's now find when the beta angle reaches a target value

solution = propagator.calculate_state_to_condition(
    initial_state,
    initial_state.get_instant() + Duration.hours(24.0),
    beta_angle_condition,
)
print(
    f"Satellite beta angle [{get_beta_angle(solution.state)}] at {solution.state.get_instant().to_string()}"
)

Satellite beta angle [21.16999999999848] at 2000-01-01 18:52:09.813.823.449 [UTC]

---