source/wham/src-HCD-5D/contact.f

   1       subroutine contact(lprint,ncont,icont,ist,ien,ipermmin)
   2       implicit none
   3       include 'DIMENSIONS'
   4       include 'DIMENSIONS.ZSCOPT'
   5       include 'COMMON.CONTROL'
   6       include 'COMMON.IOUNITS'
   7       include 'COMMON.CHAIN'
   8       include 'COMMON.INTERACT'
   9       include 'COMMON.FFIELD'
  10       include 'COMMON.NAMES'
  11       include 'COMMON.CALC'
  12       include 'COMMON.CONTPAR'
  13       include 'COMMON.LOCAL'
  14       integer ist,ien,kkk,iti,itj,itypi,itypj,i1,i2,it1,it2
  15       real*8 csc,dist
  16       real*8 cscore(maxcont),omt1(maxcont),omt2(maxcont),omt12(maxcont),
  17      & ddsc(maxcont),ddla(maxcont),ddlb(maxcont)
  18       integer ncont,icont(2,maxcont)
  19       real*8 u,v,a(3),b(3),dla,dlb
  20       logical lprint
  21       integer iperm,ipermmin,ii,jj
  22       ncont=0
  23       kkk=3
  24       if (lprint) then
  25       do i=1,nres
  26         write (iout,110) restyp(itype(i)),i,c(1,i),c(2,i),
  27      &    c(3,i),dc(1,nres+i),dc(2,nres+i),dc(3,nres+i),
  28      &    dc_norm(1,nres+i),dc_norm(2,nres+i),dc_norm(3,nres+i)
  29       enddo
  30       endif
  31   110 format (a,'(',i3,')',9f8.3)
  32       do i=ist,ien-kkk
  33         iti=iabs(itype(i))
  34         ii = iperm(i,ipermmin)
  35         if (iti.le.0 .or. iti.gt.ntyp) cycle
  36         do j=i+kkk,ien
  37           jj = iperm(j,ipermmin)
  38           itj=iabs(itype(j))
  39           if (itj.le.0 .or. itj.gt.ntyp) cycle
  40           itypi=iti
  41           itypj=itj
  42           xj = c(1,nres+jj)-c(1,nres+ii)
  43           yj = c(2,nres+jj)-c(2,nres+ii)
  44           zj = c(3,nres+jj)-c(3,nres+ii)
  45           dxi = dc_norm(1,nres+ii)
  46           dyi = dc_norm(2,nres+ii)
  47           dzi = dc_norm(3,nres+ii)
  48           dxj = dc_norm(1,nres+jj)
  49           dyj = dc_norm(2,nres+jj)
  50           dzj = dc_norm(3,nres+jj)
  51           do k=1,3
  52             a(k)=dc(k,nres+ii)
  53             b(k)=dc(k,nres+jj)
  54           enddo
  55 c          write (iout,*) (a(k),k=1,3),(b(k),k=1,3)
  56           if (icomparfunc.eq.1) then
  57             call contfunc(csc,iti,itj)
  58           else if (icomparfunc.eq.2) then
  59             call scdist(csc,iti,itj)
  60           else if (icomparfunc.eq.3 .or. icomparfunc.eq.5) then
  61             csc = dist(nres+i,nres+j)
  62           else if (icomparfunc.eq.4) then
  63             call odlodc(c(1,i),c(1,j),a,b,u,v,dla,dlb,csc)
  64           else
  65             write (*,*) "Error - Unknown sidechain contact function"
  66             write (iout,*) "Error - Unknown sidechain contact function"
  67           endif
  68           if (csc.lt.sc_cutoff(iti,itj)) then
  69 c            write(iout,*) "i",i," j",j," dla",dla,dsc(iti),
  70 c     &      " dlb",dlb,dsc(itj)," csc",csc,sc_cutoff(iti,itj),
  71 c     &      dxi,dyi,dzi,dxi**2+dyi**2+dzi**2,
  72 c     &      dxj,dyj,dzj,dxj**2+dyj**2+dzj**2,om1,om2,om12,
  73 c     &      xj,yj,zj
  74 c            write(iout,*)'egb',itypi,itypj,chi1,chi2,chip1,chip2,
  75 c     &       sig0ij,rij,rrij,om1,om2,om12,chiom1,chiom2,chiom12,
  76 c     &       chipom1,chipom2,chipom12,sig,eps2rt,rij_shift,e2,evdw,
  77 c     &       csc
  78             ncont=ncont+1
  79             cscore(ncont)=csc
  80             icont(1,ncont)=i
  81             icont(2,ncont)=j
  82             omt1(ncont)=om1
  83             omt2(ncont)=om2
  84             omt12(ncont)=om12
  85             ddsc(ncont)=1.0d0/rij
  86             ddla(ncont)=dla
  87             ddlb(ncont)=dlb
  88           endif
  89         enddo
  90       enddo
  91       if (lprint) then
  92         write (iout,'(a)') 'Contact map:'
  93         do i=1,ncont
  94           i1=icont(1,i)
  95           i2=icont(2,i)
  96           it1=itype(i1)
  97           it2=itype(i2)
  98           write (iout,'(i3,2x,a,i4,2x,a,i4,5f8.3,3f10.5)')
  99      &     i,restyp(it1),i1,restyp(it2),i2,cscore(i),
 100      &     sc_cutoff(iabs(it1),iabs(it2)),ddsc(i),ddla(i),ddlb(i),
 101      &     omt1(i),omt2(i),omt12(i)
 102         enddo
 103       endif
 104       return
 105       end
 106 c----------------------------------------------------------------------------
 107       double precision function contact_fract(ncont,ncont_ref,
 108      &                                     icont,icont_ref)
 109       implicit none
 110       include 'DIMENSIONS'
 111       include 'COMMON.IOUNITS'
 112       integer i,j,nmatch
 113       integer ncont,ncont_ref,icont(2,maxcont),icont_ref(2,maxcont)
 114       nmatch=0
 115 c     print *,'ncont=',ncont,' ncont_ref=',ncont_ref
 116 c     write (iout,'(20i4)') (icont_ref(1,i),i=1,ncont_ref)
 117 c     write (iout,'(20i4)') (icont_ref(2,i),i=1,ncont_ref)
 118 c     write (iout,'(20i4)') (icont(1,i),i=1,ncont)
 119 c     write (iout,'(20i4)') (icont(2,i),i=1,ncont)
 120       do i=1,ncont
 121         do j=1,ncont_ref
 122           if (icont(1,i).eq.icont_ref(1,j) .and.
 123      &        icont(2,i).eq.icont_ref(2,j)) nmatch=nmatch+1
 124         enddo
 125       enddo
 126 c     print *,' nmatch=',nmatch
 127 c     contact_fract=dfloat(nmatch)/dfloat(max0(ncont,ncont_ref))
 128       contact_fract=dfloat(nmatch)/dfloat(ncont_ref)
 129       return
 130       end
 131 c------------------------------------------------------------------------------
 132       subroutine pept_cont(lprint,ncont,icont)
 133       implicit none
 134       include 'DIMENSIONS'
 135       include 'DIMENSIONS.ZSCOPT'
 136       include 'COMMON.IOUNITS'
 137       include 'COMMON.CHAIN'
 138       include 'COMMON.INTERACT'
 139       include 'COMMON.FFIELD'
 140       include 'COMMON.NAMES'
 141       integer ncont,icont(2,maxcont)
 142       integer i,j,k,kkk,i1,i2,it1,it2
 143       logical lprint
 144       real*8 dist
 145       real*8 rcomp /5.5d0/
 146       ncont=0
 147       kkk=0
 148       print *,'Entering pept_cont: nnt=',nnt,' nct=',nct
 149       do i=nnt,nct-3
 150         do k=1,3
 151           c(k,2*nres+1)=0.5d0*(c(k,i)+c(k,i+1))
 152         enddo
 153         do j=i+2,nct-1
 154           do k=1,3
 155             c(k,2*nres+2)=0.5d0*(c(k,j)+c(k,j+1))
 156           enddo
 157           if (dist(2*nres+1,2*nres+2).lt.rcomp) then
 158             ncont=ncont+1
 159             icont(1,ncont)=i
 160             icont(2,ncont)=j
 161           endif
 162         enddo
 163       enddo
 164       if (lprint) then
 165         write (iout,'(a)') 'PP contact map:'
 166         do i=1,ncont
 167           i1=icont(1,i)
 168           i2=icont(2,i)
 169           it1=itype(i1)
 170           it2=itype(i2)
 171           write (iout,'(i3,2x,a,i4,2x,a,i4)')
 172      &     i,restyp(it1),i1,restyp(it2),i2
 173         enddo
 174       endif
 175       return
 176       end