问题描述
我有两个地理数据框:一个包含许多线串,另一个包含许多多边形。 这些线和多边形彼此相交。 我试图实现的输出是一个新的地理数据框,其中包含在与多边形相交的任何位置分割的链接。
简化的测试代码如下:
import geopandas
from shapely.geometry import Polygon, LineString
#Create geodataframe containing 2 links:
linkID = ['1','2']
link_geom = [LineString([(0, 0), (10, 10)]),LineString([(10, 10), (20, 10)])]
a = {'linkID':linkID,'geometry':link_geom}
gdf_links = geopandas.GeoDataFrame(a)
#Create geodataframe containing 2 polygons:
polyID = ['100','200']
poly_geom = [Polygon([(2, 1), (2, 3), (4, 3), (4, 1)]),Polygon([(15, 7), (15, 13), (18, 13), (18, 7)])]
b = {'polyID':polyID,'geometry':poly_geom}
gdf_poly = geopandas.GeoDataFrame(b)
通过运行以下代码,您可以看到这些线与多边形相交:
links = gdf_links.unary_union
polys = gdf_poly.unary_union
geopandas.GeoSeries([links,polys]).plot(cmap='tab10')
它生成下图(以红色/绿色添加注释以显示所涉及的各种元素):
将两条线在与 polgyon 相交的地方拆分将产生 6 个线段(如上图中的绿色标记)。
理想情况下,我想要实现的输出如下所示,“poly_intersect”列列出了发生分割的多边形的 ID。 请注意,项目 2 和 5 列出了两个多边形 ID,因为它们与多边形相交两次。
1楼
您可以通过将线与多边形进行“叠加”(在叠加时使用“联合”操作)来获得这些线串。 但是,它不会像您为“poly_intersect”列描述的那样准确地给出预期结果:
>> geopandas.overlay(gdf_links, gdf_poly, how='union').explode().reset_index(drop=True)
linkID polyID geometry
0 1 100 LINESTRING (2 2, 3 3)
1 2 200 LINESTRING (15 10, 18 10)
2 1 NaN LINESTRING (0 0, 2 2)
3 1 NaN LINESTRING (3 3, 10 10)
4 2 NaN LINESTRING (10 10, 15 10)
5 2 NaN LINESTRING (18 10, 20 10)
一些备注:
-
我正在使用
explode()
,因为overlay()
操作为与多边形相交之前和之后的一个线串的部分返回MultiLineStrings。 爆炸将分裂那些“多”几何。 - “polyID”列指示哪些线串源自与哪些多边形相交。 但它没有关于它们接触的多边形的其他线串的信息。
因此,如果您还需要该信息(完整的“poly_intersect”列),您可以进行一些后续分析,检查哪些线串正在接触哪些多边形。