X-Git-Url: http://mmka.chem.univ.gda.pl/gitweb/?p=unres.git;a=blobdiff_plain;f=source%2Funres%2Fsrc_CSA_DiL%2Freadpdb.F;fp=source%2Funres%2Fsrc_CSA_DiL%2Freadpdb.F;h=0000000000000000000000000000000000000000;hp=eb4ba3f5fddce94692a64359763fd292c7768ff7;hb=de4bc5453ea46e111d936cb85e1758ed21c08fcd;hpb=b75425747e3e2b448ca5e0ef8367712e1f339124 diff --git a/source/unres/src_CSA_DiL/readpdb.F b/source/unres/src_CSA_DiL/readpdb.F deleted file mode 100644 index eb4ba3f..0000000 --- a/source/unres/src_CSA_DiL/readpdb.F +++ /dev/null @@ -1,428 +0,0 @@ - subroutine readpdb -C Read the PDB file and convert the peptide geometry into virtual-chain -C geometry. - implicit real*8 (a-h,o-z) - include 'DIMENSIONS' - include 'COMMON.LOCAL' - include 'COMMON.VAR' - include 'COMMON.CHAIN' - include 'COMMON.INTERACT' - include 'COMMON.IOUNITS' - include 'COMMON.GEO' - include 'COMMON.NAMES' - include 'COMMON.CONTROL' - include 'COMMON.DISTFIT' - include 'COMMON.SETUP' - character*3 seq,atom,res - character*80 card - dimension sccor(3,20) - double precision e1(3),e2(3),e3(3) - logical fail - integer rescode - ibeg=1 - lsecondary=.false. - nhfrag=0 - nbfrag=0 - do i=1,10000 - read (ipdbin,'(a80)',end=10) card - if (card(:5).eq.'HELIX') then - nhfrag=nhfrag+1 - lsecondary=.true. - read(card(22:25),*) hfrag(1,nhfrag) - read(card(34:37),*) hfrag(2,nhfrag) - endif - if (card(:5).eq.'SHEET') then - nbfrag=nbfrag+1 - lsecondary=.true. - read(card(24:26),*) bfrag(1,nbfrag) - read(card(35:37),*) bfrag(2,nbfrag) -crc---------------------------------------- -crc to be corrected !!! - bfrag(3,nbfrag)=bfrag(1,nbfrag) - bfrag(4,nbfrag)=bfrag(2,nbfrag) -crc---------------------------------------- - endif - 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) then - if (unres_pdb) then - do j=1,3 - dc(j,ires+nres)=sccor(j,iii) - enddo - else - call sccenter(ires,iii,sccor) - endif - endif -C Start new residue. - read (card(24:26),*) ires - 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 - endif - ires=ires-ishift - 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) -c if(me.eq.king.or..not.out1file) -c & write (iout,'(2i3,2x,a,3f8.3)') -c & ires,itype(ires),res,(c(j,ires),j=1,3) - iii=1 - do j=1,3 - sccor(j,iii)=c(j,ires) - enddo - 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 if(me.eq.king.or..not.out1file) - & write (iout,'(a,i5)') ' Nres: ',ires -C Calculate the CM of the last side chain. - if (unres_pdb) then - do j=1,3 - dc(j,ires+nres)=sccor(j,iii) - enddo - else - call sccenter(ires,iii,sccor) - endif - nres=ires - nsup=nres - nstart_sup=1 - if (itype(nres).ne.10) then - nres=nres+1 - itype(nres)=21 - if (unres_pdb) then -C 2/15/2013 by Adam: corrected insertion of the last dummy residue - call refsys(nres-3,nres-2,nres-1,e1,e2,e3,fail) - if (fail) then - e2(1)=0.0d0 - e2(2)=1.0d0 - e2(3)=0.0d0 - endif - do j=1,3 - c(j,nres)=c(j,nres-1)-3.8d0*e2(j) - enddo - else - 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 - 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 - if (unres_pdb) then -C 2/15/2013 by Adam: corrected insertion of the first dummy residue - call refsys(2,3,4,e1,e2,e3,fail) - if (fail) then - e2(1)=0.0d0 - e2(2)=1.0d0 - e2(3)=0.0d0 - endif - do j=1,3 - c(j,1)=c(j,2)-3.8d0*e2(j) - enddo - else - 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 - endif -C Calculate internal coordinates. - if(me.eq.king.or..not.out1file)then - write (iout,'(a)') - & "Backbone and SC coordinates as read from the PDB" - 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,nres+ires),j=1,3) - enddo - endif - call int_from_cart(.true.,.false.) - call sc_loc_geom(.false.) - do i=1,nres - thetaref(i)=theta(i) - phiref(i)=phi(i) - enddo - do i=1,nres-1 - do j=1,3 - dc(j,i)=c(j,i+1)-c(j,i) - dc_norm(j,i)=dc(j,i)*vbld_inv(i+1) - enddo - enddo - do i=2,nres-1 - do j=1,3 - dc(j,i+nres)=c(j,i+nres)-c(j,i) - dc_norm(j,i+nres)=dc(j,i+nres)*vbld_inv(i+nres) - enddo -c write (iout,*) i,(dc(j,i+nres),j=1,3),(dc_norm(j,i+nres),j=1,3), -c & vbld_inv(i+nres) - enddo -c call chainbuild -C Copy the coordinates to reference coordinates - do i=1,2*nres - do j=1,3 - cref(j,i)=c(j,i) - enddo - enddo - - - do j=1,nbfrag - do i=1,4 - bfrag(i,j)=bfrag(i,j)-ishift - enddo - enddo - - do j=1,nhfrag - do i=1,2 - hfrag(i,j)=hfrag(i,j)-ishift - enddo - enddo - - return - end -c--------------------------------------------------------------------------- - subroutine int_from_cart(lside,lprn) - implicit real*8 (a-h,o-z) - include 'DIMENSIONS' -c include "mpif.h" - include 'COMMON.LOCAL' - include 'COMMON.VAR' - include 'COMMON.CHAIN' - include 'COMMON.INTERACT' - include 'COMMON.IOUNITS' - include 'COMMON.GEO' - include 'COMMON.NAMES' - include 'COMMON.CONTROL' - include 'COMMON.SETUP' - character*3 seq,atom,res - character*80 card - dimension sccor(3,20) - integer rescode - logical lside,lprn - if(me.eq.king.or..not.out1file)then - if (lprn) then - write (iout,'(/a)') - & 'Internal coordinates calculated from crystal structure.' - if (lside) then - write (iout,'(8a)') ' Res ',' dvb',' Theta', - & ' Gamma',' Dsc_id',' Dsc',' Alpha', - & ' Beta ' - else - write (iout,'(4a)') ' Res ',' dvb',' Theta', - & ' Gamma' - endif - endif - endif - do i=1,nres-1 - iti=itype(i) - if (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 -ctest stop - endif - vbld(i+1)=dist(i,i+1) - vbld_inv(i+1)=1.0d0/vbld(i+1) - if (i.gt.1) 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 -c if (unres_pdb) then -c if (itype(1).eq.21) then -c theta(3)=90.0d0*deg2rad -c phi(4)=180.0d0*deg2rad -c vbld(2)=3.8d0 -c vbld_inv(2)=1.0d0/vbld(2) -c endif -c if (itype(nres).eq.21) then -c theta(nres)=90.0d0*deg2rad -c phi(nres)=180.0d0*deg2rad -c vbld(nres)=3.8d0 -c vbld_inv(nres)=1.0d0/vbld(2) -c endif -c endif - if (lside) then - do i=2,nres-1 - do j=1,3 - c(j,maxres2)=0.5D0*(2*c(j,i)+(c(j,i-1)-c(j,i))*vbld_inv(i) - & +(c(j,i+1)-c(j,i))*vbld_inv(i+1)) - enddo - iti=itype(i) - di=dist(i,nres+i) -C 9/29/12 Adam: Correction for zero SC-SC bond length - if (itype(i).ne.10 .and. itype(i).ne.21. and. di.eq.0.0d0) - & di=dsc(itype(i)) - vbld(i+nres)=di - if (itype(i).ne.10) then - vbld_inv(i+nres)=1.0d0/di - else - vbld_inv(i+nres)=0.0d0 - endif - if (iti.ne.10) then - alph(i)=alpha(nres+i,i,maxres2) - omeg(i)=beta(nres+i,i,maxres2,i+1) - endif - if(me.eq.king.or..not.out1file)then - if (lprn) - & write (iout,'(a3,i4,7f10.3)') restyp(iti),i,vbld(i), - & rad2deg*theta(i),rad2deg*phi(i),dsc(iti),vbld(nres+i), - & rad2deg*alph(i),rad2deg*omeg(i) - endif - enddo - else if (lprn) then - do i=2,nres - iti=itype(i) - if(me.eq.king.or..not.out1file) - & 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 sc_loc_geom(lprn) - implicit real*8 (a-h,o-z) - include 'DIMENSIONS' -c include "mpif.h" - include 'COMMON.LOCAL' - include 'COMMON.VAR' - include 'COMMON.CHAIN' - include 'COMMON.INTERACT' - include 'COMMON.IOUNITS' - include 'COMMON.GEO' - include 'COMMON.NAMES' - include 'COMMON.CONTROL' - include 'COMMON.SETUP' - double precision x_prime(3),y_prime(3),z_prime(3) - logical lprn - do i=1,nres-1 - do j=1,3 - dc_norm(j,i)=vbld_inv(i+1)*(c(j,i+1)-c(j,i)) - enddo - enddo - do i=2,nres-1 - if (itype(i).ne.10) then - do j=1,3 - dc_norm(j,i+nres)=vbld_inv(i+nres)*(c(j,i+nres)-c(j,i)) - enddo - else - do j=1,3 - dc_norm(j,i+nres)=0.0d0 - enddo - endif - enddo - do i=2,nres-1 - costtab(i+1) =dcos(theta(i+1)) - sinttab(i+1) =dsqrt(1-costtab(i+1)*costtab(i+1)) - cost2tab(i+1)=dsqrt(0.5d0*(1.0d0+costtab(i+1))) - sint2tab(i+1)=dsqrt(0.5d0*(1.0d0-costtab(i+1))) - cosfac2=0.5d0/(1.0d0+costtab(i+1)) - cosfac=dsqrt(cosfac2) - sinfac2=0.5d0/(1.0d0-costtab(i+1)) - sinfac=dsqrt(sinfac2) - it=itype(i) - if (it.ne.10) then -c -C Compute the axes of tghe local cartesian coordinates system; store in -c x_prime, y_prime and z_prime -c - do j=1,3 - x_prime(j) = 0.00 - y_prime(j) = 0.00 - z_prime(j) = 0.00 - enddo - do j = 1,3 - x_prime(j) = (dc_norm(j,i) - dc_norm(j,i-1))*cosfac - y_prime(j) = (dc_norm(j,i) + dc_norm(j,i-1))*sinfac - enddo - call vecpr(x_prime,y_prime,z_prime) -c -C Transform the unit vector of the ith side-chain centroid, dC_norm(*,i), -C to local coordinate system. Store in xx, yy, zz. -c - xx=0.0d0 - yy=0.0d0 - zz=0.0d0 - do j = 1,3 - xx = xx + x_prime(j)*dc_norm(j,i+nres) - yy = yy + y_prime(j)*dc_norm(j,i+nres) - zz = zz + z_prime(j)*dc_norm(j,i+nres) - enddo - - xxref(i)=xx - yyref(i)=yy - zzref(i)=zz - else - xxref(i)=0.0d0 - yyref(i)=0.0d0 - zzref(i)=0.0d0 - endif - enddo - if (lprn) then - do i=2,nres - iti=itype(i) - if(me.eq.king.or..not.out1file) - & write (iout,'(a3,i4,3f10.5)') restyp(iti),i,xxref(i), - & yyref(i),zzref(i) - enddo - endif - return - end -c--------------------------------------------------------------------------- - subroutine sccenter(ires,nscat,sccor) - implicit real*8 (a-h,o-z) - include 'DIMENSIONS' - include 'COMMON.CHAIN' - dimension sccor(3,20) - 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 -c--------------------------------------------------------------------------- - subroutine bond_regular - implicit real*8 (a-h,o-z) - include 'DIMENSIONS' - include 'COMMON.VAR' - include 'COMMON.LOCAL' - include 'COMMON.CALC' - include 'COMMON.INTERACT' - include 'COMMON.CHAIN' - do i=1,nres-1 - vbld(i+1)=vbl - vbld_inv(i+1)=1.0d0/vbld(i+1) - vbld(i+1+nres)=dsc(itype(i+1)) - vbld_inv(i+1+nres)=dsc_inv(itype(i+1)) -c print *,vbld(i+1),vbld(i+1+nres) - enddo - return - end