Random point on hyperspherical cap

A point distribution is uniform if the density is everywhere proportional to the "area". If you take spherical coordinates in $n$ dimension, $(r,\phi_1,\ldots,\phi_{n-1})$, wikipedia says that the spherical volume element is $$\mathrm d^nV = r^{n-1}\sin^{n-2}(\phi_1)\sin^{n-3}(\phi_2)\cdots\sin(\phi_{n-2})\mathrm dr\mathrm d\phi_1\mathrm d\phi_2\cdots\mathrm d\phi_{n-1}$$ from which we deduce that the spherical surface element is $$\mathrm d^{n-1}S = \mathrm d^nV/\mathrm dr = r^{n-1}\prod_{i=1}^{n-1}\sin^{n-1-i}(\phi_i)\mathrm d\phi_i$$

EDIT: The above can be rewritten as $$\mathrm d^{n-1}S = \left(r^{n-2}\prod_{i=2}^{n-1}\sin^{n-1-i}(\phi_i)\mathrm d\phi_i\right) \times\left(r\sin^{n-2}(\phi_1)\mathrm d\phi_1\right) = r\sin^{n-2}(\phi_1)\mathrm d\phi_1\mathrm d^{n-2}S$$ thus as Chris Jones justly remarked in the comments, it's sufficient to sample $\phi_1$ from a distribution proportional to $\sin^{n-2}(x)$, and sample a point from a lower dimensional sphere uniformly, to reach the desired result.

From there, we thus have to pick the spherical coordinate $\phi_i$ from a pseudo-distribution whose pdf is $\sin^{n-1-i}$. Pseudo-distribution because a proper pdf must integrate to 1. To get the actual pdf, we have to normalize the function, so it'd be $$\mathrm{pdf}_i(x) = \frac{\sin^{n-1-i}(x)}{\int_{\phi_{i,\min}}^{\phi_{i,\max}}\sin^{n-1-i}(t)\mathrm dt}$$ For a spherical cap that admits $Ox_1$ as its axis, $\phi_1$ should range from $0$ to the angular opening of the cap, $\arccos\left(\frac{\vec c.\vec e}{r^2}\right)$, $\phi_{n-1}$ ranges from $0$ to $2\pi$, and all the remaining $\phi_i$ range from $0$ to $\pi$.

If you sample every $\phi_i$ with the pdf's above, you will end up with a point in the spherical cap, and the distribution will be uniform.

---

As of currently, I don't see a simple way to sample these random variables like for the $2$-sphere. The "trick" of uniformly sampling a random variable $z$ on the range $[r\cos\theta,r]$ only works because it is basically an application of inverse transform sampling, where everything is easy to compute.

I don't know what language or software you are using, but in all likelihood, you should be able to find tools that can sample an arbitrary random variable if you provide it with the pdf. I have no clue how fast/efficient these methods are.