Curve Reconstruction from Point Cloud

Reconstruct a polygonal boundary from a cloud of points by using a Delaunay triangulation. The reconstruction is based on the elegant Crust algorithm.

Create a set of points representing the point cloud.

numpts = 192;
t = linspace(-pi,pi,numpts + 1)';
t(end) = [];
r = 0.1 + 5*sqrt(cos(6*t).^2 + 0.7.^2);
x = r.*cos(t);
y = r.*sin(t);
ri = randperm(numpts);
x = x(ri);
y = y(ri);

Construct a Delaunay Triangulation of the point set.

dt = delaunayTriangulation(x,y);
tri = dt(:,:);

Insert the location of the Voronoi vertices into the existing triangulation.

V = voronoiDiagram(dt);

Remove the infinite vertex and filter out duplicate points using unique.

V(1,:) = [];
numv = size(V,1);
dt.Points(end+(1:numv),:) = unique(V,"rows");

The Delaunay edges that connect pairs of sample points represent the boundary.

delEdges = edges(dt);
validx = delEdges(:,1) <= numpts;
validy = delEdges(:,2) <= numpts;
boundaryEdges = delEdges((validx & validy),:)';
xb = x(boundaryEdges);
yb = y(boundaryEdges);
clf
triplot(tri,x,y)
axis equal
hold on
plot(x,y,"*r")
plot(xb,yb,"r")
xlabel("Curve reconstruction from point cloud",FontWeight="b")
hold off

References

N. Amenta, M. Bern, and D. Eppstein. The crust and the beta-skeleton: combinatorial curve reconstruction. Graphical Models and Image Processing, 60:125-135, 1998.