Package ffx.numerics.quickhull

Class QuickHull3D

java.lang.Object

ffx.numerics.quickhull.QuickHull3D

public class QuickHull3D
extends Object

Computes the convex hull of a set of three dimensional points.

The algorithm is a three dimensional implementation of Quickhull, as described in Barber, Dobkin, and Huhdanpaa, ``The Quickhull Algorithm for Convex Hulls'' (ACM Transactions on Mathematical Software, Vol. 22, No. 4, December 1996), and has a complexity of O(n log(n)) with respect to the number of points. A well-known C implementation of Quickhull that works for arbitrary dimensions is provided by qhull.

A hull is constructed by providing a set of points to either a constructor or a build method. After the hull is built, its vertices and faces can be retrieved using getVertices and getFaces. A typical usage might look like this:

 // x y z coordinates of 6 points
 Point3d[] points = new Point3d[]{
     new Point3d(0.0, 0.0, 0.0),
     new Point3d(1.0, 0.5, 0.0),
     new Point3d(2.0, 0.0, 0.0),
     new Point3d(0.5, 0.5, 0.5),
     new Point3d(0.0, 0.0, 2.0),
     new Point3d(0.1, 0.2, 0.3),
     new Point3d(0.0, 2.0, 0.0),
 };

 QuickHull3D hull = new QuickHull3D();
 hull.build(points);

 System.out.println("Vertices:");
 Point3d[] vertices = hull.getVertices();
 for (int i = 0; i < vertices.length; i++) {
     Point3d pnt = vertices[i];
     System.out.println(pnt.x + " " + pnt.y + " " + pnt.z);
 }

 System.out.println("Faces:");
 int[][] faceIndices = hull.getFaces();
 for (int i = 0; i < faceIndices.length; i++) {
     for (int k = 0; k < faceIndices[i].length; k++) {
         System.out.print(faceIndices[i][k] + " ");
     }
     System.out.println("");
 }
 

As a convenience, there are also build and getVertex methods which pass point information using an array of doubles.

Robustness

Because this algorithm uses floating point arithmetic, it is potentially vulnerable to errors arising from numerical imprecision. We address this problem in the same way as href=http://www.qhull.org>qhull, by merging faces whose edges are not clearly convex. A face is convex if its edges are convex, and an edge is convex if the centroid of each adjacent plane is clearly below the plane of the other face. The centroid is considered below a plane if its distance to the plane is less than the negative of a distance tolerance. This tolerance represents the smallest distance that can be reliably computed within the available numeric precision. It is normally computed automatically from the point data, although an application may set this tolerance explicitly.

Numerical problems are more likely to arise in situations where data points lie on or within the faces or edges of the convex hull. We have tested QuickHull3D for such situations by computing the convex hull of a random point set, then adding additional randomly chosen points which lie very close to the hull vertices and edges, and computing the convex hull again. The hull is deemed correct if check returns true. These tests have been successful for a large number of trials and so we are confident that QuickHull3D is reasonably robust.

Merged Faces

The merging of faces means that the faces returned by QuickHull3D may be convex polygons instead of triangles. If triangles are desired, the application may triangulate the faces, but it should be noted that this may result in triangles which are very small or thin and hence difficult to perform reliable convexity tests on. In other words, triangulating a merged face is likely to restore the numerical problems which the merging process removed. Hence is it possible that, after triangulation, check will fail (the same behavior is observed with triangulated output from qhull).

Degenerate Input

It is assumed that the input points are non-degenerate in that they are not coincident, colinear, or colplanar, and thus the convex hull has a non-zero volume. If the input points are detected to be degenerate within the distance tolerance, an IllegalArgumentException will be thrown.

Since:

1.0

Author:

John E. Lloyd, Fall 2004, Michael J. Schnieders

Field Summary

Fields

