subroutine match_contact(ishif1,ishif2,nc_match,nc_match1_max,
     &   ncont_ref,icont_ref,ncont,icont,jfrag,n_shif1,n_shif2,
     &   nc_frac,nc_req_set,istr,llocal,lprn)
      implicit real*8 (a-h,o-z)
      include 'DIMENSIONS'
      include 'COMMON.IOUNITS'
      integer ncont_ref,icont_ref(2,maxcont),ncont,icont(2,maxcont),
     &   ishift,ishif2,nc_match
      double precision nc_frac
      logical llocal,lprn
      nc_match_max=0
      do i=1,ncont_ref
        nc_match_max=nc_match_max+
     &   min0(icont_ref(2,i)-icont_ref(1,i)-1,3)
      enddo
      if (istr.eq.3) then
        nc_req=0
      else if (nc_req_set.eq.0) then
        nc_req=nc_match_max*nc_frac
      else
        nc_req = dmin1(nc_match_max*nc_frac+0.5d0,
     &    dfloat(nc_req_set)+1.0d-7)
      endif
c      write (iout,*) "match_contact: nc_req:",nc_req
c      write (iout,*) "nc_match_max",nc_match_max
c      write (iout,*) "jfrag",jfrag," n_shif1",n_shif1,
c     &   " n_shif2",n_shif2
C Match current contact map against reference contact map; exit, if at least
C half of the contacts match
      call ncont_match(nc_match,nc_match1,0,0,ncont_ref,icont_ref,
     &    ncont,icont,jfrag,llocal,lprn)
      nc_match1_max=nc_match1
      if (lprn .and. nc_match.gt.0) write (iout,*) 
     &  "Shift:",0,0," nc_match1",nc_match1,
     &  " nc_match=",nc_match," req'd",nc_req
      if (nc_req.gt.0 .and. nc_match.ge.nc_req .or.
     &    nc_req.eq.0 .and. nc_match.eq.1) then
         ishif1=0
         ishif2=0
         return
      endif
C If sufficient matches are not found, try to shift contact maps up to three
C positions.
      if (n_shif1.gt.0) then
      do is=1,n_shif1
C The following four tries help to find shifted beta-sheet patterns
C Shift "left" strand backward
        call ncont_match(nc_match,nc_match1,-is,0,ncont_ref,
     &    icont_ref,ncont,icont,jfrag,llocal,lprn)
        if (nc_match1.gt.nc_match1_max) nc_match1_max=nc_match1
        if (lprn .and. nc_match.gt.0) write (iout,*) 
     &    "Shift:",-is,0," nc_match1",nc_match1,
     &    " nc_match=",nc_match," req'd",nc_req
        if (nc_req.gt.0 .and. nc_match.ge.nc_req .or.
     &     nc_req.eq.0 .and. nc_match.eq.1) then
          ishif1=-is
          ishif2=0
          return
        endif
C Shift "left" strand forward
        call ncont_match(nc_match,nc_match1,is,0,ncont_ref,
     &      icont_ref,ncont,icont,jfrag,llocal,lprn)
        if (nc_match1.gt.nc_match1_max) nc_match1_max=nc_match1
        if (lprn .and. nc_match.gt.0) write (iout,*) 
     &   "Shift:",is,0," nc_match1",nc_match1,
     &   " nc_match=",nc_match," req'd",nc_req
        if (nc_req.gt.0 .and. nc_match.ge.nc_req .or.
     &     nc_req.eq.0 .and. nc_match.eq.1) then
          ishif1=is
          ishif2=0
          return
        endif
      enddo
      if (nc_req.eq.0) return
C Shift "right" strand backward
      do is=1,n_shif1
        call ncont_match(nc_match,nc_match1,0,-is,ncont_ref,
     &     icont_ref,ncont,icont,jfrag,llocal,lprn)
        if (nc_match1.gt.nc_match1_max) nc_match1_max=nc_match1
        if (lprn .and. nc_match.gt.0) write (iout,*) 
     &    "Shift:",0,-is," nc_match1",nc_match1,
     &    " nc_match=",nc_match," req'd",nc_req
        if (nc_match.ge.nc_req) then
          ishif1=0
          ishif2=-is
          return
        endif
