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