To make a break from the more theoretical posts, I’ll describe how to simulate a homogeneous Poisson point process on the surface of a sphere. I have already simulated this point process on a rectangle, disk and triangle. As with those shapes, we need to randomly position points in a uniform manner. This is not too complicated for the surface of a sphere, as it’s just an extension of randomly positioning points on a disk. I recommend reading that post first, as a lot of the material presented here builds off it.

## Steps

Simulating a Poisson point process requires two steps: simulating the random number of points and then randomly positioning each point.

**Number of points**

The number of points of a Poisson point process on the surface of sphere circle of radius \(r>0\) is a Poisson random variable with mean \(\lambda A\), where \(A=4\pi r^2\) is the area of the sphere. As long as your preferred programming language can produce (pseudo-)random numbers according to a Poisson distribution, you can simulate a homogeneous Poisson point process on a sphere.

To do this in MATLAB, use the poissrnd function with the argument \(\lambda A\). In Python, we can use either the scipy.stats.poisson or numpy.random.poisson function from the SciPy or NumPy libraries.

(If you’re curious about how Poisson simulation works, I suggest seeing this post for details on simulating or sampling Poisson random variables or, more accurately, *variates*.)

**Locations of points**

To position all the points randomly and uniformly on a sphere, we need to work in spherical coordinates. I’ll denote the two angular coordinates with \(\theta\) and \(\theta\), which respectively represent the polar and azimuth angles. These variable are angles on a sphere, meaning \(0\leq \theta \leq \pi\) and \(0 \leq \phi\leq 2\pi\). The radial distance is just \(\rho\), which is a fixed constance for positioning points on the surface.

To generate the random angular values, we simply produce uniform random variables between zero and one. We then multiply all these two sets of numbers respectively by \(\pi\), and \(2\pi\).

(If you wanted to sample the points uniformly **in** the sphere, you would need to generate uniform variables on the unit interval \((0,1)\), take their cubic roots, and then, multiply them by the radius \(r\). This is akin to the step of taking the square root in the disk setting.)

##### From polar to Cartesian coordinates

We have generated spherical coordinates for points randomly and uniformly located on the surface of a disk. But to plot the points, we have to convert coordinates back to Cartesian form. This is easily done in MATLAB with the sph2cart function.

In languages without such a function, we just use the standard change of coordinates: \(x=\rho\sin(\theta)\cos(\phi)\), \(y=\rho\sin(\theta)\sin(\phi)\) and \(y=\rho\cos(\theta)\).

##### Equal x and y axes

Sometimes the plotted points more resemble points on an ellipse than a disk due to the different scaling of the x and y axes. To fix this in MATLAB, run the command: axis square. In Python, it’s quite difficult to do this, as it seems to lack an automatic function.

## Results

I have presented the results produced by code written in MATLAB and Python. The blue points are the Poisson points on the sphere. I have used a surface plot (with clear faces) to illustrate the underling sphere in black.

##### MATLAB

##### Python

Note: The aspect ratio in 3-D Python plots tends to squash the sphere slightly, but it is a sphere.

## Code

The code for all my posts is located online here. For this post, the code in MATLAB and Python is here. In Python I used the library mpl_toolkits for doing 3-D plots.

## Further reading

The material in this post is just a natural extension of simulating points uniformly on a disk, which is covered in the third edition of the classic book *Stochastic Geometry and its Applications* by Chiu, Stoyan, Kendall and Mecke details on page 54 how to uniformly place points on a disk, which they call the *radial way*. But I wouldn’t consult this book for simulation methods in general.

I found this post on a similar topic.

Here is some sample Python code for creating a 3-D scatter plot.