C Shift "right" strand upward
        call ncont_match(nc_match,nc_match1,0,is,ncont_ref,
     &    icont_ref,ncont,icont,jfrag,llocal,lprn)
        if (nc_match1.gt.nc_match1_max) nc_match1_max=nc_match1
        if (lprn .and. nc_match.gt.0) write (iout,*) 
     &    "Shift:",0,is," nc_match1",nc_match1,
     &    " nc_match=",nc_match," req'd",nc_req
        if (nc_match.ge.nc_req) then
          ishif1=0
          ishif2=is
          return
        endif
      enddo ! is
C Now try to shift both residues in contacts.
      do is=1,n_shif1
        do js=1,is
          if (js.ne.is) then
            call ncont_match(nc_match,nc_match1,-is,-js,ncont_ref,
     &        icont_ref,ncont,icont,jfrag,llocal,lprn)
            if (nc_match1.gt.nc_match1_max) nc_match1_max=nc_match1
            if (lprn .and. nc_match.gt.0) write (iout,*) 
     &         "Shift:",-is,-js," nc_match1",nc_match1,
     &         " nc_match=",nc_match," req'd",nc_req
            if (nc_match.ge.nc_req) then
              ishif1=-is
              ishif2=-js
              return
            endif
            call ncont_match(nc_match,nc_match1,is,js,ncont_ref,
     &        icont_ref,ncont,icont,jfrag,llocal,lprn)
            if (nc_match1.gt.nc_match1_max) nc_match1_max=nc_match1
            if (lprn .and. nc_match.gt.0) write (iout,*) 
     &        "Shift:",is,js," nc_match1",nc_match1,
     &        " nc_match=",nc_match," req'd",nc_req
            if (nc_match.ge.nc_req) then
              ishif1=is
              ishif2=js
              return
            endif
c
            call ncont_match(nc_match,nc_match1,-js,-is,ncont_ref,
     &        icont_ref,ncont,icont,jfrag,llocal,lprn)
            if (nc_match1.gt.nc_match1_max) nc_match1_max=nc_match1
            if (lprn .and. nc_match.gt.0) write (iout,*) 
     &        "Shift:",-js,-is," nc_match1",nc_match1,
     &        " nc_match=",nc_match," req'd",nc_req
            if (nc_match.ge.nc_req) then
              ishif1=-js
              ishif2=-is
              return
            endif
c
            call ncont_match(nc_match,nc_match1,js,is,ncont_ref,
     &        icont_ref,ncont,icont,jfrag,llocal,lprn)
            if (nc_match1.gt.nc_match1_max) nc_match1_max=nc_match1
            if (lprn .and. nc_match.gt.0) write (iout,*) 
     &        "Shift:",js,is," nc_match1",nc_match1,
     &        " nc_match=",nc_match," req'd",nc_req
            if (nc_match.ge.nc_req) then
              ishif1=js
              ishif2=is
              return
            endif
          endif
c
          if (is+js.le.n_shif1) then
          call ncont_match(nc_match,nc_match1,-is,js,ncont_ref,
     &      icont_ref,ncont,icont,jfrag,llocal,lprn)
          if (nc_match1.gt.nc_match1_max) nc_match1_max=nc_match1
          if (lprn .and. nc_match.gt.0) write (iout,*) 
     &     "Shift:",-is,js," nc_match1",nc_match1,
     &     " nc_match=",nc_match," req'd",nc_req
          if (nc_match.ge.nc_req) then
            ishif1=-is
            ishif2=js
            return
          endif
