Secure your code as it's written. Use Snyk Code to scan source code in minutes - no build needed - and fix issues immediately.
class LinearGaussian(MeasurementModel, LinearModel, GaussianModel):
r"""This is a class implementation of a time-invariant 1D
Linear-Gaussian Measurement Model.
The model is described by the following equations:
.. math::
y_t = H_k*x_t + v_k,\ \ \ \ v(k)\sim \mathcal{N}(0,R)
where ``H_k`` is a (:py:attr:`~ndim_meas`, :py:attr:`~ndim_state`) \
matrix and ``v_k`` is Gaussian distributed.
"""
noise_covar = Property(CovarianceMatrix, doc="Noise covariance")
@property
def ndim_meas(self):
"""ndim_meas getter method
Returns
-------
:class:`int`
The number of measurement dimensions
"""
return len(self.mapping)
def matrix(self, **kwargs):
"""Model matrix :math:`H(t)`
def covar(self, time_interval=timedelta(seconds=1)):
r"""Return the transition covariance matrix
Parameters
----------
time_interval: :math:`\delta t` :attr:`datetime.timedelta`, optional
The time interval over which to test the new state (default is 1
second)
Returns
-------
: CovarianceMatrix
The transition covariance matrix
"""
return CovarianceMatrix(self.noise * time_interval.total_seconds())
time_interval_sec = time_interval.total_seconds()
# Only leading terms get calculated for speed.
covar = sp.array(
[[sp.power(time_interval_sec, 5) / 20,
sp.power(time_interval_sec, 4) / 8,
sp.power(time_interval_sec, 3) / 6],
[sp.power(time_interval_sec, 4) / 8,
sp.power(time_interval_sec, 3) / 3,
sp.power(time_interval_sec, 2) / 2],
[sp.power(time_interval_sec, 3) / 6,
sp.power(time_interval_sec, 2) / 2,
time_interval_sec]]
) * self.noise_diff_coeff
return CovarianceMatrix(covar)
k = self.damping_coeff
q = self.noise_diff_coeff
dt = time_interval.total_seconds()
exp_kdt = sp.exp(-k*dt)
exp_2kdt = sp.exp(-2*k*dt)
q11 = q*(dt - 2/k*(1 - exp_kdt) + 1/(2*k)*(1 - exp_2kdt))/(k**2)
q12 = q*((1 - exp_kdt)/k - 1/(2*k)*(1 - exp_2kdt))/k
q22 = q*(1 - exp_2kdt)/(2*k)
covar = sp.array([[q11, q12],
[q12, q22]])
return CovarianceMatrix(covar)
dt = time_interval_sec
N = self.constant_derivative
if N == 1:
covar = sp.array([[dt**3 / 3, dt**2 / 2],
[dt**2 / 2, dt]])
else:
Fmat = self.matrix(time_interval, **kwargs)
Q = sp.zeros((N + 1, N + 1))
Q[N, N] = 1
igrand = Fmat @ Q @ Fmat.T
covar = sp.zeros((N + 1, N + 1))
for l in range(0, N + 1):
for k in range(0, N + 1):
covar[l, k] = (igrand[l, k]*dt / (1 + N**2 - l - k))
covar *= self.noise_diff_coeff
return CovarianceMatrix(covar)
-------
:class:`stonesoup.types.state.CovarianceMatrix` of shape\
(:py:attr:`~ndim_state`, :py:attr:`~ndim_state`)
The process noise covariance.
"""
time_interval_sec = time_interval.total_seconds()
base_covar = sp.array([[sp.power(time_interval_sec, 3) / 3,
sp.power(time_interval_sec, 2) / 2],
[sp.power(time_interval_sec, 2) / 2,
time_interval_sec]])
covar_list = [base_covar*sp.power(self.noise_diff_coeffs[0], 2),
base_covar*sp.power(self.noise_diff_coeffs[1], 2)]
covar = sp.linalg.block_diag(*covar_list)
return CovarianceMatrix(covar)
class GaussianStatePrediction(Prediction, GaussianState):
""" GaussianStatePrediction type
This is a simple Gaussian state prediction object, which, as the name
suggests, is described by a Gaussian distribution.
"""
class GaussianMeasurementPrediction(MeasurementPrediction, GaussianState):
""" GaussianMeasurementPrediction type
This is a simple Gaussian measurement prediction object, which, as the name
suggests, is described by a Gaussian distribution.
"""
cross_covar = Property(CovarianceMatrix,
doc="The state-measurement cross covariance matrix",
default=None)
def __init__(self, state_vector, covar, timestamp=None,
cross_covar=None, *args, **kwargs):
if(cross_covar is not None
and cross_covar.shape[1] != state_vector.shape[0]):
raise ValueError("cross_covar should have the same number of \
columns as the number of rows in state_vector")
super().__init__(state_vector, covar, timestamp,
cross_covar, *args, **kwargs)
class ParticleStatePrediction(Prediction, ParticleState):
"""ParticleStatePrediction type
time_interval_sec = time_interval.total_seconds()
# Only leading terms get calculated for speed.
covar = sp.array(
[[sp.power(time_interval_sec, 5) / 20,
sp.power(time_interval_sec, 4) / 8,
sp.power(time_interval_sec, 3) / 6],
[sp.power(time_interval_sec, 4) / 8,
sp.power(time_interval_sec, 3) / 3,
sp.power(time_interval_sec, 2) / 2],
[sp.power(time_interval_sec, 3) / 6,
sp.power(time_interval_sec, 2) / 2,
time_interval_sec]]
) * self.noise_diff_coeff
return CovarianceMatrix(covar)
"""
coordinates = Property(
strg, default="cartesian",
doc="The parameterisation used on initiation. Acceptable values "
"are 'Cartesian', 'Keplerian', 'TLE', or 'Equinoctial'. All"
"other inputs will return errors"
)
grav_parameter = Property(
float, default=3.986004418e14,
doc=r"Standard gravitational parameter :math:`\mu = G M` in units of "
r":math:`\mathrm{m}^3 \mathrm{s}^{-2}`")
covar = Property(
CovarianceMatrix, default=None,
doc="The covariance matrix. Care should be exercised in that its coordinate"
"frame isn't defined, and output will be highly dependant on which"
"parameterisation is chosen."
)
_eanom_precision = Property(
float, default=1e-8,
doc="Precision used for the stopping point in determining eccentric "
"anomaly from mean anomaly"
)
metadata = Property(
{}, default={}, doc="Dictionary containing metadata about orbit")
def __init__(self, state_vector, *args, **kwargs):
r"""Can be initialised in a number of different ways according to
def __init__(self, state_vector, covar, *args, **kwargs):
covar = CovarianceMatrix(covar)
super().__init__(state_vector, covar, *args, **kwargs)
if self.state_vector.shape[0] != self.covar.shape[0]:
raise ValueError(
"state vector and covar should have same dimensions")