scatter randomized points with circles

Hello,

How we can draw circles of radius R on randomized points in a area100*100? Knowing that its points are located in the following way:

The minimum distance between all points> 6 meters.

do you have an idea ?

Thank you

[Jan]

What is the difference to your question https://www.mathworks.com/matlabcentral/answers/322431-randompoints-condition-distance? Do you have the coordinates of the points already and only want to draw the circles?

Answer 1

[John BG]

Hi Marwen

Here are 2 functions that solve the question, not randomly guessing points from the entire area, but excluding gradually the area sized by already allocated points:

    function [X,Y,Nmax,Dmatrix,delay]=scatter_points7(void)

    % given H and W size of 2D rectangle scatter_points7.m does the following:
    % 1. calculate Nmax, the maximum amount of circles radius R0 that orderly fit in rectangle HxW
    % 2. scatter Ap random points with spatial resolution dx =- 0.1 and dy=0.1
    %           and spaced at least distance R0
    %           if requested amount Ap>Nmax rand2Dpoints breaks because there
    %           is no space to fit in so many random points complying with min
    %           distance R0
    % 3. calculate distance matrix Dmatrix among all points
    % 4. return the coordinates of all generated points in X and Y, along with Dmatrix and Nmax
    % 5. plot points
    % 6. plot safety circles to visually verify 
    % 
    % the amount of generated points is numel(X) and cannot be larger than Nmax 
    % because of the request for the points to be randomly generated.
    % 
    % in this initial version, only manual input through message box.
    % 
    % call examples:
    % 1.
    % [X,Y,Dmatrix,Nmax]=scatter_points7
    % 
    % 2.
    % [X,Y,Nmax]=scatter_points7
    % 2 examples how to verify distance values and check distances meet requirement > R0
    % 1.
    % L=combinator(Ap,2,'c'); 
    % relD2=((X(L(:,2))-X(L(:,1))).^2+(Y(L(:,2))-Y(L(:,1))).^2).^.5
    % find(relD2<R0)
    % 
    % 2.
    % L2=combinator(Ap,2);  
    % Dmatrix=reshape(((X(L2(:,2))-X(L2(:,1))).^2+(Y(L2(:,2))-Y(L2(:,1))).^2).^.5,[Ap Ap]);
    % Dmatrix([1:21:end])=NaN ;
    % Dmatrix(Dmatrix<6)
    % 
    % February 10th 2017
    % version: 1.0
    % author: John Bofarull Guix, any feedback to build next version is welcome at
    % [e-mail] or through the Mathworks website
    % this script was inspired by Mr Marwen Tarhoumi [e-mail]
    % and Matt Fig's [e-mail] mighty function combinator.m available
    % from Mathworks File exchange

    clear all;clc;close all
    format bank
    rng('Shuffle')

    prompt = {'rectangle width W: ','rectange height H: ','amount points to scatter: ','safety distance Radius: '};
    dlg_title = 'Input';
    num_lines = 1;
    defaultans = {'100','100','20','3'};
    input_answer = inputdlg(prompt,dlg_title,num_lines,defaultans);

    Diam0=2*str2double(input_answer{4});
    W=str2double(input_answer{1})+Diam0;
    H=str2double(input_answer{2})+Diam0;
    Ap=str2double(input_answer{3});
    Ap=floor(Ap);
    Nmax=calc_amount_circles(W-Diam0,H-Diam0,Diam0/2);
    fprintf('\nRectangle %f x %f has\n max capacity: %i circles radius %f\n',W-Diam0,H-Diam0,Nmax,Diam0/2);
    if Nmax<Ap 
        fprintf('\nCannot fit in more than %i circles.\n',Nmax);
        X=0;Y=0;Nmax=0;Dmatrix=0;
        return; 
    end                                                         % error message in case requested amount points above Nmax

    As=21;                                                   % amount sides polygon approximating circles
    a=linspace(0,2*pi,As);                       % angle for  circles

    dx=1;dy=1;refine=0;
    X=zeros(1,Ap);Y=zeros(1,Ap);
    if W*H>=900^2
    dx=1;dy=1;
    refine=1;
    else
        dx=.1;dy=.1;
        refine=0;
    end

    x_grid=[(-W+Diam0)/2:dx:(W-Diam0)/2];y_grid=[(-H+Diam0)/2:dy:(H-Diam0)/2];       % avoid circles hitting frame
    [X_grid,Y_grid]=meshgrid(x_grid,y_grid);
    P=[X_grid(:)';Y_grid(:)'];[sz1P,sz2P]=size(P);

    if refine==0
    xc2_base=Diam0/2*cos(a);yc2_base=Diam0/2*sin(a); 
    elseif refine==1
        xc2_base=(Diam0/2+1)*cos(a);yc2_base=(Diam0/2+1)*sin(a); 
    end

    figure;ax=gca;ax.DataAspectRatio=[1 1 1];
    ax.XLim=[-W/2 W/2];ax.YLim=[-H/2 H/2];
    ax.XTick=[(-W+Diam0)/2:10:(W-Diam0)/2];ax.YTick=[(-H+Diam0)/2:10:(H-Diam0)/2];grid on;hold all;
    perimeter=[(-W+Diam0)/2-.25 (W-Diam0)/2+.25 (W-Diam0)/2+.25 (-W+Diam0)/2-.25 (-W+Diam0)/2-.25;
                         (-H+Diam0)/2-.25 (-H+Diam0)/2-.25 (H-Diam0)/2+.25 (H-Diam0)/2+.25 (-H+Diam0)/2-.25];
    plot(perimeter(1,:),perimeter(2,:),'Color',[.3 .3 .3]);
    tic;
    for k=1:1:Ap
       [sz1P,sz2P]=size(P);
       if sz2P>0
           nP = randi(sz2P,1,1);
           else
                break
       end

       if refine==0
           X(k)=P(1,nP);Y(k)=P(2,nP);
           elseif refine==1
                dec_xP=randi([0 499],1,1);dec_xP=dec_xP/1e3;  % worst case added jitter bringing 2 points on crash course
                dec_yP=randi([0 499],1,1);dec_yP=dec_yP/1e3;
                X(k)=P(1,nP)+dec_xP;Y(k)=P(2,nP);Y(k)=P(2,nP)+dec_yP;
       end

       xc2=xc2_base+X(k);yc2=yc2_base+Y(k);
       in=inpolygon(P(1,:),P(2,:),xc2,yc2);
       P(:,find(in>0))=[];                                                                          % exclude area already busy
        figure(1);plot(X(k),Y(k),'r*');                                                         % centre circles
        figure(1);plot(xc2,yc2,'Color',[0.8 0.8 1]);                                % circles radius Diam0/2
    end

    L2=combinator(Ap,2);  % test 2
    Dmatrix=reshape(((X(L2(:,2))-X(L2(:,1))).^2+(Y(L2(:,2))-Y(L2(:,1))).^2).^.5,[Ap Ap]);
    Dmatrix([1:Ap+1:end])=NaN ;
    hold off
    delay=toc;
    end

    function Anxny=calc_amount_circles(H_,L_,D_)
    % calculates 1. amount of circles in hex pattern that fit within 2D rectangle L (columns) x H (tall, lines or rows)
    % used graph from http://www.engineeringtoolbox.com/circles-within-rectangle-d_1905.html to calibrate

    R_=D_/2;

    if L_<(2*R_) 
        Nx=0; 
        else Nx=floor(L_/(2*R_));
    end;

    if H_<(2*R_) 
        Ny=0;
    else
        s=1;
        while H_/(2*R_+s*R_*3^.5)>1;
            s=s+1; 
        end
        Ny=s;
    end;

    if rem(L_,2*R_)>=R_
       min1_evenlines=0;
       else min1_evenlines=1; 
    end;

    Anxny=Nx*Ny-floor(Ny/2)*min1_evenlines;
    end

    % 
    % function handling_input_errors
    % % input checks
    % narginchk(5,5);nargoutchk(4,4);
    % 
    % err_message={'error input type W';'error input type H';'error input type R0';'error input type R0';'error input saturate'};
    % 
    % if(~isreal(H) || ~isscalar(W) || W<=0 )    
    %     error(err_message{1}); 
    % end
    % if(~isreal(H) || ~isscalar(H) || H<=0 )    
    %     error(err_message{2}); 
    % end
    % if(~isreal(R0) || ~isscalar(R0) || R0<=0 )    
    %     error(err_message{3});
    % end
    % 
    % if(~isreal(Ap) || ~isscalar(Ap) || Ap<=0 )    
    %     error(err_message{3}); 
    % end
    % if(~isreal(saturate) || ~isscalar(saturate) || saturate<0 || saturate>1)    
    %     error(err_message{3}); 
    % end
    % end