c
          call ncont_match(nc_match,nc_match1,js,-is,ncont_ref,
     &      icont_ref,ncont,icont,jfrag,llocal,lprn)
          if (nc_match1.gt.nc_match1_max) nc_match1_max=nc_match1
          if (lprn .and. nc_match.gt.0) write (iout,*) 
     &     "Shift:",js,-is," nc_match1",nc_match1,
     &     " nc_match=",nc_match," req'd",nc_req
          if (nc_match.ge.nc_req) then
            ishif1=js
            ishif2=-is
            return
          endif
          endif
c
        enddo !js
      enddo !is
      endif

      if (n_shif2.gt.0) then
      do is=1,n_shif2
        call ncont_match(nc_match,nc_match1,-is,-is,ncont_ref,
     &    icont_ref,ncont,icont,jfrag,llocal,lprn)
        if (nc_match1.gt.nc_match1_max) nc_match1_max=nc_match1
        if (lprn .and. nc_match.gt.0) write (iout,*) 
     &     "Shift:",-is,-is," nc_match1",nc_match1,
     &     " nc_match=",nc_match," req'd",nc_req
        if (nc_match.ge.nc_req) then
          ishif1=-is
          ishif2=-is
          return
        endif
        call ncont_match(nc_match,nc_match1,is,is,ncont_ref,
     &    icont_ref,ncont,icont,jfrag,llocal,lprn)
        if (nc_match1.gt.nc_match1_max) nc_match1_max=nc_match1
        if (lprn .and. nc_match.gt.0) write (iout,*) 
     &    "Shift:",is,is," nc_match1",nc_match1,
     &    " nc_match=",nc_match," req'd",nc_req
        if (nc_match.ge.nc_req) then
          ishif1=is
          ishif2=is
          return
        endif
      enddo
      endif
C If this point is reached, the contact maps are different. 
      nc_match=0
      ishif1=0
      ishif2=0
      return
      end
c-------------------------------------------------------------------------
      subroutine ncont_match(nc_match,nc_match1,ishif1,ishif2,
     &   ncont_ref,icont_ref,ncont,icont,jfrag,llocal,lprn)
      implicit real*8 (a-h,o-z)
      include 'DIMENSIONS'
      include 'DIMENSIONS.ZSCOPT'
      include 'DIMENSIONS.COMPAR'
      include 'COMMON.IOUNITS'
      include 'COMMON.INTERACT'
      include 'COMMON.GEO'
      include 'COMMON.COMPAR'
      logical llocal,lprn
      integer ncont_ref,icont_ref(2,maxcont),ncont,icont(2,maxcont),
     &   icont_match(2,maxcont),ishift,ishif2,nang_pair,
     &   iang_pair(2,maxres)
C Compare the contact map against the reference contact map; they're stored
C in ICONT and ICONT_REF, respectively. The current contact map can be shifted.
      if (lprn) write (iout,'(80(1h*))')
      nc_match=0
      nc_match1=0
c Check the local structure by comparing dihedral angles.
c      write (iout,*) "ncont_match: ncont_ref",ncont_ref," llocal",llocal
      if (llocal .and. ncont_ref.eq.0) then
c If there are no contacts just compare the dihedral angles and exit.
        call angnorm(jfrag,ishif1,ishif2,ang_cut1(jfrag),diffang,fract,
     &    lprn)
        if (lprn) write (iout,*) "diffang:",diffang*rad2deg,
     &   " ang_cut:",ang_cut(jfrag)*rad2deg," fract",fract
        if (diffang.le.ang_cut(jfrag) .and. fract.ge.frac_min(jfrag)) 
     &  then
          nc_match=1
        else
          nc_match=0
        endif
        return
      endif
      nang_pair=0
      do i=1,ncont
        ic1=icont(1,i)+ishif1
        ic2=icont(2,i)+ishif2
