--- /dev/null
+ subroutine readpdb
+C Read the PDB file and convert the peptide geometry into virtual-chain
+C geometry.
+ implicit none
+ include 'DIMENSIONS'
+ include 'DIMENSIONS.ZSCOPT'
+ include 'COMMON.CONTROL'
+ include 'COMMON.LOCAL'
+ include 'COMMON.VAR'
+ include 'COMMON.CHAIN'
+ include 'COMMON.INTERACT'
+ include 'COMMON.IOUNITS'
+ include 'COMMON.GEO'
+ include 'COMMON.NAMES'
+ character*3 seq,atom,res
+ character*80 card
+ double precision sccor(3,20)
+ integer i,j,iii,ibeg,ishift,ishift1,ity,ires,ires_old
+ double precision dcj
+ integer rescode
+ ibeg=1
+ ishift1=0
+ do i=1,10000
+ read (ipdbin,'(a80)',end=10) card
+ if (card(:3).eq.'END' .or. card(:3).eq.'TER') goto 10
+C Fish out the ATOM cards.
+ if (index(card(1:4),'ATOM').gt.0) then
+ read (card(14:16),'(a3)') atom
+ if (atom.eq.'CA' .or. atom.eq.'CH3') then
+C Calculate the CM of the preceding residue.
+ if (ibeg.eq.0) call sccenter(ires,iii,sccor)
+C Start new residue.
+ ires_old=ires+ishift-ishift1
+ read (card(23:26),*) ires
+c print *,"ires_old",ires_old," ires",ires
+ if (card(27:27).eq."A" .or. card(27:27).eq."B") then
+c ishift1=ishift1+1
+ endif
+ read (card(18:20),'(a3)') res
+ if (ibeg.eq.1) then
+ ishift=ires-1
+ if (res.ne.'GLY' .and. res.ne. 'ACE') then
+ ishift=ishift-1
+ itype(1)=21
+ endif
+ ibeg=0
+ else
+ ishift=ishift+ires-ires_old-1
+ endif
+ ires=ires-ishift+ishift1
+ if (res.eq.'ACE') then
+ ity=10
+ else
+ itype(ires)=rescode(ires,res,0)
+ endif
+ read(card(31:54),'(3f8.3)') (c(j,ires),j=1,3)
+ write (iout,'(2i3,2x,a,3f8.3)')
+ & ires,itype(ires),res,(c(j,ires),j=1,3)
+ iii=1
+ do j=1,3
+ sccor(j,iii)=c(j,ires)
+ enddo
+c write (*,*) card(23:27),ires,itype(ires)
+ else if (atom.ne.'O '.and.atom(1:1).ne.'H' .and.
+ & atom.ne.'N ' .and. atom.ne.'C ') then
+ iii=iii+1
+ read(card(31:54),'(3f8.3)') (sccor(j,iii),j=1,3)
+ endif
+ endif
+ enddo
+ 10 write (iout,'(a,i5)') ' Nres: ',ires
+C Calculate the CM of the last side chain.
+ call sccenter(ires,iii,sccor)
+ nres=ires
+ nsup=nres
+ nstart_sup=1
+ if (itype(nres).ne.10) then
+ nres=nres+1
+ itype(nres)=21
+ do j=1,3
+ dcj=c(j,nres-2)-c(j,nres-3)
+ c(j,nres)=c(j,nres-1)+dcj
+ c(j,2*nres)=c(j,nres)
+ enddo
+ endif
+ do i=2,nres-1
+ do j=1,3
+ c(j,i+nres)=dc(j,i)
+ enddo
+ enddo
+ do j=1,3
+ c(j,nres+1)=c(j,1)
+ c(j,2*nres)=c(j,nres)
+ enddo
+ if (itype(1).eq.21) then
+ nsup=nsup-1
+ nstart_sup=2
+ do j=1,3
+ dcj=c(j,4)-c(j,3)
+ c(j,1)=c(j,2)-dcj
+ c(j,nres+1)=c(j,1)
+ enddo
+ endif
+C Copy the coordinates to reference coordinates
+ do i=1,2*nres
+ do j=1,3
+ cref(j,i)=c(j,i)
+ enddo
+ enddo
+C Calculate internal coordinates.
+ do ires=1,nres
+ write (iout,'(2i3,2x,a,3f8.3,5x,3f8.3)')
+ & ires,itype(ires),restyp(itype(ires)),(c(j,ires),j=1,3),
+ & (c(j,ires+nres),j=1,3)
+ enddo
+ call flush(iout)
+ call int_from_cart(.true.,.true.)
+ do i=1,nres
+ phi_ref(i)=phi(i)
+ theta_ref(i)=theta(i)
+ alph_ref(i)=alph(i)
+ omeg_ref(i)=omeg(i)
+ enddo
+ ishift_pdb=ishift
+ return
+ end
+c---------------------------------------------------------------------------
+ subroutine int_from_cart(lside,lprn)
+ implicit none
+ include 'DIMENSIONS'
+ include 'DIMENSIONS.ZSCOPT'
+ include 'COMMON.LOCAL'
+ include 'COMMON.VAR'
+ include 'COMMON.CHAIN'
+ include 'COMMON.INTERACT'
+ include 'COMMON.IOUNITS'
+ include 'COMMON.GEO'
+ include 'COMMON.NAMES'
+ character*3 seq,atom,res
+ character*80 card
+ double precision sccor(3,20)
+ integer rescode
+ double precision dist,alpha,beta,di
+ integer i,j,iti
+ logical lside,lprn
+ if (lprn) then
+ write (iout,'(/a)')
+ & 'Internal coordinates calculated from crystal structure.'
+ if (lside) then
+ write (iout,'(8a)') ' Res ',' dvb',' Theta',
+ & ' Phi',' Dsc_id',' Dsc',' Alpha',
+ & ' Omega'
+ else
+ write (iout,'(4a)') ' Res ',' dvb',' Theta',
+ & ' Phi'
+ endif
+ endif
+ do i=2,nres
+ iti=itype(i)
+ write (iout,*) i,i-1,(c(j,i),j=1,3),(c(j,i-1),j=1,3),dist(i,i-1)
+ if (itype(i-1).ne.21 .and. itype(i).ne.21 .and.
+ & (dist(i,i-1).lt.2.0D0 .or. dist(i,i-1).gt.5.0D0)) then
+ write (iout,'(a,i4)') 'Bad Cartesians for residue',i
+ stop
+ endif
+ theta(i+1)=alpha(i-1,i,i+1)
+ if (i.gt.2) phi(i+1)=beta(i-2,i-1,i,i+1)
+ enddo
+ if (itype(1).eq.21) then
+ do j=1,3
+ c(j,1)=c(j,2)+(c(j,3)-c(j,4))
+ enddo
+ endif
+ if (itype(nres).eq.21) then
+ do j=1,3
+ c(j,nres)=c(j,nres-1)+(c(j,nres-2)-c(j,nres-3))
+ enddo
+ endif
+ if (lside) then
+ do i=2,nres-1
+ do j=1,3
+ c(j,maxres2)=0.5D0*(c(j,i-1)+c(j,i+1))
+ enddo
+ iti=itype(i)
+ di=dist(i,nres+i)
+ if (iti.ne.10) then
+ alph(i)=alpha(nres+i,i,maxres2)
+ omeg(i)=beta(nres+i,i,maxres2,i+1)
+ endif
+ if (lprn)
+ & write (iout,'(a3,i4,7f10.3)') restyp(iti),i,dist(i,i-1),
+ & rad2deg*theta(i),rad2deg*phi(i),dsc(iti),di,
+ & rad2deg*alph(i),rad2deg*omeg(i)
+ enddo
+ else if (lprn) then
+ do i=2,nres
+ iti=itype(i)
+ write (iout,'(a3,i4,7f10.3)') restyp(iti),i,dist(i,i-1),
+ & rad2deg*theta(i),rad2deg*phi(i)
+ enddo
+ endif
+ return
+ end
+c---------------------------------------------------------------------------
+ subroutine sccenter(ires,nscat,sccor)
+ implicit none
+ include 'DIMENSIONS'
+ include 'COMMON.CHAIN'
+ integer ires,nscat,i,j
+ double precision sccor(3,20),sccmj
+ do j=1,3
+ sccmj=0.0D0
+ do i=1,nscat
+ sccmj=sccmj+sccor(j,i)
+ enddo
+ dc(j,ires)=sccmj/nscat
+ enddo
+ return
+ end
projects / unres.git / blobdiff