C
C $Id: getpep.F,v 1.2 1998/07/16 16:39:56 jjv5 Exp arjan $
C
C------------------------------------------------------------------------
      SUBROUTINE GETPEP(IERROR)
C
C     LOCATES AND STORES ATOMS INVOLVED IN EACH PEPTIDE BOND.  THE
C     ARRANGEMENT H-N-C-O IS LOCATED ON THE BASIS OF INTERATOMIC
C     DISTANCES.  THE ATOMS H-N-C-O ARE STORED IN IABC, AS WELL AS
C     X-N-C-O, WHERE X IS THE OTHER ATOM BONDED TO THE NITROGEN:
C
C                            X
C                           /
C                       H--N
C                           \
C                            C=O
C                           /
C
C     THE TOTAL NUMBER OF QUADRUPLES STORED IN IABC IS N2PEP.
C     IABC AND N2PEP ARE RETURNED VIA COMMON BLOCK /PEPTID/
C
C     ERROR FLAG IS SET TO 1 IF STORAGE LIMITS ARE EXCEEDED.
C
      IMPLICIT DOUBLE PRECISION (A-H,O-Z)
#include "divcon.dim"
#include "divcon.h"
C
C     LOCAL:
C
      DIMENSION IHNCO(MAXATM),ISTART(5),NAT(4),IAT(4)
C
      IERROR = 0
      N2PEP = 0
C
C     H,N,C,O ATOMIC NUMBERS:
C
      IAT(1) = 1
      IAT(2) = 7
      IAT(3) = 6
      IAT(4) = 8
C
C     MAKE FOUR PASSES THROUGH THE ATOMS AND STORE H, N, C, AND, O
C     ATOMS IN SEPARATE BLOCKS OF IHNCO.  ISTART WILL MARK THE BEGINNING
C     OF EACH BLOCK.
C
      DO 10 I=1,4
        ISTART(I) = 1
 10   CONTINUE
      IPOINT = 1
      DO 100 IBLOCK=1,4
        DO 50 IATM=1,NATOMS
          IF(IATNUM(IATM).EQ.IAT(IBLOCK))THEN
            IHNCO(IPOINT) = IATM
            IPOINT = IPOINT + 1
          ENDIF
 50     CONTINUE
        ISTART(IBLOCK+1) = IPOINT
 100  CONTINUE
      DO 120 I=1,4
        NAT(I) = ISTART(I+1) - ISTART(I)
 120  CONTINUE
C
C     DON'T BOTHER LOOKING IF WE'RE MISSING ANY OF THE FOUR ATOM TYPES.
C
      NMIN = MIN(NAT(1),NAT(2),NAT(3),NAT(4))
      IF(NMIN.EQ.0) RETURN
C
C     NOW INTERSECT THE FOUR BLOCKS AND CHECK INTERATOMIC DISTANCES TO
C     SEE IF WE HAVE A PEPTIDE FUNCTION.
C
      DO 300 I1=1,ISTART(2)-1
        IH = IHNCO(I1)
        XH = XYZ(1,IH)
        YH = XYZ(2,IH)
        ZH = XYZ(3,IH)
        DO 280 I2=ISTART(2),ISTART(3)-1
          IN = IHNCO(I2)
          XN = XYZ(1,IN)
          YN = XYZ(2,IN)
          ZN = XYZ(3,IN)
          DHN = DSQRT((XH-XN)**2 + (YH-YN)**2 + (ZH-ZN)**2)
          IF(DHN.GT.1.3D0) GO TO 280
          DO 260 I3=ISTART(3),ISTART(4)-1
            IC = IHNCO(I3)
            XC = XYZ(1,IC)
            YC = XYZ(2,IC)
            ZC = XYZ(3,IC)
            DNC = DSQRT((XN-XC)**2 + (YN-YC)**2 + (ZN-ZC)**2)
            IF(DNC.GT.1.6D0) GO TO 260
            DO 240 I4=ISTART(4),ISTART(5)-1
              IO = IHNCO(I4)
              XO = XYZ(1,IO)
              YO = XYZ(2,IO)
              ZO = XYZ(3,IO)
              DCO = DSQRT((XC-XO)**2 + (YC-YO)**2 + (ZC-ZO)**2)
              IF(DCO.GT.1.3D0) GO TO 240
C
C             IF WE GET HERE THEN WE'VE FOUND THE H-N-C-O FUNCTIONALITY.
C             MAKE SURE THERE'S SOMETHING ELSE BONDED TO NITROGEN, THEN
C             STORE H-N-C-O AND X-N-C-O.
C
              DO 220 IX=1,NATOMS
                IF(IX.EQ.IH.OR.
     .             IX.EQ.IN.OR.
     .             IX.EQ.IC.OR.
     .             IX.EQ.IO.OR.
     .             IATNUM(IX).EQ.0) GO TO 220
                XX = XYZ(1,IX)
                YX = XYZ(2,IX)
                ZX = XYZ(3,IX)
                DXN = DSQRT((XX-XN)**2 + (YX-YN)**2 + (ZX-ZN)**2)
                IF(DXN.GT.1.7D0) GO TO 220
C
C               STORE H-N-C-O AND X-N-C-O.
C
                N2PEP = N2PEP + 1
                IF(N2PEP.GT.MAXPEP)THEN
                  IERROR = 1
C-RDC                   WRITE(IOUT,
C-RDC      .                 '(/" STORAGE LIMIT EXCEEDED FOR PEPTIDE ",
C-RDC      .                  " BONDS"/" -- INCREASE MAXPEP PARAMETER IN",
C-RDC      .                  " divcon.dim AND RECOMPILE")')
                  RETURN
                ENDIF
                IABC(1,N2PEP) = IH
                IABC(2,N2PEP) = IN
                IABC(3,N2PEP) = IC
                IABC(4,N2PEP) = IO
                N2PEP = N2PEP + 1
                IABC(1,N2PEP) = IX
                IABC(2,N2PEP) = IN
                IABC(3,N2PEP) = IC
                IABC(4,N2PEP) = IO
                GO TO 240
 220          CONTINUE
 240        CONTINUE
 260      CONTINUE
 280    CONTINUE
 300  CONTINUE
      RETURN
      END