source: getmol.f@ e40e335

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

Initial import to BerliOS corresponding to 3.0.4

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

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