12 include
'geodesic.inc'
14 double precision a, f, lat1, lon1, azi1, lat2, lon2, azi2, s12,
15 + dummy1, dummy2, dummy3, dummy4, dummy5
25 read(*, *, end=90, err=90) lat1, lon1, lat2, lon2
26 if (abs(lat1) .gt. 90 .or. abs(lat2) .gt. 90)
then
28 15
format(1x,
'lat1 and lat2 must be in [-90,90]')
31 call invers(a, f, lat1, lon1, lat2, lon2,
32 + s12, azi1, azi2, omask,
33 + dummy1, dummy2, dummy3, dummy4 , dummy5)
34 print 20, azi1, azi2, s12
35 20
format(1x, f20.15, 1x, f20.15, 1x, f19.10)
program geodinverse
A simple program to solve the inverse geodesic problem.