QHULL

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/

Examples

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.

PRO ex_qhull

 

   ; 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, $

      XSTYLE=4, YSTYLE=4

   PLOTS, x, y, PSYM=4

 

   ; Plot the Voronoi diagram.

   FOR i=0,N_ELEMENTS(vdiagram[2,*])-1 DO BEGIN

      vdiag = vdiagram[*, i]

      j = vdiag[2] + 1

      ; Bounded or unbounded?

      IF (j gt 0) THEN BEGIN   ; Bounded.

         PLOTS, vvert[*, vdiag[2:3]], PSYM=-5

         CONTINUE

      ENDIF

 

      ; Unbounded, retrieve starting vertex.

      xystart = vvert[*, vdiag[3]]

      ; Determine the line equation.

      ; Vnorm[0]*x + Vnorm[1]*y + Vnorm[2] = 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[0] + 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.

      IF (idx ne vdiag[0] && idx ne vdiag[1]) THEN BEGIN

         xunbound = xystart[0] - 5

         yunbound = slope*xunbound + intercept

      ENDIF

 

      PLOTS, [[xystart], [xunbound, yunbound]]

 

   ENDFOR

 

END

For some other examples using the QHULL procedure, see the QGRID3 function.

Syntax

QHULL, V,  Tr

or,

QHULL, V0 , V1, [, V2 ... [, V6] ] , Tr  

Keywords: [, BOUNDS=variable ] [, CONNECTIVITY=variable ] [, /DELAUNAY ] [, SPHERE=variable ] [, VDIAGRAM=variable ] [, VNORMALS=variable ] [, VVERTICES=variable ]

Arguments

V

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.

Tr

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.

Keywords

BOUNDS

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.

CONNECTIVITY

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].

DELAUNAY

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.

SPHERE

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.

VDIAGRAM

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].

VNORMALS

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.

VVERTICES

Set this keyword equal to a named variable that will contain the Voronoi vertices.

Version History

5.5

Introduced

See Also

QGRID3