source: getmol.f@ cb47b9c

Last change on this file since cb47b9c was cb47b9c, checked in by baerbaer <baerbaer@…>, 14 years ago

Explicitly declare variables.

All variables should be declared so that we can remove the implicit statements
from the beginning of the INCL.H file.

git-svn-id: svn+ssh://svn.berlios.de/svnroot/repos/smmp/trunk@25 26dc1dd8-5c4e-0410-9ffe-d298b4865968

  • Property mode set to 100644
File size: 13.2 KB
Line 
1! **************************
2! **************************************************************
3!
4! This file contains the subroutines: getmol,redres
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! **************************************************************
12
13
14 subroutine getmol(nml)
15
16! ...................................................................
17! PURPOSE: assemble data for molecule 'nml' according to
18! its sequence using residue library 'reslib'
19!
20! ! Molecules must be assembled in sequential order (1 -> ntlml)
21! (or number of atoms & variables must remain the same)
22!
23! INPUT: irsml1(nml),irsml2(nml),seq(irsml1()...irsml2())
24! nml>1: irsml1(nml-1),iatrs2(irsml2(nml-1))
25! ivrrs1(irsml2(nml-1)),nvrrs(irsml2(nml-1))
26!
27! OUTPUT: molecule - ivrml1,nvrml
28! residues - iatrs1,ixatrs,iatrs2,ivrrs1,nvrrs
29! atoms - nmat,ityat,cgat,blat,baat,csbaat,snbaat,
30! toat,cstoat,sntoat
31! bonds - nbdat,iowat,iyowat,ibdat(1-mxbd,),iybdat(1-mxbd,)
32! ! 1st atom of 'nml': iowat indicates 1st bond
33! to a FOLLOWING atom (not previous) !
34! variables - ityvr,iclvr,iatvr,nmvr
35!
36! CALLS: iopfil,redres,iendst
37! ...................................................................
38
39 include 'INCL.H'
40! arguments
41 integer nml
42! functions
43 integer iopfil, iendst
44
45 integer i, iat, iow, ifirs, ibd, ilars, j, mvr, nbd, n, nbd1, nj
46 integer nh, nrs, ntlat, ntlvr, nat, nxt
47 double precision ba, t, to
48 character res*4
49
50 if (iopfil(lunlib,reslib,'old','formatted').le.izero) then
51 write (*,'(a,/,a,i3,2a)')
52 & ' getmol> ERROR opening library of residues:',
53 & ' LUN=',lunlib,' FILE=',reslib(1:iendst(reslib))
54 stop
55 endif
56
57 if (nml.eq.1) then
58 ntlat=0
59 ntlvr=0
60 else
61 i=irsml2(nml-1) ! last res. of previous mol.
62 ntlat=iatrs2(i)
63 ntlvr=ivrrs1(i)+nvrrs(i)-1
64 endif
65 ivrml1(nml)=ntlvr + 1
66
67 ilars=irsml2(nml)
68 ifirs=irsml1(nml)
69
70 do nrs=ifirs,ilars ! Residues in molecule
71
72 res = seq(nrs)
73 call tolost(res) ! ensure lower case for residue name
74
75 if (res(:3).eq.'nme'.and.nrs.ne.ilars) then
76 write (*,'(3a)') ' getmol> residue >',res,
77 & '< allowed at C-terminus only !'
78 close(lunlib)
79 stop
80 elseif (res(:3).eq.'ace'.and.nrs.ne.ifirs) then
81 write (*,'(3a)') ' getmol> residue >',res,
82 & '< allowed at N-terminus only !'
83 close(lunlib)
84 stop
85 endif
86
87 call redres(res,nat,nxt,nvr)
88
89 if ((nat+ntlat).gt.mxat) then
90 write (*,'(a,i5)') ' getmol> number of atoms > ',mxat
91 close(lunlib)
92 stop
93 endif
94 if ((nvr+ntlvr).gt.mxvr) then
95 write (*,'(a,i5)') ' getmol> number of variables > ',mxvr
96 close(lunlib)
97 stop
98 endif
99
100 rewind lunlib
101
102! ___________________________________________________________ Atoms
103 do i=1,nat
104 n=i+ntlat
105 nmat(n)=nmath(i)
106 ityat(n)=ityath(i)
107 cgat(n)=cgath(i)
108 blat(n)=blath(i)
109 ba=baath(i)
110 baat(n)=ba
111 csbaat(n)=cos(ba)
112 snbaat(n)=sin(ba)
113 to=toath(i)
114 toat(n)=to
115 cstoat(n)=cos(to)
116 sntoat(n)=sin(to)
117! ______________________________ bonds to previous & following atoms
118 iow=iowath(i)
119 if (iow.eq.0) then ! 1st atom of residue
120 if (nrs.eq.ifirs) then ! 1st atom of 'nml'
121
122 iowat(n)=ibdath(1,i)+ntlat
123 iyowat(n)=iybdath(1,i)
124
125 do j = 1,mxbd-1
126 ibdath(j,i)=ibdath(j+1,i)
127 iybdath(j,i)=iybdath(j+1,i)
128 enddo
129 ibdath(mxbd,i)=0
130 iybdath(mxbd,i)=1
131
132 nbdath(i)=nbdath(i)-1
133 else ! connected with prev. res.
134 nh=ixatrs(nrs-1) ! atom to 'next' residue of prev.res.
135 iowat(n)=nh
136 iyowat(n)=1 !!! only single bonds assumed !!!
137
138! ___________________________ correct atom to 'next' res.
139 nbd=nbdat(nh)
140 if (nbd.eq.mxbd) then
141 write(*,'(a,i2,a,i4,2a,i4,a)')
142 & ' getmol> need ',(mxbd+2),
143 & 'th bond to connect residues ',
144 & nrs-1,seq(nrs-1),' and ',nrs,seq(nrs)
145 close(lunlib)
146 stop
147 else ! correct atom to 'next' res.
148! _______________________________!! dihedrals for atoms bound to 'nh'
149! are assumed to be phase angles !!
150 do j=1,nbd
151
152 nj=ibdat(j,nh)
153 t=toat(nj)
154
155 if (t.eq.0.0) then
156 write (*,'(3a,/,2a)')
157 & ' getmol> DIHEDRAL for atom ',nmat(nj),
158 & ' should be PHASE angle with respect to atom ',
159 & nmat(n),' & therefore must be not 0.0 !!'
160 close(lunlib)
161 stop
162 endif
163
164 t=t+to
165 if (abs(t).gt.pi) t=t-sign(pi2,t)
166 toat(nj)=t
167 cstoat(nj)=cos(t)
168 sntoat(nj)=sin(t)
169
170 enddo
171
172 nbd1=nbd+1
173 ibdat(nbd1,nh)=n
174 iybdat(nbd1,nh)=1 ! (only single bonds !)
175 nbdat(nh)=nbd1
176
177 endif
178
179 endif
180
181 else ! connected within res.
182
183 iowat(n)=ntlat+iow
184 iyowat(n)=iyowath(i)
185
186 endif
187
188 nbdat(n)=nbdath(i)
189
190 do j=1,mxbd
191
192 ibd=ibdath(j,i)
193
194 if (ibd.ne.0) then
195 ibdat(j,n)=ibd+ntlat
196 iybdat(j,n)=iybdath(j,i)
197 else
198 ibdat(j,n)=0
199 iybdat(j,n)=1
200 endif
201
202 enddo
203
204 enddo ! ... atoms
205
206! ________________________________________________________ Variables
207 ivrrs1(nrs)=ntlvr+1
208 mvr=0
209
210 do i=1,nvr
211
212 if (nrs.eq.ifirs) then
213
214 iat=iatvrh(i)
215! ____________________________________ Exclude all variables for 1st atom
216! & torsion for atoms bound to it
217 if ( iat.eq.1.or.
218 & (iowath(iat).eq.1.and.ityvrh(i).eq.3)) goto 1
219
220 endif
221
222 mvr=mvr+1
223 ntlvr=ntlvr+1
224 ityvr(ntlvr)=ityvrh(i)
225 iclvr(ntlvr)=iclvrh(i)
226 iatvr(ntlvr)=iatvrh(i)+ntlat
227 nmvr(ntlvr)=nmvrh(i)
228
229 1 enddo ! ... Variables
230
231 nvrrs(nrs)=mvr
232
233 iatrs1(nrs)=ntlat+1 ! first backbone atom of res.
234 ixatrs(nrs)=ntlat+nxt ! last backbone atom
235 ntlat=ntlat+nat
236 iatrs2(nrs)=ntlat ! last atom of res.
237
238 enddo ! ... residues
239
240 close(lunlib)
241
242! _______________________________ Variables
243 if (nml.eq.1) then
244 nvrml(nml)=ntlvr
245 else
246 nvrml(nml)=ntlvr-ivrml1(nml) + 1
247 endif
248
249 return
250 end
251! **************************************
252 subroutine redres(res,nat,nxt,nvrr)
253
254! .......................................................
255! PURPOSE: read atom data for residue 'res' from library
256! (file 'lunlib' 'reslib' opened in routine calling
257! this one)
258!
259! OUTPUT: nat - number of atoms in residue
260! nxt - atom which may bind to following residue
261! nvrr - number of variables in residue
262! for atoms - nmath,blath,baath(rad),toath(rad),
263! ityath,iyowath,iowath (INSIDE residue,
264! =0 if 1st atom)
265! for variables - ityvrh (1=bl/2=ba/3=to),iclvrh,iatvrh,nmvrh
266!
267! LIBRARY: residue-lines:
268! '#', res, nat, nxt; Format: a1,a4,2i4
269! atom-lines:
270! nmat,3{"fix" =' ', clvr,nmvr, blat/baat(deg)/toat(deg)},
271! cgat, ityat, iowat,ibdat1,ibdat2,ibdat3;
272! Format: a4, 3(1x,i2,a1,a3,f9.3), f7.4, i4,4i4
273!
274! CALLS: iendst,tolost
275!
276! .......................................................
277
278 include 'INCL.H'
279! arguments
280 integer nat, nxt, nvrr
281
282! functions
283 integer iendst
284
285 integer icl, ibd, ib3, i, ib1, ib2, ic, iow, iexcp, ity, jow, lg
286 integer nln, j
287 dimension icl(3),ibd(mxbd)
288 character blnk,fix(3),nm(3)*3,res*4,resl*4,line*132
289 data blnk/' '/
290
291 double precision ba, to
292
293 nln=0
294 do i=1,mxbd
295 ibd(i)=0
296 enddo
297
298 resl=res
299
300 call tolost(resl) ! ensure lower case for residue name
301
302! ________________________________ find residue 'resl'
303 1 line=blnk
304 nln=nln+1
305
306 read (lunlib,'(a)',end=2,err=3) line
307 lg=iendst(line)
308
309 if (lg.ge.13.and.line(1:1).eq.'#'.and.line(2:5).eq.resl) then
310! _____________________________________________ read atom data for 'resl'
311 read (line(6:13),'(2i4)',err=3) nat,nxt
312
313 if (nat.gt.mxath) then
314 write (*,'(a,i5)') ' redres> number of atoms > ',mxath
315 close(lunlib)
316 stop
317 endif
318
319 nvrr=0
320 do i=1,nat
321
322 nln=nln+1
323
324 read (lunlib,'(a4,3(1x,i2,a1,a3,d9.3),d7.4,i4,4i4)',
325 & end=3,err=3)
326 & nmath(i),icl(1),fix(1),nm(1),blath(i),icl(2),fix(2),nm(2),ba,
327 & icl(3),fix(3),nm(3),to,cgath(i),ity,iow,(ibd(j),j=1,mxbd)
328
329 if (ity.le.0.or.ity.gt.mxtyat) goto 6
330 ityath(i)=ity
331
332 jow=abs(iow)
333
334 if (res(:3).eq.'ace'.and.i.eq.1) then ! exception from following check
335 iexcp = 1
336 else
337 iexcp = 0
338 endif
339
340 if (iexcp.eq.0.and.i.le.jow) then
341 if (i.eq.jow) then
342 write (*,'(5a)') ' redres> atom ',nmath(i),' of ',
343 & resl,' cannot preceed itself '
344 else
345 write (*,'(5a,i4)') ' redres> atom ',nmath(i),' of ',
346 & resl,' should be placed AFTER atom #',jow
347 endif
348 goto 5
349 endif
350
351 iowath(i)=jow
352 iyowath(i)=sign(1,iow)
353! ____________________________________ check order & find number of bonds
354! (bonds closing ring must be last !)
355 ib1=abs(ibd(1))
356 ib2=abs(ibd(2))
357 ib3=abs(ibd(3))
358
359 if (ib1.eq.i.or.ib2.eq.i.or.ib3.eq.i) goto 4
360
361 if (ib1.eq.0) then ! no bond to following
362 if (ib2.ne.0.or.ib3.ne.0) goto 4
363 nbdath(i)=0
364 else
365 if (ib1.eq.jow) goto 4
366 if (ib2.eq.0) then
367 if (ib3.ne.0) goto 4
368 nbdath(i)=1
369 else
370 if ( ib2.eq.jow.or.ib2.eq.ib1.or.
371 & (ib2.gt.i.and.ib2.lt.ib1) ) goto 4
372 if (ib3.eq.0) then
373 nbdath(i)=2
374 else
375 if (ib3.eq.jow.or.ib3.eq.ib1.or.ib3.eq.ib2.or.
376 & (ib3.gt.i.and.(ib3.lt.ib1.or.ib3.lt.ib2)) ) goto 4
377 nbdath(i)=3
378 endif
379 endif
380 endif
381
382 do j=1,mxbd
383 ibdath(j,i)=abs(ibd(j))
384 iybdath(j,i)=sign(1,ibd(j))
385 enddo
386
387 baath(i)=ba*cdr ! convert angles into 'radians'
388 toath(i)=to*cdr
389
390! ______________________________ internal degrees of freedom
391 do j=1,3
392 if (fix(j).ne.blnk) then
393 nvrr=nvrr+1
394
395 if (nvrr.gt.mxvrh) then
396 write (*,'(a,i5)') ' redres> number of variables > ',
397 & mxvrh
398 close(lunlib)
399 stop
400 endif
401
402 ic=icl(j)
403
404 if ( ic.le.0
405 & .or.(j.eq.3.and.ic.gt.mxtyto) ! dihedral
406 & .or.(j.eq.2.and.ic.gt.mxtyba) ! bond angle
407 & .or.(j.eq.1.and.ic.gt.mxtybl) ) goto 7 ! b. length
408
409 ityvrh(nvrr)=j
410 iclvrh(nvrr)=ic
411 iatvrh(nvrr)=i
412 nmvrh(nvrr)=nm(j)
413
414 endif
415
416 enddo
417 enddo ! ... atoms
418
419 return
420 endif
421
422 goto 1
423
424! ____________________________________________________________ ERRORS
425 2 write (*,'(4a)') ' redres> residue >',resl,'< NOT FOUND in ',
426 &reslib(1:iendst(reslib))
427 close(lunlib)
428 stop
429
430 3 write (*,'(a,i4,2a)') ' redres> ERROR reading line No. ',nln,
431 &' in ',reslib(1:iendst(reslib))
432 close(lunlib)
433 stop
434
435 4 write (*,'(4a)') ' redres> Incorrect order of bonds for atom ',
436 & nmath(i),' of ',resl
437
438 5 write (*,'(8x,2a)') '... must correct ',
439 & reslib(1:iendst(reslib))
440 close(lunlib)
441 stop
442
443 6 write (*,'(a,i2,4a)') ' redres> unknown type :',ity,
444 & ': for atom ',nmath(i),' in residue ',resl
445 close(lunlib)
446 stop
447
448 7 write (*,'(a,i2,4a)') ' redres> unknown class :',ic,
449 & ': for variable ',nm(j),' in residue ',resl
450 close(lunlib)
451 stop
452
453 end
454
Note: See TracBrowser for help on using the repository browser.