基于双线偏振参量对北京X波段雷达非降水回波识别方法的研究*

夏 凡 张乐坚 张 林

1 山东省气象科学研究所，济南 250031 2 中国气象局气象探测中心，北京 100081

提 要：为了提高中国X波段双线偏振雷达(简称双偏振雷达)的数据质量，提升应用水平，基于美国多雷达多传感器定量降水估测业务系统质量控制(简称dpQC)方案，在此基础上增加了差分反射率阈值判别并对方案中的计算参数进行本地化调整，设计了针对X波段双偏振雷达的质量控制(简称dpXQC)方案，原理主要包括：将双偏振参量相关系数(CC)区分降水回波与非降水回波的阈值调整为0.9，小于0.9的回波判别为非降水回波；通过计算回波顶高与强回波区位置对冰雹与非均一性波束填充区域进行标识；通过统计探空资料0℃高度上下相关系数的特征，调整了融化层的识别参数；优化算法将双偏振参量差分反射率绝对值大于5 dB判别为非降水回波；加入了条状杂波识别、连续性检验与斑块杂波识别算法。利用dpXQC方案对北京五部X波段多普勒双偏振业务雷达进行试验，检验评估其对非降水回波识别效果，获得以下结果与结论：大部分非降水回波可以通过CC阈值判别被滤除；对于CC阈值判别无法滤除的条幅状、孤立点、块状非降水回波，dpXQC方案利用条幅状杂波识别、连续性检验与差分反射率阈值判别去滤除；dpXQC方案有效保护了冰雹、非均一性波束填充与融化层区域中CC<0.9的降水回波不被滤除；检验结果显示，dpXQC方案对非降水回波正确识别率为91.8%，错误识别率为20.6%，均要优于dpQC方案，但dpXQC方案的错误识别率依然很高，主要有两方面原因，一是由于X波段雷达反射率衰减严重，方案会对一些离雷达较远的非均一性填充区域出现遗漏识别，二是融化层识别算法在所有方位角的位置与厚度是固定的，与实际观测不符，这均会导致被遗漏标识区域中CC<0.9的降水回波被误判为非降水回波。总体来看，dpXQC方案在滤除非降水回波的同时，最大程度保护了降水回波不被滤除，为雷达资料定量化应用等相关业务提供保障。

引 言

我国目前正逐步将CINRAD单线偏振雷达(简称单偏振雷达)升级改造为双线偏振雷达(简称双偏振雷达)，与单偏振雷达相比，双偏振雷达可以接收水平和垂直方向的回波信号，除了反射率(ZH)、径向速度(Vr)、谱宽(SPW)，双偏振雷达还可以探测接收零滞后相关系数(CC)、差分反射率(ZDR)、差分传播相移(ФDP)等双偏振参量。双偏振雷达资料在强对流天气预警(冯晋勤等，2018)、定量降水估测(Zhang et al, 2016；陈超等，2019)、水凝物粒子分类(冯亮等，2018)领域发挥了重要的作用，我国自主研发的GRAPES(Global and Regional Assimilation and Prediction System)数值预报模式(Chen et al，2015)，也正尝试将双偏振雷达数据同化进模式初始场中，提高强对流天气的预报效果。然而这些不可避免地会受到非降水回波的影响，高质量的数据是保证雷达自动化应用的关键，因此雷达基数据在应用前需要进行质量控制，识别并滤除非降水回波。

