CN113281754B - 一种雨量站融合雷达定量估测降雨的WRF-Hydro关键参数率定方法 - Google Patents

一种雨量站融合雷达定量估测降雨的WRF-Hydro关键参数率定方法 Download PDF

Info

Publication number
CN113281754B
CN113281754B CN202110841644.3A CN202110841644A CN113281754B CN 113281754 B CN113281754 B CN 113281754B CN 202110841644 A CN202110841644 A CN 202110841644A CN 113281754 B CN113281754 B CN 113281754B
Authority
CN
China
Prior art keywords
rainfall
radar
wrf
station
hydro
Prior art date
Legal status (The legal status is an assumption and is not a legal conclusion. Google has not performed a legal analysis and makes no representation as to the accuracy of the status listed.)
Active
Application number
CN202110841644.3A
Other languages
English (en)
Other versions
CN113281754A (zh
Inventor
刘佳
李传哲
王维
于福亮
顾艳玲
Current Assignee (The listed assignees may be inaccurate. Google has not performed a legal analysis and makes no representation or warranty as to the accuracy of the list.)
China Institute of Water Resources and Hydropower Research
Original Assignee
China Institute of Water Resources and Hydropower Research
Priority date (The priority date is an assumption and is not a legal conclusion. Google has not performed a legal analysis and makes no representation as to the accuracy of the date listed.)
Filing date
Publication date
Application filed by China Institute of Water Resources and Hydropower Research filed Critical China Institute of Water Resources and Hydropower Research
Priority to CN202110841644.3A priority Critical patent/CN113281754B/zh
Publication of CN113281754A publication Critical patent/CN113281754A/zh
Application granted granted Critical
Publication of CN113281754B publication Critical patent/CN113281754B/zh
Active legal-status Critical Current
Anticipated expiration legal-status Critical

Links

Images

Classifications

    • G—PHYSICS
    • G01—MEASURING; TESTING
    • G01S—RADIO DIRECTION-FINDING; RADIO NAVIGATION; DETERMINING DISTANCE OR VELOCITY BY USE OF RADIO WAVES; LOCATING OR PRESENCE-DETECTING BY USE OF THE REFLECTION OR RERADIATION OF RADIO WAVES; ANALOGOUS ARRANGEMENTS USING OTHER WAVES
    • G01S13/00—Systems using the reflection or reradiation of radio waves, e.g. radar systems; Analogous systems using reflection or reradiation of waves whose nature or wavelength is irrelevant or unspecified
    • G01S13/88—Radar or analogous systems specially adapted for specific applications
    • G01S13/95—Radar or analogous systems specially adapted for specific applications for meteorological use
    • G—PHYSICS
    • G01—MEASURING; TESTING
    • G01S—RADIO DIRECTION-FINDING; RADIO NAVIGATION; DETERMINING DISTANCE OR VELOCITY BY USE OF RADIO WAVES; LOCATING OR PRESENCE-DETECTING BY USE OF THE REFLECTION OR RERADIATION OF RADIO WAVES; ANALOGOUS ARRANGEMENTS USING OTHER WAVES
    • G01S7/00—Details of systems according to groups G01S13/00, G01S15/00, G01S17/00
    • G01S7/02—Details of systems according to groups G01S13/00, G01S15/00, G01S17/00 of systems according to group G01S13/00
    • G01S7/41—Details of systems according to groups G01S13/00, G01S15/00, G01S17/00 of systems according to group G01S13/00 using analysis of echo signal for target characterisation; Target signature; Target cross-section
    • G—PHYSICS
    • G01—MEASURING; TESTING
    • G01W—METEOROLOGY
    • G01W1/00—Meteorology
    • G01W1/02—Instruments for indicating weather conditions by measuring two or more variables, e.g. humidity, pressure, temperature, cloud cover or wind speed
    • G—PHYSICS
    • G01—MEASURING; TESTING
    • G01W—METEOROLOGY
    • G01W1/00—Meteorology
    • G01W1/10—Devices for predicting weather conditions
    • G—PHYSICS
    • G01—MEASURING; TESTING
    • G01W—METEOROLOGY
    • G01W1/00—Meteorology
    • G01W1/14—Rainfall or precipitation gauges
    • Y—GENERAL TAGGING OF NEW TECHNOLOGICAL DEVELOPMENTS; GENERAL TAGGING OF CROSS-SECTIONAL TECHNOLOGIES SPANNING OVER SEVERAL SECTIONS OF THE IPC; TECHNICAL SUBJECTS COVERED BY FORMER USPC CROSS-REFERENCE ART COLLECTIONS [XRACs] AND DIGESTS
    • Y02—TECHNOLOGIES OR APPLICATIONS FOR MITIGATION OR ADAPTATION AGAINST CLIMATE CHANGE
    • Y02A—TECHNOLOGIES FOR ADAPTATION TO CLIMATE CHANGE
    • Y02A10/00—TECHNOLOGIES FOR ADAPTATION TO CLIMATE CHANGE at coastal zones; at river basins
    • Y02A10/40—Controlling or monitoring, e.g. of flood or hurricane; Forecasting, e.g. risk assessment or mapping
    • Y—GENERAL TAGGING OF NEW TECHNOLOGICAL DEVELOPMENTS; GENERAL TAGGING OF CROSS-SECTIONAL TECHNOLOGIES SPANNING OVER SEVERAL SECTIONS OF THE IPC; TECHNICAL SUBJECTS COVERED BY FORMER USPC CROSS-REFERENCE ART COLLECTIONS [XRACs] AND DIGESTS
    • Y02—TECHNOLOGIES OR APPLICATIONS FOR MITIGATION OR ADAPTATION AGAINST CLIMATE CHANGE
    • Y02A—TECHNOLOGIES FOR ADAPTATION TO CLIMATE CHANGE
    • Y02A90/00—Technologies having an indirect contribution to adaptation to climate change
    • Y02A90/10—Information and communication technologies [ICT] supporting adaptation to climate change, e.g. for weather forecasting or climate simulation

Landscapes

  • Engineering & Computer Science (AREA)
  • Environmental & Geological Engineering (AREA)
  • Life Sciences & Earth Sciences (AREA)
  • Remote Sensing (AREA)
  • Radar, Positioning & Navigation (AREA)
  • Environmental Sciences (AREA)
  • Ecology (AREA)
  • Physics & Mathematics (AREA)
  • Biodiversity & Conservation Biology (AREA)
  • Atmospheric Sciences (AREA)
  • Computer Networks & Wireless Communication (AREA)
  • General Physics & Mathematics (AREA)
  • Hydrology & Water Resources (AREA)
  • Electromagnetism (AREA)
  • Radar Systems Or Details Thereof (AREA)

Abstract

本发明提供一种雨量站融合雷达定量估测降雨的WRF‑Hydro关键参数率定方法，包括设置WRF中的气象驱动数据与预热时间；雨量站融合多普勒天气雷达定量估测降雨获得修正降雨数据；降雨驱动数据的适用性判断；基于并行的动态多维搜索算法（P‑DDS）的WRF‑Hydro模式参数校准；验证校准参数后的WRF‑Hydro对洪水过程模拟效果。雷达与雨量站融合降雨的能够提高输入模式的时空准确性，能够更加稳定的刻画时空分布不均匀的降雨，以其作为率定WRF‑Hydro陆面水文模式的降雨驱动数据，并配合并行动态多维搜索算法（P‑DDS），能够改进WRF‑Hydro的参数校准过程，提高场次降雨洪水的模拟精度。

Description

一种雨量站融合雷达定量估测降雨的WRF-Hydro关键参数率 定方法
技术领域
本发明涉及水文气象的技术领域，特别是一种雨量站融合雷达定量估测降雨的WRF-Hydro关键参数率定方法。
背景技术
WRF-Hydro是近几年国际上比较复杂的陆面水文模式代表，可通过WRF和WRF-Hydro的耦合提供水文过程时空要素的模拟/预报信息，但是陆面水文模式WRF-Hydro参数的率定却受制于输入数据以及率定方式的影响。通常的参数定采用WRF与WRF-Hydro耦合方式，输入的降雨并非实际的“真值”，会极大的影响模拟/预报效果。
雨量计作为一种直接观测降雨的传感器可以在相对精确的位置对降雨进行点雨量观测，但雨量站的维护检测成本相对较高，在北方地区，其站网密度往往不能满足流域空间建模的需求，较难捕捉到降雨在空间上的异质性。雷达站资料具有对流尺度分辨率高，观测准确度高等优点，单个雷达站覆盖半径200~300km，非常适合中小尺度气象观测。从1940年利用天气雷达观测雨量的空间特性开始，到近几十年间天气雷达数据的不断更新，在雨量估测精度方面不断提升，应用也越来越广泛。将多雷达观测数据转化为降雨强度的过程比较复杂，主要原因在于观测数据往往存在误差，现有的雷达站网由于布设较早，硬件设置并不完善，预报质量尚存在一定的不确定性。雨量站与雷达各有优势，融合高分辨率的降雨产品是获得稳定的高分辨率的降雨数据的有效途径，准确的降雨时空分布将进一步提高所构建的陆面水文模式系统模拟精度。
WRF-Hydro在水文、气象等多个学科交叉领域上均的广阔应用前景，近年来受到国内外学者的广泛关注。当前WRF-Hydro已经作为美国国家水预报模型（National WaterModel）之一用于指导美国西北部的水文预报业务。2017年Naabil E、Lamptey BL、ArnaultJ、Olufayo A、Kunstmann H的《利用WRF-Hydro系统进行水资源管理：以西非托诺大坝为例》用单向耦合的WRF-Hydro成功模拟了西非的一个小流域2014年的径流，但是发现部分时段内存在快速回落情况；2017年Silver M、Karnieli A、Ginat H、Meiri E、Fredj E的《干旱区WRF-Hydro水文参数率定的一种新方法》对美国干旱与半干旱区的5场洪水进行模拟，并通过土地利用类型将参数变量值进行了分类修正，该方法调参后的流量过程明显优于系统默认参数取值下的结果；2020年孙明坤、李致家、刘志雨、侯爱中、霍文博、温娅惠、孔祥意、戴金旺、梁世强的《WRF-Hydro模型与新安江模型在陈河流域的应用对比》比较了WRF-Hydro和新安江模型在陕西半湿润区的应用效果，结果表明WRF-Hydro能够更好的反映洪水细节，但对中小流域的洪水模拟精度仍需加强。当陆面水文模式用于对洪水过程进行时空动态模拟时，模式关键产汇流参数仍需依赖可靠的高分辨降雨数据，上述研究对此并没有很好的解决方案。
发明内容
为了解决上述的技术问题，本发明提出的一种雨量站融合雷达定量估测降雨的WRF-Hydro关键参数率定方法，其解决的主要技术问题是针对低密度雨量站无法满足空间建模精度的需求，采用雨量站与多普勒天气雷达定量估测降雨融合的方法提高降雨时空分布特征。
本发明的目的是提供一种雨量站融合雷达定量估测降雨的WRF-Hydro关键参数率定方法，包括以下步骤：
步骤1：设置天气预报模式WRF中的气象驱动数据与预热时间，所述气象驱动数据包括温度驱动数据、压强驱动数据、风速驱动数据和降雨驱动数据；
步骤2：雨量站融合多普勒天气雷达定量估测降雨获得修正降雨数据；
步骤3：降雨驱动数据的适用性判断；
步骤4：基于并行的动态多维搜索算法的WRF-Hydro模式参数校准；
步骤5：验证校准参数后的WRF-Hydro对洪水过程模拟效果。
优选的是，所述WRF采用三层嵌方式并取最内层区域用于为所述WRF-Hydro提供降雨、辐射、湿度、温度、气压、风速中至少一种模式的输入。
在上述任一方案中优选的是，在所述WRF的预热过程中设置2个预热时段：第一个预热时段为降雨开始前的16天，并在第15天时输出Restart文件，使陆面预热能够获得相对准确的模式下边界数据，并在预热结束生成Restart；第二个预热时段为降雨开始前的24h，第二个预热时段调用所述Restart文件中的下边界数据与WRF驱动资料中的侧边界数据。
在上述任一方案中优选的是，所述修正降雨数据的处理方法包括以下三种：雷达定量估测降雨法、基于反距离插值权重平均的雨量站插值降雨法和外部漂移克里金融合雨量站与雷达降雨法。
在上述任一方案中优选的是，所述雷达定量估测降雨法对多普勒雷达获得的降雨分别进行了电磁波衰减校正、非降雨回波抑制、不同层反射率组合以及雨强-反射率转换，采用HB方法来订正雷达回波衰减，设定路径积分衰减系数PIA的上限为10dB；采用基于三维组合反射率抑制非降雨回波，并设定仰角4°以下的回波为杂波；结合DEM数据和雷达波束标准传播方程获得了不同高度层的组合反射率选取规则；结合流域雨滴谱对不同类型的反射率因子与雨强系数进行调整。
在上述任一方案中优选的是，所述基于反距离插值权重平均的雨量站插值降雨法根据插值区域内各采样点之间的相似度，设定网格对应未知点处的估计值由局部邻域内所有数据点的距离加权平均值表示，公式为
Figure 873368DEST_PATH_IMAGE001
其中，Z 0 为对应网格插值点处的估计值，Z i 为第i个雨量站观测数据；d i 为需要插值的网格点与雨量站点间的空间距离；n为研究区雨量站的个数。
在上述任一方案中优选的是，所述外部漂移克里金融合雨量站与雷达降雨法中的外部漂移克里金插值中主变量的期望值与附加变量之间符合如下线性关系：
Figure 144860DEST_PATH_IMAGE002
其中，主变量的期望G(x)为对应x点处的雨量站的观测数据，附加变量R(x)为雷达数据得到的降雨估计值；a和b为由假设条件确立的线性方程的系数。
在上述任一方案中优选的是，对于给定位置x 0的估计由线性估计器计算得到，对应的权重因子表示为
Figure 299898DEST_PATH_IMAGE003
Figure 303626DEST_PATH_IMAGE004
其中，
Figure 761283DEST_PATH_IMAGE005
为在x 0点处与降雨估计值相关的权重因子，n为研究区雨量站的个数，R 1 (x i )为雨量站点i处的雨量值，R(x 0)表示x 0点经过KED融合后的降雨估算值。
在上述任一方案中优选的是，在所述外部漂移克里金融合雨量站与雷达降雨法中增加额外的约束条件用于计算雷达场和雨量站场残差的协方差：
Figure 748831DEST_PATH_IMAGE006
Figure 441980DEST_PATH_IMAGE007
Figure 113264DEST_PATH_IMAGE008
其中，C表示残差在不同点间的残差协方差结果，x i 为雨量站点i处，x j 为雨量站点j 处，μ 0和μ 1为拉格朗日乘数，R 2 (x j )为雷达网格点j处的雨量值，R 2 (x 0 )表示x 0 点经过权重处 理后得到的雷达降雨估算值，
Figure 866457DEST_PATH_IMAGE009
为在x 0 处与雷达估计值相关的权重因子。
在上述任一方案中优选的是，在雨量站位置处插值雷达时，上述约束使用一组权重对雨量站数据进行加权，同时这些权重也会增加估计点x 0 处的雷达值效果。
在上述任一方案中优选的是，所述步骤3包括对3种基于雷达或雨量站的空间降雨驱动数据精度进行分析，包括以下子步骤：
步骤31：针对每种降雨空间数据，均采用雨量站点的留一交叉验证方式，选择了系统偏差BIAS、均方根误差RMSE和均方根转换误差MRTE三个降雨精度评价指标；
步骤32：选择了雨量站点空间散点和降雨数据的累积降雨量空间分布图来对3种降雨数据和WRF降雨的空间降雨情况进行综合评价。
在上述任一方案中优选的是，三个所述降雨精度评价指标的计算公式为
Figure 606880DEST_PATH_IMAGE010
Figure 369299DEST_PATH_IMAGE011
Figure 895090DEST_PATH_IMAGE012
其中，P t 为雨量站测得的第t时刻的降雨量，i为第i个雨量站，n为研究区雨量站的 个数，
Figure 819183DEST_PATH_IMAGE013
为不同降雨驱动数据下对应雨量站位置处的估计值。
在上述任一方案中优选的是，所述步骤4包括采用适用于模式的并行自动调参方法对敏感参数进行率定，在所述率定的过程中选取克林-古普塔效率系数KGE作为目标函数对模拟的洪水过程进行评估，其中，敏感参数为容易引起WRF-Hydro产汇流过程发生明显变化的参数。
在上述任一方案中优选的是，所述克林-古普塔效率系数KGE的公式为
Figure 515744DEST_PATH_IMAGE014
其中，γ为模拟系列与观测系列标准偏差的比率，α为模拟系列与观测系列平均值的比率，β则为模拟系列与观测系列的相关系数。
在上述任一方案中优选的是，所述步骤4还包括采用并行动态多维搜索算法进行路面水温模式率定过程。
在上述任一方案中优选的是，所述并行动态多维搜索算法的串行计算步骤如下
步骤41：定义算法输入变量；
步骤42：模式初值分析；
步骤43：解集领域概率筛选；
步骤44：加入额外扰动因子r；
步骤45：候选解集X new 校准更新；
步骤46：终止条件检查。
在上述任一方案中优选的是，所述算法输入变量包括邻近扰动距离因子r、算法最大校准次数m、针对D维变量每一个决策变量的初始解的下限X min 和上限值X max 、算法初始解向量为X 0 =[x 1 ，x 2 ，…，x D ]和分配的模式处理节点数N。
在上述任一方案中优选的是，所述步骤42包括迭代次数i从1开始，计算当前初始解对应的目标函数值F（X 0 ），将初值解暂时赋值为最优解，即X best =X 0 ，F best =F（X 0 ）。
在上述任一方案中优选的是，所述步骤43包括随机选择D个决策变量中的J个产生新的解集邻域M，计算决策变量在当前迭代计数下在M中发生的概率：P（i）=1-ln（i）/ln（m），对于d=1，…，D个决策变量，按求得的概率P（i）将d个决策变量添加到解集邻域M中，当M出现空集时，则任选d加入M中。
在上述任一方案中优选的是，所述步骤44包括对于M中的决策变量，确定j=1，…，J上的X best ，并在j中加入新的扰动因子，并假设加入的新的扰动服从B（0，1）的标准正态分布，则加入扰动因子后j个决策变量的波动范围可表示为：
Figure 423132DEST_PATH_IMAGE015
其中，
Figure 928063DEST_PATH_IMAGE016
为新的参数值，
Figure 147692DEST_PATH_IMAGE017
为当前参数最佳值，
Figure 3652DEST_PATH_IMAGE018
为新扰动因子系数，
Figure 717662DEST_PATH_IMAGE019
为当前 参数范围最大值，
Figure 342678DEST_PATH_IMAGE020
为当前参数范围最小值。
在上述任一方案中优选的是，所述步骤45包括计算所述步骤44中生成的解集X new 的目标函数F（X new ），当F（x new ）≥F best ，F best =F（x new ），x best =x new ，更新当前最优解。
在上述任一方案中优选的是，所述步骤46包括当i<m时，返回所述步骤43，并更新迭代次数i=i+1；否则，算法结束。
在上述任一方案中优选的是，所述步骤5包括使用评价指标对洪水过程的适应性进行判断。
在上述任一方案中优选的是，所述评价指标包括针对不同降雨-参数组合的泰勒图、克林-古普塔效率系数、洪峰误差、洪量误差和纳什效率系数。
在上述任一方案中优选的是，所述洪峰误差R p 的公式表示为
Figure 733208DEST_PATH_IMAGE021
其中，
Figure 810885DEST_PATH_IMAGE022
为模拟的洪峰流量，q p 为观测的洪峰流量。
在上述任一方案中优选的是，所述洪量误差R r 的公式表示为
Figure 63006DEST_PATH_IMAGE023
其中，
Figure 870425DEST_PATH_IMAGE024
为模拟洪量，r r 为观测洪量。
在上述任一方案中优选的是，所述纳什效率系数NSE的公式表示为
Figure 307223DEST_PATH_IMAGE025
其中，i为时间步长，T为洪水工程的总时长，
Figure 747563DEST_PATH_IMAGE026
为模拟流量，q i 为观测流量，
Figure 928008DEST_PATH_IMAGE027
为平 均流量。
本发明提出了一种雨量站融合雷达定量估测降雨的WRF-Hydro关键参数率定方法，采用并行动态多维搜索算法（P-DDS）应用于WRF/WRF-Hydro陆气耦合过程参数校准，提升陆面水文模式对洪水过程的模拟精度。
WRF是指天气预报模式。
WRF-Hydro是指近年来逐渐被国际关注的陆面水文模式。
HB方法是指Hitschfeld and Borden提出的衰减订正方法。
DEM数据是指数字高程模型数据。
附图说明
图1为按照本发明的雨量站融合雷达定量估测降雨的WRF-Hydro关键参数率定方法的一优选实施例的流程图。
图2为按照本发明的雨量站融合雷达定量估测降雨的WRF-Hydro关键参数率定方法的另一优选实施例的流程图。
图3为按照本发明的雨量站融合雷达定量估测降雨的WRF-Hydro关键参数率定方法的雷达位置与覆盖范围的一实施例的示意图。
图4为按照本发明的雨量站融合雷达定量估测降雨的WRF-Hydro关键参数率定方法的WRF嵌套方案的一实施例的示意图。
图5为按照本发明的雨量站融合雷达定量估测降雨的WRF-Hydro关键参数率定方法的不同降雨数据下留一交叉验证指标的一实施例的结果示意图。
图6为按照本发明的雨量站融合雷达定量估测降雨的WRF-Hydro关键参数率定方法的雨量站观测值与4种降雨数据在对应雨量站点模拟值散点的一实施例的结果示意图。
图7为按照本发明的雨量站融合雷达定量估测降雨的WRF-Hydro关键参数率定方法的4种降雨数据的累积降雨量空间的一实施例的分布示意图。
图8为按照本发明的雨量站融合雷达定量估测降雨的WRF-Hydro关键参数率定方法的RIDW和KED降雨-参数方案下WRF-Hydro模拟的洪水过程的一实施例的指标箱型图。
图9为按照本发明的雨量站融合雷达定量估测降雨的WRF-Hydro关键参数率定方法的WRF-RIDW和WRF-KED降雨-参数方案下WRF-Hydro模拟的洪水过程的一实施例的指标箱型图。
图10为按照本发明的雨量站融合雷达定量估测降雨的WRF-Hydro关键参数率定方法的4种降雨-参数方案WRF-Hydro模拟的洪水过程的一实施例的泰勒图。
具体实施方式
下面结合附图和具体的实施例对本发明做进一步的阐述。
实施例一
如图1所示，执行步骤110，设置天气预报模式WRF中的气象驱动数据与预热时间，所述气象驱动数据包括温度驱动数据、压强驱动数据、风速驱动数据和降雨驱动数据。所述WRF采用三层嵌方式并取最内层区域用于为所述WRF-Hydro提供降雨、辐射、湿度、温度、气压、风速中至少一种模式的输入。在所述WRF的预热过程中设置2个预热时段：第一个预热时段为降雨开始前的16天，并在第15天时输出Restart文件，使陆面预热能够获得相对准确的模式下边界数据，并在预热结束生成Restart；第二个预热时段为降雨开始前的24h，第二个预热时段调用所述Restart文件中的下边界数据与WRF驱动资料中的侧边界数据。
执行步骤120，雨量站融合多普勒天气雷达定量估测降雨获得修正降雨数据。所述修正降雨数据的处理方法包括以下三种：雷达定量估测降雨法、基于反距离插值权重平均的雨量站插值降雨法和外部漂移克里金融合雨量站与雷达降雨法。
所述雷达定量估测降雨法对多普勒雷达获得的降雨分别进行了电磁波衰减校正、非降雨回波抑制、不同层反射率组合以及雨强-反射率转换，采用HB方法来订正雷达回波衰减，设定路径积分衰减系数PIA的上限为10dB；采用基于三维组合反射率抑制非降雨回波，并设定仰角4°以下的回波为杂波；结合DEM数据和雷达波束标准传播方程获得了不同高度层的组合反射率选取规则；结合流域雨滴谱对不同类型的反射率因子与雨强系数进行调整。
所述基于反距离插值权重平均的雨量站插值降雨法根据插值区域内各采样点之间的相似度，设定网格对应未知点处的估计值由局部邻域内所有数据点的距离加权平均值表示，公式为
Figure 652251DEST_PATH_IMAGE028
其中，Z 0 为对应网格插值点处的估计值，Z i 为第i个雨量站观测数据，d i 为需要插值的网格点与雨量站点间的空间距离n为研究区雨量站的个数。
所述外部漂移克里金融合雨量站与雷达降雨法中的外部漂移克里金插值中主变量的期望值与附加变量之间符合如下线性关系：
Figure 525529DEST_PATH_IMAGE002
其中，主变量的期望G(x)为对应x点处的雨量站的观测数据，附加变量R(x)为雷达数据得到的降雨估计值；a和b为由假设条件确立的线性方程的系数。
对于给定位置x 0的估计由线性估计器计算得到，对应的权重因子表示为
Figure 184656DEST_PATH_IMAGE003
Figure 168792DEST_PATH_IMAGE004
其中，
Figure 747541DEST_PATH_IMAGE029
为在x 0点处与降雨估计值相关的权重因子，n为研究区雨量站的个数，R 1 (x i )为雨量站点i处的雨量值，R(x 0)表示x 0点经过KED融合后的降雨估算值。
在所述外部漂移克里金融合雨量站与雷达降雨法中增加额外的约束条件用于计算雷达场和雨量站场残差的协方差：
Figure 604770DEST_PATH_IMAGE030
Figure 878756DEST_PATH_IMAGE031
Figure 525638DEST_PATH_IMAGE008
其中，C表示残差在不同点间的残差协方差结果，x i 为雨量站点i处，x j 为雨量站点j 处，μ 0和μ 1为拉格朗日乘数，R 2 (x j )为雷达网格点j处的雨量值，R 2 (x 0 )表示x 0 点经过权重处 理后得到的雷达降雨估算值，
Figure 834260DEST_PATH_IMAGE009
为在x 0 处与雷达估计值相关的权重因子。
在雨量站位置处插值雷达时，上述约束使用一组权重对雨量站数据进行加权，同时这些权重也会增加估计点x 0 处的雷达值效果。
执行步骤130，降雨驱动数据的适用性判断，对3种基于雷达或雨量站的空间降雨驱动数据精度进行分析，包括以下子步骤：
步骤131：针对每种降雨空间数据，均采用雨量站点的留一交叉验证方式，选择了系统偏差BIAS、均方根误差RMSE和均方根转换误差MRTE三个降雨精度评价指标。三个所述降雨精度评价指标的计算公式为
Figure 924707DEST_PATH_IMAGE010
Figure 951569DEST_PATH_IMAGE011
Figure 136562DEST_PATH_IMAGE012
其中，P t 为雨量站测得的第t时刻的降雨量，i为第i个雨量站，n为研究区雨量站的 个数，
Figure 565270DEST_PATH_IMAGE013
为不同降雨驱动数据下对应雨量站位置处的估计值。
步骤132：选择了雨量站点空间散点和降雨数据的累积降雨量空间分布图来对3种降雨数据和WRF降雨的空间降雨情况进行综合评价。
执行步骤140，基于并行的动态多维搜索算法（P-DDS）的WRF-Hydro模式参数校准，采用适用于模式的并行自动调参方法对敏感参数进行率定，在所述率定的过程中选取克林-古普塔效率系数KGE作为目标函数对模拟的洪水过程进行评估，其中，敏感参数为容易引起WRF-Hydro产汇流过程发生明显变化的参数。
克林-古普塔效率系数KGE的公式为
Figure 764301DEST_PATH_IMAGE014
其中，γ为模拟系列与观测系列标准偏差的比率，α为模拟系列与观测系列平均值的比率，β则为模拟系列与观测系列的相关系数。
步骤140还包括采用并行动态多维搜索算法进行路面水温模式率定过程，并行动态多维搜索算法的串行计算步骤如下：
步骤141：定义算法输入变量。算法输入变量包括邻近扰动距离因子r、算法最大校准次数m、针对D维变量每一个决策变量的初始解的下限X min 和上限值X max 、算法初始解向量为X 0 =[x 1 ，x 2 ，…，x D ]和分配的模式处理节点数N。
步骤142：模式初值分析；迭代次数i从1开始，计算当前初始解对应的目标函数值F （X 0 ），将初值解暂时赋值为最优解，即X best =X 0 ，F best =F（X 0 ）。
步骤143：解集领域概率；随机选择D个决策变量中的J个产生新的解集邻域M，计算决策变量在当前迭代计数下在M中发生的概率：P（i）=1-ln（i）/ln（m），对于d=1，…，D个决策变量，按求得的概率P（i）将d个决策变量添加到解集邻域M中，当M出现空集时，则任选d加入M中。
步骤144：加入额外扰动因子r；对于M中的决策变量，确定j=1，…，J上的X best ，并在j中加入新的扰动因子，并假设加入的新的扰动服从B（0，1）的标准正态分布，则加入扰动因子后j个决策变量的波动范围可表示为：
Figure 137513DEST_PATH_IMAGE015
其中，
Figure 1564DEST_PATH_IMAGE016
为新的参数值，
Figure 157215DEST_PATH_IMAGE017
为当前参数最佳值，
Figure 448519DEST_PATH_IMAGE018
为新扰动因子系数，
Figure 574606DEST_PATH_IMAGE019
为当前 参数范围最大值，
Figure 242348DEST_PATH_IMAGE020
为当前参数范围最小值。
步骤145：候选解集X new 校准更新；计算所述步骤144中生成的解集X new 的目标函数F （X new ），当F（x new ）≥F best ，F best =F（x new ），x best =x new ，更新当前最优解。
步骤146：终止条件检查；当i<m时，返回所述步骤43，并更新迭代次数i=i+1；否则，算法结束。
执行步骤150，验证校准参数后的WRF-Hydro对洪水过程模拟效果，使用评价指标对洪水过程的适应性进行判断。评价指标包括针对不同降雨-参数组合的泰勒图、克林-古普塔效率系数、洪峰误差、洪量误差和纳什效率系数。
洪峰误差R p 的公式表示为
Figure 255435DEST_PATH_IMAGE021
其中，
Figure 717640DEST_PATH_IMAGE022
为模拟的洪峰流量，q p 为观测的洪峰流量。
洪量误差R r 的公式表示为
Figure 65445DEST_PATH_IMAGE023
其中，
Figure 615506DEST_PATH_IMAGE024
为模拟洪量，r r 为观测洪量。
纳什效率系数NSE的公式表示为
Figure 607733DEST_PATH_IMAGE025
其中，i为时间步长，T为洪水工程的总时长，
Figure 365473DEST_PATH_IMAGE026
为模拟流量，q i 为观测流量，
Figure 75940DEST_PATH_IMAGE027
为平 均流量。
实施例二
如图2所示，本发明设计了一种雨量站融合雷达定量估测降雨的WRF-Hydro关键参数率定方法，其解决的主要技术问题是针对低密度雨量站无法满足空间建模精度的需求，采用雨量站与多普勒天气雷达定量估测降雨融合的方法提高降雨时空分布特征，并采用并行动态多维搜索算法（P-DDS）应用于WRF/WRF-Hydro陆气耦合过程参数校准，提升陆面水文模式对洪水过程的模拟精度。
为了解决上述存在的技术问题，本发明采用了以下方案：
一种雨量站融合雷达定量估测降雨的WRF-Hydro关键参数率定方法，包括以下几个步骤：
步骤1、设置WRF中的气象驱动数据与预热时间；
步骤2、雨量站融合多普勒天气雷达定量估测降雨（QPE）获得修正降雨数据；
步骤3、降雨驱动数据的对比验证；
步骤4、基于并行的动态多维搜索算法（P-DDS）的WRF-Hydro模式参数校准；
步骤5、验证校准参数后的WRF-Hydro对洪水过程模拟效果。
进一步，所述步骤1中WRF采用三层嵌方式（每层嵌套比例为1:3）并取最内层区域用于为WRF-Hydro提供包括降雨、辐射、湿度、温度、气压、风速等在内的模式输入。WRF的预热过程中设置2个预热时段：第一个预热时段为降雨开始前的16天，并在第15天时输出Restart文件，使陆面预热能够获得相对准确的模式下边界数据（例如土壤含水量条件）并在预热结束生成Restart；第二个预热时间为降雨开始前的24h，第二个预热时段调用Restart文件中的下边界数据与WRF驱动资料（ERA-Interim）中的侧边界数据（温度、辐射、风场等），此时的预热协调了两种边界数据，确保了初始WRF模式运行的稳定性。
进一步，所述步骤2中给出了3种雨量站融合多普勒天气雷达定量估测降雨（QPE）获得修正降雨数据处理方法。
雷达定量估测降雨（QPE），对多普勒雷达获得的降雨分别进行了电磁波衰减校正、非降雨回波抑制、不同层反射率组合以及雨强-反射率转换。用HB方法来订正雷达回波衰减，设定路径积分衰减系数（Path Integration Attenuation，PIA）PIA的上限为10dB；采用基于三维组合反射率抑制非降雨回波，并设定仰角4°以下的回波为杂波；结合DEM数据和雷达波束标准传播方程获得了不同高度层的组合反射率选取规则；结合流域雨滴谱对不同类型的反射率因子与雨强系数进行调整。
基于反距离插值权重平均（RIDW）的雨量站插值降雨，该方法根据插值区域内各采样点之间的相似度，假设网格对应未知点处的估计值可由局部邻域内所有数据点的距离加权平均值表示，该方法既考虑了常用的泰森多边形中邻点的概念也考虑了趋势面分析中的渐变过程，公式可以表达为：
Figure 960851DEST_PATH_IMAGE032
（1）
式中：Z 0 为对应网格插值点处的估计值；Z i 为第i个雨量站观测数据；d i 为需要插值的网格点与雨量站点间的空间距离；n为研究区雨量站的个数。
外部漂移克里金（KED）插值方法是地统计方法中广义克里金插值法的一种扩展，因此该方法仍将插值变量用确定性项（漂移项）和随机项（残差项）的加和表示。不同的是，带外部漂移的克里金允许合并多个附加变量作为背景信息来对主变量进行插值。研究雷达与雨量站资料融合时，雷达资料被视为唯一的附加变量。KED插值中主变量的期望值与附加变量之间符合如下线性关系：
Figure 73163DEST_PATH_IMAGE002
（2）
式中：G(x)代表主变量的期望，这里G(x)指的是对应x点处的雨量站的观测数据；R (x)为附加变量，这里表示雷达数据得到的降雨估计值；a和b为由假设条件确立的线性方程的系数。
上述公式中给定位置x 0 的估计可由线性估计器计算得到，对应的权重因子表示为：
Figure 267384DEST_PATH_IMAGE003
（3）
Figure 199568DEST_PATH_IMAGE004
（4）
式中：R(x 0)表示x 0 点经过KED融合后的降雨估算值；R 1 (x i )为雨量站点i处的雨量 值；n为研究区雨量站的个数；
Figure 885240DEST_PATH_IMAGE033
表示在x 0 点处与降雨估计值相关的权重因子，随着估算 点的改变，其值也相应发生改变。
上式与传统克里金法的描述一致，但是在具体权重因子 计算上，KED方法增加了额外的约束条件（即与雷达观测有关的数据约束）用于计算雷达场和雨量站场残差的协方差：
Figure 914376DEST_PATH_IMAGE034
（5）
Figure 154864DEST_PATH_IMAGE007
（6）
Figure 449710DEST_PATH_IMAGE008
（7）
式中：R 2 (x 0 )表示x 0 点经过权重处理后得到的雷达降雨估算值；为雷达网格点j处 的雨量值；μ为拉格朗日乘数；C表示残差在不同点间的残差协方差结果，可由雷达估测的降 雨场与估算趋势场（漂移场）的差表示；
Figure 801057DEST_PATH_IMAGE009
为在x 0 处与雷达估计值相关的权重因子，其他物 理意义同上。研究中的漂移场是通过传统克里金方式对雷达测算的雨量站位置处降雨进行 空间变异性处理获得，并通过FFT方法，将基于雨量站点的协方差矩阵C转换为整个流域网 格上的频谱矩阵。
在雨量站位置处插值雷达时，上述约束使用一组权重对雨量站数据进行加权，同时这些权重也会增加估计点x 0 处的雷达值效果。需要特别指出的是，该方法无需考虑前期的转换模型及其参数，通过FFT方式实现了完整的雷达数据在流域网格上的残差计算，该方法把雷达数据转换成具有连续、严格累积分布，同时将雷达数据的空间可变性完全引入到雨量站值的内插过程中，能更好的解译降雨的时空特征，从而更适合陆气耦合建模工作。
进一步，所述步骤3对步骤2中3种基于雷达或雨量站的空间降雨驱动数据精度进行分析，主要体现在：首先针对每种降雨空间数据，均采用雨量站点的留一交叉验证方式，选择了系统偏差BIAS、均方根误差RMSE和均方根转换误差MRTE三个降雨精度评价指标；选择了雨量站点空间散点和降雨数据的累积降雨量空间分布图来对3种降雨数据和WRF降雨的空间降雨情况进行综合评价。
3个指标的计算公式如下：
Figure 12596DEST_PATH_IMAGE010
（8）
Figure 423986DEST_PATH_IMAGE011
（9）
Figure 206128DEST_PATH_IMAGE012
（10）
式中：P t 为雨量站测得的第t时刻的降雨量；i为第i个雨量站；n为研究区雨量站的 个数；
Figure 361166DEST_PATH_IMAGE013
为不同降雨驱动数据下对应雨量站位置处的估计值。针对10场降雨，RMSE和MRTE的 范围为0到+∞，最优值为0，BIAS的结果为整个降雨统计时段和所有雨量站的均值，范围在- ∞到+∞，最优值为0。
进一步，所述步骤4中采用适用于模式的并行自动调参方法对敏感参数进行率定，其中，敏感参数为容易引起WRF-Hydro产汇流过程发生明显变化的参数，主要体现在将基于雨量站融合多普勒天气雷达定量估测的降雨、WRF模拟得到的其它气象数据和WRF-Hydro一起用于陆气耦合研究，WRF-Hydro中采用了降尺度聚解以及升尺度聚合方式对产汇流变量进行网格细分，降尺度比例设定一般在3-10倍，选择克林-古普塔效率系数（Kling-Guptaefficiency，KGE）作为目标函数，率定时采用了并行的动态多维搜索算法（P-DDS）。
率定过程中需要对模拟的洪水过程进行评估，此处选取的目标函数为克林-古普塔效率系数（Kling-Gupta efficiency，KGE）。该指标权衡了洪水过程的相关性，水量平衡以及方差误差。与常用的纳什效率系数关注大的洪水过程相比，KGE能够同时将量级较小的洪水过程也考虑进来，更有利于对洪水过程的综合评判。考虑到研究区可供率定的洪水过程不多，该目标函数对量级较小洪水的考虑有助于更好的率定模式参数，KGE可以用下面的公式表示：
Figure 364894DEST_PATH_IMAGE035
（11）
式中：γ为模拟系列与观测系列标准偏差的比率；α为模拟系列与观测系列平均值的比率；β则为模拟系列与观测系列的相关系数。当KGE等于1时，说明模拟能力达到最佳（此时表示与理想模拟系列的欧几里得距离最小），而当KGE为负值时，则认为模拟的效果较差。
为了加快计算资源获取高质量解集的收敛速度，在进行路面水温模式率定过程中采用了并行动态多维搜索算法（Parallel Dynamically Dimensioned Search，P-DDS），动态多维搜索算法（Dynamically Dimensioned Search，DDS）属于随机搜索算法的一种，适合于自动校准计算成本高、参数化程度高的气候与环境模式。DDS的串行计算步骤如下：
Step1：定义算法输入变量。邻近扰动距离因子r（算法缺省值为0.2）；算法最大校准次数m；针对D维变量每一个决策变量的初始解的下限X min 和上限值X max ；定义算法初始解向量为X 0 =[x 1 ，x 2 ，…，x D ]；
Step2：模式初值分析。迭代次数i从1开始，并计算当前初始解对应的目标函数值F （X 0 ），将初值解暂时赋值为最优解，即X best =X 0 ，F best =F（X 0 ）；
Step3：解集领域概率筛选。随机选择D个决策变量中的J个产生新的解集邻域M，计算决策变量在当前迭代计数下在M中发生的概率：P（i）=1-ln（i）/ln（m），对于d=1，…，D个决策变量，按求得的概率P（i）将d个决策变量添加到解集邻域M中，当M出现空集时，则任选d加入M中；
Step4：加入额外扰动因子r。对于M中的决策变量，确定j=1，…，J上的X best ，并在j 中加入新的扰动因子，并假设加入的新的扰动服从B（0，1）的标准正态分布，则加入扰动因 子后j个决策变量的波动范围可表示为：
Figure 822551DEST_PATH_IMAGE015
；
Step5：候候选解集X new 校准更新。计算上一步中生成的解集X new 的目标函数F （X new ），当F（x new ）≥F best ，F best =F（x new ），x best =x new ，更新当前最优解。
Step6：终止条件检查。若i<m，从Step3重新开始循环，并更新迭代次数i=i+1；否则，算法结束。
与DDS相比，P-DDS在上述步骤中加入了计算节点的分配，即在处理器原有N个节点基础上分配N-1个节点用于解集计算（计算节点），1个节点用于初始的计算任务分配与模式输出判断（处理节点），此时对应步骤可更新为：Step1和Step2：均在处理节点中产生，设置Step1中用于分配的模式处理节点数N；Step3和Step4：计算节点进行，此时j=1,2…N-1；Step5：处理节点进行，更新迭代次数i=i+N-1，将新的候选解集方案发送到计算处理器进行评估，每个计算节点评估唯一的F（Xnew）值，确定在所有计算节点中新获得的最佳候选解集。
进一步，所述步骤5中时空模拟效果好的降雨以及步骤4率定出的对应参数组对洪水过程的适用性进行了判断，评价指标选择了针对不同降雨-参数组合的泰勒图以及克林-古普塔效率系数（KGE）、洪峰误差、洪量误差、纳什效率系数（NSE）。KGE已在步骤4中给出，另外3个指标公式可以表示为：
Figure 216623DEST_PATH_IMAGE025
（12）
Figure 34407DEST_PATH_IMAGE021
（13）
Figure 830324DEST_PATH_IMAGE023
（14）
式中：i为时间步长；T为洪水过程的总时长；
Figure 190374DEST_PATH_IMAGE026
、q i 和
Figure 337322DEST_PATH_IMAGE027
分别为模拟、观测和平均流 量（m3/s）；
Figure 958796DEST_PATH_IMAGE022
和q p 分别为模拟、观测的洪峰流量（m3/s）；
Figure 609220DEST_PATH_IMAGE024
和r r 分别为模拟和观测洪量（mm）。
实施例三
下面结合实施实例，对本发明做进一步说明：步骤1、设置WRF中的气象驱动数据与预热时间；步骤2、雨量站融合多普勒天气雷达定量估测降雨（QPE）获得修正降雨数据；步骤3、降雨驱动数据的对比验证；步骤4、基于并行的动态多维搜索算法（P-DDS）的WRF-Hydro模式参数校准；步骤5、验证校准参数后的WRF-Hydro对洪水过程模拟效果。
步骤1，设置WRF中的气象驱动数据与预热时间，这里选择中国北方大清河流域南北支的两个中尺度流域阜平（2210km2）和紫荆关（1760km2），雷达站选择了流域内的石家庄雷达站（如图3所示），WRF采用三层嵌方式（每层嵌套比例为1:3，从外到内分辨率分别为9km，3km，1km，如图4所示）并取最内层区域用于为WRF-Hydro提供包括降雨、辐射、湿度、温度、气压、风速等在内的模式输入。WRF的预热过程中设置2个预热时段：第一个预热时段为降雨开始前的16天，并在第15天时输出Restart文件；第二个预热时间为降雨开始前的24h，第二个预热时段调用Restart文件中的下边界数据与WRF驱动资料（ERA-Interim）中的侧边界数据（温度、辐射、风场等）。
步骤2，雨量站融合多普勒天气雷达定量估测降雨（QPE）获得修正降雨数据，选择了流域内2005-2018年间的10场降雨场次，针对每场降雨过程均分析雷达定量估测降雨（QPE）、基于反距离插值权重平均（RIDW）的雨量站插值降雨、外部漂移克里金（KED）融合雨量站与雷达降雨。3种方式的处理方法：
雷达定量估测降雨（QPE），对多普勒雷达获得的降雨分别进行了电磁波衰减校正、非降雨回波抑制、不同层反射率组合以及雨强-反射率转换。用HB方法来订正雷达回波衰减，设定路径积分衰减系数（Path Integration Attenuation，PIA）PIA的上限为10dB；采用基于三维组合反射率抑制非降雨回波，并设定仰角4°以下的回波为杂波；结合DEM数据和雷达波束标准传播方程获得了不同高度层的组合反射率选取规则，如图5所示；结合流域雨滴谱对不同类型的反射率因子与雨强系数进行调整。
基于反距离插值权重平均（RIDW）的雨量站插值降雨，该方法根据插值区域内各采样点之间的相似度，假设网格对应未知点处的估计值可由局部邻域内所有数据点的距离加权平均值表示，该方法既考虑了常用的泰森多边形中邻点的概念也考虑了趋势面分析中的渐变过程，公式可以表达为：
Figure 408680DEST_PATH_IMAGE036
（1）
式中：Z 0 为对应网格插值点处的估计值；Z i 为第i个雨量站观测数据；d i 为需要插值的网格点与雨量站点间的空间距离；n为研究区雨量站的个数。
研究雷达与雨量站资料融合时，雷达资料被视为唯一的附加变量。KED插值中主变量的期望值与附加变量之间符合如下线性关系：
Figure 511765DEST_PATH_IMAGE002
（2）
式中：G(x)代表主变量的期望，这里G(x)指的是对应x点处的雨量站的观测数据；R (x)为附加变量，这里表示雷达数据得到的降雨估计值；a和b为由假设条件确立的线性方程的系数。
上述公式中给定位置x 0 的估计可由线性估计器计算得到，对应的权重因子表示为：
Figure 874613DEST_PATH_IMAGE003
（3）
Figure 520489DEST_PATH_IMAGE004
（4）
式中：R i (x 0)表示x 0 点经过KED融合后的降雨估算值；R 1 (x i )为雨量站点i处的雨量 值；n为研究区雨量站的个数；
Figure 615484DEST_PATH_IMAGE033
表示在x 0 点处与降雨估计值相关的权重因子，随着估算 点的改变，其值也相应发生改变。
上式与传统克里金法的描述一致，但是在具体权重因子计算上，KED方法增加了额外的约束条件（即与雷达观测有关的数据约束）用于计算雷达场和雨量站场残差的协方差：
Figure 330499DEST_PATH_IMAGE037
（5）
Figure 434722DEST_PATH_IMAGE031
（6）
Figure 935104DEST_PATH_IMAGE008
（7）
式中：R 2 (x 0 )表示x 0 点经过权重处理后得到的雷达降雨估算值；为雷达网格点j处 的雨量值；μ为拉格朗日乘数；C表示残差在不同点间的残差协方差结果，可由雷达估测的降 雨场与估算趋势场（漂移场）的差表示；
Figure 201001DEST_PATH_IMAGE009
为在x 0 处与雷达估计值相关的权重因子，其他 物理意义同上。研究中的漂移场是通过传统克里金方式对雷达测算的雨量站位置处降雨进 行空间变异性处理获得，并通过FFT方法，将基于雨量站点的协方差矩阵C转换为整个流域 网格上的频谱矩阵。
在雨量站位置处插值雷达时，上述约束使用一组权重对雨量站数据进行加权，同时这些权重也会增加估计点x 0 处的雷达值效果。
步骤3，所述步骤3对步骤2中3种基于雷达或雨量站的空间降雨驱动数据精度进行分析，主要体现在：首先针对每种降雨空间数据，均采用雨量站点的留一交叉验证方式，选择了系统偏差BIAS、均方根误差RMSE和均方根转换误差MRTE三个降雨精度评价指标；选择了雨量站点空间散点和降雨数据的累积降雨量空间分布图来对3种降雨数据和WRF降雨的空间降雨情况进行综合评价，如图5-图7所示。
3个指标的计算公式如下：
Figure 403312DEST_PATH_IMAGE010
（8）
Figure 855765DEST_PATH_IMAGE011
（9）
Figure 335288DEST_PATH_IMAGE012
（10）
式中：P t 是雨量站测得的第t时刻的降雨量；n为研究区雨量站的个数；
Figure 896720DEST_PATH_IMAGE013
是不同降 雨驱动数据下对应雨量站位置处的估计值。针对10场降雨，RMSE和MRTE的范围为0到+∞，最 优值为0，BIAS的结果为整个降雨统计时段和所有雨量站的均值，范围在-∞到+∞，最优值 为0。模拟的流域降雨结果如表1和图5所示，可以看出WRF模拟的降雨空间精度较差，单独采 用QPE方式的降雨则量级偏低，而RIDW方式与基于雨量站融合多普勒天气雷达定量估测的 降雨效果较好，但是对于2012年7月21日这种在华北地区造成较大洪水灾害的洪水场次而 言，融合方式的降雨空间效果较好。
Figure 461693DEST_PATH_IMAGE038
表1 组合反射率选取规则
步骤4，采用适用于模式的并行自动调参方法对敏感参数进行率定，其中，敏感参数为容易引起wrf-hydro产汇流过程发生明显变化的参数主要体现在将基于雨量站融合多普勒天气雷达定量估测的降雨、WRF模拟得到的其它气象数据和WRF-Hydro一起用于陆气耦合研究，率定的参数包括土壤下渗系数（R efkdt ）、地表粗糙度（O vrought ）、河道曼宁糙率系数（M ann ）和河道边坡坡度（C hslp ）。WRF-Hydro中采用了降尺度聚解以及升尺度聚合方式对产汇流变量进行网格细分，降尺度比例设定为10，即WRF-Hydro的汇流分辨率为100m，选择克林-古普塔效率系数（Kling-Gupta efficiency，KGE）作为目标函数，KGE可以用下面的公式表示：
Figure 783084DEST_PATH_IMAGE035
（11）
式中：r为模拟系列与观测系列标准偏差的比率；α为模拟系列与观测系列平均值的比率；β则为模拟系列与观测系列的相关系数。当KGE等于1时，说明模拟能力达到最佳（此时表示与理想模拟系列的欧几里得距离最小），而当KGE为负值时，则认为模拟的效果较差。率定时采用了并行的动态多维搜索算法（P-DDS）。P-DDS加入了计算节点的分配，即在处理器原有N个节点基础上分配N-1个节点用于解集计算（计算节点），1个节点用于初始的计算任务分配与模式输出判断（处理节点）。率定的参数结果如表2所示，选择2016年前的8场洪水率定，2016年的2场洪水验证。
Figure 117114DEST_PATH_IMAGE039
表 2 RIDW和KED降雨数据下率定的WRF-Hydro参数组
步骤5，所述步骤5中时空模拟效果好的降雨以及步骤4率定出的对应参数组对洪水过程的适用性进行了判断，评价指标选择了针对不同降雨-参数组合的泰勒图以及克林-古普塔效率系数（KGE）、洪峰误差、洪量误差、纳什效率系数（NSE）。KGE已在步骤4中给出，另外3个指标公式可以表示为：
Figure 849446DEST_PATH_IMAGE040
（12）
Figure 901716DEST_PATH_IMAGE021
（13）
Figure 964481DEST_PATH_IMAGE023
（14）
式中：i为时间步长；T为洪水过程的总时长；
Figure 277651DEST_PATH_IMAGE026
、q i 和
Figure 56251DEST_PATH_IMAGE027
分别为模拟、观测和平均流 量（m3/s）；
Figure 471183DEST_PATH_IMAGE022
和q p 分别为模拟、观测的洪峰流量（m3/s）；
Figure 259010DEST_PATH_IMAGE024
和r r 分别为模拟和观测洪量（mm）， 结果如图8-图10所示。采用融合方式降雨及其率定出的参数组对于洪水指标的模拟效果均 得到了提升，采用基于雨量站融合多普勒天气雷达定量估测降雨的WRF-Hydro率定的关键 参数更适用于WRF/WRF-Hydro的模拟/预报需求。
为了更好地理解本发明，以上结合本发明的具体实施例做了详细描述，但并非是对本发明的限制。凡是依据本发明的技术实质对以上实施例所做的任何简单修改，均仍属于本发明技术方案的范围。本说明书中每个实施例重点说明的都是与其它实施例的不同之处，各个实施例之间相同或相似的部分相互参见即可。对于系统实施例而言，由于其与方法实施例基本对应，所以描述的比较简单，相关之处参见方法实施例的部分说明即可。

Claims (6)

1.一种雨量站融合雷达定量估测降雨的WRF-Hydro关键参数率定方法，包括以下步骤：
步骤1：设置天气预报模式WRF中的气象驱动数据与预热时间，所述气象驱动数据包括温度驱动数据、压强驱动数据、风速驱动数据和降雨驱动数据；
步骤2：雨量站融合多普勒天气雷达定量估测降雨获得修正降雨数据；所述修正降雨数据的处理方法包括以下三种：雷达定量估测降雨法或基于反距离插值权重平均的雨量站插值降雨法或外部漂移克里金融合雨量站与雷达降雨法；
步骤3：降雨驱动数据的适用性判断；对3种基于雷达或雨量站的空间降雨驱动数据精度进行分析，包括以下子步骤：
步骤31：针对每种降雨空间数据，均采用雨量站点的留一交叉验证方式，选择了系统偏差BIAS、均方根误差RMSE和均方根转换误差MRTE三个降雨精度评价指标；
步骤32：选择了雨量站点空间散点和降雨数据的累积降雨量空间分布图来对3种降雨数据和WRF降雨的空间降雨情况进行综合评价；
步骤4：基于并行的动态多维搜索算法的WRF-Hydro模式参数校准，采用适用于模式的并行自动调参方法对敏感参数进行率定，在所述率定的过程中选取克林-古普塔效率系数KGE作为目标函数对模拟的洪水过程进行评估，其中，敏感参数为容易引起WRF-Hydro产汇流过程发生明显变化的参数，
克林-古普塔效率系数KGE的公式为
Figure 379339DEST_PATH_IMAGE001
其中，γ为模拟系列与观测系列标准偏差的比率，α为模拟系列与观测系列平均值的比率，β则为模拟系列与观测系列的相关系数；
步骤4还包括采用并行动态多维搜索算法进行路面水温模式率定过程，并行动态多维搜索算法的串行计算步骤如下：
步骤41：定义算法输入变量；算法输入变量包括邻近扰动距离因子r、算法最大校准次数m、针对D维变量每一个决策变量的初始解的下限X min 和上限值X max 、算法初始解向量为X 0 = [x 1 ，x 2 ，…，x D ]和分配的模式处理节点数N；
步骤42：模式初值分析；迭代次数i从1开始，计算当前初始解对应的目标函数值F（X 0 ），将初值解暂时赋值为最优解，即X best =X 0 ，F best =F（X 0 ）；
步骤43：解集领域概率；随机选择D个决策变量中的J个产生新的解集邻域M，计算决策变量在当前迭代计数下在M中发生的概率：P（i）=1-ln（i）/ln（m），对于d=1，…，D个决策变量，按求得的概率P（i）将d个决策变量添加到解集邻域M中，当M出现空集时，则任选d加入M中；
步骤44：加入额外扰动因子r；对于M中的决策变量，确定j=1，…，J上的X best ，并在j中加入新的扰动因子，并假设加入的新的扰动服从B（0，1）的标准正态分布，则加入扰动因子后j个决策变量的波动范围可表示为：
Figure 466243DEST_PATH_IMAGE002
其中，
Figure 189218DEST_PATH_IMAGE003
为新的参数值，
Figure 423890DEST_PATH_IMAGE004
为当前参数最佳值，
Figure 578928DEST_PATH_IMAGE005
为新扰动因子系数，
Figure 848235DEST_PATH_IMAGE006
为当前参数 范围最大值，
Figure 696105DEST_PATH_IMAGE007
为当前参数范围最小值；
步骤45：候选解集X new 校准更新；计算所述步骤44中生成的解集X new 的目标函数F（X new ），当F（x new ）≥F best ，F best =F（x new ），x best =x new ，更新当前最优解；
步骤46：终止条件检查；当i<m时，返回所述步骤43，并更新迭代次数i=i+1；否则，算法结束；
步骤5：验证校准参数后的WRF-Hydro对洪水过程模拟效果。
2.如权利要求1所述的雨量站融合雷达定量估测降雨的WRF-Hydro关键参数率定方法，其特征在于，所述雷达定量估测降雨法对多普勒雷达获得的降雨分别进行了电磁波衰减校正、非降雨回波抑制、不同层反射率组合以及雨强-反射率转换，采用HB方法来订正雷达回波衰减，设定路径积分衰减系数PIA的上限为10dB；采用基于三维组合反射率抑制非降雨回波，并设定仰角4°以下的回波为杂波；结合DEM数据和雷达波束标准传播方程获得了不同高度层的组合反射率选取规则；结合流域雨滴谱对不同类型的反射率因子与雨强系数进行调整。
3.如权利要求2所述的雨量站融合雷达定量估测降雨的WRF-Hydro关键参数率定方法，其特征在于，所述基于反距离插值权重平均的雨量站插值降雨法根据插值区域内各采样点之间的相似度，设定网格对应未知点处的估计值由局部邻域内所有数据点的距离加权平均值表示，公式为
Figure 418074DEST_PATH_IMAGE008
其中，Z 0 为对应网格插值点处的估计值，Z i 为第i个雨量站观测数据，d i 为需要插值的网格点与雨量站点间的空间距离，n为研究区雨量站的个数。
4.如权利要求3所述的雨量站融合雷达定量估测降雨的WRF-Hydro关键参数率定方法，其特征在于，所述外部漂移克里金融合雨量站与雷达降雨法中的外部漂移克里金插值中主变量的期望值与附加变量之间符合如下线性关系：
Figure 376802DEST_PATH_IMAGE009
其中，主变量的期望G(x)为对应x点处的雨量站的观测数据，附加变量R(x)为雷达数据得到的降雨估计值；a和b为由假设条件确立的线性方程的系数。
5.如权利要求4所述的雨量站融合雷达定量估测降雨的WRF-Hydro关键参数率定方法，其特征在于，对于给定位置x 0的估计由线性估计器计算得到，对应的权重因子表示为
Figure 985769DEST_PATH_IMAGE010
Figure 4541DEST_PATH_IMAGE011
其中，
Figure 213805DEST_PATH_IMAGE012
为在x 0点处与降雨估计值相关的权重因子，n为研究区雨量站的个数，R 1 (x i ) 为雨量站点i处的雨量值，R(x 0)表示x 0点经过KED融合后的降雨估算值。
6.如权利要求5所述的雨量站融合雷达定量估测降雨的WRF-Hydro关键参数率定方法，其特征在于，在所述外部漂移克里金融合雨量站与雷达降雨法中增加额外的约束条件用于计算雷达场和雨量站场残差的协方差：
Figure 772963DEST_PATH_IMAGE013
Figure 688966DEST_PATH_IMAGE014
Figure 675377DEST_PATH_IMAGE015
其中，C表示残差在不同点间的残差协方差结果，x i 为雨量站点i处，x j 为雨量站点j处，μ 0和μ 1为拉格朗日乘数，R 2 (x j )为雷达网格点j处的雨量值，R 2 (x 0 )表示x 0 点经过权重处理后 得到的雷达降雨估算值，
Figure 309620DEST_PATH_IMAGE016
为在x 0 处与雷达估计值相关的权重因子。
CN202110841644.3A 2021-07-26 2021-07-26 一种雨量站融合雷达定量估测降雨的WRF-Hydro关键参数率定方法 Active CN113281754B (zh)

Priority Applications (1)

Application Number Priority Date Filing Date Title
CN202110841644.3A CN113281754B (zh) 2021-07-26 2021-07-26 一种雨量站融合雷达定量估测降雨的WRF-Hydro关键参数率定方法

Applications Claiming Priority (1)

Application Number Priority Date Filing Date Title
CN202110841644.3A CN113281754B (zh) 2021-07-26 2021-07-26 一种雨量站融合雷达定量估测降雨的WRF-Hydro关键参数率定方法

Publications (2)

Publication Number Publication Date
CN113281754A CN113281754A (zh) 2021-08-20
CN113281754B true CN113281754B (zh) 2021-10-01

Family

ID=77281330

Family Applications (1)

Application Number Title Priority Date Filing Date
CN202110841644.3A Active CN113281754B (zh) 2021-07-26 2021-07-26 一种雨量站融合雷达定量估测降雨的WRF-Hydro关键参数率定方法

Country Status (1)

Country Link
CN (1) CN113281754B (zh)

Families Citing this family (3)

* Cited by examiner, † Cited by third party
Publication number Priority date Publication date Assignee Title
CN113640900A (zh) * 2021-08-13 2021-11-12 中国科学院西北生态环境资源研究院 雨雪量测量方法、装置、设备及存储介质
CN113887847B (zh) * 2021-12-08 2022-03-11 中国水利水电科学研究院 一种基于WRF-Hydro模型的混合产流区次洪预报方法
CN116679357B (zh) * 2023-06-05 2024-04-19 广东省水文局梅州水文分局 一种对雷达定量降水估测进行水量平衡处理的方法

Citations (4)

* Cited by examiner, † Cited by third party
Publication number Priority date Publication date Assignee Title
CN102103195A (zh) * 2009-12-18 2011-06-22 东软飞利浦医疗设备系统有限责任公司 一种宽频带数字磁共振射频接收实现装置及方法
CN106485366A (zh) * 2016-10-31 2017-03-08 武汉大学 一种复杂梯级水库群蓄水期优化调度方法
CN108051786A (zh) * 2017-10-30 2018-05-18 北京航天福道高技术股份有限公司 一种宽带目标模拟器验证平台及验证方法
CN111625993A (zh) * 2020-05-25 2020-09-04 中国水利水电科学研究院 一种基于山区地形及降雨特征预测的小流域面雨量插值方法

Family Cites Families (13)

* Cited by examiner, † Cited by third party
Publication number Priority date Publication date Assignee Title
US4649388A (en) * 1985-11-08 1987-03-10 David Atlas Radar detection of hazardous small scale weather disturbances
DE10315012B4 (de) * 2003-04-02 2005-05-12 Eads Deutschland Gmbh Verfahren zur Linearisierung von FMCW-Radargeräten
US8898436B2 (en) * 2009-04-20 2014-11-25 Oracle America, Inc. Method and structure for solving the evil-twin problem
JPWO2015182142A1 (ja) * 2014-05-28 2017-06-15 メトロウェザー株式会社 気象予報システム
CN105842692B (zh) * 2016-03-17 2018-07-03 中国科学院遥感与数字地球研究所 一种insar测量中的大气校正方法
CN107403073B (zh) * 2017-10-03 2020-07-21 中国水利水电科学研究院 一种基于数据同化改进预报降雨的集合洪水预报方法
CN107942413A (zh) * 2017-10-27 2018-04-20 中国水利水电科学研究院 一种降雨均匀度测试装置及测试方法
CN108761576B (zh) * 2018-05-28 2020-11-13 国网山西省电力公司电力科学研究院 一种x波段气象雷达与雨量站数据融合方法及系统
CN110222783B (zh) * 2019-06-13 2023-04-07 南京信息工程大学 基于小波域正则化的地基和星载雷达降水数据融合方法
CN111915158A (zh) * 2020-07-15 2020-11-10 云南电网有限责任公司带电作业分公司 一种基于Flood Area模型的暴雨灾害天气风险评估方法、装置及设备
CN112213727B (zh) * 2020-10-15 2024-01-02 国家卫星气象中心(国家空间天气监测预警中心) 一种基于主被动微波联合探测的星载雷达的降水订正方法
CN112651118B (zh) * 2020-12-21 2023-07-28 中国科学院地理科学与资源研究所 一种气候-陆面-水文过程全耦合模拟方法
CN112711083B (zh) * 2021-01-07 2022-07-05 国网福建省电力有限公司 基于自适应权重特征的多源降水数据动态融合方法及系统

Patent Citations (4)

* Cited by examiner, † Cited by third party
Publication number Priority date Publication date Assignee Title
CN102103195A (zh) * 2009-12-18 2011-06-22 东软飞利浦医疗设备系统有限责任公司 一种宽频带数字磁共振射频接收实现装置及方法
CN106485366A (zh) * 2016-10-31 2017-03-08 武汉大学 一种复杂梯级水库群蓄水期优化调度方法
CN108051786A (zh) * 2017-10-30 2018-05-18 北京航天福道高技术股份有限公司 一种宽带目标模拟器验证平台及验证方法
CN111625993A (zh) * 2020-05-25 2020-09-04 中国水利水电科学研究院 一种基于山区地形及降雨特征预测的小流域面雨量插值方法

Non-Patent Citations (4)

* Cited by examiner, † Cited by third party
Title
" Evaluation of Three Satellite-Based Precipitation Products Over the Lower Mekong River Basin Using Rain Gauge Observations and Hydrological Modeling";Yishan Li 等;《IEEE Journal of Selected Topics in Applied Earth Observations and Remote Sensing》;20191231;第2357-2373页 *
"imaging radiometer wind speed and rain rate retrieval: Part-2. Analysis of retrieval accuracy";Ruba A. Amarin 等;《2010 11th Specialist Meeting on Microwave Radiometry and Remote Sensing of the Environment》;20101231;第279-369页 *
"基于不同降水产品的WRF-Hydro模式径流模拟——以漳河流域为例";高玉芳 等;《热带气象学报》;20200630;第299-306页 *
"珠江流域分布式水文模拟及水库压咸优化调度研究 ";孙甲岚;《中国博士学位论文全文数据库 (工程科技Ⅱ辑)》;20151130;全文 *

Also Published As

Publication number Publication date
CN113281754A (zh) 2021-08-20

Similar Documents

Publication Publication Date Title
CN113281754B (zh) 一种雨量站融合雷达定量估测降雨的WRF-Hydro关键参数率定方法
Stathopoulos et al. Wind power prediction based on numerical and statistical models
Blanchet et al. Mapping snow depth return levels: smooth spatial modeling versus station interpolation
Yuan et al. Validation of China-wide interpolated daily climate variables from 1960 to 2011
Yakir et al. Hydrologic response of a semi-arid watershed to spatial and temporal characteristics of convective rain cells
CN112711899B (zh) 一种蒸发波导高度的融合预测方法
CN112785024B (zh) 一种基于流域水文模型的径流计算和预测方法
CN111027175A (zh) 基于耦合模型集成模拟的洪水对社会经济影响的评估方法
CN114741987B (zh) 考虑洪水预报模型绝对误差拟合残差分布的洪水概率预报模型
Atencia et al. Effect of radar rainfall time resolution on the predictive capability of a distributed hydrologic model
CN112800636A (zh) 一种估算无资料地区流域地表水资源量的方法及系统
CN115544785B (zh) 一种无资料梯级水库流域水文模拟方法和系统
Rendon et al. Continuous forecasting and evaluation of derived Z-R relationships in a sparse rain gauge network using NEXRAD
CN111915158A (zh) 一种基于Flood Area模型的暴雨灾害天气风险评估方法、装置及设备
Fuentes et al. Spatial–temporal mesoscale modeling of rainfall intensity using gage and radar data
CN114819322A (zh) 湖泊入湖流量的预报方法
CN110186533A (zh) 一种高精度的河口短期潮位预报方法
Larabi et al. Using functional data analysis to calibrate and evaluate hydrological model performance
CN115575914B (zh) 一种多波段双线偏振天气雷达观测量误差量化方法
Bakiş et al. Analysis and comparison of spatial rainfall distribution applying different interpolation methods in Porsuk river basin, Turkey
Fischer et al. Seasonal cycle in German daily precipitation extremes
CN113887847B (zh) 一种基于WRF-Hydro模型的混合产流区次洪预报方法
Li et al. A geostatistical method for Texas NexRad data calibration
Guo et al. Daily runoff simulation in Poyang Lake Intervening Basin based on remote sensing data
Gao et al. Radar-rainfall estimation from S-band radar and its impact on the runoff simulation of a heavy rainfall event in the Huaihe River Basin

Legal Events

Date Code Title Description
PB01 Publication
PB01 Publication
SE01 Entry into force of request for substantive examination
SE01 Entry into force of request for substantive examination
GR01 Patent grant
GR01 Patent grant