What's an efficient way to find if a point lies in the convex hull of a point cloud?

Here is an easy solution that requires only scipy:

def in_hull(p, hull):
    """
    Test if points in `p` are in `hull`

    `p` should be a `NxK` coordinates of `N` points in `K` dimensions
    `hull` is either a scipy.spatial.Delaunay object or the `MxK` array of the 
    coordinates of `M` points in `K`dimensions for which Delaunay triangulation
    will be computed
    """
    from scipy.spatial import Delaunay
    if not isinstance(hull,Delaunay):
        hull = Delaunay(hull)

    return hull.find_simplex(p)>=0

It returns a boolean array where True values indicate points that lie in the given convex hull. It can be used like this:

tested = np.random.rand(20,3)
cloud  = np.random.rand(50,3)

print in_hull(tested,cloud)

If you have matplotlib installed, you can also use the following function that calls the first one and plots the results. For 2D data only, given by Nx2 arrays:

def plot_in_hull(p, hull):
    """
    plot relative to `in_hull` for 2d data
    """
    import matplotlib.pyplot as plt
    from matplotlib.collections import PolyCollection, LineCollection

    from scipy.spatial import Delaunay
    if not isinstance(hull,Delaunay):
        hull = Delaunay(hull)

    # plot triangulation
    poly = PolyCollection(hull.points[hull.vertices], facecolors='w', edgecolors='b')
    plt.clf()
    plt.title('in hull')
    plt.gca().add_collection(poly)
    plt.plot(hull.points[:,0], hull.points[:,1], 'o', hold=1)


    # plot the convex hull
    edges = set()
    edge_points = []

    def add_edge(i, j):
        """Add a line between the i-th and j-th points, if not in the list already"""
        if (i, j) in edges or (j, i) in edges:
            # already added
            return
        edges.add( (i, j) )
        edge_points.append(hull.points[ [i, j] ])

    for ia, ib in hull.convex_hull:
        add_edge(ia, ib)

    lines = LineCollection(edge_points, color='g')
    plt.gca().add_collection(lines)
    plt.show()    

    # plot tested points `p` - black are inside hull, red outside
    inside = in_hull(p,hull)
    plt.plot(p[ inside,0],p[ inside,1],'.k')
    plt.plot(p[-inside,0],p[-inside,1],'.r')

I would not use a convex hull algorithm, because you do not need to compute the convex hull, you just want to check whether your point can be expressed as a convex combination of the set of points of whom a subset defines a convex hull. Moreover, finding the convex hull is computationally expensive, especially in higher dimensions.

In fact, the mere problem of finding out whether a point can be expressed as a convex combination of another set of points can be formulated as a linear programming problem.

import numpy as np
from scipy.optimize import linprog

def in_hull(points, x):
    n_points = len(points)
    n_dim = len(x)
    c = np.zeros(n_points)
    A = np.r_[points.T,np.ones((1,n_points))]
    b = np.r_[x, np.ones(1)]
    lp = linprog(c, A_eq=A, b_eq=b)
    return lp.success

n_points = 10000
n_dim = 10
Z = np.random.rand(n_points,n_dim)
x = np.random.rand(n_dim)
print(in_hull(Z, x))

For the example, I solved the problem for 10000 points in 10 dimensions. The executions time is in the ms range. Would not want to know how long this would take with QHull.

Hi I am not sure about how to use your program library to achieve this. But there is a simple algorithm to achieve this described in words:

1. create a point that is definitely outside your hull. Call it Y
2. produce a line segment connecting your point in question (X) to the new point Y.
3. loop around all edge segments of your convex hull. check for each of them if the segment intersects with XY.
4. If the number of intersection you counted is even (including 0), X is outside the hull. Otherwise X is inside the hull.
5. if so occurs XY pass through one of your vertexes on the hull, or directly overlap with one of your hull's edge, move Y a little bit.
6. the above works for concave hull as well. You can see in below illustration (Green dot is the X point you are trying to determine. Yellow marks the intersection points.

First, obtain the convex hull for your point cloud.

Then loop over all of the edges of the convex hull in counter-clockwise order. For each of the edges, check whether your target point lies to the "left" of that edge. When doing this, treat the edges as vectors pointing counter-clockwise around the convex hull. If the target point is to the "left" of all of the vectors, then it is contained by the polygon; otherwise, it lies outside the polygon.

This other Stack Overflow topic includes a solution to finding which "side" of a line a point is on: Determine Which Side of a Line a Point Lies

Note that this will work only for convex polygons. But you're dealing with a convex hull, so it should suit your needs.

It looks like you already have a way to get the convex hull for your point cloud. But if you find that you have to implement your own, Wikipedia has a nice list of convex hull algorithms here: Convex Hull Algorithms

Use the equations attribute of ConvexHull:

def point_in_hull(point, hull, tolerance=1e-12):
    return all(
        (np.dot(eq[:-1], point) + eq[-1] <= tolerance)
        for eq in hull.equations)

In words, a point is in the hull if and only if for every equation (describing the facets) the dot product between the point and the normal vector (eq[:-1]) plus the offset (eq[-1]) is less than or equal to zero. You may want to compare to a small, positive constant tolerance = 1e-12 rather than to zero because of issues of numerical precision (otherwise, you may find that a vertex of the convex hull is not in the convex hull).

Demonstration:

import matplotlib.pyplot as plt
import numpy as np
from scipy.spatial import ConvexHull

points = np.array([(1, 2), (3, 4), (3, 6), (2, 4.5), (2.5, 5)])
hull = ConvexHull(points)

np.random.seed(1)
random_points = np.random.uniform(0, 6, (100, 2))

for simplex in hull.simplices:
    plt.plot(points[simplex, 0], points[simplex, 1])

plt.scatter(*points.T, alpha=.5, color='k', s=200, marker='v')

for p in random_points:
    point_is_in_hull = point_in_hull(p, hull)
    marker = 'x' if point_is_in_hull else 'd'
    color = 'g' if point_is_in_hull else 'm'
    plt.scatter(p[0], p[1], marker=marker, color=color)