MATLAB Answers

Pseudo-random distribution of points with minimum distance

26 views (last 30 days)
Florian
Florian on 28 Mar 2019
Edited: Bruno Luong on 26 Sep 2019
Hi all,
I use the following code to generate an exact number of points (650 in this case) within a one hectar area:
numPoints = 650; % number of points to be generated
width = 100; % length in m of one square side
x = 0;
y = 0;
figure('Position', [300 300 900 900])
rectangle('Position', [x, y, width, width],'LineWidth',2,'LineStyle','--');
grid on;
hold on;
xRandom = 50 + (width * rand(1, numPoints) - width / 2);
yRandom = 50 + (width * rand(1, numPoints) - width / 2);
plot(xRandom, yRandom, 'b.', 'MarkerSize', 8);
hold on;
plot(xRandom, yRandom, 'ro', 'MarkerSize', 20);
title(['',num2str(numPoints),' points inside one hectar'], 'Interpreter', 'None');
xlabel('length in meters');
ylabel('length in meters');
axis equal tight;
What I would like to achieve is that the minimum distance between these points is 3 meters. This condition is fulfilled if the red circles (with 3 m diameter - only approximated here using MarkerSize) only touch each other but don't overlap as they currently do (see image below).
Does anybody know how to accomplish this?
PS: at a completely regular spacing of 650 points in one hectar there would be 4 m distance between each point. Maybe a possible solution would be to create a regular spacing first and then add or substract a smaller random increment?

  3 Comments

James Tursa
James Tursa on 28 Mar 2019
If there are only a few points, you could just push the offenders away from each other. But if you have a large number of points, this will turn into more of a packing problem than a random distribution problem. Is that acceptable?
Florian
Florian on 28 Mar 2019
I need to end up with the exact number, so eliminating offenders (as proposed in previous solutions) is not an option. Maybe 'rand' is not the best function to start with...
James Tursa
James Tursa on 28 Mar 2019
I wasn't suggesting eliminating the offenders, but pushing them away from each other until they are not offending anymore.

Sign in to comment.

Accepted Answer

Torsten
Torsten on 29 Mar 2019
See Bruno Luong's code under
https://de.mathworks.com/matlabcentral/answers/432516-model-of-a-crowd-on-concert-venue-or-how-to-distribute-random-points-according-to-the-2d-window-dist

  6 Comments

Show 3 older comments
Bruno Luong
Bruno Luong on 29 Mar 2019
For reference I put here the code with uniform distribution
L = 100; % <-- Choose length of square sides
x0 = 0; y0 = 0; % <-- Choose center of square
n = 300; % <-- Choose number of points
% Generate uniform-distribution on the square
X = rand(n,2)*L - (L/2*[1,1]+[x0,y0]);
XYR = [x0,y0]+[[-1;1;1;-1;-1],[-1;-1;1;1;-1]]*L/2;
XB = interp1((0:4)'*L,XYR,linspace(0,4*L,200));
XB(end,:) = [];
nrepulsion = 500;
% Repulsion of seeds to avoid them to be too close to each other
n = size(X,1);
Xmin = [x0-L/2,y0-L/2];
Xmax = [x0+L/2,y0+L/2];
% Point on boundary
XR = x0+[-1,1,1,-1,-1]*L/2;
YR = y0+[-1,-1,1,1,-1]*L/2;
cla;
hold on
plot(XR,YR,'r-');
h = plot(X(:,1),X(:,2),'b.');
axis equal
dmin = 3; % minima distance between 2 objects
d2min = dmin*dmin;
beta = 0.5;
for k = 1:nrepulsion
XALL = [X; XB];
DT = delaunayTriangulation(XALL);
T = DT.ConnectivityList;
containX = ismember(T,1:n);
b = any(containX,2);
TX = T(b,:);
[r,i0] = find(containX(b,:));
i = mod(i0+(-1:1),3)+1;
i = TX(r + (i-1)*size(TX,1));
T = accumarray([i(:,1);i(:,1)],[i(:,2);i(:,3)],[n 1],@(x) {x});
maxd2 = 0;
R = zeros(n,2);
move = false(n,1);
for i=1:n
Ti = T{i};
P = X(i,:) - XALL(Ti,:);
nP2 = sum(P.^2,2);
if any(nP2<4*d2min)
move(i) = true;
move(Ti(Ti<=n)) = true;
end
maxd2 = maxd2 + mean(nP2);
b = Ti > n;
nP2(b) = nP2(b)*5; % reduce repulsion from each point of the border
R(i,:) = sum(P./max((nP2-d2min),1e-3),1);
end
if ~any(move)
break
end
if k==1
v0 = (L*5e-3)/sqrt(maxd2/n);
end
R = R(move,:);
v = v0/sqrt(max(sum(R.^2,2)));
X(move,:) = X(move,:) + v*R;
% Project back if points falling outside the rectangle
X = min(max(X,Xmin),Xmax);
set(h,'XData',X(:,1),'YData',X(:,2));
pause(0.01);
end
theta = linspace(0,2*pi,65);
xc = dmin/2*sin(theta);
yc = dmin/2*cos(theta);
% plot circles f diameter dmin around random points
for i=1:n
plot(X(i,1)+xc,X(i,2)+yc,'k');
end

Sign in to comment.

More Answers (1)

Image Analyst
Image Analyst on 29 Mar 2019
See my attached demo that I've posted before.
Capture.PNG

  2 Comments

Florian
Florian on 29 Mar 2019
This script doesn't seem to work for me using MATLAB Version: 9.2.0.556344 - R2017a).
I use maxPoints = 1000 and numPoints = 650 with minClosestDistance = 0.03;
The script gives me a message that 650 points have been computed, but when I open x and y there are only NaNs in there.
Image Analyst
Image Analyst on 29 Mar 2019
You should not have set numPoints = 650. Then it will only place an additional 350 points. You should have left it at zero and let it do the count by itself.
When I ran it, I tried a million times to place 1000 points and I was only able to place 789. It looks like Bruno's algorithm works though, at least for the one time I tried it.

Sign in to comment.

Sign in to answer this question.