      SUBROUTINE XYZINT(XYZ,NUMAT,NA,NB,NC,DEGREE,GEO)
      implicit double precision (A-H,O-Z)
      dimension XYZ(3,*), NA(*), NB(*), NC(*), GEO(3,*)
      include 'MOP6FILE'
***********************************************************************
*
* XYZINT WORKS OUT THE INTERNAL COORDINATES OF A MOLECULE.
*        THE "RULES" FOR THE CONNECTIVITY ARE AS FOLLOWS:
*        ATOM I IS DEFINED AS BEING AT A DISTANCE FROM THE NEAREST
*        ATOM J, ATOM J ALREADY HAVING BEEN DEFINED.
*        ATOM I MAKES AN ANGLE WITH ATOM J AND THE ATOM K, WHICH HAS
*        ALREADY BEEN DEFINED, AND IS THE NEAREST ATOM TO J
*        ATOM I MAKES A DIHEDRAL ANGLE WITH ATOMS J, K, AND L. L HAVING
*        BEEN DEFINED AND IS THE NEAREST ATOM TO K, AND J, K AND L
*        HAVE A CONTAINED ANGLE IN THE RANGE 15 TO 165 DEGREES,
*        IF POSSIBLE.
*
*        IF(NA(2).EQ.1 THEN THE ORIGINAL CONNECTIVITY IS USED.
*
*        NOTE THAT GEO AND XYZ MUST NOT BE THE SAME IN THE call.
*
*   ON INPUT XYZ    = CARTESIAN ARRAY OF NUMAT ATOMS
*            DEGREE = 1 IF ANGLES ARE TO BE IN RADIANS
*            DEGREE = 57.29578 IF ANGLES ARE TO BE IN DEGREES
*
***********************************************************************
      common /GEOOK/ IGEOOK
      common /NUMCAL/ NUMCAL
      data ICALCN/0/
      IGEOOK=99
      NAI1=0
      NAI2=0
      if (.NOT.(ICALCN.NE.NUMCAL).AND.NA(2).EQ.-1) then
         NA(2)=1
         do 10 I=2,NUMAT
            J=NA(I)
            if (I.GT.3) call DIHED(XYZ,I,J,NB(I),NC(I),GEO(3,I))
            if (I.GT.2) call BANGLE(XYZ,I,J,NB(I),GEO(2,I))
            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)
   10    continue
      else
         if (NA(2).EQ.-1) ICALCN=NUMCAL
         do 30 I=1,NUMAT
            NA(I)=2
            NB(I)=3
            NC(I)=4
            IM1=I-1
            if (IM1.EQ.0) go to 30
            SUM=1.D30
            do 20 J=1,IM1
               R=(XYZ(1,I)-XYZ(1,J))**2+
     1          (XYZ(2,I)-XYZ(2,J))**2+
     2          (XYZ(3,I)-XYZ(3,J))**2
               if (R.LT.SUM.AND.NA(J).NE.J.AND.NB(J).NE.J) then
                  SUM=R
                  K=J
               end if
   20       continue
C
C   ATOM I IS NEAREST TO ATOM K
C
            NA(I)=K
            if (I.GT.2) NB(I)=NA(K)
            if (I.GT.3) NC(I)=NB(K)
C
C   FIND ANY ATOM TO RELATE TO NA(I)
C
   30    continue
      end if
      NA(1)=0
      NB(1)=0
      NC(1)=0
      NB(2)=0
      NC(2)=0
      NC(3)=0
      call XYZGEO(XYZ,NUMAT,NA,NB,NC,DEGREE,GEO)
      return
      end
