      SUBROUTINE GETGEO(IREAD,LABELS,GEO,LOPT,NA,NB,NC,AMS,NATOMS,INT)
      IMPLICIT DOUBLE PRECISION (A-H,O-Z)
      INCLUDE 'mopac.inc'
      DIMENSION GEO(3,*),NA(*),NB(*),NC(*),AMS(*), LOPT(3,*)
     1,LABELS(*)
      LOGICAL INT
************************************************************************
*
*   GETGEO READS IN THE GEOMETRY. THE ELEMENT IS SPECIFIED BY IT'S
*          CHEMICAL SYMBOL, OR, OPTIONALLY, BY IT'S ATOMIC NUMBER.
*
*  ON INPUT   IREAD  = CHANNEL NUMBER FOR READ, NORMALLY 5
*             AMS    = DEFAULT ATOMIC MASSES.
*
* ON OUTPUT LABELS = ATOMIC NUMBERS OF ALL ATOMS, INCLUDING DUMMIES.
*           GEO    = INTERNAL COORDINATES, IN ANGSTROMS, AND DEGREES.
*           LOPT   = INTEGER ARRAY, A '1' MEANS OPTIMIZE THIS PARAMETER,
*                    '0' MEANS DO NOT OPTIMIZE, AND A '-1' LABELS THE
*                    REACTION COORDINATE.
*           NA     = INTEGER ARRAY OF ATOMS (SEE DATA INPUT)
*           NB     = INTEGER ARRAY OF ATOMS (SEE DATA INPUT)
*           NC     = INTEGER ARRAY OF ATOMS (SEE DATA INPUT)
*           ATMASS = ATOMIC MASSES OF ATOMS.
************************************************************************
      COMMON /PATH  / DUMMY,REACT(3,66), DUMM1,DUMM2
      COMMON /ATMASS/ ATMASS(NUMATM)
      COMMON /KEYWRD/ KEYWRD
      DIMENSION ISTART(40), XYZ(3,NUMATM), VALUE(4)
      LOGICAL LEADSP, IRCDRC
      CHARACTER*80 KEYWRD
      CHARACTER ELEMNT(107)*2, LINE*80, SPACE*1, NINE*1,ZERO*1,
     1TAB*1, COMMA*1, STRING*80, ELE*2, TURN*1
      DATA (ELEMNT(I),I=1,107)/'H','HE',
     1 'LI','BE','B','C','N','O','F','NE',
     2 'NA','MG','AL','SI','P','S','CL','AR',
     3 'K','CA','SC','TI','V','CR','MN','FE','CO','NI','CU',
     4 'ZN','GA','GE','AS','SE','BR','KR',
     5 'RB','SR','Y','ZR','NB','MO','TC','RU','RH','PD','AG',
     6 'CD','IN','SN','SB','TE','I','XE',
     7 'CS','BA','LA','CE','PR','ND','PM','SM','EU','GD','TB','DY',
     8 'HO','ER','TM','YB','LU','HF','TA','W','RE','OS','IR','PT',
     9 'AU','HG','TL','PB','BI','PO','AT','RN',
     1 'FR','RA','AC','TH','PA','U','NP','PU','AM','CM','BK','CF','XX',
     2 'FM','MD','CB','++','+','--','-','TV'/
      DATA COMMA,SPACE,NINE,ZERO/',',' ','9','0'/
      TAB=CHAR(9)
      IRCDRC=(INDEX(KEYWRD,'IRC')+INDEX(KEYWRD,'DRC') .NE.0)
      ILOWA = ICHAR('a')
      ILOWZ = ICHAR('z')
      ICAPA = ICHAR('A')
      NATOMS=0
      NUMAT=0
   10 READ(IREAD,'(A)',END=90,ERR=180)LINE
      IF(LINE.EQ.' ') GO TO 90
      NATOMS=NATOMS+1
      IF(NATOMS.GT.NUMATM)THEN
         WRITE(6,'(//10X,''****  MAX. NUMBER OF ATOMS ALLOWED:'',I4)')
     1NUMATM
         CALL WRAPUP()
      ENDIF
