source/cluster/wham/src-M-SAXS-homology/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 c      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         open(ipdbin,file=pdbfile,status='old',err=33)
  38         goto 34
  39   33    write (iout,'(a)') 'Error opening PDB file.'
  40         return1
  41   34    continue
  42         do i=1,nres
  43           itype_pdb(i)=itype(i)
  44         enddo
  45         call readpdb(.true.)
  46         do i=1,2*nres
  47           do j=1,3
  48             cref_pdb(j,i)=c(j,i)
  49           enddo
  50         enddo
  51         do i=1,nres
  52           iaux=itype_pdb(i)
  53           itype_pdb(i)=itype(i)
  54           itype(i)=iaux
  55         enddo
  56         close (ipdbin)
  57         nres_pdb=nres
  58         nres=nres0
  59         nstart_seq=nnt
  60         if (nsup.le.(nct-nnt+1)) then
  61           do i=0,nct-nnt+1-nsup
  62             if (seq_comp(itype(nnt+i),itype_pdb(nstart_sup),
  63      &        nsup)) then
  64               do j=nnt+nsup-1,nnt,-1
  65                 do k=1,3
  66                   cref_pdb(k,nres+j+i)=cref_pdb(k,nres_pdb+j)
  67                 enddo
  68               enddo
  69               do j=nnt+nsup-1,nnt,-1
  70                 do k=1,3
  71                   cref_pdb(k,j+i)=cref_pdb(k,j)
  72                 enddo
  73                 phi_ref(j+i)=phi_ref(j)
  74                 theta_ref(j+i)=theta_ref(j)
  75                 alph_ref(j+i)=alph_ref(j)
  76                 omeg_ref(j+i)=omeg_ref(j)
  77               enddo
  78 #ifdef DEBUG
  79               do j=nnt,nct
  80                 write (iout,'(i5,3f10.5,5x,3f10.5)')
  81      &            j,(cref_pdb(k,j),k=1,3),(cref_pdb(k,j+nres),k=1,3)
  82               enddo
  83 #endif
  84               nstart_seq=nnt+i
  85               nstart_sup=nnt+i
  86               goto 111
  87             endif
  88           enddo
  89           write (iout,'(a)')
  90      &            'Error - sequences to be superposed do not match.'
  91           return1
  92         else
  93           do i=0,nsup-(nct-nnt+1)
  94             if (seq_comp(itype(nnt),itype_pdb(nstart_sup+i),
  95      &        nct-nnt+1))
  96      &      then
  97               nstart_sup=nstart_sup+i
  98               nsup=nct-nnt+1
  99               goto 111
 100             endif
 101           enddo
 102           write (iout,'(a)')
 103      &            'Error - sequences to be superposed do not match.'
 104         endif
 105   111   continue
 106         write (iout,'(a,i5)')
 107      &   'Experimental structure begins at residue',nstart_seq
 108       else
 109         call read_angles(inp,*38)
 110         goto 39
 111    38   write (iout,'(a)') 'Error reading reference structure.'
 112         return1
 113    39   call chainbuild
 114         nstart_sup=nnt
 115         nstart_seq=nnt
 116         nsup=nct-nnt+1
 117         do i=1,2*nres
 118           do j=1,3
 119             cref_pdb(j,i)=c(j,i)
 120           enddo
 121         enddo
 122       endif
 123       nend_sup=nstart_sup+nsup-1
 124       do i=1,2*nres
 125         do j=1,3
 126           c(j,i)=cref_pdb(j,i)
 127         enddo
 128       enddo
 129       do i=1,nres
 130         do j=1,3
 131           dc(j,nres+i)=cref_pdb(j,nres+i)-cref_pdb(j,i)
 132         enddo
 133         if (itype(i).ne.10) then
 134           ddsc = dist(i,nres+i)
 135           do j=1,3
 136             dc_norm(j,nres+i)=dc(j,nres+i)/ddsc
 137           enddo
 138         else
 139           do j=1,3
 140             dc_norm(j,nres+i)=0.0d0
 141           enddo
 142         endif
 143 c        write (iout,*) "i",i," dc_norm",(dc_norm(k,nres+i),k=1,3),
 144 c         " norm",dc_norm(1,nres+i)**2+dc_norm(2,nres+i)**2+
 145 c         dc_norm(3,nres+i)**2
 146         do j=1,3
 147           dc(j,i)=c(j,i+1)-c(j,i)
 148         enddo
 149         ddsc = dist(i,i+1)
 150         do j=1,3
 151           dc_norm(j,i)=dc(j,i)/ddsc
 152         enddo
 153       enddo
 154        write (iout,'(a,i3,a,i3,a,i3,a)')
 155      &    'Number of residues to be superposed:',nsup,
 156      &    ' (from residue',nstart_sup,' to residue',
 157      &    nend_sup,').'
 158       return
 159       end