Quaternion-based orientation: UAV Systems Engineering (M5-T4)

Learning objectives

After this lesson, you should be able to:

Define a unit quaternion from axis and angle and state the double cover.
Compose rotations and rotate a vector with quaternions.
Integrate angular rate into a quaternion and renormalise.
Convert to and from a rotation matrix and verify against the matrix exponential.

Conceptual explanation

A rotation by angle $\theta$ about the unit axis $n$ is the unit quaternion $q=[\cos (\theta /2),\,n\,\sin (\theta /2)]$ (scalar part first, Hamilton convention). Quaternions compose by the Hamilton product, and the rotation of a vector $v$ is $q$ ⊗ [0, $v$ ] ⊗ q*. The quaternions $q$ and − $q$ describe the same rotation (a double cover), which matters when you compare or average them.

Axis-angle meaning of a unit quaternion.

Convention warning. Libraries differ: scalar-first [ $w$ , $x$ , $y$ , $z$ ] (PX4 messages, Eigen’s constructor) versus scalar-last [ $x$ , $y$ , $z$ , $w$ ] (ROS geometry_msgs, Eigen’s coefficient order, ROS TF); Hamilton versus JPL multiplication order; active versus passive. Write the convention next to every quaternion you store.

Animation of quaternion attitude integration for one minute at a constant angular rate. A first-order update without renormalisation lets the quaternion norm grow steadily above one, so a unit vector rotated by it is stretched instead of just turned. A first-order update followed by renormalisation keeps the norm at exactly one, and the exact exponential-map update never leaves the unit sphere. Two strips show the norm error and the length of a rotated unit vector. — $q$ ← $q$ + ½ $q\otimes \omega$ dt is only first-order accurate and does not preserve $|q|=1$ : the norm grows by about (ω dt) $^{2}/2$ per step, and rotating a vector with a non-unit $q$ also scales it by $|q|^{2}$ . Renormalising removes the norm error; the exponential map $q$ ← $q\otimes \exp$ (ω dt/2) is exact for constant $\omega$ . $\omega =2\,\mathrm{rad}/\mathrm{s}$ , dt = 20 ms, 60 s. Computed. (Looping, about 6 seconds; drawn and computed for this course with illustrative values.)

Engineering intuition

Euler angles describe an orientation as three consecutive turns and break down when two axes align (gimbal lock). A quaternion has no such failure: any orientation is a point on the 3-sphere in four dimensions, and small rotations are small steps on that sphere. The price is one constraint $(\Vert q\Vert =1)$ that numerical integration slowly violates, so you renormalise.

Mathematical formulation

Quaternion algebra (Hamilton, scalar first)

\mathbf q = [w,\ \mathbf v], \qquad w = \cos\tfrac{\theta}{2}, \qquad \mathbf v = \mathbf n\sin\tfrac{\theta}{2}, \qquad \lVert\mathbf q\rVert = 1

Product

\mathbf q_1\otimes\mathbf q_2 = \left[\,w_1 w_2 - \mathbf v_1\cdot\mathbf v_2,\ \ w_1\mathbf v_2 + w_2\mathbf v_1 + \mathbf v_1\times\mathbf v_2\,\right]

Rotation of a vector