c        write (iout,*) "i",i," ic1",ic1," ic2",ic2
        if (ic1.lt.nnt .or. ic2.gt.nct) goto 10
        do j=1,ncont_ref
          if (ic1.eq.icont_ref(1,j).and.ic2.eq.icont_ref(2,j)) then
            nc_match=nc_match+min0(icont_ref(2,j)-icont_ref(1,j)-1,3)
            nc_match1=nc_match1+1
            icont_match(1,nc_match1)=ic1
            icont_match(2,nc_match1)=ic2
c            call add_angpair(icont(1,i),icont_ref(1,j),
c     &         nang_pair,iang_pair)
c            call add_angpair(icont(2,i),icont_ref(2,j),
c     &         nang_pair,iang_pair) 
            if (lprn) write (iout,*) "Contacts:",icont(1,i),icont(2,i),
     &       " match",icont_ref(1,j),icont_ref(2,j),
     &       " shifts",ishif1,ishif2
            goto 10
          endif
        enddo 
   10   continue
      enddo
      if (lprn) then
        write (iout,*) "nc_match",nc_match," nc_match1",nc_match1
        write (iout,*) "icont_match"
        do i=1,nc_match1
          write (iout,*) icont_match(1,i),icont_match(2,i)
        enddo
      endif
      if (llocal .and. nc_match.gt.0) then
        call angnorm2(jfrag,ishif1,ishif2,nc_match1,icont_match,lprn,
     &    ang_cut1(jfrag),diffang,fract)
        if (lprn) write (iout,*) "diffang:",diffang*rad2deg,
     &   " ang_cut:",ang_cut(jfrag)*rad2deg,
     &   " ang_cut1",ang_cut1(jfrag)*rad2deg
        if (diffang.gt.ang_cut(jfrag) 
     &    .or. fract.lt.frac_min(jfrag)) nc_match=0
      endif
c      if (nc_match.gt.0) then
c        diffang = angnorm1(nang_pair,iang_pair,lprn)
c        if (diffang.gt.ang_cut(jfrag)) nc_match=0
c      endif
      if (lprn) write (iout,*) "ishif1",ishif1," ishif2",ishif2,
     &   " diffang",rad2deg*diffang," nc_match",nc_match
      return
      end
c------------------------------------------------------------------------------
      subroutine match_secondary(jfrag,isecstr,nsec_match,lprn)
c This subroutine compares the secondary structure (isecstr) of fragment jfrag 
c conformation considered to that of the reference conformation.
c Returns the number of equivalent residues (nsec_match).
      implicit real*8 (a-h,o-z)
      include 'DIMENSIONS'
      include 'DIMENSIONS.ZSCOPT'
      include 'DIMENSIONS.COMPAR'
      include 'COMMON.IOUNITS'
      include 'COMMON.CHAIN'
      include 'COMMON.PEPTCONT'
      include 'COMMON.COMPAR'
      logical lprn
      integer isecstr(maxres)
      npart = npiece(jfrag,1)
      nsec_match=0
      if (lprn) then
        write (iout,*) "match_secondary jfrag",jfrag," ifrag",
     &        (ifrag(1,i,jfrag),ifrag(2,i,jfrag),i=1,npart)
        write (iout,'(80i1)') (isec_ref(j),j=1,nres)
        write (iout,'(80i1)') (isecstr(j),j=1,nres)
      endif
      do i=1,npart
        do j=ifrag(1,i,jfrag),ifrag(2,i,jfrag)
c The residue has equivalent conformational state to that of the reference
c structure, if:
c  a) the conformational states are equal or
c  b) the reference state is a coil and that of the conformation considered 
c     is a strand or
c  c) the conformational state of the conformation considered is a strand
c     and that of the reference conformation is a coil.
c 10/28/02 - case (b) deleted.
          if (isecstr(j).eq.isec_ref(j) .or. 
c     &        isecstr(j).eq.0 .and. isec_ref(j).eq.1 .or.
     &        isec_ref(j).eq.0 .and. isecstr(j).eq.1) 
     &      nsec_match=nsec_match+1 
        enddo
      enddo
      return
      end