模糊逻辑法(Kessinger et al，2003；刘黎平等，2007；江源等，2009)与神经网络法(Lakshmanan et al，2007)是较为成熟的雷达质量控制算法，在区分降水回波与非降水回波中取得了较好的效果。随着双偏振参量的引入，研究人员可以在原有算法的基础上增加更多判别条件。杜牧云等(2013a)在模糊逻辑算法中增加了CC、ZDR与ФDP纹理隶属成员函数，改善了算法对零速度区降水回波错误识别问题，同时提升了混合型回波的识别效果，但是仍遗漏了不少地物杂波；Lakshmanan et al(2013)在原有神经网络算法的基础上，训练数据中增加了CC、ZDR与ФDP参量来区分降水与非降水回波，结果显示识别技巧要高于原有的神经网络法，对混合型回波的识别效果较好，然而该方法质量控制效果依赖于搭建网络模型的训练数据集。Tang et al(2014)基于ZH、CC与探空温度数据，针对双偏振雷达设计了一套质量控制(dual polarization quality control,dpQC)方案，其核心是利用CC识别大部分非降水回波与降水回波，然后利用回波顶高等参数处理CC无法区分的回波。统计结果显示，dpQC方案的质量控制技巧高于神经网络法，而且计算时间远少于神经网络方法。这一方法已经在美国多雷达多传感器(MRMS)定量降水估测(QPE)系统中业务应用。在国内，张林和杨洪平(2018)借鉴了dpQC算法流程，分类识别上海南汇业务双偏振雷达的降水回波与非降水回波。还有一些研究人员单纯对一个双偏振参量进行质量控制(杜牧云等，2013b；吴林林等，2015；肖艳姣等，2012；赵川鸿等，2019；马建立等，2019)，为后续反射率质量控制做好铺垫。

目前已经业务化的双偏振雷达以S与C波段为主，很多质量控制工作也是围绕这两种波段的雷达数据展开，但是两者在局地强对流的精细观测方面有所不足。与S、C波段雷达相比，X波段雷达观测数据具有更高的时空分辨率，可以提供更细微水凝物信息，在灾害天气预警中发挥了更重要的补充作用。如北京自2015年搭建了5部X波段双偏振雷达并投入到科研与业务中，在随后的几年应用中发现X波段雷达容易受到环境噪声与信号衰减等因素影响，观测数据掺杂了晴空回波、电磁干扰杂波与生物回波等非降水回波。目前国内X波段双偏振雷达非降水回波识别的相关研究较少(郭春辉等，2019；王超等，2019)，研究人员基于CC去识别非降水回波，这可能将非均一填充与融化层区域中的降水回波错误识别为非降水回波。为了弥补以往研究的不足，本文对Tang et al(2014)设计的dpQC方案进行优化改进，设计了针对我国X波段双偏振雷达的质量控制方案来识别并滤除非降水回波，为业务化应用X波段双偏振雷达数据提供更好的质量控制方案。

1 资料与方法

1.1 数据来源

本文使用了北京昌平站、房山站、密云站、顺义站与通州站5部X波段双偏振天气雷达在2019年5—7月观测的基数据资料，雷达采用VCP21体扫模式，每个体扫有9个仰角，大约3 min完成一次体扫，距离库长为75 m。用于计算融化层高度的0℃高度来自北京探空数据。

根据探空站距离昌平、房山、密云、顺义与通州雷达站距离大致为42、33，70、39与31 km，因此，探空站0℃高度可以代表5部雷达。

1.2 质量控制方案原理

由于北京X波段雷达在生成基数据时已经做了地物回波的滤除处理，生物回波个例较难收集，本文设计的质量控制方案识别的非降水回波主要是电磁干扰回波与晴空回波。图1给出了方案的整体流程，包含了7个步骤，分别是：(1)CC阈值判别；(2)冰雹、非均一性波束填充(NBF)区域识别；(3)融化层识别；(4)ZDR阈值判别；(5)条状杂波滤除；(6)水平连续性检验；(7)斑块杂波滤除并在降水回波区域弥补小面积空洞。雷达基数据依次经过这7步自动化处理，最后保留的即为降水回波。

图1 北京X波段双偏振雷达质量控制方案流程图

与原始dpQC方案(Tang et al，2014)相比，本文的质量控制方案移除了CC纹理阈值判别，增加了ZDR阈值判别，这是由于通过大量质量控制试验发现：CC纹理阈值判别会滤除大量降水回波，ZDR阈值判别会去除面积较大的电磁干扰杂波；方案还基于统计结果对CC阈值判别与融化层区域识别的计算参数进行调整。本文将新的质量控制方案称为dpXQC(dual polarization X-band quality control)，下面详细介绍各步骤原理。

1.2.1 CC阈值判别

