如何在同一地理数据框中获取多边形与其他多边形的交集?
How to get the intersection of a polygon with others in the same geodataframe?
我有一个包含数千个多边形的 shapefile。他们中的很多人接触但不交叉。我需要获取触摸多边形的 public 线。
我尝试使用以下函数来实现我的目的,但输出显示一些 MultiLineString
线只有两个点,应该是一个整体 LineString
。
def calcu_intersect_lines(cgidf):
intersection = gpd.GeoDataFrame(columns=['geometry'], crs=cgidf.crs)
while len(cgidf) > 1:
choose = cgidf.iloc[0]
cgidf.drop(cgidf.index[0], inplace=True)
for i in range(len(cgidf.index)):
cgids = cgidf.iloc[i]
if choose.geometry.exterior.intersects(cgids.geometry.exterior):
intersects = choose.geometry.exterior.intersection(cgids.geometry.exterior)
index = len(intersection)
intersection.loc[index] = [intersects]
else:
continue
return intersection
对于 MultiLineString
,我尝试使用 shapely.geometry.LineString.union()
函数将两条短线连接到同一 MultiLineString
中(如果它们相互接触)。但结果也显示 MultiLineString
。
geopandas 本身的交集函数似乎也会产生 MultiLineString
。
是否有任何方法return获得正常结果(对于连续的 public 行,LineString
而不是 MultiLineString
)?
这是输入和输出数据的一个小例子:
a = Polygon(((0, 0), (0, 0.5), (0.5, 1), (1, 0.5), (1, 0), (0.5, -0.5), (0, 0)))
b = Polygon(((0, 0.5), (0.5, 1), (1, 0.5), (1, 2), (0, 0.5)))
c = Polygon(((1, 0.5), (1, 0), (0.5, -0.5), (1.5, -1), (1, 0.5)))
gdf = gpd.GeoDataFrame(columns=['geometry'], data = [a, b, c])
h = calcu_intersect_lines(gdf)
这里是 h
的值:
index geometry
0 MULTILINESTRING ((0 0.5, 0.5 1), (0.5 1, 1 0.5))
1 MULTILINESTRING ((1 0.5, 1 0), (1 0, 0.5 -0.5))
两个MultiLineString
中的LineString
分别有public点(0.5, 1)
和(1, 0)
。
我想要的结果是这样的:
index geometry
0 LINESTRING (0 0.5, 0.5 1, 1 0.5))
1 LINESTRING (1 0.5, 1 0, 0.5 -0.5))
可能的解决方案:
在评论中,有人建议我替换下面这行
intersection.loc[index] = [intersects]
来自
intersection.loc[index] = [LineString([*intersects[0].coords, *map(lambda x: x.coords[1], intersects[1:])])]
在我的简单示例中效果很好。然而,对于真正的 shapefile,它会比这复杂得多。可能有以下几种情况:
两个多边形超过一条 public 线。
from shapely.geometry import Polygon
a = Polygon(((0., 0.), (0., 0.5), (0.5, 1.), (1., 0.5), (1., 0.), (0.5, -0.5), (0., 0.)))
b = Polygon(((0., 0.5), (0.5, 1.), (1.2, 0.7), (1., 0.), (0.5, -0.5), (2., 0.5), (0., 2.)))
对于a
和b
,他们有两条公共线路LineString(((0., 0.5), (0.5, 1.)))
和LineString(((1., 0.), (0.5, -0.5)))
。在这种情况下,我可以简单地使用 intersects
函数来测试线条是否接触。但是接下来还有一个问题:
MultiLineString
中的行顺序不对。
from shapely.geometry import MultiLineString
ml = MultiLineString((((2, 3), (3, 4)), ((0, 2), (2, 3))))
对于ml
,这个建议会return一个错误的结果。
你对上面的第二个例子有什么想法吗?
感谢 Georgy 和其他贡献者的帮助,我已经解决了我的问题。
here中介绍的函数shapely.ops.linemerge()
是我解决的重点
我post我的解决方案在这里:
from shapely import ops
def union_multils(ml):
'''Union touched LineStrings in MultiLineString or GeometryCollection.
Parameter
---------
ml: GeometryCollection, MultiLineString or LineString
return
------
ul: MultiLineString or LineString: a MultiLineString suggest the LineStrings
in input ml is not connect entitly.
'''
# Drop Point and other geom_type(if exist) out
ml = list(ml)
ml = [l for l in ml if l.geom_type == 'LineString']
# Union
if len(ml) == 1 and ml[0].geom_type == 'LineString':
ul = ml[0]
else:
ul = ops.linemerge(ml)
return ul
我有一个包含数千个多边形的 shapefile。他们中的很多人接触但不交叉。我需要获取触摸多边形的 public 线。
我尝试使用以下函数来实现我的目的,但输出显示一些 MultiLineString
线只有两个点,应该是一个整体 LineString
。
def calcu_intersect_lines(cgidf):
intersection = gpd.GeoDataFrame(columns=['geometry'], crs=cgidf.crs)
while len(cgidf) > 1:
choose = cgidf.iloc[0]
cgidf.drop(cgidf.index[0], inplace=True)
for i in range(len(cgidf.index)):
cgids = cgidf.iloc[i]
if choose.geometry.exterior.intersects(cgids.geometry.exterior):
intersects = choose.geometry.exterior.intersection(cgids.geometry.exterior)
index = len(intersection)
intersection.loc[index] = [intersects]
else:
continue
return intersection
对于 MultiLineString
,我尝试使用 shapely.geometry.LineString.union()
函数将两条短线连接到同一 MultiLineString
中(如果它们相互接触)。但结果也显示 MultiLineString
。
geopandas 本身的交集函数似乎也会产生 MultiLineString
。
是否有任何方法return获得正常结果(对于连续的 public 行,LineString
而不是 MultiLineString
)?
这是输入和输出数据的一个小例子:
a = Polygon(((0, 0), (0, 0.5), (0.5, 1), (1, 0.5), (1, 0), (0.5, -0.5), (0, 0)))
b = Polygon(((0, 0.5), (0.5, 1), (1, 0.5), (1, 2), (0, 0.5)))
c = Polygon(((1, 0.5), (1, 0), (0.5, -0.5), (1.5, -1), (1, 0.5)))
gdf = gpd.GeoDataFrame(columns=['geometry'], data = [a, b, c])
h = calcu_intersect_lines(gdf)
这里是 h
的值:
index geometry
0 MULTILINESTRING ((0 0.5, 0.5 1), (0.5 1, 1 0.5))
1 MULTILINESTRING ((1 0.5, 1 0), (1 0, 0.5 -0.5))
两个MultiLineString
中的LineString
分别有public点(0.5, 1)
和(1, 0)
。
我想要的结果是这样的:
index geometry
0 LINESTRING (0 0.5, 0.5 1, 1 0.5))
1 LINESTRING (1 0.5, 1 0, 0.5 -0.5))
可能的解决方案:
在评论中,有人建议我替换下面这行
intersection.loc[index] = [intersects]
来自
intersection.loc[index] = [LineString([*intersects[0].coords, *map(lambda x: x.coords[1], intersects[1:])])]
在我的简单示例中效果很好。然而,对于真正的 shapefile,它会比这复杂得多。可能有以下几种情况:
两个多边形超过一条 public 线。
from shapely.geometry import Polygon a = Polygon(((0., 0.), (0., 0.5), (0.5, 1.), (1., 0.5), (1., 0.), (0.5, -0.5), (0., 0.))) b = Polygon(((0., 0.5), (0.5, 1.), (1.2, 0.7), (1., 0.), (0.5, -0.5), (2., 0.5), (0., 2.)))
对于
a
和b
,他们有两条公共线路LineString(((0., 0.5), (0.5, 1.)))
和LineString(((1., 0.), (0.5, -0.5)))
。在这种情况下,我可以简单地使用intersects
函数来测试线条是否接触。但是接下来还有一个问题:MultiLineString
中的行顺序不对。from shapely.geometry import MultiLineString ml = MultiLineString((((2, 3), (3, 4)), ((0, 2), (2, 3))))
对于
ml
,这个建议会return一个错误的结果。 你对上面的第二个例子有什么想法吗?
感谢 Georgy 和其他贡献者的帮助,我已经解决了我的问题。
here中介绍的函数shapely.ops.linemerge()
是我解决的重点
我post我的解决方案在这里:
from shapely import ops
def union_multils(ml):
'''Union touched LineStrings in MultiLineString or GeometryCollection.
Parameter
---------
ml: GeometryCollection, MultiLineString or LineString
return
------
ul: MultiLineString or LineString: a MultiLineString suggest the LineStrings
in input ml is not connect entitly.
'''
# Drop Point and other geom_type(if exist) out
ml = list(ml)
ml = [l for l in ml if l.geom_type == 'LineString']
# Union
if len(ml) == 1 and ml[0].geom_type == 'LineString':
ul = ml[0]
else:
ul = ops.linemerge(ml)
return ul