Compute Jacobian through finite difference calculation
Parameters:
fun (function handle) – A (non-linear) transition function
Must be of the form “y = fun(x)”, where y can be a scalar or numpy.ndarray of shape (Nd, 1) or (Nd,)
x (State) – A state with state vector of shape (Ns, 1)
Approximate a given distribution to a Gaussian, using a
deterministically selected set of sigma points.
Parameters:
state (State) – A state object capable of returning a StateVector of
shape (Ns, 1) representing the Gaussian mean and a
CovarianceMatrix of shape (Ns, Ns) which is the
covariance of the distribution
alpha (float, optional) – Spread of the sigma points. Typically 1e-3.
(default is 1)
beta (float, optional) – Used to incorporate prior knowledge of the distribution
2 is optimal if the state is normally distributed.
(default is 2)
kappa (float, optional) – Secondary spread scaling parameter
(default is calculated as 3-Ns)
Returns:
list of length 2*Ns+1 – An list of States containing the locations of the sigma points.
Note that only the state_vector attribute in these
States will be meaningful. Other quantities, like covar
will be inherited from the input and don’t really make sense
for a sigma point.
numpy.ndarray of shape (2*Ns+1,) – An array containing the sigma point mean weights
numpy.ndarray of shape (2*Ns+1,) – An array containing the sigma point covariance weights
Apply the Unscented Transform to a set of sigma points
Apply f to points (with secondary argument points_noise, if available),
then approximate the resulting mean and covariance. If sigma_noise is
available, treat it as additional variance due to additive noise.
Parameters:
sigma_points_states (StateVectors of shape (Ns, 2*Ns+1)) – An array containing the locations of the sigma points
mean_weights (numpy.ndarray of shape (2*Ns+1,)) – An array containing the sigma point mean weights
covar_weights (numpy.ndarray of shape (2*Ns+1,)) – An array containing the sigma point covariance weights
fun (function handle) – A (non-linear) transition function
Must be of the form “y = fun(x,w)”, where y can be a scalar or numpy.ndarray of shape (Ns, 1) or (Ns,)
covar_noise (CovarianceMatrix of shape (Ns, Ns), optional) – Additive noise covariance matrix
(default is None)
points_noise (numpy.ndarray of shape (Ns, 2*Ns+1,), optional) – points to pass into f’s second argument
(default is None)
theta (Union[float, np.ndarray]) – Rotation angle specified as a real-valued number or an np.ndarray of reals. The rotation angle is positive if the rotation is in the clockwise direction when viewed by an observer looking down the x-axis towards the origin. Angle units are in radians.
theta (Union[float, np.ndarray]) – Rotation angle specified as a real-valued number or an np.ndarray of reals. The rotation angle is positive if the rotation is in the clockwise direction when viewed by an observer looking down the y-axis towards the origin. Angle units are in radians.
theta (Union[float, np.ndarray]) – Rotation angle specified as a real-valued number or an np.ndarray of reals. The rotation angle is positive if the rotation is in the clockwise direction when viewed by an observer looking down the z-axis towards the origin. Angle units are in radians.
means (StateVector, StateVectors, np.ndarray of shape (num_dims, num_components)) – The means of GM components
covars (np.ndarray of shape (num_components, num_dims, num_dims) or list of np.ndarray of shape (num_dims, num_dims)) – Covariance matrices of the GM components
Calculates and returns the (3D) axis rotation matrix given a vector of
three angles:
[roll, pitch/elevation, yaw/azimuth]
Order of rotations is in reverse: yaw, pitch, roll (z, y, x)
This is the rotation matrix that implements the rotations that convert the input
angle_vector to match the x-axis.
Parameters:
angle_vector (numpy.ndarray of shape (3, 1): the rotations)
Calculates and returns the (3D) axis rotation matrix given a vector of
three angles:
[roll, pitch/elevation, yaw/azimuth]
Order of rotations is roll, pitch, yaw (x, y, z)
This is the rotation matrix that implements the rotations that convert a vector aligned to the
x-axis to the input angle_vector.
Parameters:
angle_vector (numpy.ndarray of shape (3, 1): the rotations)
Evaluate the cubature points for an input Gaussian state. This is done under the assumption
that the input state is \(\mathcal{N}(\mathbf{\mu}, \Sigma)\) of dimension \(n\). We
calculate the square root of the covariance (via Cholesky factorization), and find the cubature
points, \(X\), as,
\[ \begin{align}\begin{aligned}\Sigma &= S S^T\\X_i &= S \xi_i + \mathbf{\mu}\end{aligned}\end{align} \]
for \(i = 1,...,2n\), where \(\xi_i = \sqrt{ \alpha n} [\pm \mathbf{1}]_i\) and
\([\pm \mathbf{1}]_i\) are the positive and negative unit vectors in each dimension. We
include a scaling parameter \(\alpha\) to allow the selection of cubature points closer to
the mean or more in the tails, as a potentially useful free parameter.
Parameters:
state (GaussianState) – A Gaussian state with mean and covariance
alpha (float, optional) – scaling parameter allowing the selection of cubature points closer to the mean (lower
values) or further from the mean (higher values)
Returns:
Cubature points (as a StateVectors of dimension \(n \times 2n\))
where \(Q\) is an optional additive noise matrix. The scaling parameter \(\alpha\)
allow the for cubature points closer to the mean or more in the tails,
Parameters:
cubature_points (StateVectors) – Cubature points (as a StateVectors of dimension \(n \times 2n\))
covar_noise (CovarianceMatrix of shape (Ns, Ns), optional) – Additive noise covariance matrix
(default is None)
alpha (float, optional) – scaling parameter allowing the nomination of cubature points closer to the mean (lower
values) or further from the mean (higher values)
Undertakes the cubature transform as described in [3]
Given a Gaussian distribution, calculates the set of cubature points using
gauss2cubature(), then passes these through the given function and reconstructs the
Gaussian using cubature2gauss(). Returns the mean, covariance, cross covariance and
transformed cubature points. This instance includes a scaling parameter \(\alpha\), not
included in the reference detailed above, which allows for the selection of cubature points
closer to, or further from, tne mean.
Parameters:
state (GaussianState) – A Gaussian state with mean and covariance
fun (function handle) – A (non-linear) transition function
Must be of the form “y = fun(x,w)”, where y can be a scalar or numpy.ndarray of shape (Ns, 1) or (Ns,)
covar_noise (CovarianceMatrix of shape (Ns, Ns), optional) – Additive noise covariance matrix
(default is None)
points_noise (numpy.ndarray of shape (Ns, 2*Ns+1,), optional) – points to pass into f’s second argument
(default is None)
alpha (float, optional) – scaling parameter allowing the selection of cubature points closer to the mean (lower
values) or further from the mean (higher values)
Calculate the universal anomaly via Newton’s method. Algorithm 3.3 in [1].
Parameters:
o_state_vector (StateVector, StateVectors) – The orbital state vector formed as
\([r_x, r_y, r_z, \dot{r}_x, \dot{r}_y, \dot{r}_z]^T\)
delta_t (timedelta) – The time interval over which to estimate the universal anomaly
grav_parameter (float, optional) – The universal gravitational parameter. Defaults to that of the
Earth, \(3.986004418 \times 10^{14} \ \mathrm{m}^{3} \
\mathrm{s}^{-2}\)
precision (float, optional) – For Newton’s method, the difference between new and old estimates of the universal anomaly
below which the iteration stops and the answer is returned, (default = 1e-8)
max_iterations (float, optional) – Maximum number of iterations allowed in while loop (default = 1e5)
Calculate the Lagrangian coefficients, f and g, and their time derivatives, by way of the
universal anomaly and the Stumpff functions [2].
Parameters:
o_state_vector (StateVector, StateVectors) – The (Cartesian) orbital state vector,
\([r_x, r_y, r_z, \dot{r}_x, \dot{r}_y, \dot{r}_z]^T\)
delta_t (timedelta) – The time interval over which to calculate
grav_parameter (float, optional) – The universal gravitational parameter. Defaults to that of the
Earth, \(3.986004418 \times 10^{14} \ \mathrm{m}^{3} \
\mathrm{s}^{-2}\). Note that the units of time must be seconds.
precision (float, optional) – Precision to which to calculate the universalanomaly() (default = 1e-8). See the doc
section for that function
max_iterations (float, optional) – Maximum number of iterations in determining universal anomaly (default = 1e5)
Returns:
The Lagrange coefficients, \(f, g, \dot{f}, \dot{g}\), in that order.
\[\begin{split}X = [e, a, i, \Omega, \omega, \theta]^{T} \\\end{split}\]
where:
\(e\) is the orbital eccentricity (unitless),
\(a\) the semi-major axis (m),
\(i\) the inclination (rad),
\(\Omega\) is the longitude of the ascending node (rad),
\(\omega\) the argument of periapsis (rad), and
\(\theta\) the true anomaly (rad)
grav_parameter (float, optional) – Standard gravitational parameter \(\mu = G M\). The default is \(3.986004418
\times 10^{14} \mathrm{m}^3 \mathrm{s}^{-2}\)
Returns:
Orbital state vector as \([r_x, r_y, r_z, \dot{r}_x, \dot{r}_y, \dot{r}_z]\)
This function performs linear interpolation on a StateMutableSequence. The function
has two slightly different forms:
If an individual datetime is inputted for the variable times then a
State is returned corresponding to times.
If a list of datetime is inputted for the variable times then a
StateMutableSequence is returned with the states in the sequence corresponding to
times.
Note
This function does not extrapolate. Times outside the range of the time range of sms
are discarded and warning is given. If all times values are outside the time range of
sms then an IndexError is raised.
Unique states for each time are required for interpolation. If there are multiple states with
the same time in sms the later state in the sequence is used.
For Track inputs the metadatas is removed as it can’t be interpolated.