function R = RotationMatrix3(theta,omega)
% function R = RotationMatrix3(theta,omega)
% constructs a 3d rotation matrix for a rotation through an angle theta (in radians) about an axis given by the vector omega
%
% Example: 
% suppose we have a vector v1 = [3;4;5]; and we want to rotate this vector through 30degrees about the axis [1 1 1].
% Then the rotated vector is given by v2 = R * v1; where R = RotationMatrix3(pi/6,[1 1 1]).
%
% Jenny Read 2/4/2003

if prod(size(theta))~=1
    disp('theta must be a scalar')
    return
end
if prod(size(omega))~=3
    disp('omega must be a 3-entry vector')
    return
end
%Normalise omega
omega = omega ./ sqrt(sum(omega.^2));
c = cos(theta);
s = sin(theta);
w1 = omega(1);
w2 = omega(2);
w3 = omega(3);
R = [    1+(1-c)*(-w3^2-w2^2), -s*w3+(1-c)*w2*w1,  s*w2+(1-c)*w3*w1 ; ...
        s*w3+(1-c)*w2*w1,    1+(1-c)*(-w3^2-w1^2), -s*w1+(1-c)*w3*w2 ; ...
        -s*w2+(1-c)*w3*w1,  s*w1+(1-c)*w3*w2,    1+(1-c)*(-w2^2-w1^2)];

