How to get center coordinates of multiple cylinders

Here is the code that plots multiple cylinders.

% Inputs
Diameter = 18;
Height = 100;
rows = 685;
columns = 300;
GapWidth = 2;
Distance = 10;
r = Diameter/2;
scalingFactor = sqrt(5) / 2;
gridSpacing = r/scalingFactor+GapWidth;

% Cuboid
P = [columns/2, rows/2, Height/2];
L = [columns, rows, Height];
O = P-L/2;
plotcube(L,O,.8,[1 1 1]);   %plotcube(L,O,.8,[1 1 1]);   %if present  
hold on
plot3(P(1),P(2),P(3),'*k')
alpha(.01)

%Cylinder
[X,Y,Z] = cylinder(r);
Z = Z*Height;
   
for i = r : scalingFactor * gridSpacing*2 : columns
    for ii = r : gridSpacing*8+Distance*2 : rows
        surf(X+i,Y+r,Z);               
    end
end
for k = r+gridSpacing: scalingFactor * gridSpacing*2 : columns
    for kk= r+gridSpacing*2+Distance : gridSpacing*8+Distance*2 : rows
        surf(X+k,Y+kk,Z); 
    end
end
if Distance>0
    b = 1;
else
    b = 0;
end
for l = r+(scalingFactor * gridSpacing*2)*b:scalingFactor * gridSpacing*2:columns
    for ll = r+gridSpacing*4+Distance : gridSpacing*8+Distance*2 : rows
        surf(X+l,Y+ll,Z);
    end
end
for j = r+gridSpacing:scalingFactor * gridSpacing*2:columns-Distance
    for jj = r+gridSpacing*6+Distance*2 : gridSpacing*8+Distance*2 : rows
        surf(X+j,Y+jj,Z);
    end
end
for n = r : scalingFactor * gridSpacing*2 : columns-Distance
    for nn = r+gridSpacing*8+Distance*2: gridSpacing*8+Distance*2 : rows
        surf(X+n,Y+nn,Z);               
    end
end

axis image
alpha(.5)
view(3)

Function to plot the cuboid:

function plotcube(varargin)

% Default input arguments
inArgs = { ...
  [10 56 100] , ... % Default edge sizes (x,y and z)
  [10 10  10] , ... % Default coordinates of the origin point of the cube
  .7          , ... % Default alpha value for the cube's faces
  [1 0 0]       ... % Default Color for the cube
  };

% Replace default input arguments by input values
inArgs(1:nargin) = varargin;

% Create all variables
[edges,origin,alpha,clr] = deal(inArgs{:});

XYZ = { ...
  [0 0 0 0]  [0 0 1 1]  [0 1 1 0] ; ...
  [1 1 1 1]  [0 0 1 1]  [0 1 1 0] ; ...
  [0 1 1 0]  [0 0 0 0]  [0 0 1 1] ; ...
  [0 1 1 0]  [1 1 1 1]  [0 0 1 1] ; ...
  [0 1 1 0]  [0 0 1 1]  [0 0 0 0] ; ...
  [0 1 1 0]  [0 0 1 1]  [1 1 1 1]   ...
  };

XYZ = mat2cell(...
  cellfun( @(x,y,z) x*y+z , ...
    XYZ , ...
    repmat(mat2cell(edges,1,[1 1 1]),6,1) , ...
    repmat(mat2cell(origin,1,[1 1 1]),6,1) , ...
    'UniformOutput',false), ...
  6,[1 1 1]);


cellfun(@patch,XYZ{1},XYZ{2},XYZ{3},...
  repmat({clr},6,1),...
  repmat({'FaceAlpha'},6,1),...
  repmat({alpha},6,1)...
  );

view(3);
end

I need to find the x and y coordinates of the centers and also the total number of the cylinders.

Answer 1

[William Rose]

@Ege Arsan,

Counting cylinders: I added counter numCylinder to count the number of cylinders. Search the code for numCylinder. The counter is incremented by one each time a cylinder is plotted. The counter reports there are 422 cylinders. The plot shows 345-350 obvious cylinders. Where are the 70 or so cylinders that are not obvious? I suspect that some cylinders overlap or superimpose.

Finding the centers: I added array CylCtr (1000 by 3) to save the cylinder centers. Search the code for CylCtr. One row is filled in each time a cylinder is plotted. The basic cylinder which is replicated many times has its center at [0,0,Height/2]. Therefore we determine each cylider's center by the offsets that are added to (X,Y,Z) when each cylinder is plotted with surf(). I chose 1000 rows for the initial size of CylCtr, in order to be more than big enough.

Notice the console output, which was not present before:

Number of cylinders=422. Center of 1=(9.00,9.00,50.00). Center of 422=(278.67,611.39,50.00).

numCylinder=0;          %cylinder counter
CylCtr=zeros(1000,3);   %allocate array for center locations
% Inputs
Diameter = 18;
Height = 100;
rows = 685;
columns = 300;
GapWidth = 2;
Distance = 10;
r = Diameter/2;
scalingFactor = sqrt(5) / 2;
gridSpacing = r/scalingFactor+GapWidth;
% Cuboid
P = [columns/2, rows/2, Height/2];
L = [columns, rows, Height];
O = P-L/2;
plotcube(L,O,.8,[1 1 1]);   %plotcube(L,O,.8,[1 1 1]);   %if present  
hold on
plot3(P(1),P(2),P(3),'*k')
alpha(.01)
%Cylinder
[X,Y,Z] = cylinder(r);
Z = Z*Height;

