如何在GeoDataFrame中创建含MultiPolygon的共享边多边形字典？

我需要在包含Polygon与MultiPolygon的GeoDataFrame中，创建一个字典，记录所有共享边的多边形（多边形仅相交但不交叉）。

示例数据代码 #

from shapely.geometry import Polygon, MultiPolygon
import geopandas as gpd
import matplotlib.pyplot as plt

polys = gpd.GeoSeries([Polygon([(0,0), (2,0), (2, 1.5), (2,2), (0,2)]),
                     Polygon([(0,2), (2,2), (2,4), (0,4)]),
                     Polygon([(2,0), (5,0), (5,1.5), (2,1.5)]),
                     Polygon([(3,3), (5,3), (5,5), (3,5)]),
                       MultiPolygon([Polygon([(6,1), (8, 1), (8, 3), (6, 3)], 
                                     [[(6.5, 1.5), (7.5, 1.5), (7.5, 2.5), (6.5, 2.5)][::-1]]
                                            )])])

fp = gpd.GeoDataFrame({'geometry': polys, 'name': ['a', 'b', 'c', 'd', 'e'],
                      'grnd': [25, 25, 25, 25, 25],
                      'rf': [29, 35, 26, 31, 28]})

fig, ax = plt.subplots(figsize=(6, 6))
fp.plot(ax=ax, alpha=0.3, cmap='tab10', edgecolor='k',)
fp.apply(lambda x: ax.annotate(text=x['name'], xy=x.geometry.centroid.coords[0], ha='center'), axis=1)
plt.show()

尝试的代码及报错 #

我尝试了以下代码：

from itertools import chain
from shapely.geometry import LineString

# 创建所有边线段的序列
lines = fp.geometry.apply(lambda x: (list(map(
    LineString, 
    zip(x.boundary.coords[:-1], x.boundary.coords[1:])))
                         if isinstance(x, MultiPolygon)
                         else list(chain(*list(list(map(
                             LineString, 
    zip(x.boundary.coords[:-1], x.boundary.coords[1:])))
                                               for poly in x.geoms)))
                         )).explode()

result = {
        line : list(fp.loc[
            (fp.geometry.touches(line)) # 线段与多边形接触
            & (fp.geometry.intersection(line).length > 0), # 交集长度大于0（非点接触）
            'rf'].values)
        for line in lines}

出现报错：AttributeError: 'Polygon' object has no attribute 'geoms'

解决方案 #

报错原因是判断逻辑搞反了：MultiPolygon对象才有geoms属性，Polygon对象没有。修正线段提取逻辑后，即可正常运行。同时建议用线段的WKT字符串作为字典键，避免shapely对象作为键的哈希问题。

修正后的完整代码：

from itertools import chain
from shapely.geometry import Polygon, MultiPolygon, LineString
import geopandas as gpd
import matplotlib.pyplot as plt

# 示例数据
polys = gpd.GeoSeries([Polygon([(0,0), (2,0), (2, 1.5), (2,2), (0,2)]),
                     Polygon([(0,2), (2,2), (2,4), (0,4)]),
                     Polygon([(2,0), (5,0), (5,1.5), (2,1.5)]),
                     Polygon([(3,3), (5,3), (5,5), (3,5)]),
                       MultiPolygon([Polygon([(6,1), (8, 1), (8, 3), (6, 3)], 
                                     [[(6.5, 1.5), (7.5, 1.5), (7.5, 2.5), (6.5, 2.5)][::-1]]
                                            )])])

fp = gpd.GeoDataFrame({'geometry': polys, 'name': ['a', 'b', 'c', 'd', 'e'],
                      'grnd': [25, 25, 25, 25, 25],
                      'rf': [29, 35, 26, 31, 28]})

# 定义线段提取函数，区分Polygon和MultiPolygon
def extract_boundary_segments(geom):
    if isinstance(geom, Polygon):
        # Polygon直接提取边界的连续线段
        coords = geom.boundary.coords
        return [LineString(pair) for pair in zip(coords[:-1], coords[1:])]
    elif isinstance(geom, MultiPolygon):
        # MultiPolygon遍历每个子Polygon提取线段
        segments = []
        for poly in geom.geoms:
            coords = poly.boundary.coords
            segments.extend([LineString(pair) for pair in zip(coords[:-1], coords[1:])])
        return segments
    return []

# 生成所有边界线段
lines = fp.geometry.apply(extract_boundary_segments).explode()

# 生成目标字典，键为线段WKT字符串，值为对应多边形的rf值列表
result = {}
for line in lines:
    # 筛选出与线段有非点接触的多边形，提取rf值
    matching_rfs = fp.loc[
        fp.geometry.touches(line) & (fp.geometry.intersection(line).length > 0),
        'rf'
    ].tolist()
    result[line.wkt] = matching_rfs

# 打印结果示例
for key, value in list(result.items())[:10]:
    print(f"{key}: {value}")

运行后会得到类似如下的结果：

LINESTRING (0 0, 2 0): [29]
LINESTRING (2 0, 2 1.5): [29, 26]
LINESTRING (2 1.5, 2 2): [29]
LINESTRING (2 2, 0 2): [29, 35]
LINESTRING (0 2, 0 0): [29]
LINESTRING (0 2, 2 2): [35, 29]
LINESTRING (2 2, 2 4): [35]
LINESTRING (2 4, 0 4): [35]
LINESTRING (0 4, 0 2): [35]
LINESTRING (2 0, 5 0): [26]

内容的提问来源于stack exchange，提问作者arkriger