雷达站心坐标与gps坐标转换matlab代码

下面是雷达站心坐标与GPS坐标转换的MATLAB代码：

% 假设雷达站的经纬度为 (lat0, lon0) ，高度为 h0，目标点的经纬度为 (lat, lon)，高度为 h
% 先将经纬度转换为 ECEF 坐标系下的坐标
R = 6378137.0; % 地球半径，单位为米
f = 1/298.257223563; % 地球扁率
e2 = f*(2-f); % 第一偏心率的平方
lat0 = deg2rad(lat0); lon0 = deg2rad(lon0);
lat = deg2rad(lat); lon = deg2rad(lon);
N0 = R/sqrt(1-e2*sin(lat0)^2); % 卯酉圈曲率半径
X0 = (N0+h0)*cos(lat0)*cos(lon0);
Y0 = (N0+h0)*cos(lat0)*sin(lon0);
Z0 = (N0*(1-e2)+h0)*sin(lat0);
N = R/sqrt(1-e2*sin(lat)^2);
X = (N+h)*cos(lat)*cos(lon);
Y = (N+h)*cos(lat)*sin(lon);
Z = (N*(1-e2)+h)*sin(lat);

% 将 ECEF 坐标系下的点转换为雷达站心坐标系下的点
dx = X - X0; dy = Y - Y0; dz = Z - Z0;
T = [-sin(lon0)             cos(lon0)               0;
     -sin(lat0)*cos(lon0)   -sin(lat0)*sin(lon0)    cos(lat0);
     cos(lat0)*cos(lon0)    cos(lat0)*sin(lon0)     sin(lat0)];
P = T*[dx; dy; dz];
x = P(1); y = P(2); z = P(3);

其中，将经纬度转换为 ECEF 坐标系下的坐标需要用到地球半径、地球扁率等参数，可以参考 WGS 84 椭球体参数。转换后的坐标为三维向量（X，Y，Z），单位为米。

将 ECEF 坐标系下的点转换为雷达站心坐标系下的点，则需要用到旋转矩阵 T，其中的元素可以根据雷达站的经纬度计算得到。转换后的坐标为三维向量（x，y，z），单位为米。