Generating a random 3D orientation for a rigid body - for a simulated die roll, a randomized camera pose, or scattering debris after a collision - needs a precise idea of what "uniform" even means for rotations: a distribution that looks the same no matter which fixed rotation you compose it with first, the Haar measure of the rotation group \(SO(3)\). The obvious approach, drawing three Euler angles uniformly, distorts orientation space exactly the way that latitude/longitude sampling distorts a sphere, as covered in Uniform Sampling on a Sphere: naive angular sampling clusters instead of covering evenly. This chapter derives the correct construction by extending that same sphere construction one dimension up, from the ordinary sphere to the 3-sphere of unit quaternions.
Live preview of the construction derived in this chapter: the arrow is a fixed reference direction; each time it arrives somewhere new, a fresh uniformly random rotation is sampled and the arrow smoothly slerps along the great circle to where that rotation carries it. Drag to rotate the view.
Quaternions and the Double Cover
Every rotation can be represented by a unit quaternion \(q=(q_0,q_1,q_2,q_3)\) with \(q_0^2+q_1^2+q_2^2+q_3^2=1\), i.e. a point on the unit sphere \(S^3\subset\mathbb{R}^4\); the algebra of composing and applying quaternions to rotate a vector is covered there and is not repeated here. What matters for sampling is that the map from \(S^3\) onto \(SO(3)\) is exactly two-to-one: \(q\) and \(-q\) produce the same rotation. Multiplying every quaternion on \(S^3\) by a fixed unit quaternion \(h\) on the left is a rotation of \(\mathbb{R}^4\), so it preserves the uniform surface measure on \(S^3\); and because the double cover is a group homomorphism, that same left-multiplication corresponds exactly to composing with a fixed rotation on the \(SO(3)\) side. So if \(q\) is uniform on \(S^3\), its image is invariant under composition with every rotation reachable this way - and since \(SO(3)\) is compact and connected, that is already every rotation. A compact group has only one probability measure with that invariance property, so the image of the uniform measure on \(S^3\) must be exactly the Haar measure on \(SO(3)\). Sampling a uniformly random rotation therefore reduces to sampling a uniformly random point on \(S^3\).
Hopf Coordinates on \(S^3\)
\(S^2\) needed one substitution to remove its area distortion; \(S^3\) needs the same idea with one extra angle. Parametrize a point on \(S^3\) by an angle \(\theta\in[0,\tfrac{\pi}{2}]\) and two independent angles \(\xi_1,\xi_2\in[0,2\pi)\):
\[ q(\theta,\xi_1,\xi_2)=(\cos\theta\cos\xi_1,\ \cos\theta\sin\xi_1,\ \sin\theta\sin\xi_2,\ \sin\theta\cos\xi_2). \]Every such \(q\) satisfies \(q_0^2+q_1^2+q_2^2+q_3^2=\cos^2\theta+\sin^2\theta=1\), so it does sit on \(S^3\); at \(\theta=0\) it traces the circle \(q_2=q_3=0\) in the \((q_0,q_1)\)-plane, and at \(\theta=\tfrac{\pi}{2}\) the circle \(q_0=q_1=0\) in the \((q_2,q_3)\)-plane - two circles in completely orthogonal planes, joined by a continuum of intermediate tori as \(\theta\) sweeps between them.
Schematic view of \(S^3\) foliated by \(\theta\) - not a literal embedding, since \(S^3\) does not fit in three dimensions: a circle of radius \(\cos\theta\) in the \(q_0q_1\)-plane (blue, \(\theta=0\)) and a circle of radius \(\sin\theta\) in the \(q_2q_3\)-plane (orange, \(\theta=\pi/2\)), joined into a torus for each fixed \(\theta\) in between (gray). Drag to orbit the camera.
To repeat the sphere's trick here, the first step is finding how a small parameter box \(d\theta\,d\xi_1\,d\xi_2\) stretches into an actual patch of volume on \(S^3\) - the 3-dimensional analogue of the area element \(dA=\sin\theta\,d\theta\,d\varphi\) found in Uniform Sampling on a Sphere. There, the stretching factor came from the cross product of the two tangent vectors \(\partial\mathbf{p}/\partial\theta\) and \(\partial\mathbf{p}/\partial\varphi\); a cross product is specific to two vectors in three dimensions, so it does not generalize directly to three tangent vectors in four dimensions. The recipe that does generalize is the Gram matrix: the matrix \(G_{ij}=\frac{\partial q}{\partial u_i}\cdot\frac{\partial q}{\partial u_j}\) of dot products between every pair of tangent vectors, whose entries measure how much the parametrization stretches each parameter direction and how much those directions leak into one another. The square root of its determinant, \(\sqrt{\det G}\), is the local volume-scaling factor - exactly what the \(\sin\theta\) in \(dA\) was doing above, just written in a form that works in any dimension. When the parametrization happens to have no cross terms, i.e. each tangent direction is already orthogonal to the others, \(G\) is diagonal and \(\sqrt{\det G}\) collapses to the plain product of each tangent vector's own length - which is exactly what happens here, as the computation below shows. Concretely, that means computing the three tangent vectors of \(q(\theta,\xi_1,\xi_2)\) and checking whether they are mutually orthogonal.
The volume element is the square root of the determinant of the metric induced by this parametrization, i.e. the Gram matrix of the three partial derivatives of \(q\). Computing all three tangent vectors,
\[ \frac{\partial q}{\partial\theta}=(-\sin\theta\cos\xi_1,\ -\sin\theta\sin\xi_1,\ \cos\theta\sin\xi_2,\ \cos\theta\cos\xi_2), \] \[ \frac{\partial q}{\partial\xi_1}=(-\cos\theta\sin\xi_1,\ \cos\theta\cos\xi_1,\ 0,\ 0), \qquad \frac{\partial q}{\partial\xi_2}=(0,\ 0,\ \sin\theta\cos\xi_2,\ -\sin\theta\sin\xi_2), \]every cross term vanishes (each pairing splits into two equal-and-opposite products of sines and cosines), and the three diagonal entries are
\[ \left\lVert\frac{\partial q}{\partial\theta}\right\rVert^2=\sin^2\theta+\cos^2\theta=1, \qquad \left\lVert\frac{\partial q}{\partial\xi_1}\right\rVert^2=\cos^2\theta, \qquad \left\lVert\frac{\partial q}{\partial\xi_2}\right\rVert^2=\sin^2\theta. \]The Gram matrix is therefore diagonal, and the volume element is the square root of its determinant:
\[ \boxed{dV=\cos\theta\sin\theta\,d\theta\,d\xi_1\,d\xi_2}. \]As a check, integrating this over the whole parameter domain gives \(\int_0^{\pi/2}\cos\theta\sin\theta\,d\theta\int_0^{2\pi}d\xi_1\int_0^{2\pi}d\xi_2=\tfrac12\cdot2\pi\cdot2\pi=2\pi^2\), exactly the known volume of the unit 3-sphere.
Sampling \(\theta\)
Because the Gram matrix is diagonal, the volume factors into a part depending only on \(\theta\) and a uniform part in \(\xi_1,\xi_2\): the two angles are independent of each other and of \(\theta\), and each is uniform on \([0,2\pi)\) exactly as \(\varphi\) was for the sphere. What remains is sampling \(\theta\) with density proportional to \(\cos\theta\sin\theta=\tfrac12\sin2\theta\) on \([0,\tfrac{\pi}{2}]\). Its cumulative distribution is
\[ F(\theta)=\int_0^\theta \sin(2t)\,dt=\frac{1-\cos(2\theta)}{2}=\sin^2\theta. \]so the substitution \(U=\sin^2\theta\) turns \(\theta\)'s awkward trigonometric density into the identity: \(U\) itself is uniform on \([0,1]\), which is precisely the same phenomenon as \(z=\cos\theta\) turning \(dA\) into a flat measure for the sphere. Reading \(\cos\theta=\sqrt{1-U}\) and \(\sin\theta=\sqrt{U}\) off that substitution avoids computing \(\theta\) itself entirely.
The Random Quaternion
Putting the three independent uniforms \(U_1,U_2,U_3\sim\mathrm{unif}(0,1)\) together gives a uniformly random point on \(S^3\), i.e. a uniformly random rotation:
\[ \boxed{ q=\left(\sqrt{1-U_1}\cos(2\pi U_2),\ \sqrt{1-U_1}\sin(2\pi U_2),\ \sqrt{U_1}\sin(2\pi U_3),\ \sqrt{U_1}\cos(2\pi U_3)\right) }. \]Swapping sine and cosine within either pair of components only relabels which point of a uniform circle is reached at \(\xi_1=0\) or \(\xi_2=0\); it is a fixed rotation of that circle and does not change that the circle is covered uniformly. Both orderings, and the one usually quoted for this construction, describe the same distribution.
Since \(q\) and \(-q\) represent the same rotation and both occur with equal probability once \(q\) is uniform on all of \(S^3\), no correction is needed for correctness. It is still convenient to fix a sign convention - for instance always flipping \(q\) to \(-q\) whenever \(q_0<0\) - purely so that every rotation is represented by a single, canonical quaternion; since the map \(q\mapsto-q\) is measure-preserving on \(S^3\), applying this flip afterwards does not disturb the uniform distribution already established.
The Same Substitution, Twice
Both this chapter's construction and the one in Uniform Sampling on a Sphere hinge on the same move. A naive angular parametrization comes with a density built from sines and cosines that the uniform distribution does not have; picking the one substitution whose derivative reproduces that exact density turns its cumulative distribution function into the identity, so that a uniform random number is already the correctly distributed parameter with no rejection and no inverse trigonometric function needed at sampling time. For \(S^2\) that substitution is \(z=\cos\theta\); one dimension up, on \(S^3\), it is \(U=\sin^2\theta\). Both are instances of the same principle applied to progressively higher dimensional spheres.
References
- Shoemake92K. Shoemake (1992) Uniform Random Rotations, in Graphics Gems III (D. Kirk, ed.), Academic Press, 124-132
- LaValle06S. M. LaValle (2006) Planning Algorithms, Cambridge University Press, Section 5.2.4
For implementing quaternion arithmetic and applying the resulting rotations to vectors in JavaScript, Quaternion.js covers the algebra this chapter builds on.