Rotation Matrix

Andrej Sali sali at tamika.harvard.edu
Thu Nov 24 11:02:07 EST 1994


In article <1994Nov24.152711.1606 at inca.comlab.ox.ac.uk> smb at bioch.ox.ac.uk  
(Simon Brocklehurst) writes:
> Is there a kind sole out there with easy access to  Crystallography
> Table type books?
> 
> What I need is the matrix for rotation about an arbritary axis.
> 
> Pleeeeeeaaaasssssseeeeee save me from having to work it out from  
scratch!
> 
> Thanks in advance,
> 
>  
_________________________________________________________________________
> |
> |  ,_ o     Simon M. Brocklehurst,
> | /  //\,   Oxford Centre for Molecular Sciences,
> |   \>> |   Department of Biochemistry, University of Oxford,
> |    \\,    Oxford, UK.
> |           E-mail: smb at bioch.ox.ac.uk
>  
|________________________________________________________________________



c --- This subroutine rotates N points in X,Y,Z around axis through the
c     origin. Angle of rotation is in radians (ANG). Exact formula is 
c     used and the routine is
c     therefore appropriate for large angles, too. Result is in X,Y,Z.
c     Translate the points before and after calling this routine if you 
c     want rotation around an axis that does not go through origin.

      subroutine scrot1(x,y,z,n,xaxis,yaxis,zaxis,ang)
      implicit none
      integer n,i
      real x(n), y(n), z(n), xaxis, yaxis, zaxis, v(3), rtl(3)
      real sinl,cosl,rl,ang,xa,ya,za
      xa = xaxis
      ya = yaxis
      za = zaxis
      call norml1(xa,ya,za)
      sinl = sin(ang)
      cosl = 1.0-cos(ang)
      do  i = 1, n
        call vectr1(xa,ya,za,x(i),y(i),z(i),v(1),v(2),v(3))
        rl = x(i)*xa + y(i)*ya + z(i)*za
        rtl(1) = rl*xa
        rtl(2) = rl*ya
        rtl(3) = rl*za
        x(i) = x(i) + v(1)*sinl + (rtl(1)-x(i))*cosl
        y(i) = y(i) + v(2)*sinl + (rtl(2)-y(i))*cosl
        z(i) = z(i) + v(3)*sinl + (rtl(3)-z(i))*cosl
      end do
      return
      end

c --- vector (x1,x2,x3) is normalized to unit length
      subroutine norml1(x1,x2,x3)
      implicit none
      real eps
      parameter (eps=1.e-8)
      real size,x1,x2,x3
      size = sqrt(x1*x1+x2*x2+x3*x3)
      if (abs(size) .lt. eps) stop 'norml1_E> size < eps.'
      x1 = x1/size
      x2 = x2/size
      x3 = x3/size
      return
      end

c --- calculates the vector product
      subroutine vectr1(a1,a2,a3,b1,b2,b3,ab1,ab2,ab3)
      implicit none
      real a1,a2,a3,b1,b2,b3,ab1,ab2,ab3
      ab1 = a2*b3 - a3*b2
      ab2 = a3*b1 - a1*b3
      ab3 = a1*b2 - a2*b1
      return
      end





More information about the Xtal-log mailing list