最近需要实现一个计算非凸多边形面积的功能,需要输入是顺次排序的多边形顶点坐标,假设输入的多边形顶点是V={v0, v1, v2, …, vn-1},则多边形的边为E={<v0, v1>, <v1, v2>, <v2, v3>,...,<vn-2, vn-1>, <vn-1, v0>}。要求输出该多边形所围成的面积。
对于凸多边形来讲,比较容易操作,只需要选择一个点作为公共顶点,然后将这个选定的公共顶点与所有其他顶点相连,即可将这个凸多边形分割为n-2个三角形,求这些三角形的面积之和即可。这里n为多边形的顶点个数,这些三角形都有这个选定的点作为一个公共顶点。
而对于非凸多边形显然无法直接这样操作。这一篇文章给出了一种求解任意多边形当然包括非凸多边形面积的算法。他的做法是给每个三角形的面积值一个符号,其基本思想仍然是构造一个类似于凸多边形的情况(并不一定是凸多边形,只是类似),然后减去多余的面积。这个思路是不错的,但是他在接下来具体实现的时候做法并不高明。具体在计算每个三角形的面积的时候,使用的方法不太好,不仅时间复杂度高,并且当两条线比较靠近的时候,使用开方操作,可能会导致溢出而无法开根号。因此我在这里对原作者的实现方法进行了一些改进。
三角形带符号的面积,容易让人想起向量的叉乘,之前在Graham-Scan算法计算凸包的Python代码实现一文中也使用了向量的叉乘。向量的叉乘是带符号的,并且其绝对值为两个向量共起点时所构成的平行四边形的面积。而平行四边形的对角线将整个平行四边形分割成了两个全等的三角形,因而每个三角形的面积为整个平行四边形的一半。
具体的Python实现及测试代码如下:
import numpy as np
"""
计算多边形的面积,包括非凸多边形的情况,
要求输入的顶点是按照逆时针顺序顺次排列的。
所要求解面积的多边形即为这些顶点一个一个地相连而成
E = {<v0, v1>, <v1, v2>, <v2, v3>,...,<vn-2, vn-1>, <vn-1, v0>}
@author: sdu_brz
@date: 2019/02/18
"""
def coss_multi(v1, v2):
"""
计算两个向量的叉乘
:param v1:
:param v2:
:return:
"""
return v1[0]*v2[1] - v1[1]*v2[0]
def polygon_area(polygon):
"""
计算多边形的面积,支持非凸情况
:param polygon: 多边形顶点,已经进行顺次逆时针排序
:return: 该多边形的面积
"""
n = len(polygon)
if n < 3:
return 0
vectors = np.zeros((n, 2))
for i in range(0, n):
vectors[i, :] = polygon[i, :] - polygon[0, :]
area = 0
for i in range(1, n):
area = area + coss_multi(vectors[i-1, :], vectors[i, :]) / 2
return area
if __name__ == "__main__":
"""测试"""
polygon1 = np.array([[0, 0],
[1, 0],
[1, 1],
[0, 1]])
print(polygon_area(polygon1))
polygon2 = np.array([[0, 0],
[5, 0],
[5, 4],
[4, 4],
[4, 1],
[1, 1],
[1, 7],
[0, 7]])
print(polygon_area(polygon2))