python - 在 Python 3.6 中计算多边形和 shapefile 之间的重叠

Question

我想计算 shapefile 和多边形之间的重叠百分比。我正在使用 Cartopy 和 Matplotlib 并创建了此处显示的 map :

显示了欧洲的一部分(使用下载的形状文件 here )和任意矩形。假设我想计算矩形覆盖比利时的百分比。我该怎么做？下面显示了到目前为止的代码。

import matplotlib.pyplot as plt
import cartopy.crs as ccrs
import cartopy.io.shapereader as shapereader
from shapely.geometry import Polygon
from descartes import PolygonPatch

#create figure
fig1 = plt.figure(figsize=(10,10))  
PLT = plt.axes(projection=ccrs.PlateCarree())
PLT.set_extent([-10,10,45,55])
PLT.gridlines()

#import and display shapefile
fname = r'C:\Users\Me\ne_50m_admin_0_countries.shp'
adm1_shapes = list(shapereader.Reader(fname).geometries())
PLT.add_geometries(adm1_shapes, ccrs.PlateCarree(),
              edgecolor='black', facecolor='gray', alpha=0.5)

#create arbitrary polygon
x3 = 4
x4 = 5
y3 = 50
y4 = 52

poly = Polygon([(x3,y3),(x3,y4),(x4,y4),(x4,y3)])
PLT.add_patch(PolygonPatch(poly,  fc='#cc00cc', ec='#555555', alpha=0.5, 
zorder=5))

Answer 1

基于其他人发布的所有答案/评论，最终当我在最后添加这些行时它起作用了。没有响应者的帮助，我不可能做到这一点:

#find the Belgium polygon. 
for country in shapereader.Reader(fname).records(): 
    if country.attributes['SOV_A3'] == 'BEL': 
        Belgium = country.geometry 
        break 
PLT.add_patch(PolygonPatch(poly.intersection(Belgium), fc='#009900', alpha=1)) 

#calculate coverage 
x = poly.intersection(Belgium) 
print ('Coverage: ', x.area/Belgium.area*100,'%')