判断一个坐标点是否在不规则多边形内部的算法

墨蓝 2023-09-30 13:20 26阅读 0赞

参考：https://wrf.ecse.rpi.edu//Research/Short\_Notes/pnpoly.html

在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;
    }

<table> 
 <thead> 
  <tr> 
   <th>Argument</th> 
   <th>Meaning</th> 
  </tr> 
 </thead> 
 <tbody> 
  <tr> 
   <td>nvert</td> 
   <td>Number of vertices in the polygon. Whether to repeat the first vertex at the end is discussed below.</td> 
  </tr> 
  <tr> 
   <td>vertx, verty</td> 
   <td>Arrays containing the x- and y-coordinates of the polygon’s vertices.</td> 
  </tr> 
  <tr> 
   <td>testx, testy</td> 
   <td>X- and y-coordinate of the test point.</td> 
  </tr> 
 </tbody> 
</table>

额，代码就这么简单，但到底啥意思呢：

首先，参数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两点间的斜率公式，要用到一点初中解析几何和不等式的知识范畴，对广大码农来说小菜一碟。

然后每次这两个条件同时满足的时候我们把返回的布尔量取反。

可这到底是啥意思啊？

这个表达式的意思是说，随便画个多边形，随便定一个点，然后通过这个点水平划一条线，先数数看这条横线和多边形的边相交几次，（或者说先排除那些不相交的边，第一个判断条件），然后再数这条横线穿越多边形的次数是否为奇数，如果是奇数，那么该点在多边形内，如果是偶数，则在多边形外。详细的数学证明这里就不做了，不过读者可以自行画多边形进行验证。

判断一个点是否在多边形内部 - 射线法思路  
比如说，我就随便涂了一个多边形和一个点，现在我要给出一种通用的方法来判断这个点是不是在多边形内部（别告诉我用肉眼观察……）。

![这里写图片描述][7d9e2c7493e64f0ca913283ba6a559af.png]

首先想到的一个解法是从这个点做一条射线，计算它跟多边形边界的交点个数，如果交点个数为奇数，那么点在多边形内部，否则点在多边形外。

![在这里插入图片描述][8fb9367127b84b3c804b8d44ef658897.png]

这个结论很简单，那它是怎么来的？下面就简单讲解一下。

首先，对于平面内任意闭合曲线，我们都可以直观地认为，曲线把平面分割成了内、外两部分，其中“内”就是我们所谓的多边形区域。

![在这里插入图片描述][4acdf864e997426fab501a2bbe15cb3f.png]

基于这一认识，对于平面内任意一条直线，我们可以得出下面这些结论：

直线穿越多边形边界时，有且只有两种情况：进入多边形或穿出多边形。  
在不考虑非欧空间的情况下，直线不可能从内部再次进入多边形，或从外部再次穿出多边形，即连续两次穿越边界的情况必然成对。  
直线可以无限延伸，而闭合曲线包围的区域是有限的，因此最后一次穿越多边形边界，一定是穿出多边形，到达外部。  
![在这里插入图片描述][ca041e2826a64bce83c04fb5fb7417fe.png]

现在回到我们最初的题目。假如我们从一个给定的点做射线，还可以得出下面两条结论：

如果点在多边形内部，射线第一次穿越边界一定是穿出多边形。  
如果点在多边形外部，射线第一次穿越边界一定是进入多边形。  
![在这里插入图片描述][b2ef2aa132c54404966118e8d388d973.png]

把上面这些结论综合起来，我们可以归纳出：

当射线穿越多边形边界的次数为偶数时，所有第偶数次（包括最后一次）穿越都是穿出，因此所有第奇数次（包括第一次）穿越为穿入，由此可推断点在多边形外部。  
![在这里插入图片描述][3f2424231d06448a8c44eb15866de952.png]

当射线穿越多边形边界的次数为奇数时，所有第奇数次（包括第一次和最后一次）穿越都是穿出，由此可推断点在多边形内部。

![在这里插入图片描述][52cf93ea52dc4741be82c289b4578701.png]

到这里，我们已经了解了这个解法的思路，大家可以试着自己写一个实现出来。

不知道大家思考得怎么样，有没有遇到一些不好处理的特殊情况。今天就来讲讲射线法在实际应用中的一些问题和解决方案。

1点在多边形的边上

前面我们讲到，射线法的主要思路就是计算射线穿越多边形边界的次数。那么对于点在多边形的边上这种特殊情况，射线出发的这一次，是否应该算作穿越呢？

