Up: Home page for Qhull (http://www.geom.umn.edu/locate/qhull)
Up: Qhull manual

To: FAQ: Table of Contents (please wait while loading)


[4-d cube] Frequently Asked Questions about Qhull

If your question does not appear here, see:

Qhull is a general dimension code for computing convex hulls, Delaunay triangulations, halfspace intersections, 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, February 18, 1998

Copyright © 1998 The Geometry Center, Minneapolis MN


»FAQ: Table of Contents


»Startup questions

»How do I run Qhull from Windows?

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.

»How do I enter points for Qhull?

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.

»How do I learn to use Qhull?

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.

»Why does Qhull have so many options?

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.

»Delaunay and Voronoi diagram questions

»How do I construct a 3-d Delaunay triangulation?

The documentation is confusing for Delaunay triangulation and Voronoi diagrams. Qhull computes these structures by computing a convex hull in one higher dimension. In 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, option 'i' produces the vertices of the convex hull. These correspond to input sites of the Delaunay triangulation. A facet is non-simplicial if multiple input sites are cospherical. For 3-d Delaunay triangulations (a 4-d convex hull), option 'i' triangulates non-simplicial facets by adding a point to the facet. Option 'i' is not convenient in this case. As shown below, better choices are to prevent cospherical 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

»Can Qhull produce a triangular mesh for an object?

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.

»How do I get the triangles for a 2-d Delaunay triangulation and the vertices of its Voronoi diagram?

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

»Can Qhull compute the unbounded rays of the Voronoi diagram?

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.

»Can Qhull triangulate a hundred 16-d points?

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.

»Approximation questions

»How do I approximate data with a simplex

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

»Qhull library questions

»How do I list the vertices?

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

»How do I test code that uses the Qhull library?

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.

»When I compute a plane equation from a facet, I sometimes get an outward-pointing normal and sometimes an inward-pointing normal

Qhull orients simplicial facets, and prints oriented output for 'i', 'Ft', and other options. The orientation depends on both the vertex order and the flag facet->toporient.

Qhull does not orient non-simplicial facets. Instead it orients the facet's ridges. These are printed with the 'Ft' option. The facet's hyperplane is oriented.


Up: Home page for Qhull
Up: Qhull manual
To:
FAQ: Table of Contents


The Geometry Center Home Page

Comments to: qhull@geom.umn.edu
Created: Sept. 25, 1995 --- Last modified: see top