3D: avoid pinching at poles when creating sphere from polar coordinates

Question

I'm using Wikipedia's spherical coordinate system article to create a sphere made out of particles in Three.js. Based on this article, I created a small Polarizer class that takes in polar coordinates with setPolar(rho, theta, phi) and it returns its corresponding x, y, z

Here's the setPolar() function:

// Rho: radius
// theta θ: polar angle on Y axis
// phi φ: azimuthal angle on Z axis

Polarizer.prototype.setPolar = function(rho, theta, phi){
  // Limit values to zero
  this.rho = Math.max(0, rho);
  this.theta = Math.max(0, theta);
  this.phi = Math.max(0, phi);

  // Calculate x,y,z
  this.x = this.rho * Math.sin(this.theta) * Math.sin(this.phi);
  this.y = this.rho * Math.cos(this.theta);
  this.z = this.rho * Math.sin(this.theta) * Math.cos(this.phi);

  return this;
}

I'm using it to position my particles as follows:

var tempPolarizer = new Polarizer();

for(var i = 0; i < geometry.vertices.length; i++){
  tempPolarizer.setPolar(
    50,                         // Radius of 50
    Math.random() * Math.PI,    // Theta ranges from 0 - PI
    Math.random() * 2 * Math.PI // Phi ranges from 0 - 2PI
  );

  // Set new vertex positions
  geometry.vertices[i].set(
    tempPolarizer.x,
    tempPolarizer.y,
    tempPolarizer.z
  );
}

It works wonderfully, except that I'm getting high particle densities, or "pinching" at the poles:

I'm stumped as to how to avoid this from happening. I thought of passing a weighted random number to the latitude, but I'm hoping to animate the particles without the longitude also slowing down and bunching up at the poles.

Is there a different formula to generate a sphere where the poles don't get as much weight? Should I be using quaternions instead?

Answer 1

In order to avoid high density at the poles, I had to lower the likelihood of theta (latitude) landing close to 0 and PI. My input of

Math.random() * Math.PI, for theta gives an equal likelihood to all values (orange).

Math.acos((Math.random() * 2) - 1) perfectly weights the output to make 0 and PI less likely along the sphere's surface (yellow)

Now I can't even tell where the poles are!