怎么用Python解析气象雷达基数据？pyart与wradlib对比

一、气象雷达基数据的组织特征与解析痛点

气象雷达基数据并不是普通的二维栅格图像，也不是常规的三维直角坐标矩阵。现代天气雷达在进行全方位观测时，天线按照预设的一系列仰角逐层向上旋转扫描，这一整套观测被称为体扫（Volume Scan）。

在每一个固定的扫描仰角（Sweep）下，雷达天线在 360 度圆周运动过程中不断向外发射微波脉冲并接收回波，沿不同方位角采集若干条径向射线（Rays）。每条射线上又按照固定的脉冲采样时间间隔，离散化为成百上千个距离库（Range Gates）。最终记录的反射率因子（Reflectivity，dBZ）、径向速度（Radial Velocity，m/s）以及双偏振物理量（差分反射率 ZDR、差分相位常数 KDP、互相关系数 CC），都依附在由仰角、方位角和距离库构成的极坐标体空间中。

图 1：气象雷达基数据多层级立体组织结构解析示意

直接解析基数据的难点在于雷达格式高度碎片化：国际通用的 NEXRAD Level II、欧洲气象常用的 ODIM HDF5、CF-Radial NetCDF，以及中国气象局新一代天气雷达常用的 CINRAD 格式（SA/SB/CB/CC/CD 等）。每种格式都有独特的二进制文件头、字节序编码、变长数据块和浮点标定系数（Scale & Offset）。直接编写底层解析器成本极高，因此在工程与科研中必须依赖 Py-ART 或 wradlib 这种经过严格检验的专业驱动库。

二、Python 解析雷达基数据的核心工作流

无论是使用 Py-ART 还是 wradlib，从雷达原始基数据到能够开展下游科研分析或制图输出，通用的技术管线都可以划分为四个阶段：

1. I/O 驱动层：探测文件格式头信息，加载雷达天线经纬度、站点海拔、扫描模式等元数据。
2. 空间几何重构：提取各个仰角的方位角与距离序列，将纯数据数组与空间球坐标对应起来。
3. 质量控制（QC）：剔除地物杂波（Ground Clutter）、退多普勒速度模糊（De-aliasing）、消除异常射束衰减。
4. 科学分析与栅格化：生成单仰角平面显示图（PPI）、垂直剖面图（RHI）或三维空间笛卡尔网格插值。

图 2：气象雷达基数据解析端到端处理管线图

三、使用 Py-ART 解析雷达基数据

Py-ART 是美国能源部 ARM（Atmospheric Radiation Measurement）项目开源的旗舰库。它的核心思想是高层面向对象封装。在 Py-ART 中，无论是何种雷达输入文件，只要调用统一的 pyart.io.read()，就会解析为一个标准化的 Radar 对象。

这个 Radar 对象自带所有物理场字段字典、仰角索引、扫描范围、地理信息及质量控制属性，非常适合大气气象领域的快速科研探索与制图。

import pyart
import matplotlib.pyplot as plt

# 1. 载入雷达基数据文件（自动探测并适配 NEXRAD、Cf/Radial 等多类驱动）
radar = pyart.io.read("radar_data_sample.nc")

# 2. 检查雷达基础元数据与搭载场变量
print(f"雷达站点名称: {radar.metadata.get('instrument_name', 'Radar_Station')}")
print(f"体扫仰角总层数: {radar.nsweeps}")
print(f"各仰角角度列表: {radar.fixed_angle['data']}")
print(f"包含的物理观测场: {list(radar.fields.keys())}")

# 3. 提取第 0 仰角（最底层）的反射率因子数据
sweep_slice = radar.get_slice(0)
reflectivity_field = radar.fields['reflectivity']['data'][sweep_slice]
azimuths = radar.azimuth['data'][sweep_slice]
ranges = radar.range['data']

print(f"仰角0包含的径向线数量: {len(azimuths)}, 单射线距离库数量: {len(ranges)}")

# 4. 利用内置专用绘图器快速渲染标准极坐标 PPI 图像
display = pyart.graph.RadarDisplay(radar)
fig, ax = plt.subplots(figsize=(8, 6))