and the saturation:

function [X,Y,Nmax,Dmatrix]=scatter_points_saturate(W0,H0,R0)
% scatter_points_saturate.m
% given H and W size of 2D rectangle scatter_points7.m does the following:
% 1. calculate Nmax, the maximum amount of circles radius R0 that orderly fit in rectangle HxW
% 2. scatter as many random points as possible with spatial resolution dx =- 0.1 and dy=0.1
%           and spaced at least distance R0
% 3. calculate distance matrix Dmatrix among all points
% 4. return the coordinates of all generated points in X and Y, along with Dmatrix and Nmax
% 5. plot points
% 6. plot safety circles to visually verify 
% 
% the amount of generated points is numel(X) and cannot be larger than Nmax 
% 
% 
% call examples:
% 1.
% [X,Y,Dmatrix,Nmax]=scatter_points_saturate(100,100,3)
% 
% 2.
% [X,Y,Nmax]=scatter_points7
% 2 examples how to verify distance values and check distances meet requirement > R0
% 1.
% L=combinator(Ap,2,'c'); 
% relD2=((X(L(:,2))-X(L(:,1))).^2+(Y(L(:,2))-Y(L(:,1))).^2).^.5
% find(relD2<R0)
% 
% 2.
% L2=combinator(Ap,2);  
% Dmatrix=reshape(((X(L2(:,2))-X(L2(:,1))).^2+(Y(L2(:,2))-Y(L2(:,1))).^2).^.5,[Ap Ap]);
% Dmatrix([1:21:end])=NaN ;
% Dmatrix(Dmatrix<6)
% 
% February 10th 2017
% version: 1.0
% author: John Bofarull Guix, any feedback to build next version is welcome at
% [e-mail] or through the Mathworks website
% this script was inspired by Mr Marwen Tarhoumi [e-mail]
% and Matt Fig's [e-mail] mighty function combinator.m available
% from Mathworks File exchange

