source/analysis/src-M-prop/read_ref_str.F

   1       subroutine read_ref_structure(*)
   2 C
   3 C Read the reference structure from the PDB file or from a PDB file or in the form of the dihedral
   4 C angles.
   5 C
   6       implicit none
   7       include 'DIMENSIONS'
   8       include 'sizesclu.dat'
   9       include 'COMMON.IOUNITS'
  10       include 'COMMON.GEO'
  11       include 'COMMON.VAR'
  12       include 'COMMON.INTERACT'
  13       include 'COMMON.LOCAL'
  14       include 'COMMON.NAMES'
  15       include 'COMMON.CHAIN'
  16       include 'COMMON.FFIELD'
  17       include 'COMMON.SBRIDGE'
  18       include 'COMMON.HEADER'
  19       include 'COMMON.CONTROL'
  20       include 'COMMON.TIME1'
  21       character*4 sequence(maxres)
  22       integer rescode
  23       double precision x(maxvar)
  24       integer itype_pdb(maxres)
  25       logical seq_comp
  26       integer i,j,k,nres_pdb,iaux
  27       double precision ddsc,dist
  28       integer ilen
  29       external ilen
  30 C
  31       nres0=nres
  32       write (iout,*) "pdbref",pdbref
  33       if (pdbref) then
  34         read(inp,'(a)') pdbfile
  35         write (iout,'(2a,1h.)') 'PDB data will be read from file ',
  36      &    pdbfile(:ilen(pdbfile))
  37         call flush(iout)
  38         open(ipdbin,file=pdbfile,status='old',err=33)
  39         goto 34
  40   33    write (iout,'(a)') 'Error opening PDB file.'
  41         return1
  42   34    continue
  43         do i=1,nres
  44           itype_pdb(i)=itype(i)
  45         enddo
  46         call readpdb(.true.)
  47         do i=1,nres
  48           iaux=itype_pdb(i)
  49           itype_pdb(i)=itype(i)
  50           itype(i)=iaux
  51         enddo
  52         close (ipdbin)
  53         nres_pdb=nres
  54         nres=nres0
  55         nstart_seq=nnt
  56         if (nsup.le.(nct-nnt+1)) then
  57           do i=0,nct-nnt+1-nsup
  58             if (seq_comp(itype(nnt+i),itype_pdb(nstart_sup),
  59      &        nsup)) then
  60               do j=nnt+nsup-1,nnt,-1
  61                 do k=1,3
  62                   cref(k,nres+j+i,1)=cref(k,nres_pdb+j,1)
  63                 enddo
  64               enddo
  65               do j=nnt+nsup-1,nnt,-1
  66                 do k=1,3
  67                   cref(k,j+i,1)=cref(k,j,1)
  68                 enddo
  69                 phi_ref(j+i)=phi_ref(j)
  70                 theta_ref(j+i)=theta_ref(j)
  71                 alph_ref(j+i)=alph_ref(j)
  72                 omeg_ref(j+i)=omeg_ref(j)
  73               enddo
  74               nstart_seq=nnt+i
  75               nstart_sup=nnt+i
  76               goto 111
  77             endif
  78           enddo
  79 #define DEBUG
  80 #ifdef DEBUG
  81           do j=nnt,nct
  82             write (iout,'(i5,3f10.5,5x,3f10.5)')
  83      &         j,(cref(k,j,1),k=1,3),(cref(k,j+nres,1),k=1,3)
  84           enddo
  85           call flush(iout)
  86 #endif
  87 #undef DEBUG
  88           write (iout,'(a)')
  89      &            'Error - sequences to be superposed do not match.'
  90           return1
  91         else
  92           do i=0,nsup-(nct-nnt+1)
  93             if (seq_comp(itype(nnt),itype_pdb(nstart_sup+i),
  94      &        nct-nnt+1))
  95      &      then
  96               nstart_sup=nstart_sup+i
  97               nsup=nct-nnt+1
  98               goto 111
  99             endif
 100           enddo
 101           write (iout,'(a)')
 102      &            'Error - sequences to be superposed do not match.'
 103         endif
 104   111   continue
 105         write (iout,'(a,i5)')
 106      &   'Experimental structure begins at residue',nstart_seq
 107 #define DEBUG
 108 #ifdef DEBUG
 109         do j=nnt,nct
 110           write (iout,'(i5,3f10.5,5x,3f10.5)')
 111      &        j,(cref(k,j,1),k=1,3),(cref(k,j+nres,1),k=1,3)
 112         enddo
 113         call flush(iout)
 114 #endif
 115 #undef DEBUG
 116       else
 117         call read_angles(inp,*38)
 118         goto 39
 119    38   write (iout,'(a)') 'Error reading reference structure.'
 120         return1
 121    39   call chainbuild
 122         nstart_sup=nnt
 123         nstart_seq=nnt
 124         nsup=nct-nnt+1
 125         do i=1,2*nres
 126           do j=1,3
 127             cref(j,i,1)=c(j,i)
 128           enddo
 129         enddo
 130       endif
 131       nend_sup=nstart_sup+nsup-1
 132       do i=1,2*nres
 133         do j=1,3
 134           c(j,i)=cref(j,i,1)
 135         enddo
 136       enddo
 137       do i=1,nres
 138         do j=1,3
 139           dc(j,nres+i)=cref(j,nres+i,1)-cref(j,i,1)
 140         enddo
 141         if (itype(i).ne.10) then
 142           ddsc = dist(i,nres+i)
 143           do j=1,3
 144             dc_norm(j,nres+i)=dc(j,nres+i)/ddsc
 145           enddo
 146         else
 147           do j=1,3
 148             dc_norm(j,nres+i)=0.0d0
 149           enddo
 150         endif
 151 c        write (iout,*) "i",i," dc_norm",(dc_norm(k,nres+i),k=1,3),
 152 c         " norm",dc_norm(1,nres+i)**2+dc_norm(2,nres+i)**2+
 153 c         dc_norm(3,nres+i)**2
 154         do j=1,3
 155           dc(j,i)=c(j,i+1)-c(j,i)
 156         enddo
 157         ddsc = dist(i,i+1)
 158         do j=1,3
 159           dc_norm(j,i)=dc(j,i)/ddsc
 160         enddo
 161       enddo
 162       write (iout,'(a,i3,a,i3,a,i3,a)')
 163      &    'Number of residues to be superposed:',nsup,
 164      &    ' (from residue',nstart_sup,' to residue',
 165      &    nend_sup,').'
 166       call flush(iout)
 167       return
 168       end