Discover MakerZone

MATLAB and Simulink resources for Arduino, LEGO, and Raspberry Pi

Learn more

Discover what MATLAB® can do for your career.

Opportunities for recent engineering grads.

Apply Today

Thread Subject:
Random points in a circular zone

Subject: Random points in a circular zone

From: Valentina

Date: 23 May, 2012 08:42:06

Message: 1 of 13

Hi, I want to create a function that give me N random points in a circular zone.

The circular zone can be represented as
A={(x,y)=(R*cos(theta), R*t*sin(theta)), t in [-1, 1], theta in [a,b]}
where 0<a<b<pi.

I have written a function for the circular segment, that is when a=0, but this function doesn't give me a uniform distribution and I don't understand how apply it to the zone.

function [x,y]= RandSeg(N, R, w)
    angle = -w+2*w*rand(N,1); % theta is in [0, w]
    k=linspace(-1,1,N)';
    x=R.*cos(angle);
    y=R.*k.*sin(angle);

Can anyone help me?
Thank you!
Valentina

Subject: Random points in a circular zone

From: Peter Bone

Date: 23 May, 2012 14:53:06

Message: 2 of 13

"Valentina" wrote in message <jpi7su$g55$1@newscl01ah.mathworks.com>...
> Hi, I want to create a function that give me N random points in a circular zone.
>
> The circular zone can be represented as
> A={(x,y)=(R*cos(theta), R*t*sin(theta)), t in [-1, 1], theta in [a,b]}
> where 0<a<b<pi.
>
> I have written a function for the circular segment, that is when a=0, but this function doesn't give me a uniform distribution and I don't understand how apply it to the zone.
>
> function [x,y]= RandSeg(N, R, w)
> angle = -w+2*w*rand(N,1); % theta is in [0, w]
> k=linspace(-1,1,N)';
> x=R.*cos(angle);
> y=R.*k.*sin(angle);
>
> Can anyone help me?
> Thank you!
> Valentina

How about:

d = 2*R;
r2 = R*R;
x = zeros(1,N);
y = x;
count = 1;
while count <= N
  xp = d*(rand-0.5);
  yp = d*(rand-0.5);
  if (xp*xp + yp*yp) < r2
    x(count) = xp;
    y(count) = yp;
    count = count + 1;
  end
end

Not sure it can be vectorised.

Subject: Random points in a circular zone

From: Peter Bone

Date: 23 May, 2012 15:18:07

Message: 3 of 13

"Valentina" wrote in message <jpi7su$g55$1@newscl01ah.mathworks.com>...
> Hi, I want to create a function that give me N random points in a circular zone.
>
> The circular zone can be represented as
> A={(x,y)=(R*cos(theta), R*t*sin(theta)), t in [-1, 1], theta in [a,b]}
> where 0<a<b<pi.
>
> I have written a function for the circular segment, that is when a=0, but this function doesn't give me a uniform distribution and I don't understand how apply it to the zone.
>
> function [x,y]= RandSeg(N, R, w)
> angle = -w+2*w*rand(N,1); % theta is in [0, w]
> k=linspace(-1,1,N)';
> x=R.*cos(angle);
> y=R.*k.*sin(angle);
>
> Can anyone help me?
> Thank you!
> Valentina

Actually, here's a better vectorised method. The square root provides the uniform distribution.

function [x,y]= RandSeg(N, R, w)
  angle = -w+2*w*rand(N,1);
  rad = rand(N,1).^0.5;
  x=R*rad.*cos(angle);
  y=R*rad.*sin(angle);

Subject: Random points in a circular zone

From: Bruno Luong

Date: 23 May, 2012 15:24:05

Message: 4 of 13

I though similar question was answered in a previous thread:

R = 7;
a = pi/3;
b = pi;

n=1e4;
w = b-a;
r=R*sqrt(rand(1,n));
theta = a + rand(1,n)*w;
x=r.*cos(theta);
y=r.*sin(theta);

plot(x,y,'.')
axis equal

% Bruno

