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