判断点是否在多边形内的Python实现及小应用(射线法)

本文介绍了一种判断点是否在多边形内的方法——射线法,并提供了Python实现。射线法通过计算从点出发的射线与多边形边的交点个数,奇数则点在多边形内,偶数则在多边形外。适用于GIS应用,如点选功能、数据筛选等。

摘要生成于 C知道 ,由 DeepSeek-R1 满血版支持, 前往体验 >

判断点是否在多边形内的Python实现及小应用(射线法)

转 https://www.jianshu.com/p/ba03c600a557

判断一个点是否在多边形内是处理空间数据时经常面对的需求,例如GIS中的点选功能、根据多边形边界筛选出位于多边形内的点、求交集、筛选不在多边形内的点等等。
判断一个点是否在多边形内有几种不同的思路,相应的方法(感觉还谈不上算法)有:

  • 射线法:从判断点向某个统一方向作射线,依交点个数的奇偶判断;
  • 转角法:按照多边形顶点逆时针顺序,根据顶点和判断点连线的方向正负(设定角度逆时针为正)求和判断;
  • 夹角和法:求判断点与所有边的夹角和,等于360度则在多边形内部。
  • 面积和法:求判断点与多边形边组成的三角形面积和,等于多边形面积则点在多边形内部。

面积和法涉及多个面积的计算,比较复杂,夹角和法以及转角法用到角度计算,会涉及反三角函数,计算开销比较大,而射线法主要涉及循环多边形的每条边进行求交运算,但大部分边可以通过简单坐标比对直接排除,因此这是比较好的方法。

射线法的实现

射线法就是以判断点开始,向右(或向左)的水平方向作一射线,计算该射线与多边形每条边的交点个数,如果交点个数为奇数,则点位于多边形内,偶数则在多边形外。该算法对于复合多边形也能正确判断。

复合多边形的情况


射线法的关键是正确计算射线与每条边是否相交。并且规定线段与射线重叠或者射线经过线段下端点属于不相交。首先排除掉不相交的情况,下图的情况都是需要排除掉的:

求交之前可排除的情况


排除掉这些情况的函数如下:

 

def isRayIntersectsSegment(poi,s_poi,e_poi): #[x,y] [lng,lat]
    #输入:判断点,边起点,边终点,都是[lng,lat]格式数组
    if s_poi[1]==e_poi[1]: #排除与射线平行、重合,线段首尾端点重合的情况
        return False
    if s_poi[1]>poi[1] and e_poi[1]>poi[1]: #线段在射线上边
        return False
    if s_poi[1]<poi[1] and e_poi[1]<poi[1]: #线段在射线下边
        return False
    if s_poi[1]==poi[1] and e_poi[1]>poi[1]: #交点为下端点,对应spoint
        return False
    if e_poi[1]==poi[1] and s_poi[1]>poi[1]: #交点为下端点,对应epoint
        return False
    if s_poi[0]<poi[0] and e_poi[1]<poi[1]: #线段在射线左边
        return False

    xseg=e_poi[0]-(e_poi[0]-s_poi[0])*(e_poi[1]-poi[1])/(e_poi[1]-s_poi[1]) #求交
    if xseg<poi[0]: #交点在射线起点的左侧
        return False
    return True  #排除上述情况之后

排除掉上述情况真正需要求交点来判断的情况只有两种:

 

需要求交的两种情况

函数isRayIntersectsSegment()里求交的部分就是利用两个三角形的比例关系求出交点在起点的左边还是右边;用图去理解如下:

 

求交的具体过程

最后判断的代码如下:

def isPoiWithinPoly(poi,poly):
    #输入:点,多边形三维数组
    #poly=[[[x1,y1],[x2,y2],……,[xn,yn],[x1,y1]],[[w1,t1],……[wk,tk]]] 三维数组

    #可以先判断点是否在外包矩形内 
    #if not isPoiWithinBox(poi,mbr=[[0,0],[180,90]]): return False
    #但算最小外包矩形本身需要循环边,会造成开销,本处略去
    sinsc=0 #交点个数
    for epoly in poly: #循环每条边的曲线->each polygon 是二维数组[[x1,y1],…[xn,yn]]
        for i in range(len(epoly)-1): #[0,len-1]
            s_poi=epoly[i]
            e_poi=epoly[i+1]
            if isRayIntersectsSegment(poi,s_poi,e_poi):
                sinsc+=1 #有交点就加1

    return True if sinsc%2==1 else  False

