Products & Services Solutions Academia Support User Community Company

Learn more about MATLAB   

Delaunay Triangulation and Related Constructs

The Delaunay triangulation is the most widely used triangulation in scientific computing. While there are numerous algorithms for computing Delaunay triangulations, it is the favorable properties of the triangulation that has promoted widespread adoption. The fundamental property is the Delaunay criterion, or in the case of 2-D triangulations the empty circumcircle criterion. Given a set of points in 2-D, the Delaunay criterion ensures that the circumcircle associated with each triangle contains no other point in its interior. This defines a nearest-neighbor relation and leads to the creation of "well shaped" triangles. These characteristics have important implications in practice, and motivate the use of Delaunay triangulations in scattered data interpolation.

To create a 2-D Delaunay triangulation, start with a random 2-D dataset:

x=rand(10,1);
y=rand(10,1);
scatter(x,y)

The DelaunayTri constructor creates the triangulation object representing the Delaunay triangulation of a set of points.

dt=DelaunayTri(x,y);

Find the properties of the Delaunay triangulation object dt with the properties function:

properties(dt)

Properties for class DelaunayTri:

Constraints
X
Triangulation

The X property is the original set of data points.

isequal(dt.X,[x,y])

ans =

 1

The Triangulation property is a matrix representing the set of simplices that make up the triangulation. in this case, the simplices are triangles.

dt.Triangulation

ans =

 7 4 5
 7 3 6
 7 1 4
 9 5 4
 8 3 7
 1 7 6
 5 8 7
 2 8 5
 2 3 8
 2 5 9
 2 9 10

You can access the triangles by indexing directly into the dt object. For example, to return the vertices of the first triangle, use MATLAB's indexing capability.

ind = dt(1,:)

ind =

 7 4 5
tri1 = dt.X(ind,:)

tri1 =

0.2785 0.4218
0.9134 0.4854
0.6324 0.8003
scatter(tri1(:,1),tri1(:,2),'MarkerFaceColor','y')

You can remove the point from the triangulation:

dt.X(1,:)=[];
triplot(dt,'r')

You can use the methods of the DelaunayTri class to find other geometric representations of the point set. The Voronoi diagram of a discrete set of points X decomposes the space around each point X(i) into a region of influence R{i}. Locations within the region R{i} are closer to point i than any other point in X. The region of influence is the Voronoi region. The collection of all the Voronoi regions is the Voronoi diagram. [V, R] = voronoiDiagram(dt) returns the vertices V and regions R of the Voronoi diagram of the points dt.X. The region R{i} is a cell array of indices into V that represents the Voronoi vertices bounding the region. The Voronoi diagram of a set of points is the dual of the Delaunay triangulation;

hold off
dt=DelaunayTri(x,y);
voronoi(dt, '-.r')

The convex hull of a set of points is the smallest convex set that includes all the points. All triangulations of a set have the convex hull as their outer boundary. The convexHull method returns the indices of the vertices on the convex hull.

ch=convexHull(dt)

ch =

 1
 4
 9
 10
 2
 3
 6
 1

The nearestNeighbor method returns the index of nearest point in dt.X closest to a specified point. ThepointLocation method returns the triangle that encloses a specified point.

To see all methods of the DelaunayTri class, type:

methods DelaunayTri -full

While the Delaunay property is well defined, the topology of the triangulation is not unique in the presence of degenerate pointsets. In two dimensions, degeneracies arise when four or more unique points lie on the same circle. The vertices of a square, for example, have a nonunique Delaunay triangulation.

Constrained Delaunay Triangulation

Another property of the DelaunayTri object is Constraints. Given a set of points:

X = [0 0; 16 0; 16 2; 2 2; 2 3; 8 3; 8 5; 0 5];

and its Delaunay triangulation:

DT = DelaunayTri(X);

You can constrain the triangulation to use a particular set of edges.

First, plot the point set X:

subplot(2,1,1);
hold on;
scatter(X(:,1),X(:,2))
axis([-1 17 -1 6]);

Then, plot the Delaunay triangulation of X:

triplot(DT)

Given a set of constraining points C, you can make the triangulation use the edges in C. For example, the edges represented by the set C are plotted below.

subplot(2,1,2);
C = [1 2; 2 3; 3 4; 4 5; 5 6; 6 7; 7 8; 8 1];
	axis([-1 17 -1 6]);
plot(X(C'),X(C'+size(X,1)),'-r', 'LineWidth', 2);

Set the Constraints property of the DT object to C. The constraints in C require that there be an edge between points 1 and 2; between points 2 and 3; between points 3 and 4; etc.

DT.Constraints = C;

A constrained triangulation is one in which specified edges are required to be part of the triangulation. Since there is no edge between points 3 and 4 in the original triangulation, the Constraints property creates one. Plot the triangulation which is constrained by the edges represented by C:

scatter(X(:,1),X(:,2))
triplot(dt)

The top plot shows the original unconstrained Delaunay triangulation. The bottom plot shows the triangulation subject to the constraints C.

You can use the inOutStatus method to constrain a triangulation to the interior of a region. To return only the edges that are inside the region delineated by C:

hold off;
subplot(2,1,1);
hold on;
IO = inOutStatus(DT);
subplot(2,1,1);
hold on;
axis([-1 17 -1 6]);
plot(X(C'),X(C'+size(X,1)),'-r', 'LineWidth', 2);
hold off;
subplot(2,1,2);
hold on;
axis([-1 17 -1 6]);
triplot(DT(IO, :), DT.X(:,1), DT.X(:,2), 'LineWidth', 2)
hold off;

  


Recommended Products

Includes the most popular MATLAB recorded presentations with Q&A sessions led by MATLAB experts.

 © 1984-2009- The MathWorks, Inc.    -   Site Help   -   Patents   -   Trademarks   -   Privacy Policy   -   Preventing Piracy   -   RSS