% clear all;clc;close all
% format bank
% rng('Shuffle')

% R0=35
Diam0=2*R0;                                                                                      % distance requirement among adjacent points not to be closer than
% W0=50
% H0=70
W=W0+Diam0;                                                                                    % rectangle width, columns
H=H0+Diam0;                                                                                     % rectangle tall, row

   Nmax=calc_amount_circles(W0,H0,Diam0/2)
   Ap=20

As=21; a=linspace(0,2*pi,As);xc=Diam0*cos(a);yc=Diam0*sin(a);

dx=1;dy=1;
X=zeros(1,Ap);Y=zeros(1,Ap);
if W*H>=1e6
dx=1;dy=1;
else dx=.1;dy=.1;
end

X=zeros(1,Ap);Y=zeros(1,Ap);
x_grid=[-W0/2:dx:W0/2];y_grid=[-H0/2:dy:H0/2];                              % avoid circles hitting frame
[X_grid,Y_grid]=meshgrid(x_grid,y_grid);
% X_grid0=X_grid;Y_grid0=Y_grid;

P=[X_grid(:)';Y_grid(:)'];[sz1P sz2P]=size(P);

xc2_base=Diam0/2*cos(a);yc2_base=Diam0/2*sin(a); 
% xc1_base=Diam0*cos(a);yc1_base=Diam0*sin(a);
figure;ax=gca;ax.DataAspectRatio=[1 1 1]
ax.XLim=[-W/2 W/2];ax.YLim=[-H/2 H/2];
ax.XTick=[(-W+Diam0)/2:10:(W-Diam0)/2];ax.YTick=[(-H+Diam0)/2:10:(H-Diam0)/2];grid on;hold all;
perimeter=[(-W+Diam0)/2-.25 (W-Diam0)/2+.25 (W-Diam0)/2+.25 (-W+Diam0)/2-.25 (-W+Diam0)/2-.25;
                     (-H+Diam0)/2-.25 (-H+Diam0)/2-.25 (H-Diam0)/2+.25 (H-Diam0)/2+.25 (-H+Diam0)/2-.25];
plot(perimeter(1,:),perimeter(2,:),'Color',[.3 .3 .3]);

k=1
while sz2P>Diam0*.2
   k=k+1;
   [sz1P sz2P]=size(P);
   if sz2P>0
   nP = randi(sz2P,1,1);
   else
       break
   end
   X(k)=P(1,nP);Y(k)=P(2,nP);
