Well-behaved Euler Angles
Greg Abbas

Here's another technical article. The last one was about rendering; this time I'm switching to animation (specifically, rotation).

Euler angles have well-known limitations, such as gimbal lock. Another slightly-less-well-known issue is the way in which there are many different (rx,ry,rz) triplets that all encode the same orientation, and writing a program to choose the "right" triplet for a given orientation is sometimes a lot more difficult than it would seem. This problem can come up when implementing a user-interface control for instance, or a dynamics engine. Any time that a program has an non-Euler orientation (expressed as a quaternion, a matrix, and axis and an angle, or whatever) and it has to convert it to a triplet of Euler angles that must fit in smoothly with a sequence of other Euler angles, you have to be careful how you choose which triplet to use.

It's not hard to find convert a quaternion to a valid triplet triplet of Euler angles, of course. Ken Shoemake's article in Graphics Gems IV (for instance) demonstrates how to do that. But applying his recipe to a sequence of smoothly-varying orientations doesn't necessarily produce a sequence of smoothly-varying Euler triplets. For instance, because you can add or subtract any multiple of 360 degrees to any of the three Euler angles and the triplet still describes the same orientation, our valid sequence of Euler triplets might have objectionable 360-degree discontinuities in it. We could make a rule that indicates which value is preferred (like, the angle with the smallest absolute value), but that would only serve to define where the discontinuities occur, not eliminate them. For instance, imagine that we're writing a physics simulator and we run it on an object with a angular velocity around the x-axis. We would expect the x-angle to be a linear function of time, but if compute it by applying Ken's formulae to quaternions, we get a sawtooth function instead:

Sawtooth

This problem is solved easily enough by taking the previous frame's Euler angles into account when converting each frame's orientation from a quaternion to Euler angles. After finding a triplet that corresponds to the given quaternion (which probably has angles in some nominal range like -π ... π), we can add or subtract some multiple of 2π to each angle so that it's as close to the previous values as possible. That'll restore our x-angle function to the linear shape that we'd expect it to have.

There's one more problem though, and this one's subtler. If you implement a rotation manipulator and you notice that 180-degree "flips" are creeping into your results, you'll realize there's another choice (like the ±2π one described above) that our conversion has to make. To understand this one, let's start by visualizing an example. Suppose we have an object that's rotated by 180 degrees around its z axis, and -90 degrees around its x axis (assuming ZXY rotation order). That object would be facing up, with its head towards us. Now suppose we took a different object and rotated it by -90 degrees around its x axis, and 180 degrees around its y axis. The second object would be in exactly the same position: facing up, with its head towards us. But these two sets of angles (-90,0,180 and -90,180,0) aren't merely different by multiples of 360 degrees... something else is going on. With a little bit more spatial visualization, we notice that it seems we can generalize it with the following relation:

Y(φ) X(θ) Z(ψ) = Y(π + φ) X(π - θ) Z(π + ψ)

where X, Y, and Z are rotation operators, and [θ, φ, ψ] is our ZXY Euler triplet. But if we're going to use this "180 degree identity", we need more than spatial intuition... we must show that it's really correct. We can do that by converting the rotation operators to quaternions, so that we have an algebra with which to manipulate them. Because we can express rotation θ around a unit vector a using the quaternion

q = (cos(θ/2), a sin(θ/2)),

our three rotation operators become

X(θ) = cos(θ/2) + i sin(θ/2)
Y(φ) = cos(φ/2) + j sin(φ/2)
Z(ψ) = cos(ψ/2) + k sin(ψ/2)

And the two sides of the 180-degree identity can be written

(cos(φ/2) + j sin(φ/2)) (cos(θ/2) + i sin(θ/2)) (cos(ψ/2) + k sin(ψ/2))

and

(-sin(φ/2) + j cos(φ/2)) (sin(θ/2) + i cos(θ/2)) (-cos(ψ/2) + k cos(ψ/2)).

If you multiply-out both of those expressions (keeping in mind that i2 = j2 = k2 = i j k = -1) then you'll find them equal. That's more algebra than I'm going to include in this article; I verified it using Mathematica and you should feel free to convince yourself, too. The point is that you can use this identity (in combination with the 360-degree property mentioned above) to convert quaternions to Euler angles gracefully, as follows:

