CN112666559B - 一种多普勒天气雷达数据质量控制方法 - Google Patents

Abstract

本发明公开了一种多普勒天气雷达数据质量控制方法，包括如下步骤：基于多普勒天气雷达历史基数据资料，获得气象回波数据分布规律和地物杂波数据分布规律；基于气象回波数据分布规律和地物杂波数据分布规律，利用偏度方法识别地物杂波；根据数据分布规律和偏度方法设计雷达地物杂波滤波器，利用该雷达地物杂波滤波器进行滤波，获得去除地物杂波后的雷达滤波数据；基于多普勒天气雷达原始基数据和雷达滤波数据，确定雷达波束阻挡区的雷达波束阻挡率；进行雷达滤波数据修正。该方法可以有效识别、去除雷达基数据中的地物杂波，并对滤波后的雷达数据进行填补修正，实现对雷达基数据的质量控制，从而提高气象降水估测的准确率。

Description

一种多普勒天气雷达数据质量控制方法

技术领域

本发明属于气象分析技术领域，尤其涉及一种多普勒天气雷达数据质量控制方法。

背景技术

多普勒天气雷达以其高时空分辨率的特点成为气象业务中监测对流系统的重要工具。天气雷达在正常工作时不仅能接收气象回波信号，同时也能接收非气象回波信号，如高山或建筑物等产生的回波；诸类障碍物的反射率因子高，往往会导致降雨强度被高估；而被山区或其他障碍物遮挡的区域，降雨强度则会被低估。因此雷达资料的质量控制对于雷达能有效使用显得尤为重要。目前在雷达信号处理部分已应用滤波器消除非气象回波，但仍有杂波残留，从而影响雷达定量估测降水结果。若要提高雷达定量估测降水的精度，首要问题是提高雷达数据质量，合理剔除非气象回波；目前影响降水估测精度最严重的非气象回波为地物杂波。为此，人们研究了多种地物杂波的识别和去除方法。Hynek等对地面上的静止地物杂波开展了研究，为消除地物杂波，一般选用短时间内、晴朗的天气条件下进行实验，因为在有降水的天气过程中进行实验，降水回波和地物杂波会混在一起，影响杂波去除的效果。选择短期内无降水的时间段，在杂波测量过程中，需对仰角保持均匀的高度扫描，因此，精确的仰角选择非常重要。Steiner等通过雷达反射率三维结构特征，根据决策树判别法来评估雷达回波的垂直梯度和空间的差异性，自动来识别被气象回波包围的地物杂波。Tabary等发现在多普勒天气雷达扫描时，在某些特定地区经常会接收到一些奇异的地物杂波，进一步研究发现受影响的地物杂波为风力发电机，针对该特定地物杂波，开展了相应的识别和剔除方法研究工作。Kucera等分析了10m、30m和90m分辨率的地形数据，并对美国关岛雷达波束遮挡进行模拟实验，指出地形特别复杂的区域，分辨率对波束遮挡模拟结果的影响较大。

国内，气象学者围绕地物回波特征、地物阻挡情况等方面也开展了相应的研究工作。刘黎平等基于模糊逻辑的原理，通过单部雷达反射率因子强度纹理、垂直变化、径向方向变化及地物回波与地形海拔高度的关系等特征，计算径向速度区域平均值和方差、速度谱宽区域平均值等物理量，并确定了模糊逻辑识别超折射地物杂波的隶属函数关系，通过调整地物杂波判别阈值，分步识别并剔除超折射地物。王红艳等根据雷达波束对远近距离地物阻挡影响情况，进行了模拟和实测检验。文浩等利用统计和典型个例分析方法对地物回波及电磁干扰质控算法进行了初步业务试运行的效果评估，验证了算法的可靠性，但该算法仍需人工干预。也有学者着手开展雷达径向杂波的质控方法；吕博等尝试利用偏振雷达观测的层状云和对流云降水资料进行偏振参量的质量分析。

目前，国内大部分业务使用的多普勒天气雷达，都不能有效的去除地物杂波；且目前国内外去除雷达地物杂波的方法基本都是针对小样本的以成片地物杂波或是特定某一径向杂波为主要研究对象开展研发工作，并且每部雷达地物差异性较大，相同算法杂波去除率仍有限。因此，需要对长序列、大样本的雷达数据逐网格识别和剔除地物杂波的方法开展研究。

