MATLAB Examples

# Theory of convex hulls, Delaunay triangulations and Voronoi diagrams

In this section theoretical details about convex hulls, Delaunay triangulations and Voronoi Diagrams, as well as their implementation in this package are presented.

## Convex hull

help convhull_nd 
 Calculate the CONVex HULL of a set of points in N Dimensions. Description #varargout#=convhull_nd(#varargin#) gives the facets of the convex hull of given points. Definitions: (1) A set S is convex if for any two points p,q in S, the line segment pq lies into S. (2) For a subset P of Rd, the convex hull conv(P) is defined as: (2.1) The (unique) smallest convex set C in Rd containing P (P<=C) (2.2) The intersection of all convex sets containing P (2.3) The set of all convex combinations of points in P (2.4) The union of all simplices with vertices in P (3) The convex hull computation means the determination of conv(S) for a given finite set of n points (#points#) in Rd. Properties: The convex-hull operator conv() has the characteristic properties of a closure operator: (1) It is extensive, meaning that the convex hull of every set P is a superset of P. (2) It is non-decreasing, meaning that, for every two sets X and Y with X<=Y, the convex hull of X is a subset of the convex hull of Y. (3) It is idempotent, meaning that for every P, the convex hull of the convex hull of P is the same as the convex hull of P. (4) The Delaunay triangulation of a point set and its dual, the Voronoi diagram, are mathematically related to convex hulls: the Delaunay triangulation of a point set in R(#d#) can be viewed as the projection of a convex hull in R(#d#+1) (Brown, K. Q., 1979). Algorithm: To compute the convex hull (#faces#) of a set of n points (#points#) this algorithm uses the following scheme: () set #faces#, #cf# and #df# to an arbitrary simplex defined by the first #d#+1 points in #points# () for each face #k# of the convex hull created (#faces#) () define the simplex formed by the face #k# and the remaining point not belonging to the hyperplane defined by the face () calculate its volume () if the volume is negative () reverse the order of the last two vertices of the #k#th face in #faces# to change the volume sign () change the sign of the plane coefficients #cf# and #df# of the #k#th face in #faces# () end if () end for () find the relative distance of the remaining points from the center of the point set. () sort the remaining points from the highest to the lowest relative distance from the center of the point set. () set a remaining point vector (#pleft#) equal to the above order () while #pleft# is not empty () select the first point #i# from #pleft#. Delete it from #pleft#. () * find which of the faces of the current convex hull are visible from this point () if there are any visible faces (#visible#) () for each visible face #v# from #visible# () find the nonvisible faces #u# connected to (i.e. having one common edge with) the visible face #v#, if any () for each nonvisible face #k# from #u# () add the boundary between the visible face #v# and the nonvisible face #k# connected to the visible face #v# to the horizon (#horizon#) () end for () end for () delete the visible faces (#visible#) and the corresponding plane coefficients #cf# and #df# () for each boundary #j# from #horizon# () add the face defined by horizon and the selected point #i# to the convex hull (#faces#) () calculate the plane coefficients of the above face and store them appropriately in #cf# and #df# () end for () for each face #k# from the new faces added to the convex hull by the previous for-loop () choose a point non-coplaner to the face () define the simplex formed by the face #k# and the above point () calculate its volume () if the volume is negative () reverse the order of the last two vertices of the #k#th face to change the volume sign () change the sign of the plane coefficients #cf# and #df# of the #k#th face () end if () end for clear #horizon# () end if () if the number of iterations has reached the specified number #iter2del# () ** delete the points which are not visible by all of the faces of the current convex hull () end if () end while Input arguments #varargin# ({#points# #iter2del#}) contains the following arguments: #points# ([#n# x #d#]) is a matrix containing the coordinates of the points the convex hull of which is to be calculated. #n# is the number of points and #d# is the dimensionality of the problem. #iter2del# (scalar) is the number of iterations every which deletion of the points inside the current convex hull occurs. It is recommended that it is specified (i.e. it is assigned a finite value) if the number of points becomes relatively large. There is an optimum value for this variable in general. Low values lead to increased number of internal point deletions which minimizes the number of points deleted in each repetition of the deletion process, and is generally more computationally costly if the number of points is small; the same number of iterations (one for each point) would have much lower computational cost. On the other hand, if this argument takes large values, more iterations take place before each deletion, which would possibly check internal points and would not have been done if these points were deleted. The final selection of the value of this parameter is judged according to the ratio of the execution time of step ()* to the execution time of step ()** shown in the algorithm above. If this ratio is large, #iter2del# should be decreased, and in the opposite case #iter2del# should be increased. If the points belonging to the convex hull is a significant portion of the initial point set, point deletion is not large, and consequently this parameter has little or no effect. It should be ignored or set to a large value. The opposite happens if the portion of the initial point set belonging to the convex hull is small. Output arguments #varargout# ({#faces# #cf# #df# #Vol#}) contains the following arguments: #faces# ([#nfaces# x #d#]) is a matrix containing the identities of the points defining each facet of the convex hull calculated. #nfaces# is the number of facets of the convex hull. #cf# ([#nfaces# x #d#]) contains the coefficients of the planes defined by #faces#. Each row of #cf# contains the coefficients of the corresponding row of #faces#. The equation of the plane #i# is: #cf#(#i#,:)*[x1;x2;...xd]+#df#(#i#)=0 where #d# is the dimension of the problem. #df# ([#nfaces# x 1]) contains the constant terms of the planes defined by #faces#. Each row of #df# contains the constant term of the corresponding row of #faces#. The equation of the plane #i# is: #cf#(#i#,:)*[x1;x2;...xd]+#df#(#i#)=0 where #d# is the dimension of the problem. #Vol# (scalar) is the volume of the convex hull. Example: [X,Y,Z]=meshgrid(0:1); points=[X(:),Y(:),Z(:)]; points=points+0.001*rand(size(points)); [faces,~,~,V]=convhull_nd(points) trisurf(faces,points(:,1),points(:,2),points(:,3)) Parents (calling functions) delaunay_nd > convhull_nd voronoi_nd > convhull_nd Children (called functions) convhull_nd > ismembc convhull_nd > plane_nd 

