Power Cosine Weighted Hemisphere Cap

The probability of each sample is weighted by a power of its cosine with the normal. Samples are chosen from a cap of the hemisphere, defined by a maximum elevation angle.

Here is the result:

function pow_cos_hemisphere_cap(theta_max,n){

    const u_1 = Math.random();
    const u_2 = Math.random();

    const theta = Math.acos(Math.pow(
        1.0 - u_1 * (1.0 - Math.pow(Math.cos(theta_max),n+1)), 
        1.0/(n+1)
    ));
    const phi = u_2 * 2.0 * Math.PI;

    return [theta,phi];
}

Algorithm

Generate points weighted by a power cosine with power $n$ on the sphere cap defined by the maximum elevation angle $\theta_{max}$

Choose uniform numbers $u_1, u_2$ in $[0,1]$
Calculate $\theta = \cos^{-1}(\sqrt[n+1]{1 - u_1(1- \cos^{n+1}\theta_{max})} )$
Calculate $\phi = 2\pi u_2$

Optionally convert to cartesian coordinates with:

x = \sin\theta\cos\phi

y = \sin\theta\sin\phi

z = \cos\theta

PDF Used In Monte Carlo Integration

p(\omega) = \frac{n+1}{2\pi(1 - \cos^{n+1}\theta_{max}) }\cos^n\theta

p(\theta,\phi) = \frac{n+1}{2\pi(1 - \cos^{n+1}\theta_{max}) }\cos^n\theta\sin\theta

Derivation

This is a special case of the power cosine weighted hemisphere sector with $\phi_{min}=0, \phi_{max}=2\pi, \theta_{min} = 0$