一、测雨雷达定量监测降雨产品与雨量站的点对点评估流程(含python代码)
以雨量站降雨量观测数据为基准,与测雨雷达定量监测降雨数据建立匹配关系,得到雨量站累计雨量和对应点的雷达累计雨量数据,对测雨雷达定量监测降雨结果进行检验评估。良好的测雨雷达降雨算法产品是流域面雨量应用的前提和基础。
业务流程:
选取包含完整降雨过程的时段(该降雨过程持续时间至少为 6 小时,西北地区则至少为3小时),对整点测雨雷达1小时累计降水量产品(OHP)进行点对点评估。
(1)降雨过程筛选:雷达定量监测和雨量计观测都有降水的数目为有效点数为6, 一小时雨量大于10 mm的区间的点数占总有效点数比例少于15%时,则所选时段不具有代表性,不对此过程进行评定;另外,当各降雨量区间存在无样本的情况时(即1 小时雨量均小于等于10 mm),评分不生成。
(2)基于经纬度将所选时段逐小时整点的 OHP 产品与雷达范围内雨量站的整点一小时累计降水量进行匹配;
(3)计算所选时间段内的有效格点的RME和RMB,并利用计算得到的RME和RMB计算定量监测降雨产品综合评分。
说明:均方根误差(RMSE)、相对平均误差(RME)、相对平均偏差(RMB)各指标定义如下:
(1)均方根误差(RMSE)的大小反映了测雨雷达监测降水与雨量计观测降水的离散程度,值越小离散程度越小。
(2)相对平均误差(RME)的大小反映了测雨雷达监测降水与雨量计观测降水之间的平均相对差异程度,值越小差异程度越小。计算公式为:
(3)相对平均偏差(RMB)是指测雨雷达监测降水与雨量计观测降水的实际偏差占观测平均值比例,用来衡量单项测定结果对平均值的偏离程度,值越小偏离程度越低。
测雨雷达定量监测降雨产品评分方法
使用综合评分对雷达定量监测降雨产品进行评价并排序,评估OHP。其次需满足水利部发布的《水利测雨雷达系统建设应用技术要求》中的RMSE标准。综合评分Score计算公式:
= ((5) ×(1-|(5)|)
+(10) ×(1-|(10)|)
+(20) × (1-|(20)|)
+(5) ×(1-(5))
+(10) ×(1-(10))
+(20) ×(1-(20)))×100
其中,RME(5),RME(10),RME(20)分别代表的是当雨量站OHP为5 ≤ OHP < 10、10 ≤ OHP < 20、OHP ≥ 20 时,测雨雷达监测雨量与雨量站观测雨量之间的RME值,(5),(10),(20)代表的是对应的权重系数。RMB(5),RMB(10),RMB(20)分别代表的是当雨量站OHP为5 ≤ OHP < 10、10 ≤ OHP < 20、OHP ≥ 20时,测雨雷达监测雨量与雨量站观测雨量之间的RMB值。(5),(10),(20)代表的是对应的权重系数;综合评分结果换算成百分制后保留两位小数。各个权重系系数
wRME(5) 对应权重 0.05
wRME(10) 0.1
wRME(20) 0.25
wRMB(5) 0.1
wRMB(10) 0.2
wRMB(20) 0.3
针对降雨估计结果进行等级评定,划分为“优秀”、“良好”、“一般”和“差”四个等级,分别对应Score ≥ 70、60 ≤ Score < 70、50 ≤ Score< 60和Score < 50;其中,“优秀”、“良好”、“一般”三个等级代表业务可用,而等级“差”代表业务不可用。
二、淮河上游未来3天降雨、流量预测。
淮河干流洪河口以上为上游,河长约 360 km,地面落差 178 m,流域面积 3.06 万km^2。
台风降暴雨吹到西北400毫米等降雨量域外了。淮河上游降雨偏少,河道处于低水位运行。
(注:图中数据、底图来源于网络,预测模型基于天然地形地貌和气候预报降雨,仅制图作学术探讨供参考。气象水文预报见气象、水文部门、水利部门!)
三、测雨雷达定量监测降雨产品与雨量站的点对点评估python代码:
【import numpy as np
import pandas as pd
def calculate_rmse(radar_data, gauge_data):
"""计算均方根误差(RMSE)"""
return np.sqrt(np.mean((radar_data - gauge_data) ** 2))
def calculate_rme(radar_data, gauge_data):
"""计算相对平均误差(RME)"""
gauge_mean = np.mean(gauge_data)
if gauge_mean == 0:
return 0 # 避免除以零
return np.mean(radar_data - gauge_data) / gauge_mean
def calculate_rmb(radar_data, gauge_data):
"""计算相对平均偏差(RMB)"""
gauge_mean = np.mean(gauge_data)
if gauge_mean == 0:
return 0 # 避免除以零
return np.mean(np.abs(radar_data - gauge_data)) / gauge_mean
def generate_simulated_data(n=1000, precision=0.1):
"""
生成模拟数据:雷达监测降水和雨量计观测降水
数据精确到0.1毫米
"""
# 生成雨量计观测降水(OHP),符合实际降雨分布
# 大部分时间降雨量较小,偶尔有较大降雨
gauge_obs = np.random.exponential(scale=5, size=n) # 指数分布,均值为5
# 雷达监测降水,基于雨量计观测值加入一些误差
# 误差随雨量增大而增大,但保持在合理范围内
error_scale = 0.05 * gauge_obs + 0.3 # 误差尺度
radar_obs = gauge_obs + np.random.normal(loc=0, scale=error_scale, size=n)
# 确保降雨量不为负
radar_obs = np.maximum(radar_obs, 0)
gauge_obs = np.maximum(gauge_obs, 0)
# 四舍五入到指定精度(0.1毫米)
# 先乘以10,四舍五入后再除以10,确保精确到0.1毫米
radar_obs = np.round(radar_obs / precision) * precision
gauge_obs = np.round(gauge_obs / precision) * precision
return pd.DataFrame({
'radar_obs': radar_obs, # 雷达监测降水
'gauge_obs': gauge_obs # 雨量计观测降水(OHP)
})
def evaluate_rainfall_product(data):
"""评估降雨产品并返回综合评分和等级"""
# 按雨量区间分组
group1 = data[(data['gauge_obs'] >= 5) & (data['gauge_obs'] < 10)] # 5 ≤ OHP < 10
group2 = data[(data['gauge_obs'] >= 10) & (data['gauge_obs'] < 20)] # 10 ≤ OHP < 20
group3 = data[data['gauge_obs'] >= 20] # OHP ≥ 20
# 计算各区间的RME和RMB,空组则为0
rme5 = calculate_rme(group1['radar_obs'], group1['gauge_obs']) if not group1.empty else 0
rme10 = calculate_rme(group2['radar_obs'], group2['gauge_obs']) if not group2.empty else 0
rme20 = calculate_rme(group3['radar_obs'], group3['gauge_obs']) if not group3.empty else 0
rmb5 = calculate_rmb(group1['radar_obs'], group1['gauge_obs']) if not group1.empty else 0
rmb10 = calculate_rmb(group2['radar_obs'], group2['gauge_obs']) if not group2.empty else 0
rmb20 = calculate_rmb(group3['radar_obs'], group3['gauge_obs']) if not group3.empty else 0
# 权重系数(来自表15)
weights = {
'w_rme5': 0.05,
'w_rme10': 0.1,
'w_rme20': 0.25,
'w_rmb5': 0.1,
'w_rmb10': 0.2,
'w_rmb20': 0.3
}
# 计算综合评分(公式10)
score = (
weights['w_rme5'] * (1 - abs(rme5)) +
weights['w_rme10'] * (1 - abs(rme10)) +
weights['w_rme20'] * (1 - abs(rme20)) +
weights['w_rmb5'] * (1 - rmb5) +
weights['w_rmb10'] * (1 - rmb10) +
weights['w_rmb20'] * (1 - rmb20)
) * 100
# 确保评分在0-100之间
score = max(0, min(100, score))
# 确定等级
if score >= 70:
grade = "优秀"
elif score >= 60:
grade = "良好"
elif score >= 50:
grade = "一般"
else:
grade = "差"
# 计算整体RMSE(用于检查是否满足准入条件)
overall_rmse = calculate_rmse(data['radar_obs'], data['gauge_obs'])
return {
'RME(5)': round(rme5, 4),
'RME(10)': round(rme10, 4),
'RME(20)': round(rme20, 4),
'RMB(5)': round(rmb5, 4),
'RMB(10)': round(rmb10, 4),
'RMB(20)': round(rmb20, 4),
'综合评分': round(score, 2),
'等级': grade,
'整体RMSE': round(overall_rmse, 4)
}
def main():
# 生成模拟数据,精确到0.1毫米
print("正在生成模拟数据(精确到0.1毫米)...")
data = generate_simulated_data(n=2000, precision=0.1)
# 输出数据基本信息
print("\n模拟数据统计:")
print(f"总数据点数量: {len(data)}")
print(f"5 ≤ OHP < 10 数据点数量: {len(data[(data['gauge_obs'] >= 5) & (data['gauge_obs'] < 10)])}")
print(f"10 ≤ OHP < 20 数据点数量: {len(data[(data['gauge_obs'] >= 10) & (data['gauge_obs'] < 20)])}")
print(f"OHP ≥ 20 数据点数量: {len(data[data['gauge_obs'] >= 20])}")
# 展示部分数据样例,验证精度
print("\n数据样例(精确到0.1毫米):")
print(data.head(10))
# 评估降雨产品
print("\n正在评估降雨产品...")
results = evaluate_rainfall_product(data)
# 输出评估结果
print("\n=== 降雨产品评估结果 ===")
for key, value in results.items():
print(f"{key}: {value}")
# 业务可用性判断
print("\n业务可用性判断:", "可用" if results['等级'] != "差" else "不可用")
if __name__ == "__main__":
main()
】