与叠加的集合操作#
在处理多个空间数据集时——特别是多个多边形或线数据集——用户常常希望根据这些数据集重叠（或不重叠）的位置创建新形状。这些操作通常使用集合的语言来描述——交集、并集和差集。通过overlay()方法，GeoPandas库提供了这些类型的操作。
基本思想如下面的图所示，但请记住，叠加操作在数据框级别进行，而不是针对单个几何图形，并且两个数据框的属性都会保留。实际上，对于左边的每个形状GeoDataFrame，此操作都是针对右边的每个形状GeoDataFrame执行的：
来源：QGIS文档
注意
对于熟悉shapely库的用户：overlay()可以被视为提供标准shapely集合操作的版本，这些版本处理将集合操作应用于两个GeoSeries的复杂性。标准shapely集合操作同样可作为GeoSeries方法使用。
不同的叠加操作#
首先，创建一些示例数据：
In [1]: from shapely.geometry import Polygon
In [2]: polys1 = geopandas.GeoSeries([Polygon([(0,0), (2,0), (2,2), (0,2)]),
...: Polygon([(2,2), (4,2), (4,4), (2,4)])])
...:
In [3]: polys2 = geopandas.GeoSeries([Polygon([(1,1), (3,1), (3,3), (1,3)]),
...: Polygon([(3,3), (5,3), (5,5), (3,5)])])
...:
In [4]: df1 = geopandas.GeoDataFrame({'geometry': polys1, 'df1':[1,2]})
In [5]: df2 = geopandas.GeoDataFrame({'geometry': polys2, 'df2':[1,2]})
这两个GeoDataFrames有一些重叠区域：
In [6]: ax = df1.plot(color='red');
In [7]: df2.plot(ax=ax, color='green', alpha=0.5);
上述示例说明了不同的叠加模式。
overlay() 方法将确定通过叠加两个输入的 GeoDataFrames 得到的所有个体几何形状的集合。此结果涵盖两个输入 GeoDataFrames 所覆盖的区域，并且保留由两个 GeoDataFrames 的组合边界定义的所有独特区域。
注意
出于历史原因，叠加方法也作为一个顶层函数 overlay() 可用。建议使用该方法，因为该函数可能在将来被弃用。
使用 how='union' 时，所有可能的几何图形都会被返回：
In [8]: res_union = df1.overlay(df2, how='union')
In [9]: res_union
Out[9]:
df1 df2 geometry
0 1.0 1.0 POLYGON ((2 2, 2 1, 1 1, 1 2, 2 2))
1 2.0 1.0 POLYGON ((2 2, 2 3, 3 3, 3 2, 2 2))
2 2.0 2.0 POLYGON ((4 4, 4 3, 3 3, 3 4, 4 4))
3 1.0 NaN POLYGON ((2 0, 0 0, 0 2, 1 2, 1 1, 2 1, 2 0))
4 2.0 NaN MULTIPOLYGON (((3 4, 3 3, 2 3, 2 4, 3 4)), ((4...
5 NaN 1.0 MULTIPOLYGON (((2 3, 2 2, 1 2, 1 3, 2 3)), ((3...
6 NaN 2.0 POLYGON ((3 5, 5 5, 5 3, 4 3, 4 4, 3 4, 3 5))
In [10]: ax = res_union.plot(alpha=0.5, cmap='tab10')
In [11]: df1.plot(ax=ax, facecolor='none', edgecolor='k');
In [12]: df2.plot(ax=ax, facecolor='none', edgecolor='k');
其他 how 操作将返回这些几何形状的不同子集。
使用 how='intersection'，它只返回那些被两个 GeoDataFrames 包含的几何形状：
In [13]: res_intersection = df1.overlay(df2, how='intersection')
In [14]: res_intersection
Out[14]:
df1 df2 geometry
0 1 1 POLYGON ((2 2, 2 1, 1 1, 1 2, 2 2))
1 2 1 POLYGON ((2 2, 2 3, 3 3, 3 2, 2 2))
2 2 2 POLYGON ((4 4, 4 3, 3 3, 3 4, 4 4))
In [15]: ax = res_intersection.plot(cmap='tab10')
In [16]: df1.plot(ax=ax, facecolor='none', edgecolor='k');
In [17]: df2.plot(ax=ax, facecolor='none', edgecolor='k');
how='symmetric_difference' 是 'intersection' 的相反，返回仅属于一个 GeoDataFrame 但不属于两个 GeoDataFrame 的几何形状：
In [18]: res_symdiff = df1.overlay(df2, how='symmetric_difference')
In [19]: res_symdiff
Out[19]:
df1 df2 geometry
0 1.0 NaN POLYGON ((2 0, 0 0, 0 2, 1 2, 1 1, 2 1, 2 0))
1 2.0 NaN MULTIPOLYGON (((3 4, 3 3, 2 3, 2 4, 3 4)), ((4...
2 NaN 1.0 MULTIPOLYGON (((2 3, 2 2, 1 2, 1 3, 2 3)), ((3...
3 NaN 2.0 POLYGON ((3 5, 5 5, 5 3, 4 3, 4 4, 3 4, 3 5))
In [20]: ax = res_symdiff.plot(cmap='tab10')
In [21]: df1.plot(ax=ax, facecolor='none', edgecolor='k');
In [22]: df2.plot(ax=ax, facecolor='none', edgecolor='k');
要获取属于 df1 的几何图形但不包含在 df2 中的部分，可以使用 how='difference':
In [23]: res_difference = df1.overlay(df2, how='difference')
In [24]: res_difference
Out[24]:
geometry df1
0 POLYGON ((2 0, 0 0, 0 2, 1 2, 1 1, 2 1, 2 0)) 1
1 MULTIPOLYGON (((3 4, 3 3, 2 3, 2 4, 3 4)), ((4... 2
In [25]: ax = res_difference.plot(cmap='tab10')
In [26]: df1.plot(ax=ax, facecolor='none', edgecolor='k');
In [27]: df2.plot(ax=ax, facecolor='none', edgecolor='k');
最后，使用 how='identity'，结果包含了 df1 的表面，但几何体是通过将 df1 与 df2 叠加而获得的：
In [28]: res_identity = df1.overlay(df2, how='identity')
In [29]: res_identity
Out[29]:
df1 df2 geometry
0 1.0 1.0 POLYGON ((2 2, 2 1, 1 1, 1 2, 2 2))
1 2.0 1.0 POLYGON ((2 2, 2 3, 3 3, 3 2, 2 2))
2 2.0 2.0 POLYGON ((4 4, 4 3, 3 3, 3 4, 4 4))
3 1.0 NaN POLYGON ((2 0, 0 0, 0 2, 1 2, 1 1, 2 1, 2 0))
4 2.0 NaN MULTIPOLYGON (((3 4, 3 3, 2 3, 2 4, 3 4)), ((4...
In [30]: ax = res_identity.plot(cmap='tab10')
In [31]: df1.plot(ax=ax, facecolor='none', edgecolor='k');
In [32]: df2.plot(ax=ax, facecolor='none', edgecolor='k');
叠加杂货示例#
首先，加载芝加哥社区区域和杂货示例数据集并选择：
In [33]: import geodatasets
In [34]: chicago = geopandas.read_file(geodatasets.get_path("geoda.chicago_commpop"))
In [35]: groceries = geopandas.read_file(geodatasets.get_path("geoda.groceries"))
# Project to crs that uses meters as distance measure
In [36]: chicago = chicago.to_crs("ESRI:102003")
In [37]: groceries = groceries.to_crs("ESRI:102003")
为了说明overlay()方法，考虑以下情况，其中希望识别每个区域的“服务”部分——定义为距杂货店1公里内的区域——使用一个GeoDataFrame的社区区域和一个GeoDataFrame的杂货店。
# Look at Chicago:
In [38]: chicago.plot();
# Now buffer groceries to find area within 1km.
# Check CRS -- USA Contiguous Albers Equal Area, units of meters.
In [39]: groceries.crs
Out[39]:
<Projected CRS: ESRI:102003>
Name: USA_Contiguous_Albers_Equal_Area_Conic
Axis Info [cartesian]:
- E[east]: Easting (metre)
- N[north]: Northing (metre)
Area of Use:
- name: United States (USA) - CONUS onshore - Alabama; Arizona; Arkansas; California; Colorado; Connecticut; Delaware; Florida; Georgia; Idaho; Illinois; Indiana; Iowa; Kansas; Kentucky; Louisiana; Maine; Maryland; Massachusetts; Michigan; Minnesota; Mississippi; Missouri; Montana; Nebraska; Nevada; New Hampshire; New Jersey; New Mexico; New York; North Carolina; North Dakota; Ohio; Oklahoma; Oregon; Pennsylvania; Rhode Island; South Carolina; South Dakota; Tennessee; Texas; Utah; Vermont; Virginia; Washington; West Virginia; Wisconsin; Wyoming.
- bounds: (-124.79, 24.41, -66.91, 49.38)
Coordinate Operation:
- name: USA_Contiguous_Albers_Equal_Area_Conic
- method: Albers Equal Area
Datum: North American Datum 1983
- Ellipsoid: GRS 1980
- Prime Meridian: Greenwich
# make 1km buffer
In [40]: groceries['geometry']= groceries.buffer(1000)
In [41]: groceries.plot();
要选择距离杂货店1公里内的社区区域部分，请将 how 选项设置为“intersect”，这将创建一个新的多边形集合，其中这两个图层重叠：
In [42]: chicago_cores = chicago.overlay(groceries, how='intersection')
In [43]: chicago_cores.plot(alpha=0.5, edgecolor='k', cmap='tab10');
更改 how 选项允许进行不同类型的叠加操作。例如，如果您对距离杂货店 远 的芝加哥部分感兴趣（边缘地带），您可以计算这两者之间的差异。
In [44]: chicago_peripheries = chicago.overlay(groceries, how='difference')
In [45]: chicago_peripheries.plot(alpha=0.5, edgecolor='k', cmap='tab10');
keep_geom_type 关键字#
在默认设置下，overlay() 仅返回与 GeoDataFrame（左侧）具有相同几何类型的几何图形，其中多边形和多面体被视为同一类型（其他类型同样如此）。你可以使用 keep_geom_type 选项控制这一行为，该选项默认设置为 True。一旦设置为 False，overlay 将返回由所选集合操作结果产生的所有几何类型。不同的类型可以例如由于相交的几何体而产生，其中两个多边形在一条线或一个点上相交。