Damu的个人博客 - 技术分享与学习笔记

WRF输出数据降雨量来自于三个变量：

RAINC: ACCUMULATED TOTAL CUMULUS PRECIPITATION

RAINNC: ACCUMULATED TOTAL GRID SCALE PRECIPITATION

RAINSH: ACCUMULATED SHALLOW CUMULUS PRECIPITATION

RAINC: 积云深对流过程产生的累积降水量，也就是模式中的积云对流参数化方案导致的降雨( cu_physics)。对于高分辨率的模拟，比如dx<5km，通常会将积云对流参数化方案关闭，此时RAINC为0。
RANNC: 此类降雨来源于云微物理参数化方案(mp_physics)，如大尺度抬升过程产生的凝结等微物理过程降水，也就是非对流产生的降水。
RAINSH: 积云对流参数化方案主要是反映深对流的降水过程，但是一些积云对流参数化方案，能够支持浅对流导致的降水，此时总降水还需要加上RAINSH。WRF中支持浅对流的参数化方案(cu_physics)有以下几种：KF，SAS，G3，BMJ，Tiedtke。WRF中也有独立于深对流过程的浅对流方案，通过namelist中设置shcu_physics。一般情况下，浅对流产生的降水量较小。

此外，固态降水例如雪(SNOWC/SNOWNC)、霰(GRAUPELC/GRAUPELNC)等降水，它们是降水的不同相态，已经都包含在RAINC/RAINNC中，不需要额外添加。 总结为一句话：**总降水=RAINNC + RAINC + RAINSH**，其中RAINC和RAINSH根据物理参数化方案的设置可能为0。 需要注意的是，降雨是个**累积量**，而不是瞬时值，在wrfout中如果需要得到模拟起止日期内的降雨，需要将最后一个时刻的降雨量减去第一个时刻的降雨量.

import numpy as np
from netCDF4 import Dataset
from wrf import (getvar, to_np) 
import matplotlib.pyplot as plt
import cartopy.crs as ccrs
from cartopy.mpl.gridliner import LONGITUDE_FORMATTER, LATITUDE_FORMATTER
from cartopy.io.shapereader import BasicReader
import matplotlib.ticker as mticker
# 读取数据
f='E:\BaiduSyncdisk\Code\python\WRF\Data\wrfout_d01_2020-06-11_06_00_00'
ncfile = Dataset(f)
R = to_np(getvar(ncfile,"RAINC"))+to_np(getvar(ncfile,"RAINNC"))+to_np(getvar(ncfile,"RAINSH"))
lon = np.array(ncfile['XLONG'])[0,:,:]
lat = np.array(ncfile['XLAT'])[0,:,:]

# shape文件读取
countries=BasicReader(r'E:\BaiduSyncdisk\Code\Data\shape数据\国界\国家矢量.shp')
ten=BasicReader(r'E:\BaiduSyncdisk\Code\Data\shape数据\九段线\九段线.shp')
province=BasicReader(r'E:\BaiduSyncdisk\Code\Data\shape数据\2023年省级\2023年省级.shp')


# 颜色条设置
rgb = ([237, 237, 237],[209, 209, 209],[173, 173, 173],[131, 131, 131],
       [93 , 93 , 93],[151, 198, 223],[111, 176, 214],[49 , 129, 189],
       [26 , 104, 174],[8  , 79 , 153],[62 , 168, 91],[110, 193, 115],
       [154, 214, 149],[192, 230, 185],[223, 242, 217],[255, 255, 164],
       [255, 243, 0],[255, 183, 0],[255, 123, 0],[255, 62 , 0],
       [255, 2  , 0],[196, 0  , 0],[136, 0  , 0],
      )
clors = np.array(rgb)/255.
clevs = [0.1, 1,  2,  5, 7.5, 10, 13, 16, 20, 25, 30, 35, 40, 50, 60,
         70, 80, 90, 100, 125, 150, 175, 200, 250
        ]

# 画图
fig = plt.figure(figsize=(10, 10))
ax = fig.add_subplot(111, projection=ccrs.PlateCarree())
contour = ax.contourf(lon,lat, R, clevs, colors=clors, transform=ccrs.PlateCarree())
cb_cf = fig.colorbar(contour, ax=ax,ticks=clevs, shrink=.65)

ax.add_geometries(countries.geometries(),linewidth=.5,edgecolor='black',crs=ccrs.PlateCarree(),facecolor='none') 
ax.add_geometries(ten.geometries(),linewidth=.5,edgecolor='black',crs=ccrs.PlateCarree(),facecolor='none') 
ax.add_geometries(province.geometries(),linewidth=.5,edgecolor='black',crs=ccrs.PlateCarree(),facecolor='none') 
ax.set_extent([60,140,10,60],crs=ccrs.PlateCarree()) # 小范围


g1=ax.gridlines(draw_labels=True,linewidth=1,color='none',alpha=0.5,linestyle='--',x_inline=False,y_inline=False)
g1.top_labels=False
g1.right_labels=False
g1.xformatter=LONGITUDE_FORMATTER
g1.yformatter=LATITUDE_FORMATTER
g1.rotate_labels=False
g1.xlocator = mticker.FixedLocator(np.arange(60,140,10))
g1.ylocator = mticker.FixedLocator(np.arange(10,60,5))
g1.xlabel_style={'size':17}
g1.ylabel_style={'size':17}

参考资料：

https://blog.csdn.net/islandowner2017/article/details/119719854

https://www.heywhale.com/mw/project/6284b478901f1de4da8c8169