For symmetry:

procedure DistanceToLatLong(const lat0, long0, dist, bearing: double;

var lat1, long1: double);

var

R: double = 6371000.0; // earth radius in meters

d, x, y: double;

begin

d := dist/R; // angular distance

lat1 := ArcSin(Sin(lat0) * Cos(d) + Cos(lat0) * Sin(d) * Cos(bearing));

x := Cos(d) - Sin(lat0) * Sin(lat1);

y := Sin(bearing) * Sin(d) * Cos(lat0);

long1 := long0 + ArcTan2(y, x);

end;

Hi!

This needs some speedup - there are situations where it is heavily used.

* Some values are computed twice with sin or cos -

you can compute theses variables first

* To speedup you can do math.sincos where sinus and cosinus are

computed in one job.

* We should take the latest exact value for the earth radius from WGS84 = 6378137.0 meters

procedure DistanceToLatLong(const lat0, long0, dist, bearing: double;

var lat1, long1: double);

var

//R: double = 6371000.0; // earth radius in meters

R : double = 6378137.0; // WGS84

d, x, y: double;

cosd, sind, coslat0, sinlat0 : double;

begin

d := dist/R; // angular distance

math.sincos (d,sind,cosd);

math.sincos (lat0,sinlat0, coslat0);

lat1 := ArcSin(sinlat0 * cosd + coslat0 * sind * cos(bearing) );

//lat1 := ArcSin(Sin(lat0) * Cos(d) + Cos(lat0) * Sin(d) * Cos(bearing));

x := cosd - sinlat0 * sin(lat1);

//x := Cos(d) - Sin(lat0) * Sin(lat1);

y := sin(bearing) * sind * coslat0;

//y := Sin(bearing) * Sin(d) * Cos(lat0);

long1 := long0 + ArcTan2(y, x);

end;

And don't forget:

The earth is an ellipsoid - at the poles the radius is about 21 km less than at the equator.

(So it's easierer for the Trolls from Midearth to the surface ....)

This can lead to 0.5% deviation.

Winni