# Find closest sphere in a cluster with random spheres having different diameters-2

1 view (last 30 days)

Show older comments

Hi,

I posted this question before and i got good help in how to find closest sphere in a cluster with random spheres having different diameters. I am making clusters of spheres using random number generator. Some of the spheres are connected and some of them are not connected. I want to exclude the not connected ones, so I only get the clusters/connected spheres .I have random spheres having random diameters between 0.1 and 0.5.

here is my code till now:

%% Data for Agglomerates

N = 15; % number of spheres

a = 0.1 ; % lowest diameter of sphere

b = 0.5 ; % highest diameter of sphere

Diam = a + (b-a).*rand(N,1);

aaa= 0; % minimum center x and y coordinate limit for spheres

bbb= 1.5 ; % maximum center x and y coordinat limit for sphere

M=3 ;

Axes= zeros(N,M);

Axes(:,1)=aaa+(bbb-aaa)*rand(N,1);

for k=2:M

aaa=randperm(N);

Axes(:,k)=Axes(aaa,1);

end

Axes_Label ={'X axis','Y axis','Z axis'};

Data_agglo = [Diam Axes];

Data_Label ={'Diameter','X axis','Y axis','Z axis'};

R = Diam ./2;

%% Algorithm to find connected spheres and eliminate outliers

distances = zeros(size(Data_agglo,1));

active_spheres = zeros(size(Data_agglo,1));

delta = 1e-4;

for p =1:1:size(Data_agglo,1)

for q = 1:1:size(Data_agglo,1)

if size(Data_agglo,1)==1

distances(p,q) = sqrt((Data_agglo(q,2)-Data_agglo(p,2)).^2 + (Data_agglo(q,3)-Data_agglo(p,3)).^2 + (Data_agglo(q,4)-Data_agglo(p,4)).^2);

active_spheres(p,q) = (distances(p,q) <= (Data_agglo(p,1)/2)+(Data_agglo(q,1)/2) - delta );

else

if p ~= q

distances(p,q) = sqrt((Data_agglo(q,2)-Data_agglo(p,2)).^2 + (Data_agglo(q,3)-Data_agglo(p,3)).^2 + (Data_agglo(q,4)-Data_agglo(p,4)).^2);

active_spheres(p,q) = (distances(p,q) <= (Data_agglo(p,1)/2)+(Data_agglo(q,1)/2) - delta );

end

end

end

end

rem_particle = sum(active_spheres);

rem_particle = find(rem_particle == 0);

Data_agglo2 = Data_agglo;

Data_agglo(rem_particle,:) = [];

R(rem_particle) = [];

if isempty(R)

error('No spheres are connected!\n');

end

%% generate mesh

dipS = 0.01;

xmin = min(Data_agglo(:,2)-R);

xmax = max(Data_agglo(:,2)+R);

ymin = min(Data_agglo(:,3)-R);

ymax = max(Data_agglo(:,3)+R);

zmin = min(Data_agglo(:,4)-R);

zmax = max(Data_agglo(:,4)+R);

[Xgrid,Ygrid,Zgrid]= ndgrid((xmin:dipS:xmax)-(xmin+xmax)/2,(ymin:dipS:ymax)-(ymin+ymax)/2,(zmin:dipS:zmax)-(zmin+zmax)/2);

Data_agglo(:,2:4) = Data_agglo(:,2:4) - [(xmin+xmax)/2,(ymin+ymax)/2,(zmin+zmax)/2];

% plot3(Xgrid(:),Ygrid(:), Zgrid(:),'.','MarkerFaceColor','blue');

% hold all

%% get active dipoles

tic

active = false(size(Xgrid));

for i =1:1:size(Data_agglo,1)

active = active | (((Xgrid - Data_agglo(i,2)).^2 + (Ygrid - Data_agglo(i,3)).^2 + (Zgrid - Data_agglo(i,4)).^2) <= R(i).^2);

end

plot3(Xgrid(active),Ygrid(active),Zgrid(active),'o','MarkerFaceColor','red');

axis equal

grid on

toc

The problem I am getting is :

I eliminate the spheres, but I am left with two sets of connected spheres. I want to choose one of them and delete the other one. there is no pirority which one to use.

I am attaching the figure file as well. Other problem I get is that when i give the number of spheres =1, my program doesnt work. Also if it is possible to find the largest set of connected spheres as shown in the picture and keep them?

Does anyone knows how it can be done?

##### 4 Comments

Walter Roberson
on 4 May 2022

Untested.

equivalence_class_numbers = (1:number_of_spheres).';

equivalence_class_members = num2cell(equilvance_class_numbers);

