77
  • How do I generate $1000$ points $\left(x, y, z\right)$ and make sure they land on a sphere whose center is $\left(0, 0, 0\right)$ and its diameter is $20$ ?.
  • Simply, how do I manipulate a point's coordinates so that the point lies on the sphere's "surface" ?.
Felix Marin
  • 94,079
Filip
  • 779
  • 1
  • 6
  • 5

6 Answers6

103

Using Gaussian distribution for all three coordinates of your point will ensure an uniform distribution on the surface of the sphere. You should proceed as follows

  1. Generate three random numbers $x, y, z$ using Gaussian distribution
  2. Multiply each number by $1/\sqrt{x^2+y^2+z^2}$ (a.k.a. Normalise) . You should handle what happens if $x=y=z=0$.
  3. Multiply each number by the radius of your sphere.
Community
  • 1
Henricus V.
  • 19,100
  • 2
    Excellent! This works in any dimension and seems not to be widely known, although I've seen some older paper (late $'40$'s) where the authors claim that they've learned it from Harald Cramér. – orangeskid Dec 22 '15 at 22:43
  • 34
    Not widely known? Every probabilist should know that the multivariate standard normal distribution is spherically symmetric, from which this follows immediately. – Robert Israel Dec 23 '15 at 00:40
  • 2
    Does it have to be a Gaussian distribution? Does the standard deviation of the distribution matter? What if I were to use a uniform distribution? – NicNic8 Jul 06 '19 at 03:19
  • 1
    @NicNic8 It has to be spherically symmetric, so for example the density (before adjusting for the magnitude) at $(1.25,0,0,0)$ has to be the same as the density at $(0,0.75,1,0)$ and a multivariate Gaussian distribution with covariance matrix which is a constant times the identity matrix achieves this. The standard deviation does not matter so long as it is positive. A uniform distribution on say $[0,1]$ would not therefore not work – Henry Jul 06 '19 at 11:04
  • @Henry's comment is correct. Intuitively, you can reason about this in 2D. Using the above technique, any sample that you take that occurs along the same radial line will be mapped to the same point on the perimeter of the unit circle. Since the Gaussian has an identity covariance (which can be scaled by a constant), the area along any radial line from 0 to infinity is the same. Note that this is not true if the Gaussian is not symmetric. Hence, each point on the perimeter of the circle has an equal probability of being generated. Therefore, samples on the perimeter are uniformly distributed. – Vivek Subramanian Sep 21 '19 at 18:01
  • Why is a Gaussian distribution used for the 3 coordinates, instead of, say a uniform distribution? – user5965026 Jun 29 '20 at 22:05
  • @user5965026 because if you use a uniform distribution you'll get biases in the corners. Compare it to a 2D image, if you splat out coordinates based on a uniform distribution they'll neatly fill a square. If you then transform that square into a circular shape by normalizing all the vectors you will get all the points that were outside the circle squashed into the bits of the circle closest to the corners and your distribution isn't uniform within the circle (or on the circumference). – alkanen Jan 10 '25 at 09:08
93

Use the fact that if you cut a sphere of a given radius with two parallel planes, the area of the strip of spherical surface between the planes depends only on the distance between the planes, not on where they cut the sphere. Thus, you can get a uniform distribution on the surface using two uniformly distributed random variables:

  • a $z$-coordinate, which in your case should be chosen between $-10$ and $10$; and
  • an angle in $[0,2\pi)$ corresponding to a longitude.

From those it’s straightforward to generate the $x$- and $y$-coordinates.