发明内容

为解决上述技术问题，本发明的目的在于提供一种多普勒天气雷达数据质量控制方法，该方法使用大量基数据历史资料，通过雷达气候学的方法，发现气象回波与地物杂波的分布规律，并基于该分布规律可以利用数学统计方法识别、去除雷达基数据中的地物杂波，并对滤波后的雷达数据进行填补修正，实现对雷达基数据的质量控制，从而提高气象降水估测的准确率。

为实现上述技术目的，达到上述技术效果，本发明通过以下技术方案实现：

一种多普勒天气雷达数据质量控制方法，包括如下步骤：

S1，基于多普勒天气雷达历史基数据资料，获得气象回波数据分布规律和地物杂波数据分布规律；

其中，气象回波数据分布规律为：气象回波呈现较低的反射率因子分布，为正偏度；地物杂波数据分布规律为：地物杂波呈现较高的反射率因子分布，为负偏度；

S2，基于气象回波数据分布规律和地物杂波数据分布规律，利用偏度方法识别地物杂波；

在该步骤S2中，将雷达回波强度数据采用强度间隔设计方法分割为多档，根据分割后的数据，分别计算雷达回波强度数据的平均值、标准方差值和偏度值，根据偏度值识别需要进行处理的多普勒天气雷达原始基数据中的地物杂波；

S3，根据步骤S1所得的数据分布规律和步骤S2的偏度方法设计雷达地物杂波滤波器，利用该雷达地物杂波滤波器进行滤波，获得去除地物杂波后的雷达滤波数据；

S4，基于多普勒天气雷达原始基数据和雷达滤波数据，确定雷达波束阻挡区的雷达波束阻挡率；

在该步骤S4中，雷达波束阻挡率POD值的计算公式为：

式中，P为多普拉天气雷达原始基数据中的基于波束阻挡情况下的某一个径向圈上实际测得的反射率因子频数，P＇为雷达滤波数据中的基于波束阻挡情况下的对应径向圈上的反射率因子频数，P₉₉为多普拉天气雷达原始基数据中的基于波束阻挡情况下的99百分位的地物杂波频数，其仅次于地物杂波最大极值情况下的频数；

S5，雷达滤波数据修正

分别针对低层多级仰角的多普勒天气雷达原始基数据，按照步骤S3进行地物杂波识别，并按照步骤S4和步骤S5的方法获得对应的雷达滤波数据和雷达波束阻挡率；根据地物杂波识别过程中的偏度值以及雷达波束阻挡率，利用上级仰角偏度值为正的回波强度替代下级仰角偏度值为负的对应格点的回波强度的方法，对部分雷达波束阻挡区域进行填补，获得修正后的多普勒天气雷达数据。

进一步的，步骤S2中，将雷达回波强度数据0-60dBZ按照每2.5dBZ间隔一档，共24档。

更进一步的，步骤S2中，雷达回波强度数据的平均值的计算公式为式(1)：

式(1)中，N代表24档，X_i为某一档的雷达回波强度值。

更进一步的，步骤S2中，雷达回波强度数据的标准方差值的计算公式为式(2)：

式(2)中，N为24档，X_i为某一档的雷达回波强度值，μ为式(1)求得的平均值。

更进一步的，步骤S2中，雷达回波强度数据的偏度值的计算公式为式(3)：

式(3)中X_i为某一档的雷达回波强度值，N代表24档，μ为式(1)求得的平均值，σ为式(2)求得的标准方差值。

进一步的，步骤S5中，低层多级仰角分别包括0.5°仰角、1.5°仰角和2.5°仰角。

进一步的，该控制方法还包括对流风暴数据控制，具体为：

根据对流风暴回波类型的空间分布特征，识别修正后的多普勒天气雷达数据中的对流风暴回波，并提取各高度层对流风暴回波中回波强度的最大值，分层统计其50百分位和90百分位的回波强度值，以获得对流风暴的垂直结构分布特征；然后提取对流风暴回波1km×1km格点上的所有回波强度，获得每个格点上的50百分位和90百分位的回波强度值。

更进一步的，对流风暴回波类型的空间分布特征是基于多普勒天气雷达历史基数据资料获得。

本发明的有益效果是：

