function [x,y,z] = wgs2xyz(lam,phi,h)

% WGS2XYZ	converts lam(longitude) phi(latitude) ellipsoidal coordinates
%		from WGS-84 to ECEF cartesian coordinates
%		[x,y,z] = wgs2xyz(longitude,latitude,h)

% semimajor and semiminor axis for WGS-84
a = 6378137.0;
b = 6356752.314;
f = 1.0/298.257222101;
ee = 2*f - f^2;

% degrees to radians
lam = lam*pi/180;
phi = phi*pi/180;

% radius of curvature in prime vertical
N = a / sqrt(1-(sin(phi)).^2*ee);
%N = a^2 / sqrt((cos(phi)).^2*a^2 + (sin(phi)).^2*b^2);

x = cos(phi).*cos(lam).*(N+h);
y = cos(phi).*sin(lam).*(N+h);
z = sin(phi).*(N*(b^2/a^2) + h);