! include 'COMMON.CHAIN'
! include 'COMMON.INTERACT'
! include 'COMMON.SBRIDGE'
- use geometry_data, only: nres,c,nfrag
- use energy_data, only: itype,nnt,nct,nss,ihpb,jhpb,uconst_back
+ use geometry_data, only: nres,c,nfrag,boxxsize,boxysize,boxzsize
+ use energy_data, only: itype,nnt,nct,nss,ihpb,jhpb,uconst_back,molnum
use io_base!, only: maxres,ucase,iblnk,rescode,pdbout
-
+! use geometry, only: returnbox
implicit none
real(kind=4),allocatable,dimension(:,:) :: coord !(3,2*maxres)
real(kind=4) :: prec,time,potE,uconst,t_bath,qfrag(100)
real(kind=8) :: etot
- character(len=80) :: arg,seqfile,pdbfileX
+ character(len=80) :: arg,seqfile,pdbfileX,boxtraj
character(len=3),allocatable,dimension(:) :: sequenc !(maxres)
character(len=50) :: tytul
character(len=8) :: onethree,cfreq,cntraj,citraj
! external rescode
integer :: i,ii,j,k,kk,is,ie,ifreq,mnum,molec,iext,isize
integer :: iret,ixdrf,itraj,ntraj
-
+ print *,"begin program"
+! if(.not.allocated(molnum)) allocate(molnum(maxres))
ifreq=1
nres=0
mnum=1
molec=1
allocate(sequenc(maxres))
+ allocate(molnum(maxres))
allocate(itype(maxres,5))
- if (iargc().lt.5) then
+ if (iargc().lt.6) then
print '(2a)',&
"Usage: xdrf2pdb-m one/three seqfile cxfile ntraj itraj",&
- " [pdbfile] [freq]"
+ "boxsize [pdbfile] [freq]"
stop
endif
call getarg(1,onethree)
itype(i,mnum)=rescode(i,sequenc(i),1,molec)
enddo
else
- read(1,'(20(a3,1x))',end=11,err=11) (sequenc(i),i=1,maxres)
+ read(1,'(20(a3,i1))',end=11,err=11) (sequenc(i),molnum(i),i=1,maxres)
11 continue
nres=i
i=0
enddo
nres=i
do i=1,nres
- itype(i,mnum)=rescode(i,sequenc(i),0,molec)
+ itype(i,molnum(i))=rescode(i,sequenc(i),0,molnum(i))
enddo
print *,nres
print '(20(a3,1x))',(sequenc(i),i=1,nres)
read (cntraj,*) ntraj
call getarg(5,citraj)
read (citraj,*) itraj
- if (iargc().gt.5) then
- call getarg(6,pdbfileX)
+ call getarg(6,boxtraj)
+ read (boxtraj,*) boxxsize
+ boxysize=boxxsize
+ boxzsize=boxxsize
+ if (iargc().gt.6) then
+ call getarg(7,pdbfileX)
else
write(licz,'(bz,i3.3)') itraj
pdbfileX=arg(:iext)//'_'//licz//'.pdb'
endif
- if (iargc().gt.6) then
- call getarg(7,cfreq)
+ if (iargc().gt.7) then
+ call getarg(8,cfreq)
read (cfreq,*) ifreq
endif
! print *,"ifreq",ifreq," ntraj",ntraj," itraj",itraj
open(9,file=pdbfileX)
nnt = 1
+ mnum=molnum(1)
if (itype(1,mnum).eq.ntyp1) nnt = 2
nct=nres
+ mnum=molnum(nct)
if (itype(nres,mnum).eq.ntyp1) nct = nres-1
print *,"nnt",nnt," nct",nct
call xdrfopen(ixdrf,arg, "r", iret)
! write (*,'(e15.8,2e15.5,f12.5,$)') time,potE,uconst,t_bath
! write (*,'(i4,$)') nss,(ihpb(j),jhpb(j),j=1,nss)
! write (*,'(i4,20f7.4)') nfrag,(qfrag(i),i=1,nfrag)
-! write (*,'(8f10.5)') ((coord(k,j),k=1,3),j=1,isize)
+ do j=1,isize
+ write (*,'(3f10.5,i3)') (coord(k,j),k=1,3),j
+ enddo
if (mod(kk/ntraj,ifreq).eq.0 .and. mod(kk,ntraj).eq.itraj) then
if (isize .ne. nres+nct-nnt+1) then
print *,"Error: inconsistent sizes",isize,nres+nct-nnt+1
c(j,i+nres)=coord(j,ii+nres)
enddo
enddo
+ call returnbox
etot=potE
write (tytul,'(a,i6)') "Structure",kk
call pdbout(etot,tytul,9)
enddo
end
+ subroutine returnbox
+ use geometry_data, only: nres,c,nfrag,boxxsize,boxysize,boxzsize
+ use energy_data, only: itype,nnt,nct,nss,ihpb,jhpb,uconst_back,molnum
+ use io_base!, only: maxres,ucase,iblnk,rescode,pdbout
+ integer :: allareout,i,j,k,nojumpval,chain_beg,mnum
+ integer :: chain_end,ireturnval
+ real*8 :: difference,xi,boxsize,x,xtemp,box2shift
+
+ do i=1,nres
+ if (molnum(i).eq.5) then
+ c(1,i)=dmod(c(1,i),boxxsize)
+! if (c(1,i).lt.0) c(1,i)=c(1,i)+boxxsize
+ x=c(1,i)-c(1,2)
+ print *,"return box,before wrap",c(1,i)
+ boxsize=boxxsize
+ xtemp=dmod(x,boxsize)
+ if (dabs(xtemp-boxsize).lt.dabs(xtemp)) then
+ box2shift=xtemp-boxsize
+ else if (dabs(xtemp+boxsize).lt.dabs(xtemp)) then
+ box2shift=xtemp+boxsize
+ else
+ box2shift=xtemp
+ endif
+ xi=box2shift
+! print *,c(1,2),xi,xtemp,box2shift
+ c(1,i)=c(1,2)+xi
+ print *,"after",c(1,2),xi,xtemp,box2shift,c(1,i)
+
+ c(2,i)=dmod(c(2,i),boxysize)
+! if (c(2,i).lt.0) c(2,i)=c(2,i)+boxysize
+ x=c(2,i)-c(2,2)
+ boxsize=boxysize
+ xtemp=dmod(x,boxsize)
+ if (dabs(xtemp-boxsize).lt.dabs(xtemp)) then
+ box2shift=xtemp-boxsize
+ else if (dabs(xtemp+boxsize).lt.dabs(xtemp)) then
+ box2shift=xtemp+boxsize
+ else
+ box2shift=xtemp
+ endif
+ xi=box2shift
+ c(2,i)=c(2,2)+xi
+! c(3,i)=dmod(c(3,i),boxzsize)
+ if (c(3,i).lt.0) c(3,i)=c(3,i)+boxzsize
+ x=c(3,i)-c(3,2)
+ boxsize=boxzsize
+ xtemp=dmod(x,boxsize)
+ if (dabs(xtemp-boxsize).lt.dabs(xtemp)) then
+ box2shift=xtemp-boxsize
+ else if (dabs(xtemp+boxsize).lt.dabs(xtemp)) then
+ box2shift=xtemp+boxsize
+ else
+ box2shift=xtemp
+ endif
+ xi=box2shift
+ c(3,i)=c(3,2)+xi
+
+ c(1,i+nres)=dmod(c(1,i+nres),boxxsize)
+ c(2,i+nres)=dmod(c(2,i+nres),boxysize)
+ c(3,i+nres)=dmod(c(3,i+nres),boxzsize)
+ endif
+ enddo
+ return
+ end subroutine returnbox
+
projects / unres4.git / blobdiff