[0,\ \mathbf v'] = \mathbf q\otimes[0,\ \mathbf v]\otimes\mathbf q^{*}, \qquad \mathbf q^{*} = [w,\ -\mathbf v]

Integration of angular rate

With body angular rate $\omega$ (rad/s, in the body frame) and $q$ the rotation from body to world:

Quaternion kinematics and its exact discrete step

\dot{\mathbf q} = \tfrac{1}{2}\,\mathbf q\otimes[0,\ \boldsymbol\omega]

Exact step

\mathbf q_{k+1} = \mathbf q_k\otimes\left[\cos\frac{\lVert\boldsymbol\omega\rVert\Delta t}{2},\ \ \frac{\boldsymbol\omega}{\lVert\boldsymbol\omega\rVert}\sin\frac{\lVert\boldsymbol\omega\rVert\Delta t}{2}\right]

Exact over a step Δt with constant ω.

The first-order form $q$ + ½ $q\otimes [0,\,\omega ]\Delta t$ followed by normalisation is adequate at high rate; the exact step is as cheap and avoids the integration error for large $\omega \Delta t$ . Renormalise after each step: $q$ ← $q/\Vert q\Vert$ .

To a rotation matrix

Rotation matrix of a unit quaternion

\mathbf R(\mathbf q) = \begin{bmatrix}1 - 2(y^{2}+z^{2}) & 2(xy - wz) & 2(xz + wy)\\ 2(xy + wz) & 1 - 2(x^{2}+z^{2}) & 2(yz - wx)\\ 2(xz - wy) & 2(yz + wx) & 1 - 2(x^{2}+y^{2})\end{bmatrix}

Attitude error for feedback is best formed as a rotation, $q_{e}=q_{\text{ref}}$ * ⊗ $q$ , with the sign chosen so that $w\ge 0$ (take the short way round), and the small-error vector part ≈ half the rotation vector.

Architecture and data flow

The IMU gives $\omega$ ; the attitude filter integrates it into $q$ and corrects it with accelerometer (gravity), magnetometer or vision. The control stack receives $q$ (and its derived $R$ ), forms the attitude error, and outputs body rates or torques. Conventions change at three boundaries: PX4 (FRD to NED, scalar first), ROS (FLU to ENU, scalar last in messages), and your own code. Convert once, at the boundary (lessons M1-T4, M6-T10).

Real-world UAV applications

Attitude estimation and control in every autopilot.
Error-state filters that estimate a small rotation correction (lesson M3-T6).
Camera pose and IMU preintegration in VIO (lesson M3-T8).
Smooth interpolation between orientations in trajectory generation.

Implementation example

This runs as shown: quaternion rotation and composition are checked against rotation matrices, the exact integration step is checked against the matrix exponential, and renormalisation is shown to matter.

import numpy as np

def qmul(a, b):
    w1, x1, y1, z1 = a; w2, x2, y2, z2 = b
    return np.array([w1*w2 - x1*x2 - y1*y2 - z1*z2, w1*x2 + x1*w2 + y1*z2 - z1*y2,
                     w1*y2 - x1*z2 + y1*w2 + z1*x2, w1*z2 + x1*y2 - y1*x2 + z1*w2])
def qconj(q): return np.array([q[0], -q[1], -q[2], -q[3]])
def axis_angle(n, th): n = np.asarray(n, float) / np.linalg.norm(n); return np.concatenate([[np.cos(th / 2)], n * np.sin(th / 2)])
def to_R(q):
    w, x, y, z = q
    return np.array([[1-2*(y*y+z*z), 2*(x*y-w*z), 2*(x*z+w*y)], [2*(x*y+w*z), 1-2*(x*x+z*z), 2*(y*z-w*x)], [2*(x*z-w*y), 2*(y*z+w*x), 1-2*(x*x+y*y)]])
def rotate(q, v): return qmul(qmul(q, np.concatenate([[0.0], v])), qconj(q))[1:]
def expm_so3(w):
    th = np.linalg.norm(w)
    if th < 1e-12: return np.eye(3)
    k = w / th; K = np.array([[0, -k[2], k[1]], [k[2], 0, -k[0]], [-k[1], k[0], 0]])
    return np.eye(3) + np.sin(th) * K + (1 - np.cos(th)) * K @ K

# 1) rotating a vector: quaternion and matrix agree; 90 deg about z maps x to y
q = axis_angle([0, 0, 1], np.pi / 2)
assert np.allclose(rotate(q, [1, 0, 0]), [0, 1, 0]) and np.allclose(to_R(q) @ [1, 0, 0], [0, 1, 0])
# 2) composition: q1 then q2 (q2 applied last, in the world frame) equals R2 @ R1
q1, q2 = axis_angle([1, 0, 0], 0.4), axis_angle([0, 1, 1], -0.9)
assert np.allclose(to_R(qmul(q2, q1)), to_R(q2) @ to_R(q1))
# 3) double cover: q and -q are the same rotation
assert np.allclose(to_R(q1), to_R(-q1))
# 4) integrate a constant body rate: exact step vs matrix exponential
w = np.array([0.3, -0.5, 0.8]); dt, steps = 0.01, 500
qe = np.array([1.0, 0, 0, 0]); qf = qe.copy()
for _ in range(steps):
    th = np.linalg.norm(w) * dt
    qe = qmul(qe, np.concatenate([[np.cos(th / 2)], w / np.linalg.norm(w) * np.sin(th / 2)]))       # exact step
    qf = qf + 0.5 * qmul(qf, np.concatenate([[0.0], w])) * dt                                       # first-order, NOT renormalised
R_exact = expm_so3(w * dt * steps)
assert np.allclose(to_R(qe), R_exact, atol=1e-12)
print("exact step matches expm; first-order norm drift:", round(abs(np.linalg.norm(qf) - 1), 5))
assert abs(np.linalg.norm(qf) - 1) > 1e-4 and abs(np.linalg.norm(qe) - 1) < 1e-12                     # un-normalised first-order step drifts off the unit sphere
assert np.allclose(to_R(qf / np.linalg.norm(qf)), R_exact, atol=5e-3)                                  # after normalising it is still accurate

Quaternion algebra checked against rotation matrices and the matrix exponential.

Common engineering mistakes

Mixing scalar-first and scalar-last quaternions

The rotation silently becomes a different one; label the convention and convert at the boundary.

Not renormalising

The norm drifts and the rotation matrix stops being orthogonal.

Multiplying in the wrong order

Hamilton product is not commutative; verify against rotation matrices in a test.

Comparing q and −q directly

They are the same rotation; compare with the sign fixed or use the rotation angle of $q_{1}$ * ⊗ $q_{2}$ .

Averaging quaternions component-wise with mixed signs

Align signs first, or use a proper rotation average.

Failure modes

Failure modes, symptoms and causes
Failure mode	How it shows up	Likely cause
Vehicle flips on frame change	Attitude appears inverted after a conversion	Scalar order or inverse (body-to-world vs world-to-body) confusion
Attitude drift in the estimator	Slowly rotating estimate with a still vehicle	Missing normalisation or integration error at large rates
Spin the long way round	Controller commands a $350^{\circ}$ turn instead of − $10^{\circ}$	Error quaternion not sign-corrected
Jump in interpolated orientation	Sudden flip between samples	Opposite-sign quaternions interpolated linearly

Debugging approach

Write unit tests like the ones above for every new convention.
Log the quaternion norm: it should stay within 1e-6 of 1.
Check the rotation of a known vector (forward, down) at known attitudes against expectation.
Compare the attitude from your code with the autopilot’s on the same log.

Practical exercises

Show that $q$ ⊗ q* = [1, 0, 0, 0] for a unit quaternion.
Convert the ZYX Euler angles (roll $10^{\circ}$ , pitch $20^{\circ}$ , yaw $30^{\circ}$ ) to a quaternion and verify with the matrix of lesson M1-T4.
Implement the error quaternion with the short-way correction and the small-angle rotation-vector estimate.
Integrate a 360°/s rate for one second with first-order and exact steps at $\Delta t=0.05\,\mathrm{s}$ and compare.

Interview questions

Write an answer in your own words, then reveal the model answer. Nothing is sent anywhere.

Question 1

Why use quaternions instead of Euler angles?

Model answer

Euler angles have singularities (gimbal lock) and non-uniform interpolation; quaternions represent any attitude with four numbers and one constraint, compose cheaply and integrate without singularities.

Question 2

Why must you renormalise a quaternion after integration?

Model answer

Numerical integration violates the unit-norm constraint; without renormalisation the quaternion stops representing a pure rotation and the derived matrix is not orthogonal.

Question 3

What is the double cover?

Model answer

$q$ and − $q$ represent the same rotation, so comparisons and averaging must account for the sign.

References and further reading

Primary papers, official documentation and standards. Each claim above that is specific to a tool or version points here.

J. Solà, “Quaternion kinematics for the error-state Kalman filter”, arXiv:1711.02508
$T$ . $D$ . Barfoot, State Estimation for Robotics, Cambridge University Press, 2017, chapter 6.
REP 103 (ROS conventions)

What’s next

Mark this lesson

Saved in this browser and synced to your account when you are signed in.

Next lesson

M5-T5: Attitude control

Describe the cascaded attitude controller (angle loop outside a rate loop) and the limits that set its bandwidth.