--- /dev/null
+ subroutine add_angpair(ici,icj,nang_pair,iang_pair)
+ implicit real*8 (a-h,o-z)
+ include 'DIMENSIONS'
+ include 'COMMON.IOUNITS'
+ include 'COMMON.CHAIN'
+ integer ici,icj,nang_pair,iang_pair(2,maxres)
+ integer i,ian1,ian2
+c write (iout,*) "add_angpair: ici",ici," icj",icj,
+c & " nang_pair",nang_pair
+ ian1=ici+2
+ if (ian1.lt.4 .or. ian1.gt.nres) return
+ ian2=icj+2
+c write (iout,*) "ian1",ian1," ian2",ian2
+ if (ian2.lt.4 .or. ian2.gt.nres) return
+ do i=1,nang_pair
+ if (ian1.eq.iang_pair(1,i) .and. ian2.eq.iang_pair(2,i)) return
+ enddo
+ nang_pair=nang_pair+1
+ iang_pair(1,nang_pair)=ian1
+ iang_pair(2,nang_pair)=ian2
+ return
+ end
+c-------------------------------------------------------------------------
+ subroutine angnorm(jfrag,ishif1,ishif2,diffang_max,angn,fract,
+ & lprn)
+ implicit real*8 (a-h,o-z)
+ include 'DIMENSIONS'
+ include 'DIMENSIONS.ZSCOPT'
+ include 'DIMENSIONS.COMPAR'
+ include 'COMMON.IOUNITS'
+ include 'COMMON.VAR'
+ include 'COMMON.COMPAR'
+ include 'COMMON.CHAIN'
+ include 'COMMON.GEO'
+ double precision pinorm,deltang
+ logical lprn
+ if (lprn) write (iout,'(80(1h*))')
+ angn=0.0d0
+ nn = 0
+ fract = 1.0d0
+ npart = npiece(jfrag,1)
+ nn4 = nstart_sup+3
+ nne = min0(nend_sup,nres)
+ if (lprn) write (iout,*) "nn4",nn4," nne",nne
+ do i=1,npart
+ nbeg = ifrag(1,i,jfrag) + 3 - ishif1
+ if (nbeg.lt.nn4) nbeg=nn4
+ nend = ifrag(2,i,jfrag) + 1 - ishif2
+ if (nend.gt.nne) nend=nne
+ if (nend.ge.nbeg) then
+ nn = nn + nend - nbeg + 1
+ if (lprn) write (iout,*) "i=",i," nbeg",nbeg," nend",nend,
+ & " nn",nn," ishift1",ishif1," ishift2",ishif2
+ if (lprn) write (iout,*) "angles"
+ longest=0
+ ll = 0
+ do j=nbeg,nend
+c deltang = pinorm(phi(j)-phi_ref(j+ishif1))
+ deltang=spherang(phi_ref(j+ishif1),theta_ref(j-1+ishif1),
+ & theta_ref(j+ishif1),phi(j),theta(j-1),theta(j))
+ if (dabs(deltang).gt.diffang_max) then
+ if (ll.gt.longest) longest = ll
+ ll = 0
+ else
+ ll=ll+1
+ endif
+ if (ll.gt.longest) longest = ll
+ if (lprn) write (iout,'(i5,3f10.5)')j,rad2deg*phi(j),
+ & rad2deg*phi_ref(j+ishif1),rad2deg*deltang
+ angn=angn+dabs(deltang)
+ enddo
+ longest=longest+3
+ ff = dfloat(longest)/dfloat(nend - nbeg + 4)
+ if (lprn) write (iout,*)"segment",i," longest fragment within",
+ & diffang_max*rad2deg,":",longest," fraction",ff
+ if (ff.lt.fract) fract = ff
+ endif
+ enddo
+ if (nn.gt.0) then
+ angn = angn/nn
+ else
+ angn = dwapi
+ endif
+ if (lprn) write (iout,*) "nn",nn," norm",rad2deg*angn,
+ & " fract",fract
+ return
+ end
+c-------------------------------------------------------------------------
+ subroutine angnorm2(jfrag,ishif1,ishif2,ncont,icont,lprn,
+ & diffang_max,anorm,fract)
+ implicit real*8 (a-h,o-z)
+ include 'DIMENSIONS'
+ include 'DIMENSIONS.ZSCOPT'
+ include 'DIMENSIONS.COMPAR'
+ include 'COMMON.IOUNITS'
+ include 'COMMON.VAR'
+ include 'COMMON.COMPAR'
+ include 'COMMON.CHAIN'
+ include 'COMMON.GEO'
+ integer ncont,icont(2,ncont),longest
+ double precision anorm,diffang_max,fract
+ integer npiece_c,ifrag_c(2,maxpiece),ishift_c(maxpiece)
+ double precision pinorm
+ logical lprn
+ if (lprn) write (iout,'(80(1h*))')
+c
+c Determine the segments for which angles will be compared
+c
+ nn4 = nstart_sup+3
+ nne = min0(nend_sup,nres)
+ if (lprn) write (iout,*) "nn4",nn4," nne",nne
+ npart=npiece(jfrag,1)
+ npiece_c=0
+ do i=1,npart
+c write (iout,*) "i",i," ifrag",ifrag(1,i,jfrag),ifrag(2,i,jfrag)
+ if (icont(1,ncont).lt.ifrag(1,i,jfrag) .or.
+ & icont(1,1).gt.ifrag(2,i,jfrag)) goto 11
+ jstart=1
+ do while (jstart.lt.ncont .and.
+ & icont(1,jstart).lt.ifrag(1,i,jfrag))
+c write (iout,*) "jstart",jstart," icont",icont(1,jstart),
+c & " ifrag",ifrag(1,i,jfrag)
+ jstart=jstart+1
+ enddo
+c write (iout,*) "jstart",jstart," icont",icont(1,jstart),
+c & " ifrag",ifrag(1,i,jfrag)
+ if (icont(1,jstart).lt.ifrag(1,i,jfrag)) goto 11
+ npiece_c=npiece_c+1
+ ic1=icont(1,jstart)
+ ifrag_c(1,npiece_c)=icont(1,jstart)
+ jend=ncont
+ do while (jend.gt.1 .and. icont(1,jend).gt.ifrag(2,i,jfrag))
+c write (iout,*) "jend",jend," icont",icont(1,jend),
+c & " ifrag",ifrag(2,i,jfrag)
+ jend=jend-1
+ enddo
+c write (iout,*) "jend",jend," icont",icont(1,jend),
+c & " ifrag",ifrag(2,i,jfrag)
+ ic2=icont(1,jend)
+ ifrag_c(2,npiece_c)=icont(1,jend)+1
+ ishift_c(npiece_c)=ishif1
+c write (iout,*) "1: i",i," jstart:",jstart," jend",jend,
+c & " ic1",ic1," ic2",ic2,
+c & " ifrag",ifrag(1,i,jfrag),ifrag(2,i,jfrag)
+ 11 continue
+ if (ncont.eq.1 .or. icont(2,ncont).gt.icont(2,1)) then
+ idi=1
+ else
+ idi=-1
+ endif
+c write (iout,*) "idi",idi
+ if (idi.eq.1) then
+ if (icont(2,1).gt.ifrag(2,i,jfrag) .or.
+ & icont(2,ncont).lt.ifrag(1,i,jfrag)) goto 12
+ jstart=1
+ do while (jstart.lt.ncont .and.
+ & icont(2,jstart).lt.ifrag(1,i,jfrag))
+c write (iout,*) "jstart",jstart," icont",icont(2,jstart),
+c & " ifrag",ifrag(1,i,jfrag)
+ jstart=jstart+1
+ enddo
+c write (iout,*) "jstart",jstart," icont",icont(2,jstart),
+c & " ifrag",ifrag(1,i,jfrag)
+ if (icont(2,jstart).lt.ifrag(1,i,jfrag)) goto 12
+ npiece_c=npiece_c+1
+ ic1=icont(2,jstart)
+ ifrag_c(2,npiece_c)=icont(2,jstart)+1
+ jend=ncont
+ do while (jend.gt.1 .and. icont(2,jend).gt.ifrag(2,i,jfrag))
+c write (iout,*) "jend",jend," icont",icont(2,jend),
+c & " ifrag",ifrag(2,i,jfrag)
+ jend=jend-1
+ enddo
+c write (iout,*) "jend",jend," icont",icont(2,jend),
+c & " ifrag",ifrag(2,i,jfrag)
+ else if (idi.eq.-1) then
+ if (icont(2,ncont).gt.ifrag(2,i,jfrag) .or.
+ & icont(2,1).lt.ifrag(1,i,jfrag)) goto 12
+ jstart=ncont
+ do while (jstart.gt.ncont .and.
+ & icont(2,jstart).lt.ifrag(1,i,jfrag))
+c write (iout,*) "jstart",jstart," icont",icont(2,jstart),
+c & " ifrag",ifrag(1,i,jfrag)
+ jstart=jstart-1
+ enddo
+c write (iout,*) "jstart",jstart," icont",icont(2,jstart),
+c & " ifrag",ifrag(1,i,jfrag)
+ if (icont(2,jstart).lt.ifrag(1,i,jfrag)) goto 12
+ npiece_c=npiece_c+1
+ ic1=icont(2,jstart)
+ ifrag_c(2,npiece_c)=icont(2,jstart)+1
+ jend=1
+ do while (jend.lt.ncont .and.
+ & icont(2,jend).gt.ifrag(2,i,jfrag))
+c write (iout,*) "jend",jend," icont",icont(2,jend),
+c & " ifrag",ifrag(2,i,jfrag)
+ jend=jend+1
+ enddo
+c write (iout,*) "jend",jend," icont",icont(2,jend),
+c & " ifrag",ifrag(2,i,jfrag)
+ endif
+ ic2=icont(2,jend)
+ if (ic2.lt.ic1) then
+ iic = ic1
+ ic1 = ic2
+ ic2 = iic
+ endif
+c write (iout,*) "2: i",i," ic1",ic1," ic2",ic2,
+c & " jstart:",jstart," jend",jend,
+c & " ifrag",ifrag(1,i,jfrag),ifrag(2,i,jfrag)
+ ifrag_c(1,npiece_c)=ic1
+ ifrag_c(2,npiece_c)=ic2+1
+ ishift_c(npiece_c)=ishif2
+ 12 continue
+ enddo
+ if (lprn) then
+ write (iout,*) "Before merge: npiece_c",npiece_c
+ do i=1,npiece_c
+ write (iout,*) ifrag_c(1,i),ifrag_c(2,i),ishift_c(i)
+ enddo
+ endif
+c
+c Merge overlapping segments (e.g., avoid splitting helices)
+c
+ i=1
+ do while (i .lt. npiece_c)
+ if (ishift_c(i).eq.ishift_c(i+1) .and.
+ & ifrag_c(2,i).gt.ifrag_c(1,i+1)) then
+ ifrag_c(2,i)=ifrag_c(2,i+1)
+ do j=i+1,npiece_c
+ ishift_c(j)=ishift_c(j+1)
+ ifrag_c(1,j)=ifrag_c(1,j+1)
+ ifrag_c(2,j)=ifrag_c(2,j+1)
+ enddo
+ npiece_c=npiece_c-1
+ else
+ i=i+1
+ endif
+ enddo
+ if (lprn) then
+ write (iout,*) "After merge: npiece_c",npiece_c
+ do i=1,npiece_c
+ write (iout,*) ifrag_c(1,i),ifrag_c(2,i),ishift_c(i)
+ enddo
+ endif
+c
+c Compare angles
+c
+ angn=0.0d0
+ anorm=0
+ nn = 0
+ fract = 1.0d0
+ npart = npiece_c
+ do i=1,npart
+ ishifc=ishift_c(i)
+ nbeg = ifrag_c(1,i) + 3 - ishifc
+ if (nbeg.lt.nn4) nbeg=nn4
+ nend = ifrag_c(2,i) - ishifc + 1
+ if (nend.gt.nne) nend=nne
+ if (nend.ge.nbeg) then
+ nn = nn + nend - nbeg + 1
+ if (lprn) write (iout,*) "i=",i," nbeg",nbeg," nend",nend,
+ & " nn",nn," ishifc",ishifc
+ if (lprn) write (iout,*) "angles"
+ longest=0
+ ll = 0
+ do j=nbeg,nend
+c deltang = pinorm(phi(j)-phi_ref(j+ishifc))
+ deltang=spherang(phi_ref(j+ishifc),theta_ref(j-1+ishifc),
+ & theta_ref(j+ishifc),phi(j),theta(j-1),theta(j))
+ if (dabs(deltang).gt.diffang_max) then
+ if (ll.gt.longest) longest = ll
+ ll = 0
+ else
+ ll=ll+1
+ endif
+ if (ll.gt.longest) longest = ll
+ if (lprn) write (iout,'(i5,3f10.5)')j,rad2deg*phi(j),
+ & rad2deg*phi_ref(j+ishifc),rad2deg*deltang
+ angn=angn+dabs(deltang)
+ enddo
+ longest=longest+3
+ ff = dfloat(longest)/dfloat(nend - nbeg + 4)
+ if (lprn) write (iout,*)"segment",i," longest fragment within",
+ & diffang_max*rad2deg,":",longest," fraction",ff
+ if (ff.lt.fract) fract = ff
+ endif
+ enddo
+ if (nn.gt.0) anorm = angn/nn
+ if (lprn) write (iout,*) "nn",nn," norm",anorm," fract:",fract
+ return
+ end
+c-------------------------------------------------------------------------
+ double precision function angnorm1(nang_pair,iang_pair,lprn)
+ implicit real*8 (a-h,o-z)
+ include 'DIMENSIONS'
+ include 'DIMENSIONS.ZSCOPT'
+ include 'DIMENSIONS.COMPAR'
+ include 'COMMON.IOUNITS'
+ include 'COMMON.VAR'
+ include 'COMMON.COMPAR'
+ include 'COMMON.CHAIN'
+ include 'COMMON.GEO'
+ logical lprn
+ integer nang_pair,iang_pair(2,maxres)
+ double precision pinorm
+ angn=0.0d0
+ if (lprn) write (iout,'(80(1h*))')
+ if (lprn) write (iout,*) "nang_pair",nang_pair
+ if (lprn) write (iout,*) "angles"
+ do j=1,nang_pair
+ ia1 = iang_pair(1,j)
+ ia2 = iang_pair(2,j)
+c deltang = pinorm(phi(ia1)-phi_ref(ia2))
+ deltang=spherang(phi_ref(ia2),theta_ref(ia2-1),
+ & theta_ref(ia2),phi(ia2),theta(ia2-1),theta(ia2))
+ if (lprn) write (iout,'(3i5,3f10.5)')j,ia1,ia2,rad2deg*phi(ia1),
+ & rad2deg*phi_ref(ia2),rad2deg*deltang
+ angn=angn+dabs(deltang)
+ enddo
+ if (lprn)
+ &write (iout,*)"nang_pair",nang_pair," angn",rad2deg*angn/nang_pair
+ angnorm1 = angn/nang_pair
+ return
+ end
+c------------------------------------------------------------------------------
+ subroutine angnorm12(diff)
+ implicit real*8 (a-h,o-z)
+ include 'DIMENSIONS'
+ include 'DIMENSIONS.ZSCOPT'
+ include 'DIMENSIONS.COMPAR'
+ include 'COMMON.IOUNITS'
+ include 'COMMON.VAR'
+ include 'COMMON.COMPAR'
+ include 'COMMON.CHAIN'
+ include 'COMMON.GEO'
+ double precision pinorm
+ diff=0.0d0
+ nn4 = nstart_sup+3
+ nne = min0(nend_sup,nres)
+c do j=nn4-1,nne
+c diff = diff+rad2deg*dabs(pinorm(theta(j)-theta_ref(j)))
+c enddo
+ do j=nn4,nne
+c diff = diff+rad2deg*dabs(pinorm(phi(j)-phi_ref(j)))
+ diff=diff+spherang(phi_ref(j),theta_ref(j-1),
+ & theta_ref(j),phi(j),theta(j-1),theta(j))
+ enddo
+ return
+ end
+c--------------------------------------------------------------------------------
+ double precision function spherang(gam1,theta11,theta12,
+ & gam2,theta21,theta22)
+ implicit none
+ double precision gam1,theta11,theta12,gam2,theta21,theta22,
+ & x1,x2,xmed,f1,f2,fmed
+ double precision tolx /1.0d-4/, tolf /1.0d-4/
+ double precision sumcos
+ double precision arcos,pinorm,sumangp
+ integer it,maxit /100/
+c Calculate the difference of the angles of two superposed 4-redidue fragments
+c
+c O P
+c \ /
+c O'--C--C
+c \
+c P'
+c
+c The fragment O'-C-C-P' is rotated by angle fi about the C-C axis
+c to achieve the minimum difference between the O'-C-O and P-C-P angles;
+c the sum of these angles is the difference returned by the function.
+c
+c 4/28/04 AL
+c If thetas match, take the difference of gamma and exit.
+ if (dabs(theta11-theta12).lt.tolx
+ & .and. dabs(theta21-theta22).lt.tolx) then
+ spherang=dabs(pinorm(gam2-gam1))
+ return
+ endif
+c If the gammas are the same, take the difference of thetas and exit.
+ x1=0.0d0
+ x2=0.5d0*pinorm(gam2-gam1)
+ if (dabs(x2) .lt. tolx) then
+ spherang=dabs(theta11-theta21)+dabs(theta12-theta22)
+ return
+ else if (x2.lt.0.0d0) then
+ x1=x2
+ x2=0.0d0
+ endif
+c Else apply regula falsi method to compute optimum overlap of the terminal Calphas
+ f1=sumangp(gam1,theta11,theta12,gam2,theta21,theta22,x1)
+ f2=sumangp(gam1,theta11,theta12,gam2,theta21,theta22,x2)
+ do it=1,maxit
+ xmed=x1-f1*(x2-x1)/(f2-f1)
+ fmed=sumangp(gam1,theta11,theta12,gam2,theta21,theta22,xmed)
+c write (*,*) 'it',it,' xmed ',xmed,' fmed ',fmed
+ if ( (dabs(xmed-x1).lt.tolx .or. dabs(x2-xmed).lt.tolx)
+ & .and. dabs(fmed).lt.tolf ) then
+ x1=xmed
+ f1=fmed
+ goto 10
+ else if ( fmed*f1.lt.0.0d0 ) then
+ x2=xmed
+ f2=fmed
+ else
+ x1=xmed
+ f1=fmed
+ endif
+ enddo
+ 10 continue
+ spherang=arcos(dcos(theta11)*dcos(theta12)
+ & +dsin(theta11)*dsin(theta12)*dcos(x1))+
+ & arcos(dcos(theta21)*dcos(theta22)+
+ & dsin(theta21)*dsin(theta22)*dcos(gam2-gam1+x1))
+ return
+ end
+c--------------------------------------------------------------------------------
+ double precision function sumangp(gam1,theta11,theta12,gam2,
+ & theta21,theta22,fi)
+ implicit none
+ double precision gam1,theta11,theta12,gam2,theta21,theta22,fi,
+ & cost11,cost12,cost21,cost22,sint11,sint12,sint21,sint22,cosd1,
+ & cosd2
+c derivarive of the sum of the difference of the angles of a 4-residue fragment.
+ double precision arcos
+ cost11=dcos(theta11)
+ cost12=dcos(theta12)
+ cost21=dcos(theta21)
+ cost22=dcos(theta22)
+ sint11=dsin(theta11)
+ sint12=dsin(theta12)
+ sint21=dsin(theta21)
+ sint22=dsin(theta22)
+ cosd1=cost11*cost12+sint11*sint12*dcos(fi)
+ cosd2=cost21*cost22+sint21*sint22*dcos(gam2-gam1+fi)
+ sumangp=sint11*sint12*dsin(fi)/dsqrt(1.0d0-cosd1*cosd1)
+ & +sint21*sint22*dsin(gam2-gam1+fi)/dsqrt(1.0d0-cosd2*cosd2)
+ return
+ end
projects / unres.git / blobdiff