      SUBROUTINE XYZGEO(XYZ,NUMAT,NA,NB,NC,DEGREE,GEO)
      implicit double precision (A-H,O-Z)
      include 'MOP6FILE'
      dimension XYZ(3,*), NA(*), NB(*), NC(*), GEO(3,*)
***********************************************************************
*
*   XYZGEO CONVERTS COORDINATES FROM CARTESIAN TO INTERNAL.
*
*     ON INPUT XYZ  = ARRAY OF CARTESIAN COORDINATES
*              NUMAT= NUMBER OF ATOMS
*              NA   = NUMBERS OF ATOM TO WHICH ATOMS ARE RELATED
*                     BY DISTANCE
*              NB   = NUMBERS OF ATOM TO WHICH ATOMS ARE RELATED
*                     BY ANGLE
*              NC   = NUMBERS OF ATOM TO WHICH ATOMS ARE RELATED
*                     BY DIHEDRAL
*
*    ON OUTPUT GEO  = INTERNAL COORDINATES IN ANGSTROMS, RADIANS,
*                     AND RADIANS
*
***********************************************************************
      do 30 I=2,NUMAT
         J=NA(I)
         K=NB(I)
         L=NC(I)
         if (I.LT.3) go to 30
         II=I
         call BANGLE(XYZ,II,J,K,GEO(2,I))
         GEO(2,I)=GEO(2,I)*DEGREE
         if (I.LT.4) go to 30
C
C   MAKE SURE DIHEDRAL IS MEANINGLFUL
C
         call BANGLE(XYZ,J,K,L,ANGL)
         TOL=0.2617994D0
         if (ANGL.GT.3.1415926D0-TOL.OR.ANGL.LT.TOL) then
C
C  ANGLE IS UNSATISFACTORY, LET'S SEARCH FOR ANOTHER ATOM FOR
C  DEFINING THE DIHEDRAL.
   10       SUM=100.D0
            do 20 I1=1,II-1
               R=(XYZ(1,I1)-XYZ(1,K))**2+
     1          (XYZ(2,I1)-XYZ(2,K))**2+
     2          (XYZ(3,I1)-XYZ(3,K))**2
               if (R.LT.SUM.AND.I1.NE.J.AND.I1.NE.K) then
                  call BANGLE(XYZ,J,K,I1,ANGL)
                  if (ANGL.LT.3.1415926D0-TOL.AND.ANGL.GT.TOL) then
                     SUM=R
                     L=I1
                     NC(II)=L
                  end if
               end if
   20       continue
            if (SUM.GT.99.D0.AND.TOL.GT.0.1D0) then
C
C ANYTHING WITHIN 5 DEGREES?
C
               TOL=0.087266D0
               go to 10
            end if
         end if
         call DIHED(XYZ,II,J,K,L,GEO(3,I))
         GEO(3,I)=GEO(3,I)*DEGREE
   30 GEO(1,I)= SQRT((XYZ(1,I)-XYZ(1,J))**2+
     1                   (XYZ(2,I)-XYZ(2,J))**2+
     2                   (XYZ(3,I)-XYZ(3,J))**2)
      GEO(1,1)=0.D0
      GEO(2,1)=0.D0
      GEO(3,1)=0.D0
      GEO(2,2)=0.D0
      GEO(3,2)=0.D0
      GEO(3,3)=0.D0
      return
      end