CC表示回波的水平偏振分量与垂直偏振分量的相关程度。研究表明，单一类型水凝物降水回波，其对应的CC较大，通常接近于1.0，而非降水回波的CC较小(Kumjian,2013)。基于CIMISS数据库中国地面逐小时资料的降水量与总云量、地面分钟降水资料的降水量数据选取了非降水回波与降水回波的个例(表1)，统计两者CC的概率分布(图2)。由图可见，降水回波CC主要分布在0.8～1.0区间内，非降水回波分布范围较广，在0.2～1.0。

表1 2019年5—7月降水过程时间及体扫数量

CC阈值的选取直接影响方案对非降水回波与降水回波的识别效果，由图2可知，阈值选取较低可减少对降水回波的错误识别，但会使少量非降水回波无法被滤除；反之则可以有效地识别非降水回波，但降水回波可能被错误识别。综合考虑对非降水回波滤除效果与降水回波错误识别因素，本文以0.05为步长，依次计算CC区间[0.5,1]对应的降水回波概率分布阈值右侧的累计概率值与非降水回波概率分布阈值左侧累计概率值之和，从选出最大值，最终将0.9定为阈值。如果雷达回波对应CC<0.9，将其判定为非降水回波。

1.2.2 冰雹或者非均一性波束填充区域判别

从图2可以看出，小部分降水回波对应的CC<0.9。Ryzhkov et al(2005)研究发现，冰雹区经常出现部分回波对应的CC较小，Ryzhkov(2007)还发现NBF区域也会出现较小的CC，这主要是由于随着雷达探测波束的延长,单位距离库内所探测的体积不断增大，强降水回波中心占据单位探测体积的一部分或者不同相态的水凝物会出现在同一单位探测体积中，使得ФDP梯度变大导致CC降低。由于冰雹或者NBF的低CC区域通常出现在强降水天气中，其对应的回波较强，回波顶较高，因此利用公式对其进行识别：

图2 北京X波段雷达降水回波与非降水回波对应的CC的概率分布

式中:ETOP18 dBz与ETOP0 dBz分别表示被检测的距离库之上18 dBz与0 dBz反射率的高度，rstormcore表示为强回波区到雷达站的径向距离，定义强回波区为某一径向大于45 dBz回波累计距离大于1 km的区域。被标识为冰雹与NBF区域中CC<0.9的回波不被滤除。

1.2.3 融化层亮带识别

融化层包含了融化的雪、霰与雹等不同相态的粒子，因此雷达探测到部分回波对应的CC较小。为了使融化层融区域中CC<0.9的降水回波不被滤除，借鉴了dpQC方案，即基于探空资料的0℃高度与CC来对其进行识别，步骤如下:

(1)首先将每个径向的雷达数据分为三组，融化层之下(0℃高度其下1 km至其下2 km)，融化层(0℃高度至其下1 km)与融化层之上(0℃高度至其上1 km)，计算每组CC平均数(分别为CCb、CCin、CCa)，公式如下：

(3)

式中:H0℃表示探空资料0℃高度;h(i)表示第i个距离库的高度;n表示每个分组中距离库的总数。

(2)当某个径向上三组数据CC的平均值满足下式时，这个径向距离库高度在[H0℃-1 km,H0℃)区间内被判定为融化层

CCin≥0.9∩{[CCin<(CCb-0.02)∩CCin<

(CCa-0.02)]∪CCin<(CCb-0.05)}

(4)

(3)在判定的融化层区域中，小于0.7的回波被标识为非降水回波并滤除，CC在[0.7,0.9)剩余的回波暂不处理，需要在接下来的步骤中继续判别回波属性。

本文统计了177个出现融化层的仰角，计算了融化层之下、融化层之中与融化层之上的CC的平均值(图3)，可以发现，融化层组的CC平均值在0.85～0.97波动，而且部分仰角融化层之上、之下CC平均值与融化层差较小，如果基于dpQC方案对融化层进行识别，会出现遗漏，基于统计结果，dpXQC修改了式(4)中的阈值，新的判别条件如下：

图3 177个仰角融化层之上、融化层与融化层之下CC的平均值

CCin≥0.85∩{[CCin<(CCb-0.01)∩CCin<

(CCa-0.01)]∪CCin<(CCb-0.03)}

(5)

1.2.4ZDR阈值判别

