% ECEF直角坐标系转LLA(WGS-84)坐标系
% Longitude 经度 （单位 度）
% Latitude 纬度
% H 海拔高度（ 相对地球表面 ） (单位 km)
% LLA = [ Longitude ; Latitude ; H ]
function LLA = Ecef2Lla( E )
a = 6378.137; % 地球赤道半径
e = 0.0818191908426214957; % 地球偏心率
b = a*sqrt(1-e^2); % 地球极半径
for i = 1:size(E,2)
p = sqrt( E(1,i)^2 + E(2,i)^2 );
theta = atan( E(3,i)*a/p/b );
e_2 = (a^2 - b^2)/b^2;
latitude = atan( ( E(3,i) + e_2*b*sin(theta)^3 )/( p - e^2*a*cos(theta)^3 ) );
longitude = atan2( E(2,i),E(1,i) );
H = p/cos(latitude) - a/sqrt(1-e^2*sin(latitude)^2);
Longitude = longitude/pi*180;
Latitude = latitude/pi*180;
LLA(:,i) = [ Longitude ; Latitude ; H ];
end
