--- /dev/null
+ subroutine molread(*)
+C
+C Read molecular data.
+C
+ implicit real*8 (a-h,o-z)
+ include 'DIMENSIONS'
+ include 'DIMENSIONS.ZSCOPT'
+ include 'COMMON.IOUNITS'
+ include 'COMMON.GEO'
+ include 'COMMON.VAR'
+ include 'COMMON.INTERACT'
+ include 'COMMON.LOCAL'
+ include 'COMMON.NAMES'
+ include 'COMMON.CHAIN'
+ include 'COMMON.FFIELD'
+ include 'COMMON.SBRIDGE'
+ include 'COMMON.TORCNSTR'
+ include 'COMMON.CONTROL'
+ character*4 sequence(maxres)
+ integer rescode
+ double precision x(maxvar)
+ character*320 controlcard,ucase
+ dimension itype_pdb(maxres)
+ logical seq_comp
+ call card_concat(controlcard,.true.)
+ call reada(controlcard,'SCAL14',scal14,0.4d0)
+ call reada(controlcard,'SCALSCP',scalscp,1.0d0)
+ call reada(controlcard,'CUTOFF',cutoff_corr,7.0d0)
+ call reada(controlcard,'DELT_CORR',delt_corr,0.5d0)
+ r0_corr=cutoff_corr-delt_corr
+ call readi(controlcard,"NRES",nres,0)
+ iscode=index(controlcard,"ONE_LETTER")
+ if (nres.le.0) then
+ write (iout,*) "Error: no residues in molecule"
+ return1
+ endif
+ if (nres.gt.maxres) then
+ write (iout,*) "Error: too many residues",nres,maxres
+ endif
+ write(iout,*) 'nres=',nres
+C Read sequence of the protein
+ if (iscode.gt.0) then
+ read (inp,'(80a1)') (sequence(i)(1:1),i=1,nres)
+ else
+ read (inp,'(20(1x,a3))') (sequence(i),i=1,nres)
+ endif
+C Convert sequence to numeric code
+ do i=1,nres
+ itype(i)=rescode(i,sequence(i),iscode)
+ enddo
+ write (iout,*) "Numeric code:"
+ write (iout,'(20i4)') (itype(i),i=1,nres)
+ do i=1,nres-1
+#ifdef PROCOR
+ if (itype(i).eq.21 .or. itype(i+1).eq.21) then
+#else
+ if (itype(i).eq.21) then
+#endif
+ itel(i)=0
+#ifdef PROCOR
+ else if (itype(i+1).ne.20) then
+#else
+ else if (itype(i).ne.20) then
+#endif
+ itel(i)=1
+ else
+ itel(i)=2
+ endif
+ enddo
+ call read_bridge
+
+ if (with_dihed_constr) then
+
+ read (inp,*) ndih_constr
+ if (ndih_constr.gt.0) then
+ read (inp,*) ftors
+ write (iout,*) 'FTORS',ftors
+ read (inp,*) (idih_constr(i),phi0(i),drange(i),i=1,ndih_constr)
+ write (iout,*)
+ & 'There are',ndih_constr,' constraints on phi angles.'
+ do i=1,ndih_constr
+ write (iout,'(i5,2f8.3)') idih_constr(i),phi0(i),drange(i)
+ enddo
+ do i=1,ndih_constr
+ phi0(i)=deg2rad*phi0(i)
+ drange(i)=deg2rad*drange(i)
+ enddo
+ endif
+
+ endif
+
+ nnt=1
+ nct=nres
+ if (itype(1).eq.21) nnt=2
+ if (itype(nres).eq.21) nct=nct-1
+ write(iout,*) 'NNT=',NNT,' NCT=',NCT
+c Read distance restraints
+ if (constr_dist.gt.0) then
+ if (refstr) call read_ref_structure(*11)
+ call read_dist_constr
+ call hpb_partition
+ endif
+
+ call setup_var
+ call init_int_table
+ if (ns.gt.0) then
+ write (iout,'(/a,i3,a)') 'The chain contains',ns,
+ & ' disulfide-bridging cysteines.'
+ write (iout,'(20i4)') (iss(i),i=1,ns)
+ write (iout,'(/a/)') 'Pre-formed links are:'
+ do i=1,nss
+ i1=ihpb(i)-nres
+ i2=jhpb(i)-nres
+ it1=itype(i1)
+ it2=itype(i2)
+ write (iout,'(2a,i3,3a,i3,a,3f10.3)')
+ & restyp(it1),'(',i1,') -- ',restyp(it2),'(',i2,')',
+ & dhpb(i),ebr,forcon(i)
+ enddo
+ endif
+ write (iout,'(a)')
+ return
+ 11 stop "Error reading reference structure"
+ end
+c-----------------------------------------------------------------------------
+ logical function seq_comp(itypea,itypeb,length)
+ implicit none
+ integer length,itypea(length),itypeb(length)
+ integer i
+ do i=1,length
+ if (itypea(i).ne.itypeb(i)) then
+ seq_comp=.false.
+ return
+ endif
+ enddo
+ seq_comp=.true.
+ return
+ end
+c-----------------------------------------------------------------------------
+ subroutine read_bridge
+C Read information about disulfide bridges.
+ implicit real*8 (a-h,o-z)
+ include 'DIMENSIONS'
+ include 'DIMENSIONS.ZSCOPT'
+ include 'COMMON.IOUNITS'
+ include 'COMMON.GEO'
+ include 'COMMON.VAR'
+ include 'COMMON.INTERACT'
+ include 'COMMON.NAMES'
+ include 'COMMON.CHAIN'
+ include 'COMMON.FFIELD'
+ include 'COMMON.SBRIDGE'
+C Read bridging residues.
+ read (inp,*) ns,(iss(i),i=1,ns)
+ print *,'ns=',ns
+ write (iout,*) 'ns=',ns,' iss:',(iss(i),i=1,ns)
+C Check whether the specified bridging residues are cystines.
+ do i=1,ns
+ if (itype(iss(i)).ne.1) then
+ write (iout,'(2a,i3,a)')
+ & 'Do you REALLY think that the residue ',restyp(iss(i)),i,
+ & ' can form a disulfide bridge?!!!'
+ write (*,'(2a,i3,a)')
+ & 'Do you REALLY think that the residue ',restyp(iss(i)),i,
+ & ' can form a disulfide bridge?!!!'
+ stop
+ endif
+ enddo
+C Read preformed bridges.
+ if (ns.gt.0) then
+ read (inp,*) nss,(ihpb(i),jhpb(i),i=1,nss)
+ write (iout,*) 'nss=',nss,' ihpb,jhpb: ',(ihpb(i),jhpb(i),i=1,nss)
+ if (nss.gt.0) then
+ nhpb=nss
+C Check if the residues involved in bridges are in the specified list of
+C bridging residues.
+ do i=1,nss
+ do j=1,i-1
+ if (ihpb(i).eq.ihpb(j).or.ihpb(i).eq.jhpb(j)
+ & .or.jhpb(i).eq.ihpb(j).or.jhpb(i).eq.jhpb(j)) then
+ write (iout,'(a,i3,a)') 'Disulfide pair',i,
+ & ' contains residues present in other pairs.'
+ write (*,'(a,i3,a)') 'Disulfide pair',i,
+ & ' contains residues present in other pairs.'
+ stop
+ endif
+ enddo
+ do j=1,ns
+ if (ihpb(i).eq.iss(j)) goto 10
+ enddo
+ write (iout,'(a,i3,a)') 'Pair',i,' contains unknown cystine.'
+ 10 continue
+ do j=1,ns
+ if (jhpb(i).eq.iss(j)) goto 20
+ enddo
+ write (iout,'(a,i3,a)') 'Pair',i,' contains unknown cystine.'
+ 20 continue
+ dhpb(i)=dbr
+ forcon(i)=fbr
+ enddo
+ do i=1,nss
+ ihpb(i)=ihpb(i)+nres
+ jhpb(i)=jhpb(i)+nres
+ enddo
+ endif
+ endif
+ return
+ end
+c------------------------------------------------------------------------------
+ subroutine read_angles(kanal,iscor,energ,iprot,*)
+ implicit real*8 (a-h,o-z)
+ include 'DIMENSIONS'
+ include 'DIMENSIONS.ZSCOPT'
+ include 'COMMON.INTERACT'
+ include 'COMMON.SBRIDGE'
+ include 'COMMON.GEO'
+ include 'COMMON.VAR'
+ include 'COMMON.CHAIN'
+ include 'COMMON.IOUNITS'
+ character*80 lineh
+ read(kanal,'(a80)',end=10,err=10) lineh
+ read(lineh(:5),*,err=8) ic
+ read(lineh(6:),*,err=8) energ
+ goto 9
+ 8 ic=1
+ print *,'error, assuming e=1d10',lineh
+ energ=1d10
+ nss=0
+ 9 continue
+ read(lineh(18:),*,end=10,err=10) nss
+ IF (NSS.LT.9) THEN
+ read (lineh(20:),*,end=10,err=10)
+ & (IHPB(I),JHPB(I),I=1,NSS),iscor
+ ELSE
+ read (lineh(20:),*,end=10,err=10) (IHPB(I),JHPB(I),I=1,8)
+ read (kanal,*,end=10,err=10) (IHPB(I),JHPB(I),
+ & I=9,NSS),iscor
+ ENDIF
+c print *,"energy",energ," iscor",iscor
+ read (kanal,*,err=10,end=10) (theta(i),i=3,nres)
+ read (kanal,*,err=10,end=10) (phi(i),i=4,nres)
+ read (kanal,*,err=10,end=10) (alph(i),i=2,nres-1)
+ read (kanal,*,err=10,end=10) (omeg(i),i=2,nres-1)
+ do i=1,nres
+ theta(i)=deg2rad*theta(i)
+ phi(i)=deg2rad*phi(i)
+ alph(i)=deg2rad*alph(i)
+ omeg(i)=deg2rad*omeg(i)
+ enddo
+ return
+ 10 return1
+ end
+c-------------------------------------------------------------------------------
+ subroutine read_dist_constr
+ implicit real*8 (a-h,o-z)
+ include 'DIMENSIONS'
+ include 'COMMON.CONTROL'
+ include 'COMMON.CHAIN'
+ include 'COMMON.IOUNITS'
+ include 'COMMON.SBRIDGE'
+ integer ifrag_(2,100),ipair_(2,100)
+ double precision wfrag_(100),wpair_(100)
+ character*500 controlcard
+c write (iout,*) "Calling read_dist_constr"
+c write (iout,*) "nres",nres," nstart_sup",nstart_sup," nsup",nsup
+c call flush(iout)
+ call card_concat(controlcard,.true.)
+ call readi(controlcard,"NFRAG",nfrag_,0)
+ call readi(controlcard,"NPAIR",npair_,0)
+ call readi(controlcard,"NDIST",ndist_,0)
+ call reada(controlcard,'DIST_CUT',dist_cut,5.0d0)
+ call multreadi(controlcard,"IFRAG",ifrag_(1,1),2*nfrag_,0)
+ call multreadi(controlcard,"IPAIR",ipair_(1,1),2*npair_,0)
+ call multreada(controlcard,"WFRAG",wfrag_(1),nfrag_,0.0d0)
+ call multreada(controlcard,"WPAIR",wpair_(1),npair_,0.0d0)
+ write (iout,*) "NFRAG",nfrag_," NPAIR",npair_," NDIST",ndist_
+ write (iout,*) "IFRAG"
+ do i=1,nfrag_
+ write (iout,*) i,ifrag_(1,i),ifrag_(2,i),wfrag_(i)
+ enddo
+ write (iout,*) "IPAIR"
+ do i=1,npair_
+ write (iout,*) i,ipair_(1,i),ipair_(2,i),wpair_(i)
+ enddo
+ call flush(iout)
+ if (.not.refstr .and. nfrag_.gt.0) then
+ write (iout,*)
+ & "ERROR: no reference structure to compute distance restraints"
+ write (iout,*)
+ & "Restraints must be specified explicitly (NDIST=number)"
+ stop
+ endif
+ if (nfrag_.lt.2 .and. npair_.gt.0) then
+ write (iout,*) "ERROR: Less than 2 fragments specified",
+ & " but distance restraints between pairs requested"
+ stop
+ endif
+ call flush(iout)
+ do i=1,nfrag_
+ if (ifrag_(1,i).lt.nstart_sup) ifrag_(1,i)=nstart_sup
+ if (ifrag_(2,i).gt.nstart_sup+nsup-1)
+ & ifrag_(2,i)=nstart_sup+nsup-1
+c write (iout,*) i,ifrag_(1,i),ifrag_(2,i),wfrag_(i)
+ call flush(iout)
+ if (wfrag_(i).gt.0.0d0) then
+ do j=ifrag_(1,i),ifrag_(2,i)-1
+ do k=j+1,ifrag_(2,i)
+ write (iout,*) "j",j," k",k
+ ddjk=dist(j,k)
+ if (constr_dist.eq.1) then
+ nhpb=nhpb+1
+ ihpb(nhpb)=j
+ jhpb(nhpb)=k
+ dhpb(nhpb)=ddjk
+ forcon(nhpb)=wfrag_(i)
+ else if (constr_dist.eq.2) then
+ if (ddjk.le.dist_cut) then
+ nhpb=nhpb+1
+ ihpb(nhpb)=j
+ jhpb(nhpb)=k
+ dhpb(nhpb)=ddjk
+ forcon(nhpb)=wfrag_(i)
+ endif
+ else
+ nhpb=nhpb+1
+ ihpb(nhpb)=j
+ jhpb(nhpb)=k
+ dhpb(nhpb)=ddjk
+ forcon(nhpb)=wfrag_(i)*dexp(-0.5d0*(ddjk/dist_cut)**2)
+ endif
+ write (iout,'(a,3i5,f8.2,1pe12.2)') "+dist.constr ",
+ & nhpb,ihpb(nhpb),jhpb(nhpb),dhpb(nhpb),forcon(nhpb)
+ enddo
+ enddo
+ endif
+ enddo
+ do i=1,npair_
+ if (wpair_(i).gt.0.0d0) then
+ ii = ipair_(1,i)
+ jj = ipair_(2,i)
+ if (ii.gt.jj) then
+ itemp=ii
+ ii=jj
+ jj=itemp
+ endif
+ do j=ifrag_(1,ii),ifrag_(2,ii)
+ do k=ifrag_(1,jj),ifrag_(2,jj)
+ nhpb=nhpb+1
+ ihpb(nhpb)=j
+ jhpb(nhpb)=k
+ forcon(nhpb)=wpair_(i)
+ dhpb(nhpb)=dist(j,k)
+ write (iout,'(a,3i5,f8.2,f10.1)') "+dist.constr ",
+ & nhpb,ihpb(nhpb),jhpb(nhpb),dhpb(nhpb),forcon(nhpb)
+ enddo
+ enddo
+ endif
+ enddo
+ do i=1,ndist_
+ read (inp,*) ihpb(nhpb+1),jhpb(nhpb+1),dhpb(i),dhpb1(i),
+ & ibecarb(i),forcon(nhpb+1)
+ if (forcon(nhpb+1).gt.0.0d0) then
+ nhpb=nhpb+1
+ if (ibecarb(i).gt.0) then
+ ihpb(i)=ihpb(i)+nres
+ jhpb(i)=jhpb(i)+nres
+ endif
+ if (dhpb(nhpb).eq.0.0d0)
+ & dhpb(nhpb)=dist(ihpb(nhpb),jhpb(nhpb))
+ endif
+ enddo
+ do i=1,nhpb
+ write (iout,'(a,3i5,2f8.2,i2,f10.1)') "+dist.constr ",
+ & i,ihpb(i),jhpb(i),dhpb(i),dhpb1(i),ibecarb(i),forcon(i)
+ enddo
+ call flush(iout)
+ return
+ end
projects / unres.git / blobdiff