Mike.Riddle@f27.n285.z1.fidonet.org (Mike Riddle) (06/27/90)
The following routines were mentioned by Jim Riddle in a recent
Digest article. The response has already been substantial; however,
Jim has no easy way of sending files through the Internet.
Since several people have already responded, perhaps you could run
them in the Digest or mention they are available wherever you put them
in the archives?
C TO CONVERT EITHER EARTH-CENTERED TO LAT/LON OR VICE VERSA
SUBROUTINE M2CONV (OPT,NUMPTS,LATS,LONS,XLIST,YLIST,ZLIST)
INTEGER*4 OPT,NUMPTS,I
REAL*8 LATS(NUMPTS),LONS(NUMPTS),THETA,PHI,DTR
REAL*8 XLIST(NUMPTS),YLIST(NUMPTS),ZLIST(NUMPTS)
PI = 3.141592653589793238462643
DTR = PI/180.
IF (OPT .EQ. 1) THEN
DO FOR I = 1,NUMPTS
THETA = LONS(I) * DTR
PHI = (90.0D+0 - LATS(I)) * DTR
XLIST(I) = DCOS(THETA) * DSIN(PHI)
YLIST(I) = DSIN(THETA) * DSIN(PHI)
ZLIST(I) = DCOS(PHI)
ENDDO
ELSEIF (OPT .EQ. 2) THEN
DO FOR I = 1,NUMPTS
LATS(I) = 90.0D+0 - DACOS(ZLIST(I)) / DTR
IF (XLIST(I) .EQ. 0.0D+0) THEN
IF (YLIST(I) .LT. 0.0D+0) THEN
LONS(I) = -90.0D+0
ELSE
LONS(I) = 90.0D+0
ENDIF
ELSE
LONS(I) = DATAN(YLIST(I)/XLIST(I)) / DTR
IF (XLIST(I) .LT. 0.0D+0) THEN
IF (YLIST(I) .LT. 0.0D+0) THEN
LONS(I) = LONS(I) - 180.0D+0
ELSE
LONS(I) = LONS(I) + 180.0D+0
ENDIF
ENDIF
ENDIF
ENDDO
ENDIF
END
SUBROUTINE LLAXYZ(LAT,LON,ALT,X,Y,Z)
C SUBROUTINE TO CALCULATE ECI FROM LAT/LON/ALT
C
C RE IS RADIUS OF THE EARTH
C RE = 3437.75 NAUTICAL MILES, MORE OR LESS
C 1 NAUTICAL MILE = 1852 METERS (EXACT)
C 1 STATUTE MILE = 1609.344 METERS (EXACT)
C SUBROUTINE _DOES_ REQUIRE ALTITUDE -- CHECK THE LOCAL AIRPORT (!)
C INPUT REQUIRED IS LAT, LON, ALT OF POINT AND RETURNS X, Y, Z
C
REAL LAT,LON,ALT
RE = 3437.75
C RE = 3437.75 * 1.852 FOR KILOMETERS
C RE = 3437.75 * 1852/1609.344 FOR STATUTE MILES
C ON MOST MACHINES YOU MAY WANT TO GO DOUBLE PRECISION, BY THE WAY
RANGE = RE + ALT
CLAT = COS(LAT)
SLAT = SIN(LAT)
CLON = COS(LON)
SLON = SIN(LON)
Z = RANGE * SLAT
X = RANGE*CLAT * CLON
Y = RANGE*CLAT * SLON
RETURN
END
C CONVERTS LAT AND LONG TO EUCLIDEAN COORDINATES
SUBROUTINE M2EUCL
REAL*4 DEGRAD, PI, ECLX(N), ECLY(N), ECLZ(N)
REAL*4 LOND(N), LATD(N)
INTEGER*4 NUMPTS,I
PI = 3.141592653589793238462643
DEGRAD = PI/180.
DO FOR I = 1,NUMPTS
ECLX(I) = COS(LOND(I)*DEGRAD)*COS(LATD(I)*DEGRAD)
ECLY(I) = SIN(LOND(I)*DEGRAD)*COS(LATD(I)*DEGRAD)
ECLZ(I) = SIN(LATD(I)*DEGRAD)
ENDDO
RETURN
END
SUBROUTINE M2DIST(XPT,YPT,ZPT,NUMLST)
C
INTEGER*4 NUMLST,N,INDEX
C CALCULATE DISTANCES (GREAT CIRCLE) FROM A LIST OF POINTS TO A FOCAL
C POINT
C XPT, YPT, ZPT ARE X, Y, Z COORDINATES OF THE FOCAL POINT
C XLST, YLST, ZLST ARE LISTS OF X, Y AND Z COORDINATES
C RADIUS IS RADIUS OF THE EARTH
REAL*4 XPT,YPT,ZPT,CSQRED,CTHETA,XLST(N),YLST(N),ZLST(N),DLST(N)
DO FOR INDEX = 1, NUMLST
CSQRED = (XLST(INDEX) - XPT)**2 + (YLST(INDEX) - YPT)
1 **2 + (ZLST(INDEX) - ZPT)**2
CTHETA = 1.0 - CSQRED/2.0
DLST(INDEX) = ACOS(CTHETA) * RADIUS
ENDDO
RETURN
END
Ybbat (DRBBS) 8.9 v. 3.11 r.3
[1:285/27@fidonet] The Inns of Court 402/593-1192 (1:285/27.0)
--- Through FidoNet gateway node 1:16/390
Mike.Riddle@f27.n285.z1.fidonet.org