# Python+Cartopy实战：用WRFout数据绘制专业级台风路径与降水分布图（附完整代码） 如果你正在处理WRF模式输出的数据，尤其是像台风“贝碧嘉”这类复杂的天气系统，那么你很可能已经体会过那种面对NetCDF文件时无从下手的迷茫。气象数据本身维度多、体量大，而要将这些冰冷的格点数据转化为一张张能清晰讲述天气故事的图表，更是需要跨越数据提取、坐标转换、地图投影和美学设计等多重关卡。我最初接触WRFout文件时，也曾在各种论坛和代码片段中摸索，试图拼凑出一个能用的绘图脚本，结果往往是图表勉强能看，但距离发表级的标准相去甚远，更别提代码的复用性和可读性了。 这篇文章正是为了解决这些痛点而生。我们不打算重复那些基础的“Hello World”式教程，而是直接切入核心，分享一套经过实战检验的、模块化的Python绘图方案。我们将以2024年强台风“贝碧嘉”的模拟数据为例，手把手带你使用Cartopy和Matplotlib，从原始的`wrfout_d01`文件开始，一步步绘制出专业级的**500hPa位势高度叠加风矢量时序组图**和**累计降水热力图**。更重要的是，我会把过程中踩过的坑、调试过的参数以及提升图表专业度的技巧毫无保留地分享出来，并提供**可直接嵌入你项目**的完整函数模块。无论你是气象专业的研究生，还是从事气候或环境数据分析的开发者，这套代码都能让你在数据可视化环节节省大量时间，把精力更多地投入到科学问题本身。 ## 1. 环境准备与核心库深度解析 在开始绘图之前，搭建一个稳定且高效的环境至关重要。不同于简单的`pip install`，针对WRF数据处理和科学绘图，我们需要对几个核心库的版本和协作方式有更清晰的认识。 ### 1.1 创建专属虚拟环境与库安装 我强烈建议为气象数据处理项目创建独立的虚拟环境，避免与其他项目的库版本冲突。使用conda管理环境能很好地处理一些具有非Python依赖（如NetCDF库）的包。 ```bash # 创建并激活一个名为wrf_plot的新环境，指定Python版本 conda create -n wrf_plot python=3.9 conda activate wrf_plot # 通过conda安装具有复杂依赖的核心科学计算库 conda install -c conda-forge numpy scipy xarray pandas # 安装NetCDF文件处理及WRF数据读取的关键库 conda install -c conda-forge netcdf4 wrf-python # 安装绘图核心库：Matplotlib和Cartopy。Cartopy的安装务必通过conda-forge conda install -c conda-forge matplotlib cartopy ``` > **注意**：`wrf-python`是NCAR官方维护的用于便捷读取和计算WRF数据的工具包，它封装了许多复杂的插值和计算函数（如`getvar`, `interplevel`），能极大简化我们的代码。`cartopy`用于地理信息可视化，其地图投影和海岸线数据通过conda-forge渠道安装最为可靠。 ### 1.2 关键库的角色与版本协同 仅仅安装成功还不够，理解这些库在流程中的作用，才能在使用时得心应手。下面这个表格梳理了它们的分工： | 库名称 | 核心用途 | 在本项目中的关键功能/类 | 推荐版本 | | :--- | :--- | :--- | :--- | | **netCDF4** | 读取/写入NetCDF格式文件 | `Dataset` 类，用于打开`wrfout`文件 | >=1.6.0 | | **wrf-python** | 简化WRF数据提取与计算 | `getvar`, `interplevel`, `to_np` | >=1.3.2 | | **xarray** | 处理带标签的多维数组（可选但推荐） | 高级数据筛选、切片、运算 | >=0.20.0 | | **cartopy** | 地理空间数据可视化 | `ccrs.PlateCarree`, `cfeat.COASTLINE` | >=0.20.0 | | **matplotlib** | 基础绘图与图形定制 | `plt.subplots`, `colors.ListedColormap` | >=3.5.0 | **版本协同要点**：`cartopy` 0.18.0之后，其默认的地图数据源有所变化，可能导致之前代码中的海岸线加载失败。确保使用较新版本并配合`conda-forge`的更新，可以避免许多奇怪的问题。`wrf-python`的API相对稳定，但建议关注其更新日志，以利用性能优化和新功能。 ### 1.3 数据准备：理解你的WRFout文件 在运行代码前，你需要将WRF模式输出的NetCDF文件（通常命名为`wrfout_d0*_*`）放在项目目录中。使用`ncdump -h`命令可以快速查看文件结构，但用Python查看更直接： ```python # 快速查看wrfout文件变量列表 from netCDF4 import Dataset ncfile = Dataset('your_wrfout_file.nc') print(ncfile.variables.keys()) # 查看所有变量名 print(ncfile.variables['Times'][:]) # 查看时间维度的具体值 ncfile.close() ``` 通常我们需要关注的变量包括： - **`Times`**: 时间维度（字符串格式） - **`XLAT`, `XLONG`**: 网格点的纬度和经度（二维数组） - **`z`**: 位势高度 - **`ua`, `va`**: U/V方向的风速分量 - **`pressure`**: 气压层 - **`RAINC`, `RAINNC`**: 对流性和非对流性累积降水 ## 2. 500hPa位势高度与风场时序组图绘制 绘制高空形势场是分析天气系统动力结构的基础。我们将创建一个函数，一次性生成多个时次的子图，便于对比系统随时间演变。 ### 2.1 函数框架与数据提取 首先，我们构建一个名为`plot_500hpa_timeseries`的函数。它的核心任务是：**循环读取指定时次的数据，在统一的地图投影和色标下，绘制位势高度填色图、等高线并叠加风矢量**。 ```python import numpy as np import matplotlib.pyplot as plt import cartopy.crs as ccrs import cartopy.feature as cfeat from cartopy.mpl.ticker import LongitudeFormatter, LatitudeFormatter from matplotlib.colors import ListedColormap from netCDF4 import Dataset from wrf import getvar, interplevel, to_np def plot_500hpa_timeseries(ncfile_path, time_indices, output_filename='500hpa_wind_series.png'): """ 绘制指定时次的500hPa位势高度与风场组图。 参数 ---------- ncfile_path : str WRFout NetCDF文件路径。 time_indices : list of int 需要绘图的时次索引列表，例如 [0, 3, 6, 9]。 output_filename : str, 可选 输出图片的文件名。 """ # 打开文件 ncfile = Dataset(ncfile_path) times = ncfile.variables['Times'][:] # 所有时间步 # 预设地图范围（根据你的模拟区域调整） map_extent = [110, 140, 20, 37] # [lon_min, lon_max, lat_min, lat_max] # 创建图形和子图网格 n_times = len(time_indices) # 设计布局：如果时次多，可以考虑2行N列 ncols = min(3, n_times) # 每行最多3列 nrows = (n_times + ncols - 1) // ncols # 计算所需行数 fig, axes = plt.subplots(nrows=nrows, ncols=ncols, subplot_kw={'projection': ccrs.PlateCarree()}, figsize=(6*ncols, 4.5*nrows), dpi=300) # 如果只有一个子图，确保axes是数组形式以便统一循环 if n_times == 1: axes = np.array([axes]) axes = axes.flatten() # --- 为所有时次计算统一的色标范围，确保颜色对比一致 --- all_ht_500 = [] for idx in time_indices: ht = getvar(ncfile, "z", timeidx=idx, units="dm") # 单位：位势什米 p = getvar(ncfile, "pressure", timeidx=idx) ht_500 = to_np(interplevel(ht, p, 500)) # 插值到500hPa all_ht_500.append(ht_500) all_data = np.concatenate([arr.flatten() for arr in all_ht_500]) vmin, vmax = np.percentile(all_data, [2, 98]) # 使用2%和98%分位数，避免极端值影响色标 # 定义专业气象配色（类似NCL的‘BlAqGrYeOrReVi200’） colors = ['#2b83ba', '#64abb0', '#9cb886', '#c7c266', '#e4a940', '#f48d43', '#d35f4d'] custom_cmap = ListedColormap(colors) levels = np.arange(np.floor(vmin/4)*4, np.ceil(vmax/4)*4 + 1, 4) # 生成以4为间隔的等值线 ``` **代码解析**： - **统一色标**：在循环绘图前，先预读所有指定时次500hPa位势高度数据，计算全局的`vmin`和`vmax`。这是制作专业组图的关键，它保证了不同时次之间的颜色代表相同的数值范围，便于对比。 - **动态布局**：根据输入的时次数量，自动计算子图的行列布局，使函数更具通用性。 - **配色选择**：我们使用了从蓝到红的渐变色系，这是位势高度场的常用配色方案，能清晰反映高压（暖色）和低压（冷色）系统。 ### 2.2 循环绘图与地图要素叠加 接下来，在循环中为每个时次绘制子图。 ```python # 循环绘制每个时次的子图 for ax_idx, (ax, time_idx) in enumerate(zip(axes[:n_times], time_indices)): # 1. 提取当前时次数据 time_str = ''.join([x.decode('utf-8') for x in times[time_idx]]) # 解码时间字符串 ht = getvar(ncfile, "z", timeidx=time_idx, units="dm") u = getvar(ncfile, "ua", timeidx=time_idx) v = getvar(ncfile, "va", timeidx=time_idx) p = getvar(ncfile, "pressure", timeidx=time_idx) ht_500 = to_np(interplevel(ht, p, 500)) u_500 = to_np(interplevel(u, p, 500)) v_500 = to_np(interplevel(v, p, 500)) lats = to_np(getvar(ncfile, "XLAT", timeidx=time_idx)) lons = to_np(getvar(ncfile, "XLONG", timeidx=time_idx)) # 2. 设置地图范围 ax.set_extent(map_extent, crs=ccrs.PlateCarree()) # 3. 绘制位势高度填色图 cf = ax.contourf(lons, lats, ht_500, levels=levels, cmap=custom_cmap, vmin=vmin, vmax=vmax, extend='both', transform=ccrs.PlateCarree()) # 4. 叠加等高线（黑色细线） cs = ax.contour(lons, lats, ht_500, levels=levels, colors='k', linewidths=0.6, transform=ccrs.PlateCarree()) ax.clabel(cs, inline=True, fontsize=8, fmt='%1.0f') # 等高线标注 # 5. 叠加风矢量（适当稀疏化，避免过于密集） skip = max(1, lons.shape[0] // 25) # 每25个网格点取一个风矢量 qv = ax.quiver(lons[::skip, ::skip], lats[::skip, ::skip], u_500[::skip, ::skip], v_500[::skip, ::skip], transform=ccrs.PlateCarree(), scale=400, width=0.0025, headwidth=3, color='#333333') # 6. 添加地理要素 ax.add_feature(cfeat.COASTLINE.with_scale('50m'), linewidth=0.8) ax.add_feature(cfeat.BORDERS.with_scale('50m'), linestyle=':', linewidth=0.5) # 可以在此处添加省界、河流等，需要相应的shapefile # 7. 设置坐标轴和标题 ax.set_xticks(np.arange(map_extent[0], map_extent[1]+1, 5)) ax.set_yticks(np.arange(map_extent[2], map_extent[3]+1, 2)) ax.xaxis.set_major_formatter(LongitudeFormatter(zero_direction_label=True)) ax.yaxis.set_major_formatter(LatitudeFormatter()) ax.tick_params(labelsize=9) ax.set_title(f'Time: {time_str}', fontsize=11, pad=6) # 8. 仅在第一个子图添加风矢量图例 if ax_idx == 0: # 将图例放在子图内的一个固定位置（例如右上角） ax.quiverkey(qv, 0.85, 0.95, 20, '20 m/s', labelpos='E', coordinates='axes', fontproperties={'size': 9}) # 隐藏多余的子图（如果时次数不是行列数的整数倍） for ax in axes[n_times:]: ax.axis('off') ``` **关键技巧**： - **风矢量稀疏化**：WRF模式分辨率通常很高，直接绘制所有格点的风矢量会导致图形杂乱。通过`skip`参数进行等间隔采样，能在保留风场整体特征的同时保持图面清晰。 - **等高线标注**：`clabel`函数自动在等高线上添加数值标注。`fmt='%1.0f'`表示格式化为整数，这对于位势高度（单位：位势什米）是合适的。 - **地图细节**：使用`cfeat.COASTLINE.with_scale('50m')`指定中等分辨率的海岸线数据，在保证清晰度的同时减少图形文件大小。 ### 2.3 添加全局色标与保存输出 最后，为整个图组添加一个共用的色标，并调整布局保存图片。 ```python # 添加全局色标 # 调整色标位置，使其位于所有子图的右侧 cbar_ax = fig.add_axes([0.92, 0.15, 0.015, 0.7]) # [左, 下, 宽, 高] cbar = fig.colorbar(cf, cax=cbar_ax, orientation='vertical') cbar.set_label('500hPa Geopotential Height (dam)', fontsize=12) cbar.ax.tick_params(labelsize=10) # 添加总标题 fig.suptitle('WRF Simulation: 500hPa Geopotential Height & Wind Field', fontsize=14, y=0.98) # 调整子图间距并保存 plt.tight_layout(rect=[0, 0, 0.9, 0.96]) # 为色标和总标题留出空间 plt.savefig(output_filename, dpi=300, bbox_inches='tight') plt.close(fig) ncfile.close() print(f"图表已保存至: {output_filename}") ``` 现在，你只需要调用这个函数，传入文件路径和时次列表，即可一键生成专业的组图： ```python # 示例：绘制第0, 6, 12, 18时次（假设时间间隔为6小时） plot_500hpa_timeseries('wrfout_d01_2024-09-15_12_00_00', time_indices=[0, 1, 2, 3], output_filename='typhoon_500hpa_evolution.png') ``` ## 3. 累计降水热力图绘制与高级配色方案 降水是台风影响评估的核心要素。WRF模式输出的`RAINC`（对流性降水）和`RAINNC`（非对流性/网格尺度降水）是累积量，我们需要将其合并并绘制成空间分布图。 ### 3.1 计算累计降水与数据维度对齐 首先创建一个专门绘制降水的函数。这里有一个常见的“坑”：WRF中降水变量的维度有时可能与经纬度变量顺序不一致，需要进行转置。 ```python def plot_accumulated_precipitation(ncfile_path, time_idx, output_filename='accumulated_rain.png'): """ 绘制指定时次的累计降水空间分布图。 参数 ---------- ncfile_path : str WRFout NetCDF文件路径。 time_idx : int 需要绘图的时次索引。 output_filename : str, 可选 输出图片的文件名。 """ ncfile = Dataset(ncfile_path) times = ncfile.variables['Times'][:] # 获取当前时次的时间字符串 time_str = ''.join([x.decode('utf-8') for x in times[time_idx]]) # 提取经纬度 lats = to_np(getvar(ncfile, "XLAT", timeidx=time_idx)) lons = to_np(getvar(ncfile, "XLONG", timeidx=time_idx)) # 计算累计降水：RAINC + RAINNC （注意：RAINSH在某些参数化方案下可能为0） rainc = to_np(getvar(ncfile, "RAINC", timeidx=time_idx)) rainnc = to_np(getvar(ncfile, "RAINNC", timeidx=time_idx)) total_rain = rainc + rainnc # **关键步骤：检查并处理维度不匹配** if lats.shape != total_rain.shape: # 常见情况：经纬度是二维数组 (y, x)，但降水数据可能是 (x, y) if lats.shape == total_rain.shape[::-1]: total_rain = total_rain.T # 进行转置 print(f"提示：降水数据维度已自动转置以匹配经纬度。") else: raise ValueError(f"数据维度不匹配且无法自动处理。纬度: {lats.shape}, 降水: {total_rain.shape}") # 设置地图范围（可与高度场图保持一致，也可根据降水范围调整） map_extent = [110, 140, 20, 37] ``` ### 3.2 设计专业降水配色与分级 降水图的配色直接影响到对雨带强度和分布的解读。气象上常用从蓝到红，再到紫的色系表示降水量级递增。 ```python # 定义专业降水配色（24个颜色，对应25个分级界限） precip_colors = [ '#FFFFFF', '#C6E1FF', '#8CB3FF', '#528BFF', '#0064FF', # 0-10mm: 白到蓝 '#00FFFF', '#00E1C6', '#00C38B', '#00A550', '#008000', # 10-50mm: 青到绿 '#FFFF00', '#E1E100', '#C3C300', '#A5A500', '#878700', # 50-100mm: 黄 '#FFA500', '#FF8000', '#FF5A00', '#FF3500', '#FF0000', # 100-200mm: 橙到红 '#D900D9', '#B300B3', '#8C008C', '#660066', '#400040' # >200mm: 紫 ] # 定义降水等级界限（单位：mm） levels = [0.1, 2, 5, 10, 15, 20, 25, 30, 40, 50, 75, 100, 150, 200, 300, 400] # 创建自定义颜色映射和标准化对象 cmap = ListedColormap(precip_colors[:len(levels)]) # 颜色数比等级数多1 norm = BoundaryNorm(levels, cmap.N) # 创建图形 fig = plt.figure(figsize=(12, 8), dpi=300) ax = plt.subplot(projection=ccrs.PlateCarree()) ax.set_extent(map_extent, crs=ccrs.PlateCarree()) ``` **配色设计思路**： - **小雨（<10mm）**：使用白色到蓝色的柔和渐变，避免在无雨或小雨区域产生视觉压迫感。 - **中雨（10-50mm）**：过渡到青色和绿色，自然醒目。 - **大雨到暴雨（50-200mm）**：采用黄色、橙色到红色的暖色调，警示性强。 - **大暴雨以上（>200mm）**：使用深紫色，突出极端降水区域。 ### 3.3 绘制降水填色与叠加关键等值线 使用`pcolormesh`绘制填色图，并叠加主要量级的等值线以突出雨带轮廓。 ```python # 使用pcolormesh绘制降水填色（比contourf边界更清晰） rain_plot = ax.pcolormesh(lons, lats, total_rain, cmap=cmap, norm=norm, shading='auto', # 自动选择 shading 方式（‘nearest’或‘gouraud’） transform=ccrs.PlateCarree()) # 叠加重要量级的等值线（如50mm, 100mm, 200mm），避免图面过于杂乱 heavy_rain_levels = [50, 100, 200] cs = ax.contour(lons, lats, total_rain, levels=heavy_rain_levels, colors='k', linewidths=1.5, transform=ccrs.PlateCarree()) ax.clabel(cs, inline=True, fontsize=10, fmt='%1.0f mm') # 添加地理要素 ax.add_feature(cfeat.COASTLINE.with_scale('50m'), linewidth=1.2) ax.add_feature(cfeat.BORDERS.with_scale('50m'), linestyle='-', linewidth=0.8, alpha=0.7) # 可在此处加载河流、湖泊等自然地理要素 # 设置美观的坐标轴 ax.set_xticks(np.arange(map_extent[0], map_extent[1]+1, 5)) ax.set_yticks(np.arange(map_extent[2], map_extent[3]+1, 2)) ax.xaxis.set_major_formatter(LongitudeFormatter(zero_direction_label=True)) ax.yaxis.set_major_formatter(LatitudeFormatter()) ax.tick_params(labelsize=11) # 添加色标 cbar = fig.colorbar(rain_plot, ax=ax, orientation='horizontal', fraction=0.05, pad=0.08, aspect=30) cbar.set_label('Accumulated Precipitation (mm)', fontsize=13) cbar.ax.tick_params(labelsize=11) # 设置色标刻度为定义的等级 cbar.set_ticks(levels) cbar.set_ticklabels([str(l) for l in levels]) ``` ### 3.4 添加自定义标题与台风路径标记（进阶） 为了使图表信息更完整，我们可以添加模拟的台风路径。这需要从多个时次中提取台风中心位置（例如最小海平面气压点）。 ```python # --- 进阶：叠加台风路径（需要循环读取多个时次）--- # 假设我们已经有一个包含各时次台风中心经纬度的列表：track_lons, track_lats # track_lons = [125.1, 124.8, 124.3, ...] # track_lats = [22.5, 23.1, 23.8, ...] # 绘制台风路径线 # ax.plot(track_lons, track_lats, color='magenta', linewidth=2.5, # transform=ccrs.PlateCarree(), marker='o', markersize=4, # label='Simulated Typhoon Track') # 在路径上标记特定时次（如登陆时刻） # ax.scatter(landfall_lon, landfall_lat, s=120, color='red', # edgecolors='k', linewidth=1.5, zorder=5, # transform=ccrs.PlateCarree(), label='Landfall') # 添加图例（如果绘制了路径） # ax.legend(loc='upper right', fontsize=10, framealpha=0.9) # 添加主标题和副标题 forecast_hour = int(getvar(ncfile, 'XTIME', timeidx=time_idx) / 60) # 计算预报小时数 plt.title(f'WRF Simulated Accumulated Precipitation: {forecast_hour}h Forecast\n' f'Valid Time: {time_str}', fontsize=14, pad=12) # 保存图形 plt.tight_layout() plt.savefig(output_filename, dpi=300, bbox_inches='tight') plt.close(fig) ncfile.close() print(f"降水图已保存至: {output_filename}") ``` 调用这个函数，你就能得到一张信息丰富的累计降水分布图： ```python # 绘制第24个时次（假设是模拟结束时刻）的累计降水 plot_accumulated_precipitation('wrfout_d01_2024-09-15_12_00_00', time_idx=24, output_filename='typhoon_total_rainfall.png') ``` ## 4. 代码整合、优化与批量处理实战 将上述函数整合到一个模块中，并添加一些实用功能，如批量处理多个文件、生成动画等，可以极大提升工作效率。 ### 4.1 创建可复用的绘图模块 建议创建一个名为`wrf_plotter.py`的Python模块，将上述函数封装起来，并添加一些辅助函数。 ```python """ wrf_plotter.py 专业WRFout数据可视化模块。 包含绘制位势高度场、风场、降水场等函数。 """ import numpy as np import matplotlib.pyplot as plt import cartopy.crs as ccrs import cartopy.feature as cfeat from cartopy.mpl.ticker import LongitudeFormatter, LatitudeFormatter from matplotlib.colors import ListedColormap, BoundaryNorm from netCDF4 import Dataset from wrf import getvar, interplevel, to_np import os # 设置全局绘图样式（可选） plt.style.use('seaborn-v0_8-whitegrid') # 使用seaborn的白色网格样式，更清爽 plt.rcParams['font.family'] = 'DejaVu Sans' # 确保字体支持中文（如有需要） plt.rcParams['axes.titlesize'] = 14 plt.rcParams['axes.labelsize'] = 12 # 在此处粘贴之前定义的 plot_500hpa_timeseries 和 plot_accumulated_precipitation 函数 # ... def plot_multiple_files_rainfall(file_pattern, output_dir='./rainfall_plots/'): """ 批量处理多个wrfout文件，为每个文件生成累计降水图。 适用于分析不同参数化方案或不同起报时间的模拟结果。 参数 ---------- file_pattern : str 用于匹配文件的glob模式，例如 'wrfout_d01_*'。 output_dir : str 输出图片的目录。 """ import glob if not os.path.exists(output_dir): os.makedirs(output_dir) file_list = sorted(glob.glob(file_pattern)) print(f"找到 {len(file_list)} 个文件待处理。") for file_path in file_list: # 从文件名中提取标识信息（如日期时间） basename = os.path.basename(file_path) # 假设文件名格式为 wrfout_d01_2024-09-15_12_00_00 time_str_from_file = basename.split('_')[2] + '_' + basename.split('_')[3].replace('-', ':') # 为每个文件绘制最后一个时次的降水 ncfile = Dataset(file_path) last_time_idx = len(ncfile.variables['Times'][:]) - 1 ncfile.close() output_name = os.path.join(output_dir, f'rainfall_{time_str_from_file}.png') try: plot_accumulated_precipitation(file_path, last_time_idx, output_name) print(f"成功处理: {basename}") except Exception as e: print(f"处理 {basename} 时出错: {e}") ``` ### 4.2 性能优化与常见问题排查 处理高分辨率WRF数据时，可能会遇到内存占用大、绘图速度慢的问题。以下是一些优化建议： 1. **数据切片**：如果区域很大，可以只提取感兴趣的子区域进行绘图，减少数据量。 ```python # 在提取变量后，可以对数组进行切片 lat_slice = slice(50, 150) # 行索引范围 lon_slice = slice(30, 130) # 列索引范围 lats_sub = lats[lat_slice, lon_slice] lons_sub = lons[lat_slice, lon_slice] data_sub = total_rain[lat_slice, lon_slice] ``` 2. **使用更快的绘图方法**：对于非常大的网格，`pcolormesh`可能比`contourf`更快，但`contourf`在平滑性上更好。可以尝试调整`shading`参数。 3. **关闭交互模式**：在脚本运行时，使用`plt.ioff()`关闭交互模式，并在最后统一`plt.show()`或`plt.savefig()`，可以避免内存累积。 **常见问题排查清单**： | 问题现象 | 可能原因 | 解决方案 | | :--- | :--- | :--- | | 地图空白或变形 | 地图范围设置错误或投影不匹配 | 检查`set_extent`参数单位是否为度，确保与`PlateCarree`投影匹配。 | | 风矢量方向/大小异常 | `ua`, `va`单位可能不是m/s，或维度顺序错误 | 使用`wrf.getvar`时确认单位，用`to_np`检查数组形状。 | | 色标显示异常（全图一个颜色） | 数据值范围异常（如全为0或包含NaN） | 打印`data.min(), data.max()`检查，使用`np.nan_to_num`处理NaN。 | | 绘图速度极慢 | 网格分辨率过高，或同时打开过多图形窗口 | 对风矢量进行更大幅度的稀疏化(`skip`)，使用子图时关闭未使用的轴。 | | 中文标签显示为方框 | 系统缺少中文字体 | 指定一个系统中存在的中文字体，或使用英文标签。 | ### 4.3 进阶应用：制作降水演变动画 利用`matplotlib.animation`模块，你可以轻松地将一系列时次的降水图制作成GIF或MP4动画，直观展示雨带的移动和发展。 ```python # 示例：创建降水演变动画（需单独保存每个时次的图或直接在内存中处理） import matplotlib.animation as animation from PIL import Image # 用于后期处理GIF def create_rainfall_animation(ncfile_path, output_gif='rainfall_evolution.gif'): """创建累计降水随时间演变的动画。""" ncfile = Dataset(ncfile_path) n_times = len(ncfile.variables['Times'][:]) fig, ax = plt.subplots(subplot_kw={'projection': ccrs.PlateCarree()}, figsize=(10, 7)) map_extent = [110, 140, 20, 37] ax.set_extent(map_extent) # 初始化一个空的绘图对象 im = ax.pcolormesh([], [], [], cmap=plt.cm.viridis) # 先用一个临时colormap ax.add_feature(cfeat.COASTLINE, linewidth=1.0) def update(frame): """更新函数，为每一帧绘制数据。""" ax.clear() ax.set_extent(map_extent) # 提取第frame时次的数据并绘图（此处省略具体提取代码） # lats, lons, total_rain = extract_data_for_frame(ncfile, frame) # rain_plot = ax.pcolormesh(lons, lats, total_rain, ...) # ax.set_title(f'Time: {time_str}') return [im] # 返回更新后的图形对象列表 ani = animation.FuncAnimation(fig, update, frames=n_times, interval=200, blit=False) ani.save(output_gif, writer='pillow', fps=5) print(f"动画已保存为: {output_gif}") ncfile.close() ``` 将这套代码应用到你的WRF后处理流程中，无论是用于论文图表制作、项目报告还是业务产品展示，都能显著提升成果的可视化质量和你的工作效率。记住，好的数据可视化本身就是一种强大的分析工具和沟通语言。