Subject: Random points in a circular zone

From: Valentina

Date: 23 May, 2012 16:26:05

Message: 5 of 13

"Bruno Luong" <b.luong@fogale.findmycountry> wrote in message <jpivel$ra5$1@newscl01ah.mathworks.com>...
> I though similar question was answered in a previous thread:
>
> R = 7;
> a = pi/3;
> b = pi;
>
> n=1e4;
> w = b-a;
> r=R*sqrt(rand(1,n));
> theta = a + rand(1,n)*w;
> x=r.*cos(theta);
> y=r.*sin(theta);
>
> plot(x,y,'.')
> axis equal
>
> % Bruno

Hi Bruno, your function generates n random points in a circular sector, but unfortunately I want a function that generates them in a circular zone and I don't understand how to change it.
Valentina

Subject: Random points in a circular zone

From: Bruno Luong

Date: 23 May, 2012 16:52:05

Message: 6 of 13

"Valentina" wrote in message <jpj32t$el4$1@newscl01ah.mathworks.com>...
> Hi Bruno, your function generates n random points in a circular sector, but unfortunately I want a function that generates them in a circular zone and I don't understand how to change it.

How do you define a circular sector (I do not understand you first post)? In the mean time try this

% Radius min and max
Rmin = 2;
Rmax = 5;
% Angle min and max
a = pi/3;
b = pi;

n=1e4;
w = b-a;
r = sqrt(rand(1,n)*(Rmax^2-Rmin^2)+Rmin^2);
theta = a + rand(1,n)*w;
x = r.*cos(theta);
y = r.*sin(theta);

plot(x,y,'.')
axis equal

% Bruno

Subject: Random points in a circular zone

From: Valentina

Date: 23 May, 2012 17:17:09

Message: 7 of 13

"Bruno Luong" <b.luong@fogale.findmycountry> wrote in message <jpj4jl$lb3$1@newscl01ah.mathworks.com>...
> How do you define a circular sector (I do not understand you first post)? In the mean time try this
>
> % Radius min and max
> Rmin = 2;
> Rmax = 5;
> % Angle min and max
> a = pi/3;
> b = pi;
>
> n=1e4;
> w = b-a;
> r = sqrt(rand(1,n)*(Rmax^2-Rmin^2)+Rmin^2);
> theta = a + rand(1,n)*w;
> x = r.*cos(theta);
> y = r.*sin(theta);
>
> plot(x,y,'.')
> axis equal
>
> % Bruno

Your last function gives random points in an anular sector, but I want a circular zone. A circular zone is the portion of a disk included between any two parallel chords and their intercepted arcs and I found that it can be represented as I write in my first post.
Valentina

Subject: Random points in a circular zone

From: Bruno Luong

Date: 23 May, 2012 17:39:06

Message: 8 of 13

I see.

I have the feeling that there is no close form formula that can generate unform distribution on a sector.

Bruno

Subject: Random points in a circular zone

From: Alan_Weiss

Date: 23 May, 2012 19:12:37

Message: 9 of 13

On 5/23/2012 1:39 PM, Bruno Luong wrote:
> I see.
>
> I have the feeling that there is no close form formula that can
> generate unform distribution on a sector.
>
> Bruno

An inefficient but theoretically correct scheme is known as rejection.
You generate a bunch of points that are uniformly distributed on a
region that is larger than the region you have in mind, and then throw
away all points that are outside the region you have in mind.

So, for example, you could use Bruno's method of generating points in a
disk, and then remove the points that are not in your zone. The result
will have points uniformly distributed within your zone.

The problem, of course, is you are never sure beforehand how many points
to generate so that you keep the number you are aiming for. In other
words, when you use rejection, you have to generate more points than you
are aiming for, and you don't know how many more.

Alan Weiss
MATLAB mathematical toolbox documentation

Subject: Random points in a circular zone

From: ImageAnalyst

Date: 23 May, 2012 11:54:01

Message: 10 of 13

FYI: Original thread: http://www.mathworks.com/matlabcentral/newsreader/view_thread/319832#876000

