1 | ! **************************************************************
|
---|
2 | !
|
---|
3 | !
|
---|
4 | ! This file contains the subroutines: addend, redchg, rplgrp
|
---|
5 | !
|
---|
6 | ! Copyright 2003-2005 Frank Eisenmenger, U.H.E. Hansmann,
|
---|
7 | ! Shura Hayryan, Chin-Ku
|
---|
8 | ! Copyright 2007 Frank Eisenmenger, U.H.E. Hansmann,
|
---|
9 | ! Jan H. Meinke, Sandipan Mohanty
|
---|
10 | !
|
---|
11 | ! $Id: addend.f 334 2007-08-07 09:23:59Z meinke $
|
---|
12 | ! **************************************************************
|
---|
13 | subroutine addend(nml,grpn,grpc)
|
---|
14 |
|
---|
15 | ! ..............................................................
|
---|
16 | ! PURPOSE: modify terminal residues to complete bonding scheme
|
---|
17 | ! with residue 'grpn' at N- and residue 'grpc' at C-terminus
|
---|
18 | ! ! need initial co-ordinates for residues to modify
|
---|
19 | ! ! for N-terminus: may add only simple groups
|
---|
20 | !
|
---|
21 | ! CALLS: rplgrp,tolost,redchg
|
---|
22 | ! ..............................................................
|
---|
23 |
|
---|
24 | include 'INCL.H'
|
---|
25 |
|
---|
26 | character grpn*4,grpc*4
|
---|
27 |
|
---|
28 | character res*4,rpat*4,sbrs*4,grn*4,grc*4
|
---|
29 |
|
---|
30 |
|
---|
31 | grn = grpn
|
---|
32 | call tolost(grn)
|
---|
33 | grc = grpc
|
---|
34 | call tolost(grc)
|
---|
35 |
|
---|
36 | if (grn(:3).eq.'ace'.or.grc(:3).eq.'ace'
|
---|
37 | &.or.grn(:3).eq.'nme'.or.grc(:3).eq.'nme') then
|
---|
38 |
|
---|
39 | write(*,'(2a)') ' addend> N-Acetyl (ace) or N-Methylamide (nme)'
|
---|
40 | & ,' should be put in SEQUENCE file, not added as end groups'
|
---|
41 |
|
---|
42 | stop
|
---|
43 | endif
|
---|
44 |
|
---|
45 | ! __________________________________________ N-terminus
|
---|
46 | ifirs=irsml1(nml)
|
---|
47 | rpat='n '
|
---|
48 | res=seq(ifirs)
|
---|
49 | call tolost(res)
|
---|
50 | if (res(:3).ne.'ace') then
|
---|
51 | if (grn(:3).eq.'nh2') then
|
---|
52 | if (res(:3).eq.'hyp'.and..not.flex) then
|
---|
53 | if(sh2) then
|
---|
54 | sbrs='nh2+'
|
---|
55 | else
|
---|
56 | write (*,'(2a)') ' addend> ',
|
---|
57 | & ' No N-terminal Hyp possible with ECEPP/3 dataset'
|
---|
58 | stop
|
---|
59 | endif
|
---|
60 |
|
---|
61 | elseif (res(:3).eq.'pro') then
|
---|
62 | sbrs='nh1 '
|
---|
63 | elseif (res(:3).eq.'gly'.and..not.flex) then
|
---|
64 | if (sh2) then
|
---|
65 | sbrs='nh2 '
|
---|
66 | else
|
---|
67 | sbrs='nh2g'
|
---|
68 | endif
|
---|
69 | else
|
---|
70 | sbrs='nh2 '
|
---|
71 | endif
|
---|
72 |
|
---|
73 | elseif (grn.eq.'nh3+') then
|
---|
74 | if (res(:3).eq.'pro'.and..not.flex) then
|
---|
75 | sbrs='nh2+'
|
---|
76 | elseif (res(:3).eq.'gly'.and..not.flex) then
|
---|
77 | sbrs='nh3g'
|
---|
78 | else
|
---|
79 | sbrs='nh3+'
|
---|
80 | endif
|
---|
81 |
|
---|
82 | else
|
---|
83 |
|
---|
84 | write(*,'(2a)') ' addend> Can add only ',
|
---|
85 | & 'nh2 or nh3+ to N-terminus'
|
---|
86 | stop
|
---|
87 |
|
---|
88 | endif
|
---|
89 |
|
---|
90 | call rplgrp(nml,ifirs,rpat,sbrs)
|
---|
91 | if (flex) call redchg(nml,ifirs,rpat,sbrs) ! Flex dataset
|
---|
92 |
|
---|
93 | else ! ace
|
---|
94 |
|
---|
95 | write(*,'(2a)') ' addend> Acetyl group',
|
---|
96 | & ' at N-terminus not modified'
|
---|
97 | endif
|
---|
98 |
|
---|
99 | ! __________________________________________ C-terminus
|
---|
100 | ilars=irsml2(nml)
|
---|
101 | rpat='c '
|
---|
102 | res=seq(ilars)
|
---|
103 | call tolost(res)
|
---|
104 |
|
---|
105 | if (res(:3).ne.'nme') then
|
---|
106 |
|
---|
107 | if (grc.eq.'cooh') then
|
---|
108 |
|
---|
109 | if (res(:3).eq.'gly'.and..not.sh2) then
|
---|
110 | sbrs='gooh'
|
---|
111 | else
|
---|
112 | sbrs='cooh'
|
---|
113 | endif
|
---|
114 |
|
---|
115 | elseif (grc.eq.'coo-') then
|
---|
116 |
|
---|
117 | if (res(:3).eq.'gly'.and..not.sh2) then
|
---|
118 | sbrs='goo-'
|
---|
119 | else
|
---|
120 | sbrs='coo-'
|
---|
121 | endif
|
---|
122 |
|
---|
123 | else
|
---|
124 |
|
---|
125 | write(*,'(2a)') ' addend> Can add only ',
|
---|
126 | & 'cooh or coo- to C-terminus'
|
---|
127 | stop
|
---|
128 |
|
---|
129 | endif
|
---|
130 |
|
---|
131 | call rplgrp(nml,ilars,rpat,sbrs)
|
---|
132 | if (flex) call redchg(nml,ilars,rpat,sbrs) ! Flex dataset
|
---|
133 |
|
---|
134 | else ! N'-methylamide
|
---|
135 |
|
---|
136 | write(*,'(2a)') ' addend> N-Methylamide',
|
---|
137 | & ' at C-terminus not modified'
|
---|
138 |
|
---|
139 | endif
|
---|
140 |
|
---|
141 | ! ----------------------------- net charge of molecule
|
---|
142 | cg = 0.d0
|
---|
143 | do i=iatrs1(irsml1(nml)),iatrs2(irsml2(nml))
|
---|
144 | cg = cg + cgat(i)
|
---|
145 | enddo
|
---|
146 | if (abs(cg).gt.1.d-5) write(*,'(a,i2,a,f7.3,/)')
|
---|
147 | & ' addend> Net charge of molecule #'
|
---|
148 | & ,nml,': ',cg
|
---|
149 |
|
---|
150 | return
|
---|
151 | end
|
---|
152 | ! ****************************************
|
---|
153 | subroutine rplgrp(nml,nrs,rpat,sbrs)
|
---|
154 |
|
---|
155 | ! ...............................................................
|
---|
156 | ! PURPOSE: replace atom(s) rooted at atom 'rpat' in residue
|
---|
157 | ! 'nrs' of molecule 'nml' by atom(s) rooted at
|
---|
158 | ! 'rpat' of residue 'sbrs' (same name of root
|
---|
159 | ! atom 'rpat' maintains bonding geometry for
|
---|
160 | ! preceeding atoms in 'nrs')
|
---|
161 | !
|
---|
162 | ! is NOT performed if 'rpat' is within mainchain,
|
---|
163 | ! except it is first/last mainchain atom of 'nml'
|
---|
164 | !
|
---|
165 | ! CALLS: dihedr,iopfil,iendst,eyring,fndbrn,redres,setsys,valang
|
---|
166 | ! ...............................................................
|
---|
167 |
|
---|
168 | include 'INCL.H'
|
---|
169 |
|
---|
170 | character rpat*4,sbrs*4
|
---|
171 | logical ntbb,bb
|
---|
172 |
|
---|
173 | dimension ibd(mxbd+1),iybd(mxbd+1)
|
---|
174 |
|
---|
175 |
|
---|
176 | ifirs=irsml1(nml)
|
---|
177 | ilars=irsml2(nml)
|
---|
178 |
|
---|
179 | nxt=ixatrs(nrs)
|
---|
180 | nfi=iatrs1(nrs)
|
---|
181 | nla=iatrs2(nrs)
|
---|
182 | ! __________________________ indices of atoms to be replaced
|
---|
183 | do i=nfi,nla
|
---|
184 | if (rpat.eq.nmat(i)) then
|
---|
185 | nfirp=i
|
---|
186 | goto 1
|
---|
187 | endif
|
---|
188 | enddo
|
---|
189 | write (*,'(4a,i4,a,i4)') ' rplgrp> cannot find atom >',rpat,
|
---|
190 | &'< to be replaced in residue ',seq(nrs),nrs,' of molecule ',nml
|
---|
191 | stop
|
---|
192 |
|
---|
193 | 1 call fndbrn(nml,nrs,nfirp,nlarp,irng1,irng2,bb)
|
---|
194 | if (irng1.gt.0) goto 10
|
---|
195 | ntbb=.false.
|
---|
196 |
|
---|
197 | if (bb) then ! ..... backbone
|
---|
198 |
|
---|
199 | if (nfirp.eq.nxt.and.nrs.eq.ilars) goto 2 ! last mainchain atom
|
---|
200 |
|
---|
201 | if (nfirp.eq.nfi.and.nrs.eq.ifirs) then ! 1st MAINCHAIN ATOM !
|
---|
202 |
|
---|
203 | ntbb=.true.
|
---|
204 |
|
---|
205 | ibd(1)=iowat(nfirp)
|
---|
206 | iybd(1)=iyowat(nfirp)
|
---|
207 |
|
---|
208 | do i=1,mxbd
|
---|
209 | ibd(i+1)=ibdat(i,nfirp)
|
---|
210 | iybd(i+1)=iybdat(i,nfirp)
|
---|
211 | enddo
|
---|
212 |
|
---|
213 | ibdrg=0 ! __________________ check ring
|
---|
214 | iybdrg=1
|
---|
215 |
|
---|
216 | do i=1,nbdat(nfirp)+1
|
---|
217 | if (iowat(ibd(i)).ne.nfirp) then
|
---|
218 | if (ibdrg.ne.0) then
|
---|
219 | write (*,'(2a,i3)')
|
---|
220 | & ' rplgrp> Can handle only simple ring at 1st',
|
---|
221 | & ' atom of molecule #',nml
|
---|
222 | stop
|
---|
223 | endif
|
---|
224 | ibdrg=ibd(i)
|
---|
225 | iybdrg=iybd(i)
|
---|
226 | endif
|
---|
227 | enddo
|
---|
228 | nxtbb1=nlarp+1 ! _________________ next backbone atoms
|
---|
229 | isgbb1=iyowat(nxtbb1)
|
---|
230 |
|
---|
231 | if (nxtbb1.eq.nxt) then
|
---|
232 | if (nrs.lt.ilars) then
|
---|
233 | nxtbb2=iatrs1(nrs+1)
|
---|
234 | goto 3
|
---|
235 | endif
|
---|
236 | else
|
---|
237 | do i=nxt,nxtbb1+1,-1
|
---|
238 | if (iowat(i).eq.nxtbb1) then
|
---|
239 | nxtbb2=i
|
---|
240 | goto 3
|
---|
241 | endif
|
---|
242 | enddo
|
---|
243 | endif
|
---|
244 | goto 11
|
---|
245 | else
|
---|
246 | write (*,'(4a,i4,a,i4)')
|
---|
247 | & ' rplgrp> Cannot replace BACKBONE atom ',rpat,
|
---|
248 | & ' of residue ',seq(nrs),nrs,' in molecule #',nml
|
---|
249 | stop
|
---|
250 | endif
|
---|
251 |
|
---|
252 | endif ! N-terminus
|
---|
253 | ! _________________________________ previous atoms
|
---|
254 | 2 if (nfirp.eq.nfi.and.nrs.eq.ifirs) goto 11
|
---|
255 | nxtbb1=iowat(nfirp)
|
---|
256 | if (nxtbb1.eq.nfi.and.nrs.eq.ifirs) goto 11
|
---|
257 | nxtbb2=iowat(nxtbb1)
|
---|
258 | ! _______________________________ get data for substituent atoms
|
---|
259 | 3 if (iopfil(lunlib,reslib,'old','formatted').le.izero) then
|
---|
260 | write (*,'(a,/,a,i3,2a)')
|
---|
261 | & ' rplgrp> ERROR opening library of residues:',
|
---|
262 | & ' LUN=',lunlib,' FILE=',reslib(1:iendst(reslib))
|
---|
263 | stop
|
---|
264 | endif
|
---|
265 | call redres(sbrs,natsb,nxtsb,nvrsb)
|
---|
266 | close (lunlib)
|
---|
267 | ! __________________________ indices of substituent atoms
|
---|
268 | do i=1,natsb
|
---|
269 | if (rpat.eq.nmath(i)) then
|
---|
270 | nfisb=i
|
---|
271 | goto 4
|
---|
272 | endif
|
---|
273 | enddo
|
---|
274 | write (*,'(4a)') ' rplgrp> Cannot find atom >',rpat,
|
---|
275 | &'< in substituent residue ',sbrs
|
---|
276 | stop
|
---|
277 |
|
---|
278 | 4 nlasb=nfisb
|
---|
279 | do i=1,nbdath(nfisb)
|
---|
280 | ib=ibdath(i,nfisb)
|
---|
281 | if (iowath(ib).lt.nfisb) goto 10
|
---|
282 | do j=ib,natsb
|
---|
283 | if (j.gt.ib.and.iowath(j).lt.ib) goto 5
|
---|
284 | do k=1,nbdath(j)
|
---|
285 | if (ibdath(k,j).lt.nfisb) goto 10
|
---|
286 | enddo
|
---|
287 | nlasb=j
|
---|
288 | enddo ! ... branch atoms
|
---|
289 | 5 enddo ! ... branches
|
---|
290 | ! _________________________________________________ local axes at 'nfirp'
|
---|
291 | call setsys(nxtbb1,nfirp,nxtbb2,x1,x2,x3,y1,y2,y3,z1,z2,z3)
|
---|
292 |
|
---|
293 | xtoat(nfirp)=x1
|
---|
294 | ytoat(nfirp)=x2
|
---|
295 | ztoat(nfirp)=x3
|
---|
296 | xbaat(nfirp)=z1
|
---|
297 | ybaat(nfirp)=z2
|
---|
298 | zbaat(nfirp)=z3
|
---|
299 |
|
---|
300 | ! _____________________ add virtual atoms
|
---|
301 | if (ntbb) then
|
---|
302 |
|
---|
303 | ct=cstoat(nxtbb2) ! t.angle_(+2)
|
---|
304 | st=sntoat(nxtbb2)
|
---|
305 | ca=csbaat(nxtbb1) ! b.angle_(+1)
|
---|
306 | sa=snbaat(nxtbb1)
|
---|
307 |
|
---|
308 | ! ------------------- Eyring
|
---|
309 | h2=-sa*ct
|
---|
310 | h3=-sa*st
|
---|
311 | x1=-ca*x1+h2*y1+h3*z1
|
---|
312 | x2=-ca*x2+h2*y2+h3*z2
|
---|
313 | x3=-ca*x3+h2*y3+h3*z3
|
---|
314 | dx=one/sqrt(x1*x1+x2*x2+x3*x3)
|
---|
315 | x1=x1*dx
|
---|
316 | x2=x2*dx
|
---|
317 | x3=x3*dx
|
---|
318 |
|
---|
319 | xat(izero)=xat(nfirp)+x1
|
---|
320 | yat(izero)=yat(nfirp)+x2
|
---|
321 | zat(izero)=zat(nfirp)+x3
|
---|
322 | z1=-st*y1+ct*z1
|
---|
323 | z2=-st*y2+ct*z2
|
---|
324 | z3=-st*y3+ct*z3
|
---|
325 | dz=one/sqrt(z1*z1+z2*z2+z3*z3)
|
---|
326 | z1=z1*dz
|
---|
327 | z2=z2*dz
|
---|
328 | z3=z3*dz
|
---|
329 |
|
---|
330 | ct=cstoat(nxtbb1) ! t.angle_(+1)
|
---|
331 | st=sntoat(nxtbb1)
|
---|
332 |
|
---|
333 | ! -------------------- Eyring with b.angle = 90 deg.
|
---|
334 | xat(-ione)=xat(izero)-ct*(z2*x3-z3*x2)-st*z1
|
---|
335 | yat(-ione)=yat(izero)-ct*(z3*x1-z1*x3)-st*z2
|
---|
336 | zat(-ione)=zat(izero)-ct*(z1*x2-z2*x1)-st*z3
|
---|
337 |
|
---|
338 | endif
|
---|
339 | ! _____________________________________________ Shift atom data
|
---|
340 | nrp=nlarp-nfirp
|
---|
341 | nsb=nlasb-nfisb
|
---|
342 | if (nrp.ne.nsb) then
|
---|
343 | nsh=nsb-nrp
|
---|
344 | do i=nfi,nfirp-1 ! bonds to atoms after repl. group
|
---|
345 | do j=1,mxbd
|
---|
346 | ib=ibdat(j,i)
|
---|
347 | if (ib.gt.nlarp.and.ib.le.nla) ibdat(j,i)=ib+nsh
|
---|
348 | enddo
|
---|
349 | enddo
|
---|
350 | ilaat=iatrs2(irsml2(ntlml))
|
---|
351 | if (nrp.gt.nsb) then ! less atoms
|
---|
352 | i1=nlarp+1
|
---|
353 | i2=ilaat
|
---|
354 | i3=1
|
---|
355 | else ! more atoms
|
---|
356 | if ((ilaat+nsh).gt.mxat) then
|
---|
357 | write (*,'(a,i5)') ' rplgrp> number of atoms > ',mxat
|
---|
358 | stop
|
---|
359 | endif
|
---|
360 |
|
---|
361 | i1=ilaat
|
---|
362 | i2=nlarp+1
|
---|
363 | i3=-1
|
---|
364 | endif
|
---|
365 |
|
---|
366 | do i=i1,i2,i3
|
---|
367 | ii=i+nsh
|
---|
368 | nbdat(ii)=nbdat(i)
|
---|
369 |
|
---|
370 | ibd(1)=iowat(i)
|
---|
371 | iybd(1)=iyowat(i)
|
---|
372 | do j = 1,mxbd
|
---|
373 | ibd(j+1)=ibdat(j,i)
|
---|
374 | iybd(j+1)=iybdat(j,i)
|
---|
375 | enddo
|
---|
376 |
|
---|
377 | do j=1,mxbd+1
|
---|
378 | if (ibd(j).gt.nfirp) ibd(j)=ibd(j)+nsh
|
---|
379 | enddo
|
---|
380 |
|
---|
381 | iowat(ii)=ibd(1)
|
---|
382 | iyowat(ii)=iybd(1)
|
---|
383 | do j = 1,mxbd
|
---|
384 | ibdat(j,ii)=ibd(j+1)
|
---|
385 | iybdat(j,ii)=iybd(j+1)
|
---|
386 | enddo
|
---|
387 |
|
---|
388 | ityat(ii)=ityat(i)
|
---|
389 | nmat(ii)=nmat(i)
|
---|
390 | cgat(ii)=cgat(i)
|
---|
391 | blat(ii)=blat(i)
|
---|
392 | baat(ii)=baat(i)
|
---|
393 | snbaat(ii)=snbaat(i)
|
---|
394 | csbaat(ii)=csbaat(i)
|
---|
395 | toat(ii)=toat(i)
|
---|
396 | sntoat(ii)=sntoat(i)
|
---|
397 | cstoat(ii)=cstoat(i)
|
---|
398 | xat(ii)=xat(i)
|
---|
399 | yat(ii)=yat(i)
|
---|
400 | zat(ii)=zat(i)
|
---|
401 | xtoat(ii)=xtoat(i)
|
---|
402 | ytoat(ii)=ytoat(i)
|
---|
403 | ztoat(ii)=ztoat(i)
|
---|
404 | xbaat(ii)=xbaat(i)
|
---|
405 | ybaat(ii)=ybaat(i)
|
---|
406 | zbaat(ii)=zbaat(i)
|
---|
407 |
|
---|
408 | enddo
|
---|
409 | ! ____________________________________________ Shift residue data
|
---|
410 | do i=nrs+1,irsml2(ntlml)
|
---|
411 | iatrs1(i)=iatrs1(i)+nsh
|
---|
412 | iatrs2(i)=iatrs2(i)+nsh
|
---|
413 | ixatrs(i)=ixatrs(i)+nsh
|
---|
414 | enddo
|
---|
415 | iatrs2(nrs)=nla+nsh
|
---|
416 | if (nxt.gt.nlarp) ixatrs(nrs)=nxt+nsh
|
---|
417 | else
|
---|
418 | nsh=0
|
---|
419 | endif
|
---|
420 | ! _________________________________________ Correct data of 'nfirp'
|
---|
421 | ish=nfirp-nfisb
|
---|
422 | ityat(nfirp)=ityath(nfisb)
|
---|
423 |
|
---|
424 | if(.not.sh2) cgat(nfirp)=cgath(nfisb) ! NOT for 'ECEPP/2'
|
---|
425 |
|
---|
426 | nb=nbdath(nfisb) ! _______________ Bonds
|
---|
427 |
|
---|
428 | do i=1,mxbd
|
---|
429 | ib=ibdath(i,nfisb)
|
---|
430 | if (ib.ne.0) then
|
---|
431 | ibd(i)=ib+ish
|
---|
432 | iybd(i)=iybdath(i,nfisb)
|
---|
433 | else
|
---|
434 | ibd(i)=0
|
---|
435 | iybd(i)=1
|
---|
436 | endif
|
---|
437 |
|
---|
438 | enddo
|
---|
439 |
|
---|
440 | if (ntbb) then
|
---|
441 |
|
---|
442 | ibd(mxbd+1)=0
|
---|
443 | iybd(mxbd+1)=1
|
---|
444 |
|
---|
445 | ibd(nb+1)=nxtbb1+nsh ! bond to next backbone atom
|
---|
446 | iybd(nb+1)=isgbb1
|
---|
447 |
|
---|
448 | if (ibdrg.ne.0) then
|
---|
449 | nb=nb+1
|
---|
450 | if (nb.gt.mxbd) then
|
---|
451 | write (*,'(6a,/,2a,3(i4,a))')
|
---|
452 | & ' rplgrp> Cannot add atoms following ',rpat,
|
---|
453 | & ' from group ',sbrs,' to atom ',rpat,
|
---|
454 | & ' of residue ',seq(nrs),nrs,' in molecule #',nml,
|
---|
455 | & ' because need >',(mxbd+1),' bonds'
|
---|
456 | stop
|
---|
457 | endif
|
---|
458 | ibd(nb+1)=ibdrg+nsh
|
---|
459 | iybd(nb+1)=iybdrg
|
---|
460 | endif
|
---|
461 | iowat(nfirp)=ibd(1)
|
---|
462 | iyowat(nfirp)=iybd(1)
|
---|
463 | do i=1,mxbd
|
---|
464 | ibdat(i,nfirp)=ibd(i+1)
|
---|
465 | iybdat(i,nfirp)=iybd(i+1)
|
---|
466 | enddo
|
---|
467 | else ! not 'ntbb'
|
---|
468 |
|
---|
469 | do i=1,mxbd
|
---|
470 | ibdat(i,nfirp)=ibd(i)
|
---|
471 | iybdat(i,nfirp)=iybd(i)
|
---|
472 | enddo
|
---|
473 |
|
---|
474 | endif
|
---|
475 | nbdat(nfirp)=nb
|
---|
476 | ! _________________________________________ Add data for substituent
|
---|
477 | ii=nfirp
|
---|
478 | do i=nfisb+1,nlasb
|
---|
479 | ii=ii+1
|
---|
480 | nbdat(ii)=nbdath(i)
|
---|
481 |
|
---|
482 | iow=iowath(i)+ish
|
---|
483 | iowat(ii)=iow
|
---|
484 | iyowat(ii)=iyowath(i)
|
---|
485 |
|
---|
486 | do j=1,mxbd
|
---|
487 | ib=ibdath(j,i)
|
---|
488 | if (ib.ge.nfisb) then
|
---|
489 | ibdat(j,ii)=ib+ish
|
---|
490 | else
|
---|
491 | ibdat(j,ii)=ib
|
---|
492 | endif
|
---|
493 | iybdat(j,ii)=iybdath(j,i)
|
---|
494 | enddo
|
---|
495 |
|
---|
496 | ityat(ii)=ityath(i)
|
---|
497 | nmat(ii)=nmath(i)
|
---|
498 | cgat(ii)=cgath(i)
|
---|
499 | blat(ii)=blath(i)
|
---|
500 | ba=baath(i)
|
---|
501 | baat(ii)=ba
|
---|
502 | csbaat(ii)=cos(ba)
|
---|
503 | snbaat(ii)=sin(ba)
|
---|
504 | to=toath(i)
|
---|
505 | toat(ii)=to
|
---|
506 | cstoat(ii)=cos(to)
|
---|
507 | sntoat(ii)=sin(to)
|
---|
508 | call eyring(ii,iow)
|
---|
509 | if (ntbb) then ! reset some internal coordinates
|
---|
510 | if (iow.eq.nfirp) then
|
---|
511 | ba=valang(izero,nfirp,ii)
|
---|
512 | baat(ii)=ba
|
---|
513 | csbaat(ii)=cos(ba)
|
---|
514 | snbaat(ii)=sin(ba)
|
---|
515 | to=dihedr(-ione,izero,nfirp,ii)
|
---|
516 | toat(ii)=to
|
---|
517 | cstoat(ii)=cos(to)
|
---|
518 | sntoat(ii)=sin(to)
|
---|
519 | elseif (iowat(iow).eq.nfirp) then
|
---|
520 | to=dihedr(izero,nfirp,iow,ii)
|
---|
521 | toat(ii)=to
|
---|
522 | cstoat(ii)=cos(to)
|
---|
523 | sntoat(ii)=sin(to)
|
---|
524 | endif
|
---|
525 |
|
---|
526 | endif ! ntbb
|
---|
527 |
|
---|
528 | enddo ! substituent atoms
|
---|
529 | ! ___________________________________________________ Take care of Variables
|
---|
530 | ! (assume variables of replaced group/substituent to be stored CONSECUTIVELY)
|
---|
531 |
|
---|
532 | ilavr=ivrml1(ntlml)+nvrml(ntlml)-1
|
---|
533 | ifivr=ivrrs1(nrs) ! variables to be replaced (#=ivrrp,last=lvrrp)
|
---|
534 | ivrrp=0
|
---|
535 | do i=ifivr,ifivr+nvrrs(nrs)-1
|
---|
536 | iat=iatvr(i)
|
---|
537 | if (nfirp.lt.iat.and.iat.le.nlarp) then
|
---|
538 | ivrrp=ivrrp+1
|
---|
539 | lvrrp=i
|
---|
540 | endif
|
---|
541 | enddo
|
---|
542 | if (ivrrp.eq.0) then ! No variables to replace
|
---|
543 | do i=ifivr,ilavr
|
---|
544 | if (iatvr(i).gt.nlarp) then
|
---|
545 | lvrrp=i-1
|
---|
546 | goto 6
|
---|
547 | endif
|
---|
548 | enddo
|
---|
549 | lvrrp=ilavr
|
---|
550 | endif
|
---|
551 | 6 ivrsb=0 ! variables from substituent (#=ivrsb)
|
---|
552 | do i=1,nvrsb
|
---|
553 | iat=iatvrh(i)
|
---|
554 | if (nfisb.lt.iat.and.iat.le.nlasb) then
|
---|
555 | ity=ityvrh(i)
|
---|
556 | if (ntbb.and.iowath(iat).eq.nfisb.and.ity.gt.2) goto 7
|
---|
557 | ivrsb=ivrsb+1
|
---|
558 | nmvrh(ivrsb)=nmvrh(i)
|
---|
559 | ityvrh(ivrsb)=ityvrh(i)
|
---|
560 | iclvrh(ivrsb)=iclvrh(i)
|
---|
561 | iatvrh(ivrsb)=iat
|
---|
562 | endif
|
---|
563 | 7 enddo
|
---|
564 | if (nsh.ne.0) then ! if # of atoms changed
|
---|
565 | do i=lvrrp+1,ilavr
|
---|
566 | iatvr(i)=iatvr(i)+nsh
|
---|
567 | enddo
|
---|
568 | endif
|
---|
569 |
|
---|
570 | if (ivrrp.ne.ivrsb) then ! shift data for variables
|
---|
571 | jsh=ivrsb-ivrrp
|
---|
572 | if (ivrrp.gt.ivrsb) then
|
---|
573 | i1=lvrrp+1
|
---|
574 | i2=ilavr
|
---|
575 | i3=1
|
---|
576 | else
|
---|
577 | if ((ilavr+jsh).gt.mxvr) then
|
---|
578 | write (*,'(a,i5)') ' rplgrp> number of variables > ',mxvr
|
---|
579 | stop
|
---|
580 | endif
|
---|
581 | i1=ilavr
|
---|
582 | i2=lvrrp+1
|
---|
583 | i3=-1
|
---|
584 | endif
|
---|
585 | do i=i1,i2,i3
|
---|
586 | ii=i+jsh
|
---|
587 | ityvr(ii)=ityvr(i)
|
---|
588 | iclvr(ii)=iclvr(i)
|
---|
589 | nmvr(ii)=nmvr(i)
|
---|
590 | iatvr(ii)=iatvr(i)
|
---|
591 | enddo
|
---|
592 |
|
---|
593 | do i=nrs+1,irsml2(ntlml)
|
---|
594 | ivrrs1(i)=ivrrs1(i)+jsh
|
---|
595 | enddo
|
---|
596 | nvrrs(nrs)=nvrrs(nrs)+jsh
|
---|
597 | nvrml(nml)=nvrml(nml)+jsh
|
---|
598 | do i=nml+1,ntlml
|
---|
599 | ivrml1(i)=ivrml1(i)+jsh
|
---|
600 | enddo
|
---|
601 | endif
|
---|
602 | ii=lvrrp-ivrrp ! Add variables for substitutent
|
---|
603 | do i=1,ivrsb
|
---|
604 | ii=ii+1
|
---|
605 | nmvr(ii)=nmvrh(i)
|
---|
606 | ityvr(ii)=ityvrh(i)
|
---|
607 | iclvr(ii)=iclvrh(i)
|
---|
608 | iatvr(ii)=iatvrh(i)+ish
|
---|
609 | enddo
|
---|
610 |
|
---|
611 | return
|
---|
612 | ! __________________________________________ Errors
|
---|
613 | 10 write (*,'(3a,/,2a,i4,a,i4,/,2a)')
|
---|
614 | & ' rplgrp> Cannot replace atom(s) following ',rpat,
|
---|
615 | & ' from INSIDE a ring',' in residue: ',seq(nrs),nrs,
|
---|
616 | & ' in molecule #',nml,' or in substitute: ',sbrs
|
---|
617 | stop
|
---|
618 | 11 write (*,'(4a,i4,a,i4,/,a)')
|
---|
619 | & ' rplgrp> Cannot replace atom(s) following ',rpat,
|
---|
620 | & ' of residue ',seq(nrs),nrs,' in molecule #',nml,
|
---|
621 | & ' since necessary 2 previous atoms are not available'
|
---|
622 | stop
|
---|
623 |
|
---|
624 | end
|
---|
625 | ! ****************************************
|
---|
626 | subroutine redchg(nml,nrs,rpat,sbrs)
|
---|
627 |
|
---|
628 | ! .........................................................
|
---|
629 | ! PURPOSE: read and place atomic point charges from residue
|
---|
630 | ! 'sbrs' to residue 'nrs' of molecule 'nml'
|
---|
631 | ! from library 'chglib' with LUN=lunchg, if ilib=1
|
---|
632 | ! 'reslib' with LUN=lunlib, if ilib=2
|
---|
633 | !
|
---|
634 | ! CALLS: iopfil,iendst,tolost
|
---|
635 | ! ........................................................
|
---|
636 |
|
---|
637 | include 'INCL.H'
|
---|
638 |
|
---|
639 | character rpat*4,sbrs*4,atnm*4,res*4,cgty*5,line*100
|
---|
640 |
|
---|
641 | ifirs=irsml1(nml)
|
---|
642 | ilars=irsml2(nml)
|
---|
643 | res=seq(nrs)
|
---|
644 | call tolost(res)
|
---|
645 |
|
---|
646 | if (nrs.eq.ifirs.or.nrs.eq.ilars) then
|
---|
647 |
|
---|
648 | if (rpat.eq.'n '.or.rpat.eq.'c ') then
|
---|
649 | if (nrs.eq.ifirs.and.nrs.eq.ilars) goto 10 ! Dont have this yet
|
---|
650 |
|
---|
651 | if (rpat.eq.'n '.and.nrs.eq.ifirs) then
|
---|
652 | if (res(1:3).eq.'pro'.or.res(1:3).eq.'hyp') then
|
---|
653 | if (sbrs.eq.'nh1 ') cgty='n'//res
|
---|
654 | if (sbrs.eq.'nh2+') cgty='+'//res
|
---|
655 | else
|
---|
656 | if (sbrs.eq.'nh2 ') cgty='n'//res
|
---|
657 | if (sbrs.eq.'nh3+') cgty='+'//res
|
---|
658 | endif
|
---|
659 | elseif (rpat.eq.'c '.and.nrs.eq.ilars) then
|
---|
660 | if (sbrs.eq.'cooh') cgty='c'//res
|
---|
661 | if (sbrs.eq.'coo-') cgty='-'//res
|
---|
662 | else
|
---|
663 | goto 10
|
---|
664 | endif
|
---|
665 | else
|
---|
666 | write (*,*) ' redchg> dont know which end goup is present'
|
---|
667 | stop
|
---|
668 | endif
|
---|
669 | ilib=1
|
---|
670 | else
|
---|
671 | cgty(1:)=sbrs
|
---|
672 | ilib=2
|
---|
673 | endif
|
---|
674 |
|
---|
675 | if (ilib.eq.1) then
|
---|
676 |
|
---|
677 | if (iopfil(lunchg,chgfil,'old','formatted').le.izero) then
|
---|
678 | write (*,'(a,/,a,i3,2a)')
|
---|
679 | & ' redchg> ERROR opening library of charges:',
|
---|
680 | & ' LUN=',lunchg,' FILE=',chgfil(1:iendst(chgfil))
|
---|
681 | stop
|
---|
682 | endif
|
---|
683 |
|
---|
684 | 1 line=' '
|
---|
685 | read (lunchg,'(a)',end=3) line
|
---|
686 | l=iendst(line)
|
---|
687 |
|
---|
688 | if (l.ge.10.and.line(1:1).eq.'#'.and.line(2:6).eq.cgty) then
|
---|
689 | read (line(7:10),'(i4)') nchg
|
---|
690 | if (nchg.le.mxath) then
|
---|
691 | read (lunchg,'(10(2x,a4,1x))') (nmath(i),i=1,nchg)
|
---|
692 | read (lunchg,'(10(f6.4,1x))') (cgath(i),i=1,nchg)
|
---|
693 | close(lunchg)
|
---|
694 | do i=iatrs1(nrs),iatrs2(nrs)
|
---|
695 | atnm=nmat(i)
|
---|
696 | do j=1,nchg
|
---|
697 | if (nmath(j).eq.atnm) then
|
---|
698 | cgat(i)=cgath(j)
|
---|
699 | goto 2
|
---|
700 | endif
|
---|
701 | enddo
|
---|
702 | write (*,'(6a)') ' redchg> Cannot find atom: ',atnm,
|
---|
703 | & ' for entry: ',cgty,' in library: ',
|
---|
704 | & chgfil(1:iendst(chgfil))
|
---|
705 | stop
|
---|
706 | 2 enddo
|
---|
707 | return
|
---|
708 | else
|
---|
709 | write (*,'(4a)')
|
---|
710 | & ' redchg> must increase MXATH to read data for entry: ',
|
---|
711 | & cgty,' in library: ',chgfil(1:iendst(chgfil))
|
---|
712 | close(lunchg)
|
---|
713 | stop
|
---|
714 | endif
|
---|
715 | endif
|
---|
716 | goto 1
|
---|
717 | 3 write (*,'(4a)')
|
---|
718 | & ' redchg> Cannot find entry: ',cgty,' in library: ',
|
---|
719 | & chgfil(1:iendst(chgfil))
|
---|
720 | close(lunchg)
|
---|
721 | stop
|
---|
722 |
|
---|
723 | elseif (ilib.eq.2) then
|
---|
724 |
|
---|
725 | if (iopfil(lunlib,reslib,'old','formatted').le.izero) then
|
---|
726 | write (*,'(a,/,a,i3,2a)')
|
---|
727 | & ' redchg> ERROR opening library of residues:',
|
---|
728 | & ' LUN=',lunlib,' FILE=',reslib(1:iendst(reslib))
|
---|
729 | stop
|
---|
730 | endif
|
---|
731 |
|
---|
732 | 4 line=' '
|
---|
733 | read (lunlib,'(a)',end=6) line
|
---|
734 | l=iendst(line)
|
---|
735 |
|
---|
736 | if (l.ge.9.and.line(1:1).eq.'#'.and.line(2:5).eq.cgty(1:4)) then
|
---|
737 | read (line(6:9),'(i4)') nchg
|
---|
738 | if (nchg.le.mxath) then
|
---|
739 | read (lunlib,'(a4,42x,d7.4)') (nmath(i),cgath(i),i=1,nchg)
|
---|
740 | close(lunlib)
|
---|
741 | do i=iatrs1(nrs),iatrs2(nrs)
|
---|
742 | atnm=nmat(i)
|
---|
743 | do j=1,nchg
|
---|
744 | if (nmath(j).eq.atnm) then
|
---|
745 | cgat(i)=cgath(j)
|
---|
746 | goto 5
|
---|
747 | endif
|
---|
748 | enddo
|
---|
749 | write (*,'(6a)') ' redchg> Cannot find atom: ',atnm,
|
---|
750 | & ' for entry: ',cgty,' in library: ',
|
---|
751 | & reslib(1:iendst(reslib))
|
---|
752 | stop
|
---|
753 | 5 enddo
|
---|
754 | return
|
---|
755 | else
|
---|
756 | write (*,'(4a)')
|
---|
757 | & ' redchg> must increase MXATH to read data for entry: ',
|
---|
758 | & cgty,' in library: ',reslib(1:iendst(reslib))
|
---|
759 | close(lunchg)
|
---|
760 | stop
|
---|
761 | endif
|
---|
762 | endif
|
---|
763 | goto 4
|
---|
764 | 6 write (*,'(4a)')
|
---|
765 | & ' redchg> Cannot find entry: ',cgty,' in library: ',
|
---|
766 | & reslib(1:iendst(reslib))
|
---|
767 | close(lunchg)
|
---|
768 | stop
|
---|
769 |
|
---|
770 | endif
|
---|
771 |
|
---|
772 | 10 write (*,'(4a)')
|
---|
773 | & ' redchg> Do not have charges for N/C-terminal residue ',
|
---|
774 | & res,' modified with group :',sbrs
|
---|
775 | stop
|
---|
776 |
|
---|
777 | end
|
---|
778 |
|
---|