source/maxlik/src_FPy.org/readrtns_compar.F.org

   1       subroutine read_ref_structure(iprot,*)
   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 'DIMENSIONS.ZSCOPT'
   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       include 'COMMON.COMPAR'
  22       include 'COMMON.VMCPAR'
  23       character*4 sequence(maxres)
  24       integer rescode
  25       double precision x(maxvar)
  26       integer itype_pdb(maxres)
  27       logical seq_comp
  28       integer i,j,k,ib,iref,iprot,nres_pdb
  29       double precision ddsc,dist
  30       integer ilen
  31       external ilen
  32 C
  33       write (2,*) "PDBREF ",pdbref
  34       nres0=nres
  35       DO IB=1,NBETA(1,iprot)
  36       if (pdbref) then
  37         iref=0
  38         write (iout,'(2a,1h.)') 'PDB data will be read from file ',
  39      &    pdbfile(ib,iprot)(:ilen(pdbfile(ib,iprot)))
  40         open(ipdbin,file=pdbfile(ib,iprot),status='old',err=33)
  41         goto 34
  42   33    write (iout,'(a)') 'Error opening PDB file.'
  43         return1
  44   34    continue
  45         nres_pdb=nres
  46         DO IREF=1,MAXREF
  47         write (2,*) "Reading reference structure",iref
  48         call readpdb(print_refstr,iprot)
  49         write (2,*) "After readpdb"
  50         call flush(iout)
  51         if (nres_pdb.eq.0) then
  52           exit
  53         endif
  54 C Copy the coordinates to reference coordinates
  55         do i=1,2*nres
  56           do j=1,3
  57             cref(j,i,iref,ib,iprot)=c(j,i)
  58           enddo
  59         enddo
  60         do i=1,nres
  61           itype_pdb(i)=itype(i)
  62         enddo
  63         nres_pdb=nres
  64         nres=nres0
  65         nstart_seq=nnt
  66         if (nsup(iprot).le.(nct-nnt+1)) then
  67           do i=0,nct-nnt+1-nsup(iprot)
  68             if (seq_comp(itype(nnt+i),itype_pdb(nstart_sup(iprot)),
  69      &        nsup(iprot))) then
  70               do j=nnt+nsup(iprot)-1,nnt,-1
  71                 do k=1,3
  72                   cref(k,j+i,iref,ib,iprot)=cref(k,j,iref,ib,iprot)
  73                   cref(k,nres+j+i,iref,ib,iprot)=
  74      &                             cref(k,nres_pdb+j,iref,ib,iprot)
  75                 enddo
  76 c                phi_ref(j+i,iprot)=phi_ref(j,iprot)
  77 c                theta_ref(j+i,iprot)=theta_ref(j,iprot)
  78 c                alph_ref(j+i,iprot)=alph_ref(j,iprot)
  79 c                omeg_ref(j+i,iprot)=omeg_ref(j,iprot)
  80               enddo
  81 #define DEBUG
  82 #ifdef DEBUG
  83               do j=nnt,nct
  84                 write (iout,'(i5,3f10.5,5x,3f10.5)')
  85      &            j,(cref(k,j,iref,ib,iprot),k=1,3),
  86      &              (cref(k,j+nres,iref,ib,iprot),k=1,3)
  87               enddo
  88 #endif
  89 #undef DEBUG
  90               nstart_seq=nnt+i
  91               nstart_sup(iprot)=nnt+i
  92               goto 111
  93             endif
  94           enddo
  95           write (iout,'(a)')
  96      &            'Error - sequences to be superposed do not match.'
  97           return1
  98         else
  99           do i=0,nsup(iprot)-(nct-nnt+1)
 100             if (seq_comp(itype(nnt),itype_pdb(nstart_sup(iprot)+i),
 101      &        nct-nnt+1))
 102      &      then
 103               nstart_sup(iprot)=nstart_sup(iprot)+i
 104               nsup(iprot)=nct-nnt+1
 105               goto 111
 106             endif
 107           enddo
 108           write (iout,'(a)')
 109      &            'Error - sequences to be superposed do not match.'
 110         endif
 111   111   continue
 112         write (iout,'(a,i5)')
 113      &   'Experimental structure begins at residue',nstart_seq
 114         ENDDO ! IREF
 115         nref(ib,iprot)=iref-1
 116         write (iout,*) "NREF",ib,iprot,nref(ib,iprot)
 117         call flush(iout)
 118       else
 119         DO IREF=1,MAXREF
 120         call read_angles(inp,*38)
 121         goto 39
 122    38   nref(ib,iprot)=iref-1
 123         exit
 124    39   call chainbuild
 125         nstart_sup(iprot)=nnt
 126         nstart_seq=nnt
 127         nsup(iprot)=nct-nnt+1
 128         do i=1,2*nres
 129           do j=1,3
 130             cref(j,i,iref,ib,iprot)=c(j,i)
 131           enddo
 132         enddo
 133         ENDDO ! IREF
 134         nref(ib,iprot)=iref-1
 135       endif
 136       nend_sup(iprot)=nstart_sup(iprot)+nsup(iprot)-1
 137 c      do i=1,2*nres
 138 c        do j=1,3
 139 c         c(j,i)=cref(j,i,iprot)
 140 c        enddo
 141 c      enddo
 142       do i=1,nres
 143 c        do j=1,3
 144 c          dc(j,nres+i)=cref(j,nres+i,iprot)-cref(j,i,iprot)
 145 c        enddo
 146 c        if (itype(i).ne.10) then
 147 c          ddsc = dist(i,nres+i)
 148 c          do j=1,3
 149 c            dc_norm(j,nres+i)=dc(j,nres+i)/ddsc
 150 c          enddo
 151 c        else
 152 c          do j=1,3
 153 c            dc_norm(j,nres+i)=0.0d0
 154 c          enddo
 155 c        endif
 156 c        write (iout,*) "i",i," dc_norm",(dc_norm(k,nres+i),k=1,3),
 157 c         " norm",dc_norm(1,nres+i)**2+dc_norm(2,nres+i)**2+
 158 c         dc_norm(3,nres+i)**2
 159         do j=1,3
 160           dc(j,i)=c(j,i+1)-c(j,i)
 161         enddo
 162         ddsc = dist(i,i+1)
 163         do j=1,3
 164           dc_norm(j,i)=dc(j,i)/ddsc
 165         enddo
 166       enddo
 167 c      print *,"Calling contact"
 168 c      call contact(print_contact,ncont_ref(iprot),icont_ref(1,1,iprot),
 169 c     &  nstart_sup(iprot),nend_sup(iprot))
 170 c      print *,"Calling elecont"
 171 c      call elecont(print_contact,ncont_pept_ref(iprot),
 172 c     &   icont_pept_ref(1,1,iprot),
 173 c     &   nstart_sup(iprot),nend_sup(iprot))
 174 c       write (iout,'(a,i3,a,i3,a,i3,a)')
 175 c     &    'Number of residues to be superposed:',nsup(iprot),
 176 c     &    ' (from residue',nstart_sup(iprot),' to residue',
 177 c     &    nend_sup(iprot),').'
 178       ENDDO ! IB
 179       close (ipdbin)
 180       return
 181       end