请用matlab语言编写一个轨道六根数与速度位置转换的程序

以下是一个简单的程序，实现了轨道六根数与速度位置之间的转换：

%% 轨道六根数与速度位置转换程序

% 输入轨道六根数
a = input('请输入半长轴a（km）：');
e = input('请输入离心率e：');
i = input('请输入轨道倾角i（度）：');
Omega = input('请输入升交点赤经Ω（度）：');
w = input('请输入近地点幅角ω（度）：');
M = input('请输入平近点角M（度）：');

% 常数定义
mu = 3.986004418e5; % 地球引力常数（km^3/s^2）

% 计算轨道周期T
T = 2*pi*sqrt(a^3/mu);

% 计算平近点角速度n
n = 2*pi/T;

% 计算偏近点角E
E0 = M;
E = E0;
for i = 1:10
E = M + e*sin(E);
end

% 计算真近点角v
v = 2*atan(sqrt((1+e)/(1-e))*tan(E/2));

% 计算轨道半径r
r = a*(1 - e*cos(E));

% 计算速度大小V
V = sqrt(mu*(2/r - 1/a));

% 计算速度分量
Vr = V*sin(v);
Vt = V*cos(v);

% 计算位置分量
x = r*cos(Omega)*cos(w+v) - r*sin(Omega)*sin(w+v)*cos(i);
y = r*sin(Omega)*cos(w+v) + r*cos(Omega)*sin(w+v)*cos(i);
z = r*sin(i)*sin(w+v);

% 输出结果
fprintf('速度大小V = %f km/s\n', V);
fprintf('速度分量Vr = %f km/s\n', Vr);
fprintf('速度分量Vt = %f km/s\n', Vt);
fprintf('位置分量x = %f km\n', x);
fprintf('位置分量y = %f km\n', y);
fprintf('位置分量z = %f km\n', z);