1      logical function argos_prep_wrtpdb(lfnout,lfnpdb,lrgpdb,
2     + filpdb,iopt,box,num,amass,mat,
3     + csa,isat,isgm,imol,ifra,xs,vs,msa,nsa,
4     + cwa,iwat,xw,vw,mwm,mwa,nwm,nwa,
5     + xwc,vwc,mwmc,nwmc,slvnam,iropt,irrand,nxrep,nyrep,nzrep,drep,
6     + msb,nsb,idsb,rdist,nskip,iskip,lang,lfnmrg,nmerge,xmerge,filmrg,
7     + irenum,invert,ihop,ips)
8c
9c $Id$
10c
11      implicit none
12c
13#include "util.fh"
14#include "mafdecls.fh"
15c
16      character*2 argos_prep_atnam
17      real*8 argos_prep_atsiz
18      integer argos_prep_merge
19      external argos_prep_atnam,argos_prep_atsiz,argos_prep_merge
20c
21      integer lfnout,lfnpdb,iopt,msa,nsa,mwm,nwm,mwa,nwa,mat,mwmc,nwmc
22      integer lfnmrg,lrgpdb
23      character*255 filpdb,filmrg(100)
24      integer num(mat),isat(msa),isgm(msa),imol(msa),ifra(msa),iwat(mwa)
25      integer ihop(msa),ips(msa)
26      character*16 cwa(mwa),csa(msa)
27      character*3 slvnam
28      real*8 amass(mat),xs(3,msa),vs(3,msa),xw(3,mwa,mwm),vw(3,mwa,mwm)
29      real*8 xwc(3,mwa,mwmc),vwc(3,mwa,mwmc),box(3)
30      integer iropt,irrand,nxrep,nyrep,nzrep,nskip
31      real*8 rdist,drep(3),xmerge(3,100)
32      integer msb,nsb,nmols,irenum,invert,nmerge
33      integer idsb(2,msb),iskip(3,nskip),lang(*)
34c
35      integer length,ioff
36c
37      integer i,j,k,l,nats,nums,numslv,isoff,iaoff,natom,isegm,nsegm
38      integer ixrep,iyrep,izrep
39      real*8 temp,dbx,dby,dbz,dist,dbz0,angle
40      real*8 o(3),px(3),py(3),pz(3),x(3),y(3),v(3)
41      character*3 sname
42      character*1 charr
43      integer icharr,ilang,inum
44c
45      real*8 dumdst
46c
47      dumdst=0.025d0
48c
49      isoff=0
50      iaoff=0
51      inum=-1
52      if(irenum.ne.0) inum=0
53      charr=' '
54      icharr=ichar(charr)
55      dbz0=0.0d0
56c
57      o(1)=0.0d0
58      o(2)=0.0d0
59      o(3)=0.0d0
60      px(1)=1.0d0
61      px(2)=0.0d0
62      px(3)=0.0d0
63      py(1)=0.0d0
64      py(2)=1.0d0
65      py(3)=0.0d0
66      pz(1)=0.0d0
67      pz(2)=0.0d0
68      pz(3)=1.0d0
69c
70c     open PDB file
71c
72      length=index(filpdb,' ')-1
73      open(unit=lfnpdb,file=filpdb(1:length),form='formatted',
74     + status='unknown',err=9999)
75c
76      if(nxrep.lt.1) nxrep=1
77      if(nyrep.lt.1) nyrep=1
78      if(nzrep.lt.-2) nzrep=1
79      nums=0
80c
81      if(nxrep.le.1.and.nyrep.le.1.and.iabs(nzrep).le.1) then
82c
83      if(lrgpdb.eq.0) then
84      write(lfnpdb,1000) (1.0d1*box(i),i=1,3),90.0,90.0,90.0
85 1000 format(
86     + 'HEADER',/,
87     + 'TITLE     ',/,
88     + 'TITLE    2',/,
89     + 'TITLE    3',/,
90     + 'REMARK   4 XXXX COMPLIES WITH FORMAT V. 2.1, 25-OCT-1996',/,
91     + 'CRYST1',3f9.3,3f7.2)
92      else
93      write(lfnpdb,2000) (1.0d1*box(i),i=1,3),90.0,90.0,90.0
94 2000 format(
95     + 'HEADER',/,
96     + 'TITLE     ',/,
97     + 'TITLE    2',/,
98     + 'TITLE    3',/,
99     + 'REMARK   4 YYYY DOES NOT COMPLY WITH FORMAT V. 2.1,',
100     + ' 25-OCT-1996',/,
101     + 'CRYST1',3f9.3,3f7.2,/,
102     + 'LRGPDB')
103      endif
104c
105      ioff=0
106c      if(filmrg(1:1).ne.' ') then
107      if(nmerge.gt.0) then
108      do 6060 k=1,nmerge
109      ioff=ioff+argos_prep_merge(lfnmrg,xmerge(1,k),filmrg(k),
110     + lfnout,lfnpdb,slvnam,inum,lrgpdb)
111 6060 continue
112      ioff=inum
113c      ioff=argos_prep_merge(lfnmrg,filmrg,lfnout,lfnpdb,slvnam,inum)
114      if(ioff.lt.0) call md_abort('Error in argos_prep_merge',ioff)
115      if(inum.ge.0) ioff=ioff+1-isgm(1)
116      else
117      if(irenum.ne.0) ioff=-isgm(1)+1
118      endif
119c
120      nums=0
121      do 11 k=1,nsb
122      i=idsb(1,k)
123      j=idsb(2,k)
124c      write(*,'(5i5)') i,j,k,isgm(i),isgm(j)
125      if(isgm(i).ne.isgm(j)) then
126      dist=sqrt((xs(1,j)-xs(1,i))*(xs(1,j)-xs(1,i))+
127     + (xs(2,j)-xs(2,i))*(xs(2,j)-xs(2,i))+
128     + (xs(3,j)-xs(3,i))*(xs(3,j)-xs(3,i)))
129      if(dist.gt.argos_prep_atsiz(num(isat(j)))+
130     + argos_prep_atsiz(num(isat(i)))) then
131      write(lfnpdb,1101) csa(i)(11:14),csa(i)(1:3),' ',isgm(i)+ioff,
132     + csa(j)(11:14),csa(j)(1:3),isgm(j)+ioff
133 1101 format('LINK',t13,a4,1x,a3,1x,a1,i4,16x,a4,1x,a3,2x,i4)
134      endif
135      endif
136   11 continue
137      do 1 i=1,nsa
138      if(ips(i).ge.0.or.ihop(i).eq.0) then
139      x(1)=xs(1,i)
140      x(2)=xs(2,i)
141      x(3)=xs(3,i)
142      else
143      x(1)=xs(1,i)-xs(1,i+ihop(i))
144      x(2)=xs(2,i)-xs(2,i+ihop(i))
145      x(3)=xs(3,i)-xs(3,i+ihop(i))
146      dist=sqrt(x(1)*x(1)+x(2)*x(2)+x(3)*x(3))
147      dist=dumdst/dist
148      x(1)=xs(1,i+ihop(i))+dist*x(1)
149      x(2)=xs(2,i+ihop(i))+dist*x(2)
150      x(3)=xs(3,i+ihop(i))+dist*x(3)
151      endif
152      temp=4.009103873d1*amass(isat(i))*
153     + (vs(1,i)*vs(1,i)+vs(2,i)*vs(2,i)+vs(3,i)*vs(3,i))
154      if(temp.gt.999.0) temp=999.989
155      sname=csa(i)(1:3)
156      if(sname(1:1).eq.'_') sname(1:1)=' '
157      if(sname(2:2).eq.'_') sname(2:2)=' '
158      if(sname(3:3).eq.'_') sname(3:3)=' '
159      if(lrgpdb.eq.0) then
160      write(lfnpdb,1001) i,csa(i)(11:14),sname,' ',
161     + isgm(i)+ioff,(1.0d1*x(k),k=1,3),temp,
162     + argos_prep_atnam(num(isat(i)))
163      else
164      write(lfnpdb,2001) i,csa(i)(11:14),sname,
165     + isgm(i)+ioff,(1.0d1*x(k),k=1,3),temp,
166     + argos_prep_atnam(num(isat(i)))
167      endif
168 1001 format('ATOM',i7,1x,a4,1x,a3,1x,a1,i4,4x,3f8.3,6x,f6.2,4x,a2)
169 2001 format('ATOM',i7,1x,a4,1x,a3,i6,4x,3f8.3,6x,f6.2,4x,a2)
170      if(isgm(i).gt.nums) nums=isgm(i)
171    1 continue
172      else
173      dbz0=0.0d0
174      charr='A'
175      icharr=ichar(charr)
176      if(nzrep.lt.0) then
177      if(nsa.gt.0.and.nwm.eq.0) dbz0=xs(3,1)
178      if(nsa.eq.0.and.nwm.gt.0) dbz0=xw(3,1,1)
179      if(nsa.gt.0.and.nwm.gt.0) dbz=min(xs(3,1),xw(3,1,1))
180      if(invert.eq.0) then
181      do 31 i=1,nsa
182      dbz0=min(dbz0,xs(3,i))
183   31 continue
184      do 32 i=1,nwm
185      do 33 j=1,nwa
186      dbz0=min(dbz0,xw(3,j,i))
187   33 continue
188   32 continue
189      dbz0=dbz0-0.5d0*rdist
190      else
191      do 41 i=1,nsa
192      dbz0=max(dbz0,xs(3,i))
193   41 continue
194      do 42 i=1,nwm
195      do 43 j=1,nwa
196      dbz0=max(dbz0,xw(3,j,i))
197   43 continue
198   42 continue
199      dbz0=dbz0+0.5d0*rdist
200      endif
201      endif
202c
203      nmols=0
204      do 34 i=1,nsa
205      if(iropt.eq.3) then
206      nmols=max(nmols,ifra(i))
207      else
208      nmols=max(nmols,imol(i))
209      endif
210   34 continue
211c
212      if(lrgpdb.eq.0) then
213      write(lfnpdb,1000) 1.0d1*dble(nxrep)*drep(1),
214     + 1.0d1*dble(nyrep)*drep(2),1.0d1*dble(nzrep)*drep(3),
215     + 90.0,90.0,90.0
216      else
217      write(lfnpdb,2000) 1.0d1*dble(nxrep)*drep(1),
218     + 1.0d1*dble(nyrep)*drep(2),1.0d1*dble(nzrep)*drep(3),
219     + 90.0,90.0,90.0
220      endif
221c
222      ioff=0
223      inum=0
224      if(nmerge.gt.0) then
225      do 6061 k=1,nmerge
226      ioff=ioff+argos_prep_merge(lfnmrg,xmerge(1,k),filmrg(k),
227     + lfnout,lfnpdb,slvnam,inum,lrgpdb)
228 6061 continue
229      ioff=inum
230c      if(filmrg(1:1).ne.' ') then
231c      ioff=argos_prep_merge(lfnmrg,filmrg,lfnout,lfnpdb,slvnam,inum)
232      if(ioff.lt.0) call md_abort('Error 2 in argos_prep_merge',0)
233      endif
234c
235      if(iropt.le.1) then
236      ilang=0
237      do 4 izrep=1,iabs(nzrep)
238      if(nzrep.gt.0) then
239      dbz=(dble(izrep)-0.5d0*dble(nzrep+1))*drep(3)
240      else
241      dbz=dbz0
242      endif
243      do 5 iyrep=1,nyrep
244      dby=(dble(iyrep)-0.5d0*dble(nyrep+1))*drep(2)
245      do 6 ixrep=1,nxrep
246      dbx=(dble(ixrep)-0.5d0*dble(nxrep+1))*drep(1)
247c
248      do 66 k=1,nskip
249      if(ixrep.eq.iskip(1,k).and.iyrep.eq.iskip(2,k).and.
250     + izrep.eq.iskip(3,k)) goto 6
251   66 continue
252c
253      angle=8.0d0*atan(1.0d0)*util_random(0)
254      ilang=ilang+1
255      lang(ilang)=angle
256c
257      if(iropt.eq.0) then
258      icharr=ichar(' ')
259      else
260      isoff=0
261      endif
262c
263      do 111 k=1,nsb
264      i=idsb(1,k)
265      j=idsb(2,k)
266      if(isgm(i).ne.isgm(j)) then
267      dist=sqrt((xs(1,j)-xs(1,i))*(xs(1,j)-xs(1,i))+
268     + (xs(2,j)-xs(2,i))*(xs(2,j)-xs(2,i))+
269     + (xs(3,j)-xs(3,i))*(xs(3,j)-xs(3,i)))
270      if(dist.gt.argos_prep_atsiz(num(isat(j)))+
271     + argos_prep_atsiz(num(isat(i)))) then
272      write(lfnpdb,1101) csa(i)(11:14),csa(i)(1:3),char(icharr),
273     + isgm(i)+ioff+isoff,csa(j)(11:14),csa(j)(1:3),isgm(j)+ioff+isoff
274      endif
275      endif
276  111 continue
277c
278      nums=0
279c
280      do 7 i=1,nsa
281c
282      x(1)=xs(1,i)
283      x(2)=xs(2,i)
284      x(3)=xs(3,i)
285      if(ips(i).ge.0.or.ihop(i).eq.0) then
286      x(1)=xs(1,i)
287      x(2)=xs(2,i)
288      x(3)=xs(3,i)
289      else
290      x(1)=xs(1,i)-xs(1,i+ihop(i))
291      x(2)=xs(2,i)-xs(2,i+ihop(i))
292      x(3)=xs(3,i)-xs(3,i+ihop(i))
293      dist=sqrt(x(1)*x(1)+x(2)*x(2)+x(3)*x(3))
294      dist=dumdst/dist
295      x(1)=xs(1,i+ihop(i))+dist*x(1)
296      x(2)=xs(2,i+ihop(i))+dist*x(2)
297      x(3)=xs(3,i+ihop(i))+dist*x(3)
298      endif
299      y(1)=x(1)
300      y(2)=x(2)
301      y(3)=x(3)
302      if(irrand.eq.1.or.irrand.eq.4) call rotate(o,px,angle,x,y)
303      x(1)=y(1)
304      x(2)=y(2)
305      x(3)=y(3)
306      if(irrand.eq.2.or.irrand.eq.4) call rotate(o,py,angle,x,y)
307      x(1)=y(1)
308      x(2)=y(2)
309      x(3)=y(3)
310      if(irrand.eq.3.or.irrand.eq.4) call rotate(o,pz,angle,x,y)
311      x(1)=y(1)
312      x(2)=y(2)
313      x(3)=y(3)
314c
315      temp=4.009103873d1*amass(isat(i))*
316     + (vs(1,i)*vs(1,i)+vs(2,i)*vs(2,i)+vs(3,i)*vs(3,i))
317      if(temp.gt.999.0) temp=999.989
318      if(lrgpdb.eq.0) then
319      if(nzrep.gt.0) then
320      write(lfnpdb,1001) i+iaoff,csa(i)(11:14),csa(i)(1:3),
321     + char(icharr),isgm(i)+ioff+isoff,1.0d1*(y(1)+dbx),
322     + 1.0d1*(y(2)+dby),1.0d1*(y(3)+dbz),
323     + temp,argos_prep_atnam(num(isat(i)))
324      else
325      if(izrep.eq.1) then
326      write(lfnpdb,1001) i+iaoff,csa(i)(11:14),csa(i)(1:3),
327     + char(icharr),isgm(i)+ioff+isoff,1.0d1*(y(1)+dbx),
328     + 1.0d1*(y(2)+dby),1.0d1*(y(3)-dbz),
329     + temp,argos_prep_atnam(num(isat(i)))
330      else
331      write(lfnpdb,1001) i+iaoff,csa(i)(11:14),csa(i)(1:3),
332     + char(icharr),isgm(i)+ioff+isoff,-1.0d1*(y(1)+dbx),
333     + -1.0d1*(y(2)+dby),-1.0d1*(y(3)-dbz),
334     + temp,argos_prep_atnam(num(isat(i)))
335      endif
336      endif
337      else
338      if(char(icharr).ne.' ') call md_abort('Error in LRGPDB ',icharr)
339      if(nzrep.gt.0) then
340      write(lfnpdb,2001) i+iaoff,csa(i)(11:14),csa(i)(1:3),
341     + isgm(i)+ioff+isoff,1.0d1*(y(1)+dbx),
342     + 1.0d1*(y(2)+dby),1.0d1*(y(3)+dbz),
343     + temp,argos_prep_atnam(num(isat(i)))
344      else
345      if(izrep.eq.1) then
346      write(lfnpdb,2001) i+iaoff,csa(i)(11:14),csa(i)(1:3),
347     + isgm(i)+ioff+isoff,1.0d1*(y(1)+dbx),
348     + 1.0d1*(y(2)+dby),1.0d1*(y(3)-dbz),
349     + temp,argos_prep_atnam(num(isat(i)))
350      else
351      write(lfnpdb,2001) i+iaoff,csa(i)(11:14),csa(i)(1:3),
352     + isgm(i)+ioff+isoff,-1.0d1*(y(1)+dbx),
353     + -1.0d1*(y(2)+dby),-1.0d1*(y(3)-dbz),
354     + temp,argos_prep_atnam(num(isat(i)))
355      endif
356      endif
357      endif
358      if(isgm(i).gt.nums) nums=isgm(i)
359    7 continue
360      isoff=isoff+nums
361      iaoff=iaoff+nsa
362      icharr=icharr+1
363    6 continue
364    5 continue
365    4 continue
366c
367      else
368      natom=0
369      nsegm=0
370      isegm=0
371      isoff=0
372      do 20 l=1,nmols
373      ilang=0
374      do 14 izrep=1,iabs(nzrep)
375      if(nzrep.gt.0) then
376      dbz=(dble(izrep)-0.5d0*dble(nzrep+1))*drep(3)
377      else
378      dbz=dbz0
379      endif
380      do 15 iyrep=1,nyrep
381      dby=(dble(iyrep)-0.5d0*dble(nyrep+1))*drep(2)
382      do 16 ixrep=1,nxrep
383      dbx=(dble(ixrep)-0.5d0*dble(nxrep+1))*drep(1)
384c
385      do 67 k=1,nskip
386      if(ixrep.eq.iskip(1,k).and.iyrep.eq.iskip(2,k).and.
387     + izrep.eq.iskip(3,k)) goto 16
388   67 continue
389c
390      ilang=ilang+1
391c
392      if(l.eq.1) then
393      angle=8.0d0*atan(1.0d0)*util_random(0)
394      lang(ilang)=angle
395      else
396      angle=lang(ilang)
397      endif
398c
399      icharr=ichar(' ')
400c
401      do 212 i=1,nsa
402      if(iropt.eq.3) then
403      if(ifra(i).eq.l) then
404      isoff=nsegm-isgm(i)+1
405      goto 213
406      endif
407      else
408      if(imol(i).eq.l) then
409      isoff=nsegm-isgm(i)+1
410      goto 213
411      endif
412      endif
413  212 continue
414  213 continue
415c
416      do 211 k=1,nsb
417      i=idsb(1,k)
418      j=idsb(2,k)
419      if(iropt.eq.3) then
420      if(ifra(i).eq.l.and.ifra(j).eq.l.and.isgm(i).ne.isgm(j)) then
421      dist=sqrt((xs(1,j)-xs(1,i))*(xs(1,j)-xs(1,i))+
422     + (xs(2,j)-xs(2,i))*(xs(2,j)-xs(2,i))+
423     + (xs(3,j)-xs(3,i))*(xs(3,j)-xs(3,i)))
424      if(dist.gt.argos_prep_atsiz(num(isat(j)))+
425     + argos_prep_atsiz(num(isat(i)))) then
426      write(lfnpdb,1101) csa(i)(11:14),csa(i)(1:3),char(icharr),
427     + isgm(i)+ioff+isoff,csa(j)(11:14),csa(j)(1:3),isgm(j)+ioff+isoff
428      endif
429      endif
430      else
431      if(imol(i).eq.l.and.imol(j).eq.l.and.isgm(i).ne.isgm(j)) then
432      dist=sqrt((xs(1,j)-xs(1,i))*(xs(1,j)-xs(1,i))+
433     + (xs(2,j)-xs(2,i))*(xs(2,j)-xs(2,i))+
434     + (xs(3,j)-xs(3,i))*(xs(3,j)-xs(3,i)))
435      if(dist.gt.argos_prep_atsiz(num(isat(j)))+
436     + argos_prep_atsiz(num(isat(i)))) then
437      write(lfnpdb,1101) csa(i)(11:14),csa(i)(1:3),char(icharr),
438     + isgm(i)+ioff+isoff,csa(j)(11:14),csa(j)(1:3),isgm(j)+ioff+isoff
439      endif
440      endif
441      endif
442  211 continue
443c
444      nums=0
445c
446      do 17 i=1,nsa
447c
448      x(1)=xs(1,i)
449      x(2)=xs(2,i)
450      x(3)=xs(3,i)
451      if(ips(i).ge.0.or.ihop(i).eq.0) then
452      x(1)=xs(1,i)
453      x(2)=xs(2,i)
454      x(3)=xs(3,i)
455      else
456      x(1)=xs(1,i)-xs(1,i+ihop(i))
457      x(2)=xs(2,i)-xs(2,i+ihop(i))
458      x(3)=xs(3,i)-xs(3,i+ihop(i))
459      dist=sqrt(x(1)*x(1)+x(2)*x(2)+x(3)*x(3))
460      dist=dumdst/dist
461      x(1)=xs(1,i+ihop(i))+dist*x(1)
462      x(2)=xs(2,i+ihop(i))+dist*x(2)
463      x(3)=xs(3,i+ihop(i))+dist*x(3)
464      endif
465      y(1)=x(1)
466      y(2)=x(2)
467      y(3)=x(3)
468      if(irrand.eq.1.or.irrand.eq.4) call rotate(o,px,angle,x,y)
469      x(1)=y(1)
470      x(2)=y(2)
471      x(3)=y(3)
472      if(irrand.eq.2.or.irrand.eq.4) call rotate(o,py,angle,x,y)
473      x(1)=y(1)
474      x(2)=y(2)
475      x(3)=y(3)
476      if(irrand.eq.3.or.irrand.eq.4) call rotate(o,pz,angle,x,y)
477      x(1)=y(1)
478      x(2)=y(2)
479      x(3)=y(3)
480c
481      if(iropt.eq.3) then
482      if(ifra(i).eq.l) then
483      temp=4.009103873d1*amass(isat(i))*
484     + (vs(1,i)*vs(1,i)+vs(2,i)*vs(2,i)+vs(3,i)*vs(3,i))
485      if(temp.gt.999.0) temp=999.989
486      natom=natom+1
487      if(lrgpdb.eq.0) then
488      if(nzrep.gt.0) then
489      write(lfnpdb,1001) natom,csa(i)(11:14),csa(i)(1:3),
490     + char(icharr),isgm(i)+ioff+isoff,1.0d1*(y(1)+dbx),
491     + 1.0d1*(y(2)+dby),1.0d1*(y(3)+dbz),
492     + temp,argos_prep_atnam(num(isat(i)))
493      else
494      if(izrep.eq.1) then
495      write(lfnpdb,1001) natom,csa(i)(11:14),csa(i)(1:3),
496     + char(icharr),isgm(i)+ioff+isoff,1.0d1*(y(1)+dbx),
497     + 1.0d1*(y(2)+dby),1.0d1*(y(3)-dbz),
498     + temp,argos_prep_atnam(num(isat(i)))
499      else
500      write(lfnpdb,1001) natom,csa(i)(11:14),csa(i)(1:3),
501     + char(icharr),isgm(i)+ioff+isoff,-1.0d1*(y(1)+dbx),
502     + -1.0d1*(y(2)+dby),-1.0d1*(y(3)-dbz),
503     + temp,argos_prep_atnam(num(isat(i)))
504      endif
505      endif
506      else
507      if(char(icharr).ne.' ') call md_abort('Error in LRGPDB ',icharr)
508      if(nzrep.gt.0) then
509      write(lfnpdb,2001) natom,csa(i)(11:14),csa(i)(1:3),
510     + isgm(i)+ioff+isoff,1.0d1*(y(1)+dbx),
511     + 1.0d1*(y(2)+dby),1.0d1*(y(3)+dbz),
512     + temp,argos_prep_atnam(num(isat(i)))
513      else
514      if(izrep.eq.1) then
515      write(lfnpdb,2001) natom,csa(i)(11:14),csa(i)(1:3),
516     + isgm(i)+ioff+isoff,1.0d1*(y(1)+dbx),
517     + 1.0d1*(y(2)+dby),1.0d1*(y(3)-dbz),
518     + temp,argos_prep_atnam(num(isat(i)))
519      else
520      write(lfnpdb,2001) natom,csa(i)(11:14),csa(i)(1:3),
521     + isgm(i)+ioff+isoff,-1.0d1*(y(1)+dbx),
522     + -1.0d1*(y(2)+dby),-1.0d1*(y(3)-dbz),
523     + temp,argos_prep_atnam(num(isat(i)))
524      endif
525      endif
526      endif
527      nsegm=isgm(i)+isoff
528      endif
529      else
530      if(imol(i).eq.l) then
531      temp=4.009103873d1*amass(isat(i))*
532     + (vs(1,i)*vs(1,i)+vs(2,i)*vs(2,i)+vs(3,i)*vs(3,i))
533      if(temp.gt.999.0) temp=999.989
534      natom=natom+1
535      if(lrgpdb.eq.0) then
536      if(nzrep.gt.0) then
537      write(lfnpdb,1001) natom,csa(i)(11:14),csa(i)(1:3),
538     + char(icharr),isgm(i)+ioff+isoff,1.0d1*(y(1)+dbx),
539     + 1.0d1*(y(2)+dby),1.0d1*(y(3)+dbz),
540     + temp,argos_prep_atnam(num(isat(i)))
541      else
542      if(izrep.eq.1) then
543      write(lfnpdb,1001) natom,csa(i)(11:14),csa(i)(1:3),
544     + char(icharr),isgm(i)+ioff+isoff,1.0d1*(y(1)+dbx),
545     + 1.0d1*(y(2)+dby),1.0d1*(y(3)-dbz),
546     + temp,argos_prep_atnam(num(isat(i)))
547      else
548      write(lfnpdb,1001) natom,csa(i)(11:14),csa(i)(1:3),
549     + char(icharr),isgm(i)+ioff+isoff,-1.0d1*(y(1)+dbx),
550     + -1.0d1*(y(2)+dby),-1.0d1*(y(3)-dbz),
551     + temp,argos_prep_atnam(num(isat(i)))
552      endif
553      endif
554      else
555      if(char(icharr).ne.' ') call md_abort('Error in LRGPDB ',icharr)
556      if(nzrep.gt.0) then
557      write(lfnpdb,2001) natom,csa(i)(11:14),csa(i)(1:3),
558     + isgm(i)+ioff+isoff,1.0d1*(y(1)+dbx),
559     + 1.0d1*(y(2)+dby),1.0d1*(y(3)+dbz),
560     + temp,argos_prep_atnam(num(isat(i)))
561      else
562      if(izrep.eq.1) then
563      write(lfnpdb,2001) natom,csa(i)(11:14),csa(i)(1:3),
564     + isgm(i)+ioff+isoff,1.0d1*(y(1)+dbx),
565     + 1.0d1*(y(2)+dby),1.0d1*(y(3)-dbz),
566     + temp,argos_prep_atnam(num(isat(i)))
567      else
568      write(lfnpdb,2001) natom,csa(i)(11:14),csa(i)(1:3),
569     + isgm(i)+ioff+isoff,-1.0d1*(y(1)+dbx),
570     + -1.0d1*(y(2)+dby),-1.0d1*(y(3)-dbz),
571     + temp,argos_prep_atnam(num(isat(i)))
572      endif
573      endif
574      endif
575      nsegm=isgm(i)+isoff
576      endif
577      endif
578   17 continue
579   16 continue
580   15 continue
581   14 continue
582   20 continue
583      endif
584c
585      endif
586c
587      write(lfnpdb,1002)
588 1002 format('TER')
589c
590      numslv=iopt
591      if(iopt.lt.0) numslv=nwmc
592      if(iopt.gt.nwm+nwmc) numslv=nwm+nwmc
593c
594      if(numslv.gt.0) then
595c
596      if(nxrep.le.1.and.nyrep.le.1.and.iabs(nzrep).le.1) then
597c
598      nats=nsa
599      do 2 i=1,numslv
600      do 3 j=1,nwa
601      nats=nats+1
602      if(i.le.nwmc) then
603      temp=4.009103873d1*amass(iwat(j))*(vwc(1,j,i)*vwc(1,j,i)+
604     + vwc(2,j,i)*vwc(2,j,i)+vwc(3,j,i)*vwc(3,j,i))
605      if(temp.gt.999.0) temp=999.989
606      write(lfnpdb,1003) nats,cwa(j)(11:14),slvnam,
607     + i+nums,(1.0d1*xwc(k,j,i),k=1,3),temp,
608     + argos_prep_atnam(num(iwat(j)))
609      else
610      temp=4.009103873d1*amass(iwat(j))*(vw(1,j,i-nwmc)*vw(1,j,i-nwmc)+
611     + vw(2,j,i-nwmc)*vw(2,j,i-nwmc)+vw(3,j,i-nwmc)*vw(3,j,i-nwmc))
612      if(temp.gt.999.0) temp=999.989
613      write(lfnpdb,1003) nats,cwa(j)(11:14),slvnam,
614     + i+nums,(1.0d1*xw(k,j,i-nwmc),k=1,3),temp,
615     + argos_prep_atnam(num(iwat(j)))
616      endif
617 1003 format('ATOM',i7,1x,a4,1x,a3,i6,4x,3f8.3,6x,f6.2,4x,a2)
618    3 continue
619    2 continue
620c
621      else
622c
623      nats=nxrep*nyrep*iabs(nzrep)*nsa
624      ilang=0
625      do 44 izrep=1,iabs(nzrep)
626      if(nzrep.gt.0) then
627      dbz=(dble(izrep)-0.5d0*dble(nzrep+1))*drep(3)
628      else
629      dbz=dbz0
630      endif
631      do 45 iyrep=1,nyrep
632      dby=(dble(iyrep)-0.5d0*dble(nyrep+1))*drep(2)
633      do 46 ixrep=1,nxrep
634      dbx=(dble(ixrep)-0.5d0*dble(nxrep+1))*drep(1)
635c
636      ilang=ilang+1
637      angle=lang(ilang)
638c
639      nums=isoff
640      do 47 i=1,numslv
641      do 48 j=1,nwa
642      nats=nats+1
643c
644      if(i.le.nwmc) then
645      x(1)=xwc(1,j,i)
646      x(2)=xwc(2,j,i)
647      x(3)=xwc(3,j,i)
648      y(1)=xwc(1,j,i)
649      y(2)=xwc(2,j,i)
650      y(3)=xwc(3,j,i)
651      v(1)=vwc(1,j,i)
652      v(2)=vwc(2,j,i)
653      v(3)=vwc(3,j,i)
654      else
655      x(1)=xw(1,j,i-nwmc)
656      x(2)=xw(2,j,i-nwmc)
657      x(3)=xw(3,j,i-nwmc)
658      y(1)=xw(1,j,i-nwmc)
659      y(2)=xw(2,j,i-nwmc)
660      y(3)=xw(3,j,i-nwmc)
661      v(1)=vw(1,j,i-nwmc)
662      v(2)=vw(2,j,i-nwmc)
663      v(3)=vw(3,j,i-nwmc)
664      endif
665      if(irrand.eq.1.or.irrand.eq.4) call rotate(o,px,angle,x,y)
666      x(1)=y(1)
667      x(2)=y(2)
668      x(3)=y(3)
669      if(irrand.eq.2.or.irrand.eq.4) call rotate(o,py,angle,x,y)
670      x(1)=y(1)
671      x(2)=y(2)
672      x(3)=y(3)
673      if(irrand.eq.3.or.irrand.eq.4) call rotate(o,pz,angle,x,y)
674      x(1)=y(1)
675      x(2)=y(2)
676      x(3)=y(3)
677c
678      if(nzrep.gt.0) then
679      temp=4.009103873d1*amass(iwat(j))*(v(1)*v(1)+
680     + v(2)*v(2)+v(3)*v(3))
681      if(temp.gt.999.0) temp=999.989
682      write(lfnpdb,1003) nats,cwa(j)(11:14),slvnam,
683     + i+nums,1.0d1*(y(1)+dbx),1.0d1*(y(2)+dby),
684     + 1.0d1*(y(3)+dbz),temp,argos_prep_atnam(num(iwat(j)))
685      else
686      if(izrep.eq.1) then
687      temp=4.009103873d1*amass(iwat(j))*(v(1)*v(1)+
688     + v(2)*v(2)+v(3)*v(3))
689      if(temp.gt.999.0) temp=999.989
690      write(lfnpdb,1003) nats,cwa(j)(11:14),slvnam,
691     + i+nums,1.0d1*(y(1)+dbx),1.0d1*(y(2)+dby),
692     + 1.0d1*(y(3)-dbz),temp,argos_prep_atnam(num(iwat(j)))
693      else
694      temp=4.009103873d1*amass(iwat(j))*(v(1)*v(1)+
695     + v(2)*v(2)+v(3)*v(3))
696      if(temp.gt.999.0) temp=999.989
697      write(lfnpdb,1003) nats,cwa(j)(11:14),slvnam,
698     + i+nums,1.0d1*(y(1)+dbx),1.0d1*(y(2)+dby),
699     + 1.0d1*(y(3)-dbz),temp,argos_prep_atnam(num(iwat(j)))
700      endif
701      endif
702   48 continue
703   47 continue
704      nums=nums+numslv
705   46 continue
706   45 continue
707   44 continue
708c
709      endif
710c
711      endif
712c
713      write(lfnpdb,1004)
714 1004 format('END')
715c
716c
717      close(unit=lfnpdb)
718c
719      if(util_print('files',print_default)) then
720      write(lfnout,3000) filpdb(1:length)
721 3000 format(' Created pdb',t40,a,/)
722      endif
723c
724      argos_prep_wrtpdb=.true.
725      return
726c
727 9999 continue
728      argos_prep_wrtpdb=.false.
729      return
730      end
731