source/xdrfpdb/src/xdrf2pdb.F

   1       implicit real*8 (a-h,o-z)
   2       include 'DIMENSIONS'
   3       include 'COMMON.CHAIN'
   4       include 'COMMON.INTERACT'
   5       include 'COMMON.SBRIDGE'
   6       include 'COMMON.DISTFIT'
   7       real*4 coord(3,2*maxres)
   8       real*4 prec,time,potE,uconst,t_bath,qfrag(100)
   9       real*8 etot
  10       character*80 arg,seqfile,pdbfile
  11       character*3 sequenc(maxres)
  12       character*50 tytul
  13       character*8 onethree,cfreq
  14       character*8 ucase
  15       external ucase
  16       logical oneletter
  17       integer rescode
  18       external rescode
  19       logical iblnk
  20       external iblnk
  21       double precision cm(3)
  22
  23       ifreq=1
  24       is=1
  25       ie=1000000000
  26       if (iargc().lt.3) then
  27         print '(2a)',
  28      &    "Usage: xdrf2pdb one/three seqfile cxfile [freq]",
  29      &    " [start] [end] [pdbfile]"
  30         stop
  31       endif
  32       call getarg(1,onethree)
  33       onethree = ucase(onethree)
  34       if (onethree.eq.'ONE') then
  35         oneletter = .true.
  36       else if (onethree.eq.'THREE') then
  37         oneletter = .false.
  38       else
  39         print *,"ONE or THREE must be specified"
  40       endif
  41       call getarg(2,seqfile)
  42       open (1,file=seqfile,status='old')
  43       if (oneletter) then
  44         read(1,'(80a1)',end=10,err=10) (sequenc(i)(1:1),i=1,maxres)
  45    10   continue
  46         nres=i
  47         i=0
  48         do while (.not.iblnk(sequenc(i+1)(1:1)))
  49           i=i+1
  50         enddo
  51         nres=i
  52         do i=1,nres
  53           itype(i)=rescode(i,sequenc(i),1)
  54         enddo
  55       else
  56         read(1,'(20(a3,1x))',end=11,err=11) (sequenc(i),i=1,maxres)
  57    11   continue
  58         nres=i
  59         i=0
  60         do while (.not.iblnk(sequenc(i+1)(1:1)))
  61           i=i+1
  62         enddo
  63         nres=i
  64         do i=1,nres
  65           itype(i)=rescode(i,sequenc(i),0)
  66         enddo
  67       endif
  68       call getarg(3,arg)
  69       iext = index(arg,'.cx') - 1
  70       if (iext.lt.0) then
  71         print *,"Error - not a cx file"
  72         stop
  73       endif
  74       if (iargc().gt.3) then
  75         call getarg(4,cfreq)
  76         read (cfreq,*) ifreq
  77       endif
  78       if (iargc().gt.4) then
  79         call getarg(5,cfreq)
  80         read (cfreq,*) is
  81       endif
  82       if (iargc().gt.5) then
  83         call getarg(6,cfreq)
  84         read (cfreq,*) ie
  85       endif
  86       if (iargc().gt.6) then
  87         call getarg(7,pdbfile)
  88       else
  89         pdbfile=arg(:iext)//'.pdb'
  90       endif
  91       open(9,file=pdbfile)
  92       nnt = 1
  93       if (itype(1).eq.21) nnt = 2
  94       nct=nres
  95       if (itype(nres).eq.21) nct = nres-1
  96       print *,"nnt",nnt," nct",nct
  97       print *,"file",arg
  98       do i=nnt,nct
  99         if (itype(i).ne.20) then
 100           itel(i)=1
 101         else
 102           itel(i)=2
 103         endif
 104       enddo
 105       call xdrfopen(ixdrf,arg, "r", iret)
 106       if (iret.eq.0) stop
 107 c      print *,"iret",iret
 108 c      print *,"is",is," ie",ie
 109       kk = 0
 110       do while(is.eq.0 .or. kk.lt.ie)
 111        call xdrffloat(ixdrf, time, iret)
 112 c       print *,"time",time," iret",iret
 113        if(iret.eq.0) exit
 114        kk = kk + 1
 115        call xdrffloat(ixdrf, potE, iret)
 116        call xdrffloat(ixdrf, uconst, iret)
 117 c       print *,"potE",potE," uconst",uconst
 118 #ifdef NEWUNRES
 119        call xdrffloat(ixdrf, uconst_back, iret)
 120 c       print *,"uconst_back",uconst_back
 121 #endif
 122        call xdrffloat(ixdrf, t_bath, iret)
 123 c       print *,"t_bath",t_bath
 124        call xdrfint(ixdrf, nss, iret)
 125        do j=1,nss
 126         call xdrfint(ixdrf, ihpb(j), iret)
 127         call xdrfint(ixdrf, jhpb(j), iret)
 128        enddo
 129 c       print *,"nss",nss
 130 c       print *,"   ",(ihpb(j),"-",jhpb(j),j=1,nss)
 131        call xdrfint(ixdrf, nfrag, iret)
 132        do i=1,nfrag
 133         call xdrffloat(ixdrf, qfrag(i), iret)
 134        enddo
 135 c       print *,"nfrag",nfrag
 136        prec=10000.0
 137
 138        isize=0
 139        call xdrf3dfcoord(ixdrf, coord, isize, prec, iret)
 140
 141
 142 c       write (*,'(e15.8,2e15.5,f12.5,$)') time,potE,uconst,t_bath
 143 c       write (*,'(i4,$)') nss,(ihpb(j),jhpb(j),j=1,nss)
 144 c       write (*,'(i4,20f7.4)') nfrag,(qfrag(i),i=1,nfrag)
 145 c       write (*,'(8f10.5)') ((coord(k,j),k=1,3),j=1,isize)
 146        if (kk.ge.is .and. mod(kk,ifreq).eq.0) then
 147          if (isize .ne. nres+nct-nnt+1) then
 148            print *,"Error: inconsistent sizes",isize,nres+nct-nnt+1
 149          endif
 150          do i=1,nres
 151            do j=1,3
 152              c(j,i)=coord(j,i)
 153            enddo
 154          enddo
 155          ii = 0
 156          do i=nnt,nct
 157            ii = ii + 1
 158            do j=1,3
 159              c(j,i+nres)=coord(j,ii+nres)
 160            enddo
 161          enddo
 162 c Calculate the CM
 163          do j=1,3
 164            cm(j)=0.0d0
 165            do i=1,nres
 166              cm(j)=cm(j)+c(j,i)
 167            enddo
 168            cm(j)=cm(j)/nres
 169          enddo
 170          do i=1,nres
 171            do j=1,3
 172              c(j,i)=c(j,i)-cm(j)
 173              c(j,i+nres)=c(j,i+nres)-cm(j)
 174            enddo
 175          enddo
 176          etot=potE
 177          write (tytul,'(a,i6)') "Structure",kk
 178 #ifdef SECONDARY
 179          call secondary2(.false.)
 180 #endif
 181          call pdbout(etot,tytul,9)
 182        endif
 183       enddo
 184
 185       end