Similarity Transformation - from Cartesian to TLE Covariance#
One obvious application of dsgp4, where the use of its differentiability can be useful, is for covariance transformation and propagation (see tutorial covariance_propagation.ipynb).
In this notebook, we will discuss how to use dsgp4 for the first case.
import torch
import dsgp4
import numpy as np
import warnings
warnings.filterwarnings('ignore', category=DeprecationWarning)
Covariance Transformation#
We start with the state corresponding to a given TLE at a certain time: this is in Cartesian TEME reference frame.
Then, we assume that we have the covariance matrix associated to this state in the RTN frame (see, for instance, “Fundamental of Astrodynamics and Applications” by Vallado, for a discussion and theoretical background on these frames). We want to:
first convert the RTN covariance into TEME
then leverage
dsgp4to transform the covariance matrix from position and velocity coordinates to TLE elements, leveraging the similarity transformation:
(1)#\[\begin{equation}
P_{y} = m P_x m^T\text{,}
\end{equation}\]
where:
(2)#\[\begin{equation}
m_{ij}=\dfrac{\partial y_i}{\partial x_j}
\end{equation}\]
once the above is done, we are left with a covariance directly in TLE parameters, which can be used for several applications. An example, is to directly use the covariance to generate perturbed TLEs (from their mean) and feed them to \(\partial\textrm{SGP4}\) to propagate the noisy observations
As always, let’s start by constructing the TLE:
lines=[]
lines.append("0 COSMOS 2251 DEB")
lines.append("1 34427U 93036RU 22068.94647328 .00008100 00000-0 11455-2 0 9999")
lines.append("2 34427 74.0145 306.8269 0033346 13.0723 347.1308 14.76870515693886")
my_tle=dsgp4.tle.TLE(lines)
tle_elements=dsgp4.initialize_tle(my_tle,with_grad=True);
Let’s now construct the initial RTN covariance matrix, and the initial TEME state:
def matrix_rtn(cov):
cov_matrix_rtn=np.array([cov['sigma_rr'],
cov['sigma_rt'],
cov['sigma_rn'],
cov['sigma_vr'],
cov['sigma_vt'],
cov['sigma_vn']])
cov_matrix_rtn=np.diag(cov_matrix_rtn)
return cov_matrix_rtn
covariance_diagonal={"sigma_rr":10,
"sigma_rt":100,
"sigma_rn":10,
"sigma_vr":1e-3,
"sigma_vt":1e-1,
"sigma_vn":1e-3}
cov_matrix_rtn=matrix_rtn(covariance_diagonal)
print(f"RTN covariance matrix is: {cov_matrix_rtn}")
state_teme=dsgp4.propagate(my_tle,0.)
print(f"Initial state in Cartesian TEME is: {state_teme}")
RTN covariance matrix is: [[1.e+01 0.e+00 0.e+00 0.e+00 0.e+00 0.e+00]
[0.e+00 1.e+02 0.e+00 0.e+00 0.e+00 0.e+00]
[0.e+00 0.e+00 1.e+01 0.e+00 0.e+00 0.e+00]
[0.e+00 0.e+00 0.e+00 1.e-03 0.e+00 0.e+00]
[0.e+00 0.e+00 0.e+00 0.e+00 1.e-01 0.e+00]
[0.e+00 0.e+00 0.e+00 0.e+00 0.e+00 1.e-03]]
Initial state in Cartesian TEME is: tensor([[ 4.1941e+03, -5.6009e+03, 4.2616e-02],
[ 1.6578e+00, 1.2582e+00, 7.2697e+00]], grad_fn=<TransposeBackward0>)
Let’s now transform the covariance from RTN to Cartesian TEME:
def rotation_matrix(state):
"""
Computes the UVW rotation matrix.
Args:
state (`numpy.array`): numpy array of 2 rows and 3 columns, where
the first row represents position, and the second velocity.
Returns:
`numpy.array`: numpy array of the rotation matrix from the cartesian state.
"""
r, v = state[0], state[1]
u = r / np.linalg.norm(r)
w = np.cross(r, v)
w = w / np.linalg.norm(w)
v = np.cross(w, u)
return np.vstack((u, v, w))
def from_cartesian_to_rtn(state, cartesian_to_rtn_rotation_matrix=None):
"""
Converts a cartesian state to the RTN frame.
Args:
state (`numpy.array`): numpy array of 2 rows and 3 columns, where
the first row represents position, and the second velocity.
cartesian_to_rtn_rotation_matrix (`numpy.array`): numpy array of the rotation matrix from the cartesian state. If None, it is computed.
Returns:
`numpy.array`: numpy array of the RTN state.
"""
# Use the supplied rotation matrix if available, otherwise compute it
if cartesian_to_rtn_rotation_matrix is None:
cartesian_to_rtn_rotation_matrix = rotation_matrix(state)
r, v = state[0], state[1]
r_rtn = np.dot(cartesian_to_rtn_rotation_matrix, r)
v_rtn = np.dot(cartesian_to_rtn_rotation_matrix, v)
return np.stack([r_rtn, v_rtn]), cartesian_to_rtn_rotation_matrix
###### RTN->Cartesian Rotation ######
state_rtn, cartesian_to_rtn_rotation_matrix = from_cartesian_to_rtn(state_teme.detach().numpy())
transformation_matrix_cartesian_to_rtn = np.zeros((6,6))
transformation_matrix_cartesian_to_rtn[0:3, 0:3] = cartesian_to_rtn_rotation_matrix
transformation_matrix_cartesian_to_rtn[3:,3:] = cartesian_to_rtn_rotation_matrix
C_teme = np.matmul(np.matmul(transformation_matrix_cartesian_to_rtn.T, cov_matrix_rtn),transformation_matrix_cartesian_to_rtn)
print(f"Cartesian TEME covariance matrix is: {C_teme}")
Cartesian TEME covariance matrix is: [[1.43678346e+01 3.27091412e+00 1.90611714e+01 0.00000000e+00
0.00000000e+00 0.00000000e+00]
[3.27091412e+00 1.24494699e+01 1.42742252e+01 0.00000000e+00
0.00000000e+00 0.00000000e+00]
[1.90611714e+01 1.42742252e+01 9.31826956e+01 0.00000000e+00
0.00000000e+00 0.00000000e+00]
[0.00000000e+00 0.00000000e+00 0.00000000e+00 5.80461804e-03
3.59800553e-03 2.09672885e-02]
[0.00000000e+00 0.00000000e+00 0.00000000e+00 3.59800553e-03
3.69441685e-03 1.57016478e-02]
[0.00000000e+00 0.00000000e+00 0.00000000e+00 2.09672885e-02
1.57016478e-02 9.25009651e-02]]
#quick check to confirm the transformation was correct:
C_rtn2 = np.matmul(np.matmul(transformation_matrix_cartesian_to_rtn, C_teme),transformation_matrix_cartesian_to_rtn.T)
print(f"Was the RTN->Cartesian transformation correct? {np.allclose(C_rtn2,cov_matrix_rtn)}")
Was the RTN->Cartesian transformation correct? True
C_teme
array([[1.43678346e+01, 3.27091412e+00, 1.90611714e+01, 0.00000000e+00,
0.00000000e+00, 0.00000000e+00],
[3.27091412e+00, 1.24494699e+01, 1.42742252e+01, 0.00000000e+00,
0.00000000e+00, 0.00000000e+00],
[1.90611714e+01, 1.42742252e+01, 9.31826956e+01, 0.00000000e+00,
0.00000000e+00, 0.00000000e+00],
[0.00000000e+00, 0.00000000e+00, 0.00000000e+00, 5.80461804e-03,
3.59800553e-03, 2.09672885e-02],
[0.00000000e+00, 0.00000000e+00, 0.00000000e+00, 3.59800553e-03,
3.69441685e-03, 1.57016478e-02],
[0.00000000e+00, 0.00000000e+00, 0.00000000e+00, 2.09672885e-02,
1.57016478e-02, 9.25009651e-02]])
Now, we do the similarity transformation, leveraging the partials of the TLE parameters, w.r.t. the state in Cartesian TEME. For a detailed discussion (with examples) on how to compute partials with dsgp4, check out the tutorial tle_propagation.ipynb
#let's construct the 6x9 matrix of partials:
dx_dtle = torch.zeros((6,9))
for i in range(6):
tle_elements.grad=None
state_teme.flatten()[i].backward(retain_graph=True)
dx_dtle[i,:] = tle_elements.grad
#let's print them to screen:
print(dx_dtle)
tensor([[ 0.0000e+00, 0.0000e+00, 0.0000e+00, -4.7924e+03, 1.5394e+03,
1.0479e+00, 1.5463e+03, -4.3390e+04, 5.6009e+03],
[ 0.0000e+00, 0.0000e+00, 0.0000e+00, 4.9523e+03, 1.1618e+03,
-3.3683e+00, 1.1736e+03, 5.7875e+04, 4.1941e+03],
[ 0.0000e+00, 0.0000e+00, 0.0000e+00, -3.0267e+03, 6.7277e+03,
-3.9559e+00, 6.7716e+03, -1.4301e+02, 0.0000e+00],
[ 0.0000e+00, 0.0000e+00, 0.0000e+00, 2.6362e+00, -4.5385e+00,
-5.8152e+00, -4.5521e+00, 8.6044e+00, -1.2582e+00],
[ 0.0000e+00, 0.0000e+00, 0.0000e+00, -1.3649e-01, 6.0581e+00,
-4.3574e+00, 6.0789e+00, 6.4300e+00, 1.6578e+00],
[ 0.0000e+00, 0.0000e+00, 0.0000e+00, 7.0851e+00, -5.4891e-03,
2.0766e+00, -4.7732e-05, 3.7622e+01, 0.0000e+00]])
#now, we can apply the similarity transformation, and obtain the 9x9 covariance matrix in TLE elements:
Cov_tle=np.matmul(np.matmul(np.linalg.pinv(dx_dtle),C_teme),np.linalg.pinv(dx_dtle.T))
#quick check to confirm we are doing things correctly:
C_xyz_2=np.matmul(np.matmul(dx_dtle,Cov_tle),dx_dtle.T)
print(f"Was the TLE->Cartesian transformation correct? {np.allclose(C_xyz_2, C_teme)}")
Was the TLE->Cartesian transformation correct? True