3. 雷达数据预处理：数据格式解析、质量控制、插值与网格化
各位好，我是老张。今天咱们聊聊雷达数据预处理。说实话，这步是整条链路里最枯燥、但也是最容易翻车的地方。我见过太多项目，模型跑得飞起，结果一查原始数据全是噪声——嗯，那感觉就像你辛辛苦苦盖了栋楼，最后发现地基是豆腐渣。
雷达数据从天线接收下来，到能喂给风电功率预测模型，中间要过三关：格式解析、质量控制、插值网格化。咱们一个一个说。
3.1 数据格式解析：别让格式卡住你
雷达厂商各有各的脾气。有的用二进制，有的用HDF5，还有的用NetCDF。我个人习惯，拿到数据第一件事就是看头文件。头文件里藏着雷达型号、扫描模式、时间戳、仰角信息——这些是后续处理的命根子。
常见雷达数据格式一览
| 格式 | 特点 | 常见厂商 |
|---|---|---|
| 二进制（RAW） | 体积小，解析麻烦 | 国产雷达居多 |
| HDF5 | 自描述，带元数据 | Vaisala、EEC |
| NetCDF | 气象行业标准 | 欧洲雷达网 |
| BUFR | WMO标准，压缩率高 | 气象局交换数据 |
举个例子，我处理过一批某国产雷达的二进制数据。头文件里仰角信息居然藏在第37-40字节，用16位整数存着，单位是0.1度。当时我对着文档看了半天才反应过来——说白了，就是厂商没按套路出牌。
# Python示例：解析二进制雷达头文件
import struct
def parse_radar_header(raw_bytes):
# 假设头文件长度128字节
header = raw_bytes[:128]
# 解析仰角（第37-40字节，16位整数，单位0.1度）
elevation_raw = struct.unpack('<h', header[36:38])[0]
elevation = elevation_raw * 0.1
# 解析方位角（第41-44字节，16位整数，单位0.1度）
azimuth_raw = struct.unpack('<h', header[40:42])[0]
azimuth = azimuth_raw * 0.1
return {'elevation': elevation, 'azimuth': azimuth}
我的小技巧： 拿到新格式的数据，先写个dump函数，把头文件所有字节按16进制打印出来。对照文档逐字节核对，比啥都靠谱。
3.2 质量控制：去伪存真
雷达数据里啥妖魔鬼怪都有。地物杂波、超折射、衰减、噪声——你想想看，风电场周围有山、有风机、有建筑物，这些都会产生回波。如果不做质量控制，模型学到的全是假信号。
质量控制我一般分三步走：
- 地物杂波剔除——用多普勒速度信息，静止目标的速度接近0，直接滤掉
- 超折射识别——看反射率因子的空间梯度，异常高的区域多半是超折射
- 衰减订正——用雨区后的回波衰减量，反推真实反射率
我曾经在福建一个沿海风电场遇到过怪事：雷达显示海面上有一大片强回波，但现场观测明明是大晴天。后来排查发现，是海面蒸发波导导致的超折射。那次之后，我养成了一个习惯——质量控制完一定要跟地面测风塔数据做交叉验证。
注意： 质量控制不是越狠越好。过度滤波会丢失弱回波信息，尤其是低云和雾的回波，这些对风电功率预测其实很重要。我建议保留原始数据，质量控制结果单独存一列，方便回溯。
# 简单的杂波剔除逻辑
def clutter_filter(reflectivity, velocity, threshold=0.5):
"""
reflectivity: 反射率因子 (dBZ)
velocity: 多普勒速度 (m/s)
threshold: 速度阈值，低于此值视为地物杂波
"""
clutter_mask = np.abs(velocity) < threshold
reflectivity_clean = np.where(clutter_mask, np.nan, reflectivity)
return reflectivity_clean
3.3 插值与网格化：从极坐标到笛卡尔坐标
雷达原始数据是极坐标的——每个库（range bin）按距离和方位角排列。但风电模型需要的是规则网格上的数据。所以，插值是绕不开的一步。
常用的插值方法有三种：
- 最近邻插值——快，但锯齿效应明显
- 双线性插值——平滑，但会模糊小尺度特征
- 克里金插值——考虑空间相关性，效果好但慢
我个人在风电场景里偏爱双线性插值。为什么？因为风电场尺度通常几十公里，雷达分辨率足够，双线性插值的精度完全够用，而且计算速度快。你想想看，实时预测系统每6分钟要跑一次，克里金那计算量扛不住。
下面这张图展示了极坐标到笛卡尔网格的转换逻辑：
# 双线性插值示例
from scipy.interpolate import RegularGridInterpolator
def polar_to_cartesian(r, az, data, x_grid, y_grid):
"""
r: 距离数组 (km)
az: 方位角数组 (度)
data: 极坐标数据 (r x az)
x_grid, y_grid: 目标笛卡尔网格
"""
# 将极坐标转为笛卡尔坐标
R, AZ = np.meshgrid(r, np.deg2rad(az))
X = R * np.cos(AZ)
Y = R * np.sin(AZ)
# 构建插值函数
points = np.column_stack([X.ravel(), Y.ravel()])
values = data.ravel()
# 使用griddata进行插值
from scipy.interpolate import griddata
cartesian_data = griddata(points, values, (x_grid, y_grid), method='linear')
return cartesian_data
避坑指南： 我曾经在插值时忘了处理雷达静锥区（雷达正上方无数据区域）。结果网格化后，风电场位置正好落在静锥区里，模型直接报错。后来我加了个距离权重掩膜，静锥区用最近的有效数据填充——虽然不完美，但至少不会崩。
3.4 预处理流程总结
好了，咱们捋一遍。雷达数据预处理，说白了就是三步：
- 格式解析——把二进制或HDF5文件里的数据读出来，转成numpy数组
- 质量控制——剔除杂波、识别超折射、订正衰减
- 插值网格化——从极坐标转到笛卡尔网格，方便后续模型使用
这三步做完，数据才能算是“干净”的。我见过不少团队，花大把时间调模型参数，结果预处理这步草草了事。嗯，最后模型精度上不去，还以为是算法问题——其实根子在数据上。
记住一句话：垃圾进，垃圾出。雷达数据预处理，值得你花时间。