Generating Random Points in a Sphere

PUBLISHED ON AUG 4, 2018

While working through Peter Shirley’s Ray Tracing in One Weekend, I encountered this:

We also need a way to pick a random point in a unit radius sphere centered at the origin. We’ll use what is usually the easiest algorithm: a rejection method. First, we pick a random point in the unit cube where ​x​, ​y​, and ​z​ all range from -1 to +1. We reject this point and try again if the point is outside the sphere. A do/while construct is perfect for that.

I thought this was fairly intriguing and wondered about other methods to generate such points. Two naive methods came to mind, one involving generating a random unit vector and scaling it randomly, and another where a random polar coordinate was converted to a cartesian coordinate. Surprisingly, these methods gave very unusual results. I investigated the causes of these.

Discarding values outside the sphere

This is the method described in the aforementioned book. x, y, and z coordinates are picked uniformly. If the point is outside the bounds of the sphere, it is discarded and picked again. This is done again and again till enough points are available.

function getPoint() {
    var d, x, y, z;
    do {
        x = Math.random() * 2.0 - 1.0;
        y = Math.random() * 2.0 - 1.0;
        z = Math.random() * 2.0 - 1.0;
        d = x*x + y*y + z*z;
    } while(d > 1.0);
    return {x: x, y: y, z: z};
}

This gives a well-distributed collection of points, but about 48% of the points chosen are discarded1. This may not be a huge problem if the number of points that need to be generated is small, but I figured there must be a more nicer way.

Picking a random direction and length

x, y and z coordinates are picked uniformly. A vector is formed with these numbers, and is then normalized. The normalized vector is then scaled by a uniform random number to get the position of the point.

function getPoint() {
    var x = Math.random() - 0.5;
    var y = Math.random() - 0.5;
    var z = Math.random() - 0.5;

    var mag = Math.sqrt(x*x + y*y + z*z);
    x /= mag; y /= mag; z /= mag;

    var d = Math.random();
    return {x: x*d, y: y*d, z: z*d};
}

While simple to reason about, this method causes a lot of points to cluster around the center. This clustering happens primarily because of our choice of d. This value is used to decide by how much the point should be displaced from the center in the direction chosen by (x,y,z). We can see that this value of d is distributed over all the axes. Therefore, it tends to become weak and displace the point less for smaller values. In a later method, we explore how taking the cube root of d instead allows us to spread the points better.

Picking a random spherical coordinate

Polar coordinates are an alternative way to describe locations on a round surfaces. For example, a point on a 2D plane can be described using two numbers:

  • The angle made by a line segment between the point and the origin
  • The distance from the point to the origin

These numbers can be converted to cartesian coordinates(x and y coordinates) and vice versa. A similar analogue exists for spheres, called spherical coordinates. A version of this frequently used is latitudes and longitudes on the earth’s surface.

function getPoint() {
    var theta = Math.random() * 2.0*Math.PI;
    var phi = Math.random() * Math.PI;
    var r = Math.random();
    var sinTheta = Math.sin(theta); var cosTheta = Math.cos(theta);
    var sinPhi = Math.sin(phi); var cosPhi = Math.cos(phi);
    var x = r * sinPhi * cosTheta;
    var y = r * sinPhi * sinTheta;
    var z = r * cosPhi;

    return {x: x, y: y, z: z};
}

This method results in a clearly visible streak through an axis of the sphere. This streak appears because the angle doesn’t linearly map to the curve. Curves are modeled by trigonometric functions, which have distinct curved areas which vary differently with respect to the input. In a later method, we tackle this problem.

Using normally distributed random numbers

This uses the method described in this post. It works by picking three normally distributed numbers and normalizing the vector resulting from these numbers. It is then scaled by the cube root of a uniformly chosen random number for the radius.

function getPoint() {
    var u = Math.random();
    var x1 = randn();
    var x2 = randn();
    var x3 = randn();

    var mag = Math.sqrt(x1*x1 + x2*x2 + x3*x3);
    x1 /= mag; x2 /= mag; x3 /= mag;

    // Math.cbrt is cube root
    var c = Math.cbrt(u);

    return {x: x1*c, y:x2*c, z:x3*c};
}

This is a well distributed, simple algorithm to place the points.

Better choice of spherical coordinates

This method chooses spherical coordinates randomly. But for the azimuthal angle, the uniform random number is chosen over the cosine of the value, instead of the angle itself. Then, the inverse of this is taken to find the azimuthal angle. This gets red of the axial streak. This works because azimuthal angle determines a point on the curved part of a semi-circle: choosing a uniformly random point on the curve is what we want. Taking the cos-inverse of this random number allows us use a uniform normal distribution to get the angle. This is explained more elaborately here. Additionally, we use the cube root of the chosen radius to prevent clumping in the center, like the previous method.

function getPoint() {
    var u = Math.random();
    var v = Math.random();
    var theta = u * 2.0 * Math.PI;
    var phi = Math.acos(2.0 * v - 1.0);
    var r = Math.cbrt(Math.random());
    var sinTheta = Math.sin(theta);
    var cosTheta = Math.cos(theta);
    var sinPhi = Math.sin(phi);
    var cosPhi = Math.cos(phi);
    var x = r * sinPhi * cosTheta;
    var y = r * sinPhi * sinTheta;
    var z = r * cosPhi;
    return {x: x, y: y, z: z};
}

  1. A sphere with radius 1 occupies a volume of (4/3)*π, which is about 4.18. A cube whose sides touch this sphere has each side of length 2, to give a volume of 8. The probability of a point chosen from a uniform random distribution in the cube being outside the sphere is (8-4.18)/8, which is 0.48. [return]