[GRASSLIST:5809] Re: geodesic distance

H Bowman hamish_nospam at yahoo.com
Wed Mar 19 20:41:18 EST 2003


> is there a module available in GRASS to measure geodesic distances in 
> a projected location (eqc in my case)?

FYI, there's:
d.geodesic (lat/lon only)
d.measure
d.rhumbline

but I don't think any of them are what you are really looking for.


For an algorithm, see:
http://mathworld.wolfram.com/GreatCircle.html

this might get you started:
data = [ lon1 lat1;  lon2 lat2 ];
%lon (western hemisphere is 360 - lon!)
lam1 = data(1,1) * (2*pi/360);          % convert to radians
lam2 = data(2,1) * (2*pi/360);

%lat
del1 = data(1,2) * (2*pi/360);
del2 = data(2,2) * (2*pi/360);

LatDist = ((data(1,2) - data(2,2) ) *111.191);
          % convert to km (111.191km/degree lat)

r = 6367.444;   % mean radius of the earth (km)

R1dotR2 = cos(del1) * cos(del2) * cos(lam1-lam2) + sin(del1) * sin(del2);

d = r * acos( R1dotR2 );

% mercator heading!!!!! A trip from 0 lon, 80N to 180 lon, 80N will
% tell you to go East not North!!! Only good locally!
theta = acos(LatDist/d);



Hamish




More information about the grass-user mailing list