Last time I described how to derive a geometric derivation of the the famous Rodrigues' rotation formula. While being quite easy to understand there is another way how to get the rotation matrix from the angle-axis representation. The approach I'm gonna go through today derives the rotation matrix from a constant angular velocity of a rigid body. The goal is to get rotation matrix from axis of rotation and angle of rotation I hope I won't spoil your excitement from the discovery if I show you how it's done right now.
The is 3x3 rotation matrix, is a skew symmetric matrix of the unit vector representing rotation axis and the is the rotation angle in radians. By the right hand rule the rotation is counter clockwise. The inverse rotation is described by flipping the rotation axis.
Here is how this equation can be derived. Lets define some basic variables first.
Angular velocity in . It is a pseudovector whose direction represents the axis of rotation and the magnitude specifies the angular speed. The angular speed is a signed quantity which is reflected by the direction of this vector.
Vector from the axis of rotation to the point at the time
Position of the vector at the time
Rotation matrix that moves the vector to the
Instantaneous velocity of the the point at the time
Lets start with the idea of rotating the vector around the axis from it's starting position to new position at the time .
We know that we can rotate the by pre-multiplying by an orthogonal rotation matrix . This matrix moves the vector to the new position for time so we express it as
This means rotates the vector to the new position . In the figure above I labelled the angular displacement as angle . This means the rotation matrix corresponds to the rotation by angle around axis . We can imagine that our asix is aligned with the axis and the point is rotating in the plane.
Now we need to calculate instantaneous velocity of the vector . Why we are doing this will become apparent soon. Instantaneous velocity of the point is a perpendicular vector to the vector which we get by differentiating the above equation with respect to the time.
Note is a constant. This is all nice and well but still quite far from our goal which is to calculate the rotation matrix from the axis and angle of rotation. Ok lets carry on and see whether we can do something about it.
We know the instantaneous velocity of the vector can be also calculated as a cross product of the angular velocity and the rotating point.
And now we have something interesting, two equations that calculate the instantaneous velocity of the same vector. Lets substitute to the
From the equation (1) we can substitute the with
To simplify the calculation we convert the cross product to multiplication with skew symmetric matrix
Where the is the skew-symmetric matrix of the
Now we have to solve this differential equation. Let's rearange terms into the standard form and see how we can solve it.
We can see it is a first order, linear, time invariant equation which has a known solution I'm gonna be more explicit here and show how to solve this equation by the homogeneous and particular solution.
This technique is based on finding two solutions, the particular and homogeneous. The result is then the sum of those two.
The particular solution is found by guessing a constant variable that would satisfy the equation. Since derivation of constant is zero the differential term drops out. In our case the particular solution is zero.
Homogeneous solutions is found by replacing the right side of the differential equation with zero (removing the driving force) and guessing what solutions would satisfy this equation. In our case the equation already equals zero. As the solution we guess The reason for this guess is that the derivative of the which nicely removes the derivation without any mathematical complications.
Since the particular solution is zero the complete solution is only the homogeneous one.
What we need to do now is to get rid of the time parameter somehow and relate the angular velocity to the rotation angle While this may look like a complicated task in reality it is quite simple. Since we can choose any time interval we want we chose Now the angular speed (the magnitude of the angular velocity ) is so in one second the angle displacement will be just .
And there you go we can calculate the rotation matrix from a given rotation axis and the angle of rotation. The direction of the vector is the axis of rotation and the magnitude is the angle of rotation. We can split it into the angle and axis components.
For unit vector representing the axis of rotation we get
We could stop here and be happy with this formula. In the Octave we can get the rotation matrix from angle axis straight this way.
p = pi/3; % Rotation angle in anti-clockwise direction n = [1; 2; 3]; % Axis of rotation % Normalization of the rotation axis vector n_norm = n/norm(n); % Create the skew-symmetric matrix n_skew = [ 0 -n_norm(3) n_norm(2); n_norm(3) 0 -n_norm(1); -n_norm(2) n_norm(1) 0]; % Calculate the rotation matrix R = expm(n_skew*p)
Or in the Python with SciPy library
from numpy import * from scipy.linalg import * # Rotation axis n = mat('[1. 2. 3.]'); # Normalize the rotation axis n = n/numpy.linalg.norm(n) # Build skew-symmetric matrix n_skew = mat([[0, -n[:,2], n[:,1]], [n[:,2], 0, -n[:,0]], [-n[:,1], n[:,0], 0]]) #Angle of rotation alpha = math.pi / 3. R = expm3(n_skew*alpha)
The Octave and SciPy in Python nicely hide the calculation of the matrix exponential which is great but in many practical applications we need to do it on our own, many times even with very limited resources. Next time I will dive into how to calculate the matrix exponential and what everything in under the bonnet.