--- /dev/null
+ 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
projects / unres.git / blobdiff