Prerequisites

Definition of terms:

• Euler Angles
• Axis Angle

Equations

heading = atan((x * y * (1 - cos(angle)) + z * sin(angle)) / (1 - (y² + z² ) * (1 - cos (angle)))
bank = atan((y * z * (1 - cos(angle)) + x * sin(angle)) / (1 - (x² + y²)* (1 - cos (angle)))
attitude = asin(- (x * z * (1 - cos(angle)) - y * sin(angle))

Code

x

Derivation of Equations

Euler angles represent 3 rotations about the x,y and z axis in some given order. We can replace any sequence of rotations by one single rotation about some axis.

We need to be very careful about using Euler Angles and it is best to work in terms of quaternions or matricies whenever we can. As explained in euler section there are different types of Euler angles and the result will depend on the Euler definition used. Also the following calculations use a lot of trig functions and therefor will use a lot of CPU time. The accuracy may be low, specially near the singularity of the Euler system being used.

With these warnings in mind here is the calculation.

Euler can be defined in terms of a quaternion as shown here.

heading = atan(2.0 * (qx*qy + qz*qw)/(qx² - qy² - qz² + qw²))
bank = atan(2.0 * (qy*qz + qx*qw)/(-qx² - qy² + qz² + qw²))
attitude = asin(-2.0 * (qx*qz - qy*qw))

This quaternion can be defined in terms of axis-angle as shown here.

qx = x * sin(angle/2);
qy = y * sin(angle/2);
qz = z * sin(angle/2);
qw = cos(angle/2)

So working out some products gives:

qx*qw = x * sin(angle/2) * cos(angle/2) = 0.5 * x * sin(angle) // using trig identity: sin(angle) = 2 sin(angle/2) cos(angle/2) from here

qy*qw = y * sin(angle/2) * cos(angle/2) = 0.5 * y * sin(angle)

qz*qw = z * sin(angle/2) * cos(angle/2) = 0.5 * z * sin(angle)

qx*qy = x * y * sin²(angle/2) = x * y * (0.5 - 0.5*cos(angle)) // using trig identity cos(angle) = 1 - 2 sin2(angle/2) from here

qy*qz = y * z * sin²(angle/2) = y * z * (0.5 - 0.5*cos(angle))

qx*qz = x * z * sin²(angle/2) = x * z * (0.5 - 0.5*cos(angle))

qx² = x² * sin²(angle/2) = x² * (0.5 - 0.5*cos(angle))

qy² = y² * sin²(angle/2) = y² * (0.5 - 0.5*cos(angle))

qz² = z² * sin²(angle/2) = z² * (0.5 - 0.5*cos(angle))

qw² =cos²(angle/2) = 0.5 * (cos(angle)+1) // using trig identity cos(2A) = 2 cos2(A) - 1 from here

qx² - qy² - qz² + qw²
= qx² + qy² + qz² + qw² - 2qy² - 2qz²
= 1 - 2qy² - 2qz²
= 1 - y² * (1 - cos(angle)) - z² * (1 - cos(angle))
= 1 - (y² + z² ) * (1 - cos (angle))

-qx² - qy² + qz² + qw²
= qx² + qy² + qz² + qw² -2qx² - 2qy²
= 1 - 2qx² - 2qy²
= 1 - x² * (1 - cos (angle)) - y² * (1 - cos (angle))
= 1 - (x² + y²)* (1 - cos (angle))

So substituting in euler equations gives:

heading = atan(2.0 * (qx*qy + qz*qw)/(qx² - qy² - qz² + qw²))
bank = atan(2.0 * (qy*qz + qx*qw)/(-qx² - qy² + qz² + qw²))
attitude = asin(-2.0 * (qx*qz - qy*qw))

heading = atan(2.0 * (x * y * (0.5 - 0.5*cos(angle)) + 0.5 * z * sin(angle)) / (1 - (y² + z² ) * (1 - cos (angle)))
bank = atan(2.0 * (y * z * (0.5 - 0.5*cos(angle)) + 0.5 * x * sin(angle)) / (1 - (x² + y²)* (1 - cos (angle)))
attitude = asin(-2.0 * (x * z * (0.5 - 0.5*cos(angle)) - 0.5 * y * sin(angle))

heading = atan((x * y * (1 - cos(angle)) + z * sin(angle)) / (1 - (y² + z² ) * (1 - cos (angle)))
bank = atan((y * z * (1 - cos(angle)) + x * sin(angle)) / (1 - (x² + y²)* (1 - cos (angle)))
attitude = asin(- (x * z * (1 - cos(angle)) - y * sin(angle))

Issues

• It is better to use atan2 function than atan function as described here.
• For issues concerning asin and atan see here.
• Most trig functions (except openGL) use radians not degrees.

Example

we take the 90 degree rotation from this:

to this:

As shown here the axis angle for this rotation is:

angle = 90 degrees
axis = 1,0,0

So using the above result:

heading = atan((x * y * (1 - cos(angle)) + z * sin(angle)) / (1 - (y² + z² ) * (1 - cos (angle)))
bank = atan((y * z * (1 - cos(angle)) + x * sin(angle)) / (1 - (x² + y²)* (1 - cos (angle)))
attitude = asin(- (x * z * (1 - cos(angle)) - y * sin(angle))

we substitute our axis angle: angle=90 degrees, x=1, y=0, z=0

heading = atan(0 / 1 )
bank = atan(1 / 0)
attitude = asin(0 )

which gives:

heading = 0 degrees
bank = 90 degrees (since the tangent of 90 degrees is infinity)
attitude = 0 degrees

So this gives the correct result, it is banking by 90 degrees, but we have to be very careful about the following issues:

• One of the intermediate steps involves infinity (division by 0) therefore most maths libraries have an atan2 function to avoid this step.
• Most maths libraries use radians instead of degrees (apart from OpenGL).