How to specify a $k$-periodic rotation matrix in $N$ dimensions that, when applied to any vector $k$ times, returns the original vector?

Question

I am interested in applying a $k$-periodic rotation matrix $M \in \mathbb{R}^{n \times n}$ to a vector $\vec{x} \in \mathbb{R}^n$, such that $M^k \vec{x} = M\vec{x}$ for an integer $k>1$. I know how to uniformly sample a random rotation matrix from the special orthogonal group (SO($n$)), but instead, I would like to sample a rotation matrix that rotates vectors by $\theta$ such that $k\theta = 2 \pi$.

How do I characterize these matrices in order to generate (or uniformly sample from the full set of) them?