"ConvexHull", "alpha_shape" in python generate different area values why?

Question

I am trying to compute the area of an irregular polygon.

I used "ConvexHull" and "alpha_shape" but the output of the surface(area) is not the same.

any suggestions please?

you will find the data, and the code used to compute the area below:

###################################################################

data = [(1.866200,4.379100),
(1.128700,4.814700),
(2.036000,3.074600),
(1.077500,6.100000),
(1.833300,7.126600),
(2.061200,8.399900),
(2.710500,9.418000),
(3.426000,9.141000),
(2.795400,2.644300),
(2.773300,1.400700),
(3.554100,1.084300),
(4.338900,2.025600),
(4.365600,9.981400),
(5.284800,9.657300),
(6.305400,10.194400),
(6.829000,9.124800),
(5.889000,4.120900),
(4.884500,1.854300),
(5.694300,2.821400),
(6.825100,6.487400),
(6.698000,5.189600),
(7.676500,8.884200),
(7.670500,7.610700)]

####################################################

import numpy as np

from shapely.geometry import Polygon

from scipy.spatial import ConvexHull, convex_hull_plot_2d

####################################################

hull = ConvexHull(data)

print('area of the surface = ', hull.volume)

####################################################

polygon = Polygon(data)

print('area of the surface = ', polygon.area)

A convex hull surrounding a set of points is not the same as a polygon using the points as the vertices — Mad Physicist, Mar 24 '22 at 13:04
I am not sure what this question has to do with `alphashape` but [I found this library](https://stackoverflow.com/a/70928696/8881141) not well-maintained in the past — Mr. T, Mar 24 '22 at 13:19
@MadPhysicist The idea is to compute the area formed by set of the vertices of the polygon, that should be a little less than the convex hull surrounding the same set of vertices. the huge difference between the two values makes me suspect the answer. — Mohamad, Mar 24 '22 at 13:38
@Mr.T alphashape is used to visualize the area formed by the polygon vertices (to be calculated). so I used it to see what I am computing — Mohamad, Mar 24 '22 at 13:40
Unless the points are specified in the right order, the areas won't generally be close. You want a cw or ccw ordering about the center, not approximately sorted by x. — Mad Physicist, Mar 24 '22 at 13:41
I added an answer to show how to reorder the points correctly — Mad Physicist, Mar 24 '22 at 15:45
@MadPhysicist Thank a lot for the clarification, I thought that sorting should be along x (but never thought of ordering around the center). — Mohamad, Mar 25 '22 at 08:08

Mr. T · Answer 1 · 2022-03-24T14:30:03.280

Your shapely polygon is not a hull figure (it is not even convex), so unsurprisingly its area is lower:

import numpy as np    
from shapely.geometry import Polygon    
from scipy.spatial import ConvexHull, convex_hull_plot_2d
from matplotlib import pyplot as plt
   
hull = ConvexHull(data)    
print('area of the surface = ', hull.volume)
    
polygon = Polygon(data)    
print('area of the surface = ', polygon.area)

fig, ax = plt.subplots()
convex_hull_plot_2d(hull, ax=ax)
x, y = polygon.exterior.xy
ax.plot(x, y, c="red", ls="--", label="shapely")
ax.legend()
plt.show()

If we sort the data points according to the convex hull,

data_shapely = [
    (2.710500,9.418000),
    (4.365600,9.981400),
    (6.305400,10.194400),
    (7.676500,8.884200),
    (7.670500,7.610700),
    (6.698000,5.189600),
    (5.694300,2.821400),
    (4.884500,1.854300),
    (3.554100,1.084300),
    (2.773300,1.400700),
    (1.128700,4.814700),
    (1.077500,6.100000),
    (2.061200,8.399900)    
    ]

polygon = Polygon(data_shapely)

then the areas are the same:

area of the surface =  41.115790855
area of the surface =  41.11579085500001

You also misunderstand how shapely calculates areas of non-convex polygons.

data_shapely = [
    (0,0),
    (0,1),
    (1,1),
    (1,0)    
    ]

has unsurprisingly an area of 1. One might assume that

data_shapely = [
    (0,0),
    (1,1),
    (0,1),
    (1,0)    
    ]

has an area of 0.5, given its form. The result, however, is 0, as seemingly "flipped" areas in non-convex polygons are counted as negative values.

Exactly it should be a little lower, but not by a factor 1/5, would you think sorting data in an increasing order(x-axis) could solve the problem? — Mohamad, Mar 24 '22 at 13:44
No. See the comments by MadPhysicist. I have added also explanations about some misunderstanding on your side of how shapely works. — Mr. T, Mar 24 '22 at 14:23

Mad Physicist · Accepted Answer · 2022-03-24T16:07:43.700

If your points are well behaved and you want to get a convex hull that has approximately the same area as a polygon through the same points, you need to reorder the points. For example, you could order them by angle about the barycenter:

data = np.array(data)
center = data.mean(0)
angle = np.arctan2(*(data - center).T[::-1])
index = np.argsort(angle)
polygon2 = Polygon(data[index])

The plot looks much better now:

fig, ax = plt.subplots(1, 2)

ax[0].set_title("Original")
convex_hull_plot_2d(hull, ax=ax[0])
ax[0].plot(*polygon.exterior.xy, 'r--')

ax[1].set_title("Untangled")
convex_hull_plot_2d(hull, ax=ax[1])
ax[1].plot(*polygon2.exterior.xy, 'r--')
ax[1].plot(*center, 'kx')

plt.show()

And indeed you now get that the hull is only about 10% larger than the polygon:

>>> polygon2.area
37.25662339499999
>>> hull.volume
41.115790855
>>> polygon2.area / hull.volume * 100
90.61390434247064
>>> (hull.volume / polygon2.area - 1) * 100
10.35833929200876

@Mr.T. Of course not, I forgot to subtract off the barycenter! Fixed now, and thanks for the catch. — Mad Physicist, Mar 24 '22 at 15:36

"ConvexHull", "alpha_shape" in python generate different area values why?

2 Answers2