function point2=km2lonlat(point1,d,angle)

% KM2LONLAT : gives a [lon,lat] arrival point as a function of the initial point, a distance and an angle (orientation)
% from http://www.movable-type.co.uk/scripts/latlong.html
%
% point2=km2lonlat(point1,d,angle)
% point1, point2 as [lon,lat] ou [lon,lat] columns, d in km, angle in deg :
% 0 if going full north, 90 eastward, 180 southward, -90 westward (ie angles counted clockwise starting north)
%
%  Monique, mai 2007, MBARI

R = 6378137/1000;		% Earth radius in km
nb_pts=size(point1,1);
point2=ones(nb_pts,2)*NaN;

for i=1:nb_pts

	point_encours=deg2rad(point1(i,:));
	angle=deg2rad(angle);
	lon1=point_encours(1);
	lat1=point_encours(2);

	alpha=d*1000/R;
	lat2 = lat1 + alpha*cos(angle);
	delta_phi=log(tan(lat2/2+pi/4)/tan(lat1/2+pi/4));
	if (angle==deg2rad(90)) | (angle==deg2rad(-90)), q=cos(lat1);
		else, q=(lat2-lat1)/delta_phi;, end
	delta_lon=alpha*sin(angle)/q;
	lon2=mod(lon1+delta_lon+pi,2*pi)-pi;

	%lat2 = asin(sin(lat1)*cos(d/R) + cos(lat1)*sin(d/R)*cos(angle));
	%lon2 = lon1 + atan2(sin(angle)*sin(d/R)*cos(lat1), cos(d/R)-sin(lat1)*sin(lat2));

	point2(i,:)=rad2deg([lon2,lat2]);

end

return
