Uniformly Distributed Random Points Inside a Regular Polygon

I will end this series of posts on generating uniformly random points inside geometric shapes by considering regular polygons. Regular polygons are equilateral and equiangular. When the number of sides is equal to 4 we have a square (it is fairly easy to generate uniformly distributed random points inside a square). As the number of sides goes to infinity we have a circle (I already showed how to generate uniformly distributed random points inside a circle here).

For all the other regular polygons with \(n\) sides you could use the technique presented in this paper to generate uniformly distributed random points. The paper is by Myron Hlynka and Deborah Loach from the Department of Mathematics and Statistics at University of Windsor in Canada.

Uniformly Distributed Random Points Inside a Circular Ring

I will continue my series of posts on generating uniformly distributed points within different shapes by considering a circular ring (a.k.a. annulus or doughnut). Just like my previous post we have the joint probability density function (PDF) of the \(x\) and \(y\) coordinate of the random points as:

$$ f_{X,Y}(x,y)=
\begin{cases}
\frac{1}{A}=\frac{1}{\pi (R_{c_2}^2-R_{c_1}^2)} & R_{c_1}^2 \leq x^2+y^2 \leq R_{c_2}^2 \\
0 & \text{otherwise}
\end{cases},$$
where \(R_{c_2}\) is the radius of the outer circle and \(R_{c_1}\) is the radius of the inner circle, and \(R_{c_1}<R_{c_2}\)

Using a similar technique as in the previous post we have:
$$ f_{R,\Theta}(r,\theta) = \frac{1}{2\pi} \times \frac{2r}{R_{c_2}^2-R_{c_1}^2} = f_\Theta(\theta)f_R(r) $$
where
$$f_\Theta(\theta) = \frac{1}{2\pi} \text{ for } 0 \leq \theta \leq 2\pi$$
and
$$f_R(r)=\frac{2r}{R_{c_2}^2-R_{c_1}^2} \text{ for } R_{c_1} \leq r \leq R_{c_2}.$$
  
Therefore, \(\Theta\) is uniformly distributed between \(0\) and \(2\pi\). The random variable \(R\) can be generated by first calculating its cumulative distribution function as
$$ F_R(r) = \int_{R_{c_1}}^r \frac{2\alpha}{R_{c_2}^2-R_{c_1}^2}d\alpha = \frac{r^2-R_{c_1}^2}{R_{c_2}^2-R_{c_1}^2}, $$
and then using a uniformly distributed random variable \(U\) over the interval \([0,1]\) to get
$$ U = \frac{R^2-R_{c_1}^2}{R_{c_2}^2-R_{c_1}^2} \implies R = \sqrt{(R_{c_2}^2-R_{c_1}^2)U + R_{c_1}^2}. $$

The following MATLAB code generates 1000 random numbers inside a circular ring with outer radius \(20\), and inner radius \(10\) centered at \(-30\), \(-40\).

%********************************************** 

n = 10000;
Rc2 = 20;
Rc1 = 10;
Xc = -30;
Yc = -40;

theta = rand(1,n)*(2*pi);
r = sqrt((Rc2^2-Rc1^2)*rand(1,n)+Rc1^2);
x = Xc + r.*cos(theta);
y = Yc + r.*sin(theta);

plot(x,y,'.'); axis square

%********************************************** 

Uniformly Distributed Random Points Inside a Circle (2)

So much for posting the answer in a few days! Sorry, I have been Very Very busy.

So here is the answer. Lets assume we have the Cartesian coordinate system with (\(x,y\)) representing the points. Then joint probability distribution function (PDF) of our random points inside the circle (i.e. the joint distribution of random variable \(X\) and random variable \(Y\) representing the \(x\) and \(y\) coordinate of the random point) is given by:

$$ f_{X,Y}(x,y)=
\begin{cases}
\frac{1}{A}=\frac{1}{\pi R_c^2} & x^2+y^2 \leq R_c^2 \\
0 & \text{otherwise}
\end{cases},$$
where \(R_c\) is the radius of the circle.

Now lets consider the polar coordinate system where \(r=g_1(x,y)= \sqrt{x^2+y^2}\) and \(\theta=g_2(x,y)=\arctan (y/x)\). Also, from inverting the functions \(g_1\) and \(g_2\) we have \(x=h_1(r,\theta)=r\cos \theta\) and \(y=h_2(r,\theta)=r\sin\theta\). Now calculating the Jacobian of the transformation we have:

$$ J(x,y)= \det
\begin{bmatrix}
\frac{\partial r}{\partial x}=\frac{x}{\sqrt{x^2+y^2}} & \frac{\partial r}{\partial y}=\frac{y}{\sqrt{x^2+y^2}} \\
\frac{\partial \theta}{\partial x}=\frac{-y}{x^2+y^2} & \frac{\partial \theta}{\partial y}=\frac{x}{x^2+y^2}
\end{bmatrix}= \frac{1}{\sqrt{x^2+y^2}}=\frac{1}{r}.$$
 Also we have
 $$ |J(r,\theta)|=\frac{1}{|J(x,y)|}$$

The joint PDF of random variables \(R\) and \(\Theta\) can then be calculated using the joint PDF of \(X\) and \(Y\) as:

$$ f_{R,\Theta}(r,\theta)= f_{X,Y}\left( h_1(r,\theta),h_2(r,\theta) \right) =
\begin{cases}
\frac{r}{\pi R_c^2} & 0 \leq \theta \leq 2\pi, 0 \leq r \leq R_c \\
0 & \text{otherwise}
\end{cases}.$$

We can rewrite the non zero part as:
$$ f_{R,\Theta}(r,\theta) = \frac{1}{2\pi} \times \frac{2r}{R_c^2} = f_\Theta(\theta)f_R(r) $$
where
$$f_\Theta(\theta) = \frac{1}{2\pi}$$
and
$$f_R(r)=\frac{2r}{R_c^2}.$$
Therefore, \(\Theta\) is uniformly distributed between \(0\) and \(2\pi\). The random variable \(R\) can be generated by first calculating its cumulative distribution function as
$$ F_R(r) = \int_0^r \frac{2\alpha}{R_c^2}d\alpha = \frac{r^2}{R_c^2}, $$
and then using a uniformly distributed random variable \(U\) over the interval \([0,1]\) to get
$$ U = \frac{R^2}{R_c^2} \implies R=R_c \sqrt{U}.$$

The following MATLAB code generates 1000 random numbers inside a circle with radius \(20\) centered at \(-30\), \(-40\).

%********************************************** 

n = 10000;
Rc = 20;
Xc = -30;
Yc = -40;

theta = rand(1,n)*(2*pi);
r = Rc*sqrt(rand(1,n));
x = Xc + r.*cos(theta);
y = Yc + r.*sin(theta);

plot(x,y,'.'); axis square;

%********************************************** 

If you don't understand how I did the PDF transformation see

Probability, Statistics, and Random Processes For Electrical Engineering (3rd Edition)

Chapter 6 (6.2.3 to be exact) or

Probability and Random Processes for Electrical and Computer Engineers

Chapter 8 (8.3 to be exact).

I hope this helps someone out there.