我们取一个比较复杂的多边形进行测试,多边形和一些点如图:

 

测试用的有孔洞多边形

 

用isPoiWithinPoly()的测试结果如下:

 

测试结果

点在多边形内的应用

上面第一段已经描述了一些应用场景,下面给出一个应用的例子:有一堆点数据存在csv文件里,如何检索位于某个城市的点出来,检索出来之后的分析(例如加标签、改属性、做统计还是其他)这里不讨论,检索的结果统一写到新文件里。点输入的格式如下:

id,name,wgslng,wgslat,score,adds
1,沃美,116.3309,40.0706,4.3,昌平回龙观同成街华联购物中心4楼
2,星美国际,116.446,39.916,5,金汇路8号世界城E座
3,……

城市边界为geojson格式,就是加了一些限定条件的json格式数据,如果需要详细了解geojson格式,可以参考本人之前的文章:GEOJSON标准格式学习。形如:

{
  "type": "FeatureCollection",
  "features": [{
      "type": "Feature",
      "properties": {},
      "geometry": {
        "type": "Polygon",
        "coordinates":
         [
            [
              [108.71658325195312,34.231106222010531],
              [108.96240234375,34.168635904722734],
              [109.00222778320313,34.354774165387568],
              [108.80172729492186,34.35023911062779],
              [108.71658325195312,34.231106222010531]
            ]
          ]
        }
      }
  ]
}

下面的代码只考虑了Polygon的情况,对于MultiPolygon也是比较容易改的,要改为处理kml保存的边界数据也不难改。文中代码同步于本人GitHub

import json
import csv
def pointInPolygon():
    gfile = './beijing_poly_wgs84.geojson' #utf-8编码
    cin_path = './poi_cinema_wgs84.csv'
    out_path = './beijing_poi_cinema_wgs84.csv' #输出文件

    pindex = [2, 3]  # wgslng,wgslat 在的位置

    with open(out_path, 'w', newline='') as cout_file:
        fin = open(cin_path, 'r', encoding='gbk') #出现编码错误就改编码 utf-8
        gfn = open(gfile, 'r', encoding='utf-8')
        gjson = json.load(gfn)
        polygon = gjson["features"][0]["geometry"]['coordinates'] #提取多边形,如果是4维数组需要相应的处理
        filewriter = csv.writer(cout_file, delimiter=',')
        w = 0
        for line in csv.reader(fin, delimiter=','):
            if w == 0: #写入表头 id,name,… 如果没有就去掉if语句
                filewriter.writerow(line)
                w = 1
                continue
            point = [float(line[pindex[0]]), float(line[pindex][1])]
            if isPoiWithinPoly(point, polygon): #在多边形内,写入新表
                filewriter.writerow(line)
            else:
                continue
        fin.close()
        gfn.close()
    print('done')