Algorithm

Here's a summary of this technique, then. If you have a set of "reference" angles [θ, φ, ψ] and a quaternion q′ that you want to convert to "similar" angles [θ′, φ′, ψ′], then first convert the quaternion q′ to a triplet [θ0, φ0, ψ0]. Then add/subtract multiples of 360 degrees:

θ1 = θ0 + nθ1
φ1 = φ0 + nφ1
ψ1 = ψ0 + nψ1

(where nθ1, nφ1, and nψ1 are integers) in order to minimize the distance from the reference angles,

δ1 = | θ1 - θ | + | φ1 - φ | + | ψ1 - ψ |,

do the same thing for the angles produced by the 180-degree identity:

θ2 = π + θ1 + nθ2
φ2 = π - φ1 + nφ2
ψ2 = π + ψ1 + nψ2
δ2 = | θ2 - θ | + | φ2 - φ | + | ψ2 - ψ |,

and pick whichever triplet ([&theta1, &phi1, &psi1] or [&theta2, &phi2, &psi2]) is closer to the reference triplet.

Further reading

By the way, if you're interested in this kind of thing, check out James Diebel's "quick reference guide". (Just kidding, it's anything but quick... it's exhaustive!). Also, Andrew J. Hanson is the biggest fan of quaternions that I've ever met and he has a book on visualizing them. I haven't read it, but if it's half as enthusiastic as his lectures then it's probably an insightful read.