%    xc1=xc1_base+X(k);yc1=yc1_base+Y(k);
   xc2=xc2_base+X(k);yc2=yc2_base+Y(k);
   in=inpolygon(P(1,:),P(2,:),xc2,yc2);
   in=in(:)';P(:,find(in>0))=[];                                                               % exclude area already busy
   figure(1);plot(X(k),Y(k),'g*');                                                             % centre circles
   figure(1);plot(xc2,yc2,'Color',[0.7 .9 0.7]);                                    % circles radius Diam0/2
end

L2=combinator(Ap,2);  % test 2
Dmatrix=reshape(((X(L2(:,2))-X(L2(:,1))).^2+(Y(L2(:,2))-Y(L2(:,1))).^2).^.5,[Ap Ap]);
Dmatrix([1:Ap+1:end])=NaN ;

end

function Anxny=calc_amount_circles(H_,L_,D_)
% calculates 1. amount of circles in hex pattern that fit within 2D rectangle L (columns) x H (tall, lines or rows)
% used graph from http://www.engineeringtoolbox.com/circles-within-rectangle-d_1905.html to calibrate

R_=D_/2;

if L_<(2*R_) 
    Nx=0; 
    else Nx=floor(L_/(2*R_)); 
end;

if H_<(2*R_) Ny=0;
else
    s=1;
    while H_/(2*R_+s*R_*3^.5)>1;
        s=s+1;  
    end
    Ny=s;
end;

if rem(L_,2*R_)>=R_
   min1_evenlines=0;
   else min1_evenlines=1; 
end;

Anxny=Nx*Ny-floor(Ny/2)*min1_evenlines;
end

function handling_input_errors

narginchk(5,5);nargoutchk(4,4);

err_message={'error input type W';'error input type H';'error input type R0';'error input type R0';'error input saturate'}

if(~isreal(H) || ~isscalar(W) || W<=0 )    error(err_message{1}); end
if(~isreal(H) || ~isscalar(H) || H<=0 )    error(err_message{2}); end
if(~isreal(R0) || ~isscalar(R0) || R0<=0 )    error(err_message{3}); end
Ap=floor(Ap); if(~isreal(Ap) || ~isscalar(Ap) || Ap<=0 )    error(err_message{3}); end
if(~isreal(saturate) || ~isscalar(saturate) || saturate<0 || saturate>1)    error(err_message{3}); end
end

the function combinator is available from here

https://uk.mathworks.com/matlabcentral/fileexchange/24325-combinator-combinations-and-permutations?s_tid=srchtitle

attached scatter_points7.m and scatter_points_saturate.m

Appreciating time and attention

John Bofarull Guix

<mailto:[e-mail] [e-mail]>

Answer 2

[Jan]

You can obtain the coordinates from the solution posted here: https://www.mathworks.com/matlabcentral/answers/142322-find-n-random-points-with-a-minimum-distance-r-inside-a-2d-rectangular-box#answer_254413. Then to draw the circles:

nWant  = 20;
Dist   = 6;
[X, Y] = GetPointsRandom(nWant, 100, 100, Dist);

