Code covered by the BSD License  

Highlights from
Fast Alpha Hulls (alpha shapes in 3d; parfor enabled)

5.0

5.0 | 1 rating Rate this file 44 Downloads (last 30 days) File Size: 3.97 KB File ID: #32725
image thumbnail

Fast Alpha Hulls (alpha shapes in 3d; parfor enabled)

by

 

29 Aug 2011 (Updated )

Compute the alpha hulls (exterior and interior) of a set of points.

| Watch this File

File Information
Description

See also http://www.dylan-muir.com/articles/alpha_hulls/

Usage: [triHull, vbOutside, vbInside] = AlphaHull(mfPoints, fAlphaRadius , triDelaunay)

This function computes the alpha shape / alpha hulls of a set of points; both the external hull as well as interior voids. The Matlab parfor construct is used by the function, so that this code will run quickly on a machine running several instances of the Matlab parallel computing toolbox.

This algorithm is based on qhull and the delaunay tetrahedralisation of the set of points. It will return a hull triangulation, and ignore points connected only by a line.

'mfPoints' is an Nx3 ma­trix, where each row de­fines a point in 3-space. AlphaHull will find the hull of the set of points in 'mfPoints'.

'fAlphaRadius' is a scalar dis­tance which de­fines the pa­ra­me­ter alpha of the alpha hull. This dis­tance is in­ter­preted as the ra­dius of a sphere that will "roll around" the sur­face, with the bound­ary of the sphere touch­ing one to three of the points in 'mfPoints'. The tri­an­gles of the De­lau­nay tetra­he­dral­i­sa­tion where the spere can fit with­out in­ter­sect­ing any other points will form part of the alpha hull.

The op­tional pa­ra­me­ter 'triDelaunay' can be used to pro­vide the De­lau­nary tetra­he­dral­i­sa­tion of the set of points, if it is known in ad­vance.

'triHull' will be a tri­an­gu­la­tion con­tain­ing tri­an­gles that fall ei­ther on the alpha hull sur­face, or on the in­side sur­face of an alpha void (a hole) within the point set. The boolean vec­tors 'vbOutside' and 'vbInside' de­fine which rows of 'triHull' de­fine "out­side" and "in­side" hulls. The sur­faces re­turned by AlphaHull will be con­vex to the space pa­ra­me­ter alpha.

'triHull' will be a Tx3 ma­trix, where each row ['p1' 'p2' 'p3'] de­fines a tri­an­gle on an alpha sur­face. The in­dices 'pn' refer to rows in 'mfPoints', and so de­fine tri­an­gles in­clud­ing three of the orig­i­nal points.

* Caveats and room for improvement *
The method for la­belling "in­side" and "out­side" tri­an­gles is not ideal. It works by de­cid­ing whether the nor­mal of a tri­an­gle, in the di­rec­tion away from the rest of the point cloud, points in the di­rec­tion of the point set cen­troid. A bet­ter tech­nique might be to it­er­ate along the sur­face, la­belling tri­an­gles con­sis­tently as you go. If you im­prove on this, I'd love to hear about it.

MATLAB release MATLAB 7.7 (R2008b)
Tags for This File   Please login to tag files.
Please login to add a comment or rating.
Comments and Ratings (3)
27 Feb 2012 Dylan Muir

Hi Borc,

Your point cloud contains many points and triangles that are co-planar; this causes the determination of whether or not the triangles are on the alpha surface to perform badly. If you add a tiny jitter to each point:

cloud = cloud + (rand(size(cloud))-0.5)*0.001;

then everything works as it should. If you have a suggestion for sorting out this corner case, I'd love to hear it.

Best regards,
Dylan.

27 Feb 2012 Borc

Hi,

it seems that the code isn't working well for the Data below. Please check the AlphaHull function.

kind regards

%3D Example - Ring
[x,y,z] = sphere;
ii = abs(z) < 0.4;
X = [x(ii),y(ii),z(ii)];
cloud = [X; 0.8*X];

% loop over different values for r
radiusVals = [0.01 0.1 0.5 1 2];
for radiusIx = 1:length(radiusVals)
figure(radiusIx);
clf;
radius = radiusVals(radiusIx);
fprintf('Radius is %.2g\n',radius);
tic
[triHull, vbOutside, vbInside] = AlphaHull(cloud,radius);
fprintf(' * Got hull in %.1f seconds\n',toc);
plot3(cloud(:,1),cloud(:,2),cloud(:,3),'b.');
hold on;
fprintf(' * Found %d points for the outer surface\n',sum(vbOutside));
if (sum(vbOutside) > 0)
trisurf(triHull(vbOutside,:),cloud(:,1),cloud(:,2),cloud(:,3),...
'FaceColor','cyan','FaceAlpha',0.5)
end
fprintf(' * Found %d points for the inner surface\n',sum(vbInside));
if (sum(vbInside)>0)
trisurf(triHull(vbInside,:),cloud(:,1),cloud(:,2),cloud(:,3),...
'FaceColor','green','FaceAlpha',0.5)
end
grid on;
axis equal
title(sprintf('Radius = %.1g',radius))
end

01 Feb 2012 Patrik Eschle

Worked nicely out of the box - thank you very much.

I've made a little demo that you may include.

function alphaHull_demo
% Demo for the alpha hull function
%

% define a noisy torus
N=500;
t=linspace(0,1,N)'*2*pi;
cloud=[sin(t),cos(t),ones(size(t))]+(randn(N,3)-0.5)/10;

% loop over different values for r
radiusVals = [0.01 0.1 0.5 1 2];
for radiusIx = 1:length(radiusVals)
figure(radiusIx);
clf;
radius = radiusVals(radiusIx);
fprintf('Radius is %.2g\n',radius);
tic
[triHull, vbOutside, vbInside] = AlphaHull(cloud,radius);
fprintf(' * Got hull in %.1f seconds\n',toc);
plot3(cloud(:,1),cloud(:,2),cloud(:,3),'b.');
hold on;
fprintf(' * Found %d points for the outer surface\n',sum(vbOutside));
if (sum(vbOutside) > 0)
trisurf(triHull(vbOutside,:),cloud(:,1),cloud(:,2),cloud(:,3),...
'FaceColor','cyan','FaceAlpha',0.5)
end
fprintf(' * Found %d points for the inner surface\n',sum(vbInside));
if (sum(vbInside)>0)
trisurf(triHull(vbInside,:),cloud(:,1),cloud(:,2),cloud(:,3),...
'FaceColor','green','FaceAlpha',0.5)
end
grid on;
axis equal
title(sprintf('Radius = %.1g',radius))
end

Updates
30 Sep 2011

Updated image

07 Nov 2011

Updated description

17 Jan 2014

Updated description

28 Mar 2014

Updated description

Contact us