<think>好的,我现在需要帮助用户找到用Python实现判断是否多边形内的方。用户提到了射线,但可能对具体实现不太清楚。首先,我得回忆一下射线的原理。射线的基本思路是从待测向一个方向发射射线,统计射线多边形边的交数量。如果是奇数,在内部;偶数则在外部。这个方的关键在于如何处理边界情况,比如射线穿过顶或者与边重合。 接下来,我需要查找相关的Python代码实现。根据用户提供的引用,特别是引用[1]、[2]、[3]、[4],这些资料应该提供了不同的实现。比如引用[3]中的PNPoly算,这可能是一个常用的方。我需要确保代码正确处理各种情况,比如水平射线与边的交计算,以及在边上的情况。 然后,我需要考虑代码的结构。首先定义函数,参数包括的坐标和多边形的顶列表。遍历多边形的每条边,检查射线与边的交。这里需要注意的是,射线的方向通常选择水平向右,这样计算交时比较方便。对于每条边,需要判断是否在边的两侧,然后计算交的x坐标是否在右侧。 还要处理特殊情况,比如多边形的顶上,或者射线恰好穿过顶。这时候需要根据相邻边的位置来决定是否计数。例如,如果射线穿过顶,而该顶多边形的一个凸顶,可能需要调整计数策略,避免重复或遗漏交。 另外,代码的效率也很重要。多边形可能有很多顶,所以需要循环遍历所有边,但时间复杂度是O(n),这在大多数情况下是可接受的。用户可能希望代码简洁高效,所以需要优化条件判断,避免不必要的计算。 测试代码的部分,用户可能需要示例来验证正确性。比如用一个正方形多边形,测试内部、外部和边界上的。确保代码能正确处理这些情况,并返回正确的结果。 最后,总结射线的优缺,比如对凸多边形和凹多边形都适用,但需要注意精度问题,比如浮数运算可能导致的误差。同时,可以提到其他方如转角或面积,但射线较为常用且实现简单。 现在将这些思路整理成步骤清晰的代码和解释,确保符合用户的要求,使用正确的中文表达,并引用相关的参考资料。</think>### Python实现射线判断是否多边形内 #### 算原理 射线(Ray Casting Algorithm)通过从待测$P$向水平右侧发射射线,统计与多边形边的交数量: - **奇数个交**:多边形内 - **偶数个交**:多边形外[^1] $$ \text{交数量} \mod 2 = \begin{cases} 1 & \text{内部} \\ 0 & \text{外部} \end{cases} $$ #### 实现步骤 1. **输入参数**:待测$(x,y)$,多边形列表$polygon$ 2. **遍历各边**:依次取相邻顶构成边$V_iV_{i+1}$ 3. **边界判断**:若恰好在边上则直接返回True 4. **射线交叉检测**:判断射线是否与边相交 5. **奇偶统计**:最终根据交数量奇偶性返回结果 #### Python代码实现 ```python def is_point_in_polygon(point, polygon): x, y = point n = len(polygon) inside = False for i in range(n): p1x, p1y = polygon[i] p2x, p2y = polygon[(i+1) % n] # 检查是否多边形边上 if ((p1y == y == p2y) and (min(p1x, p2x) <= x <= max(p1x, p2x))) or \ ((p1x == x == p2x) and (min(p1y, p2y) <= y <= max(p1y, p2y))): return True # 射线与边相交检测 if (p1y > y) != (p2y > y): # 边跨越射线高度 xinters = (y - p1y) * (p2x - p1x) / (p2y - p1y) + p1x # 求交X坐标 if x <= xinters: # 交射线右侧 inside = not inside return inside ``` #### 测试案例 ```python # 定义一个正方形多边形 polygon = [(0,0), (0,10), (10,10), (10,0)] # 测试 print(is_point_in_polygon((5,5), polygon)) # True(内部) print(is_point_in_polygon((15,5), polygon)) # False(外部) print(is_point_in_polygon((0,5), polygon)) # True(边界) ``` #### 关键说明 1. **射线方向**:默认使用水平向右射线,避免浮精度问题 2. **边界处理**:在边上直接返回True(根据需求可调整) 3. **顶处理**:通过`(i+1) % n`自动连接首尾顶 4. **浮运算**:使用精确除避免整数运算错误[^3] #### 应用场景 - GIS空间分析(如地图选) - 计算机图形学中的裁剪算 - 游戏开发中的碰撞检测 - CAD/CAM系统中的区域判定
评论
添加红包

请填写红包祝福语或标题

红包个数最小为10个

红包金额最低5元

当前余额3.43前往充值 >
需支付:10.00
成就一亿技术人!
领取后你会自动成为博主和红包主的粉丝 规则
hope_wisdom
发出的红包
实付
使用余额支付
点击重新获取
扫码支付
钱包余额 0

抵扣说明:

1.余额是钱包充值的虚拟货币,按照1:1的比例进行支付金额的抵扣。
2.余额无法直接购买下载,可以购买VIP、付费专栏及课程。

余额充值