source/maxlik/src_MD_T_maxlik-NEWCORR-PMF-5/molread_zs.F

   1       subroutine molread_zs(molnum)
   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.ALLPROT'
  18       include 'COMMON.TORCNSTR'
  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",molnum
  35         stop
  36       endif
  37       write(iout,*) 'molnum:',molnum,' nres=',nres
  38 C Read sequence of the protein
  39       if (iscode.gt.0) then
  40         read (inp,'(80a1)') (sequence(i)(1:1),i=1,nres)
  41       else
  42         read (inp,'(20(1x,a3))') (sequence(i),i=1,nres)
  43       endif
  44 C Convert sequence to numeric code
  45       do i=1,nres
  46         itype(i)=rescode(i,sequence(i),iscode,.true.)
  47         if (itype(i).eq.0) stop
  48       enddo
  49       write (iout,*) "Numeric code:"
  50       write (iout,'(20i4)') (itype(i),i=1,nres)
  51       do i=1,nres-1
  52 #ifdef PROCOR
  53         if (itype(i).eq.ntyp1 .or. itype(i+1).eq.ntyp1) then
  54 #else
  55         if (itype(i).eq.ntyp1) then
  56 #endif
  57           itel(i)=0
  58 #ifdef PROCOR
  59         else if (itype(i+1).ne.20) then
  60 #else
  61         else if (itype(i).ne.20) then
  62 #endif
  63           itel(i)=1
  64         else
  65           itel(i)=2
  66         endif
  67       enddo
  68       call read_bridge
  69       nnt=1
  70       nct=nres
  71       write(iout,*) 'NNT=',NNT,' NCT=',NCT
  72       if (itype(1).eq.ntyp1) nnt=2
  73       if (itype(nres).eq.ntyp1) nct=nct-1
  74 C 8/1/20 AL - splitting into chains and determining symmetry relations
  75       call setup_var
  76       call init_int_table
  77       call seq2chains(nres,itype,nchain,chain_length,chain_border)
  78       write(iout,*) "Protein",molnum," nres",nres," nchain",nchain
  79       do i=1,nchain
  80         print *,"chain",i,chain_length(i),chain_border(1,i),
  81      &    chain_border(2,i)
  82       enddo
  83       call chain_symmetry(nchain,nres,itype,chain_border,
  84      &    chain_length,npermchain,tabpermchain)
  85       call store_molinfo(molnum)
  86       if (ns_zs(molnum).gt.0) then
  87         write (iout,'(/a,i3,a)') 'The chain contains',ns_zs(molnum),
  88      &  ' disulfide-bridging cysteines.'
  89         write (iout,'(20i4)') (iss_zs(i,molnum),i=1,ns_zs(molnum))
  90         write (iout,'(/a/)') 'Pre-formed links are:'
  91         do i=1,nss_zs(1,molnum)
  92           i1=ihpb_zs(i,1,molnum)-nres_zs(molnum)
  93           i2=jhpb_zs(i,1,molnum)-nres_zs(molnum)
  94           it1=itype_zs(i1,molnum)
  95           it2=itype_zs(i2,molnum)
  96          write (iout,'(2a,i3,3a,i3,a,3f10.3)')
  97      &    restyp(it1),'(',i1,') -- ',restyp(it2),'(',i2,')',
  98      &    dhpb_zs(i,molnum),ebr,forcon_zs(i,molnum)
  99         enddo
 100       endif
 101       write (iout,*) "Protein:",molnum," leaving MOLREAD_ZS"
 102       call flush(iout)
 103       return
 104       end
 105 c-----------------------------------------------------------------------------
 106       logical function seq_comp(itypea,itypeb,length)
 107       implicit none
 108       integer length,itypea(length),itypeb(length)
 109       integer i
 110       do i=1,length
 111         if (itypea(i).ne.itypeb(i)) then
 112           seq_comp=.false.
 113           return
 114         endif
 115       enddo
 116       seq_comp=.true.
 117       return
 118       end
 119 c-----------------------------------------------------------------------------
 120       subroutine read_bridge
 121 C Read information about disulfide bridges.
 122       implicit real*8 (a-h,o-z)
 123       include 'DIMENSIONS'
 124       include 'DIMENSIONS.ZSCOPT'
 125       include 'COMMON.IOUNITS'
 126       include 'COMMON.GEO'
 127       include 'COMMON.VAR'
 128       include 'COMMON.INTERACT'
 129       include 'COMMON.NAMES'
 130       include 'COMMON.CHAIN'
 131       include 'COMMON.FFIELD'
 132       include 'COMMON.SBRIDGE'
 133       include 'COMMON.ALLPROT'
 134 C Read bridging residues.
 135       read (inp,*) ns,(iss(i),i=1,ns)
 136 c      print *,'ns=',ns
 137       write (iout,*) 'ns=',ns,' iss:',(iss(i),i=1,ns)
 138 C Check whether the specified bridging residues are cystines.
 139       do i=1,ns
 140         if (itype(iss(i)).ne.1) then
 141           write (iout,'(2a,i3,a)')
 142      &   'Do you REALLY think that the residue ',restyp(iss(i)),i,
 143      &   ' can form a disulfide bridge?!!!'
 144           write (*,'(2a,i3,a)')
 145      &   'Do you REALLY think that the residue ',restyp(iss(i)),i,
 146      &   ' can form a disulfide bridge?!!!'
 147          stop
 148         endif
 149       enddo
 150 C Read preformed bridges.
 151       if (ns.gt.0) then
 152       read (inp,*) nss,(ihpb(i),jhpb(i),i=1,nss)
 153       write (iout,*) 'nss=',nss,' ihpb,jhpb: ',(ihpb(i),jhpb(i),i=1,nss)
 154       if (nss.gt.0) then
 155         nhpb=nss
 156 C Check if the residues involved in bridges are in the specified list of
 157 C bridging residues.
 158         do i=1,nss
 159           do j=1,i-1
 160             if (ihpb(i).eq.ihpb(j).or.ihpb(i).eq.jhpb(j)
 161      &      .or.jhpb(i).eq.ihpb(j).or.jhpb(i).eq.jhpb(j)) then
 162               write (iout,'(a,i3,a)') 'Disulfide pair',i,
 163      &      ' contains residues present in other pairs.'
 164               write (*,'(a,i3,a)') 'Disulfide pair',i,
 165      &      ' contains residues present in other pairs.'
 166               stop
 167             endif
 168           enddo
 169           do j=1,ns
 170             if (ihpb(i).eq.iss(j)) goto 10
 171           enddo
 172           write (iout,'(a,i3,a)') 'Pair',i,' contains unknown cystine.'
 173    10     continue
 174           do j=1,ns
 175             if (jhpb(i).eq.iss(j)) goto 20
 176           enddo
 177           write (iout,'(a,i3,a)') 'Pair',i,' contains unknown cystine.'
 178    20     continue
 179           dhpb(i)=dbr
 180           forcon(i)=fbr
 181         enddo
 182         do i=1,nss
 183           ihpb(i)=ihpb(i)+nres
 184           jhpb(i)=jhpb(i)+nres
 185         enddo
 186       endif
 187       endif
 188       return
 189       end
 190 c------------------------------------------------------------------------------
 191       subroutine store_molinfo(molnum)
 192       implicit real*8 (a-h,o-z)
 193       include "DIMENSIONS"
 194       include 'DIMENSIONS.ZSCOPT'
 195       include "COMMON.ALLPROT"
 196       include 'COMMON.CHAIN'
 197       include 'COMMON.SBRIDGE'
 198       include "COMMON.VAR"
 199       include "COMMON.INTERACT"
 200       include "COMMON.LOCAL"
 201       include "COMMON.IOUNITS"
 202       nres_zs(molnum)=nres
 203       nnt_zs(molnum)=nnt
 204       nct_zs(molnum)=nct
 205       iatsc_s_zs(molnum)=iatsc_s
 206       iatsc_e_zs(molnum)=iatsc_e
 207       iatel_s_zs(molnum)=iatel_s
 208       iatel_e_zs(molnum)=iatel_e
 209       iturn3_start_zs(molnum)=iturn3_start
 210       iturn3_end_zs(molnum)=iturn3_end
 211       iturn4_start_zs(molnum)=iturn4_start
 212       iturn4_end_zs(molnum)=iturn4_end
 213       iatscp_s_zs(molnum)=iatscp_s
 214       iatscp_e_zs(molnum)=iatscp_e
 215       loc_start_zs(molnum)=loc_start
 216       loc_end_zs(molnum)=loc_end
 217       ithet_start_zs(molnum)=ithet_start
 218       ithet_end_zs(molnum)=ithet_end
 219       iphi_start_zs(molnum)=iphi_start
 220       iphi_end_zs(molnum)=iphi_end
 221       itau_start_zs(molnum)=itau_start
 222       itau_end_zs(molnum)=itau_end
 223       nvar_zs(molnum)=nvar
 224       nphi_zs(molnum)=nphi
 225       ntheta_zs(molnum)=ntheta
 226       nside_zs(molnum)=nside
 227       do i=1,nres
 228         ialph_zs(i,1,molnum)=ialph(i,1)
 229         ialph_zs(i,2,molnum)=ialph(i,2)
 230       enddo
 231       do i=1,nres
 232         nint_gr_zs(i,molnum)=nint_gr(i)
 233         nscp_gr_zs(i,molnum)=nscp_gr(i)
 234         do j=1,nint_gr(i)
 235           istart_zs(i,j,molnum)=istart(i,j)
 236           iend_zs(i,j,molnum)=iend(i,j)
 237         enddo
 238         itype_zs(i,molnum)=itype(i)
 239         itel_zs(i,molnum)=itel(i)
 240         ielstart_zs(i,molnum)=ielstart(i)
 241         ielend_zs(i,molnum)=ielend(i)
 242 c        write (iout,*) "i",i," nscp_gr",nscp_gr(i),
 243 c     &   " iscp",(iscpstart(i,j),iscpend(i,j),j=1,nscp_gr(i))
 244         do j=1,nscp_gr(i)
 245           iscpstart_zs(i,j,molnum)=iscpstart(i,j)
 246           iscpend_zs(i,j,molnum)=iscpend(i,j)
 247         enddo
 248       enddo
 249       ns_zs(molnum)=ns
 250 c      nss_zs(molnum)=nss
 251 c      nfree_zs(molnum)=nfree
 252 c      do i=1,ns
 253 c        iss_zs(i,molnum)=iss(i)
 254 c      enddo
 255 c      do i=1,nss
 256 c        ihpb_zs(i,molnum)=ihpb(i)
 257 c        jhpb_zs(i,molnum)=jhpb(i)
 258 c        dhpb_zs(i,molnum)=dhpb(i)
 259 c        forcon_zs(i,molnum)=forcon(i)
 260 c      enddo
 261       link_start_zs(molnum)=link_start
 262       link_end_zs(molnum)=link_end
 263       nchain_zs(molnum)=nchain
 264       npermchain_zs(molnum)=npermchain
 265       do i=1,nchain
 266         chain_border_zs(1,i,molnum)= chain_border(1,i)
 267         chain_border_zs(2,i,molnum)= chain_border(2,i)
 268         chain_length_zs(i,molnum)=chain_length(i)
 269       enddo
 270       do i=1,npermchain
 271         do j=1,nchain
 272           tabpermchain_zs(j,i,molnum)=tabpermchain(j,i)
 273         enddo
 274       enddo
 275       return
 276       end
 277 c------------------------------------------------------------------------------
 278       subroutine restore_molinfo(molnum)
 279       implicit real*8 (a-h,o-z)
 280       include "DIMENSIONS"
 281       include 'DIMENSIONS.ZSCOPT'
 282       include "COMMON.ALLPROT"
 283       include "COMMON.SBRIDGE"
 284       include "COMMON.CHAIN"
 285       include "COMMON.VAR"
 286       include "COMMON.INTERACT"
 287       include "COMMON.TORCNSTR"
 288       include "COMMON.LOCAL"
 289       ndih_constr=0
 290       nres=nres_zs(molnum)
 291       nnt=nnt_zs(molnum)
 292       nct=nct_zs(molnum)
 293       iatsc_s=iatsc_s_zs(molnum)
 294       iatsc_e=iatsc_e_zs(molnum)
 295       iatel_s=iatel_s_zs(molnum)
 296       iatel_e=iatel_e_zs(molnum)
 297       iturn3_start=iturn3_start_zs(molnum)
 298       iturn3_end=iturn3_end_zs(molnum)
 299       iturn4_start=iturn4_start_zs(molnum)
 300       iturn4_end=iturn4_end_zs(molnum)
 301       iatscp_s=iatscp_s_zs(molnum)
 302       iatscp_e=iatscp_e_zs(molnum)
 303       ithet_start=ithet_start_zs(molnum)
 304       ithet_end=ithet_end_zs(molnum)
 305       iphi_start=iphi_start_zs(molnum)
 306       iphi_end=iphi_end_zs(molnum)
 307       itau_start=itau_start_zs(molnum)
 308       itau_end=itau_end_zs(molnum)
 309       loc_start=loc_start_zs(molnum)
 310       loc_end=loc_end_zs(molnum)
 311       nvar=nvar_zs(molnum)
 312       nphi=nphi_zs(molnum)
 313       ntheta=ntheta_zs(molnum)
 314       nside=nside_zs(molnum)
 315       do i=1,nres
 316         ialph(i,1)=ialph_zs(i,1,molnum)
 317         ialph(i,2)=ialph_zs(i,2,molnum)
 318       enddo
 319       do i=1,nres
 320         nint_gr(i)=nint_gr_zs(i,molnum)
 321         nscp_gr(i)=nscp_gr_zs(i,molnum)
 322         do j=1,nint_gr(i)
 323           istart(i,j)=istart_zs(i,j,molnum)
 324           iend(i,j)=iend_zs(i,j,molnum)
 325         enddo
 326         itype(i)=itype_zs(i,molnum)
 327         itel(i)=itel_zs(i,molnum)
 328         ielstart(i)=ielstart_zs(i,molnum)
 329         ielend(i)=ielend_zs(i,molnum)
 330         do j=1,nscp_gr(i)
 331           iscpstart(i,j)=iscpstart_zs(i,j,molnum)
 332           iscpend(i,j)=iscpend_zs(i,j,molnum)
 333         enddo
 334       enddo
 335       ns=ns_zs(molnum)
 336 c      nss=nss_zs(molnum)
 337 c      nfree=nfree_zs(molnum)
 338 c      do i=1,ns
 339 c        iss(i)=iss_zs(i,molnum)
 340 c      enddo
 341 c      do i=1,nss
 342 c        ihpb(i)=ihpb_zs(i,molnum)
 343 c        jhpb(i)=jhpb_zs(i,molnum)
 344 c        dhpb(i)=dhpb_zs(i,molnum)
 345 c        forcon(i)=forcon_zs(i,molnum)
 346 c      enddo
 347       link_start=link_start_zs(molnum)
 348       link_end=link_end_zs(molnum)
 349       nchain=nchain_zs(molnum)
 350       npermchain=npermchain_zs(molnum)
 351       do i=1,nchain
 352         chain_border(1,i)=chain_border_zs(1,i,molnum)
 353         chain_border(2,i)=chain_border_zs(2,i,molnum)
 354         chain_length(i)=chain_length_zs(i,molnum)
 355       enddo
 356       do i=1,npermchain
 357         do j=1,nchain
 358           tabpermchain(j,i)=tabpermchain_zs(j,i,molnum)
 359         enddo
 360       enddo
 361       return
 362       end
 363 c------------------------------------------------------------------------------
 364       subroutine read_angles(kanal,iscor,energ,iprot,*)
 365       implicit real*8 (a-h,o-z)
 366       include 'DIMENSIONS'
 367       include 'DIMENSIONS.ZSCOPT'
 368       include 'COMMON.INTERACT'
 369       include 'COMMON.SBRIDGE'
 370       include 'COMMON.GEO'
 371       include 'COMMON.VAR'
 372       include 'COMMON.CHAIN'
 373       include 'COMMON.IOUNITS'
 374       character*80 lineh
 375       read(kanal,'(a80)',end=10,err=10) lineh
 376       read(lineh(:5),*,err=8) ic
 377       read(lineh(6:),*,err=8) energ
 378       goto 9
 379     8 ic=1
 380       print *,'error, assuming e=1d10',lineh
 381       energ=1d10
 382       nss=0
 383     9 continue
 384       read(lineh(18:),*,end=10,err=10) nss
 385       IF (NSS.LT.9) THEN
 386         read (lineh(20:),*,end=10,err=10)
 387      &  (IHPB(I),JHPB(I),I=1,NSS),iscor
 388       ELSE
 389         read (lineh(20:),*,end=10,err=10) (IHPB(I),JHPB(I),I=1,8)
 390         read (kanal,*,end=10,err=10) (IHPB(I),JHPB(I),
 391      &    I=9,NSS),iscor
 392       ENDIF
 393 c      print *,"energy",energ," iscor",iscor
 394       read (kanal,*,err=10,end=10) (theta(i),i=3,nres)
 395       read (kanal,*,err=10,end=10) (phi(i),i=4,nres)
 396       read (kanal,*,err=10,end=10) (alph(i),i=2,nres-1)
 397       read (kanal,*,err=10,end=10) (omeg(i),i=2,nres-1)
 398       do i=1,nres
 399         theta(i)=deg2rad*theta(i)
 400         phi(i)=deg2rad*phi(i)
 401         alph(i)=deg2rad*alph(i)
 402         omeg(i)=deg2rad*omeg(i)
 403       enddo
 404       return
 405    10 return1
 406       end