Modifier and Type	Field	Description
static final double	AUTOMATIC_TOLERANCE	Specifies that the distance tolerance should be computed automatically from the input point data.
protected double	charLength
static final int	CLOCKWISE	Specifies that (on output) vertex indices for a face should be listed in clockwise order.
protected double	explicitTolerance
protected Vector<Face>	faces
protected int	findIndex
protected Vector<HalfEdge>	horizon
static final int	INDEXED_FROM_ONE	Specifies that (on output) the vertex indices for a face should be numbered starting from 1.
static final int	INDEXED_FROM_ZERO	Specifies that (on output) the vertex indices for a face should be numbered starting from 0.
protected int	numFaces
protected int	numPoints
protected int	numVertices
static final int	POINT_RELATIVE	Specifies that (on output) the vertex indices for a face should be numbered with respect to the original input points.
protected Vertex[]	pointBuffer
protected double	tolerance
protected int[]	vertexPointIndices

Constructor Summary

Constructors

Constructor	Description
QuickHull3D()	Creates an empty convex hull object.
QuickHull3D(double[] coords)	Creates a convex hull object and initializes it to the convex hull of a set of points whose coordinates are given by an array of doubles.
QuickHull3D(Point3d[] points)	Creates a convex hull object and initializes it to the convex hull of a set of points.

Method Summary

Modifier and Type	Method	Description
protected void	addNewFaces(FaceList newFaces, Vertex eyeVtx, Vector<HalfEdge> horizon)	Builds and links the ring of new faces around the horizon using the eye vertex, recording them in the provided newFaces list.
protected void	addPointToHull(Vertex eyeVtx)	Incorporates the specified eye vertex into the hull: computes the horizon, creates new faces, merges non-convex adjacencies, and resolves/reassigns outside points.
void	build(double[] coords)	Constructs the convex hull of a set of points whose coordinates are given by an array of doubles.
void	build(double[] coords, int nump)	Constructs the convex hull of a set of points whose coordinates are given by an array of doubles.
void	build(Point3d[] points)	Constructs the convex hull of a set of points.
void	build(Point3d[] points, int nump)	Constructs the convex hull of a set of points.
protected void	buildHull()	Builds the full convex hull by repeatedly selecting and adding extreme points until no outside points remain; then reindexes faces and vertices.
protected void	calculateHorizon(Point3d eyePnt, HalfEdge edge0, Face face, Vector<HalfEdge> horizon)	Recursively computes the horizon (boundary) edges visible from a given eye point, marking visited faces deleted and collecting the bordering edges.
boolean	check(PrintStream ps)	Checks the correctness of the hull using the distance tolerance returned by getDistanceTolerance; see check(PrintStream,double) for details.
boolean	check(PrintStream ps, double tol)	Checks the correctness of the hull.
protected boolean	checkFaceConvexity(Face face, double tol, PrintStream ps)	Verifies that a face is locally convex by checking distances between opposite face centroids and face planes, and ensuring no redundant vertices exist.
protected boolean	checkFaces(double tol, PrintStream ps)	Performs convexity validation across all visible faces using a specified tolerance.
protected void	computeMaxAndMin()	Computes per-axis min/max vertices and sets characteristic length and tolerance.
protected void	createInitialSimplex()	Creates the initial tetrahedral simplex from extremal points, throwing if points are coincident/colinear/coplanar within the current tolerance.
protected void	deleteFacePoints(Face face, Face absorbingFace)	Removes all outside vertices from a face and either discards them as unclaimed or reassigns them to an absorbing face if they lie above it.
double	getDistanceTolerance()	Returns the distance tolerance that was used for the most recently computed hull.
double	getExplicitDistanceTolerance()	Returns the explicit distance tolerance.
int[][]	getFaces()	Returns the faces associated with this hull.
int[][]	getFaces(int indexFlags)	Returns the faces associated with this hull.
int	getNumFaces()	Returns the number of faces in this hull.
int	getNumVertices()	Returns the number of vertices in this hull.
int[]	getVertexPointIndices()	Returns an array specifing the index of each hull vertex with respect to the original input points.
Point3d[]	getVertices()	Returns the vertex points in this hull.
int	getVertices(double[] coords)	Returns the coordinates of the vertex points of this hull.
protected void	initBuffers(int nump)	(Re)initializes internal buffers and lists for a hull build of nump points.
protected Vertex	nextPointToAdd()	Selects the next vertex to add: the farthest point above some face that currently has outside points claimed.
protected double	oppFaceDistance(HalfEdge he)	Computes the signed distance from the centroid of the opposite face across this edge to the plane of this edge's face.
void	print(PrintStream ps)	Prints the vertices and faces of this hull to the stream ps.
void	print(PrintStream ps, int indexFlags)	Prints the vertices and faces of this hull to the stream ps.
void	printPoints(PrintStream ps)	print all points to the print stream (very point a line)
protected void	reindexFacesAndVertices()	Removes deleted faces, marks active vertices, and assigns new contiguous indices to both faces and vertices for compact output.
protected void	resolveUnclaimedPoints(FaceList newFaces)	Assigns previously unclaimed vertices to the newly created faces most distant above them.
void	setExplicitDistanceTolerance(double tol)	Sets an explicit distance tolerance for convexity tests.
protected void	setHull(double[] coords, int nump, int[][] faceIndices, int numf)	Initializes the hull from precomputed face indices and point coordinates.
protected void	setPoints(double[] coords, int nump)	Populates the internal vertex buffer from a flat coordinate array.
protected void	setPoints(Point3d[] pnts, int nump)	Populates the internal vertex buffer from an array of Point3d.
void	triangulate()	Triangulates any non-triangular hull faces.

