Equations
first method:
w = c1 c2 c3 - s1 s2 s3
x = s1 s2 c3 +c1 c2 s3
y = s1 c2 c3 + c1 s2 s3
z = c1 s2 c3 - s1 c2 s3
where:
- c1 = cos(heading / 2)
- c2 = cos(attitude / 2)
- c3 = cos(bank / 2)
- s1 = sin(heading / 2)
- s2 = sin(attitude / 2)
- s3 = sin(bank / 2)
An alternative form is:
w = Math.sqrt(1.0 + C1 * C2 + C1*C3
- S1 * S2 * S3 + C2*C3)
/ 2
x = (C2 * S3 + C1 * S3 + S1
* S2 * C3) / (4.0 * w)
y = (S1 * C2 + S1 * C3 + C1
* S2 * S3) / (4.0 * w)
z = (-S1 * S3 + C1 * S2 * C3
+ S2) /(4.0 * w)
where:
- C1 = cos(heading)
- C2 = cos(attitude)
- C3 = cos(bank)
- S1 = sin(heading)
- S2 = sin(attitude)
- S3 = sin(bank)
note: in the second form the angles are not divided by 2. I don't know which
of these forms is most stable? However, as William
points out the first is better because it requires the same number of trig operations,
no square root, no worry about dividing by zero, uses familiar formulae, and
is fairly clearly normalised.
Code
First method
public final void rotate(double heading, double attitude, double bank) {
// Assuming the angles are in radians.
double c1 = Math.cos(heading/2);
double s1 = Math.sin(heading/2);
double c2 = Math.cos(attitude/2);
double s2 = Math.sin(attitude/2);
double c3 = Math.cos(bank/2);
double s3 = Math.sin(bank/2);
double c1c2 = c1*c2;
double s1s2 = s1*s2;
w =c1c2*c3 - s1s2*s3;
x =c1c2*s3 + s1s2*c3;
y =s1*c2*c3 + c1*s2*s3;
z =c1*s2*c3 - s1*c2*s3;
}
Second method
public final void rotate(double heading, double attitude, double bank) {
// Assuming the angles are in radians.
double c1 = Math.cos(heading);
double s1 = Math.sin(heading);
double c2 = Math.cos(attitude);
double s2 = Math.sin(attitude);
double c3 = Math.cos(bank);
double s3 = Math.sin(bank);
w = Math.sqrt(1.0 + c1 * c2 + c1*c3 - s1 * s2 * s3 + c2*c3) / 2.0;
double w4 = (4.0 * w);
x = (c2 * s3 + c1 * s3 + s1 * s2 * c3) / w4 ;
y = (s1 * c2 + s1 * c3 + c1 * s2 * s3) / w4 ;
z = (-s1 * s3 + c1 * s2 * c3 +s2) / w4 ;
}
Derivation of Equations
This is derived from a method originally sent to me by Andy:
The quaternion for the rotation by angle a about unit vector (x1,y1,z1) is
given by:
cos(angle/2) + i ( x1 * sin(angle/2)) + j (y1 * sin(angle/2)) + k ( z1 * sin(angle/2))
Therefore:
if h = heading angle (rotation about y) |
then Qh = quaternion for pure heading rotation = cos(h/2) + j sin(h/2) = c1
+ j s1 |
if a = attitude angle (rotation about z) |
then Qa = quaternion for pure attitude rotation = cos(a/2) + k sin(a/2) = c2
+ k s2 |
if b = bank angle (rotation about x) |
then Qb = quaternion for pure bank rotation = cos(b/2) + i sin(b/2) = c3 +
i s3 |
where:
- c1 = cos(h / 2)
- c2 = cos(a / 2)
- c3 = cos(b / 2)
- s1 = sin(h / 2)
- s2 = sin(a / 2)
- s3 = sin(b / 2)
The required quaternion can be calculated by multiplying these individual quaternions
From our definitions the order of applying these rotations is heading,attitude then bank (about y,z then x). As we saw on this page the rotation applied first goes on the right hand side of the equation but since we are working in the frame of reference of the moving object the first rotation goes on the left. Applying heading then attitude gives: (Qh * Qa), then applying bank gives:
Q = (Qh * Qa)* Qb
Q = ((c1 + j s1)) * (c2 + k s2))* (c3) + i s3))
Q = (c1 c2 + j s1 c2 + k c1 s2 + j k s1 s2)* (c3 + i s3)
but jk=i which gives:
Q = (c1 c2 + i s1 s2 + j s1 c2 + k c1 s2) * (c3 + i s3)
Q = c1 c2 c3 + i s1 s2 c3 + j s1 c2 c3 + k c1 s2 c3 + i s3 c1 c2 + i i s1 s2
s3 + j i s1 c2 s3 + k i c1 s2 s3
but ii=-1 and j i = -k and k i = j which gives:
Q = c1 c2 c3 + i s1 s2 c3 + j s1 c2 c3 + k c1 s2 c3 + i c1 c2 s3 + i i s1 s2
s3 + j i s1 c2 s3 + k i c1 s2 s3
Q = (c1 c2 c3 - s1 s2 s3) + i (s1 s2 c3 +c1 c2 s3) + j (s1 c2 c3 + c1 s2 s3)
+ k (c1 s2 c3 - s1 c2 s3)
if we define Q = w + i x + j y + k z
then
w = c1 c2 c3 - s1 s2 s3
x = s1 s2 c3 +c1 c2 s3
y = s1 c2 c3 + c1 s2 s3
z = c1 s2 c3 - s1 c2 s3
William has proved here that the second method
is equivalent to the first.
Issues
It is easier to convert from euler angles to quaternions than the reverse direction,
so once you have converted to quaternions it is best to stay in that form.
If you have a different result from that shown on this page it may be that
you are using different standards, I have tried to keep the standards consistent
across this site and I have tried to define the standards that I am using here.
One of these standards is that the order that the euler rotations are applied
is 'NASA standard aeroplane' it uses a slightly different coordinate definition
from VRML (z and y axis swapped). Also when working with aeroplanes we often
work in terms of the position of external objects relative to the aircraft (i.e.
the inverse of its position transform as explained
here). therefore to get the expression normally used with NASA aeroplane
we invert all inputs (change sign of every term with odd number of sine terms)
invert output (conjugate quaternion) swap z and y for inputs (swap z and y columns)
swap z and y for outputs (swap z and y rows). In the case of aircraft it can
make sense to work in terms of the local frame of reference of the aircraft
looking out, I have a version of this page that does that here,
but be careful as this will no longer be compatible with the rest of this site.
However if you have checked these things and you still have a discrepancy then
I have probably made an error so please let me know
here.
Example
we take the 90 degree rotation from this: |
|
to this: |
|
As shown here the axis angle for this rotation
is:
heading = 0 degrees
bank = 90 degrees
attitude = 0 degrees
so substituting this in the above formula gives:
- c1 = cos(heading / 2) = 1
- c2 = cos(attitude / 2) = 1
- c3 = cos(bank / 2) = 0.7071
- s1 = sin(heading / 2) = 0
- s2 = sin(attitude / 2) = 0
- s3 = sin(bank / 2) = 0.7071
w = c1 c2 c3 - s1 s2 s3 = 0.7071
x = s1 s2 c3 +c1 c2 s3 = 0.7071
y = s1 c2 c3 + c1 s2 s3 = 0
z = c1 s2 c3 - s1 c2 s3 = 0
which gives the quaternion 0.7071 + i 0.7071
This agrees with the quaternions here.
Angle Calculator and Further examples
Enter the heading, attitude and bank, in degrees, below. Then click calculate
Be careful: very low numbers close to zero are sometimes shown in exponential form like: 6.123031769111886e-17 (I don't know enough about javascript to stop it doing this). If you see a number greater than one then check for the exponent.
If you need to have the option to show angles in radians and a graphical representation I have put a java applet here which allows the values to be entered and the converted values shown along with a graphical representation of the orientation.
Also further examples in 90 degree steps here
This site may have errors. Don't use for critical systems.