*   CLEAN THE INPUT DATA
************************************************************************
      DO 20 I=1,80
         ILINE=ICHAR(LINE(I:I))
         IF(ILINE.GE.ILOWA.AND.ILINE.LE.ILOWZ) THEN
            LINE(I:I)=CHAR(ILINE+ICAPA-ILOWA)
         ENDIF
   20 CONTINUE
************************************************************************
      ICOMMA=ICHAR(COMMA)
      ITAB=ICHAR(TAB)
      DO 30 I=1,80
      KHAR=ICHAR(LINE(I:I))
      IF(KHAR.EQ.ICOMMA.OR.KHAR.EQ.ITAB)LINE(I:I)=SPACE
  30  CONTINUE
*
*   INITIALIZE ISTART TO INTERPRET BLANKS AS ZERO'S
      DO 40 I=1,10
   40 ISTART(I)=80
* FIND INITIAL DIGIT OF ALL NUMBERS, CHECK FOR LEADING SPACES FOLLOWED
*     BY A CHARACTER AND STORE IN ISTART
      LEADSP=.TRUE.
      NVALUE=0
      DO 50 I=1,80
         IF (LEADSP.AND.LINE(I:I).NE.SPACE) THEN
            NVALUE=NVALUE+1
            ISTART(NVALUE)=I
         END IF
         LEADSP=(LINE(I:I).EQ.SPACE)
   50 CONTINUE
*
* ESTABLISH THE ELEMENT'S NAME AND ISOTOPE, CHECK FOR ERRORS OR E.O.DATA
*
      WEIGHT=0.D0
      STRING=LINE(ISTART(1):ISTART(2)-1)
      IF( STRING(1:1) .GE. ZERO .AND. STRING(1:1) .LE. NINE) THEN
*  ATOMIC NUMBER USED: NO ISOTOPE ALLOWED
         LABEL=READA(STRING,1)
         IF (LABEL.EQ.0) GO TO 80
         IF (LABEL.LT.0.OR.LABEL.GT.107) THEN
            WRITE(6,'(''  ILLEGAL ATOMIC NUMBER'')')
            GO TO 190
         END IF
         GO TO 70
      END IF
*  ATOMIC SYMBOL USED
      REAL=ABS(READA(STRING,1))
      IF (REAL.LT.1.D-15) THEN
*   NO ISOTOPE
         ELE=STRING(1:2)
      ELSE
         WEIGHT=REAL
         IF( STRING(2:2) .GE. ZERO .AND. STRING(2:2) .LE. NINE) THEN
            ELE=STRING(1:1)
         ELSE
            ELE=STRING(1:2)
         END IF
      END IF
*   CHECK FOR ERROR IN ATOMIC SYMBOL
      IF(ELE(1:1).EQ.'-'.AND.ELE(2:2).NE.'-')ELE(2:2)=' '
      DO 60 I=1,107
         IF(ELE.EQ.ELEMNT(I)) THEN
            LABEL=I
            GO TO 70
         END IF
   60 CONTINUE
      WRITE(6,'(''  UNRECOGNIZED ELEMENT NAME: ('',A,'')'')')ELE
      GOTO 190
*
* ALL O.K.
*
   70 IF (LABEL.NE.99) NUMAT=NUMAT+1
      IF(WEIGHT.NE.0.D0)THEN
         WRITE(6,'('' FOR ATOM'',I4,''  ISOTOPIC MASS:''
     1    ,F15.5)')NATOMS, WEIGHT
         ATMASS(NUMAT)=WEIGHT
      ELSE
         IF(LABEL .NE. 99)  ATMASS(NUMAT)=AMS(LABEL)
      ENDIF
      LABELS(NATOMS)   =LABEL
      GEO(1,NATOMS)    =READA(LINE,ISTART(2))
      GEO(2,NATOMS)    =READA(LINE,ISTART(4))
      GEO(3,NATOMS)    =READA(LINE,ISTART(6))
      IF(IRCDRC)THEN
         TURN=LINE(ISTART(3):ISTART(3))
         IF(TURN.EQ.'T')THEN
            LOPT(1,NATOMS)=1
            IF(NATOMS.EQ.1)WRITE(6,'(A)')' IN DRC MONITOR POTENTIAL ENER
     1GY'//' TURNING POINTS'
         ELSE
            LOPT(1,NATOMS)=0
         ENDIF
         TURN=LINE(ISTART(5):ISTART(5))
         IF(TURN.EQ.'T')THEN
            LOPT(2,NATOMS)=1
         ELSE
            LOPT(2,NATOMS)=0
         ENDIF
         TURN=LINE(ISTART(7):ISTART(7))
         IF(TURN.EQ.'T')THEN
            LOPT(3,NATOMS)=1
         ELSE
            LOPT(3,NATOMS)=0
         ENDIF
      ELSE
         LOPT(1,NATOMS)   =READA(LINE,ISTART(3))
         LOPT(2,NATOMS)   =READA(LINE,ISTART(5))
         LOPT(3,NATOMS)   =READA(LINE,ISTART(7))
      ENDIF
      NA(NATOMS)       =READA(LINE,ISTART(8))
      NB(NATOMS)       =READA(LINE,ISTART(9))
      NC(NATOMS)       =READA(LINE,ISTART(10))
      GOTO 10
