r9177 - /branches/frame_order/maths_fns/rotation_matrix.py -- July 02, 2009 - 14:03

Author: bugman
Date: Thu Jul  2 14:03:17 2009
New Revision: 9177

URL: http://svn.gna.org/viewcvs/relax?rev=9177&view=rev
Log:
Added the new function R_to_axis_angle().


Modified:
    branches/frame_order/maths_fns/rotation_matrix.py

Modified: branches/frame_order/maths_fns/rotation_matrix.py
URL: 
http://svn.gna.org/viewcvs/relax/branches/frame_order/maths_fns/rotation_matrix.py?rev=9177&r1=9176&r2=9177&view=diff
==============================================================================
--- branches/frame_order/maths_fns/rotation_matrix.py (original)
+++ branches/frame_order/maths_fns/rotation_matrix.py Thu Jul  2 14:03:17 2009
@@ -21,8 +21,8 @@
 
###############################################################################
 
 # Python module imports.
-from math import acos, cos, pi, sin
-from numpy import array, cross, dot, float64
+from math import acos, atan2, cos, pi, sin
+from numpy import array, cross, dot, float64, hypot, zeros
 from numpy.linalg import norm
 from random import gauss, uniform
 
@@ -175,6 +175,42 @@
     matrix[2, 2] = z*zC + ca
 
 
+def R_to_axis_angle(matrix):
+    """Convert the rotation matrix into the axis-angle notation.
+
+    Conversion equations
+    ====================
+
+    From Wikipedia (http://en.wikipedia.org/wiki/Rotation_matrix), the 
conversion is given by::
+
+        x = Qzy-Qyz
+        y = Qxz-Qzx
+        z = Qyx-Qxy
+        r = hypot(x,hypot(y,z))
+        t = Qxx+Qyy+Qzz
+        theta = atan2(r,t-1)
+
+    @param matrix:  The 3x3 rotation matrix to update.
+    @type matrix:   3x3 numpy array
+    @return:    The 3D rotation axis and angle.
+    @rtype:     numpy 3D rank-1 array, float
+    """
+
+    # Axes.
+    axis = zeros(3, float64)
+    axis[0] = matrix[2,1] - matrix[1,2]
+    axis[1] = matrix[0,2] - matrix[2,0]
+    axis[2] = matrix[1,0] - matrix[0,1]
+
+    # Angle.
+    r = hypot(axis[0], hypot(axis[1], axis[2]))
+    t = matrix[0,0] + matrix[1,1] + matrix[2,2]
+    theta = atan2(r, t-1)
+
+    # Return the data.
+    return axis, theta
+
+
 def R_euler_zyz(matrix, alpha, beta, gamma):
     """Function for calculating the z-y-z Euler angle convention rotation 
matrix.