Euler angles and gimbal lock Can someone show mathematically how gimbal lock happens when doing matrix rotation with Euler angles for yaw, pitch, roll?  I'm having a hard time understanding what is going on even after reading several articles on Google.
Is the only work-around to use quaternions?
 A: If you only care about the beginning and end, you can express any rotation with the Euler angles. But if you want a smooth transition between the beginning and the end, then with some starting orientations, you have a problem, and this is called the gimbal lock.
For example, here we start from (x,y,z) = (1,0,0) and simultaneously turn yaw left and lift pitch up at 5 degree steps, until both reach 45 degrees. The black dots move on the surface of the unit sphere. (Roll could be whatever, it cannot be seen as I am only plotting dots.)

But say, we start from a position pointing directly upwards, so our object is pointing to (x,y,z) = (0,0,1) and we want to rotate it down around the x-axis by 45 degrees. Now we have to turn yaw (or heading) 90 degrees, and lower the pitch 45 degrees (and perhaps roll -90 degrees, depending what orientation we want, but my dots cannot illustrate this).
Now if we do the yaw and pitch transformations simultaneously:  
(yaw=0,pitch=0), (yaw=5,pitch=-5), (yaw=10,pitch=-10), ...  
we get the black dots in the figure below, whereas a straight line by 5 degree turns around the x-axis is plotted with the white dots. We can still reach the same endpoint, but we cannot start to turn directly into the "correct" or "straight" direction.

You could redefine the Euler angles so that the first rotation is not around the z-axis, but around the x-axis, and then this situation could be remedied. But then a similar problem could happen with an object that is initially pointing to the x-direction. However you reorder the three Euler rotations, there is always some orientation from which you cannot start turning directly to a certain direction.
Here is also a Youtube-video on the topic.
A: "Euler Angles" you can think of as a function $(S^1)^3 \to SO_3$ or $\mathbb R^3 \to SO_3$.  The derivative of this function does not always have rank 3, so you have degenerate submanifolds where the function is many-to-one.  In this special case that's called "gimbal lock".   
One formalism that avoids this is quaternions.  You can of course use other formalisms, and many other formalisms are naturally related to the quaternion version, so people tend to gravitate to the quaternion version.  One version that's closely related to quaternions would be to use the exponential map for the unit quaternion group.  But this also has "gymbal lock" but of a different kind.  But it does have the rather appealing interpretation as rotations about arbitrary axis -- this is perhaps more useful if you're only interested in rotations that differ from the identity matrix (or some given matrix) by a small amount, they're very natural coordinates on "small scales" in $SO_3$. 
Are there any special properties you'd like for coordinates on $SO_3$? That might give a sense for where you want to go with this. 
Edit in response to 1st comment: There's various conventions for Euler angles. Let's use what Wikipedia calls "Proper Euler Angles"   http://en.wikipedia.org/wiki/Euler_angles
Given a unit vector $\vec v \in \mathbb R^3$ and a number $\theta \in \mathbb R$, let $R_{\vec v, \theta}$ be the rotation matrix which fixes $\vec v$ and which rotates counter-clockwise by an angle $\theta$ in the plane orthogonal to $\vec v$.  One version of "Proper Euler angles" would be:
$$R_{(1,0,0), \theta_1} R_{(0,1,0), \theta_2} R_{(1,0,0), \theta_3}$$
The "Tait-Bryan" variant would use the three distinct vectors, notice above we only used the standard vectors for the xy-plane.   The Tait-Bryan variant is the one appropriate to pitch/roll/yaw of an aeroplane.  Wikipedia has an excellent picture:

In this case you'd have 
$$R_{(1,0,0),\theta_1} R_{(0,1,0), \theta_2} R_{(0,0,1),\theta_3} $$
$\theta_3$ would be the roll, $\theta_2$ the pitch and $\theta_1$ the yaw. To relate my coordinates to the picture, $(1,0,0)$ is the direction the plane is pointing. $(0,1,0)$ is the direction of the left wing. $(0,0,1)$ is the direction of the yellow axis sticking out of the top of the plane. 
You can write out the matrices quite explicitly:
$$ R_{(1,0,0), \theta} = 
 \pmatrix{ 1 & 0 & 0 \cr 
           0 & \cos \theta & -\sin \theta \cr
           0 & \sin \theta & \cos \theta }$$
