1 subroutine add_angpair(ici,icj,nang_pair,iang_pair)
4 include 'DIMENSIONS.ZSCOPT'
5 include 'COMMON.IOUNITS'
7 integer ici,icj,nang_pair,iang_pair(2,maxres)
9 c write (iout,*) "add_angpair: ici",ici," icj",icj,
10 c & " nang_pair",nang_pair
12 if (ian1.lt.4 .or. ian1.gt.nres) return
14 c write (iout,*) "ian1",ian1," ian2",ian2
15 if (ian2.lt.4 .or. ian2.gt.nres) return
17 if (ian1.eq.iang_pair(1,i) .and. ian2.eq.iang_pair(2,i)) return
20 iang_pair(1,nang_pair)=ian1
21 iang_pair(2,nang_pair)=ian2
24 c-------------------------------------------------------------------------
25 subroutine angnorm(jfrag,ishif1,ishif2,diffang_max,angn,fract,
29 include 'DIMENSIONS.ZSCOPT'
30 include 'DIMENSIONS.COMPAR'
31 include 'COMMON.IOUNITS'
33 include 'COMMON.COMPAR'
34 include 'COMMON.CHAIN'
36 integer jfrag,ishif1,ishif2,ib,iprot
37 integer nn,npart,nn4,nne,i,nbeg,nend,longest,ll,j
38 double precision deltang,diffang_max,angn,fract
40 double precision pinorm,spherang
42 if (lprn) write (iout,'(80(1h*))')
46 npart = npiece(jfrag,1,ib,iprot)
47 nn4 = nstart_sup(iprot)+3
48 nne = min0(nend_sup(iprot),nres)
49 if (lprn) write (iout,*) "nn4",nn4," nne",nne
51 nbeg = ifrag(1,i,jfrag,ib,iprot) + 3 - ishif1
52 if (nbeg.lt.nn4) nbeg=nn4
53 nend = ifrag(2,i,jfrag,ib,iprot) + 1 - ishif2
54 if (nend.gt.nne) nend=nne
55 if (nend.ge.nbeg) then
56 nn = nn + nend - nbeg + 1
57 if (lprn) write (iout,*) "i=",i," nbeg",nbeg," nend",nend,
58 & " nn",nn," ishift1",ishif1," ishift2",ishif2
59 if (lprn) write (iout,*) "angles"
63 c deltang = pinorm(phi(j)-phi_ref(j+ishif1,iprot))
64 deltang=spherang(phi_ref(j+ishif1,iprot),
65 & theta_ref(j-1+ishif1,iprot),theta_ref(j+ishif1,iprot),
66 & phi(j),theta(j-1),theta(j))
67 if (dabs(deltang).gt.diffang_max) then
68 if (ll.gt.longest) longest = ll
73 if (ll.gt.longest) longest = ll
74 if (lprn) write (iout,'(i5,3f10.5)')j,rad2deg*phi(j),
75 & rad2deg*phi_ref(j+ishif1,iprot),rad2deg*deltang
76 angn=angn+dabs(deltang)
79 ff = dfloat(longest)/dfloat(nend - nbeg + 4)
80 if (lprn) write (iout,*)"segment",i," longest fragment within",
81 & diffang_max*rad2deg,":",longest," fraction",ff
82 if (ff.lt.fract) fract = ff
90 if (lprn) write (iout,*) "nn",nn," norm",rad2deg*angn,
94 c-------------------------------------------------------------------------
95 subroutine angnorm2(jfrag,ishif1,ishif2,ncont,icont,lprn,
96 & diffang_max,anorm,fract,ib,iprot)
99 include 'DIMENSIONS.ZSCOPT'
100 include 'DIMENSIONS.COMPAR'
101 include 'COMMON.IOUNITS'
103 include 'COMMON.COMPAR'
104 include 'COMMON.CHAIN'
107 integer ncont,icont(2,ncont),longest,jfrag,ishif1,ishif2
108 integer npiece_c,ifrag_c(2,maxpiece),ishift_c(maxpiece)
109 integer nn4,nne,npart,i,jstart,ic1,jend,ic2,idi,iic,j,
110 & nn,ishifc,nbeg,nend,ll
111 double precision diffang_max,anorm,fract
112 double precision angn, deltang,ff
113 double precision pinorm,spherang
115 if (lprn) write (iout,'(80(1h*))')
117 c Determine the segments for which angles will be compared
119 nn4 = nstart_sup(iprot)+3
120 nne = min0(nend_sup(iprot),nres)
121 if (lprn) write (iout,*) "nn4",nn4," nne",nne
122 npart=npiece(jfrag,1,ib,iprot)
125 c write (iout,*) "i",i," ifrag",ifrag(1,i,jfrag,ib,iprot),
126 c & ifrag(2,i,jfrag,ib,iprot)
127 if (icont(1,ncont).lt.ifrag(1,i,jfrag,ib,iprot) .or.
128 & icont(1,1).gt.ifrag(2,i,jfrag,ib,iprot)) goto 11
130 do while (jstart.lt.ncont .and.
131 & icont(1,jstart).lt.ifrag(1,i,jfrag,ib,iprot))
132 c write (iout,*) "jstart",jstart," icont",icont(1,jstart),
133 c & " ifrag",ifrag(1,i,jfrag,ib,iprot)
136 c write (iout,*) "jstart",jstart," icont",icont(1,jstart),
137 c & " ifrag",ifrag(1,i,jfrag,ib,iprot)
138 if (icont(1,jstart).lt.ifrag(1,i,jfrag,ib,iprot)) goto 11
141 ifrag_c(1,npiece_c)=icont(1,jstart)
143 do while(jend.gt.1.and.icont(1,jend).gt.
144 & ifrag(2,i,jfrag,ib,iprot))
145 c write (iout,*) "jend",jend," icont",icont(1,jend),
146 c & " ifrag",ifrag(2,i,jfrag,ib,iprot)
149 c write (iout,*) "jend",jend," icont",icont(1,jend),
150 c & " ifrag",ifrag(2,i,jfrag,ib,iprot)
152 ifrag_c(2,npiece_c)=icont(1,jend)+1
153 ishift_c(npiece_c)=ishif1
154 c write (iout,*) "1: i",i," jstart:",jstart," jend",jend,
155 c & " ic1",ic1," ic2",ic2,
156 c & " ifrag",ifrag(1,i,jfrag,ib,iprot),ifrag(2,i,jfrag,ib,iprot)
158 if (ncont.eq.1 .or. icont(2,ncont).gt.icont(2,1)) then
163 c write (iout,*) "idi",idi
165 if (icont(2,1).gt.ifrag(2,i,jfrag,ib,iprot) .or.
166 & icont(2,ncont).lt.ifrag(1,i,jfrag,ib,iprot)) goto 12
168 do while (jstart.lt.ncont .and.
169 & icont(2,jstart).lt.ifrag(1,i,jfrag,ib,iprot))
170 c write (iout,*) "jstart",jstart," icont",icont(2,jstart),
171 c & " ifrag",ifrag(1,i,jfrag,ib,iprot)
174 c write (iout,*) "jstart",jstart," icont",icont(2,jstart),
175 c & " ifrag",ifrag(1,i,jfrag,ib,iprot)
176 if (icont(2,jstart).lt.ifrag(1,i,jfrag,ib,iprot)) goto 12
179 ifrag_c(2,npiece_c)=icont(2,jstart)+1
181 do while(jend.gt.1.and.icont(2,jend).gt.
182 & ifrag(2,i,jfrag,ib,iprot))
183 c write (iout,*) "jend",jend," icont",icont(2,jend),
184 c & " ifrag",ifrag(2,i,jfrag,ib,iprot)
187 c write (iout,*) "jend",jend," icont",icont(2,jend),
188 c & " ifrag",ifrag(2,i,jfrag,ib,iprot)
189 else if (idi.eq.-1) then
190 if (icont(2,ncont).gt.ifrag(2,i,jfrag,ib,iprot) .or.
191 & icont(2,1).lt.ifrag(1,i,jfrag,ib,iprot)) goto 12
193 do while (jstart.gt.ncont .and.
194 & icont(2,jstart).lt.ifrag(1,i,jfrag,ib,iprot))
195 c write (iout,*) "jstart",jstart," icont",icont(2,jstart),
196 c & " ifrag",ifrag(1,i,jfrag,ib,iprot)
199 c write (iout,*) "jstart",jstart," icont",icont(2,jstart),
200 c & " ifrag",ifrag(1,i,jfrag,ib,iprot)
201 if (icont(2,jstart).lt.ifrag(1,i,jfrag,ib,iprot)) goto 12
204 ifrag_c(2,npiece_c)=icont(2,jstart)+1
206 do while (jend.lt.ncont .and.
207 & icont(2,jend).gt.ifrag(2,i,jfrag,ib,iprot))
208 c write (iout,*) "jend",jend," icont",icont(2,jend),
209 c & " ifrag",ifrag(2,i,jfrag,ib,iprot)
212 c write (iout,*) "jend",jend," icont",icont(2,jend),
213 c & " ifrag",ifrag(2,i,jfrag,ib,iprot)
221 c write (iout,*) "2: i",i," ic1",ic1," ic2",ic2,
222 c & " jstart:",jstart," jend",jend,
223 c & " ifrag",ifrag(1,i,jfrag,ib,iprot),ifrag(2,i,jfrag,ib,iprot)
224 ifrag_c(1,npiece_c)=ic1
225 ifrag_c(2,npiece_c)=ic2+1
226 ishift_c(npiece_c)=ishif2
230 write (iout,*) "Before merge: npiece_c",npiece_c
232 write (iout,*) ifrag_c(1,i),ifrag_c(2,i),ishift_c(i)
236 c Merge overlapping segments (e.g., avoid splitting helices)
239 do while (i .lt. npiece_c)
240 if (ishift_c(i).eq.ishift_c(i+1) .and.
241 & ifrag_c(2,i).gt.ifrag_c(1,i+1)) then
242 ifrag_c(2,i)=ifrag_c(2,i+1)
244 ishift_c(j)=ishift_c(j+1)
245 ifrag_c(1,j)=ifrag_c(1,j+1)
246 ifrag_c(2,j)=ifrag_c(2,j+1)
254 write (iout,*) "After merge: npiece_c",npiece_c
256 write (iout,*) ifrag_c(1,i),ifrag_c(2,i),ishift_c(i)
269 nbeg = ifrag_c(1,i) + 3 - ishifc
270 if (nbeg.lt.nn4) nbeg=nn4
271 nend = ifrag_c(2,i) - ishifc + 1
272 if (nend.gt.nne) nend=nne
273 if (nend.ge.nbeg) then
274 nn = nn + nend - nbeg + 1
275 if (lprn) write (iout,*) "i=",i," nbeg",nbeg," nend",nend,
276 & " nn",nn," ishifc",ishifc
277 if (lprn) write (iout,*) "angles"
281 c deltang = pinorm(phi(j)-phi_ref(j+ishifc,iprot))
282 deltang=spherang(phi_ref(j+ishifc,iprot),
283 & theta_ref(j-1+ishifc,iprot),theta_ref(j+ishifc,iprot),
284 & phi(j),theta(j-1),theta(j))
285 if (dabs(deltang).gt.diffang_max) then
286 if (ll.gt.longest) longest = ll
291 if (ll.gt.longest) longest = ll
292 if (lprn) write (iout,'(i5,3f10.5)')j,rad2deg*phi(j),
293 & rad2deg*phi_ref(j+ishifc,iprot),rad2deg*deltang
294 angn=angn+dabs(deltang)
297 ff = dfloat(longest)/dfloat(nend - nbeg + 4)
298 if (lprn) write (iout,*)"segment",i," longest fragment within",
299 & diffang_max*rad2deg,":",longest," fraction",ff
300 if (ff.lt.fract) fract = ff
303 if (nn.gt.0) anorm = angn/nn
304 if (lprn) write (iout,*) "nn",nn," norm",anorm," fract:",fract
307 c-------------------------------------------------------------------------
308 double precision function angnorm1(nang_pair,iang_pair,iprot,lprn)
311 include 'DIMENSIONS.ZSCOPT'
312 include 'DIMENSIONS.COMPAR'
313 include 'COMMON.IOUNITS'
315 include 'COMMON.COMPAR'
316 include 'COMMON.CHAIN'
319 integer nang_pair,iang_pair(2,maxres)
321 double precision angn,deltang
323 double precision pinorm,spherang
325 if (lprn) write (iout,'(80(1h*))')
326 if (lprn) write (iout,*) "nang_pair",nang_pair
327 if (lprn) write (iout,*) "angles"
331 c deltang = pinorm(phi(ia1)-phi_ref(ia2,iprot))
332 deltang=spherang(phi_ref(ia2,iprot),theta_ref(ia2-1,iprot),
333 & theta_ref(ia2,iprot),phi(ia2),theta(ia2-1),theta(ia2))
334 if (lprn) write (iout,'(3i5,3f10.5)')j,ia1,ia2,rad2deg*phi(ia1),
335 & rad2deg*phi_ref(ia2,iprot),rad2deg*deltang
336 angn=angn+dabs(deltang)
339 &write (iout,*)"nang_pair",nang_pair," angn",rad2deg*angn/nang_pair
340 angnorm1 = angn/nang_pair
343 c------------------------------------------------------------------------------
344 subroutine angnorm12(diff,iprot)
347 include 'DIMENSIONS.ZSCOPT'
348 include 'DIMENSIONS.COMPAR'
349 include 'COMMON.IOUNITS'
351 include 'COMMON.COMPAR'
352 include 'COMMON.CHAIN'
354 integer iprot,nn4,nne,j
355 double precision diff
356 double precision pinorm,spherang
358 nn4 = nstart_sup(iprot)+3
359 nne = min0(nend_sup(iprot),nres)
361 c diff = diff+dabs(pinorm(theta(j)-theta_ref(j,iprot)))
364 c diff = diff+dabs(pinorm(phi(j)-phi_ref(j,iprot)))
365 diff=diff+spherang(phi_ref(j,iprot),theta_ref(j-1,iprot),
366 & theta_ref(j,iprot),phi(j),theta(j-1),theta(j))
368 diff=diff/(2*(nne-nn4)+3)
371 c--------------------------------------------------------------------------------
372 double precision function spherang(gam1,theta11,theta12,
373 & gam2,theta21,theta22)
375 double precision gam1,theta11,theta12,gam2,theta21,theta22,
376 & x1,x2,xmed,f1,f2,fmed
377 double precision tolx /1.0d-4/, tolf /1.0d-4/
378 double precision sumcos
379 double precision arcos,pinorm,sumangp
380 integer it,maxit /100/
381 c Calculate the difference of the angles of two superposed 4-redidue fragments
389 c The fragment O'-C-C-P' is rotated by angle fi about the C-C axis
390 c to achieve the minimum difference between the O'-C-O and P-C-P angles;
391 c the sum of these angles is the difference returned by the function.
394 c If thetas match, take the difference of gamma and exit.
395 if (dabs(theta11-theta12).lt.tolx
396 & .and. dabs(theta21-theta22).lt.tolx) then
397 spherang=dabs(pinorm(gam2-gam1))
400 c If the gammas are the same, take the difference of thetas and exit.
402 x2=0.5d0*pinorm(gam2-gam1)
403 if (dabs(x2) .lt. tolx) then
404 spherang=dabs(theta11-theta21)+dabs(theta12-theta22)
406 else if (x2.lt.0.0d0) then
410 c Else apply regula falsi method to compute optimum overlap of the terminal Calphas
411 f1=sumangp(gam1,theta11,theta12,gam2,theta21,theta22,x1)
412 f2=sumangp(gam1,theta11,theta12,gam2,theta21,theta22,x2)
414 xmed=x1-f1*(x2-x1)/(f2-f1)
415 fmed=sumangp(gam1,theta11,theta12,gam2,theta21,theta22,xmed)
416 c write (*,*) 'it',it,' xmed ',xmed,' fmed ',fmed
417 if ( (dabs(xmed-x1).lt.tolx .or. dabs(x2-xmed).lt.tolx)
418 & .and. dabs(fmed).lt.tolf ) then
422 else if ( fmed*f1.lt.0.0d0 ) then
431 spherang=arcos(dcos(theta11)*dcos(theta12)
432 & +dsin(theta11)*dsin(theta12)*dcos(x1))+
433 & arcos(dcos(theta21)*dcos(theta22)+
434 & dsin(theta21)*dsin(theta22)*dcos(gam2-gam1+x1))
437 c--------------------------------------------------------------------------------
438 double precision function sumangp(gam1,theta11,theta12,gam2,
439 & theta21,theta22,fi)
441 double precision gam1,theta11,theta12,gam2,theta21,theta22,fi,
442 & cost11,cost12,cost21,cost22,sint11,sint12,sint21,sint22,cosd1,
444 c derivarive of the sum of the difference of the angles of a 4-residue fragment.
445 double precision arcos
454 cosd1=cost11*cost12+sint11*sint12*dcos(fi)
455 cosd2=cost21*cost22+sint21*sint22*dcos(gam2-gam1+fi)
456 sumangp=sint11*sint12*dsin(fi)/dsqrt(1.0d0-cosd1*cosd1)
457 & +sint21*sint22*dsin(gam2-gam1+fi)/dsqrt(1.0d0-cosd2*cosd2)