![在这里插入图片描述][aacd723551944d359d35c3fa89bbf133.png]

看了上面的图就会发现，不管算不算穿越，都会陷入两难的境地——同样落在多边形边上的点，可能会得到相反的结果。这显然是不正确的，因此对这种特殊情况需要特殊处理。

2点和多边形的顶点重合

![在这里插入图片描述][93c429915c314d6f9c1c1358e937d496.png]

这其实是第一种情况的一个特例。

3射线经过多边形顶点

射线刚好经过多边形顶点的时候，应该算一次还是两次穿越？这种情况比前两种复杂，也是实现中的难点，后面会讲解它的解决方案。

![在这里插入图片描述][bd2acf79c4f34b47b6ed80ae3f0e48f7.png]

4射线刚好经过多边形的一条边

这是上一种情况的特例，也就是说，射线连续经过了多边形的两个相邻顶点。

![在这里插入图片描述][c8278431fd4d4ad3bde56a790d53811d.png]

解决方案：

1判断点是否在线上的方法有很多，比较简单直接的就是计算点与两个多边形顶点的连线斜率是否相等，中学数学都学过。

2点和多边形顶点重合的情况更简单，直接比较点的坐标就行了。

3顶点穿越看似棘手，其实我们换一个角度，思路会大不相同。先来回答一个问题，射线穿越一条线段需要什么前提条件？没错，就是线段两个端点分别在射线两侧。只要想通这一点，顶点穿越就迎刃而解了。这样一来，我们只需要规定被射线穿越的点都算作其中一侧。

![在这里插入图片描述][d8417d39470a45efb868fbe251ed1221.png]

如上图，假如我们规定射线经过的点都属于射线以上的一侧，显然点D和发生顶点穿越的点C都位于射线Y的同一侧，所以射线Y其实并没有穿越CD这条边。而点C和点B则分别位于射线Y的两侧，所以射线Y和BC发生了穿越，由此我们可以断定点Y在多边形内。同理，射线X分别与AD和CD都发生了穿越，因此点X在多边形外，而射线Z没有和多边形发生穿越，点Z位于多边形外。

解决了第三点，这一点就毫无难度了。根据上面的假设，射线连续经过的两个顶点显然都位于射线以上的一侧，因此这种情况看作没有发生穿越就可以了。由于第三点的解决方案实际上已经覆盖到这种特例，因此不需要再做特别的处理。

[7d9e2c7493e64f0ca913283ba6a559af.png]: https://img-blog.csdnimg.cn/7d9e2c7493e64f0ca913283ba6a559af.png
[8fb9367127b84b3c804b8d44ef658897.png]: https://img-blog.csdnimg.cn/8fb9367127b84b3c804b8d44ef658897.png
[4acdf864e997426fab501a2bbe15cb3f.png]: https://img-blog.csdnimg.cn/4acdf864e997426fab501a2bbe15cb3f.png
[ca041e2826a64bce83c04fb5fb7417fe.png]: https://img-blog.csdnimg.cn/ca041e2826a64bce83c04fb5fb7417fe.png
[b2ef2aa132c54404966118e8d388d973.png]: https://img-blog.csdnimg.cn/b2ef2aa132c54404966118e8d388d973.png
[3f2424231d06448a8c44eb15866de952.png]: https://img-blog.csdnimg.cn/3f2424231d06448a8c44eb15866de952.png
[52cf93ea52dc4741be82c289b4578701.png]: https://img-blog.csdnimg.cn/52cf93ea52dc4741be82c289b4578701.png
[aacd723551944d359d35c3fa89bbf133.png]: https://img-blog.csdnimg.cn/aacd723551944d359d35c3fa89bbf133.png
[93c429915c314d6f9c1c1358e937d496.png]: https://img-blog.csdnimg.cn/93c429915c314d6f9c1c1358e937d496.png
[bd2acf79c4f34b47b6ed80ae3f0e48f7.png]: https://img-blog.csdnimg.cn/bd2acf79c4f34b47b6ed80ae3f0e48f7.png
[c8278431fd4d4ad3bde56a790d53811d.png]: https://img-blog.csdnimg.cn/c8278431fd4d4ad3bde56a790d53811d.png
[d8417d39470a45efb868fbe251ed1221.png]: https://img-blog.csdnimg.cn/d8417d39470a45efb868fbe251ed1221.png