这一步是在dpQC方案基础上新增的，因为在实际观测中，会出现面积较大并且CC较大的电磁干扰杂波，这些回波无法利用CC阈值判别去滤除，也无法利用水平连续性检验去滤除。通过大量观测个例发现，ZDR可以较为有效识别这类非降水回波，ZDR反映了粒子偏离球形的情况，这些电磁干扰杂波的ZDR通常偏大或者偏小。本文统计了这类电磁干扰回波与降水回波二者的ZDR概率分布(图4)发现，这类回波的ZDR分布在-7～8 dB，大概率分布在绝对值较大的ZDR区间，降水回波则分布在-4～5 dB，为了尽可能保留降水回波，本文将︱ZDR|>5 dB做为区分此类电磁干扰回波与降水回波的判别阈值。

图4 北京X波段雷达降水回波与非降水回波对应ZDR的概率分布

1.2.5 条状杂波滤除

条状杂波通常由太阳射线或者电磁干扰造成，有的可以利用CC阈值判别进行滤除，但是也存在整条径向回波的CC>0.9的情况。条状回波的主要特征是径向上连续性较好，并且条状杂波通常出现在低层仰角，在垂直方向缺乏连续性。基于这些特征，本文质量控制算法中利用下式来识别条状杂波：

(6)

1.2.6 连续性检验

实际观测中常出现散点状回波，这些回波可能是晴空回波，也可能是电磁干扰回波，散点状回波的特点主要是具有孤立性，在方位与径向上都不连续，这些回波有时CC>0.9或者对应CC为缺测。算法利用连续检验可以滤除这类非降水回波，其检验原理是，以被检测的距离库为中心，设定一个窗口并统计范围内的ZH，如果多于一半是缺测值，或者有效回波ZH的平均值小于被检测回波ZH的25%，那么该距离库的回波被判定为非降水回波，本文将这个窗口大小设置为0.75 km×2.0°。

1.2.7 斑块杂波滤除

经过前面几步质量控制后，计算剩余的回波的面积，如果小于10 km2，那么这个区域被判别为斑块杂波并滤除。这一步可以滤除利用CC、ZDR阈值判别与连续性检验无法滤除的非降水回波。

2 质量控制效果分析

下面利用上节介绍dpXQC方案中各算法进行试验，为了验证效果，每一步算法单独进行试验，不与其他算法协同试验。

2.1 CC阈值判别

图5给出了2019年5月26日04：42密云雷达站回波分布，雷达站周围存在大量散点状电磁干扰回波与晴空回波(图5a)，对应的CC主要分布在0.1～0.5(图5b)，图5c显示了用CC阈值判别后的回波分布，可以发现雷达站周围大部分非降水回波被滤除。但是仍然存在很多散点状或小面积非降水回波保存下来，这些回波对应的CC为缺测或者大于0.9，后期可以利用连续性检验或者斑块回波判别进行滤除。

图5 2019年5月26日04:42密云雷达站0.5°仰角观测的(a)质量控制前ZH，(b)CC，(c)CC阈值判别后ZH

2.2 冰雹或者非均一性波束填充区域判别

图6a给出了2019年5月17日22:29北京顺义站雷达0.5°仰角质量控制前的ZH分布，在距离雷达100 km西南偏南方向有些回波的反射率达到45～60 dBz，实际观测中这一区域出现冰雹，这些区域有的回波对应的CC较低(图6b)，只利用CC阈值判别，可以滤除雷达站周围大部分非降水回波(图6c上部椭圆)，但是部分冰雹区与其附近区域(图6c下部椭圆)的降水回波也会被滤除；对冰雹区域进行识别并标识后(图6d)，使得区域中的降水回波(圆圈区域)不被滤除，但是冰雹区附近区域降水回波没有保留下来，这里主要由于不同相态的降水粒子出现在在同一距离库中形成NBF区域，加入NBF区域识别(图6e)，有效保护了这一区域的降水回波。

图6 2019年5月17日22:29顺义雷达站0.5°仰角观测的(a)质量控制前ZH，(b)CC，(c)CC阈值判别后ZH(上部椭圆内为经过CC阈值判别滤除的非降水回波，下部椭圆内为错误滤除的降水回波)，(d)冰雹区识别后保留的ZH(圆圈内为经过冰雹识别后保留的降水回波)，(e)冰雹与NBF区域识别保留的ZH

