将美国州边界添加到绘图中

Question

有人可以告诉我如何修改下面的脚本以使用 Geopandas 添加美国国家边界吗？我已经尝试过，但似乎无法让它发挥作用。另外，如果有人可以推荐特定的州边界形状文件，请告诉我。否则，这个情节看起来就像我想要的那样。请注意，我不想使用map.drawstates()，因为它也绘制其他国家/地区的州。我只想要美国的州界。

from netCDF4 import Dataset as NetCDFFile
import matplotlib.pyplot as plt
import numpy as np
from mpl_toolkits.basemap import Basemap
import matplotlib.colors as mcolors
import geopandas as gpd

# Read the U.S. state borders shapefile
us_states = gpd.read_file('C:/Users/scapa/Downloads/tl_2012_us_state/tl_2012_us_state.shp')

nc = NetCDFFile('C:/Users/scapa/Downloads/nclwAXbfFo6UE.nc')
lat = nc.variables['lat'][:]
lon = nc.variables['lon'][:]
prate = nc.variables['VAR'][:]

# Determine the extent of the shapefile
state_bounds = us_states.total_bounds  # Get the bounding box of the state borders

# Calculate the map extent
map_extent = [state_bounds[0], state_bounds[2], state_bounds[1], state_bounds[3]]

fig, ax = plt.subplots(figsize=(10, 8))
map = Basemap(llcrnrlon=230., llcrnrlat=10., urcrnrlon=305., urcrnrlat=55.)
lons, lats = np.meshgrid(lon, lat)
x, y = map(lons, lats)

# Set the levels and colors
levels = np.arange(-3.6, 3.6, 0.4)
cmap = plt.get_cmap('BrBG', len(levels) - 1)

# Create a custom colormap with white for values between -0.4 and 0.4
cmaplist = [(1, 1, 1, 1) if -0.4 <= val <= 0.4 else cmap((val - levels.min()) / (levels.max() - levels.min())) for val in levels]

custom_cmap = mcolors.LinearSegmentedColormap.from_list('custom', cmaplist, len(levels))

# Add U.S. state borders
us_states.boundary.plot(ax=map, linewidth=0.5)

cs = map.contourf(x, y, prate, levels, cmap=custom_cmap)
cbar = plt.colorbar(cs, orientation='horizontal', ticks=levels, pad=0.07)
cbar.set_label('Surface Precipitation Rate Anomalies (mm/day)')
cbar.ax.tick_params(labelsize=8)
map.drawcoastlines()
map.drawparallels(np.arange(10, 60, 10), labels=[1, 1, 1, 1], linewidth=0.5)
map.drawmeridians(np.arange(240, 310, 20), labels=[1, 1, 0, 1], linewidth=0.5)
plt.title('El Niño Winter Precipitation Anomalies', fontsize=15, y=1.05)
plt.text(0.5, 1.02, 'Winters of 1957-58, 1965-66, 1972-73, 1982-83, 1987-88, 1991-92, 1997-98, 2009-10, 2015-16', horizontalalignment='center', fontsize=9, fontstyle='italic',transform=ax.transAxes)

plt.savefig('era5.precip.anom.states.png', dpi=600)


plt.show()

Answer 1

事实证明 Cartopy 是解决我的问题的最佳方案。以下是对上述代码的一些关键添加/更改：

import cartopy.crs as ccrs
import cartopy.feature as cfeature

# Create a Cartopy projection
projection = ccrs.PlateCarree()

# Create a figure and axis with Cartopy projection
fig, ax = plt.subplots(figsize=(10, 8), subplot_kw={'projection': projection})

# Add coastlines
ax.coastlines(linewidth=0.5)

# Add U.S. state borders
ax.add_feature(cfeature.STATES, linestyle=':', edgecolor='black', linewidth=0.5)