display.plot('reflectivity', 0, vmin=-10, vmax=65, cmap='pyart_HomeyerRainbow', ax=ax)
display.plot_range_rings([50, 100, 150], ax=ax)  # 添加 50km, 100km, 150km 距离圈
display.set_limits(xlim=(-150, 150), ylim=(-150, 150), ax=ax)

plt.title("Py-ART PPI Reflectivity (Sweep 0)")
plt.savefig("pyart_ppi_demo.png", dpi=300, bbox_inches='tight')
plt.close()

Py-ART 的代码风格高度精炼。开发者无需手动计算极坐标到 Cartesian 直角坐标的转换公式，RadarDisplay 会自动提取雷达中心经纬度和极坐标库距完成绘图。

四、使用 wradlib 解析雷达基数据

wradlib 起源于德国波恩大学的水文气象团队。它的设计理念更偏向低耦合、模块化、强算法工具箱。早期的 wradlib 直接返回原始 NumPy 数组与元数据字典；最新版本则深度融合了 xarray 与 DataTree 体系。

wradlib 在欧洲标准的 ODIM HDF5、GRIB2 以及定制雷达格式支持上表现卓越，在定量降水估算（QPE）、衰减订正和与 GIS 空间图层互操作方面具备极其灵活的底层控制力。

import wradlib as wrl
import numpy as np
import matplotlib.pyplot as plt

# 1. 使用通用层级读取接口载入数据（以常用 HDF5/ODIM 格式为例）
filename = "radar_data_sample.h5"
raw_data, metadata = wrl.io.read_generic_hdf5(filename)

# 2. 提取指定扫描层（Sweep 0）的距离与方位角序列
nbins = metadata['dataset1/where']['nbins']
rscale = metadata['dataset1/where']['rscale']
r = np.arange(nbins) * rscale

nrays = metadata['dataset1/where']['nrays']
startaz = metadata['dataset1/where']['startaz']
azres = metadata['dataset1/where']['azres']
az = startaz + np.arange(nrays) * azres

# 3. 提取反射率因子数值并应用增益与偏移标定
raw_matrix = raw_data['dataset1/data1/data']
gain = metadata['dataset1/data1/what']['gain']
offset = metadata['dataset1/data1/what']['offset']
dbz = raw_matrix * gain + offset

# 4. 使用 wradlib 的 PPI 专业可视化接口绘制图表
fig, ax = plt.subplots(figsize=(8, 6))
plot_mesh = wrl.vis.plot_ppi(dbz, r=r, az=az, vmin=-10, vmax=65, cmap="turbo", ax=ax)
plt.title("wradlib PPI Reflectivity")
plt.colorbar(plot_mesh, ax=ax, label="Reflectivity (dBZ)")

plt.savefig("wradlib_ppi_demo.png", dpi=300, bbox_inches='tight')
plt.close()

从代码可以看出，wradlib 暴露出更多的底层参数（例如距离刻度 rscale、方位角步长 azres、物理增益 gain 等）。虽然代码行数略多，但这赋予了工程师精细调整极坐标参数和介入每一步数据变换的极大自由度。

五、Py-ART 与 wradlib 核心特性横向对比

为了让读者在实际工程立项时有清晰的研判依据，我们将 Py-ART 和 wradlib 在数据模型、支持格式、算法体系与地理处理等维度进行全面对比。

图 3：Py-ART 与 wradlib 核心架构与技术能力对比图

六、如何根据业务场景选择最合适的工具？

在实际项目落地时，选用 Py-ART 还是 wradlib 不应该凭个人习惯，而应该由具体需求指标决定：

图 4：气象雷达基数据解析技术选型决策图

特别说明：国内 CINRAD 雷达数据的处理建议

国内气象业务中广泛部署的是新一代天气雷达（CINRAD，如 SA、SB、CB、CC、SC 等型号）。这些雷达的原始基数据通常是按照中国气象局业务标准封装的二进制私有流。

针对 CINRAD 数据，通常的做法是：

• 使用专用转换桥梁：借助社区成熟的针对国内雷达的开源解析器（例如 PyCINRAD）直接读取二进制基数据。
• 转换为标准结构：将解析后的极坐标矩阵直接导出为标准 Py-ART Radar 对象或通用 Cf/Radial (NetCDF) 文件。
• 统一进行科学计算：一旦完成标准格式接入，即可顺畅调用 Py-ART 的多普勒退模糊与网格化，或者调用 wradlib 的降水反演算法。