Calculate area of polygon given (x,y) coordinates
Implementation of Shoelace formula could be done in Numpy
. Assuming these vertices:
import numpy as npx = np.arange(0,1,0.001)y = np.sqrt(1-x**2)
We can redefine the function in numpy to find the area:
def PolyArea(x,y): return 0.5*np.abs(np.dot(x,np.roll(y,1))-np.dot(y,np.roll(x,1)))
And getting results:
print PolyArea(x,y)# 0.26353377782163534
Avoiding for
loop makes this function ~50X faster than PolygonArea
:
%timeit PolyArea(x,y)# 10000 loops, best of 3: 42 µs per loop%timeit PolygonArea(zip(x,y))# 100 loops, best of 3: 2.09 ms per loop.
Timing is done in Jupyter notebook.
The most optimized solution that covers all possible cases, would be to use a geometry package, like shapely, scikit-geometry or pygeos. All of them use C++ geometry packages under the hood. The first one is easy to install via pip:
pip install shapely
and simple to use:
from shapely.geometry import Polygonpgon = Polygon(zip(x, y)) # Assuming the OP's x,y coordinatesprint(pgon.area)
To build it from scratch or understand how the underlying algorithm works, check the shoelace formula:
# e.g. corners = [(2.0, 1.0), (4.0, 5.0), (7.0, 8.0)]def Area(corners): n = len(corners) # of corners area = 0.0 for i in range(n): j = (i + 1) % n area += corners[i][0] * corners[j][1] area -= corners[j][0] * corners[i][1] area = abs(area) / 2.0 return area
Since this works for simple polygons:
If you have a polygon with holes : Calculate the area of the outer ring and subtrack the areas of the inner rings
If you have self-intersecting rings : You have to decompose them into simple sectors
By analysis of Mahdi's answer, I concluded that the majority of time was spent doing np.roll()
. By removing the need of the roll, and still using numpy, I got the execution time down to 4-5µs per loop compared to Mahdi's 41µs (for comparison Mahdi's function took an average of 37µs on my machine).
def polygon_area(x,y): correction = x[-1] * y[0] - y[-1]* x[0] main_area = np.dot(x[:-1], y[1:]) - np.dot(y[:-1], x[1:]) return 0.5*np.abs(main_area + correction)
By calculating the correctional term, and then slicing the arrays, there is no need to roll or create a new array.
Benchmarks:
10000 iterationsPolyArea(x,y): 37.075µs per looppolygon_area(x,y): 4.665µs per loop
Timing was done using the time
module and time.clock()