So in this case, one occurance of "gimbal lock" is $\theta_2 = 0$.  In the Tait-Bryan variant it would be when the aeroplane is either pointing straight up or down, which is $\theta_2 = \pm \pi/2$. 
A: This refers to Ryan's answer, but was too long for a comment.
I was getting confused for a while too, but I think I understand now.
To avoid misunderstandings, the Tait–Bryan angles should be called
heading/elevation/bank instead of yaw/pitch/roll. The heading angle
refers to rotations about the $z$ axis in space, not about the yaw
axis which by definition is perpendicular to the "body-wing-plane" of
the airplane, so that "excellent picture" is in fact a bit misleading
in this context! The gimbal lock illustrated on
Wikipedia
occurs when the plane is pointing straight up or down (elevation = 90
degrees), so that changing the bank is the same thing as changing the
heading. (That is, there is momentarily no way to change the yaw by
manipulating the heading/elevation/bank angles!)
A: Many answers have been given that are mathematically very sound and correct. Here is another way, using only high school mathematics:
Suppose, we use the ZXY Euler rotations: first about the Z-axis (yaw), then about the X-axis (roll) and finally about the y-axis(pitch).
$R_z(\psi) = \left(\begin{array}{ccc}
\cos \left(\psi \right) & -\sin \left(\psi \right) & 0\\
\sin \left(\psi \right) & \cos \left(\psi \right) & 0\\
0 & 0 & 1
\end{array}\right)$
$R_x(\phi) = \left(\begin{array}{ccc}
1 & 0 & 0\\
0 & \cos \left(\phi \right) & -\sin \left(\phi \right)\\
0 & \sin \left(\phi \right) & \cos \left(\phi \right)
\end{array}\right)$
$R_y(\theta) = \left(\begin{array}{ccc}
\cos \left(\theta \right) & 0 & \sin \left(\theta \right)\\
0 & 1 & 0\\
-\sin \left(\theta \right) & 0 & \cos \left(\theta \right)
\end{array}\right)$
The rotation matrix in ZXY configuration would be:
$R = R_z(\psi)R_x(\phi)R_y(\theta)$
$R = \left(\begin{array}{ccc}
\cos \left(\psi \right)\,\cos \left(\theta \right)-\sin \left(\phi \right)\,\sin \left(\psi \right)\,\sin \left(\theta \right) & -\cos \left(\phi \right)\,\sin \left(\psi \right) & \cos \left(\psi \right)\,\sin \left(\theta \right)+\cos \left(\theta \right)\,\sin \left(\phi \right)\,\sin \left(\psi \right)\\
\cos \left(\theta \right)\,\sin \left(\psi \right)+\cos \left(\psi \right)\,\sin \left(\phi \right)\,\sin \left(\theta \right) & \cos \left(\phi \right)\,\cos \left(\psi \right) & \sin \left(\psi \right)\,\sin \left(\theta \right)-\cos \left(\psi \right)\,\cos \left(\theta \right)\,\sin \left(\phi \right)\\
-\cos \left(\phi \right)\,\sin \left(\theta \right) & \sin \left(\phi \right) & \cos \left(\phi \right)\,\cos \left(\theta \right)
\end{array}\right)$
Now, let's replace $\phi =\frac{\pi }{2}$
You will get R as:
$R = \left(\begin{array}{ccc}
\cos \left(\psi \right)\,\cos \left(\theta \right)-\sin \left(\psi \right)\,\sin \left(\theta \right) & 0 & \cos \left(\psi \right)\,\sin \left(\theta \right)+\cos \left(\theta \right)\,\sin \left(\psi \right)\\
\cos \left(\psi \right)\,\sin \left(\theta \right)+\cos \left(\theta \right)\,\sin \left(\psi \right) & 0 & \sin \left(\psi \right)\,\sin \left(\theta \right)-\cos \left(\psi \right)\,\cos \left(\theta \right)\\
0 & 1 & 0
\end{array}\right)$
Using trigonometric identities, we get:
$R = \left(\begin{array}{ccc}
\cos \left(\psi +\theta \right) & 0 & \sin \left(\psi +\theta \right)\\
\sin \left(\psi +\theta \right) & 0 & -\cos \left(\psi +\theta \right)\\
0 & 1 & 0
\end{array}\right)$
There is infinite $\psi$ and $\theta$ that create the same rotation matrix.
For $\phi =-\frac{\pi }{2}$:
$R = \left(\begin{array}{ccc}
\cos \left(\psi \right)\,\cos \left(\theta \right)+\sin \left(\psi \right)\,\sin \left(\theta \right) & 0 & \cos \left(\psi \right)\,\sin \left(\theta \right)-\cos \left(\theta \right)\,\sin \left(\psi \right)\\
\cos \left(\theta \right)\,\sin \left(\psi \right)-\cos \left(\psi \right)\,\sin \left(\theta \right) & 0 & \cos \left(\psi \right)\,\cos \left(\theta \right)+\sin \left(\psi \right)\,\sin \left(\theta \right)\\
0 & -1 & 0
\end{array}\right)$
$R = \left(\begin{array}{ccc}
\cos \left(\psi -\theta \right) & 0 & -\sin \left(\psi -\theta \right)\\
\sin \left(\psi -\theta \right) & 0 & \cos \left(\psi -\theta \right)\\
0 & -1 & 0
\end{array}\right)$
You can use Symbolic math tools to create these simplifications. I have used MATLAB:
syms psi

R_z = [ cos(psi), -sin(psi), 0; 
        sin(psi),  cos(psi), 0;
        0       ,  0       , 1
        ];
syms phi 

R_x = [ 1, 0,           0; 
        0, cos(phi), -sin(phi);
        0, sin(phi),  cos(phi)
        ];

syms theta 

R_y = [  cos(theta), 0, sin(theta); 
         0         , 1, 0;
        -sin(theta), 0, cos(theta)
        ];
R = R_z*R_x*R_y;
R_new_1_expr = subs(R, phi, pi/2)
R_new_1_simplified = simplify(R_new_1_expr);
display(R_new_1_simplified)


$\left(\begin{array}{ccc}
\cos \left(\psi +\theta \right) & 0 & \sin \left(\psi +\theta \right)\\
\sin \left(\psi +\theta \right) & 0 & -\cos \left(\psi +\theta \right)\\
0 & 1 & 0
\end{array}\right)$

