The geometry of $\operatorname{PSO}(4)$ and the quaternions Question: Given a twist of the projective space, how do I find unit quaternions that represent it?
Backgroud and what do I mean: Following Conway & Smith's On Quaternions and Octonions, every pair of unitary quaternions $l,r$ gives a map $[l,r]:x \mapsto \bar{l}xr$, and a map $\ast[l,r]:x \mapsto \bar{l} \bar{x} r$, and every element of $\operatorname{GO}(4)$ can be seen as one of these two mappings (and $\operatorname{SO}(4)$ consists of those of the first type). When working projectively, the four expressions $[\pm l, \pm r]$ are all the same projective map.
On the other hand, the elements of $\operatorname{PSO}(4)$ can be thought, in a geometric way, as twists, that is, the product of two rotations on polar lines with different angles (in particular a rotation happens when one of the two angles is zero).
I would like to see, geometrically, what is the twist obtained by the product of two twists that I know (in terms of their lines). Using quaternions that would be straightforward, if I knew, given a (geometric) twist (ie. its two polar lines and the angles) how can we find the unit quaternions $l,r$ that represent it?
 A: 
I would like to see, geometrically, what is the twist obtained by the product of two twists that I know.

Even the 3D version of this is very complicated to actually visualize: given rotations by two angles around two axes, can you see what axis and angle their product has without working anything out on paper? I imagine there is a way to do this with spherical trigonometry.

[...] given a (geometric) twist (ie. its two polar lines and the angles) how can we find the unit quaternions l,r that represent it?

Suppose $R$ is a rotation in two planes, with ordered orthonormal bases $\{w,x\}$ and $\{y,z\}$, by angles $\alpha$ and $\beta$, for which $\{w,x,y,z\}$ is ordered to match the orientation of 4D space.
If we multiply by $\overline{w}$ we get the bases $\{1,\overline{w}x\}$ and $\{\overline{w}y,\overline{w}z\}$. Note all non-$1$s in these bases are pure imaginary, so are square roots of negative one, so left multiplying by $\exp(\theta\overline{w}x)$ is a rotation by $\theta$ in both planes in the direction afforded by the ordered bases, and right multiplying by $\exp(\theta\overline{w}x)$ does the same except the opposite direction in the second plane.
Thus, we can left-multiply by $\overline{w}$ to turn $\{w,x,y,z\}$ into $\{1,\overline{w}x,\overline{w}y,\overline{w}z\}$, then left-multiply by $\exp(\frac{\alpha+\beta}{2}\overline{w}x)$ and right-multiply by $\exp(\frac{\alpha-\beta}{2}\overline{w}x)$ to get a rotation by $\alpha$ in the $\{1,\overline{w}x\}$-plane and by $\beta$ in the $\{\overline{w}y,\overline{w}z\}$-plane, then right-multiply by $w$ to turn $\{1,\overline{w}x,\overline{w}y,\overline{w}z\}$ into $\{w,x,y,z\}$. The net effect is rotation by $\alpha$ and $\beta$ respectively in the $\{w,x\}$- and $\{y,z\}$-planes. Since $w\exp(\frac{\alpha+\beta}{2}\overline{w}x)\overline{w}=\exp(\frac{\alpha+\beta}{2}x\overline{w})$. Thus $R$ is
$$ R(p)=\exp\Big(\frac{\alpha+\beta}{2}x\overline{w}\Big)\,p\,\exp\Big(\frac{\alpha-\beta}{2}\overline{w}x\Big). $$