Methods inherited from class java.lang.Object

clone, equals, finalize, getClass, hashCode, notify, notifyAll, toString, wait, wait, wait

Field Details

CLOCKWISE

public static final int CLOCKWISE

Specifies that (on output) vertex indices for a face should be listed in clockwise order.

See Also:

• Constant Field Values

INDEXED_FROM_ONE

public static final int INDEXED_FROM_ONE

Specifies that (on output) the vertex indices for a face should be numbered starting from 1.

See Also:

• Constant Field Values

INDEXED_FROM_ZERO

public static final int INDEXED_FROM_ZERO

Specifies that (on output) the vertex indices for a face should be numbered starting from 0.

See Also:

• Constant Field Values

POINT_RELATIVE

public static final int POINT_RELATIVE

Specifies that (on output) the vertex indices for a face should be numbered with respect to the original input points.

See Also:

• Constant Field Values

AUTOMATIC_TOLERANCE

public static final double AUTOMATIC_TOLERANCE

Specifies that the distance tolerance should be computed automatically from the input point data.

See Also:

• Constant Field Values

findIndex

protected int findIndex

charLength

protected double charLength

pointBuffer

protected Vertex[] pointBuffer

vertexPointIndices

protected int[] vertexPointIndices

faces

protected Vector<Face> faces

horizon

protected Vector<HalfEdge> horizon

numVertices

protected int numVertices

numFaces

protected int numFaces

numPoints

protected int numPoints

explicitTolerance

protected double explicitTolerance

tolerance

protected double tolerance

Constructor Details

QuickHull3D

public QuickHull3D()

Creates an empty convex hull object.

QuickHull3D

public QuickHull3D(double[] coords)
            throws IllegalArgumentException

Creates a convex hull object and initializes it to the convex hull of a set of points whose coordinates are given by an array of doubles.

Parameters:

coords - x, y, and z coordinates of each input point. The length of this array will be three times the the number of input points.

Throws:

IllegalArgumentException - the number of input points is less than four, or the points appear to be coincident, colinear, or coplanar.

QuickHull3D

public QuickHull3D(Point3d[] points)
            throws IllegalArgumentException

Creates a convex hull object and initializes it to the convex hull of a set of points.

Parameters:

points - input points.

