核心内容摘要
为什么92%的自研低代码平台卡在V2.0?Python内核必须攻克的5个硬核关卡:Schema演化、版本快照、跨租户隔离、插件热插拔、回滚一致性
在车联网、POI 分析、地理围栏、行政区划判断等业务中我们经常会遇到一个非常基础、但极其高频的问题给定一个经纬度点它属于哪个或哪些多边形比如一个 GPS 点属于哪个行政区一个车辆轨迹点是否落在某个物流园围栏内一个 POI 是否在规划的服务区域中这个问题在 GIS 领域被称为Point in PolygonPIP问题。
本文将结合一份完整可运行的 Python 代码从算法原理工程实现性能瓶颈加速手段包围盒、R 树几个层面系统性地讲清楚点在哪个多边形里到底该怎么算、怎么选、怎么快。
最原始的方法逐个多边形判断1️⃣ 思路最直观的方式是对每一个多边形判断点是否在这个多边形内部数学上常见的实现有射线法Ray Casting它计算从点P开始的射线穿过多边形边界的次数。
当交叉数是偶数时点在外面当交叉数是奇数时点在里面。
环绕数法Winding Number他计算多边形围着点P旋转地次数只有当圈数wn0时点在外面否则点在里面。
在工程中我们通常不会自己手写几何算法而是直接使用成熟库。
2️⃣ 代码实现defpoint_in_polygon(point_lon,point_lat,polygon_coords): 判断一个点是否在多边形内 参数: point_lon: 点的经度float point_lat: 点的纬度float polygon_coords: 多边形顶点列表格式为 [(lon1, lat
, (lon2, lat
, ...] 返回: True: 点在多边形内 False: 点不在多边形内 # 方式一基于射线法importmatplotlib.pathasmplPathimportnumpyasnp polygonmplPath.Path(np.array(polygon_coords))returnpolygon.contains_point(point)# 方式二基于射线法importcv2 cv
pointPolygonTest(np.array(polygon_coords),point,measureDistFalse)# 方式三基于环绕数法fromshapely.geometryimportPoint,Polygon pointPoint(point_lon,point_lat)polygonPolygon(polygon_coords)returnpolygon.contains(point)如果有多个多边形deffind_point_in_polygons(point_lon,point_lat,polygons_dict): 查找点所在的多边形不使用索引适合小数据集 参数: point_lon: 点的经度 point_lat: 点的纬度 polygons_dict: 多边形字典 返回: 包含该点的多边形名称列表 pointPoint(point_lon,point_lat)result[]forname,coordsinpolygons_dict.items():polygonPolygon(coords)ifpolygon.contains(point):result.append(name)returnresult3️⃣ 性能分析时间复杂度O(N)N 多边形数量如果多边形数量 100查询点数量不多完全没问题简单、清晰、好维护。
但一旦进入真实业务场景上千 / 上万围栏每秒上百 / 上千点性能会迅速成为瓶颈。
第一步加速包围盒Bounding Box过滤1️⃣ 核心思想绝大多数点连多边形的“外接矩形”都不在里面。
那我们就可以先判断点是否在多边形的包围盒内极快再对“可能命中”的少量多边形做精确判断2️⃣ 什么是包围盒对一个多边形(minx, miny, maxx, maxy)只要minx ≤ lon ≤ maxx miny ≤ lat ≤ maxy才有可能在多边形内。
3️⃣ 代码结构classPolygonIndexWithBounds: 使用包围盒快速过滤的空间索引 def__init__(self,polygons_dict):初始化包围盒索引self.polygons{}self.bounds_map{}forname,coordsinpolygons_dict.items():polygonPolygon(coords)self.polygons[name]polygon self.bounds_map[name]polygon.boundsdeffind_point_in_polygons(self,point_lon,point_lat):查找点所在的多边形pointPoint(point_lon,point_lat)result[]# 通过包围盒过滤candidates[]forname,boundsinself.bounds_map.items():minx,miny,maxx,maxyboundsifminxpoint_lonmaxxandminypoint_latmaxy:candidates.append(name)# 精确检查fornameincandidates:ifself.polygons[name].contains(point):result.append(name)returnresultdeffind_points_in_polygons_batch(self,points_list):批量查询results[]forlon,latinpoints_list:polygonsself.find_point_in_polygons(lon,lat)results.append({point:(lon,lat),polygons:polygons})returnresults4️⃣ 性能效果包围盒判断纯数值比较极快精确几何判断次数大幅减少在测试中通常可提速 10 ~ 100 倍5️⃣ 适用场景多边形数量100 ~ 1000实现成本低不依赖额外空间索引库性价比极高的工程方案
终极方案R 树Spatial Index1️⃣ 问题再升级当你遇到上万 / 数十万多边形实时或准实时查询包围盒 O(N) 扫描也开始吃力。
这时候就需要真正的空间索引。
2️⃣ 什么是 R 树R 树是一种为多维空间查询设计的树结构步骤对每个空间对象点、线、面构建最小外接矩形MBR。
如图D1|R
D2|R
D3|R10建树过程自底向上分组把若干个对象的 MBR→ 组合成一个“节点”对一个节点再算一个更大的 MBR→ 覆盖所有子节点节点数量超过阈值→ 进行 节点分裂split查询过程从根开始逐层过滤向下递归到达叶子节点特点多层包围盒类似 B 树的分层结构专为“空间相交 / 包含”优化查询复杂度接近O(log N)你可以把它理解为空间世界里的索引像数据库索引一样3️⃣ Shapely 中的 R 树Shapely
x 提供了STRtreefromshapely.strtreeimportSTRtree treeSTRtree(geometry_list)查询时candidatestree.query(point,predicateintersects)再进行一次精确判断即可。
classPolygonIndexWithRtree: 使用R树空间索引的多边形查询 性能最快
x Rtree使用动态边界体积树专门为空间查询优化 def__init__(self,polygons_dict):初始化R树索引try:fromshapely.strtreeimportSTRtreeexceptImportError:print(⚠️ 需要 shapely
0使用 pip install --upgrade shapely)raiseself.polygons{}self.geometry_list[]self.name_list[]# 构建几何体列表用于R树forname,coordsinpolygons_dict.items():polygonPolygon(coords)self.polygons[name]polygon self.geometry_list.append(polygon)self.name_list.append(name)# 创建R树索引self.treeSTRtree(self.geometry_list)deffind_point_in_polygons(self,point_lon,point_lat): 使用R树查询点所在的多边形 pointPoint(point_lon,point_lat)result[]# R树查询获取可能包含该点的多边形索引candidates_idxself.tree.query(point,predicateintersects)# 精确检查因为intersects只是粗略检查foridxincandidates_idx:ifself.geometry_list[idx].contains(point):result.append(self.name_list[idx])returnresultdeffind_points_in_polygons_batch(self,points_list):批量查询results[]forlon,latinpoints_list:polygonsself.find_point_in_polygons(lon,lat)results.append({point:(lon,lat),polygons:polygons})returnresults4️⃣ 性能表现在示例代码中使用1000 个多边形500 个查询点代码链接点击查看测试结果多边形越多R 树优势越明显。
三种方案的对比
总结方法时间复杂度构建成本查询速度适用场景直接遍历O(N)无慢小数据、一次性分析包围盒过滤O(N)极低中中等规模围栏R 树索引O(log N)中极快大规模、实时场景
工程选型建议多边形 100 个 直接 contains 判断多边形 100 ~ 1000 包围盒过滤多边形 1000 R 树索引强烈推荐 关注我 ·
获取更多内容 南墨的技术小栈这里是我的个人知识分享空间。
我会定期整理和分享工作与学习中积累的经验与资源内容涵盖算法分享 —— 深入讲解算法原理、实现思路与代码示例。
工具分享 —— 推荐实用工具与脚本包括我个人开发的小工具和精选开源工具。
开源项目 —— 精选 GitHub 上高星项目拆解原理、
使用方法和最佳实践。
数据分享 —— 工作学习中收集整理的数据资源。
无论你是技术爱好者、算法研究者还是对数据与开源感兴趣的朋友这里都希望能成为你学习、探索和实践的参考空间。
若在阅读或使用过程中有任何疑问欢迎在公众号私信我我会尽快与您交流。