0

I thought scipy.spatial.Delaunay.convex_hull is returning an array where every point/index is used twice, because one point belongs to two edges. But in my case there are several indices only one time:

hull = [[5053 6943]
        [6219 5797]
        [3441 5797]
        [7547 1405]
        [3441 6547]
        [8144 9215]
        [  47  444]
        [ 444 6219]
        [6547 5053]
        [9945 6943]
        [2695 1405]]

For example the "47" is only used once. What does this mean (geometrically) ? How can one point of a convex hull only be used for one edge?

Munchkin
  • 4,528
  • 7
  • 45
  • 93
  • `Delaunay.convex_hull` contains the triangle edges which Qhull reports as not having a neighbor. If the triangulation contains degenerate triangles, this set can contain artefacts. If you want to obtain numerically stable convex hull, use `scipy.spatial.ConvexHull`. – pv. Dec 02 '13 at 22:28

1 Answers1

0

As mentioned in the comment above, you should use the newer scipy.spatial.ConvexHull. If you look at the indices of the vertices returned from this method in my IPython examples below, you will see that they are not normally redundant for 2D or 3D data sets.

%pylab inline
import numpy
#use a simple trianglular data set:
triangle_points = numpy.array([[0,0],[-1,1],[1,1]])

#plot the triangle:
fig = plt.figure()
ax = fig.add_subplot(111,aspect='equal')
ax.scatter(triangle_points[...,0],triangle_points[...,1])

enter image description here

#now calculate the convex hull of the triangular point set:
import scipy, scipy.spatial
convex_hull_2D = scipy.spatial.ConvexHull(triangle_points)
#determine the indices of the vertices forming the convex hull (i.e., the output you get with the older scipy version):
convex_hull_vertex_indices = convex_hull_2D.vertices
#print the array of convex hull vertex indices:
convex_hull_vertex_indices
array([2, 1, 0], dtype=int32) #output

#For this simple 2D data set it is clear that scipy can define a convex hull using the index of each input point only once, so it is not doing AB, BC, CA as you might initially guess

#Let's see what happens if we add a point outside the plane of the triangle and make the data set 3D:
tetrahedral_points = numpy.array([[0,0,0],[-1,1,0],[1,1,0],[0,0.5,3]]) #the last element is the 'out of plane' new point
convex_hull = scipy.spatial.ConvexHull(tetrahedral_points)
convex_hull_vertex_indices = convex_hull.vertices
convex_hull_vertex_indices
array([0, 1, 2, 3], dtype=int32) #output  

#again, for a simple shape, even in 3D, scipy specifies the indices of the vertices of the convex hull in a non-redundant manner

#take a look at the simplices of the 2D ConvexHull object:
convex_hull_simplices = convex_hull_2D.simplices
convex_hull_simplices
array([[2, 1],
       [0, 1],
       [0, 2]], dtype=int32) #output

#so the simplices do contain duplicate indices to connect points to form edges/1-simplices
fig = plt.figure()
ax = fig.add_subplot(111,aspect='equal')
ax.set_ylim(0,1.2)
edge_colors = ['blue','red','green'] #color each 1-simplex differently for clarity
for simplex, edge_color in zip(convex_hull_simplices,edge_colors):
    plt.plot(triangle_points[simplex,0], triangle_points[simplex,1], c=edge_color)

enter image description here

treddy
  • 2,771
  • 2
  • 18
  • 31
  • sounds good, but I get an error and don't know why: `from scipy.spatial import ConvexHull import numpy working_points = numpy.array([[0,0],[-1,1],[1,1]]) conv_hull = ConvexHull(working_points) conv_hull_vert = conv_hull.vertices` I get: `AttributeError: 'ConvexHull' object has no attribute 'vertices'` – Munchkin Dec 03 '13 at 12:09
  • It works for me--you may have to use a more recent version of scipy as this version of ConvexHull is rather recent (I use 0.13.0 when the code in your comment works). I find that using pip is the easiest way to adjust the scipy version installed with: `pip install 'scipy==0.13.0'` Then `import scipy` `scipy.__version__` to make sure Python is accessing the correct version. – treddy Dec 03 '13 at 13:50
  • thank you. I know I have version 0.12.0, but I can't test it atm, will try it later. But I also have one last question: what's the difference between `ConvexHull().vertices` and `ConvexHull().simplices` in 2D ? – Munchkin Dec 03 '13 at 16:23
  • I've included an example with simplices now as well. You can see that there are duplicate indices here because these are the edges proper. This is also described to some extent in the standard documentation here: http://docs.scipy.org/doc/scipy-dev/reference/generated/scipy.spatial.ConvexHull.html – treddy Dec 17 '13 at 22:57