Throws:

IllegalArgumentException - the number of input points is less than four, or the points appear to be coincident, colinear, or coplanar.

Method Details

getDistanceTolerance

public double getDistanceTolerance()

Returns the distance tolerance that was used for the most recently computed hull. The distance tolerance is used to determine when faces are unambiguously convex with respect to each other, and when points are unambiguously above or below a face plane, in the presence of numerical imprecision. Normally, this tolerance is computed automatically for each set of input points, but it can be set explicitly by the application.

Returns:

distance tolerance

See Also:

• setExplicitDistanceTolerance(double)

setExplicitDistanceTolerance

public void setExplicitDistanceTolerance(double tol)

Sets an explicit distance tolerance for convexity tests. If AUTOMATIC_TOLERANCE is specified (the default), then the tolerance will be computed automatically from the point data.

Parameters:

tol - explicit tolerance

See Also:

• getDistanceTolerance()

getExplicitDistanceTolerance

public double getExplicitDistanceTolerance()

Returns the explicit distance tolerance.

Returns:

explicit tolerance

See Also:

• setExplicitDistanceTolerance(double)

setHull

protected void setHull(double[] coords,
 int nump,
 int[][] faceIndices,
 int numf)

Initializes the hull from precomputed face indices and point coordinates. This is primarily used for testing or reconstructing a hull state.

Parameters:

coords - flat array of x,y,z coordinates (length >= 3*nump)

nump - number of points in coords

faceIndices - face vertex indices (each row lists vertices of a face in CCW order)

numf - number of faces to add

printPoints

public void printPoints(PrintStream ps)

print all points to the print stream (very point a line)

Parameters:

ps - the print stream to write to

build

public void build(double[] coords)
           throws IllegalArgumentException

Constructs the convex hull of a set of points whose coordinates are given by an array of doubles.

Parameters:

coords - x, y, and z coordinates of each input point. The length of this array will be three times the number of input points.

Throws:

IllegalArgumentException - the number of input points is less than four, or the points appear to be coincident, colinear, or coplanar.

build

public void build(double[] coords,
 int nump)
           throws IllegalArgumentException

Constructs the convex hull of a set of points whose coordinates are given by an array of doubles.

Parameters:

coords - x, y, and z coordinates of each input point. The length of this array must be at least three times nump.

nump - number of input points

Throws:

IllegalArgumentException - the number of input points is less than four or greater than 1/3 the length of coords, or the points appear to be coincident, colinear, or coplanar.

build

public void build(Point3d[] points)
           throws IllegalArgumentException

Constructs the convex hull of a set of points.

Parameters:

points - input points

Throws:

IllegalArgumentException - the number of input points is less than four, or the points appear to be coincident, colinear, or coplanar.

build

public void build(Point3d[] points,
 int nump)
           throws IllegalArgumentException

Constructs the convex hull of a set of points.

Parameters:

points - input points

nump - number of input points

Throws:

IllegalArgumentException - the number of input points is less than four or greater then the length of points, or the points appear to be coincident, colinear, or coplanar.

triangulate

public void triangulate()

Triangulates any non-triangular hull faces. In some cases, due to precision issues, the resulting triangles may be very thin or small, and hence appear to be non-convex (this same limitation is present in qhull).

initBuffers

protected void initBuffers(int nump)

(Re)initializes internal buffers and lists for a hull build of nump points. Ensures pointBuffer capacity, clears face/claim lists, and sets counters.

Parameters:

nump - number of input points to prepare for

setPoints

protected void setPoints(double[] coords,
 int nump)

Populates the internal vertex buffer from a flat coordinate array.

Parameters:

coords - flat array of x,y,z coordinates (length >= 3*nump)

nump - number of points to load

setPoints

protected void setPoints(Point3d[] pnts,
 int nump)

Populates the internal vertex buffer from an array of Point3d.

Parameters:

pnts - array of Point3d

nump - number of points to load

