function ENUVel = xyzvel2enu( Lat, Lon, XVel, YVel, ZVel )

dRAD2DEG = 180 / pi;

Lat = Lat ./ dRAD2DEG;
Lon = Lon ./ dRAD2DEG;

sinlat = sin( Lat ) ;
coslat = cos( Lat ) ;
sinlon = sin( Lon ) ;
coslon = cos( Lon ) ;

% East/North/Up transforms  
utz =         sinlat ;
ce  =         coslat ;
et1 = -1.0 .* sinlon ;
et2 =         coslon ;

nt1 = -1.0 .* utz .* et2 ;
nt2 =         utz .* et1 ;
utx =         ce  .* et2 ;
uty = -1.0 .* ce  .* et1 ;

%Compute ENU velocity.
ENUVel(:,1) = (et1 .* XVel) + (et2 .* YVel);
ENUVel(:,2) = (nt1 .* XVel) + (nt2 .* YVel) + (ce  .* ZVel);
ENUVel(:,3) = (utx .* XVel) + (uty .* YVel) + (utz .* ZVel);


