Secure your code as it's written. Use Snyk Code to scan source code in minutes - no build needed - and fix issues immediately.
time_interval: :math:`\delta t` :attr:`datetime.timedelta`, optional
The time interval over which to test the new state (default
is 0)
Returns
-------
: float
Value of the pdf at :attr:`test_state` and :math:`t + \delta t`
Note
----
Units of mean motion must be in :math:`\mathrm{rad} s^{-1}`
"""
# First transit the prior to the current
trans_state = OrbitalState(self.transition(orbital_state,
time_interval=time_interval
),
timestamp=orbital_state.timestamp +
time_interval,
grav_parameter=orbital_state.grav_parameter)
return norm.pdf(test_state.mean_anomaly,
loc=trans_state.mean_anomaly,
scale=self.transition_noise)
:math:'lambda' is the mean longitude (rad)
Reference
---------
Broucke, R. A. & Cefola, P. J. 1972, Celestial Mechanics, Volume 5, Issue 3, pp.303-310
"""
return StateVector(np.array([[self.semimajor_axis],
[self.equinocital_h],
[self.equinocital_k],
[self.equinocital_p],
[self.equinocital_q],
[self.mean_longitude]]))
class KeplerianOrbitalState(OrbitalState):
r"""Merely a shell for the OrbitalState(coordinates='Keplerian') class, but
includes some boundary checking. As a reminder:
.. math::
X_{t_0} = [e, a, i, \Omega, \omega, \\theta]^{T} \\
where:
:math:`e` is the orbital eccentricity (unitless),
:math:`a` the semi-major axis (m),
:math:`i` the inclination (rad),
:math:`\Omega` is the longitude of the ascending node (rad),
:math:`\omega` the argument of periapsis (rad), and
:math:`\\theta` the true anomaly (rad).
"""
(:attr:`rvs()`) function instead
"""
mean_anomaly = orbital_state.mean_anomaly
mean_motion = orbital_state.mean_motion
tle_state = orbital_state.two_line_element
# TODO: Enforce the units of mean_motion are rad/s
new_mean_anomaly = np.remainder(
mean_anomaly + mean_motion * time_interval.total_seconds(),
2 * np.pi)
new_tle_state = np.insert(np.delete(tle_state, 4, 0), 4,
new_mean_anomaly, axis=0)
return OrbitalState(new_tle_state, coordinates='TLE',
timestamp=orbital_state.timestamp + time_interval,
grav_parameter=orbital_state.grav_parameter). \
cartesian_state_vector
.format(state_vector[2][0]))
if np.less(state_vector[3][0], 0.0) | np.greater(state_vector[3][0], 2 * np.pi):
raise ValueError("Longitude of Ascending Node should be between 0 and 2*pi: got {}"
.format(state_vector[3][0]))
if np.less(state_vector[4][0], 0.0) | np.greater(state_vector[4][0], 2 * np.pi):
raise ValueError("Argument of Periapsis should be between 0 and 2*pi: got {}"
.format(state_vector[4][0]))
if np.less(state_vector[5][0], 0.0) | np.greater(state_vector[5][0], 2 * np.pi):
raise ValueError("True Anomaly should be between 0 and 2*pi: got {}"
.format(state_vector[5][0]))
# And go ahead and initialise as previously
super().__init__(state_vector, coordinates='keplerian', *args, **kwargs)
class TLEOrbitalState(OrbitalState):
r"""
For the TLE state vector:
.. math::
X_{t_0} = [i, \Omega, e, \omega, n, M_0]^{T} \\
where :math:`i` the inclination (rad),
:math:`\Omega` is the longitude of the ascending node (rad),
:math:`e` is the orbital eccentricity (unitless),
:math:`\omega` the argument of perigee (rad),
:math:'M_0' the mean anomaly (rad) and
:math:`n` the mean motion (rad/[time]).
"""
def __init__(self, state_vector, *args, **kwargs):
N random samples of the state vector drawn from a normal
distribution defined by the transited mean anomaly,
:math:`M_{t_1}` and the standard deviation :math:`\epsilon`
Note
----
Units of mean motion must be in :math:`\mathrm{rad} s^{-1}`
"""
# Generate the samples
mean_anomalies = np.random.normal(0, self.transition_noise,
num_samples)
# Use the passed state, or generate a 0-state, as the mean state
if orbital_state is None:
meanstate = OrbitalState(np.zeros((6, 1)), coordinates="TLE",
timestamp=datetime(1970, 1, 1))
else:
meanstate = orbital_state
meantlestate = meanstate.two_line_element
out = np.zeros((6, 0))
for mean_anomaly in mean_anomalies:
currenttlestatev = np.remainder(meantlestate +
np.array([[0], [0], [0], [0],
[mean_anomaly], [0]]),
2*np.pi)
currentstate = TLEOrbitalState(currenttlestatev,
timestamp=orbital_state.timestamp,
grav_parameter=meanstate.
grav_parameter)
raise ValueError("Inclination should be between 0 and pi: got {}"
.format(state_vector[1][0]))
if np.less(state_vector[1][0], 0.0) | np.greater(state_vector[1][0], 2*np.pi):
raise ValueError("Longitude of Ascending Node should be between 0 and 2*pi: got {}"
.format(state_vector[2][0]))
if np.less(state_vector[3][0], 0.0) | np.greater(state_vector[3][0], 2*np.pi):
raise ValueError("Argument of Periapsis should be between 0 and 2*pi: got {}"
.format(state_vector[3][0]))
if np.less(state_vector[4][0], 0.0) | np.greater(state_vector[4][0], 2*np.pi):
raise ValueError("Mean Anomaly should be between 0 and 2*pi: got {}"
.format(state_vector[5][0]))
super().__init__(state_vector, coordinates='TLE', *args, **kwargs)
class EquinoctialOrbitalState(OrbitalState):
r"""
For the Equinoctial state vector:
.. math::
X_{t_0} = [a, h, k, p, q, \lambda]^{T} \\
where :math:`a` the semi-major axis (m),
:math:`h` is the horizontal component of the eccentricity :math:`e`,
:math:`k` is the vertical component of the eccentricity :math:`e`,
:math:`q` is the horizontal component of the inclination :math:`i`,
:math:`k` is the vertical component of the inclination :math:`i` and
:math:'lambda' is the mean longitude
"""
def __init__(self, state_vector, *args, **kwargs):
if np.less(state_vector[1][0], -1.0) | np.greater(state_vector[1][0], 1.0):