本发明基于大量的多普勒天气雷达历史基数据资料，通过雷达气候学的方法，发现气象回波与地物杂波的分布规律，利用数学统计方法获得的偏度值识别地物杂波，并根据偏度方法设计雷达地物杂波滤波器，以有效去除雷达原始基数据中的地物杂波，然后根据滤波后的雷达数据和雷达原始基数据，获得地物阻挡情况下的雷达波束阻挡率；根据雷达波束阻挡率和地物杂波识别过程中的偏度值，并结合各层仰角雷达回波强度，利用替代方法对部分雷达波束阻挡区域进行填补，获得修正后的多普勒天气雷达数据。

此外，本发明还可对对流风暴雷达数据进行控制处理，具体提取各高度层对流风暴回波中回波强度的最大值，分层统计其50百分位和90百分位的回波强度值，以获得对流风暴回波的垂直结构分布特征；然后提取对流风暴回波1km×1km格点上的所有回波强度，获得50百分位和90百分位的回波强度值。从而，根据对流风暴回波的垂直结构分布特征以及1km×1km格点上的回波强度分布，可以分析和估测强对流天气。

附图说明

图1为本发明一种多普勒天气雷达数据质量控制方法的流程示意图。

图2为本发明实施例中的2013年常州雷达站多普勒雷达0.5°仰角(格点为(2，31°))气象回波分布图。

图3为本发明实施例中的2013年常州雷达站多普勒雷达0.5°仰角(格点为(4，3°))地物杂波分布图。

具体实施方式

下面将结合具体实施例对本发明中的技术方案进行清楚、完整地描述，显然，所描述的实施例仅仅是本发明一部分实施例，而不是全部的实施例。基于本发明中的实施例，本领域普通技术人员在没有做出创造性劳动前提下所获得的所有其他实施例，都属于本发明保护的范围。

实施例

如图1所示的流程示意图，本发明实施例的多普勒天气雷达数据质量控制方法，包括如下步骤：

S1，基于多普勒天气雷达历史基数据资料，获得气象回波数据分布规律和地物杂波数据分布规律；

该实施例的多普勒天气雷达历史基数据资料采用的是2013年常州站多普勒天气雷达基数据，其存储量为552GB。其中主汛期(3-8月)雷达24h开机，存储量基本为70GB/月，非汛期每天开机时间为9:00-15:00，部分时间段因雷达停机维修等原因，存在少量资料的缺失。研究雷达扫描范围内(径向距离：0-230km，方位：0-360°，仰角：0.5°-19.4°)所有格点(745 200)的雷达资料分布特征。

结合晴空下一次体扫地物杂波分布，对上述常州站多普勒雷达每个格点、每个仰角2013年反射率因子强度进行分析，可知地物杂波基本集中在雷达100km范围内的高反射率因子值中，并且获得气象回波和地物回波的分布特征，即气象回波呈现较低的反射率因子分布，为正偏度，地物杂波具有较高的反射率因子分布，为负偏度。如图2所示为2013年常州雷达站多普勒雷达0.5°仰角(格点为(2，31°))气象回波分布图；如图3所示为2013年常州雷达站多普勒雷达0.5°仰角(格点为(4，3°))地物杂波分布图。

S2，基于气象回波数据分布规律和地物杂波数据分布规律，利用偏度方法识别地物杂波；

具体的，可以根据步骤S1获得的气象回波数据分布规律和地物杂波数据分布规律，应用数学统计方法获得的偏度值识别所需要处理的多普勒天气雷达原始基数据中的地物杂波；

具体过程为：将雷达回波强度数据采用强度间隔设计方法分割为多档；在本实施例中，将雷达回波强度数据0-60dBZ按照每2.5dBZ间隔一档，共24档；根据分割后的数据，分别计算雷达回波强度数据的平均值、标准方差值和偏度值，根据偏度值识别需要进行处理的多普勒天气雷达原始基数据中的地物杂波；

雷达回波强度数据的平均值的计算公式为式(1)：

式(1)中，N代表24档，X_i为某一档的雷达回波强度值。

雷达回波强度数据的标准方差值的计算公式为式(2)：

式(2)中，N为24档，X_i为某一档的雷达回波强度值，μ为式(1)求得的平均值。

雷达回波强度数据的偏度值的计算公式为式(3)：

式(3)中X_i为某一档的雷达回波强度值，N代表24档，μ为式(1)求得的平均值，σ为式(2)求得的标准方差值；

