Random 3D points uniformly distributed on an ellipse shaped window of a sphere

How can I generate random points uniformly distributed on the surface of a sphere such that a line that originates at the center of the sphere, and passes through one of the points, will intersect a plane within a circle. Following are more illustrations and details to make this clearer. The white dots are the points, the red is the sphere, the black circle is the circle in the plane, the green is the elliptical cone enclosing the points:

Why not Point-rejection?

To generate random points uniformly distributed on the surface of a sphere, this works (pseudo Matlab code; where $n$ is the number of the points):

z = rand(n)*2-1
θ = rand(n)*2π
x = sqrt(1. - z.*z).*cos(θ)
y = sqrt(1. - z.*z).*sin(θ)

To generate random points uniformly distributed on the surface of a sphere within a circle shaped window, the following works:

α = .9
cosα = cos(α)
z = rand(n)*(cosα - 1.) - cosα
θ = rand(n)*2π
x = sqrt(1. - z.*z).*cos(θ)
y = sqrt(1. - z.*z).*sin(θ)

Notice that this algorithm is much more efficient than first generating points uniformly distributed on a sphere (using the previous algorithm) and then rejecting all the points that are outside the circle shaped window. I'm hoping for some mathematical-wizardry to answer my question, and in the same way it does not make any sense to use a point-rejection method for the circle shaped window problem, I hope to avoid it here too.

You might try a point-rejection approach that rejects a relatively small percentage of the points.

As observed in comments, you can improve on sampling over the entire sphere by finding a circular cap that barely covers the elliptical window, and using the method you have already devised to sample uniformly within that cap. This works well for elliptical windows that are nearly circular. For elliptical windows whose "length" is much greater than their "width," you could refine the method in various ways, two of which I will try to describe.

Alternative Method 1

Find the line through the center of the sphere and the center of the elliptical "window," and construct several circles on the sphere with centers on that line. The smallest circle could have a radius slightly larger than the smaller semi-diameter of the ellipse, and the largest circle could have a radius equal to the larger semi-diameter of the ellipse.

Between each adjacent pair of circles is a spherical zone. If you use the line through the center of the elliptical window as the axis of a set of spherical coordinates, the zone lies between two lines of latitude. Now divide the zone into smaller regions along lines of longitude so that two of the smaller regions cover the intersection of the elliptical window and the zone. Once you have done this for all zones, the elliptical window will be completely contained within the union of the spherical cap bounded by the smallest circle and the selected parts of all the zones. The union of the cap and the selected parts of the zones will be your sampling region. By changing the number of circles and their radii and the longitudes of the selected regions within each zone, you can reduce the ratio between the area of the sampling region and the area of the elliptical window to be as close to $1$ as you want.

Now select a point uniformly distributed within the sampling region. To do this, you can take transformed Cartesian coordinates $x',y',z'$ with origin at the center of the sphere such that the $z'$ axis passes through the center of the elliptical window. Randomly select an appropriately distributed $z'$ coordinate using the inverse cumulative distribution of $z'$ within the sampling region, which is a piecewise linear function. Then, depending on which zone that puts the randomly generated point, you randomly select a longitude within the applicable ranges of longitude and find the $x'$ and $y'$ coordinates. Determine whether the random point is within the elliptical window; if not, reject it and try again.

In practice there is probably some number of zones that balances the extra effort of evaluating a piecewise function with that many pieces against the extra effort of generating random coordinates that will be rejected.

Alternative Method 2

A simpler approach is to orient the transformed coordinates $x',y',z'$ so that the $z'$ axis is parallel to the semi-minor axis of the elliptical window. Using the $z'$ axis as the axis of a set of spherical coordinates, the elliptical window lies between two lines of latitude and two lines of longitude. You can generate the $z'$ coordinate uniformly between its minimum and maximum values, and the $x'$ and $y'$ coordinates by generating a longitude uniformly between the two bounding longitudes. I think this means that more than $3/4$ of the random points generated will fall within the elliptical window, so fewer than $1/4$ are rejected.

You can improve the rejection ratio even further by dividing the sampling region into strips by latitude and including just enough difference in longitude in each strip to cover the elliptical window, but this may not be worthwhile.

Conclusion

I would probably implement Alternative Method 2 (without any additional subdivision of the sampling region) unless the center of the circle in the plane was so close to the coordinates $(0,0)$ (from the figures in the question) that the rejection rate would be better for a spherical cap that contained the elliptical window, and in that case I would use the cap.