11 response to "Well-behaved Euler Angles"

  1. 282 Jesse July 18, 2007 at 2:10 p.m.

    Dear Greg,
    Hmm, this is a problem I've been having in implementing a navigator control, as you mention in the article. Unfortunatly the camera is being passed to me as a matrix. If I convert my matrix to a quaternion(say using the method listed here http://web.archive.org/web/2004102900385...), then converted to euler angles using your algorithm do you think it would work?
    Regards,
    J.


  2. 283 Jesse July 18, 2007 at 2:52 p.m.

    Also could you explain a bit more about trying various integeters to minmize the distance? It seems like this would be hard to implement in code, how do you know when to stop?


  3. 312 Robert July 19, 2007 at 5:53 p.m.

    I too am VERY interested in this. My situation is that bone rotations are always represented as Euler (I have no choice in this and they can be in any one of six orders per bone: XYZ, XZY, YXZ, YZX, ZXY, ZYX) and retaining their relation is paramount for both deformations, animation, and other dependent values. What I have to start with though is a matrix that must be converted back into exact Euler representations (not equivalents like from Shoemake as you mention). Also, as you mention, a heuristic technique using the previous Euler angles can be beneficial to stabilizing the variety of discontinuous results from matrix to Euler angle conversion. I have scoured the ends of the universe seeking such information (including Google-a-plenty, ACM, IEEE, and any other source that I could lay my hands upon). The matrix can be converted to a quaternion but I'm not following your 'Algorithm' very well - on first inspection.

    I will read this two or ten more times and see if I am capable of implementation and verification of the heuristics - hopefully you will not mind answering questions if I get stuck (?).

    Thank you very much!


  4. 317 Robert July 19, 2007 at 9:42 p.m.

    This is why I hate mathematicians! ;) You provide a couple of formulae and some vague ideas of methodology but don't understand the complexity of implementing it in computer code!

    So far I can't get consistent results. You must remember that the original and new angles can be +/+, -/-, +/-, -/+. Additionally, you aren't specifying more details about the second part (pi or 180d). Is it iterative like the other (add/substract 180d to minimize the distance) or a single action. Plus I don't understand the addition of integer*2pi in that part of the equation at all (since you don't explain it thoroughly). I'm also uncertain about the change of sign in the second angle in this (pi - angle2) instead of (pi + angle2).

    I understand the basic concept - first put the angles in the closest relation with respect to full revolutions (360d) and then with respect to 180d inversions. Some form of pseudo-code would be extremely valuable - can you even do it! I've said this many times - mathematical formulae and code implementations are two different beasts. For the latter, the algorithmic breakdown must be precise enough to avoid just what I am experiencing - ambiguity.

    Thank you very much!


  5. 331 Robert July 20, 2007 at 3:01 p.m.

    Now I'm liking you much more! ;)

    Here is code that takes the current rotation matrix, the previous Euler angles, and rotation order and extracts very well-behaved Euler angles from the matrix. Note that the first section of code is using David Eberly's Euler Rotation from Matrix algorithm. This is NOT a general solution (i.e.: rotations > 360d are not considered). This may be generalized later, but since this involves constrained IK on skeletal hierarchies it is a safe bet that rotations will not exceed these values.

    http://www.kuroyumes-developmentzone.com...

    Thank you extremely much!


  6. 575 Greg Abbas July 30, 2007 at 11:45 p.m.

    Hi guys,

    Sorry, I didn't your messages until now because I've been slammed with comment spam... I need to improve my counter-measures I guess. :-/

    Anyway... Jesse, yes you can convert from a matrix to a quaternion and then to Euler angles, but you can also go directly from a matrix to Euler angles (as described on the page you mention). As for minimizing the distance, you can find the nearest [θ1, φ1, ψ1] by setting θ1=θ, φ1=φ, and ψ1=ψ, solving for n1's, and rounding to the nearest integers. Do the same for θ2 and friends, and then pick whichever triplet is closer. HTH.

    Robert, you're welcome. I'm glad it eventually made sense, even if I don't understand the complexity of computer code. :-)


  7. 777 Robert August 9, 2007 at 10:16 a.m.

    Here's an interesting problem that still seems to introduce flips - at least it is my best suspect: asin(). asin() is used in Shoemake's and Eberly's code to extract the independent rotation in the matrix (assuming these rotation orders: XYZ, XZY, YXZ, YZX, ZXY, ZYX). The fundamental problem here is that asin()'s range is limited to -pi/2 to pi/2 (-90d to 90d). That means that rotations beyond this (in the II and III quadrants) are returned incorrectly in the I or IV quadrants. I don't see any code rectification for or consideration of this in the aformentioned algorithms. Why? Even Jack Crenshaw in "Math Toolkit for Real-time Programming" devotes a section to this problem by using atan2() instead. But I can't see an easy way to use atan2() here since we have a sin but no independent cos for the same rotation in the matrix. Suggestions?

    Thank you very much!


  8. 778 Greg Abbas August 10, 2007 at 9:33 a.m.

    Yup, interesting question. If we were in two dimensions, and we just wanted to extract one angle from a 2x2 matrix, then we could certainly use atan2. In three dimensions, although intuitively it seems like we should be able to use atan2 for all three angles, we can't. We can only really use it for two of them. Any algorithm that takes just a matrix and produces just one triplet of Euler angles has to "fold" at least one of the angles (e.g. throwing away two of the four quadrants). That's a side-effect of the "180-degree identity" that I mention: because there are two ways to express an orientation (not counting all of the ±2π variants), a conversion function has to pick one. There's no way for it to tell whether you're expecting φ=π or θ=ψ=π (for instance), because they're equivalent. So maybe φ gets clamped to a range like (-π/2,π/2). The best you can do is have more control over what range the extracted angles should live in (hence this article).


  9. 786 praful April 16, 2008 at 9:44 a.m.

    Hi i am having three vectors of three axis i.e. x,y,z. i.e. i am having total 9vectors, 3 for each axis. now i want to calculate angle from this vectors.
    so can anybody tell me how should i proceed? for any angle even for negative angles also. please if possible give me arithmatic solution on my email address i.e.praful.bari@ranal.com


  10. 788 Alejandro July 17, 2008 at 5:09 p.m.

    Hi...I have a little problem...Im using the Mathematica 5.0 and I dont know how can I make the Graphics with Quaternions, I want to
    rotate a vector...if somebody knows the way to do it...please give me arithmatic solution on my email address alebyte9@hotmail.com


  11. 791 Hulud February 22, 2010 at 2:39 p.m.

    Greg,

    I've browsed half of the internet to find this answer.

    Thanks man !

    Hulud


Post a comment

Your name:

Email address:   (optional. used for gravatar but not displayed.)

Website:   (optional)

Comment:   (Limited HTML markup is allowed, including a, abbr, acronym, b, blockquote, br, em, i, li, ol, p, strong, sub, super, and ul.)