blh2xyz.m

这是国外关于卫星导航方面一书的源代码

function [x,y,z]=blh2xyz(p,l,h,ell)
% Umwandlung: Ellipsoidische in Kartesische Koordinaten
%
% Wie:     [x,y,z]=blh2xyz(p,l,h)
%          [x,y,z]=blh2xyz(p,l,h,ellips)
% Eingabe: B,L in Degree.Decimale, H in Meter
%          ELLIPS ist Textstring: 'wgs84' (default), 'grs80', oder 'bessel'
% Ausgabe: X,Y,Z in Meter
%
% Formeln aus Hofmann-Wellenhof u.a.(1994):"GPS in der Praxis"
% Keine negativen Hoehen!
%
% (c)iapg - 11/96. Keine Gewaehr, da keine Abfragen.

rho=180/pi;
p=p/rho;
l=l/rho;

if nargin < 4, ell = 'wgs84'; end	% WGS84 default
ell = lower(deblank(ell));
if strcmp(ell,'bessel')			% Bessel-Ellipsoid 
 a=6377397.155;
 f=1/299.1528128;
elseif strcmp(ell,'wgs84') 		% WGS84-Ellipsoid
 a=6378137;
 f=1/298.257223563;
elseif strcmp(ell,'grs80') 		% GRS80-Ellipsoid
 a=6378137;
 f=1/298.257222101;
else
 error('Nicht-definiertes Ellipsoid')
end

b=(1-f)*a;
e1q=((a^2)-(b^2))/(a^2);
n=a./sqrt(1-e1q.*(sin(p)).^2);

x=(n+h).*cos(p).*cos(l);
y=(n+h).*cos(p).*sin(l);
z=(n.*(1-e1q)+h).*sin(p);