python - Cartopy:绘制删除国家边界的海岸线

Question

我想用一种颜色画出中国的轮廓，同时用另一种颜色显示全局海岸线。我第一次这样做的尝试如下:

import matplotlib.pyplot as plt
import cartopy.crs as ccrs
import cartopy.feature as feature
import cartopy.io.shapereader as shapereader


countries = shapereader.natural_earth(resolution='110m',
                                      category='cultural',
                                      name='admin_0_countries')

# Find the China boundary polygon.
for country in shapereader.Reader(countries).records():
    if country.attributes['su_a3'] == 'CHN':
        china = country.geometry
        break
else:
    raise ValueError('Unable to find the CHN boundary.')


plt.figure(figsize=(8, 4))
ax = plt.axes(projection=ccrs.PlateCarree())

ax.set_extent([50, 164, 5, 60], ccrs.PlateCarree())

ax.add_feature(feature.LAND)
ax.add_feature(feature.OCEAN)
ax.add_feature(feature.COASTLINE, linewidth=4)

ax.add_geometries([china], ccrs.Geodetic(), edgecolor='red',
                  facecolor='none')

plt.show()

我把海岸线弄得很厚，这样你就可以看到它们与国家边界重叠的事实。

我的问题是:有没有办法移除国家轮廓旁边的海岸线，这样我就不会在视觉上看到两条线相互影响？

Note: This question was asked to me directly via email, and I chose to post my response here so that others may learn/benefit from a solution.

Answer 1

Natural Earth 集合中没有名为“没有中国边界的海岸线”的数据集，因此我们将不得不自己制作它。为此，我们需要使用整形操作，尤其是 difference 方法。

差异方法如下图所示(取自Shapely's docs)。下面突出显示了两个示例圆圈(a 和 b)的区别:

然后，我们的目标是达到编写 coaSTLine.difference(china) 的目的，并将这个结果可视化为我们的海岸线。强>

有很多方法可以做到这一点。 GeoPandas 和 Fiona 是两种可以提供非常可读的结果的技术。不过在这种情况下，让我们使用 cartopy 提供的工具:

首先，我们掌握了中国边界(另见:cartopy shapereader docs)。

import cartopy.io.shapereader as shapereader


countries = shapereader.natural_earth(resolution='110m',
                                      category='cultural',
                                      name='admin_0_countries')

# Find the China boundary polygon.
for country in shapereader.Reader(countries).records():
    if country.attributes['su_a3'] == 'CHN':
        china = country.geometry
        break
else:
    raise ValueError('Unable to find the CHN boundary.')

接下来，我们掌握海岸线几何图形:

coast = shapereader.natural_earth(resolution='110m',
                                  category='physical',
                                  name='coastline')

coastlines = shapereader.Reader(coast).geometries()

现在，把中国带出海岸线:

coastlines_m_china = [geom.difference(china)
                      for geom in coastlines]

当我们将其可视化时，我们发现差异并不十分完美:

我们不想要黑线的原因是 Natural Earth 海岸线数据集与国家数据集的派生方式不同，因此它们不是完全重合的坐标。

为了绕过这个事实，可以对中国边界应用一个小的“hack”来扩大边界以达到这个交叉点的目的。 buffer方法非常适合此目的。

# Buffer the Chinese border by a tiny amount to help the coordinate
# alignment with the coastlines.
coastlines_m_china = [geom.difference(china.buffer(0.001))
                      for geom in coastlines]

有了这个“hack”，我得到了以下结果(包含完整代码以确保完整性):

import matplotlib.pyplot as plt
import cartopy.crs as ccrs
import cartopy.feature as feature
import cartopy.io.shapereader as shapereader


coast = shapereader.natural_earth(resolution='110m',
                                  category='physical',
                                  name='coastline')

countries = shapereader.natural_earth(resolution='110m',
                                      category='cultural',
                                      name='admin_0_countries')

# Find the China boundary polygon.
for country in shapereader.Reader(countries).records():
    if country.attributes['su_a3'] == 'CHN':
        china = country.geometry
        break
else:
    raise ValueError('Unable to find the CHN boundary.')

coastlines = shapereader.Reader(coast).geometries() 

# Hack to get the border to intersect cleanly with the coastline.
coastlines_m_china = [geom.difference(china.buffer(0.001))
                      for geom in coastlines]

ax = plt.axes(projection=ccrs.PlateCarree())

ax.set_extent([50, 164, 5, 60], ccrs.PlateCarree())
ax.add_feature(feature.LAND)
ax.add_feature(feature.OCEAN)

ax.add_geometries(coastlines_m_china, ccrs.Geodetic(), edgecolor='black', facecolor='none', lw=4)
ax.add_geometries([china], ccrs.Geodetic(), edgecolor='red', facecolor='none')

plt.show()