效果图:
主要步骤:
1. 数据准备:下载ERA5月平均地表温度数据就好啦(不会的可以看我之前的帖子【数据下载】ERA5地表月平均数据下载-CSDN博客)
2. 数据处理:计算温度距平,掩膜四川省数据
3. 图像绘制:主要包括等值线填充,经纬网格线绘制,坐标轴标签的设置
注:找到了一版正确的地图,链接在这里国家出品审图号:GS(2024)0650号省市县shp格式地图无条件下载 (qq.com)
详细代码:着急的直接拖到最后有完整代码
步骤一:导入库包及图片存储路径并设置中文字体为宋体,西文为新罗马(没有的库包要先下好奥)
##########################################################
# 导入库包并设置字体
import numpy as np
import cartopy.crs as ccrs
import cartopy.feature as cfeature
from cartopy.mpl.gridliner import LONGITUDE_FORMATTER, LATITUDE_FORMATTER
from cartopy.mpl.ticker import LongitudeFormatter, LatitudeFormatter
import matplotlib.pyplot as plt
import matplotlib.ticker as mticker
import xarray as xr
from cnmaps import get_adm_maps
# import pygrib as pg
# 设置西文字体为新罗马字体,中文为宋体,字号为12
from matplotlib import rcParams
config = {
"font.family": 'serif',
"font.size": 12,
"mathtext.fontset": 'stix',
"font.serif": ['SimSun'],
}
rcParams.update(config)
rcParams['axes.unicode_minus']=False
步骤二:读入数据并计算平均温度及温度距平
##########################################################
datapath = r'H:/00.csdn/01data/'
figpath = r'H:/00.csdn/02fig/'
shppath = r'H:/00.csdn/04shp/ASEF:天地图_审图号:GS(2024)0650号_中国省市县地图/中国_省/中国_省1.shp'
# 读入数据
f = xr.open_dataset(datapath + 't2m_cn_jja_1950-2023.nc')
# print(f)
t = f['t2m'].loc[:, 55.:0., 70.:140.]
time = f['time']
lons = t['longitude']
lats = t['latitude']
lon, lat = np.meshgrid(lons, lats)
# 计算时间维平均,及距平
t_mean = np.mean(t, axis=0)
dt = t[-1, :, :] - t_mean
# print(np.max(dt))
# print(np.min(dt))
步骤三:掩膜数据
#掩膜数据
cn = get_adm_maps(level='省',province='四川省', only_polygon=True, record='first')
dt = cn.maskout(lon, lat, dt)
步骤四:定义绘图函数并调用函数(这里所有的东西你都可以选择要与不要,不要注释掉就好了)
##########################################################
# 画图
t_level = np.arange(0, 2, 0.2)
fig = plt.figure(figsize=(15, 6))
def plot_cont(data, location, level, extent, show):
ax = fig.add_axes(location, projection=ccrs.PlateCarree())
c1 = ax.contourf(lon, lat, data, levels=level,
extend='both', cmap='OrRd', transform=ccrs.PlateCarree())
# 添加地理信息
ax.add_feature(cfeature.OCEAN.with_scale('50m')) # 海洋填充为蓝色
ax.add_feature(cfeature.LAND.with_scale('50m'), facecolor='w') # 陆地填充为白色
# help(ax.add_feature)
# ax.add_feature(cfeature.RIVERS.with_scale('50m')) # 添加河流
# ax.add_feature(cfeature.LAKES.with_scale('50m')) # 添加湖泊
# ax.add_feature(cfeature.BORDERS.with_scale('50m'),lw=0.5) # 添加国界线
ax.add_feature(cfeature.COASTLINE.with_scale('50m'),lw=0.5) # 添加海岸线
# ax.add_geometries(Reader(shppath).geometries(), ccrs.PlateCarree(),
# facecolor='none', edgecolor='k', linewidth=0.7) #添加中国国界
# 设置显示范围
ax.set_extent(extent, crs=ccrs.PlateCarree())
# 设置是否显示经纬网格线及坐标刻度
if show == True:
# 画经纬度网格
gl = ax.gridlines(crs=ccrs.PlateCarree(), draw_labels=False, linewidth=1.2,
color='k', alpha=0.5, linestyle='--')
gl.top_labels = False # 关闭顶端的经纬度标签
gl.right_labels = False # 关闭右侧的经纬度标签
gl.xformatter = LONGITUDE_FORMATTER # x轴设为经度的格式
gl.yformatter = LATITUDE_FORMATTER # y轴设为纬度的格式
# 设置经纬度网格的间隔及坐标标签
gl.xlocator = mticker.FixedLocator(np.arange(extent[0], extent[1]+1, 2))
gl.ylocator = mticker.FixedLocator(np.arange(extent[2], extent[3]+1, 2))
ax.set_xticks(list(range(extent[0], extent[1]+1, 2)), crs=ccrs.PlateCarree())
ax.set_yticks(list(range(extent[2], extent[3]+1, 2)), crs=ccrs.PlateCarree())
lon_formatter = LongitudeFormatter(zero_direction_label=True) # zero_direction_label用来设置经度的0度加不加E和W
lat_formatter = LatitudeFormatter()
ax.xaxis.set_major_formatter(lon_formatter)
ax.yaxis.set_major_formatter(lat_formatter)
return c1
SC = plot_cont(dt, [0.1, 0.1, 1, 1], t_level, [97, 109, 26, 35], show=True)
plt.text(104, 34.5, '审图号:GS(2024)0650号', fontsize=20)
步骤五:色标设置
###############################################################################
# 色标设置
rect = [0.3, 0.015, 0.6, 0.04]
cbar_ax = fig.add_axes(rect)
cb = fig.colorbar(SC, drawedges=True, ticks=t_level,
cax=cbar_ax, orientation='horizontal')
cb.set_label('℃', fontsize=12)
cb.ax.tick_params(length=0)
步骤六:保存图像
###############################################################################
# 存图
plt.savefig(figpath+'203 单省填充', dpi=600, bbox_inches = 'tight')
plt.show()
完整代码在这里:
##########################################################
# 导入库包并设置字体
import numpy as np
import cartopy.crs as ccrs
import cartopy.feature as cfeature
from cartopy.mpl.gridliner import LONGITUDE_FORMATTER, LATITUDE_FORMATTER
from cartopy.mpl.ticker import LongitudeFormatter, LatitudeFormatter
import matplotlib.pyplot as plt
import matplotlib.ticker as mticker
import xarray as xr
from cnmaps import get_adm_maps
# import pygrib as pg
# 设置西文字体为新罗马字体,中文为宋体,字号为12
from matplotlib import rcParams
config = {
"font.family": 'serif',
"font.size": 12,
"mathtext.fontset": 'stix',
"font.serif": ['SimSun'],
}
rcParams.update(config)
rcParams['axes.unicode_minus']=False
##########################################################
datapath = r'H:/00.csdn/01data/'
figpath = r'H:/00.csdn/02fig/'
shppath = r'H:/00.csdn/04shp/ASEF:天地图_审图号:GS(2024)0650号_中国省市县地图/中国_省/中国_省1.shp'
# 读入数据
f = xr.open_dataset(datapath + 't2m_cn_jja_1950-2023.nc')
# print(f)
t = f['t2m'].loc[:, 55.:0., 70.:140.]
time = f['time']
lons = t['longitude']
lats = t['latitude']
lon, lat = np.meshgrid(lons, lats)
# 计算时间维平均,及距平
t_mean = np.mean(t, axis=0)
dt = t[-1, :, :] - t_mean
# print(np.max(dt))
# print(np.min(dt))
#掩膜数据
cn = get_adm_maps(level='省',province='四川省', only_polygon=True, record='first')
dt = cn.maskout(lon, lat, dt)
##########################################################
# 画图
t_level = np.arange(0, 2, 0.2)
fig = plt.figure(figsize=(15, 6))
def plot_cont(data, location, level, extent, show):
ax = fig.add_axes(location, projection=ccrs.PlateCarree())
c1 = ax.contourf(lon, lat, data, levels=level,
extend='both', cmap='OrRd', transform=ccrs.PlateCarree())
# 添加地理信息
ax.add_feature(cfeature.OCEAN.with_scale('50m')) # 海洋填充为蓝色
ax.add_feature(cfeature.LAND.with_scale('50m'), facecolor='w') # 陆地填充为白色
# help(ax.add_feature)
# ax.add_feature(cfeature.RIVERS.with_scale('50m')) # 添加河流
# ax.add_feature(cfeature.LAKES.with_scale('50m')) # 添加湖泊
# ax.add_feature(cfeature.BORDERS.with_scale('50m'),lw=0.5) # 添加国界线
ax.add_feature(cfeature.COASTLINE.with_scale('50m'),lw=0.5) # 添加海岸线
# ax.add_geometries(Reader(shppath).geometries(), ccrs.PlateCarree(),
# facecolor='none', edgecolor='k', linewidth=0.7) #添加中国国界
# 设置显示范围
ax.set_extent(extent, crs=ccrs.PlateCarree())
# 设置是否显示经纬网格线及坐标刻度
if show == True:
# 画经纬度网格
gl = ax.gridlines(crs=ccrs.PlateCarree(), draw_labels=False, linewidth=1.2,
color='k', alpha=0.5, linestyle='--')
gl.top_labels = False # 关闭顶端的经纬度标签
gl.right_labels = False # 关闭右侧的经纬度标签
gl.xformatter = LONGITUDE_FORMATTER # x轴设为经度的格式
gl.yformatter = LATITUDE_FORMATTER # y轴设为纬度的格式
# 设置经纬度网格的间隔及坐标标签
gl.xlocator = mticker.FixedLocator(np.arange(extent[0], extent[1]+1, 2))
gl.ylocator = mticker.FixedLocator(np.arange(extent[2], extent[3]+1, 2))
ax.set_xticks(list(range(extent[0], extent[1]+1, 2)), crs=ccrs.PlateCarree())
ax.set_yticks(list(range(extent[2], extent[3]+1, 2)), crs=ccrs.PlateCarree())
lon_formatter = LongitudeFormatter(zero_direction_label=True) # zero_direction_label用来设置经度的0度加不加E和W
lat_formatter = LatitudeFormatter()
ax.xaxis.set_major_formatter(lon_formatter)
ax.yaxis.set_major_formatter(lat_formatter)
return c1
SC = plot_cont(dt, [0.1, 0.1, 1, 1], t_level, [97, 109, 26, 35], show=True)
plt.text(104, 34.5, '审图号:GS(2024)0650号', fontsize=20)
###############################################################################
# 色标设置
rect = [0.3, 0.015, 0.6, 0.04]
cbar_ax = fig.add_axes(rect)
cb = fig.colorbar(SC, drawedges=True, ticks=t_level,
cax=cbar_ax, orientation='horizontal')
cb.set_label('℃', fontsize=12)
cb.ax.tick_params(length=0)
###############################################################################
# 存图
plt.savefig(figpath+'203 单省填充', dpi=600, bbox_inches = 'tight')
plt.show()