Adam's changes

[unres.git] / source / cluster / wham / src-HCD-5D / readpdb.F

diff --git a/source/cluster/wham/src-HCD-5D/readpdb.F b/source/cluster/wham/src-HCD-5D/readpdb.F
index 0167c00..f7cfb86 100644 (file)
--- a/source/cluster/wham/src-HCD-5D/readpdb.F
+++ b/source/cluster/wham/src-HCD-5D/readpdb.F
@@ -13,13 +13,13 @@ C geometry.
       include 'COMMON.NAMES'
       include 'COMMON.CONTROL'
       integer i,j,ibeg,ishift1,ires,iii,ires_old,ishift,ity
-      logical lprn /.false./,fail
+      logical lprn /.false./,fail,sccalc
       double precision e1(3),e2(3),e3(3)
       double precision dcj,efree_temp
       character*3 seq,res
       character*5 atom
       character*80 card
-      double precision sccor(3,20)
+      double precision sccor(3,50)
       integer rescode
       integer iterter(maxres)
       efree_temp=0.0d0
@@ -30,6 +30,7 @@ c      write (2,*) "UNRES_PDB",unres_pdb
       ires=0
       ires_old=0
       iii=0
+      sccalc=.false.
       lsecondary=.false.
       nhfrag=0
       nbfrag=0
@@ -61,6 +62,7 @@ c      write (2,*) "UNRES_PDB",unres_pdb
           itype(ires_old-1)=ntyp1
           iterter(ires_old-1)=1
           itype(ires_old)=ntyp1
+          iterter(ires_old)=1
           ishift1=ishift1+1
           ibeg=2
 !          write (iout,*) "Chain ended",ires,ishift,ires_old
@@ -71,12 +73,14 @@ c      write (2,*) "UNRES_PDB",unres_pdb
           else
             call sccenter(ires,iii,sccor)
           endif
-c          iii=0
+          iii=0
+          sccalc=.true.
         endif
 ! Read free energy
-        if (index(card,"FREE ENERGY").gt.0) read(card(35:),*) efree_temp
+c        if (index(card,"FREE ENERGY").gt.0) read(card(35:),*) efree_temp
 ! Fish out the ATOM cards.
         if (index(card(1:4),'ATOM').gt.0) then  
+          sccalc=.false.
           read (card(12:16),*) atom
 c          write (2,'(a)') card
 !          write (iout,*) "! ",atom," !",ires
@@ -94,12 +98,13 @@ c          write (2,'(a)') card
 !              write (iout,*) "Calculating sidechain center iii",iii
               if (unres_pdb) then
                 do j=1,3
-                  dc(j,ires+nres)=sccor(j,iii)
+                  dc(j,ires_old)=sccor(j,iii)
                 enddo
               else
                 call sccenter(ires_old,iii,sccor)
               endif
               iii=0
+              sccalc=.true.
             endif
 ! Start new residue.
             if (res.eq.'Cl-' .or. res.eq.'Na+') then
@@ -220,7 +225,7 @@ C          endif !unres_pdb
         endif  !itype.eq.ntyp1
       enddo
 C Calculate the CM of the last side chain.
-      call sccenter(ires,iii,sccor)
+      if (.not.sccalc) call sccenter(ires,iii,sccor)
       nsup=nres
       nstart_sup=1
       if (itype(nres).ne.10) then
@@ -327,7 +332,7 @@ c---------------------------------------------------------------------------
       character*3 seq,res
 c      character*5 atom
       character*80 card
-      dimension sccor(3,20)
+      dimension sccor(3,50)
       integer rescode
       logical lside,lprn
        if (lprn) then 
@@ -342,17 +347,20 @@ c      character*5 atom
      & '     Gamma'
         endif
        endif
-      do i=1,nres-1
+      do i=2,nres
         iti=itype(i)
-        if (dist(i,i+1).lt.2.0D0 .or. dist(i,i+1).gt.5.0D0) then
+c        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.ntyp1 .and. itype(i).ne.ntyp1 .and.
+     &    (dist(i,i-1).lt.1.0D0 .or. dist(i,i-1).gt.6.0D0)) then
           write (iout,'(a,i4)') 'Bad Cartesians for residue',i
-ctest          stop
+c          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)
+        vbld(i)=dist(i-1,i)
+        vbld_inv(i)=1.0d0/vbld(i)
+        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.ntyp1) then
 c          theta(3)=90.0d0*deg2rad
@@ -493,7 +501,7 @@ c---------------------------------------------------------------------------
       implicit real*8 (a-h,o-z)
       include 'DIMENSIONS'
       include 'COMMON.CHAIN'
-      dimension sccor(3,20)
+      dimension sccor(3,50)
       do j=1,3
         sccmj=0.0D0
         do i=1,nscat
@@ -555,7 +563,7 @@ C and convert the peptide geometry into virtual-chain geometry.
       character*3 seq,res
       character*5 atom
       character*80 card
-      double precision sccor(3,20)
+      double precision sccor(3,50)
       integer rescode,iterter(maxres)
       do i=1,maxres
          iterter(i)=0
@@ -801,7 +809,7 @@ C 2/15/2013 by Adam: corrected insertion of the first dummy residue
         endif
       endif
 C Calculate internal coordinates.
-      if (lprn) then
+      if (out_template_coord) then
       write (iout,'(/a)') 
      &  "Cartesian coordinates of the reference structure"
       write (iout,'(a,3(3x,a5),5x,3(3x,a5))') 
@@ -813,6 +821,7 @@ C Calculate internal coordinates.
       enddo
       endif
 C Calculate internal coordinates.
+#ifdef DEBUG
        write (iout,'(a)') 
      &   "Backbone and SC coordinates as read from the PDB"
        do ires=1,nres
@@ -820,7 +829,8 @@ C Calculate internal coordinates.
      &    ires,itype(ires),restyp(itype(ires)),(c(j,ires),j=1,3),
      &    (c(j,nres+ires),j=1,3)
        enddo
-      call int_from_cart(.true.,.false.)
+#endif
+      call int_from_cart(.true.,out_template_coord)
       call sc_loc_geom(.false.)
       do i=1,nres
         thetaref(i)=theta(i)
@@ -854,5 +864,3 @@ c     &   vbld_inv(i+nres)
 
       return
       end
-      
-