1 view (last 30 days)

Show older comments

Hello! I'm trying to generate a step function that will check to see if a point is inside the spherocylinder cell. After that is generated then I need to use normrnd to pick a random point within the cell. However, I'm stuck on this, therefore any help would be appreciated.

This is what I have so far:

function [X,Y,Z]=spherocylinder(pos, vec, totalLength, AR);

%% preliminary

points = 20;

vec = vec/norm(vec);

D = (0.65);

radius = D/2;

bottom = -0.675; % where bottom of cylinder is

top = 0.675; % where top of cylinder is

%% Generate cylinder with center at origin

[xcyl,ycyl,zcyl]=cylinder(radius,points);

zcyl(1,:) = bottom;

zcyl(2,:) = top;

%% Generate points for sphere

[x,y,z] = sphere(points);

pdist = points/2+1;

% split sphere into top and bottom hemispheres

xtop = x(pdist:end,:)*radius;

ytop = y(pdist:end,:)*radius;

ztop = z(pdist:end,:)*radius+top;

xbot = x(1:pdist,:)*radius;

ybot = y(1:pdist,:)*radius;

zbot = z(1:pdist,:)*radius+bottom;

%% Combine data for spherocylinder oriented along z-direction

X = [xbot; xcyl; xtop];

Y = [ybot; ycyl; ytop];

Z = [zbot; zcyl; ztop];

%% Plot spherocylinder

surf(X, Y, Z, 'EdgeColor','none',...

'LineStyle','none','FaceColor',[0.85 0.85 0.85]);

axis equal;

% length, width and height of box

L = 0.65; W = 0.65; H = 2;

% pick plenty of random points

x = L*rand(1,1); %x-coordinate of a point

y = W*rand(1,1); %y-coordinate of a point

z = H*rand(1,1); %z-coordinate of a point

scatter3(x, y, z, 'MarkerFaceColor','b','MarkerEdgeColor','b');

box

end

Thank you for any help!

DGM
on 15 Apr 2021

I would consider the line segment connecting the centers of the two hemispheres. All points within the object should be within distance=radius from that line segment.

As an aside, it's good to include some form of synopsis in your function to describe what the input (and output) arguments are.

William Rose
on 15 Apr 2021

@Mari Teli, I do not understand what spherocylinder() is supposed to do, because your question does not seem to match the output from the function. I ran the function. It produced a 3D plot of what appears to be a single point.

That plot was not enlightening, but the function also returned values of arrays x,y,z. I plotted those arrays with "surf(x,y,z)" and got a 3D spherocylinder, as shown on left below. To get the correctly scaled object, as shown on the right, in which the ends look like hemispheres, add the command "axis equal". You could add these to your function, if you want such a plot.

>> clear;

>> [x,y,z]=spherocylinder(1,1,1,1);

>> figure;

>> surf(x,y,z)

>> axis equal

becomes

Three of the arguments to the function are not used at all, an the other argument, vec, is scaled, but it is otherwise not used, so it does not affect the plot made or the output variables. What did you want those four function arguments to do?

You said in your quesiton "I'm trying to generate a step function that will check to see if a point is inside the spherocylinder cell. After that is generated then I need to use normrnd to pick a random point within the cell." That seems backwards. It would make more sense to pick a random point first, then test it to see it it is inside the cell.

The parameters of the cyliner are its radius, top height, and bottom height. Here is a function that retrns 1 if point is inside and 0 if not. I tested it for many combinations and it worked, but you should do more testing, to be sure.

function b = inSpherocylinder(p,radius,th,bh)

% WCRose 20210415

% Return 1 if p is in spherocylinder, else return 0

% Cylinder's axis is assumed to be along z and centered on x,y=0,0.

% th=level where cylinder meets top hemisphere

% bh=level where cylinder meets bottom hemisphere

% p should be a 3-vector: p=[x,y,z].

% th should be greater than or equal to bh.

if length(p)~=3

b=-99;

disp('Error in inCylinder: p should be a 3-element vector.');

return;

elseif bh>th

b=-99;

disp('Error in inCylinder: th should be >= bh.');

end

x=p(1); y=p(2); z=p(3);

if x^2+y^2>radius^2

b=0; %point is outside cylinder or its extension

return;

elseif (z>th) & (x^2+y^2+(z-th)^2>radius^2)

b=0; %point is above the top hemi

return;

elseif (z<bh) & (x^2+y^2+(z-bh)^2>radius^2)

b=0; %point is below the bottom hemi

return;

else

b=1; %point is inside

end

end

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

Start Hunting!