| C GCIRCLE: Great Circle Navigation January 1982
C
C This program calculates the great circle distance between
C a place of departure and a destination. It also calculates
C the initial course for sailing the great circle. It is based
C on the method and equations in Budlong.
C
REAL*8 DEG,DIST,HC,LAT,LATD,LONG,LONGD,M,MA,MIN,PI,X,Y,Z
C
TYPE 5
5 FORMAT (1X)
C
C Enter departure latitude in degrees and minutes
C
30 TYPE 95
95 FORMAT (1X,'DEPARTURE LATITUDE ',$)
ACCEPT *, DEG, MIN
LAT=DEG+MIN/60.0D00
C
C Enter departure longitude (assumed west) in degrees and minutes
C
TYPE 90
90 FORMAT (1X,'DEPARTURE LONGITUDE ',$)
ACCEPT *, DEG, MIN
LONG=DEG+MIN/60.0D00
TYPE 5
C
C Enter destination latitude in degrees and minutes
C
TYPE 10
10 FORMAT (1X,'DESTINATION LATITUDE ',$)
ACCEPT *, DEG, MIN
LATD=DEG+MIN/60.0D00
C
C Enter destination longitude in degrees and minutes
C
TYPE 20
20 FORMAT (1X,'DESTINATION LONGITUDE ',$)
ACCEPT *, DEG, MIN
LONGD=DEG+MIN/60.0D00
TYPE 5
C
C Find meridian angle (MA)
C If (MA .LE. 0) then course between 0 and 180
C
MA=LONG-LONGD
C
C Convert degrees to radians
C
PI=3.141592653589793D00
LAT=LAT*PI/180.0D00
LATD=LATD*PI/180.0D00
M=MA*PI/180.0D00
C
C Hc computation
C
X=SIN(LAT)*SIN(LATD)
Y=COS(LAT)*COS(LATD)*COS(M)
HC=ASIN(X+Y)
C
C Great circle distance computation
C
HC=HC*180.0D00/PI
DIST=(90.0D00-HC)*60.0D00
C
C Z computation
C
X=COS(LAT)*TAN(LATD)
Y=SIN(LAT)*COS(M)
Z=ATAN(SIN(M)/(X-Y))
Z=Z*180.0D00/PI
C
C Determine actual course -- necessary because of
C implementation of ATAN
C
IF (MA .EQ. 0.0 .AND. LATD .GT. LAT) Z=0.0D00
IF (MA .EQ. 0.0 .AND. LATD .LT. LAT) Z=180.0D00
IF (MA .GT. 0.0 .AND. Z .GE. 0.0) Z=Z
IF (MA .GT. 0.0 .AND. Z .LT. 0.0) Z=Z+180.0D00
IF (MA .LT. 0.0 .AND. Z .GE. 0.0) Z=Z+180.0D00
IF (MA .LT. 0.0 .AND. Z .LT. 0.0) Z=Z+360.0D00
C
C Results
C
TYPE *, 'DISTANCE ',DIST,'NM'
TYPE *, 'COURSE ',Z,'DEGREES'
TYPE 5
TYPE 5
TYPE 5
GOTO 30
END
|
| A good authority is Bowditch - H.O.9 from the U.S. Government
Printing Office. Any fair sized library will have a copy. The
formula has just enough special cases to preclude a full
explanation (I'm working over TYMNET which has a connection
half-life of about three minutes) and it depends on the Earth
being a sphere, which it isn't. If you'll go for the tables of
Meridional parts (derivable from an approximation of the Earth as
an oblate spheroid, which it also isn't) you can do somewhat
better.
Unless you've got a lot of money riding on this, the formula
should be enough. Bowditch is the sort of reference that
discusses corrections needed because gravity isn't straight up
and down (it usually isn't) - makes for interesting reading if
you're a Phineas Fogg type ...
|