for idx = 2 : number_of_spheres

sp = Index_of_conflicting_spheres(idx, sphere_coordinates, sphere_radii);

if isempty(sp); continue; end

sp = [sp; idx];

new_class = min(sp);

new_members = unique(cell2mat(equivalence_class_members(sp)));

equivalance_class_members(new_class) = new_members;

equivalance_class_numbers(sp) = new_class;

equivalance_class_members(setdiff(sp, new_class)) = {};

end

At the end of this, equivalance_class_numbers will be a vector of length number_of_spheres that gives the cluster number each sphere belongs to. The numbers will likely have gaps -- for example if sphere #3 is connected to sphere #1 then the class number (1) and at (3) would both be 1, and 3 would not appear in any class number. For some purposes you would want to renumber without gaps (findgroups() would help for that), but for the purpose of your particular need, you do not care about the cluster number in absolute terms, you just want to know what is joined to what.

At the end of the code, equivalance_class_members is a cell, and for any cluster number that does not exist (like the 3 in the above example, because sphere 3 is part of cluster 1), the entry in equivalance_class_members will be empty. For each cluster number that is populated, equivalance_class_members will be a list of all the sphere numbers that belong to the same cluster.

If you see a equivalance_class_members entry that is a scalar (but not empty), that means that the sphere is isolated, and in such a case equivalance_class_numbers(K) == K .

You can cellfun(@length, equivalance_class_members) and use max() to find the largest cluster (at least in terms of numbers of members -- not in terms of volume !)

### Accepted Answer

Jon
on 4 May 2022

Edited: Jon
on 4 May 2022

Here is another approach, in which we consider the spheres to be nodes in an undirected graph with edges connecting the nodes whose corresponding spheres touch. We can then use some graph theory algorithms to find the connected components, which I think are just what you are looking for. Here is some example code, untested except to see that it runs (I didn't follow your notation, or scaling of the locations etc, just wanted to get the main idea across. I think you can adapt what I did here to the specifics of your problem)

% parameters

N = 15; % number of spheres

% box dimensions that centers will fall in

xMax = 1

yMax = 2

zMax = 1

% minimum and maximum diameters

a = 0.1;

b = 0.5;

% Generate random center points

x = rand(N,1)*xMax;

y = rand(N,1)*yMax;

z = rand(N,1)*zMax;

% Generate random diameters,

d = a + (b-a)*rand(N,1);

% Find the minimum separation between sphere centers at which they will

% touch

r = d/2; % radius

Stouch = r + r'; % use scalar expansion between column and row to get every pair

% Find the distance between every pair of points, making use of scalar

% expansion of difference between a column and row to get every pair

S = sqrt((x-x').^2 + (y-y').^2 + (z -z').^2)

% Generate logical matrix whose elements are set to 1 (true) if the ith and

% jth sphere touch

Touch = S <= Stouch;

% Consider as a graph problem with a node for each sphere and edges

% connecting the vertices if the corresponding spheres touch. So our Touch

% matrix is exactly the adjacency matrix of this graph

% Obviously each sphere touches itself

% so don't include self loops (edges that connect the node to itself)

G = graph(Touch,'omitselfloops');

% Plot the graph just to visualize it,

plot(G)

% Find the connected components, there is a path between any two nodes in a

% connected component, this means the connected components are exactly the

% clusters of connected spheres you are looking for

components = conncomp(G); % components(i) gives the cluster number the ith sphere belongs to

% Find the connected component with the most spheres in it using groupcounts

% gc is the number of spheres in each group

[gcnt,grp] = groupcounts(components(:)); % use (:) to make sure it is a column

[~,idx] = max(gcnt); % index of group with the most spheres

grpMax = grp(idx); % group number of the group with the most spheres

% Get a indices of the spheres that are in the largest cluster that we wish

% to keep

idxKeep = find(components == grpMax);

>> idxKeep

idxKeep =

2 3 5 7 11 12

##### 11 Comments

Walter Roberson
on 7 Jun 2022

Jon
on 7 Jun 2022

@Walter Roberson, thanks for catching that. It has been long enough since looking at this, that I forgot that the spheres could overlap. Also I was unaware that a new question had been opened on this topic.

@hamzah khan looks like you have a lot of suggestions to work with already in https://www.mathworks.com/matlabcentral/answers/1733060-find-the-center-of-the-cluster-made-from-random-spheres-having-different-diameters I would suggest not continuing here and just moving to your new question

### More Answers (0)

### See Also

### Categories

### Community Treasure Hunt

Find the treasures in MATLAB Central and discover how the community can help you!

Start Hunting!