source/wham/src/molread_zs.F

   1       subroutine molread(*)
   2 C
   3 C Read molecular data.
   4 C
   5       implicit real*8 (a-h,o-z)
   6       include 'DIMENSIONS'
   7       include 'DIMENSIONS.ZSCOPT'
   8       include 'COMMON.IOUNITS'
   9       include 'COMMON.GEO'
  10       include 'COMMON.VAR'
  11       include 'COMMON.INTERACT'
  12       include 'COMMON.LOCAL'
  13       include 'COMMON.NAMES'
  14       include 'COMMON.CHAIN'
  15       include 'COMMON.FFIELD'
  16       include 'COMMON.SBRIDGE'
  17       include 'COMMON.TORCNSTR'
  18       include 'COMMON.CONTROL'
  19       character*4 sequence(maxres)
  20       integer rescode
  21       double precision x(maxvar)
  22       character*320 controlcard,ucase
  23       dimension itype_pdb(maxres)
  24       logical seq_comp
  25       call card_concat(controlcard,.true.)
  26       call reada(controlcard,'SCAL14',scal14,0.4d0)
  27       call reada(controlcard,'SCALSCP',scalscp,1.0d0)
  28       call reada(controlcard,'CUTOFF',cutoff_corr,7.0d0)
  29       call reada(controlcard,'DELT_CORR',delt_corr,0.5d0)
  30       r0_corr=cutoff_corr-delt_corr
  31       call readi(controlcard,"NRES",nres,0)
  32       iscode=index(controlcard,"ONE_LETTER")
  33       if (nres.le.0) then
  34         write (iout,*) "Error: no residues in molecule"
  35         return1
  36       endif
  37       if (nres.gt.maxres) then
  38         write (iout,*) "Error: too many residues",nres,maxres
  39       endif
  40       write(iout,*) 'nres=',nres
  41 C Read sequence of the protein
  42       if (iscode.gt.0) then
  43         read (inp,'(80a1)') (sequence(i)(1:1),i=1,nres)
  44       else
  45         read (inp,'(20(1x,a3))') (sequence(i),i=1,nres)
  46       endif
  47 C Convert sequence to numeric code
  48       do i=1,nres
  49         itype(i)=rescode(i,sequence(i),iscode)
  50       enddo
  51       write (iout,*) "Numeric code:"
  52       write (iout,'(20i4)') (itype(i),i=1,nres)
  53       do i=1,nres-1
  54 #ifdef PROCOR
  55         if (itype(i).eq.21 .or. itype(i+1).eq.21) then
  56 #else
  57         if (itype(i).eq.21) then
  58 #endif
  59           itel(i)=0
  60 #ifdef PROCOR
  61         else if (itype(i+1).ne.20) then
  62 #else
  63         else if (itype(i).ne.20) then
  64 #endif
  65           itel(i)=1
  66         else
  67           itel(i)=2
  68         endif
  69       enddo
  70       call read_bridge
  71
  72       if (with_dihed_constr) then
  73
  74       read (inp,*) ndih_constr
  75       if (ndih_constr.gt.0) then
  76         read (inp,*) ftors
  77         write (iout,*) 'FTORS',ftors
  78         read (inp,*) (idih_constr(i),phi0(i),drange(i),i=1,ndih_constr)
  79         write (iout,*)
  80      &   'There are',ndih_constr,' constraints on phi angles.'
  81         do i=1,ndih_constr
  82           write (iout,'(i5,2f8.3)') idih_constr(i),phi0(i),drange(i)
  83         enddo
  84         do i=1,ndih_constr
  85           phi0(i)=deg2rad*phi0(i)
  86           drange(i)=deg2rad*drange(i)
  87         enddo
  88       endif
  89
  90       endif
  91
  92       nnt=1
  93       nct=nres
  94       if (itype(1).eq.21) nnt=2
  95       if (itype(nres).eq.21) nct=nct-1
  96       write(iout,*) 'NNT=',NNT,' NCT=',NCT
  97       call setup_var
  98       call init_int_table
  99       if (ns.gt.0) then
 100         write (iout,'(/a,i3,a)') 'The chain contains',ns,
 101      &  ' disulfide-bridging cysteines.'
 102         write (iout,'(20i4)') (iss(i),i=1,ns)
 103         write (iout,'(/a/)') 'Pre-formed links are:'
 104         do i=1,nss
 105           i1=ihpb(i)-nres
 106           i2=jhpb(i)-nres
 107           it1=itype(i1)
 108           it2=itype(i2)
 109          write (iout,'(2a,i3,3a,i3,a,3f10.3)')
 110      &    restyp(it1),'(',i1,') -- ',restyp(it2),'(',i2,')',
 111      &    dhpb(i),ebr,forcon(i)
 112         enddo
 113       endif
 114       write (iout,'(a)')
 115       return
 116       end
 117 c-----------------------------------------------------------------------------
 118       logical function seq_comp(itypea,itypeb,length)
 119       implicit none
 120       integer length,itypea(length),itypeb(length)
 121       integer i
 122       do i=1,length
 123         if (itypea(i).ne.itypeb(i)) then
 124           seq_comp=.false.
 125           return
 126         endif
 127       enddo
 128       seq_comp=.true.
 129       return
 130       end
 131 c-----------------------------------------------------------------------------
 132       subroutine read_bridge
 133 C Read information about disulfide bridges.
 134       implicit real*8 (a-h,o-z)
 135       include 'DIMENSIONS'
 136       include 'DIMENSIONS.ZSCOPT'
 137       include 'COMMON.IOUNITS'
 138       include 'COMMON.GEO'
 139       include 'COMMON.VAR'
 140       include 'COMMON.INTERACT'
 141       include 'COMMON.NAMES'
 142       include 'COMMON.CHAIN'
 143       include 'COMMON.FFIELD'
 144       include 'COMMON.SBRIDGE'
 145 C Read bridging residues.
 146       read (inp,*) ns,(iss(i),i=1,ns)
 147       print *,'ns=',ns
 148       write (iout,*) 'ns=',ns,' iss:',(iss(i),i=1,ns)
 149 C Check whether the specified bridging residues are cystines.
 150       do i=1,ns
 151         if (itype(iss(i)).ne.1) then
 152           write (iout,'(2a,i3,a)')
 153      &   'Do you REALLY think that the residue ',restyp(iss(i)),i,
 154      &   ' can form a disulfide bridge?!!!'
 155           write (*,'(2a,i3,a)')
 156      &   'Do you REALLY think that the residue ',restyp(iss(i)),i,
 157      &   ' can form a disulfide bridge?!!!'
 158          stop
 159         endif
 160       enddo
 161 C Read preformed bridges.
 162       if (ns.gt.0) then
 163       read (inp,*) nss,(ihpb(i),jhpb(i),i=1,nss)
 164       write (iout,*) 'nss=',nss,' ihpb,jhpb: ',(ihpb(i),jhpb(i),i=1,nss)
 165       if (nss.gt.0) then
 166         nhpb=nss
 167 C Check if the residues involved in bridges are in the specified list of
 168 C bridging residues.
 169         do i=1,nss
 170           do j=1,i-1
 171             if (ihpb(i).eq.ihpb(j).or.ihpb(i).eq.jhpb(j)
 172      &      .or.jhpb(i).eq.ihpb(j).or.jhpb(i).eq.jhpb(j)) then
 173               write (iout,'(a,i3,a)') 'Disulfide pair',i,
 174      &      ' contains residues present in other pairs.'
 175               write (*,'(a,i3,a)') 'Disulfide pair',i,
 176      &      ' contains residues present in other pairs.'
 177               stop
 178             endif
 179           enddo
 180           do j=1,ns
 181             if (ihpb(i).eq.iss(j)) goto 10
 182           enddo
 183           write (iout,'(a,i3,a)') 'Pair',i,' contains unknown cystine.'
 184    10     continue
 185           do j=1,ns
 186             if (jhpb(i).eq.iss(j)) goto 20
 187           enddo
 188           write (iout,'(a,i3,a)') 'Pair',i,' contains unknown cystine.'
 189    20     continue
 190           dhpb(i)=dbr
 191           forcon(i)=fbr
 192         enddo
 193         do i=1,nss
 194           ihpb(i)=ihpb(i)+nres
 195           jhpb(i)=jhpb(i)+nres
 196         enddo
 197       endif
 198       endif
 199       return
 200       end
 201 c------------------------------------------------------------------------------
 202       subroutine read_angles(kanal,iscor,energ,iprot,*)
 203       implicit real*8 (a-h,o-z)
 204       include 'DIMENSIONS'
 205       include 'DIMENSIONS.ZSCOPT'
 206       include 'COMMON.INTERACT'
 207       include 'COMMON.SBRIDGE'
 208       include 'COMMON.GEO'
 209       include 'COMMON.VAR'
 210       include 'COMMON.CHAIN'
 211       include 'COMMON.IOUNITS'
 212       character*80 lineh
 213       read(kanal,'(a80)',end=10,err=10) lineh
 214       read(lineh(:5),*,err=8) ic
 215       read(lineh(6:),*,err=8) energ
 216       goto 9
 217     8 ic=1
 218       print *,'error, assuming e=1d10',lineh
 219       energ=1d10
 220       nss=0
 221     9 continue
 222       read(lineh(18:),*,end=10,err=10) nss
 223       IF (NSS.LT.9) THEN
 224         read (lineh(20:),*,end=10,err=10)
 225      &  (IHPB(I),JHPB(I),I=1,NSS),iscor
 226       ELSE
 227         read (lineh(20:),*,end=10,err=10) (IHPB(I),JHPB(I),I=1,8)
 228         read (kanal,*,end=10,err=10) (IHPB(I),JHPB(I),
 229      &    I=9,NSS),iscor
 230       ENDIF
 231 c      print *,"energy",energ," iscor",iscor
 232       read (kanal,*,err=10,end=10) (theta(i),i=3,nres)
 233       read (kanal,*,err=10,end=10) (phi(i),i=4,nres)
 234       read (kanal,*,err=10,end=10) (alph(i),i=2,nres-1)
 235       read (kanal,*,err=10,end=10) (omeg(i),i=2,nres-1)
 236       do i=1,nres
 237         theta(i)=deg2rad*theta(i)
 238         phi(i)=deg2rad*phi(i)
 239         alph(i)=deg2rad*alph(i)
 240         omeg(i)=deg2rad*omeg(i)
 241       enddo
 242       return
 243    10 return1
 244       end