"Random" generation of rotation matrices

I had this question and I solved it via QR decomposition of 3x3 matrices with normally distributed random numbers. If you have matlab, use: [Q,R] = qr(randn(3));

Q is a uniformly random rotation matrix. For more info, refer to: M. Ozols, “How to generate a random unitary matrix,” 2009.

Rotation matrices can be uniquely defined by a vector and a rotation angle. To generate the vector, you can use grandom spherical coordinates $\phi$ and $\theta$. Thus you should first generate random angles by using:

$$\theta = \arccos(2u_1 - 1)$$ $$\phi = 2\pi u_2$$ Where $u_1,u_2$ are uniformly distributed in $[0,1]$. This will give you a vector around which to rotate. Then, randomly decide the amount of rotation $\psi\in[0,2\pi]$.

I have had the same exact problem myself a while ago so I point you to this which says it very succinctly with plenty of references.