## Delaunay triangulation

help delaunay_nd 
 Calculate the DELAUNAY triangulation of a set of points in N Dimensions. Description #DT#=delaunay_nd(#varargin#) gives the point identities defining the Delaunay triangulation of #points#. Definitions: (1) A triangulation of a finite point set S is a subdivision of the #d#-dimensional domain of the convex hull of S, whose bounded faces are #d#-1 entities (each defined by #d# points) and whose vertices are the points of S. (2) For a point set P in the (#d#-dimensional) Euclidean space, a Delaunay triangulation is a triangulation DT(P) such that no point in P is inside the circum-hypersphere of any simplex in DT(P). There exists a unique Delaunay triangulation for P, if P is a set of points in general position; that is, there exists no k-flat containing k + 2 points nor a k-sphere containing k + 3 points, for 1<=k<=#d#-1 (e.g., for a point set in R3; no three points are on a line, no four on a plane, no four are on a circle, and no five on a sphere). Properties: (1) The union of all simplices in the triangulation is the convex hull of the points. (2) The Delaunay triangulation contains O(#n#^(#d#/2)) simplices. (3) In the plane (#d#=2), if there are b vertices on the convex hull, then any triangulation of the points has at most 2#n#-2-b triangles, plus one exterior face. (4) In the plane, each vertex has on average six surrounding triangles. (5) In the plane, the Delaunay triangulation maximizes the minimum angle. Compared to any other triangulation of the points, the smallest angle in the Delaunay triangulation is at least as large as the smallest angle in any other. However, the Delaunay triangulation does not necessarily minimize the maximum angle. The Delaunay triangulation also does not necessarily minimize the length of the edges. (6) A circle circumscribing any Delaunay triangle does not contain any other input points in its interior. (7) If a circle passing through two of the input points does not contain any other of them in its interior, then the segment connecting the two points is an edge of a Delaunay triangulation of the given points. (8) Each simplex of the Delaunay triangulation of a set of points in #d#-dimensional space corresponds to a facet of convex hull of the projection of the points onto a (#d#+1)-dimensional paraboloid, and vice versa (Brown, K. Q., 1979). (9) The closest neighbor b to any point p is on an edge bp in the Delaunay triangulation. (10) The shortest path between two vertices, along Delaunay edges, is known to be no longer than $\frac{4\pi}{3\sqrt{3}} \approx 2.418$ times the Euclidean distance between them. The plane tangent to the point (x0,y0,z0) of the paraboloid is used to find the point where this tangent plane crosses the (#d#+1)th axis (w). The facets of the convex hull which are visible from this point are projected appropriately to give the Delaunay triangulation. Input arguments #varargin# ({#points# #iter2del#}) contains the following arguments: #points# ([#n# x #d#]) is a matrix containing the coordinates of the points the Delaunay triangulation of which is to be calculated. #n# is the number of points and #d# is the dimensionality of the problem. #iter2del# (scalar) is a parameter controlling internal point deletion at the convex hull algorithm. Type |help convhull_nd| for more details. Output arguments #DT# ([#u# x #d#]) is a matrix containing the identities of the points defining each simplex of the Delaunay triangulation calculated. #u# is the number of simplices of the Delaunay triangulation. Example: [X,Y,Z]=meshgrid(0:1); points=[X(:),Y(:),Z(:)]; points=points+0.001*rand(size(points)); DT=delaunay_nd(points); % Plot the mesh tetramesh(DT,points) % Calculate the volume of the mesh V=0; for i=1:size(DT,1) V=V+1/6*abs(det([points(DT(i,:),:) ones(4,1)])); end Parents (calling functions) Children (called functions) delaunay_nd > convhull_nd delaunay_nd > convhull_nd > ismembc delaunay_nd > convhull_nd > plane_nd 

