Up: Home page for Qhull
(http://www.geom.uiuc.edu/locate/qhull)
Up: Qhull manual
To: FAQ: Table of Contents (please
wait while loading)
If your question does not appear here, see:
Qhull is a general dimension code for computing convex hulls, Delaunay triangulations, halfspace intersections about a point, Voronoi diagrams, furthest-site Delaunay triangulations, and furthest-site Voronoi diagrams. These structures have applications in science, engineering, statistics, and mathematics. For a detailed introduction, see O'Rourke ['94], Computational Geometry in C.
Brad Barber, Cambridge MA, April 19, 1999
Copyright © 1998-1999 The Geometry Center, Minneapolis MN
Within each category, the most recently asked questions are first.
Qhull is a console program. You will first need a DOS window (i.e., a "DOS prompt"). You can double click on 'qhull-go.bat' in the Qhull directory. It loads 'doskey' to simplify rerunning qhull and rbox.
Qhull takes its data from standard input. For example, create a file named 'data.txt' with the following contents:
2 #sample 2-d input 5 #number of points 1 2 #coordinates of points -1.1 3 3 2.2 4 5 -10 -10
Then call qhull with 'qhull < data.txt'. It will print a summary of the convex hull. Use 'qhull < data.txt o' to print the vertices and edges. See also input format.
You can generate sample data with rbox, e.g., 'rbox 10' generates 10 random points in 3-d. Use a pipe ('|') to run rbox and qhull together, e.g.,
rbox c | qhull o
computes the convex hull of a cube.
First read:
Look at Qhull's on-line documentation:
Then try out Qhull on small examples.
Install Geomview if you are running SGI Irix, Solaris, SunOS, Linux, HP, IBM RS/6000, DEC Alpha, or Next. You can then visualize the output of Qhull. Qhull comes with Geomview examples.
Then try Qhull with a small example of your application. Work out the results by hand. Then experiment with Qhull's options to find the ones that you need.
You will need to decide how Qhull should handle precision problems. It can either joggle the input ('QJ') or merge facets ('C-0' or 'Qx'). Qhull will merge facets by default.
Qhull contains many options. A simpler choice would be a customized version for each application, but this is impractical to maintain. Each option has a purpose. Most of the options extract useful information from Qhull's data structures. The remaining options configure Qhull for an application, adjust Qhull's performance for experimentation, or help locate errors.
Delaunay triangulations and Voronoi diagrams cause the most trouble. Qhull computes these structures as convex hulls in one higher dimension. The terminology is confusing since the same options may be used for convex hulls, Delaunay triangulations, Voronoi diagrams, and halfspace intersections.
The best approach is to try examples from Synopsis and examples and Main options. Also try small inputs that you have computed by hand.
You may see extra points if you use option 'i' or 'Ft'. The extra points occur when a facet is non-simplicial (i.e., a facet with more than d vertices). For example, Qhull reports the following for the convex hull of a hypercube. It uses option 'Pd0:0.5' to return the facet along the positive-x axis:
rbox c D4 | qhull i Pd0:0.5 12 17 13 14 15 17 13 12 14 17 11 13 15 17 14 11 15 17 10 11 14 17 14 12 8 17 12 13 8 17 10 14 8 17 11 10 8 17 13 9 8 17 9 11 8 17 11 9 13
The 4-d hypercube has 16 vertices; so point "17" was added by qhull. Qhull adds the point in order to report a simplicial decomposition of the facet. The point corresponds to the "centrum" which Qhull computes to test for convexity.
Use the 'Fv' option to print the vertices of simplicial and non-simplicial facets. For example, here is the same hypercube facet with option 'Fv' instead of 'i':
C:\qhull>rbox c D4 | qhull Pd0:0.5 Fv 1 8 9 10 12 11 13 14 15 8
The coordinates of the extra point are printed with the 'Ft' option.
rbox c D4 | qhull Pd0:0.5 Ft 4 17 12 3 -0.5 -0.5 -0.5 -0.5 -0.5 -0.5 -0.5 0.5 -0.5 -0.5 0.5 -0.5 -0.5 -0.5 0.5 0.5 -0.5 0.5 -0.5 -0.5 -0.5 0.5 -0.5 0.5 -0.5 0.5 0.5 -0.5 -0.5 0.5 0.5 0.5 0.5 -0.5 -0.5 -0.5 0.5 -0.5 -0.5 0.5 0.5 -0.5 0.5 -0.5 0.5 -0.5 0.5 0.5 0.5 0.5 -0.5 -0.5 0.5 0.5 -0.5 0.5 0.5 0.5 0.5 -0.5 0.5 0.5 0.5 0.5 0.5 0 0 0 4 16 13 14 15 4 16 13 12 14 4 16 11 13 15 4 16 14 11 15 4 16 10 11 14 4 16 14 12 8 4 16 12 13 8 4 16 10 14 8 4 16 11 10 8 4 16 13 9 8 4 16 9 11 8 4 16 11 9 13
There's no direct way. You can use options 'Qc FP' to report the distance to the nearest vertex for coplanar input points. Select the minimum distance for a duplicated vertex, and locate all input sites less than this distance.
For Delaunay triangulations, all coplanar points are nearly incident to a vertex. If you want a report of duplicate input sites, do not use option 'QJ'. By adding a small random quantity to each input coordinate, it prevents duplicate input sites.
For points in general position, a 3-d Delaunay triangulation generates tetrahedron. Each face of a tetrahedron is a triangle. For example, the 3-d Delaunay triangulation of random points on the surface of a cube, is a cellular structure of tetrahedron.
Use 'qhull d QJ i' to generate the Delaunay triangulation ('d'). This uses joggled input ('QJ') to guarantee tetrahedral output. Option 'i' reports each tetrahedron. The triangles are every combination of 3 vertices. Each triangle is a "ridge" of the Delaunay triangulation.
For example,
rbox 10 | qhull QJ d i 14 9 5 8 7 0 9 8 7 5 3 8 7 3 0 8 7 5 4 8 1 4 6 8 1 2 9 5 8 4 2 5 8 4 2 9 5 6 2 4 8 9 2 0 8 2 6 0 8 2 4 9 1 2 6 4 1
is the Delaunay triangulation of 10 random points. Ridge 9-5-8 occurs twice. Once for tetrahedron 9 5 8 7 and the other for tetrahedron 2 9 5 8.
You can also use the Qhull library to generate the triangles. See "How do I visit the ridges of a Delaunay triangulation?"
Three-d terrain data can be approximated with cospherical points. The Delaunay triangulation of cospherical points is the same as their convex hull. If the points lie on the unit sphere, the facet normals are the Voronoi vertices [via S. Fortune].
For example, consider the points {[1,0,0], [-1,0,0], [0,1,0], ...}. Their convex hull is:
rbox d G1 | qhull o 3 6 8 12 0 0 -1 0 0 1 0 -1 0 0 1 0 -1 0 0 1 0 0 3 3 1 4 3 1 3 5 3 0 3 4 3 3 0 5 3 2 1 5 3 1 2 4 3 2 0 4 3 0 2 5
The facet normals are:
rbox d G1 | qhull n 4 8 -0.5773502691896258 0.5773502691896258 0.5773502691896258 -0.5773502691896258 0.5773502691896258 0.5773502691896258 0.5773502691896258 -0.5773502691896258 -0.5773502691896258 0.5773502691896258 -0.5773502691896258 -0.5773502691896258 0.5773502691896258 0.5773502691896258 -0.5773502691896258 -0.5773502691896258 0.5773502691896258 -0.5773502691896258 0.5773502691896258 -0.5773502691896258 -0.5773502691896258 -0.5773502691896258 0.5773502691896258 -0.5773502691896258 -0.5773502691896258 -0.5773502691896258 -0.5773502691896258 -0.5773502691896258 0.5773502691896258 -0.5773502691896258 -0.5773502691896258 -0.5773502691896258
If you drop the offset from each line (the last number), each line is the Voronoi vertex for the corresponding facet. The neighboring facets for each point define the Voronoi region for each point. For example:
rbox d G1 | qhull FN 6 4 7 3 2 6 4 5 0 1 4 4 7 4 5 6 4 3 1 0 2 4 6 2 0 5 4 7 3 1 4
The Voronoi vertices {7, 3, 2, 6} define the Voronoi region for point 0. Point 0 is [0,0,-1]. Its Voronoi vertices are
-0.5773502691896258 0.5773502691896258 -0.5773502691896258 0.5773502691896258 0.5773502691896258 -0.5773502691896258 -0.5773502691896258 -0.5773502691896258 -0.5773502691896258 0.5773502691896258 -0.5773502691896258 -0.5773502691896258
In this case, the Voronoi vertices are oriented, but in general they are unordered.
By taking the dual of the Delaunay triangulation, you can construct the Voronoi diagram. For cospherical points, the convex hull vertices for each facet, define the input sites for each Voronoi vertex. In 3-d, the input sites are oriented. For example:
rbox d G1 | qhull i 8 3 1 4 1 3 5 0 3 4 3 0 5 2 1 5 1 2 4 2 0 4 0 2 5
The convex hull vertices for facet 0 are {3, 1, 4}. So Voronoi vertex 0 (i.e., [-0.577, 0.577, 0.577]) is the Voronoi vertex for input sites {3, 1, 4} (i.e., {[0,1,0], [0,0,1], [-1,0,0]}).
The documentation can be confusing for Delaunay triangulation and Voronoi diagrams. Qhull computes these structures by computing a convex hull in one higher dimension. When reading the documentation, you will need to translate from triangulations to the corresponding convex hull. O'Rourke's Computational Geometry in C contains a good introduction.
For example, options 'i', 'Fv', 'o', and 'Ft' produce the vertices of the convex hull. These correspond to input sites of the Delaunay triangulation. The facets of the convex hull correspond to triangles of the Delaunay triangulation. A facet is non-simplicial if multiple input sites are cospherical.
For 3-d Delaunay triangulations with cospherical input sites (a 4-d convex hull), options 'i' and 'Ft' are inconvenient. They triangulate non-simplicial facets by adding a point to the facet. As shown below, better choices are to prevent non-simplicial input sites with option 'QJ', or to report all cospherical input sites with options 'Fv' or 'o'.
If you want non-simplicial output for cospherical sites, use options 'd Qbb Fv' or 'd Qbb o'. For option 'o', ignore the last coordinate. It is the lifted coordinate for the corresponding convex hull in 4-d. Option 'Qbb' normalizes the lifted coordinate. It is especially important for coordinates greater than 100.0. The following example is a cube inside a tetrahedron. The 8-vertex facet is the cube. Ignore the last coordinates.
C:\qhull>rbox r y c G0.1 | qhull d Qbb Fv 4 12 20 44 0.5 0 0 0.3055555555555555 0 0.5 0 0.3055555555555555 0 0 0.5 0.3055555555555555 -0.5 -0.5 -0.5 0.9999999999999999 -0.1 -0.1 -0.1 -6.938893903907228e-018 -0.1 -0.1 0.1 -6.938893903907228e-018 -0.1 0.1 -0.1 -6.938893903907228e-018 -0.1 0.1 0.1 -6.938893903907228e-018 0.1 -0.1 -0.1 -6.938893903907228e-018 0.1 -0.1 0.1 -6.938893903907228e-018 0.1 0.1 -0.1 -6.938893903907228e-018 0.1 0.1 0.1 -6.938893903907228e-018 4 2 11 1 0 4 10 1 0 3 4 11 10 1 0 4 2 9 0 3 4 9 11 2 0 4 7 2 1 3 4 11 7 2 1 4 8 10 0 3 4 9 8 0 3 5 8 9 10 11 0 4 10 6 1 3 4 6 7 1 3 5 6 8 10 4 3 5 6 7 10 11 1 4 5 9 2 3 4 7 5 2 3 5 5 8 9 4 3 5 5 6 7 4 3 8 5 6 8 7 9 10 11 4 5 5 7 9 11 2
If you want simplicial output use options 'd QJ i' or 'd QJ Ft'. Option 'QJ' adds a small random quantity to the input points, e.g.,
C:\qhull>rbox r y c G0.1 | qhull d QJ Ft 3 12 32 64 0.499999999939516 -4.457283719588488e-011 3.092055050303195e-011 3.963849488044223e-012 0.4999999999660025 -5.479396782402126e-011 2.168946845951116e-011 5.258475653706346e-011 0.4999999999859073 -0.4999999999599631 -0.5000000000563027 -0.5000000000540177 -0.09999999997929608 -0.1000000000595537 -0.1000000000141032 -0.1000000000099817 -0.09999999997740613 0.1000000000107635 -0.09999999995812423 0.1000000000032576 -0.10000000004936 -0.1000000000101616 0.100000000024338 0.1000000000496365 0.09999999997126396 -0.1000000000547432 -0.0999999999714412 0.1000000000160453 -0.09999999996898205 0.1000000000594008 0.09999999996939933 0.1000000000583741 -0.0999999999730648 0.1000000000183292 0.09999999994830786 0.1000000000159239 4 2 11 1 0 4 11 7 2 1 4 10 1 0 3 4 11 10 1 0 4 10 7 11 1 4 2 9 0 3 4 9 11 2 0 4 7 9 11 2 4 8 10 11 0 4 9 8 11 0 4 8 10 0 3 4 9 8 0 3 4 5 9 7 2 4 7 5 2 3 4 4 5 7 3 4 5 9 2 3 4 9 5 7 11 4 5 8 9 3 4 8 5 4 3 4 8 5 9 11 4 5 6 4 7 4 10 6 7 1 4 6 10 7 11 4 5 6 7 11 4 5 6 8 4 4 6 7 1 3 4 6 4 7 3 4 10 6 1 3 4 6 8 10 11 4 6 5 8 11 4 8 6 10 3 4 6 8 4 3
Yes for convex objects, no for non-convex objects. For non-convex objects, it triangulates the concavities. Unless the object has many points on its surface, triangles may cross the surface.
To compute the Delaunay triangles indexed by the indices of the input sites, use
rbox 10 D2 | qhull d QJ i
To compute the Voronoi vertices and the Voronoi region for each input site, use
rbox 10 D2 | qhull v QJ o
To compute each edge ("ridge") of the Voronoi diagram for each pair of adjacent input sites, use
rbox 10 D2 | qhull v QJ Fv
To compute the area of the Voronoi region for input site 5, use
rbox 10 D2 | qhull v QJ QV5 Pg p | qhull s FS
To compute the lines ("hyperplanes") that define the Voronoi region for input site 5, use
rbox 10 D2 | qhull v QJ QV5 Pg p | qhull n
To list the extreme points of the input sites use
rbox 10 D2 | qhull d Fx
You will get the same point ids with
rbox 10 D2 | qhull Fx
Use 'Fo to compute the separating hyperplanes for unbounded Voronoi regions. The corresponding ray goes to infinity from the Voronoi vertices. If you enclose the input sites in a large enough box, the outermost bounded regions will represent the unbounded regions of the original points.
If you do not box the input sites, you can identify the unbounded regions. They list '0' as a vertex. Vertex 0 represents "infinity". Each unbounded ray includes vertex 0 in option 'Fv. See Voronoi graphics and Voronoi notes.
No. This is an immense structure. A triangulation of 19, 16-d points has 43 simplicies. If you add one point at a time, the triangulation increased as follows: 43, 189, 523, 1289, 2830, 6071, 11410, 20487. The last triangulation for 26 points used 13 megabytes of memory. When Qhull uses virtual memory, it becomes too slow to use.
Qhull may be used to help select a simplex that approximates a data set. It will take experimentation. Geomview will help to visualize the results. This task may be difficult to do in 5-d and higher. Use rbox options 'x' and 'y' to produce random distributions within a simplex. Your methods work if you can recover the simplex.
Use Qhull's precision options to get a first approximation to the hull, say with 10 to 50 facets. For example, try 'C0.05' to remove small facets after constructing the hull. Use 'W0.05' to ignore points within 0.05 of a facet. Use 'PA5' to print the five largest facets by area.
Then use other methods to fit a simplex to this data. Remove outlying vertices with few nearby points. Look for large facets in different quadrants. You can use option 'Pd0d1d2' to print all the facets in a quadrant.
In 4-d and higher, use the outer planes (option 'Fo' or 'facet->maxoutside') since the hyperplane of an approximate facet may be below many of the input points.
For example, consider fitting a cube to 1000 uniformly random points in the unit cube. In this case, the first try was good:
rbox 1000 | qhull W0.05 C0.05 PA6 Fo 4 6 0.35715408374381 0.08706467018177928 -0.9299788727015564 -0.5985514741284483 0.995841591359023 -0.02512604712761577 0.08756829720435189 -0.5258834069202866 0.02448099521570909 -0.02685210459017302 0.9993396046151313 -0.5158104982631999 -0.9990223929415094 -0.01261133513150079 0.04236994958247349 -0.509218270408407 -0.0128069014364698 -0.9998380680115362 0.01264203427283151 -0.5002512653670584 0.01120895057872914 0.01803671994177704 -0.9997744926535512 -0.5056824072956361
The Qhull library may be used to construct convex hulls and Delaunay triangulations one point at a time. It may not be used for deleting points or moving points.
Qhull is designed for batch processing. Neither Clarkson's randomized incremental algorithm nor Qhull is designed for on-line operation. For many applications, it is better to reconstruct the convex hull or Delaunay triangulation from scratch for each new point.
With random point sets and on-line processing, Clarkson's algorithm should run faster than Qhull. Clarkson uses the intermediate facets to reject new, interior points, while Qhull, when used on-line, visits every facet to reject such points. If used on-line for n points, Clarkson may take O(n) times as much memory as the average off-line case, while Qhull's space requirement does not change.
To visit the ridges of a Delaunay triangulation, visit each facet. Each ridge will appear twice since it belongs to two facets. In pseudo-code:
for each facet of the triangulation if the facet is Delaunay (i.e., part of the lower convex hull) for each ridge of the facet if the ridge's neighboring facet has not been visited ... process a ridge of the Delaunay triangulation ...
In undebugged, C code:
qh visit_id++; FORALLfacets_(facetlist) if (!facet->upperdelaunay) { facet->visitid= qh visit_id; qh_makeridges(facet); FOREACHridge_(facet->ridges) { neighbor= otherfacet_(ridge, facet); if (neighbor->visitid != qh visit_id) { /* Print ridge here with facet-id and neighbor-id */ /*fprintf(fp, "f%d\tf%d\t",facet->id,neighbor->id);*/ FOREACHvertex_(ridge->vertices) fprintf(fp,"%d ",qh_pointid (vertex->point) ); qh_printfacetNvertex_simplicial (fp, facet, format); fprintf(fp," "); if(neighbor->upperdelaunay) fprintf(fp," -1 -1 -1 -1 "); else qh_printfacetNvertex_simplicial (fp, neighbor, format); fprintf(fp,"\n"); } } } }
Qhull should be redesigned as a class library, or at least as an API. It currently provides everything needed, but the programmer has to do a lot of work. Hopefully someone will write C++ wrapper classes or a Python module for Qhull.
Qhull constructs a Delaunay triangulation but lifting the input sites to a paraboloid. The Delaunay triangulation corresponds to the lower convex hull of the lifted points. To visit each facet of the lower convex hull, use:
facetT *facet; ... FORALLfacets { if (!facet->upperdelaunay) { ... only Delaunay facets ... } }
A point is outside of a facet if it is clearly outside the facet's outer plane. The outer plane is defined by an offset (facet->maxoutside) from the facet's hyperplane.
facetT *facet; pointT *point; realT dist; ... qh_distplane(point, facet, &dist); if (dist > facet->maxoutside + 2 * qh DISTround) { /* point is clearly outside of facet */ }
A point is inside of a facet if it is clearly inside the facet's inner plane. The inner plane is computed as the maximum distance of a vertex to the facet. It may be computed for an individual facet, or you may use the maximum over all facets. For example:
facetT *facet; pointT *point; realT dist; ... qh_distplane(point, facet, &dist); if (dist < qh min_vertex - 2 * qh DISTround) { /* point is clearly inside of facet */ }
Both tests include two qh.DISTrounds because the computation of the furthest point from a facet may be off by qh.DISTround and the computation of the current distance to the facet may be off by qh.DISTround.
Use qh_findbestfacet(). For example,
coordT point[ DIM ]; boolT isoutside; realT bestdist; facetT *facet; ... set coordinates for point ... facet= qh_findbestfacet (point, qh_ALL, &bestdist, &isoutside); /* 'facet' is the closest facet to 'point' */
qh_findbestfacet() performs a directed search for the facet furthest below the point. If the point lies inside this facet, qh_findbestfacet() performs an exhaustive search of all facets. An exhaustive search may be needed because a facet on the far side of a lens-shaped distribution may be closer to a point than all of the facet's neighbors. The exhaustive search may be skipped for spherical distributions.
Also see, "How do I find the Delaunay triangle or Voronoi region that is closest to a point?"
Lift the point to the paraboloid by summing the squares of the coordinates. Use qh_findbestfacet() to find the closest Delaunay triangle. Determine the closest vertex to find the corresponding Voronoi region.
You first need to lift the point to the paraboloid. The routine, qh_setdelaunay(), lifts an array of points to the paraboloid. The following excerpt is from findclosest() in user_eg.c.
coordT point[ DIM + 1]; /* one extra coordinate for lifting the point */ boolT isoutside; realT bestdist; facetT *facet; ... set coordinates for point[] ... qh_setdelaunay (DIM+1, 1, point); facet= qh_findbestfacet (point, qh_ALL, &bestdist, &isoutside); /* 'facet' is the closest Delaunay triangle to 'point' */
The returned facet either contains the point or it is the closest Delaunay triangle along the convex hull of the input set. For a similar approach, see Mucke, et al, "Fast randomized point location without preprocessing in two- and three-dimensional Delaunay triangulations," Computational Geometry '96, p. 274-283, May 1996. Also see, "How do I find the facet that is closest to a point?"
To locate the closest Voronoi region, determine the closest vertex of the closest Delaunay triangle.
realT dist, bestdist= REALmax; vertexT *bestvertex= NULL, *vertex, **vertexp; /* 'facet' is the closest Delaunay triangle to 'point' */ FOREACHvertex_( facet->vertices ) { dist= qh_pointdist( point, vertex->point, DIM ); if (dist < bestdist) { bestdist= dist; bestvertex= vertex; } } /* 'bestvertex' represents the Voronoi region closest to 'point'. The corresponding input site is 'bestvertex->point' */
To list the vertices (i.e., extreme points) of the convex hull use
vertexT *vertex; FORALLvertices { ... // vertex->point is the coordinates of the vertex // qh_pointid (vertex->point) is the point id of the vertex ... }
Compare the output from your program with the output from the Qhull program. Use option 'T1' or 'T4' to trace what Qhull is doing. Prepare a small example for which you know the output. Run the example through the Qhull program and your code. Compare the trace outputs. If you do everything right, the two trace outputs should be almost the same. The trace output will also guide you to the functions that you need to review.
Up: Home page for Qhull
Up: Qhull manual
To: FAQ: Table of Contents
Comments to: qhull@geom.umn.edu
Created: Sept. 25, 1995 --- Last modified: see top