用matlab写一个获取卫星位置的程序

好的，下面是一个简单的 MATLAB 程序，可以获取指定时间下卫星的位置：

% 卫星编号
prn = 1;

% GPS 周数
gps_week = 2120;

% GPS 时间（秒）
gps_time = 3600;

% 计算卫星位置
[x, y, z] = satpos(prn, gps_week, gps_time);

% 显示卫星位置
fprintf('卫星 %d 在 GPS 时间 %d 周 %d 秒的位置为：\n', prn, gps_week, gps_time);
fprintf('X = %.2f km\n', x);
fprintf('Y = %.2f km\n', y);
fprintf('Z = %.2f km\n', z);

function [x, y, z] = satpos(prn, gps_week, gps_time)
    % 卫星轨道参数
    a = 26559700.0; % 卫星轨道长半轴（m）
    e = 0.0010335; % 卫星轨道偏心率
    i = deg2rad(54.7355); % 卫星轨道倾角（rad）
    omega = deg2rad(34.4228); % 升交点赤经（rad）
    w = deg2rad(32.2883); % 近地点角距（rad）
    M0 = deg2rad(0); % 平近点角（rad）
    n = sqrt(398600.4418/(a^3)); % 平均角速度（rad/s）

    % 计算时间差
    tk = gps_time - 0;
    if tk < -302400
        tk = tk + 604800;
    elseif tk > 302400
        tk = tk - 604800;
    end

    % 计算平近点角
    Mk = M0 + n*tk;

    % 迭代计算偏近点角
    E0 = Mk;
    for i = 1:10
        E1 = Mk + e*sin(E0);
        if abs(E1 - E0) < 1e-12
            break;
        end
        E0 = E1;
    end

    % 计算真近点角
    vk = atan2(sqrt(1-e^2)*sin(E0), cos(E0)-e);

    % 计算轨道平面上的卫星位置
    xk = a*(cos(E0)-e);
    yk = a*sqrt(1-e^2)*sin(E0);

    % 计算升交点到卫星的距离
    rk = sqrt(xk^2 + yk^2);
    uk = atan2(yk, xk);
    ok = omega + uk;

    % 计算卫星在地心惯性系下的位置
    x = rk*(cos(ok)*cos(vk)-sin(ok)*sin(vk)*cos(i));
    y = rk*(sin(ok)*cos(vk)+cos(ok)*sin(vk)*cos(i));
    z = rk*sin(vk)*sin(i);
end

这个程序假设卫星的轨道是一个 Kepler 椭圆，根据 Kepler 定律计算卫星的位置。程序中给出了一个例子，可以自己修改卫星编号、GPS 周数和 GPS 时间来获取不同卫星在不同时间下的位置。注意，这个程序计算的卫星位置是在地心惯性系下的，需要将其转换到地心固定系下才能和接收机位置计算距离。