computeMaxAndMin

protected void computeMaxAndMin()

Computes per-axis min/max vertices and sets characteristic length and tolerance. Used to seed the initial simplex and numerical thresholds.

createInitialSimplex

protected void createInitialSimplex()
                             throws IllegalArgumentException

Creates the initial tetrahedral simplex from extremal points, throwing if points are coincident/colinear/coplanar within the current tolerance.

Throws:

IllegalArgumentException - if a non-degenerate simplex cannot be formed

getNumVertices

public int getNumVertices()

Returns the number of vertices in this hull.

Returns:

number of vertices

getVertices

public Point3d[] getVertices()

Returns the vertex points in this hull.

Returns:

array of vertex points

See Also:

• getVertices(double[])
• getFaces()

getVertices

public int getVertices(double[] coords)

Returns the coordinates of the vertex points of this hull.

Parameters:

coords - returns the x, y, z coordinates of each vertex. This length of this array must be at least three times the number of vertices.

Returns:

the number of vertices

See Also:

• getVertices()
• getFaces()

getVertexPointIndices

public int[] getVertexPointIndices()

Returns an array specifing the index of each hull vertex with respect to the original input points.

Returns:

vertex indices with respect to the original points

getNumFaces

public int getNumFaces()

Returns the number of faces in this hull.

Returns:

number of faces

getFaces

public int[][] getFaces()

Returns the faces associated with this hull.

Each face is represented by an integer array which gives the indices of the vertices. These indices are numbered relative to the hull vertices, are zero-based, and are arranged counter-clockwise. More control over the index format can be obtained using getFaces(indexFlags).

Returns:

array of integer arrays, giving the vertex indices for each face.

See Also:

• getVertices()
• getFaces(int)

getFaces

public int[][] getFaces(int indexFlags)

Returns the faces associated with this hull.

Each face is represented by an integer array which gives the indices of the vertices. By default, these indices are numbered with respect to the hull vertices (as opposed to the input points), are zero-based, and are arranged counter-clockwise. However, this can be changed by setting POINT_RELATIVE, INDEXED_FROM_ONE, or CLOCKWISE in the indexFlags parameter.

Parameters:

indexFlags - specifies index characteristics (0 results in the default)

Returns:

array of integer arrays, giving the vertex indices for each face.

See Also:

• getVertices()

print

public void print(PrintStream ps)

Prints the vertices and faces of this hull to the stream ps.

This is done using the Alias Wavefront .obj file format, with the vertices printed first (each preceding by the letter v), followed by the vertex indices for each face (each preceded by the letter f).

The face indices are numbered with respect to the hull vertices (as opposed to the input points), with a lowest index of 1, and are arranged counter-clockwise. More control over the index format can be obtained using print(ps,indexFlags).

Parameters:

ps - stream used for printing

See Also:

• print(PrintStream, int)
• getVertices()
• getFaces()

print

public void print(PrintStream ps,
 int indexFlags)

Prints the vertices and faces of this hull to the stream ps.

This is done using the Alias Wavefront .obj file format, with the vertices printed first (each preceding by the letter v), followed by the vertex indices for each face (each preceded by the letter f).

By default, the face indices are numbered with respect to the hull vertices (as opposed to the input points), with a lowest index of 1, and are arranged counter-clockwise. However, this can be changed by setting POINT_RELATIVE, INDEXED_FROM_ZERO, or CLOCKWISE in the indexFlags parameter.

Parameters:

ps - stream used for printing

indexFlags - specifies index characteristics (0 results in the default).

See Also:

• getVertices()
• getFaces()

resolveUnclaimedPoints

protected void resolveUnclaimedPoints(FaceList newFaces)

Assigns previously unclaimed vertices to the newly created faces most distant above them.

Parameters:

newFaces - list of newly created faces during the current expansion

deleteFacePoints

protected void deleteFacePoints(Face face,
 Face absorbingFace)

