前言
1.使用os库循环读取文件夹下的wrf数据，并用nc库的dataset读取,可使用wrf_list = [Dataset(f) for f in wrf_files] ，wrf_files是os读取形成的文件列表
2.使用wrfpython的getvar读取多个wrf文件的RAINC，RAINNC，RAINSH，利用cat将多时次数据合并
例如，RAINC = getvar(wrf_list, 'RAINC', timeidx=ALL_TIMES, method='cat')
3.利用total_rain.diff(dim='Time')语句计算逐小时的降水量。或者通过for循环计算然后将数组叠加也可。
关于三个降水变量的区别可以参考
WRF后处理：降雨量的说明以及降雨的绘制_wrf模拟降水量偏小-CSDN博客 https://blog.csdn.net/islandowner2017/article/details/119719854
另外，RAINC和RAINSH根据物理参数化方案的设置，可能为0。
# 导入数据读取模块 import numpy as np import pandas as pd from netCDF4 import Dataset import xarray as xr导入可视化模块
import matplotlib.pyplot as plt
from matplotlib.cm import get_cmap
import cartopy.crs as ccrs
import cartopy.feature as cfeature
import cartopy.io.shapereader as shpreader
from cartopy.mpl.ticker import LongitudeFormatter, LatitudeFormatter
from cartopy.mpl.gridliner import LATITUDE_FORMATTER, LONGITUDE_FORMATTER
import shapely.geometry as sgeom
import cmaps
简单计算与绘图(累计值)
from netCDF4 import Dataset
from wrf import getvarWRF数据目录
wrfout_path = '/home/mw/input/typhoon9537/'
读取 WRF 模拟数据
wrf_file = Dataset(wrfout_path + 'wrfout_d01_2019-08-08_20_00_00')
提取降雨量
RAINC = getvar(wrf_file, 'RAINC')
RAINNC = getvar(wrf_file, 'RAINNC')
RAINSH = getvar(wrf_file, 'RAINSH')计算累计降雨量
total_rain00 = RAINC + RAINSH + RAINNC
total_rain00.plot.contourf(cmap=cmaps.radar)
细化绘图（累计值）
from wrf import get_cartopy,cartopy_xlim,cartopy_ylim,to_np,latlon_coords读取 WRF 模拟数据
wrf_file = Dataset(wrfout_path + 'wrfout_d01_2019-08-09_06_00_00')
提取降雨量
RAINC = getvar(wrf_file, 'RAINC')
RAINNC = getvar(wrf_file, 'RAINNC')
RAINSH = getvar(wrf_file, 'RAINSH')计算累计降雨量
total_rain = RAINC + RAINSH + RAINNC
#投影与经纬度
wrf_proj = get_cartopy(RAINC)
lats, lons = latlon_coords(RAINC)通过等值线填充绘制降雨量分布
fig = plt.figure(figsize=(10,8))
设置地图投影
ax = plt.axes(projection=wrf_proj)
设置地图范围
ax.set_xlim(cartopy_xlim(RAINC))
ax.set_ylim(cartopy_ylim(RAINC))绘制降雨量分布（pcolormesh方法进行格点绘制）
#im = ax.pcolormesh(to_np(lons),
to_np(lats),
to_np(total_rain),
vmin=0,
vmax=200,
cmap=cmaps.WhiteBlueGreenYellowRed,
transform=ccrs.PlateCarree())
im = ax.contourf(to_np(lons),
to_np(lats),
to_np(total_rain),
levels=range(0, 201, 10),
cmap=cmaps.radar,
transform=ccrs.PlateCarree())为降雨量添加colorbar
cbar = plt.colorbar(im, ax=ax, extend='max', shrink=1)
cbar.set_label('Rainfall (mm)', fontdict={'size':20})
cbar.ax.tick_params(labelsize=20)添加经纬度网格线
ax.gridlines(color='black', linestyle='dotted')
设置标签大小
plt.tick_params(labelsize=15)
添加标题
plt.title('08-08 1800 - 08-09 0600(UTC)', loc='left', fontsize=20)
plt.show()
每小时降水量组图绘制
为了代码不繁琐直接利用xarray的plot作图，更多细致的作图敬请自己实现，以下示例小时降水量的组图绘制
此处使用了xarray的data.diff计算每小时的降水量
wrfout中的降水变量都是累计降水量，因此需要根据用后一时次减去前一时次才能得出这小时下了多少。
#小练习：绘制小时降雨量与累积降雨量（用组图形式展示）
import os
import numpy as np
import matplotlib.pyplot as pltfrom netCDF4 import Dataset
from wrf import getvar, ALL_TIMES定义 WRF 文件夹路径和文件名前缀
wrfout_path = "/home/mw/input/typhoon9537/"
filename_prefix = "wrfout_d01_"获取 WRF 文件列表，并按照文件名排序
wrf_files = sorted([os.path.join(wrfout_path, f) for f in os.listdir(wrfout_path) if f.startswith(filename_prefix)])
wrf_list = [Dataset(f) for f in wrf_files]提取降雨量
RAINC = getvar(wrf_list, 'RAINC', timeidx=ALL_TIMES, method='cat')
RAINNC = getvar(wrf_list, 'RAINNC', timeidx=ALL_TIMES, method='cat')
RAINSH = getvar(wrf_list, 'RAINSH', timeidx=ALL_TIMES, method='cat')计算累计降雨量
total_rain = RAINC + RAINSH + RAINNC
#计算逐小时降水量
hourly_rain = total_rain.diff(dim='Time')绘制逐小时降水图像
fig, axes = plt.subplots(nrows=4, ncols=3, figsize=(15, 12))
for i in range(12):
row, col = divmod(i, 3)hourly_rain[i].plot(ax=axes[row, col], cmap=cmaps.radar, add_colorbar=True,levels=np.arange(0,200,10))
fig.suptitle('Hourly Rainfall from 2019-08-08 19:00:00 to 2019-08-09 06:00:00')
plt.show()
累计降水量组图绘制
# 绘制累计降水图像
fig, axes = plt.subplots(nrows=4, ncols=3, figsize=(15, 12))
for i in range(12):
row, col = divmod(i, 3)total_rain[i].plot(ax=axes[row, col], cmap=cmaps.radar, add_colorbar=True,levels=np.arange(0,200,10))axes[row, col].set_title(f'Hour')
fig.suptitle('total_rain from 2019-08-08 19:00:00 to 2019-08-09 06:00:00')
plt.show()
谢谢观看
