MATLAB: Intersection volume of two 3d alphashapes

alphashapeintersectionMATLABvolume

I have two alphaShapes in 3d that overlap each other and I would like to know the volume of the overlapsing. I already did it in 2d with polygons (but it was the overlaping area, not volume of course) using the functions intersect and polyarea but I can't find something similar for 3d shapes. I use matlab r2018a.

Best Answer

Let's say you have two shapes that were developed as follows using column vectors as inputs:

shp1=alphaShape(x1,y1,z1);
shp2=alphaShape(x2,y2,z2);

You could determine the points of shp1 that are in shp2 and vica versa the points of shp2 that are in shp1 via the inShape function, then use those to get a new shape, shp3, that represents the intersection:

id1=inShape(shp2,x1,y1,z1);
id2=inShape(shp1,x2,y2,z2);
shp3=alphaShape([x1(id1); x2(id2)], [y1(id1); y2(id2)], [z1(id1); z2(id2)]);

You can then get the volume by:

V=volume(shp3);

Related Solutions

MATLAB: Point or multiple points is/are in a triangle?

tf = inpolygon(Points(:,1),Points(:,2),Triangle(:,1),Triangle(:,2))

MATLAB: How to Generate 100 Random points inside the figure

Create the 3 shapes at first:

warning('off', 'MATLAB:alphaShape:DupPointsBasicWarnId');
Shape = cell(1, 3);
% Top sphere:
[x, y, z] = sphere;
Alpha     = 2;
Shape{1}  = alphaShape(x(:), y(:), z(:) + 9, Alpha);
% Cone:
r = 1.0;
h = 9.0;
m = h/r;
[R,A] = meshgrid(linspace(0,r,11), linspace(0,2*pi,41));
X2 = R .* cos(A);
Y2 = R .* sin(A);
Z2 = m * R;
Shape{2} = alphaShape(X2(:), Y2(:), Z2(:), Alpha);
% Bottom sphere:
Shape{3} = alphaShape(x(:) * 0.6, y(:) * 0.6, z(:) * 0.6, Alpha);
warning('on', 'MATLAB:alphaShape:DupPointsBasicWarnId');
% Surrounding cuboid:
coorMin   = [-1, -1, -0.6];  % Minimal x,y,z value
coorRange = [2, 2, 10.4];    % Max - Min of x,y,z value

Now choose random points inside the surrounding cuboid and collect the points, which are in any of the shapes:

% Get points inside the shapes:
maxIter = 1e6;  % Security limit
iter    = 0;
n       = 100;  % Number of output points
P       = zeros(n, 3);
iP      = 0;
while iP < n && iter <= maxIter
   Q = rand(1, 3) .* coorRange + coorMin;
   if inShape(Shape{1}, Q) || inShape(Shape{2}, Q) || inShape(Shape{3}, Q)
      iP       = iP + 1;
      P(iP, :) = Q;
   end
   iter = iter + 1;
end
if iter == maxIter  % Increase maxIter on demand:
   error('Cannot find enough points.');
end

This approach is simpler than creating one shape only, because you do not have to consider the overlap or know anything about the geometry or relative position of the objects. Even the sharp kink between the lower sphere and the cone is not problem here.

The smaller the volume of the shape containing the accepted points compared to the surrounding cuboid, the higher is the number of rejected points. Set maxIter accordingly.

Display the shape and the points:

% Show the shape and the points:
figure;
subplot(1,2,1, 'NextPlot', 'add');
for k = 1:3
   plot(Shape{k});
end
axis([-1,1, -1,1, -1,10], 'equal')
view(3)
subplot(1,2,2, 'NextPlot', 'add', 'XLim', [-1,1], 'YLim', [-1,1]);
plot3(P(:,1), P(:,2), P(:,3), '.');
axis([-1,1, -1,1, -1,10], 'equal')
view(3)

Here with n=1000 points:

Best Answer

Related Solutions

MATLAB: Point or multiple points is/are in a triangle?

MATLAB: How to Generate 100 Random points inside the figure

Related Question