如何绘制WRF-Chem模拟前FINN排放wrffirechemi文件的PM2.5分布图？

不管用Python还是NCL，核心都是先读取wrffirechemi文件、确认PM2.5相关变量（或计算组分）、结合地图投影绘图，以下是具体实现：

Python 实现步骤 #

用xarray + matplotlib + cartopy组合处理，这是目前Python生态中处理WRF数据的常用方案：

• 导入依赖库：

import xarray as xr
import matplotlib.pyplot as plt
import cartopy.crs as ccrs

• 读取wrffirechemi文件：

ds = xr.open_dataset('wrffirechemi_d01_2020-03-15_10:00:00')

• 确认PM2.5相关变量：执行print(ds.variables)查看所有变量，FINN排放数据的wrffirechemi可能直接提供PM25_FINE变量，或者需要将OC、EC、SO4等组分相加得到PM2.5：

# 直接读取现成变量
pm25 = ds['PM25_FINE'][0]  # 取第0个时间步

# 若需组分计算（示例）
# pm25 = ds['OC'][0] + ds['EC'][0] + ds['SO4'][0] + ds['NO3'][0] + ds['NH4'][0]

• 绘制带地图的PM2.5分布图：

plt.figure(figsize=(10, 8))
ax = plt.axes(projection=ccrs.PlateCarree())
ax.coastlines(linewidth=0.8)
# 绘制数据，自动匹配WRF投影到地理投影
pm25.plot(ax=ax, transform=ccrs.PlateCarree(), cmap='viridis', cbar_kwargs={'label': 'PM2.5 Emission (μg/m²/s)'})
plt.title('WRF-Chem FINN Fire Emissions PM2.5\n2020-03-15 10:00')
plt.show()

• 若自动投影失效，可手动提取经纬度：

lat = ds['XLAT'][0].values
lon = ds['XLONG'][0].values
plt.pcolormesh(lon, lat, pm25.values, transform=ccrs.PlateCarree(), cmap='viridis')

NCL 实现步骤 #

NCL是WRF官方推荐的绘图工具，对WRF投影的适配更原生：

• 读取文件与变量：

begin
  f = addfile("wrffirechemi_d01_2020-03-15_10:00:00.nc", "r")
  lat = f->XLAT(0, :, :)  ; 提取第0个时间步的纬度
  lon = f->XLONG(0, :, :)  ; 提取第0个时间步的经度
  
  ; 获取PM2.5变量（或组分计算）
  pm25 = f->PM25_FINE(0, :, :)
  ; pm25 = f->OC(0,:,:) + f->EC(0,:,:) + f->SO4(0,:,:)  ; 组分计算示例

• 配置绘图参数并绘图：

wks = gsn_open_wks("png", "fire_pm25_plot")  ; 输出PNG文件
  gsn_define_colormap(wks, "viridis")  ; 设置配色
  
  res = True
  res@cnFillOn = True
  res@cnLinesOn = False
  res@cnLineLabelsOn = False
  res@lbLabelBarOn = True
  res@lbLabelFontHeightF = 0.015
  
  ; 匹配WRF的Lambert投影
  res@mpProjection = "LambertConformal"
  res@mpLambertParallel1 = f->CEN_LAT
  res@mpLambertParallel2 = f->CEN_LAT
  res@mpLambertMeridian = f->CEN_LON
  res@mpLimitMode = "LatLon"
  res@mpMinLonF = min(lon)
  res@mpMaxLonF = max(lon)
  res@mpMinLatF = min(lat)
  res@mpMaxLatF = max(lat)
  
  res@tiMainString = "WRF-Chem FINN Fire Emissions PM2.5\n2020-03-15 10:00"
  res@tiMainFontHeightF = 0.02
  
  plot = gsn_csm_contour_map(wks, pm25, res)
  draw(plot)
  frame(wks)
end

关键注意事项 #

• 确认变量名：部分wrffirechemi文件可能不直接提供PM2.5变量，需要根据FINN排放清单的组分定义自行计算，具体组分可参考FINN官方文档
• 时间维度：wrffirechemi变量的维度通常是[Time, south_north, west_east]，绘图时需指定具体时间步（一般取第0个）

内容的提问来源于stack exchange，提问作者Nafb