此外，雷达回波强度数据的平均值的计算公式还可以为式(4)：

式(4)中，count为某一档的雷达回波强度出现的频数；

为0-60DBZ雷达回波强度之间的24档出现的频数之和。

雷达回波强度数据的标准方差值的计算公式还可以为式(5)：

式(5)中，count为某一档的雷达回波强度出现的频数，

为0-60DBZ雷达回波强度之间的24档出现的频数之和，mean为式(4)获得的平均值。

S3，根据步骤S1所得的数据分布规律和步骤S2的偏度方法设计雷达地物杂波滤波器，利用该雷达地物杂波滤波器对所要处理的多普勒天气雷达原始基数据进行滤波，获得去除地物杂波后的雷达滤波数据；

S4，基于多普勒天气雷达原始基数据和雷达滤波数据，确定雷达波束阻挡区的雷达波束阻挡率；

统计长时间段内雷达连续观测的数据，并经过上述步骤滤波后，得到雷达探测范围内不同仰角上回波随方位-斜距的二维概率分布，即回波探测概率分布。对于波束遮挡的区域，气象回波的回波探测概率分布沿方位方向突然减小，雷达探测到气象回波的概率相对越小，波束被遮挡的比率越大，在遮挡导致的观测盲区，气向回波的回波探测概率为零。根据雷达观测数据得到的波束遮挡分布更能体现地形和大气条件的本地特征。

以下利用POD值表征雷达波束阻挡率，以更准确的了解地形和大气条件的本地特征。

在该步骤S4中，雷达波束阻挡率POD值的计算公式为：

式中，P为多普拉天气雷达原始基数据中的基于波束阻挡情况下的某一个径向圈上实际测得的反射率因子频数，P＇为雷达滤波数据中的基于波束阻挡情况下的对应径向圈上的反射率因子频数，P₉₉为多普拉天气雷达原始基数据中的基于波束阻挡情况下的99百分位的地物杂波频数，其仅次于地物杂波最大极值情况下的频数；

当POD值(雷达波束阻挡率)越接近1，说明该波束阻挡区域基本为地物杂波，如果出现的数值接近0，则说明波束阻挡较少，基本为气象回波。

S5，雷达滤波数据修正

分别针对低层多级仰角的多普勒天气雷达原始基数据，按照步骤S3进行地物杂波识别，并按照步骤S4和步骤S5的方法获得对应的雷达滤波数据和雷达波束阻挡率；根据地物杂波识别过程中的偏度值以及雷达波束阻挡率，利用上级仰角偏度值为正的回波强度替代下级仰角偏度值为负的对应格点的回波强度的方法，对部分雷达波束阻挡区域进行填补，获得修正后的多普勒天气雷达数据。

在本实施例中，分别对低层三个仰角(0.5°、1.5°、2.5°)的雷达数据用偏度值识别地物杂波，并进行比较：首先比较0.5°仰角和1.5°仰角，相对于同一格点，如果0.5°仰角偏度值为负值，同时1.5°偏度值为正值，则用1.5°格点上的回波强度数据代替0.5°仰角对应格点上的回波强度数据值；如果0.5°和1.5°仰角对应的偏度值都为负值，则用2.5°仰角格点上的回波强度数据代替0.5°仰角对应格点上的回波强度数据。根据上述方法，对于波束阻挡率大于某一值的部分区域进行填补修正。

该控制方法还包括对对流风暴的数据控制，具体为：

对多普勒天气雷达历史基数据资料按照步骤S1-S5进行处理，获得了CA PPI格点全省雷达拼图资料，不同雷达之间可以相互补充，获得了对流风暴频数在江苏省内地区的空间分布月变化及其日变化特征；在此基础上，将对流风暴按回波顶高和面积大小各分为3类，获得了各类型对流风暴回波的空间分布月变化及其日变化特征；

根据对流风暴回波类型的空间分布特征，识别修正后的多普勒天气雷达数据中的对流风暴回波，并提取各高度层对流风暴回波中回波强度的最大值，分层统计其50百分位和90百分位的回波强度值，以获得对流风暴的垂直结构分布特征；然后提取对流风暴回波1km×1km格点上的所有回波强度，获得每个格点上的50百分位和90百分位的回波强度值。