*
* ALL DATA READ IN, CLEAN UP AND RETURN
*
   80 NATOMS=NATOMS-1
   90 NA(2)=1
      IF(NATOMS.GT.3)THEN
         INT=(NA(4).NE.0)
      ELSE
         IF(GEO(2,3).LT.10.AND.NATOMS.EQ.3)
     1WRITE(6,'(//10X,'' WARNING: INTERNAL COORDINATES ARE ASSUMED -'',/
     210X,'' FOR THREE-ATOM SYSTEMS '',//)')
         INT=.TRUE.
      ENDIF
C
C     READ IN VELOCITY VECTOR, IF PRESENT
C
      IF(INDEX(KEYWRD,'VELO').GT.0)THEN
         IF(INT)THEN
            WRITE(6,'(A)')' COORDINATES MUST BE CARTESIAN WHEN VELOCITY
     1VECTOR IS USED.'
            CALL WRAPUP()
         ENDIF
C#      WRITE(6,'(/10X,A)')'INITIAL VELOCITY VECTOR FOR DRC'
         DO 110 I=1,NATOMS
            READ(5,'(A)') LINE
            CALL NUCHAR(LINE,VALUE,NDMY)
            IF(NDMY.NE.3)THEN
               WRITE(6,'(/10X,A)')
     1'  THERE MUST BE EXACTLY THREE VELOCITY DATA PER LINE'
               CALL WRAPUP()
            ENDIF
            DO 100 J=1,3
  100       REACT(J,I+2)=VALUE(J)
C#      WRITE(6,'(2X,A2,2X,3F13.5)')ELEMNT(LABELS(I)),(VALUE(J),J=1,3)
  110    CONTINUE
         DO 120 I=1,3
            DO 120 J=1,2
  120    REACT(I,J)=GEO(I,J+1)-GEO(I,1)
C
C  NOW TO ROTATE VELOCITY VECTOR TO SUIT INTERNAL COORDINATE DEFINITION
C
C
C   ROTATE AROUND THE 1-2 X-AXIS TO AS TO ELIMINATE REACT(3,2)
C   (PUT ATOM 2 IN X-Y PLANE)
         SA=REACT(3,1)/SQRT(REACT(2,1)**2+REACT(3,1)**2+1.D-20)
         CA=SIGN(SQRT(1.D0-SA**2),REACT(2,1))
C#      LABELS(NATOMS+1)=1
C#      LABELS(NATOMS+2)=1
C#      WRITE(6,*)' FIRST ROTATION, ABOUT 1-2 X-AXIS'
         DO 130 I=1,NATOMS+2
            TEMP1= REACT(2,I)*CA+REACT(3,I)*SA
            TEMP2=-REACT(2,I)*SA+REACT(3,I)*CA
            REACT(2,I)=TEMP1
            REACT(3,I)=TEMP2
C#      WRITE(6,'(2X,A2,2X,3F13.5)')ELEMNT(LABELS(I)),(REACT(J,I),J=1,3)
  130    CONTINUE
C   ROTATE AROUND THE 1-2 Z-AXIS TO AS TO ELIMINATE REACT(2,2)
C   (PUT ATOM 2 ON X AXIS)
         CA=REACT(1,1)/SQRT(REACT(2,1)**2+REACT(1,1)**2+1.D-20)
         SA=SIGN(SQRT(1.D0-CA**2),REACT(2,1))
C#      WRITE(6,*)' SECOND ROTATION, ABOUT 1-2 Z-AXIS'
         DO 140 I=1,NATOMS+2
            TEMP1= REACT(1,I)*CA+REACT(2,I)*SA
            TEMP2=-REACT(1,I)*SA+REACT(2,I)*CA
            REACT(1,I)=TEMP1
            REACT(2,I)=TEMP2
C#      WRITE(6,'(2X,A2,2X,3F13.5)')ELEMNT(LABELS(I)),(REACT(J,I),J=1,3)
  140    CONTINUE
C   ROTATE AROUND THE 2-3 X-AXIS TO AS TO ELIMINATE REACT(3,3)
C   (PUT ATOM 3 ON X-Y PLANE)
         SA=REACT(3,2)/SQRT(REACT(2,2)**2+REACT(3,2)**2+1.D-20)
         CA=SIGN(SQRT(1.D0-SA**2),REACT(2,2))
C#      WRITE(6,*)' THIRD ROTATION, ABOUT 2-3 X-AXIS'
         DO 150 I=1,NATOMS+2
            TEMP1= REACT(2,I)*CA+REACT(3,I)*SA
            TEMP2=-REACT(2,I)*SA+REACT(3,I)*CA
            REACT(2,I)=TEMP1
            REACT(3,I)=TEMP2
C#      WRITE(6,'(2X,A2,2X,3F13.5)')ELEMNT(LABELS(I)),(REACT(J,I),J=1,3)
  150    CONTINUE
C
C  STRIP OFF FIRST TWO COORDINATES; THESE WERE THE COORDINATE AXIS
C  DEFINITIONS
C
         DO 160 I=1,NATOMS
            DO 160 J=1,3
  160    REACT(J,I)=REACT(J,I+2)
      ENDIF
      IF(  .NOT. INT ) THEN
         DO 170 I=1,NATOMS
            DO 170 J=1,3
  170    XYZ(J,I)=GEO(J,I)
         DEGREE=90.D0/ASIN(1.D0)
         CALL XYZINT(XYZ,NATOMS,NA,NB,NC,DEGREE,GEO)
         IF(ABS(GEO(2,3)-180.D0).LT.1.D-4.OR.ABS(GEO(2,3)).LT.1.D-4)
     1THEN
            WRITE(6,'(A)')' DUE TO PROGRAM BUG, THE FIRST THREE ATOMS MU
     1ST NOT LIE IN A STRAIGHT LINE.'
            CALL WRAPUP()
         ENDIF
      ELSEIF (.NOT.IRCDRC) THEN
         IF(LOPT(1,1)+LOPT(2,1)+LOPT(3,1)+LOPT(2,2)+LOPT(3,2)+
     1        LOPT(3,3) .GT. 0)THEN
            LOPT(1,1)=0
            LOPT(2,1)=0
            LOPT(3,1)=0
            LOPT(2,2)=0
            LOPT(3,2)=0
            LOPT(3,3)=0
            WRITE(6,
     +'(//10X,'' AN UNOPTIMIZABLE GEOMETRIC PARAMETER HAS'',
     1/10X,'' BEEN MARKED FOR OPTIMIZATION. THIS IS A NON-FATAL '',
     2''ERROR'')')
         ENDIF
      ENDIF
      IF(NA(3).EQ.0) THEN
         NB(3)=1
         NA(3)=2
      ENDIF
      RETURN
* ERROR CONDITIONS
  180 IF(IREAD.EQ.5) THEN
         WRITE(6,'( '' ERROR DURING READ AT ATOM NUMBER '', I3 )')NATOMS
      ELSE
         NATOMS=0
         RETURN
      ENDIF
  190 J=NATOMS-1
      WRITE(6,'('' DATA CURRENTLY READ IN ARE '')')
      DO 200 K=1,J
  200 WRITE(6,210)LABELS(K),(GEO(JJ,K),LOPT(JJ,K),JJ=1,3),
     1NA(K),NB(K),NC(K)
  210 FORMAT(I4,2X,3(F10.5,2X,I2,2X),3(I2,1X))
      CALL WRAPUP()
      END