Brian M. Scott
  • 631,399
  • 20
    This always challenges my intuition. – copper.hat Dec 22 '15 at 21:17
  • 2
    @copper.hat: Yes, I’ve always found it surprising that the flattening and stretching exactly balance the shrinking radius as the planes move towards a pole. – Brian M. Scott Dec 22 '15 at 21:22
  • 10
    So you should get a uniform distribution on the cylinder and project on the sphere, like Archimedes. Nice! – orangeskid Dec 22 '15 at 22:37
  • @orangeskid: Yep. – Brian M. Scott Dec 22 '15 at 22:39
  • @BrianM.Scott Do you have a reference for the geometrical theorem you mentioned? – Jack M Dec 22 '15 at 23:59
  • 2
    @Jack: Not off the top of my head; it’s a fairly standard second semester calculus exercise when you get to surface areas of revolution. Aha! I found one. – Brian M. Scott Dec 23 '15 at 00:00
  • 9
    This proves, by the way, that the surface area of a sphere's circumscribing cylinder (minus the endcaps) equals that of the sphere itself (and by corollary, the volume of the sphere is one-third the radius times the common surface area). – Brian Tung Dec 23 '15 at 00:11
  • 1
    Demonstrating the surface area equality is a nice little exercise in cosine cancellation. – Brian Tung Dec 23 '15 at 00:17
  • That the surface area of a segment of a given sphere is proportional to the thickness of the segment has been amazing mathematicians since Archimedes. :) – PM 2Ring Dec 23 '15 at 12:16
  • 8
    Does this generalize to higher dimensions? For a 2D circle you just pick $\theta$ uniformly. For a 3D sphere you choose a $\theta$ and $z$. Is there a similar result for 4D and higher? – Meni Rosenfeld Dec 23 '15 at 13:07
  • Isn't this actually called Archimedes's theorem? – Chappers Dec 23 '15 at 13:36
  • 17
    While this method is simple and clever, it can suffer from numerical accuracy issues near the poles (where ${\rm d}x/{\rm d}z \to \infty$). For example, let's assume that your random number generator internally generates a 32-bit integer $n$, which is then scaled to the desired interval as $z = 20n/2^{32}-10$. Thus, the minimum distance between $z$ values is $20/2^{32}\approx5\times10^{-9}$, while near the poles, the minimum distance between points along a line of longitude is $10\arccos(1-2/2^{32})\approx0.0003$. For a sphere the size of the Earth, that distance would be about 200 meters. – Ilmari Karonen Dec 23 '15 at 19:04
23

Here is a simple but less efficient way:

Generate points uniformly $x \in [-10,10]^3$ and reject if $\|x\| =0 $ (which should rarely happen) or $\|x\| > 10$ (which should happen with probability ${20^3 -{4 \over 3} \pi 10^3 \over 20^3} =1 - {\pi \over 6} \approx 48\%$). Otherwise let $y = {10 \over \|x\|} x$. Then $y$ will be distributed uniformly on the surface of the $10$-sphere.

copper.hat
  • 178,207
  • 2
    The points that are not rejected are uniformly distributed in the sphere, and projecting such points to the surface preserves uniformity. – copper.hat Dec 22 '15 at 22:48
  • This method is probably faster than the others from a computationnal point of view. The only complex operations are a division and a square root, and both are quite fast – Tryss Dec 23 '15 at 21:34
  • @Tryss Computationally, yes, but probabilistically the rejection rate is not to be ignored, especially in high dimensions (cf. curse of dimensionality). There are also quite fast implementations of standard normal random variates. It may be worth doing some simulations although I would be surprised if no one has written a paper about it. – heropup Dec 23 '15 at 22:11
  • 1
    I think for 1000 points the cost is minimal, and the coding fairly simple. – copper.hat Dec 23 '15 at 22:13
  • 3
    @heropup : yes, in high dimensions the rejection rate is big, and the gaussian method scale lineary contrary to this method. But in low dimensions, it works very well – Tryss Dec 23 '15 at 22:30
  • 4
    For those implementing this algorithm, you can improve performance by removing the square root operation and comparing the squared length of vector x with the squared radius of the circle. – adrian Aug 26 '22 at 22:50
21

In addition to Brian Scott's excellent and clever answer, here's another, more straightforward way (in case you want to approach it with a geographical intuition): From two random variables $u_1, u_2$, distributed uniformly on the interval $[0, 1]$, generate (in radians) the latitude