2.3 融化层亮带识别

图7a和7b分别给出了2019年7月5日20:03房山站雷达2.5°仰角质量控制前的ZH与CC分布，在距离雷达50～90 km、45°～180°方位角，大量回波的ZH分布在25～35 dBz，CC分布在0.6～0.9，结合当日的探空数据，这些区域出现在0℃层高度附近。单纯利用CC阈值判别会滤除大量降水回波(图7c椭圆区域内)，对融化层区域进行识别后(图7d)，保留了其中大部分CC<0.9降水回波。

图7 2019年7月5日20:03昌平雷达站2.5°仰角观测的(a)质量控制前ZH，(b)CC，(c)CC阈值判定后ZH(椭圆内白色区域表示被错误滤除的降水回波)，(d)融化层识别后保留的ZH

2.4 差分反射率判别

图8a给出了2019年7月13日12:15密云站雷达站0.5°仰角观测，在距离雷达站50～200 km、180°～240°方位角区域出现大量电磁干扰杂波，这些回波对应的CC大部分在0.9以上(图8b)，无法通过CC阈值判别来滤除，而其对应的ZDR>5 dB(图8c)，利用ZDR阈值判别后可以大部分电磁干扰回波被剔除(图8d)。

2.5 条状杂波滤除

图9a给出了2019年5月25日15:21房山站雷达站0.5°仰角回波分布，在雷达站南侧出现条状杂波，其对应的CC>0.9(图9b)，无法利用CC阈值判别进行滤除，利用了条状杂波识别算法，部分杂波可以被滤除(图9c)，剩余的回波主要由不连续的散点状或者小面积回波构成，可以利用连续性检验或者斑块回波识别进行滤除。

2.6 连续性检验

图10a和10b分别给出了图5c与图8c经过连续性检验后回波分布，通过对比可见，利用CC阈值判别或者条状杂波判别无法滤除的非降水回波，利用连续性检验可以将其滤除。

图8 2019年7月13日12:15密云站雷达0.5°仰角观测的(a)质量控制前ZH，(b)CC，(c)ZDR，(d)ZDR阈值判别后ZH

图9 2019年7月16日15:21房山雷达站0.5°仰角观测的(a)质量控制前ZH，(b)CC，(c)质量控制后ZH

图10 经过连续性检验后2019年(a)5月26日04:42北京密云雷达站0.5°仰角观测和(b)7月16日15:21房山雷达站0.5°仰角观测的ZH分布

2.7 斑块杂波滤除效果

图11给出了2019年5月26日11:54昌平站雷达0.5°仰角观测经过斑块状杂波滤除前后的ZH分布，对比图11a与11b可以发现，距离雷达站100～200 km的90°～120°方位角区域的斑块状杂波被滤除。

图11 2019年5月26日11:54昌平站雷达0.5°仰角ZH分布(a)质量控制前，(b)质量控制后

3 结果检验

为了评估质量控制方案对非降水回波的识别效果，借鉴了文浩等(2016)提出的评估方案，以一个仰角层的观测数据为单位，先人工观测进行判别，选出有、无非降水回波的仰角层数据，当质量控制方案识别非降水回波的面积超过人工判别非降水回波面积的90%，定为正确识别，当质量控制方案识别非降水回波区域在人工观测判别为降水回波，且其面积占降水回波的10%，定为错误识别，从北京5部X波段双偏振天气雷达站选取了320个体扫、共811个仰角层进行测试。

统计指标基于二分类事件评估方法，主要包括正确识别率(Hr)与错误识别率(Far)，公式如下:

(7)

(8)

式中：a为人工观测与质量控制算法均判别出现非降水回波的仰角数；b为人工观测判别没有出现非降水回波而质量控制算法判别错误的仰角数；c为人工观测判别出现非降水回波而质量控制方法判别遗漏的仰角数；d为人工观测与质量控制算法均没有判别出现非降水回波的仰角数。

