The QHULL procedure constructs convex hulls, Delaunay triangulations, and Voronoi diagrams of a set of points of 2-dimensions or higher. It uses and is based on the program QHULL, which is described in Barber, Dobkin and Huhdanpaa, “The Quickhull Algorithm for Convex Hulls,” ACM Transactions on Mathematical Software, Vol. 22, No 4, December 1996, Pages 469-483.
Note: QHULL accepts complex input but only utilizes the real part of any complex number.
For more information about QHULL see http://www.qhull.org
To run this example, copy the text of the following routine, paste it into an IDL editor window, save the window’s contents as ex_qhull.pro, and select Run ex_qhull from the Run menu.
; Create a collection of random points.
n = 20
seed = 15
x = RANDOMU(seed, n)
y = RANDOMU(seed, n)
; Construct the Delaunay triangulation
; and the Voronoi diagram.
QHULL, x, y, triangle, /DELAUNAY, $
VDIAGRAM=vdiagram, VVERTICES=vvert, VNORM=vnorm
; Plot our input points.
PLOT, [-0.1, 1.1], [-0.1, 1.1], /NODATA, $
PLOTS, x, y, PSYM=4
; Plot the Voronoi diagram.
FOR i=0,N_ELEMENTS(vdiagram[2,*])-1 DO BEGIN
vdiag = vdiagram[*, i]
j = vdiag + 1
; Bounded or unbounded?
PLOTS, vvert[*, vdiag[2:3]], PSYM=-5
; Unbounded, retrieve starting vertex.
xystart = vvert[*, vdiag]
; Determine the line equation.
; Vnorm*x + Vnorm*y + Vnorm = 0
slope = -vnorm[0,-j]/vnorm[1,-j]
intercept = -vnorm[2,-j]/vnorm[1,-j]
; Need to determine the line direction.
; Pick a point on one side along the line.
xunbound = xystart + 5
yunbound = slope*xunbound + intercept
; Find the closest original vertex.
void = MIN( (x-xunbound)^2 + (y-yunbound)^2, idx)
; By definition of Voronoi diagram, the line should
; be closest to one of the bisected points. If not,
; our line went in the wrong direction.
xunbound = xystart - 5
yunbound = slope*xunbound + intercept
PLOTS, [[xystart], [xunbound, yunbound]]
For some other examples using the QHULL procedure, see the QGRID3 function.
QHULL, V, Tr
QHULL, V0 , V1, [, V2 ... [, V6] ] , Tr
An input argument providing an nd-by-np array containing the locations of np points, in nd dimensions. The number of dimensions, nd, must be greater than or equal to 2.
V0, V1, V2, ..., V(N–1)
Input vectors of dimension np-by-1 elements each containing the i-th coordinate of np points in nd dimensions. A maximum of seven input vectors may be specified.
A named variable that will contain an nd1-by-nt array containing the indices of either the convex hull (nd1 is equal to nd), or the Delaunay triangulation (nd1 is equal to nd+1) of the input points.
Set this keyword equal to a named variable that will contain the indices of the vertices of the convex hull.
Note: The order of the vertices returned in this variable is unspecified.
Set this keyword equal to a named variable that will contain the adjacency list for each of the np nodes. The list has the following form:
Each element i, 0 ≤i < np, contains the starting index of the connectivity list (list) for node i within the list array. To obtain the adjacency list for node i, extract the list elements from list[i] to list[i+1] – 1. The adjacency list is not ordered. To obtain the connectivity list, either the DELAUNAY or SPHERE keywords must also be specified.
For example, to perform a spherical triangulation, use the following procedure call:
QHULL, lon, lat, CONNECTIVITY = list, SHPERE = sphere
which returns the adjacency list in the variable list. The subscripts of the nodes adjacent to lon[i] and lat[i] are contained in the array: list[list[i] :list[i+1] – 1].
Set this keyword to perform a Delaunay triangulation and returns the vertex indices of the resulting polyhedra; otherwise, the convex hull of the data are returned.
Set this keyword equal to a named variable that will contain the Delaunay triangulation of the points which lie on the surface of a sphere. The V0 argument contains the longitude, in degrees, and V1 contains the latitude, in degrees, of each point.
Set this keyword equal to a named variable that will contain the connectivity array for the Voronoi diagram.
For two-dimensional arrays, VDIAGRAM is a 4-by-nv integer array. For each Voronoi ridge, i, VDIAGRAM[0:1, i] contains the index of the two input points the ridge bisects. VDIAGRAM [2:3, i] contains the indices within VVERTICES of the Voronoi vertices. In the case of an unbounded half-space, VDIAGRAM[2, i] is set to a negative index, j, indicating that the corresponding Voronoi ridge is unbounded, and that the equation for the ridge is contained in VNORMAL[*, –j-1], and starts at Voronoi vertex [3, i].
For three-dimensional or higher dimensional arrays, VDIAGRAM is returned as a connectivity vector. This vector has the form [n, v0, v1, i0, i1, ..., in-3], where n is the number of points needed to describe that particular Voronoi ridge, v0 and v1 contain the indices for the two input points that the ridge bisects, and i0...in -3 contain the indices within VVERTICES of the Voronoi vertices. In the case of an unbounded half-space, VDIAGRAM[i] is set to a negative index, j, indicating that the corresponding Voronoi ridge is unbounded, and that the equation for the ridge is contained in VNORMAL[*, –j-1].
Set this keyword equal to a named variable that will contain the normals of each Voronoi ridge that is unbounded. The normals consist of a (nd+1)-by-nu array, where nd is the number of dimensions and nu is the number of unbounded vertices. Each row contains the equation for the unbounded ridge in the form:
V0X0 + V1X1 + V2X2 + ... + VndXnd + Vnd+1 = 0
where V* are the elements of the row within VNORMALS, and X* are the multidimensional coordinates. See the preceding description of VDIAGRAM for details.
Set this keyword equal to a named variable that will contain the Voronoi vertices.