Keplerian covariance propagation
[kep_t2,cov_t2]=CL_dsp_kepCovPropa(t1,kep_t1,cov_t1,t2 [,mu])
Performs the propagation of orbital elements and associated covariance matrix using a Keplerian model:
Given a set of Keplerian orbital elements kep_t1 and the associated covariance matrix (cov_t1) at time t1, the function computes the orbital elements kep_t2 and the associated covariance matrix (cov_t2) at time t2.
Note: cov_t2(:,:,i) corresponds to kep_t2(:,i).
Initial time [days] (1x1)
Keplerian elements at time t1 [sma;ecc;inc;pom;raan;anm] (6x1)
Covariance matrix at time t1 (6x6)
Final time [days] (1xN)
(optional) Gravitational constant [m^3/s^2] (default value is %CL_mu)
Keplerian elements at t2 [sma;ecc;inc;pom;raan;anm] (6xN)
Covariance matrix at time t2 (6x6xN)
CNES - DCT/SB
// Example 1 : Covariance propagation t1 = 21892.7549949074091; kep_t1 = [7029796.3;0.0036;CL_deg2rad([97.9;47.2;277.8;-59.0])]; cov_t1 = [ 8.3D+06 0.90471 -0.25296 -153.815 -0.08734 142.552 0.90471 1.3D-07 -2.6D-08 -9.9D-06 -1.5D-08 8.7D-06 -0.25296 -2.6D-08 2.7D-07 7.3D-06 2.6D-08 -6.9D-06 -153.815 -9.9D-06 7.3D-06 0.00519 1.1D-07 -0.00498 -0.08734 -1.5D-08 2.6D-08 1.1D-07 7.4D-09 5.3D-09 142.552 8.7D-06 -6.9D-06 -0.00498 5.3D-09 0.00479 ] t2 = t1 : 10/86400 :t1+10/1440; [kep_t2,cov_t2] = CL_dsp_kepCovPropa(t1,kep_t1,cov_t1,t2); // Example 2 : // Propagation of pointing errors seen from a ground station // This involves several changes of reference frame // The jacobian needed for covariance transformation is computed // Note: covariances and jacobian are modeled as hypermatrix (6,6,n) // If the previous values were in J2000 : // Transformation keplerian to cartesian (position, velocity) : [pos_car,vel_car,jacobkp2car] = CL_oe_kep2car(kep_t2); cov_car = jacobkp2car * cov_t2 * jacobkp2car'; // Frame transformation : J2000 ==> terresetrial reference frame [pos_ter,vel_ter,jacobJ20002ter] = CL_fr_J20002ter(t2,pos_car,vel_car); // Ground station definition KRU = [CL_deg2rad(-52.64);CL_deg2rad(5.1);200]; // Terrestrial to topocentric north : [pos_car,vel_car,jacobter2topo] = CL_fr_ter2topoN(KRU,pos_ter,vel_ter); jacobter2topo_t=hypermat([6 6 size(kep_t2,2)]); for j=1:size(kep_t2,2) jacobter2topo_t(:,:,j)=jacobter2topo(1:6,1:6); end // Change of coordinate(cartesian to spherical) [pos_sph,vel_sph,jacobcar2sph] = CL_co_car2sph(pos_car,vel_car); // Global Jacobian computation and covariance transformation jacobiennetot=jacobcar2sph*jacobter2topo_t*jacobJ20002ter*jacobkp2car; covariancetopoNsph=jacobiennetot*cov_t2*jacobiennetot'; // Conversions for plot err_az=matrix(covariancetopoNsph(1,1,:),1,size(covariancetopoNsph(6,6,:),3)); err_elev=matrix(covariancetopoNsph(2,2,:),1,size(covariancetopoNsph(6,6,:),3)); clear jacobter2topo jacobter2topo_t jacobter2topo covariancetopoNcar //-----------PLOTS------------------------- scf(); plot2d((t2-t2(1))*24,CL_rad2deg(sqrt(err_az)),2) plot2d((t2-t2(1))*24,CL_rad2deg(sqrt(err_elev)),5) ax=gca(); CL_g_stdaxes(); plots = CL_g_select(ax,"Polyline"); plots.thickness = 2; ax.title.text= "Pointing error"; ax.x_label.text= "Time (h)"; ax.y_label.text= "Pointing error (deg)"; ax.font_size=2; ax.title.font_size=2; ax.x_label.font_size=2; ax.y_label.font_size=2; legends(["azimuth" "elevation"],[2,5],1); | ![]() | ![]() |