在选取所有仰角层数据中，含有非降水回波的仰角层为123个，降水回波的仰角层为688个(有的仰角层也有非降水回波，但是比例非常小，忽略不计)。表2给出了dpQC方案与dpXQC方案质量控制效果统计，由表可见，非降水回波识别正确的仰角层数，dpXQC方案多于dpQC方案，dpXQC方案Hr为91.8%，dpQC方案的为85.4%，这主要是由于dpXQC方案增加了ZDR阈值判别，有效滤除了块状电磁干扰回波；降水回波被误识别为非降水回波的仰角层数，dpXQC方案少于dpQC方案，dpXQC方案的Far为20.6%，dpQC方案的为28.8%，主要是由于dpXQC方案调整了融化层识别参数，使区域中CC<0.9的降水回波更多地被正确识别，并且老方案中包含了CC纹理阈值判别，这一项有时会将降水回波判别为非降水回波滤除，dpXQC方案删除了这一判别项。

表2 dpQC方案与dpXQC方案的质量控制效果统计

然而，dpXQC方案对非降水回波错误识别率仍然很高，统计发现主要有两点原因：一是由于X波段雷达反射率随着距离的延长较S、C波段衰减更为严重，因此对于离雷达站较远的NBF区域，利用算法的原有参数去识别会出现错误，因而这些NBF区域CC<0.9的降水回波会被识别为非降水回波，目前由于缺乏个例累积，还无法对原始计算参数做出修订；二是融化层的识别算法是基于探空资料0℃层高度数据并统计其附近的CC特征来识别，探空数据一天只有三次，不能反映各时刻0℃层真实高度，同时算法中设定每个方位角的融化层位置相同且厚度固定为1 km，这与真实的融化层并不相符，因此识别区域外的真实融化层区域CC<0.9的降水回波被错误识别为非降水回波。

4 结论与讨论

本文借鉴了MRMS-QPE系统dpQC方案原理，并根据北京X波段多普勒双偏振业务雷达大量观测个例对其进行优化改进，设计了dpXQC方案，用其进行质量控制试验并且批量检验了dpXQC方案对非降水回波的识别效果，获得以下结果与结论：

(1)根据降水回波与非降水回波CC概率分布特征，dpXQC方案将区分二者的阈值调整为0.9，CC<0.9的回波判别为非降水回波，试验结果显示这可以有效去除雷达站周围的晴空回波与电磁干扰回波等非降水回波。

(2)冰雹区与NBF区域容易出现CC<0.9的降水回波，通过计算回波顶高与强回波区位置来对这些区域进行标识，标识区域内不进行CC阈值判别，可以有效保护这些区域CC<0.9的降水回波。

(3)融化层由于混合了多种相态水凝物粒子，部分降水回波的CC<0.9，dpXQC方案通过统计探空资料0℃层高度上下CC的特征来调整融化层的识别参数，以此保护融化层区域内CC<0.9的降水回波。

(4)dpXQC方案增加了ZDR阈值判别进行优化，通过统计降水回波与非降水回波ZDR概率分布特征，将|ZDR|>5 dB的回波判别为非降水回波，试验分析表明ZDR阈值判别可以去除CC>0.9的块状电磁干扰回波。

(5)对于无法利用CC与ZDR阈值判别滤除的非降水回波，如条状、散点与面积小于10 km2的非降水回波，算法通过条状杂波识别、连续性检验与斑块杂波识别分别来滤除。

(6)统计检验结果表明，dpXQC方案对非降水回波正确识别率为91.8%，错误识别率为20.6%，优于dpQC方案；其错误识别率依然较高，主要有两方面原因，一是NBF区域的识别算法的参数没有进行本地化调整，导致部分NBF区域没有被标识，这些区域内CC<0.9降水回波被误判为非降水回波；二是融化层识别算法在所有方位角的位置与厚度是固定的，与实际观测不相符，这会造成被遗漏识别的融化层区域中CC<0.9的降水回波被错误识别为非降水回波。

总体来看，本文设计的dpXQC方案，可以滤除条状、散点状与块状的电磁干扰回波和晴空回波等非降水回波，并最大程度地保护了冰雹、NBF与融化层区域中的降水回波不被滤除。同时也要看到，在部分NBF区域与融化层区域，降水回波会被错误判别为非降水回波。对于离雷达站较远的NBF区域识别，需要通过积累大量个例去修订原有的计算参数。对于融化层的识别，还要研究或借鉴更新的识别算法，这是后续改进的方向。