dfa & multichain cluster

[unres.git] / source / unres / src-HCD-5D / chainbuild.F

diff --git a/source/unres/src-HCD-5D/chainbuild.F b/source/unres/src-HCD-5D/chainbuild.F
index 51419ef..7902f15 100644 (file)
--- a/source/unres/src-HCD-5D/chainbuild.F
+++ b/source/unres/src-HCD-5D/chainbuild.F
@@ -240,6 +240,86 @@ C
       return
       end
 c-----------------------------------------------------------------------------
+      subroutine sc_coord_rebuild(i)
+C 
+C Locate the side-chain centroid i, 1 < i < NRES. Put in C(*,NRES+i).
+C
+      implicit none
+      include 'DIMENSIONS'
+      include 'COMMON.CHAIN'
+      include 'COMMON.LOCAL'
+      include 'COMMON.GEO'
+      include 'COMMON.VAR'
+      include 'COMMON.IOUNITS'
+      include 'COMMON.NAMES'
+      include 'COMMON.INTERACT'
+      integer i,j,k
+      double precision xx(3)
+      double precision dsci,dsci_inv,alphi,omegi,cosalphi,sinalphi,
+     &  cosomegi,sinomegi,rp(3),theta2,cost2,sint2,rj
+      double precision scalar
+      double precision ref(3,3),scalp,sscalp,refnorm
+      double precision alpha,beta
+c      dsci=dsc(itype(i))
+c      dsci_inv=dsc_inv(itype(i))
+      dsci=vbld(i+nres)
+      dsci_inv=vbld_inv(i+nres)
+#ifdef OSF
+      alphi=alph(i)
+      omegi=omeg(i)
+      if (alphi.ne.alphi) alphi=100.0
+      if (omegi.ne.omegi) omegi=-100.0
+#else
+      alphi=alph(i)
+      omegi=omeg(i)
+#endif
+      cosalphi=dcos(alphi)
+      sinalphi=dsin(alphi)
+      cosomegi=dcos(omegi)
+      sinomegi=dsin(omegi) 
+      rp(1)= cosalphi
+      rp(2)= sinalphi*cosomegi
+      rp(3)=-sinalphi*sinomegi
+c Build the reference system
+      do j=1,3
+        ref(j,1)=-dc_norm(j,i-1)+dc_norm(j,i)
+      enddo
+      refnorm=dsqrt(scalar(ref(1,1),ref(1,1)))
+      do j=1,3
+        ref(j,1)=ref(j,1)/refnorm
+      enddo
+      scalp=scalar(ref(1,1),dc_norm(1,i))
+      sscalp=1.0d0/dsqrt(1.0d0-scalp*scalp)
+      do j=1,3
+        ref(j,2)=(dc_norm(j,i)-scalp*ref(j,1))*sscalp
+      enddo
+      ref(1,3)= ref(2,1)*ref(3,2)-ref(3,1)*ref(2,2)
+      ref(2,3)=-ref(1,1)*ref(3,2)+ref(3,1)*ref(1,2)
+      ref(3,3)= ref(1,1)*ref(2,2)-ref(2,1)*ref(1,2)
+c      do j=1,3
+c        write (iout,*) j,scalar(ref(1,j),ref(1,1)),
+c     &  scalar(ref(1,j),ref(1,2)),scalar(ref(1,j),ref(1,3))
+c      enddo
+c Bring the coordinates to the global reference system
+      do j=1,3
+        dc_norm(j,nres+i)=0.0d0
+        do k=1,3
+          dc_norm(j,nres+i)=dc_norm(j,nres+i)+ref(j,k)*rp(k)
+        enddo
+        dc(j,nres+i)=dc_norm(j,nres+i)*dsci
+        c(j,nres+i)=c(j,i)+dc(j,nres+i)
+      enddo 
+c      write (iout,*) scalar(dc_norm(1,i+nres),dc_norm(1,i+nres)),
+c     &  dsqrt(scalar(dc(1,i+nres),dc(1,i+nres)))
+c Check the internal coordinates
+c      c(:,2*nres+1)=ref(:,1)+c(:,i)
+c      write (iout,*) "alpha",rad2deg*alphi,
+c     & rad2deg*alpha(nres+i,i,2*nres+1)
+c      write (iout,*) "omega",rad2deg*omegi,
+c     & rad2deg*beta(nres+i,i,2*nres+1,i+1)
+      return
+      end
+c-----------------------------------------------------------------------------
       subroutine locate_side_chain(i)
 C 
 C Locate the side-chain centroid i, 1 < i < NRES. Put in C(*,NRES+i).
@@ -296,7 +376,7 @@ cd   &   xp,yp,zp,(xx(k),k=1,3)
       xx(2)=xloc(2,i)*r(2,2,i-1)+xloc(3,i)*r(2,3,i-1)
       xx(3)=xloc(2,i)*r(3,2,i-1)+xloc(3,i)*r(3,3,i-1)
       do j=1,3
-       xrot(j,i)=xx(j)
+        xrot(j,i)=xx(j)
       enddo
       do j=1,3
         rj=0.0D0