SUBROUTINE GPGELA(XI,YI,ZI,XJ,YJ,ZJ,PHII,LAMI,DAZ,XXII,ETA,DIJ,AIJ        
     1,AJI,ZIJ,ZJI,DXLG,DYLG,DZLG,DXLA,DYLA,DZLA)                               
C                                                                               
C        THIS SUBROUTINE WILL COMPUTE THE INVERSE POSITIONING PROBLEM IN        
C      3-DIMENSIONS. GIVEN THE X,Y,Z COORDINATES OF 2 POINTS, IT WILL           
C      COMPUTE THE SPATIAL DISTANCE AND THE LOCAL ASTRONOMIC AZIMUTH AND        
C      ZENITH ANGLE BETWEEN THE TWO POINTS.                                     
C                                                                               
C        INPUT:                                                                 
C              XI,YI,ZI,XJ,YJ,ZJ - CARTESIAN COORDINATES OF POINTS I AND J.     
C              PHII,LAMI - GEODETIC LATITUDE AND LONGITUDE.                     
C              DAZ - ASTRO-AZIMUTH MINUS GEODETIC AZIMUTH.                      
C              XXII,ETA - DEFLECTION OF THE VERTICAL COMPONENTS OF POINT I.     
C        OUTPUT:                                                                
C              DIJ - SPATAL DISTANCE.                                           
C              AIJ,AJI - LOCAL ASTRONOMIC AZIMUTHS FROM I TO J AND FROM J TO I. 
C              ZIJ,ZJI - ZENITH ANGLE FROM I TO J AND FROM J TO I.              
C              DXLA,DYLA,DZLA - X,Y,Z LOCAL ASTRONOMIC POINT J.                 
C              DXLG,DYLG,DZLG - X,Y,Z LOCAL GEODETIC POINT J.                   
C                                                                               
C                                           WRITTEN BY G.BOWIE,NOV.1977.        
C                                                                               
      IMPLICIT REAL*8(A-Z)                                                      
      PI=DARCOS(-1.D0)                                                          
      SP=DSIN(PHII)                                                             
      SL=DSIN(LAMI)                                                             
      CP=DCOS(PHII)                                                             
      CL=DCOS(LAMI)                                                             
      DXLG=-(XJ-XI)*SP*CL-(YJ-YI)*SP*SL+(ZJ-ZI)*CP                              
      DYLG=-(XJ-XI)*SL+(YJ-YI)*CL                                               
      DZLG=(XJ-XI)*CP*CL+(YJ-YI)*CP*SL+(ZJ-ZI)*SP                               
      DXLA=DXLG-DAZ*DYLG-XXII*DZLG                                              
      DYLA=DAZ*DXLG+DYLG-ETA*DZLG                                               
      DZLA=XXII*DXLG+ETA*DYLG+DZLG                                              
      DIJ=DSQRT((XJ-XI)**2+(YJ-YI)**2+(ZJ-ZI)**2)                               
      AIJ=DATAN2(DYLA,DXLA)                                                     
      AJI=AIJ+PI                                                                
      ZIJ=DARCOS(DZLA/DIJ)                                                      
      ZJI=PI-ZIJ                                                                
      DD=DIJ                                                                    
      AA=AIJ                                                                    
      ZZ=ZIJ                                                                    
      RETURN                                                                    
      END