## Voronoi diagram

help voronoi_nd 
 Calculate the vertices of the VORONOI diagram in N Dimensions. Description #V#=voronoi_nd(#varargin#) gives the coordinates of the vertices of the Voronoi diagram of #points#. Definitions: (1) A Voronoi diagram is a partitioning of a #n#-dimensional space into regions based on 'closeness' to points in a specific subset of the space. That set of points (called seeds, sites, or generators) is specified beforehand, and for each seed there is a corresponding region consisting of all points closer to that seed than to any other. These regions are called Voronoi cells. (2) The segments of the Voronoi diagram are all the points that are equidistant to the two nearest sites. (3) The Voronoi vertices (nodes) are the points equidistant to three (or more) sites. Properties: (1) The dual graph for a Voronoi diagram (in the case of a Euclidean space with point sites) corresponds to the Delaunay triangulation for the same set of points. (2) The closest pair of points corresponds to two adjacent cells in the Voronoi diagram. (3) Assume the setting is the Euclidean plane and a group of different points are given. Then two points are adjacent on the convex hull if and only if their Voronoi cells share an infinitely long side. Similarly for higher dimensions. Input arguments #varargin# ({#points# #iter2del#}) contains the following arguments: #points# ([#n# x #d#]) is a matrix containing the coordinates of the points the Voronoi diagram of which is to be calculated. #n# is the number of points and #d# is the dimensionality of the problem. #iter2del# (scalar) is a parameter controlling internal point deletion at the convex hull algorithm. Type |help convhull_nd| for more details. Output arguments #V# ([#u# x #d#]) is a matrix containing the coordinates of the Voronoi diagram vertices. #u# is the number of vertices. Example: points=[0 0;2 1;1 2;4 0;0 4;4 4]; V=voronoi_nd(points); % Plot the initial points and the Voronoi vertices scatter(points(:,1),points(:,2)) hold on scatter(V(:,1),V(:,2)) Parents (calling functions) Children (called functions) voronoi_nd > convhull_nd voronoi_nd > convhull_nd > ismembc voronoi_nd > convhull_nd > plane_nd 

## Contact author

(c) 2014 by George Papazafeiropoulos
First Lieutenant, Infrastructure Engineer, Hellenic Air Force
Civil Engineer, M.Sc., Ph.D. candidate, NTUA

Website: http://users.ntua.gr/gpapazaf/