function [K] = calcLQR(I) Q = eye(6) * I(3, 3); Mu = 398600; r0 = [6678; 0; 0]; k1 = (I(2,2) - I(3,3)) / I(1,1); k2 = (I(1,1) - I(3,3)) / I(2,2); k3 = (I(2,2) - I(1,1)) / I(3,3); n = sqrt(Mu/(norm(r0)^3)); F = -2 * n^2 * [4 * k1, 0, 0; 0, 3 * k2, 0; 0, 0, k3]; G = n * [0, 0, (1 - k1); 0, 0, 0; (k3 - 1), 0, 0]; Iinv = [1 / I(1,1), 0, 0; 0, 1 / I(2,2), 0; 0, 0, 1 / I(3,3)]; A = [zeros(3,3), 0.5 * eye(3); F, G]; B = [zeros(3,3); Iinv]; R = eye(3)*10^-6; Rinv = R^-1; S = icare(A, B, Q, R); K = Rinv * B' * S;