As Nathan Reed and teodron exposed, the recipe for rotating a vector v by a unit-length quaternion q is:
- Create a pure quaternion p out of v. This simply means adding a fourth coordinate of 0:
$$p = (v_x, v_y, v_z, 0) \Leftrightarrow p = (v, 0)$$
- Pre-multiply it with q and post-multiply it with the conjugate q*:
$$p' = q \times p \times q*$$
- This will result in another pure quaternion which can be turned back to a vector:
$$v' = (p'_x, p'_y, p'_z)$$
This vector \$v'\$ is \$v\$ rotated by \$q\$.
This is working but far from optimal. Quaternion multiplications mean tons and tons of operations. I was curious about various implementations such as this one, and decided to find from where those came. Here are my findings.
We can also describe q as the combination of a 3-dimensional vector u and a scalar s:
$$q = (u_x, u_y, u_z, s) \Leftrightarrow q = (u, s)$$
By the rules of quaternion multiplication, and as the conjugate of a unit length quaternion is simply it's inverse, we get:
$$
\begin{align}
p' &= qpq* \\
&= (u, s)(v,0)(-u, s) \\
&= (sv + u \times v, -u \cdot v)(-u, s) \\
&= ((-u \cdot v)(-u) + s(sv + u \times\ v) + (sv + u \times v) \times (-u), \dots) \\
&= ((u \cdot v) u + s^2v + s(u \times v) + sv \times (-u) + (u \times v) \times (-u), \dots)
\end{align}
$$
The scalar part (ellipses) results in zero, as detailed here. What's interesting is the vector part, AKA our rotated vector v'. It can be simplified using some basic vector identities:
$$
\begin{align}
v' &= (u \cdot v)u + s^2v + s(u \times v) + s(u \times v) + u \times (u \times v) \\
&= (u \cdot v)u + s^2v + 2s(u \times v) + (u \cdot v)u - (u \cdot u)v \\
&= 2(u \cdot v)u + (s^2 - u \cdot u)v + 2s(u \times v)
\end{align}
$$
This is now much more optimal; two dot products, a cross product and a few extras: around half the operations. Which would give something like that in source code (assuming some generic vector math library):
void rotate_vector_by_quaternion(const Vector3& v, const Quaternion& q, Vector3& vprime)
{
// Extract the vector part of the quaternion
Vector3 u(q.x, q.y, q.z);
// Extract the scalar part of the quaternion
float s = q.w;
// Do the math
vprime = 2.0f * dot(u, v) * u
+ (s*s - dot(u, u)) * v
+ 2.0f * s * cross(u, v);
}