根据对流风暴的垂直结构分布特征以及1km×1km格点上的回波强度分布，可以分析和估测强对流天气。

Claims (4)

1.一种多普勒天气雷达数据质量控制方法，其特征在于，包括如下步骤：

S1，基于多普勒天气雷达历史基数据资料，获得气象回波数据分布规律和地物杂波数据分布规律；

其中，气象回波数据分布规律为：气象回波呈现较低的反射率因子分布，为正偏度；地物杂波数据分布规律为：地物杂波呈现较高的反射率因子分布，为负偏度；

S2，基于气象回波数据分布规律和地物杂波数据分布规律，利用偏度方法识别地物杂波；

在该步骤S2中，将雷达回波强度数据采用强度间隔设计方法分割为多档，根据分割后的数据，分别计算雷达回波强度数据的平均值、标准方差值和偏度值，根据偏度值识别需要进行处理的多普勒天气雷达原始基数据中的地物杂波；

在该步骤S2中，将雷达回波强度数据0-60dBZ按照每2.5dBZ间隔一档，共24档；

其中，雷达回波强度数据的平均值的计算公式为式(1)或式(4)：

式(1)中，N代表24档，X_i为某一档的雷达回波强度值；

式(4)中，count为某一档的雷达回波强度出现的频数；

为0-60DBZ雷达回波强度之间的24档出现的频数之和；

雷达回波强度数据的标准方差值的计算公式为式(2)或式(5)：

式(2)中，N为24档，X_i为某一档的雷达回波强度值，μ为式(1)求得的平均值；

式(5)中，count为某一档的雷达回波强度出现的频数，

count为0-60DBZ雷达回波强度之间的24档出现的频数之和，mean为式(4)获得的平均值；

雷达回波强度数据的偏度值的计算公式为式(3)：

式(3)中X_i为某一档的雷达回波强度值，N代表24档，μ为式(1)或式(4)求得的平均值，σ为式(2)或式(5)求得的标准方差值；

S3，根据步骤S1所得的数据分布规律和步骤S2的偏度方法设计雷达地物杂波滤波器，利用该雷达地物杂波滤波器进行滤波，获得去除地物杂波后的雷达滤波数据；

S4，基于多普勒天气雷达原始基数据和雷达滤波数据，确定雷达波束阻挡区的雷达波束阻挡率；

在该步骤S4中，雷达波束阻挡率POD值的计算公式为：

式中，P为多普拉天气雷达原始基数据中的基于波束阻挡情况下的某一个径向圈上实际测得的反射率因子频数，P＇为雷达滤波数据中的基于波束阻挡情况下的对应径向圈上的反射率因子频数，P₉₉为多普拉天气雷达原始基数据中的基于波束阻挡情况下的99百分位的地物杂波频数，其仅次于地物杂波最大极值情况下的频数；

S5，雷达滤波数据修正

分别针对低层多级仰角的多普勒天气雷达原始基数据，按照步骤S3进行地物杂波识别，并按照步骤S4和步骤S5的方法获得对应的雷达滤波数据和雷达波束阻挡率；根据地物杂波识别过程中的偏度值以及雷达波束阻挡率，利用上级仰角偏度值为正的回波强度替代下级仰角偏度值为负的对应格点的回波强度的方法，对部分雷达波束阻挡区域进行填补，获得修正后的多普勒天气雷达数据。

2.根据权利要求1所述的一种多普勒天气雷达数据质量控制方法，其特征在于，步骤S5中，低层多级仰角分别包括0.5°仰角、1.5°仰角和2.5°仰角。

3.根据权利要求1所述的一种多普勒天气雷达数据质量控制方法，其特征在于，该控制方法还包括对流风暴数据控制，具体为：

根据对流风暴回波类型的空间分布特征，识别修正后的多普勒天气雷达数据中的对流风暴回波，并提取各高度层对流风暴回波中回波强度的最大值，分层统计其50百分位和90百分位的回波强度值，以获得对流风暴的垂直结构分布特征；然后提取对流风暴回波1km×1km格点上的所有回波强度，获得每个格点上的50百分位和90百分位的回波强度值。

4.根据权利要求1所述的一种多普勒天气雷达数据质量控制方法，其特征在于，对流风暴回波类型的空间分布特征是基于多普勒天气雷达历史基数据资料获得。

Images