在GIS(地理信息管理系统)中,判断一个坐标是否在多边形内部是个经常要遇到的问题。乍听起来还挺复杂。根据W. Randolph Franklin 提出的PNPoly算法,只需区区几行代码就解决了这个问题。 假设多边形的坐标存放在一个数组里,首先我们需要取得该数组在横坐标和纵坐标的最大值和最小值,根据这四个点算出一个四边型,首先判断目标坐标点是否在这个四边型之内,如果在这个四边型之外,那可以跳过后面较为复杂的计算,直接返回false。 if (p.x < minX || p.x > maxX || p.y < minY || p.y > maxY) { // 这个测试都过不了。。。直接返回false; } 接下来是核心算法部分: 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分别指向储存多边形横纵坐标数组的首地址。 我们注意到,每次计算都涉及到相邻的两个点和待测试点,然后考虑两个问题: 1. 被测试点的纵坐标testy是否在本次循环所测试的两个相邻点纵坐标范围之内?即 verty[i] <testy < verty[j] 或者 verty[j] <testy < verty[i] 2. 待测点test是否在i,j两点之间的连线之下?看不懂后半短if statement的朋友请自行在纸上写下i,j两点间的斜率公式,要用到一点初中解析几何和不等式的知识范畴,对广大码农来说小菜一碟。 然后每次这两个条件同时满足的时候我们把返回的布尔量取反。 可这到底是啥意思啊? 这个表达式的意思是说,随便画个多边形,随便定一个点,然后通过这个点水平划一条射线,先数数看这条
射线和多边形的边相交几次,(或者说先排除那些不相交的边,第一个判断条件),然后再数这条射线穿越多边形的次数是否为奇数,如果是奇数,那么该点在多边形内,如果是偶数,则在多边形外。详细的数学证明这里就不做了,不过读者可以自行画多边形进行验证。
以上转载百度知道;以下转载阿凡卢,转载地址:http://www.cnblogs.com/luxiaoxun/p/3722358.html
判断点是否在多边形内部
如何判断一个点是否在多边形内部?
(1)面积和判别法:判断目标点与多边形的每条边组成的三角形面积和是否等于该多边形,相等则在多边形内部。
(2)夹角和判别法:判断目标点与所有边的夹角和是否为360度,为360度则在多边形内部。
(3)引射线法:从目标点出发引一条射线,看这条射线和多边形所有边的交点数目。如果有奇数个交点,则说明在内部,如果有偶数个交点,则说明在外部。
具体做法:将测试点的Y坐标与多边形的每一个点进行比较,会得到一个测试点所在的行与多边形边的交点的列表。在下图的这个例子中有8条边与测试点所在的行相交,而有6条边没有相交。如果测试点的两边点的个数都是奇数个则该测试点在多边形内,否则在多边形外。在这个例子中测试点的左边有5个交点,右边有三个交点,它们都是奇数,所以点在多边形内。
算法图解:
关于这个算法的具体的更多图形例子:http://alienryderflex.com/polygon/
参考代码:
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; }
来自一个polygon的内部实现:
public bool IsInside(PointLatLng p) { int count = Points.Count; if(count < 3) { return false; } bool result = false; for(int i = 0, j = count - 1; i < count; i++) { var p1 = Points[i]; var p2 = Points[j]; if(p1.Lat < p.Lat && p2.Lat >= p.Lat || p2.Lat < p.Lat && p1.Lat >= p.Lat) { if(p1.Lng + (p.Lat - p1.Lat) / (p2.Lat - p1.Lat) * (p2.Lng - p1.Lng) < p.Lng) { result = !result; } } j = i; } return result; }
特殊情况:要检测的点在多变形的一条边上,射线法判断的结果是不确定的,需要特殊处理(If the test point is on the border of the polygon, this algorithm will deliver unpredictable results)。
计算一个多边形的面积(area of a polygon):
private static double SignedPolygonArea(List<PointLatLng> points) { // Add the first point to the end. int pointsCount = points.Count; PointLatLng[] pts = new PointLatLng[pointsCount + 1]; points.CopyTo(pts, 0); pts[pointsCount] = points[0]; for (int i = 0; i < pointsCount + 1; ++i) { pts[i].Lat = pts[i].Lat * (System.Math.PI * 6378137 / 180); pts[i].Lng = pts[i].Lng * (System.Math.PI * 6378137 / 180); } // Get the areas. double area = 0; for (int i = 0; i < pointsCount; i++) { area += (pts[i + 1].Lat - pts[i].Lat) * (pts[i + 1].Lng + pts[i].Lng) / 2; } // Return the result. return area; } /// <summary> /// Get the area of a polygon /// </summary> /// <param name="points"></param> /// <returns></returns> public static double GetPolygonArea(List<PointLatLng> points) { // Return the absolute value of the signed area. // The signed area is negative if the polygon is oriented clockwise. return Math.Abs(SignedPolygonArea(points)); }
参考资料:
http://alienryderflex.com/polygon/
http://en.wikipedia.org/wiki/Point_in_polygon
http://www.codeproject.com/Tips/84226/Is-a-Point-inside-a-Polygon
作者:阿凡卢
本文版权归作者和博客园共有,欢迎转载,但未经作者同意必须保留此段声明,且在文章页面明显位置给出原文连接,否则保留追究法律责任的权利。