My goal is to uniformly sample the distribution of all unit quaternions $q$ that have angular distance $\theta$ relative to a given quaternion orientation $q_1 = [s_1, q_{1,1}, q_{1,2},q_{1,3}]$.
I calculate angular distance between two unit quaternions $q_1$ and $q_2$ as follows (this is from this post):
$$ \theta = \arccos(2\langle q_1,q_2 \rangle^2 – 1) \tag{1}$$
where
$$ \langle q_1,q_2 \rangle = s_1s_2+q_{1,1}q_{2,1}+q_{1,2}q_{2,2}+q_{1,3}q_{2,3} \tag{2} $$
In addition, I know from here that all random unit quaternions describing 3D orientations can be addressed by uniformly sampling $u_1,u_2,u_3 \in [0,1]$ and computing the quaternion
$$ q_x = (\sqrt{1-u_1}\sin2\pi u_2,\sqrt{1-u_1}\cos2\pi u_2,\sqrt{u_1}\sin2\pi u_3,\sqrt{u_1}\cos2\pi u_3) \tag{3}$$
I can combine these equations algebraically to get
$$ \sqrt{\frac{\cos(\theta)+1}{2}} = \langle q_1,q_x \rangle \tag{4}$$
Where $q_x$ is as defined in $(3)$, $q_1$ is the guiding orientation, $\theta$ is the relative angle, and $\langle q_1,q_x \rangle$ is as defined in $(2)$.
I am struggling with how to choose parameters $u_1,u_2,u_3$ to satisfy $(4)$, especially given the need for uniform sampling of the orientation space. The expanded equation $(4)$ is not trivial to simplify. Any help would be greatly appreciated.
Best Answer
Rather than trying to directly choose randomly from unit quaternions that are distributed at a fixed angle from $q_1$, why don't you try the easier task of picking those that are distributed at a fixed angle from the identity quaternion, and then displace them to be distributed around $q_1$, by acting on $q_1$ with them?
To have a unit quaternion $d$ at an angle $\theta$ from the identity $e = [1, 0, 0, 0]$, we require that $$\langle e , d \rangle = 1\cdot d_0 = d_0 = \sqrt\frac{1+\cos\theta}{2}$$
This gives no further constraints on $d_1$, $d_2$, and $d_3$, leaving only the unit constraint, which means you should take $(d_1, d_2, d_3)$ to be uniformly distributed on a sphere of radius $\sqrt\frac{1-\cos\theta}{2}$ (so that they add up to 1 when combined with the fixed $d_0$).
Generating a random vector on the unit sphere can be done by noting that any one direction is uniformly distributed (conventionally taken to be $z$), and then taking the other two in a random planar direction, appropriately normalized. In other words take two random numbers $u_1$, and $u_2$ in the unit interval [0, 1] and generate
$$z = 2u_1 - 1, x = \sqrt{1-z^2}\cos(2 \pi u_2), y = \sqrt{1-z^2} \sin(2 \pi u_2).$$
Multiply these by $\sqrt\frac{1-\cos\theta}{2}$ to get your $d$ components.
Finally, take the $d$ you've generated, and apply it to your starting quaternion $q_x$ to get a new unit quaternion at a displacement $\theta$:
$$ q_x = d q_1 $$