Removes all outside vertices from a face and either discards them as unclaimed or reassigns them to an absorbing face if they lie above it.

Parameters:

face - face whose outside vertices are being removed

absorbingFace - face to receive reassigned vertices; may be null

oppFaceDistance

protected double oppFaceDistance(HalfEdge he)

Computes the signed distance from the centroid of the opposite face across this edge to the plane of this edge's face. Positive indicates non-convexity.

Parameters:

he - the half-edge whose adjacent faces are considered

Returns:

signed distance of opposite face centroid to this face plane

calculateHorizon

protected void calculateHorizon(Point3d eyePnt,
 HalfEdge edge0,
 Face face,
 Vector<HalfEdge> horizon)

Recursively computes the horizon (boundary) edges visible from a given eye point, marking visited faces deleted and collecting the bordering edges.

Parameters:

eyePnt - the point from which visibility is tested

edge0 - starting edge within the face (may be null to start at edge 0)

face - the current face being visited

horizon - output collection of horizon half-edges

addNewFaces

protected void addNewFaces(FaceList newFaces,
 Vertex eyeVtx,
 Vector<HalfEdge> horizon)

Builds and links the ring of new faces around the horizon using the eye vertex, recording them in the provided newFaces list.

Parameters:

newFaces - list that will collect the created faces

eyeVtx - the vertex being added to the hull

horizon - ordered list of horizon half-edges

nextPointToAdd

protected Vertex nextPointToAdd()

Selects the next vertex to add: the farthest point above some face that currently has outside points claimed.

Returns:

the next Vertex to add to the hull, or null if none remain

addPointToHull

protected void addPointToHull(Vertex eyeVtx)

Incorporates the specified eye vertex into the hull: computes the horizon, creates new faces, merges non-convex adjacencies, and resolves/reassigns outside points.

Parameters:

eyeVtx - the vertex to add to the hull

buildHull

protected void buildHull()

Builds the full convex hull by repeatedly selecting and adding extreme points until no outside points remain; then reindexes faces and vertices.

reindexFacesAndVertices

protected void reindexFacesAndVertices()

Removes deleted faces, marks active vertices, and assigns new contiguous indices to both faces and vertices for compact output.

checkFaceConvexity

protected boolean checkFaceConvexity(Face face,
 double tol,
 PrintStream ps)

Verifies that a face is locally convex by checking distances between opposite face centroids and face planes, and ensuring no redundant vertices exist.

Parameters:

face - the face to validate

tol - distance tolerance for convexity checks

ps - optional PrintStream for diagnostics (null to suppress)

Returns:

true if the face passes local convexity checks

checkFaces

protected boolean checkFaces(double tol,
 PrintStream ps)

Performs convexity validation across all visible faces using a specified tolerance.

Parameters:

tol - distance tolerance for convexity checks

ps - optional PrintStream for diagnostics (null to suppress)

Returns:

true if all visible faces pass convexity checks

check

public boolean check(PrintStream ps)

Checks the correctness of the hull using the distance tolerance returned by getDistanceTolerance; see check(PrintStream,double) for details.

Parameters:

ps - print stream for diagnostic messages; may be set to null if no messages are desired.

Returns:

true if the hull is valid

See Also:

• check(PrintStream, double)

check

public boolean check(PrintStream ps,
 double tol)

Checks the correctness of the hull. This is done by making sure that no faces are non-convex and that no points are outside any face. These tests are performed using the distance tolerance tol. Faces are considered non-convex if any edge is non-convex, and an edge is non-convex if the centroid of either adjoining face is more than tol above the plane of the other face. Similarly, a point is considered outside a face if its distance to that face's plane is more than 10 times tol.

If the hull has been triangulated, then this routine may fail if some of the resulting triangles are very small or thin.

Parameters:

ps - print stream for diagnostic messages; may be set to null if no messages are desired.

tol - distance tolerance

Returns:

true if the hull is valid

See Also:

• check(PrintStream)