$$ \phi = \arccos (2u_1-1)-\frac{\pi}{2} $$

and the longitude

$$ \lambda = 2\pi u_2 $$

Then compute the rectangular coordinates accordingly:

$$ x = \cos\phi\cos\lambda $$ $$ y = \cos\phi\sin\lambda $$ $$ z = \sin\phi $$

ETA (thanks to Tanner Strunk—see comments): This will give coordinates of points on the unit sphere. To have them land on the sphere with diameter $20$ (and therefore radius $10$), simply multiply each by $10$.

ETA (2023-10-11 thanks to Silverfish): I had originally swapped the usual symbols for latitude and longitude; this is now fixed.

Brian Tung
  • 35,584
  • 1
    This is the best way. Some points: 1. arccos(-1..1) normally gives the range 0..pi, so I would be inclined to subtract pi/2 instead of pi, to get a number in the range -pi/2 .. pi/2 which is how a geographer would present the angle. Actually I think this change is required, because the current answer will always give negative z. 2. Having done that, I see nothing wrong with z = 2u -1 and cos lamba = sqrt(1-z^2) 3. phi works just as well with 2npi*u2, which may be convenient in some computer implementations. – Level River St Dec 24 '15 at 11:06
  • Oops, I think you are right. Thanks for the look out. – Brian Tung Dec 24 '15 at 15:41
  • 1
    I'm a bit late in coming to this question (ha), but you should be scaling those x,y,z expressions by 10 to land on the sphere, right? (Tiny, nit-picky detail, but I figure it should be noted for future readers.) Nice avatar by the way, Brian. – Tanner Strunk Jan 11 '18 at 17:33
  • 1
    @TannerStrunk: Yes, you're quite right, thanks! I'll add an edit. – Brian Tung Jan 11 '18 at 18:35
  • 2
    Note that the latitude generation can be simplified to $\lambda = \arcsin(2a - 1)$. – Jeff G Oct 28 '19 at 23:57
  • I think this is using a nonstandard spherical coordinate system. At least from the ones I'm familiar with the "latitude" angle is given as an angle with respect to the polar axis (from the Z axis down, see this illustration), but this is using the angle from the X/Y plane up towards the polar axis (Z axis). I'm guessing it's still correct, but it confused me since Wolfram had a slightly different formula. – Max Strater Feb 03 '23 at 08:34
  • 1
    @MaxStrater The conventional definition of latitude is measured from the fundamental (equatorial) plane... but you'd be more used to it if you ever worked with geographic coordinate systems. There are in fact many different ways to define spherical coordinates and they don't all use a polar angle. But you're right this answer is nonstandard in one respect: conventionally, latitude is $\phi$ and longitude is $\lambda$!! – Silverfish Oct 11 '23 at 10:08
  • 1
    @Silverfish: Dang, did I switch them around? I'll fix that in a bit... – Brian Tung Oct 11 '23 at 16:10
  • This is the same as Brian Scott's answer, no? – John Clements Dec 24 '23 at 20:28
  • @JohnClements: They're based on the same idea, yes (which is also useful for determining the surface area of the sphere). I don't recall if we were working on our posts at more or less the same time, and then I hastily added a reference to his answer when I realized he had posted too. Maybe...too long ago! :-) I do pass it through longitude and latitude first, but it's fair to say that's not a fundamental distinction. – Brian Tung Dec 25 '23 at 00:19
10

Same way as on a real sphere, but $(x,y,z) $ multiplied by $i.$

Narasimham
  • 42,260
3

Wolfram Mathworld provides a methodology for randomly picking a point on a sphere:

To obtain points such that any small area on the sphere is expected to contain the same number of points, choose $u$ and $ν$ to be random variates on $[0,1]$. Then: $$\begin{array}{ll}\theta=2\pi u\\ \varphi= arccos(2v - 1)\end{array}$$ gives the spherical coordinates for a set of points which are uniformly distributed over $\mathbb{S}^2$.

Daniel Fischer
  • 211,575
varunk
  • 139