Subject: Random points in a circular zone

From: Roger Stafford

Date: 24 May, 2012 07:56:09

Message: 11 of 13

"Valentina" wrote in message <jpi7su$g55$1@newscl01ah.mathworks.com>...
> Hi, I want to create a function that give me N random points in a circular zone.
> The circular zone can be represented as
> A={(x,y)=(R*cos(theta), R*t*sin(theta)), t in [-1, 1], theta in [a,b]}
> where 0<a<b<pi.
- - - - - - - - - -
  The following will give a statistically uniform distribution area-wise in your circular zone without rejection. Let R, a, and b be as defined in your original post. (I have called your t parameter s and your theta parameter t.)

 n = 8000; % Select how many random points you want to generate

 p = (2*a-sin(2*a))+(2*b-sin(2*b)-2*a+sin(2*a))*rand(n,1);
 p = mod(p+pi,2*pi)-pi;
 p = sign(p).*abs(p).^(1/3);
 c1 = 1/(6^(1/3)); % c1 and c2 used in approx. of derivative
 c2 = (c1-2/3/pi^(2/3))/pi^2;
 t = zeros(n,1);
 for k = 1:12 % Twelve trips should suffice for convergence
  p1 = t-sin(t);
  p1 = sign(p1).*abs(p1).^(1/3);
  t = t - (p1-p)./(c1-c2*t.^2); % Pseudo Newton Raphson
 end
 t = mod(t,2*pi)/2;
 s = -1+2*rand(n,1);
 x = R*cos(t); % x and y are cartesian coordinates of the random points
 y = R*s.*sin(t);

Roger Stafford

Subject: Random points in a circular zone

From: Valentina

Date: 24 May, 2012 09:36:15

Message: 12 of 13

> The following will give a statistically uniform distribution area-wise in your circular zone without rejection. Let R, a, and b be as defined in your original post. (I have called your t parameter s and your theta parameter t.)
>
> n = 8000; % Select how many random points you want to generate
>
> p = (2*a-sin(2*a))+(2*b-sin(2*b)-2*a+sin(2*a))*rand(n,1);
> p = mod(p+pi,2*pi)-pi;
> p = sign(p).*abs(p).^(1/3);
> c1 = 1/(6^(1/3)); % c1 and c2 used in approx. of derivative
> c2 = (c1-2/3/pi^(2/3))/pi^2;
> t = zeros(n,1);
> for k = 1:12 % Twelve trips should suffice for convergence
> p1 = t-sin(t);
> p1 = sign(p1).*abs(p1).^(1/3);
> t = t - (p1-p)./(c1-c2*t.^2); % Pseudo Newton Raphson
> end
> t = mod(t,2*pi)/2;
> s = -1+2*rand(n,1);
> x = R*cos(t); % x and y are cartesian coordinates of the random points
> y = R*s.*sin(t);
>
> Roger Stafford

Thank you so much Roger, it's exactly what I wanted.
Valentina

Subject: Random points in a circular zone

From: Bruno Luong

Date: 24 May, 2012 10:25:19

Message: 13 of 13

n = 8000;
a = pi/6;
b = pi/2+pi/6;
R = 7;
I think rejection method could be the fastest in average.

xmin = R*cos(b);
xmax = R*cos(a);
ymin = -R;
ymax = R;

xy = zeros(0,2);
while true
    M = 2*n - size(xy,1);
    x = xmin + (xmax-xmin)*rand(M,1);
    y = ymin + (ymax-ymin)*rand(M,1);
    in = x.^2 + y.^2 <= R^2;
    xy = cat(1, xy, [x(in) y(in)]);
    if size(xy,1) >= n
        xy(n+1:end,:) = [];
        break
    end
end

plot(xy(:,1),xy(:,2),'.');
axis equal

% Bruno

Tags for this Thread

What are tags?

A tag is like a keyword or category label associated with each thread. Tags make it easier for you to find threads of interest.

Anyone can tag a thread. Tags are public and visible to everyone.

Contact us