alpha = linspace(0, 2*pi, 32).';
% cx = X.' + 0.5 * Dist * cos(alpha);  % Matlab >= 2016b
% cy = Y.' + 0.5 * Dist * sin(alpha);  % Matlab >= 2016b
cx = bsxfun(@plus, X.', 0.5 * Dist * cos(alpha));
cy = bsxfun(@plus, Y.', 0.5 * Dist * sin(alpha));
figure;
AxesH = axes('NextPlot', 'add');
plot(cx, cy, 'b', 'Parent', AxesH);
axis equal
set(AxesH, 'XLim', [0, 100], 'YLim', [0, 100]);

[Marwen Tarhouni]

@Jan Simon With execution Of this I get nWant = 20;

Dist = 6;

[X, Y] = GetPointsRandom(nWant, 100, 100, Dist);

alpha = linspace(0, 2*pi, 15).';

cx = bsxfun(@plus, X.', 0.5 * Dist * cos(alpha));

cy = bsxfun(@plus, Y.', 0.5 * Dist * sin(alpha)); axes('XLim', [0, 100], 'YLim', [0, 100], 'NextPlot', 'add');

plot(cx, cy, 'b');

axis equal

[Jan]

@Marwen: Are you sure, that you use the code of GetPointsRandom from Answers: 142322? Here the complete code again with the bsxfun version for < R2016b:

function myTest
for nWant = [20, 200]
  Dist   = 6;
  [X, Y] = GetPointsRandom(nWant, 100, 100, Dist);
  alpha = linspace(0, 2*pi, 15).';  % 15 instead of 32, as you like
  cx    = bsxfun(@plus, X.', 0.5 * Dist * cos(alpha));
  cy    = bsxfun(@plus, Y.', 0.5 * Dist * sin(alpha));
  figure;  % Modified for nice XLimits in R2016b:
  AxesH = axes('NextPlot', 'add');
  plot(cx, cy, 'b', 'Parent', AxesH);
  axis equal
  set(AxesH, 'XLim', [0, 100], 'YLim', [0, 100]);
end
end

function [X, Y] = GetPointsRandom(nWant, XWidth, YWidth, MinDist)
X       = zeros(nWant, 1);
Y       = zeros(nWant, 1);
dist_2  = MinDist ^ 2;      % Squared once instead of SQRT each time
iLoop   = 1;                % Security break to avoid infinite loop
nValid  = 0;
while nValid < nWant && iLoop < 1e6
  newX = XWidth * rand;
  newY = YWidth * rand;
  if all(((X(1:nValid) - newX).^2 + (Y(1:nValid) - newY).^2) >= dist_2)
    % Success: The new point does not touch existing points:
    nValid    = nValid + 1;  % Append this point
    X(nValid) = newX;
    Y(nValid) = newY;
  end
  iLoop = iLoop + 1;
end

% Throw an error, if the area is filled too densely:
if nValid < nWant
  error('Cannot find wanted number of points in %d iterations.', iLoop)
end  
end

I added a figure, because in your screenshot two axes are drawn on eachother. But this could not lead to the strange overlapping circles. The result on R2016b:

This differs substantially from you screenshot and I cannot explain why. Does your code differ from this one?

You see that 200 is near to the possible limit and then all methods to find such a dense random packing can fail. In the worst case an algorithm sets the points to a regular triangular grid with the double distance minus eps, such that the area is almost empty, but this very unlikely.

Please note, that this code allows a distance greater equal 6. If only greater 6 is wanted replace ">= dist_2" with "> dist_2" .

[Marwen Tarhouni]

@Simon The first condition that the distance between center points must be> 6m

The second I want to draw circles of radii = 15. I want than the circles to overlap since the first condition is realized (distance between center of the ponits> 6meters) There is no problem when the circles overlap

The center points do not show up on your figures After this two conditions each center point must know its distance from the others (this function already exists in GetPointsRandom) The area I want to work with is the area (100 * 100) The number of randomized points = the maximum points that respect the 6 meter condition for this zone

set(AxesH, 'XLim', [-100, 100], 'YLim', [-100, 100]);

[John BG]

I completely agree with Marwen: the centre points have to be clearly indicated because each point is a mobile phone, the R circles tell each range (a fair approximation of how far each mobile phone may reach individually, am I right?), and the overlapping circles are visually helpful to verify no two points are closer than the required minimal distance.

Marwen merçi beacou und Viel Spaß

John BG

<mailto:[e-mail] [e-mail]>

[Jan]

@Marwen: The center point of the circles are not displayed in my figures, because this detail has not been mentioned in this question. If you want to draw the circles with a radius of 15, use:

alpha = linspace(0, 2*pi, 32).';
R     = 15;
cx    = bsxfun(@plus, X.', R * cos(alpha));
cy    = bsxfun(@plus, Y.', R * sin(alpha));

Add this to get the center points also:

plot(X, Y, '*');

You mention "'XLim', [-100, 100]" and "the area 100 * 100". The first is an area of 200*200.

"The maximum points that respect the 6 meter condition" is not uniquely defined, because this depends on the random locations. If you want the maximum number, and not a specified number, this might be a different problem. So I'm starting to get confused about what you need.

"Cell phones"? This detail has not been mentioned before and for the solution it does not matter.

[John BG]

the centre points are not 'a detail' but the utmost important point, it's the first thing you should plot when considering visual presentation of anything attempting to answer this and similar questions.

John BG

<mailto:[e-mail] [e-mail]>

[Jan]

@John: The question was:

 How we can draw circles of radius R on randomized points in
 a area 100*100?

Then the center is "a detail [which] has not been mentioned".

Answer 3

[KSSV]

Check this link: https://in.mathworks.com/matlabcentral/answers/158357-create-random-points-in-a-rectangular-domain-but-with-minimum-separation-distance

[Marwen Tarhouni]

@KSSV, I want to create circles surrounding random points and The minimum distance between all points> 6 meters.

[KSSV]

That link....creates pints.....around those points circles can be drawn. Isn't it?