一、PNPoly 算法
在 GIS(地理信息管理系统)中,判断一个坐标是否在多边形内部是个经常要遇到的问题。乍听起来还挺复杂。根据 W. Randolph Franklin 提出的 PNPoly 算法,只需区区几行代码就解决了这个问题。
假设多边形的坐标存放在一个数组里,首先我们需要取得该数组在横坐标和纵坐标的最大值和最小值,根据这四个点算出一个四边型,首先判断目标坐标点是否在这个四边型之内,如果在这个四边型之外,那可以跳过后面较为复杂的计算,直接返回false。
比如:在多边形的顶点中分别找出 X 坐标和 Y 坐标的最小/最大值。比如,你有点(9,1),(4,3),(2,7),(8,2),(3,6)所围成的多边形。那么 Xmin 为2,Xmax 为9,Ymin 为1,Ymax 为7。现在我们知道你的多边形中没有一个点的 X 坐标比2小或者比9大,也没有一个点的 Y 坐标比1小或者比7大。这样你就可以快速排除很多不在多边形中的点。
if (p.x < minX || p.x > maxX || p.y < minY || p.y > maxY) {
// 这个测试都过不了。。。直接返回false;
return false;
}
接下来是核心算法部分:
//nvert:多边形中的顶点数
//vertx、verty:包含多边形顶点的x和y坐标的数组。
//testx、testy:需要判断的点位置
int pnpoly (int nvert, float *vertx, float *verty, float testx, float testy) {
int i, j, c = 0;
for (i = 0, j = nvert-1; i < nvert; j = i++) {
if ( ( (verty[i]>testy) != (verty[j]>testy) ) && (testx < (vertx[j]-vertx[i]) * (testy-verty[i]) / (verty[j]-verty[i]) + vertx[i]) )
c = !c;
}
return c;
}
首先,参数nvert 代表多边形有几个点。浮点数testx, testy代表待测试点的横坐标和纵坐标,*vertx,*verty分别指向储存多边形横纵坐标数组的首地址。
我们注意到,每次计算都涉及到相邻的两个点和待测试点,然后考虑两个问题:
-
被测试点的纵坐标testy是否在本次循环所测试的两个相邻点纵坐标范围之内?即verty[i] <testy < verty[j] 或者verty[j] <testy < verty[i]
-
待测点test是否在i,j两点之间的连线之下?看不懂后半短if statement的朋友请自行在纸上写下i,j两点间的斜率公式,要用到一点初中解析几何和不等式的知识范畴。
然后每次这两个条件同时满足的时候我们把返回的布尔量取反。
这个表达式的意思是说,随便画个多边形,随便定一个点,然后通过这个点水平划一条线,先数数看这条横线和多边形的边相交几次,(或者说先排除那些不相交的边,第一个判断条件),然后再数这条横线穿越多边形的次数是否为奇数,如果是奇数,那么该点在多边形内,如果是偶数,则在多边形外。详细的数学证明这里就不做了,不过读者可以自行画多边形进行验证。
二、交点奇偶数判定
从图中可得结论:从目标点出发引一条射线,看这条射线和多边形所有边的交点数目。如果有奇数个交点,则说明在内部,如果有偶数个交点,则说明在外部。利用此结论近些代码编写:
// 功能:判断点是否在多边形内
// 方法:求解通过该点的水平线与多边形各边的交点
// 结论:单边交点为奇数,成立!
// 参数:
// POINT p 指定的某个点
// LPPOINT ptPolygon 多边形的各个顶点坐标(首末点可以不一致)
// int nCount 多边形定点的个数
public static bool PointInPolygon(Vector2 p, Vector2[] ptPolygon, int nCount)
{
int nCross = 0;
for (int i = 0; i < nCount; i++)
{
Vector2 p1 = ptPolygon[i];//当前节点
Vector2 p2 = ptPolygon[(i + 1) % nCount];//下一个节点
// 求解 y=p.y 与 p1p2 的交点
if (p1.y == p2.y) // p1p2 与 y=p0.y平行
continue;
if (p.y < Mathf.Min(p1.y, p2.y)) // 交点在p1p2延长线上
continue;
if (p.y >= Mathf.Max(p1.y, p2.y)) // 交点在p1p2延长线上
continue;
// 从P发射一条水平射线 求交点的 X 坐标 ------原理: ((p2.y-p1.y)/(p2.x-p1.x))=((y-p1.y)/(x-p1.x))
//直线k值相等 交点y=p.y
double x = (double)(p.y - p1.y) * (double)(p2.x - p1.x) / (double)(p2.y - p1.y) + p1.x;
if (x > p.x)
nCross++; // 只统计单边交点
}
// 单边交点为偶数,点在多边形之外 ---
return (nCross % 2 == 1);
}
}