matlab_gps_空间直角坐标(XYZ)转换为大地坐标(BLH)

function [B, L, H]=XYZ2BLH(X,Y,Z) % Target: % 将WGS84下的直角坐标 转换为 WGS84下的大地坐标 % Input： % X Y Z ---- WGS84下的直角坐标(m) % output： % B L H ---- WGS84下的大地坐标：大地纬度(°)、大地经度(°)、大地高(m) % 注： 北纬为正，东经为正 % 试验数据： % * SITE X Y Z || lat(B) lon(L) Height(H) % Z > 0 ---- X > 0, Y > 0 ---- OPMT: 4202777.38 171367.99 4778660.18 48.836 2.335 124.2 % X > 0, Y < 0 ---- ALGO: 918129.4 -4346071.2 4561977.8 45.959 -78.071 202.0 % X < 0, Y > 0 ---- BJFS:-2148743.8 4426641.28 4044655.99 39.609 115.89 87.4 % X < 0, Y < 0 ---- INEG:-1260435.67 -5788547.37 2360340.13 21.856 -102.284 1888.347 % Z < 0 ---- X > 0, Y > 0 ---- MAL2: 4865385.69 4110717.189 -331137.64 -2.996 40.194 -20.4 % X > 0, Y < 0 ---- COPO: 1907041 -5337379. -2916335. -27.38 -70.338 479.2 % X < 0, Y > 0 ---- ALIC:-4052051.77 4212836.22 -2545106.03 -23.67 133.89 603.3581 % X < 0, Y < 0 ---- ISPA:-1881703.51 -5359979.77 -2890599.26 -27.125 -109.34 112.49 % XYZ2BLH(4202777.38,171367.99,4778660.18) % XYZ2BLH(918129.4,-4346071.2,4561977.8) % XYZ2BLH(-2148743.8,4426641.28,4044655.99) % XYZ2BLH(-1260435.67,-5788547.37,2360340.13) % XYZ2BLH(4865385.69,4110717.189,-331137.64) % XYZ2BLH(1907041,-5337379,-2916335) % XYZ2BLH(-4052051.77,4212836.22,-2545106.03) % XYZ2BLH(-1881703.51,-5359979.77,-2890599.26) a=6378137; % WGS84椭球长半轴(m) e=sqrt(0.00669437999013); % WGS84椭球扁率 phi=Z/sqrt(X^2+Y^2); R=sqrt(X^2+Y^2+Z^2); f=inline('atan(phi*(1+a*e^2*sin(B)/Z/sqrt(1-e^2*(sin(B)^2))))-B','B','phi','a','e','Z'); B=fzero(f,0,optimset('Display','off'),phi,a,e,Z); if (Z<0) B=abs(B)*(-1); end W=sqrt(1-e^2*(sin(B)^2)); N=a/W; % 椭球的卯酉圈半径 L=atan(Y/X); if (X<0) & (Y>0) L=L+pi; end if (X<0) & (Y<0) L=L-pi; end H=sqrt(X^2+Y^2)/cos(B)-N; B=B*180/pi; L=L*180/pi; %BLH=[B L H]';