温暖的卡布奇诺 · 用python的pandas读取excel文 ...· 1 月前 · |
阳刚的小马驹 · logstash 这个怎么理解mutate ...· 6 月前 · |
孤独的鞭炮 · DeviceInformation.Crea ...· 1 年前 · |
没人理的面包 · 在 ML.NET ...· 1 年前 · |
爱听歌的铁链 · HTTP Referer 教程 - ...· 1 年前 · |
假设二维空间中的一系列不自相交的点,确定结果多边形的面积的有效方法是什么?
顺便说一句,这不是作业,我也不是在找代码。我正在寻找一种描述,我可以用来实现我自己的方法。我有从点列表中提取一系列三角形的想法,但我知道有一堆关于凸多边形和凹多边形的边缘情况,我可能无法捕捉到。
我的倾向是简单地开始切开三角形。我看不出其他任何东西可以避免可怕的毛发。
取组成多边形的三个连续点。确保角度小于180。现在你有了一个新的三角形,计算起来应该没有问题,从多边形的点列表中删除中点。重复操作,直到只剩下三个点。
一种方法是使用 decompose the polygon into triangles ,计算三角形的面积,并将总和作为多边形的面积。
或者做一个轮廓积分。斯托克斯定理允许你将面积积分表示为轮廓积分。一点高斯求积,鲍勃就是你的叔叔。
没有任何其他约束的一组点不一定是唯一定义多边形的。
因此,首先您必须决定从这些点构建哪个多边形-也许是凸包? http://en.wikipedia.org/wiki/Convex_hull
然后三角测量并计算面积。 http://www.mathopenref.com/polygonirregulararea.html
这是 the standard method ,AFAIK。基本上对每个顶点周围的叉积求和。比三角测量简单得多。
Python代码,给定一个表示为(x,y)顶点坐标列表的多边形,隐式地从最后一个顶点绕到第一个顶点:
def area(p):
return 0.5 * abs(sum(x0*y1 - x1*y0
for ((x0, y0), (x1, y1)) in segments(p)))
def segments(p):
return zip(p, p[1:] + [p[0]])
David Lehavi评论:值得一提的是该算法的工作原理:它是函数−y和x的 Green's theorem 应用程序;与 planimeter 的工作方式完全相同。更确切地说:
上面的公式=
integral_over_perimeter(-y dx + x dy) =
integral_over_area((-(-dy)/dy+dx/dx) dy dx) =
2 Area
为了扩展三角剖分和求和三角形区域,如果你碰巧有一个凸多边形,或者你碰巧选择了一个不会生成与多边形相交的其他点的线的点,那么这些区域就会起作用。
对于一般不相交的多边形,需要将矢量(参照点,a点)、(参照点,b点)的叉积求和,其中a和b彼此“相邻”。
假设您有一个按顺序定义多边形的点列表(顺序是点I和i+1形成一条多边形线):
Sum(叉积((点0,点i),(点0,点i+ 1)) for i= 1 to n-1。
取这个叉积的大小,你就得到了表面积。
这将处理凹多边形,而不必担心选择一个好的参考点;生成不在多边形内部的三角形的任何三个点都将有一个叉积,该叉积指向多边形内部的任何三角形的相反方向,因此面积可以正确求和。
在笛卡尔空间中对梯形求和比对三角形求和更好:
area = 0;
for (i = 0; i < n; i++) {
i1 = (i + 1) % n;
area += (vertex[i].y + vertex[i1].y) * (vertex[i1].x - vertex[i].x) / 2.0;
}
叉积是经典的。
如果你有大量这样的计算要做,试试下面的优化版本,它需要的乘法次数减少了一半:
area = 0;
for( i = 0; i < N; i += 2 )
area += x[i+1]*(y[i+2]-y[i]) + y[i+1]*(x[i]-x[i+2]);
area /= 2;
为了清晰起见,我使用数组下标。使用指针更有效。尽管好的编译器会帮你做到这一点。
假设多边形是“闭合”的,这意味着您将第一个点复制为下标为N的点。它还假设多边形具有偶数个点。如果N不是偶数,则附加第一个点的另一个副本。
该算法是通过展开并结合经典叉积算法的两次连续迭代而获得的。
我不太确定这两种算法在数值精度方面的比较。我的印象是,上述算法比经典算法更好,因为乘法往往会恢复减法的精度损失。当约束使用浮点数时,就像使用GPU一样,这可能会产生很大的差异。
编辑: "Area of Triangles and Polygons 2D & 3D" 描述了一种更有效的方法
// "close" polygon
x[N] = x[0];
x[N+1] = x[1];
y[N] = y[0];
y[N+1] = y[1];
// compute area
area = 0;
for( size_t i = 1; i <= N; ++i )
area += x[i]*( y[i+1] - y[i-1] );
area /= 2;
Shoelace formula 的实现可以在Numpy中完成。假设这些顶点:
import numpy as np
x = np.arange(0,1,0.001)
y = np.sqrt(1-x**2)
我们可以定义以下函数来查找面积:
def PolyArea(x,y):
return 0.5*np.abs(np.dot(x,np.roll(y,1))-np.dot(y,np.roll(x,1)))
并获得结果:
print PolyArea(x,y)
# 0.26353377782163534
避免循环使此函数比
PolygonArea
快约50倍
%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
注意:我已经为另一个 question 写了这个答案,我只是在这里提到这个,以获得一个完整的解决方案列表。
C的方法是:
float areaForPoly(const int numVerts, const Point *verts)
Point v2;
float area = 0.0f;
for (int i = 0; i<numVerts; i++){
v2 = verts[(i + 1) % numVerts];
area += verts[i].x*v2.y - verts[i].y*v2.x;
return area / 2.0f;
}
Python代码
如下所述: http://www.wikihow.com/Calculate-the-Area-of-a-Polygon
和熊猫一起
import pandas as pd
df = pd.DataFrame({'x': [10, 20, 20, 30, 20, 10, 0], 'y': [-10, -10, -10, 0, 10, 30, 20]})
df = df.append(df.loc[0])
阳刚的小马驹 · logstash 这个怎么理解mutate { split => { "genre" => "|" } remove_field => ["path", "host","@timest 6 月前 |
孤独的鞭炮 · DeviceInformation.CreateWatcher Method (Windows.Devices.Enumeration) - Windows UWP applications | Microsoft Learn 1 年前 |
爱听歌的铁链 · HTTP Referer 教程 - 阮一峰的网络日志 1 年前 |