for i = r : scalingFactor * gridSpacing*2 : columns
    for ii = r : gridSpacing*8+Distance*2 : rows
        surf(X+i,Y+r,Z);
        numCylinder=numCylinder+1;
        CylCtr(numCylinder,:)=[i,r,Height/2];
    end
end
for k = r+gridSpacing: scalingFactor * gridSpacing*2 : columns
    for kk= r+gridSpacing*2+Distance : gridSpacing*8+Distance*2 : rows
        surf(X+k,Y+kk,Z); 
        numCylinder=numCylinder+1;
        CylCtr(numCylinder,:)=[k,kk,Height/2];
    end
end
if Distance>0
    b = 1;
else
    b = 0;
end
for l = r+(scalingFactor * gridSpacing*2)*b:scalingFactor * gridSpacing*2:columns
    for ll = r+gridSpacing*4+Distance : gridSpacing*8+Distance*2 : rows
        surf(X+l,Y+ll,Z);
        numCylinder=numCylinder+1;
        CylCtr(numCylinder,:)=[l,ll,Height/2];
    end
end
for j = r+gridSpacing:scalingFactor * gridSpacing*2:columns-Distance
    for jj = r+gridSpacing*6+Distance*2 : gridSpacing*8+Distance*2 : rows
        surf(X+j,Y+jj,Z);
        numCylinder=numCylinder+1;
        CylCtr(numCylinder,:)=[j,jj,Height/2];
    end
end
for n = r : scalingFactor * gridSpacing*2 : columns-Distance
    for nn = r+gridSpacing*8+Distance*2: gridSpacing*8+Distance*2 : rows
        surf(X+n,Y+nn,Z);               
        numCylinder=numCylinder+1;
        CylCtr(numCylinder,:)=[n,nn,Height/2];
    end
end
fprintf('Number of cylinders=%d. Center of 1=(%.2f,%.2f,%.2f). Center of %d=(%.2f,%.2f,%.2f).\n',...
    numCylinder,CylCtr(1,:),numCylinder,CylCtr(numCylinder,:));

Number of cylinders=422. Center of 1=(9.00,9.00,50.00). Center of 422=(278.67,611.39,50.00).

axis image
alpha(.5)
view(3)

Cuboid-plotting function:

function plotcube(varargin)
% Default input arguments
inArgs = { ...
    [10 56 100] , ... % Default edge sizes (x,y and z)
    [10 10  10] , ... % Default coordinates of the origin point of the cube
    .7          , ... % Default alpha value for the cube's faces
    [1 0 0]       ... % Default Color for the cube
    };
% Replace default input arguments by input values
inArgs(1:nargin) = varargin;
% Create all variables
[edges,origin,alpha,clr] = deal(inArgs{:});
XYZ = { ...
    [0 0 0 0]  [0 0 1 1]  [0 1 1 0] ; ...
    [1 1 1 1]  [0 0 1 1]  [0 1 1 0] ; ...
    [0 1 1 0]  [0 0 0 0]  [0 0 1 1] ; ...
    [0 1 1 0]  [1 1 1 1]  [0 0 1 1] ; ...
    [0 1 1 0]  [0 0 1 1]  [0 0 0 0] ; ...
    [0 1 1 0]  [0 0 1 1]  [1 1 1 1]   ...
    };
XYZ = mat2cell(...
    cellfun( @(x,y,z) x*y+z , ...
    XYZ , ...
    repmat(mat2cell(edges,1,[1 1 1]),6,1) , ...
    repmat(mat2cell(origin,1,[1 1 1]),6,1) , ...
    'UniformOutput',false), ...
    6,[1 1 1]);
cellfun(@patch,XYZ{1},XYZ{2},XYZ{3},...
    repmat({clr},6,1),...
    repmat({'FaceAlpha'},6,1),...
    repmat({alpha},6,1)...
    );
view(3);
end

[Ege Arsan]

I have a follow up question. I want to iterate the plot by checking wheter the number of elements in y direction is odd and in that case i want to exclude the last created row of cylinders and set the new row legth according to the row created one before the last.

I have changed the CylCtr to 1000x2 since the z coordinates don't change and have sorted the array as the following:

    [~,idx] = sort(CylCtr(:,2));    % sort by y values
    sortCylCtr = CylCtr(idx,:);    % sort the whole matrix using the sort indices
    sortCylCtr( ~any(sortCylCtr,2), : ) = [];  %rows
    sortCylCtr( :, ~any(sortCylCtr,1) ) = [];  %columns
    
    [b,m1,n1] = unique(sortCylCtr,'first');
    [c1,d1] =sort(m1);
    b = b(d1);

I think the code i desire should look similar to something like this:

if rem(numel(b), 2) ~= 0 
     b(end,1) = 0;
     newRows = b(end,1)-r-GapWidth;
     rows = newRows;
end

But i don't know how to integrate this as an iteration to the plot.

[William Rose]

@Ege Arsan, I'm not sure. I do not follow the logic of the 5 sets of for loops inside of which the plotting happens. Those loops are over i, k, l, j, n. Once an object has been added to the plot, I don;t know how to remove it. So you might have to not plot it in the first place. Which means you have to understand the logic of those for loops.

Answer 2

[Ege Arsan]

Thank you! Works perfectly.

[William Rose]

@Ege Arsan, You're welcome. Good luck with your work.