徐舒扬 吴 翀 刘黎平
1)(成都信息工程大学电子工程学院, 成都 610225)2)(中国气象科学研究院灾害天气国家重点实验室, 北京 100081)
降水粒子群通常呈非球型且空间取向大致相似,对两个相互正交的线性极化偏振波的散射存在差异,双偏振雷达利用以上特性提取出反映两个通道信号差异的偏振参量,这些偏振参量包含了更丰富的云微物理信息[1-3],双偏振雷达因此成为水凝物相态识别和定量估测降水的重要工具。
20世纪90年代初,Straka等[4-5]将模糊逻辑法引入基于双偏振雷达数据的水凝物相态识别研究中,通过将偏振参量转化成0~1的数值再依据权重累加。其进步在于以过渡阈值取代经典布尔逻辑的固定阈值,实现信息综合分析同时简化了布尔逻辑繁琐的判断步骤。隶属函数即过渡阈值是模糊逻辑的核心。1999年Zrnic等[6]将隶属函数由二维矩阵扩充成5个参量的矩阵,将方法的分析能力从仅能区分气象和非气象目标物提高到能识别出9种水凝物的回波。过去隶属函数的确定包含很多的主观因素,1999年Vivekanandan等[7]利用绘制偏振参量的二维频率分布图确定出更客观的隶属函数。2009年Park等[8]改进了模糊逻辑相态识别算法(Hydrometeor Classification Algorithm, HCA),将隶属函数修改成非对称的梯形,为平衡数据误差引入置信矢量。2000年Liu等[9]用神经网络生成隶属函数,由于神经网络的训练需要大量真实值,而空中降水粒子相态的实测数据有限,导致方法难以广泛应用。2013年Al-Sakka等[10]发现用模式提供三维温度场辅助识别能大幅提高识别的准确度。2015年Bechini等[11]提出采用聚类方法进行预处理降低噪声对识别效果的影响。2007年曹俊武等[12]最早探讨了将方法应用到国内双偏振雷达的可行性。2017年Wu等[13]参考美国HCA方法,结合我国广州S波段、X波段双线偏振雷达数据的统计特征,初步建立了适用于广东S波段双线偏振雷达的相态识别方法。利用Zh,ZDR,ρhv参数和前一时次零度层位置确定当前零度层高度,通过零度层位置订正湿雪的识别结果,并根据湿雪和零度层的空间分布对干雪和小雨的识别结果进行订正。
HCA方法仍存在局限性,一是雷达观测量带有系统偏差和随机误差,误差如何影响识别结果尚未有系统讨论,这是方法的稳定性问题;二是如何确定识别结果是否正确,隶属函数是否为最优解,以及不同雷达适用的隶属函数是否也存在差异,这是方法的准确性问题;三是不同相态识别结果的时空分布是否符合云物理规律以及附加订正效果如何,这是方法的可靠性问题。因此,本文提出误差敏感性检验、识别效果检验和相态时空分布合理性检验方法,探讨相态识别方法的准确性、稳定性和可靠性,并提出改进方案。
目前模糊逻辑相态识别方法能识别的10种相态为地物回波、晴空回波、干雪、湿雪、冰晶、霰、大雨滴、小到中雨、大雨、冰雹(包括干雹和湿雹)。使用的6个偏振参量为雷达反射率因子Zh、差分反射率因子ZDR、协相关系数ρhv、差分相移率KDP、反射率因子标准差σZh、差分相移标准差σФDP。为使用方便,下文KDP为对数形式。HCA方法首先将经过质量控制的观测值用模糊规则计算判据,将判据和其他变量累加得到累加值A(范围为0~1),通常最大A对应相态即为识别结果。A的计算[13]见式(1):
(1)
式(1)中,P(i)(Vj)为通过隶属函数得到的判据(j指第j个偏振参量,i指第i类相态),Wij为偏振参量的权重,Qj为有关数据质量的置信矢量。
针对相态识别的局限性,本文将从以下4个方面系统地检验相态识别效果并提出改进方案:①识别效果检验,通过统计分析A的分布找到影响识别的关键因素;②方法稳定程度检验,通过引入误差检验数据质量对识别的影响;③模糊规则的改进,通过统计相态特征改进不合理的模糊规则,通过分析参量贡献改进权重矩阵;④附加步骤改进,通过空间一致性订正减少识别结果的误判。采用广州S波段双偏振雷达(简称广州雷达)数据,包括2016—2017年暖季13个不同类型降水过程,样本时段和降水云类型如表1所示。
表1 样本时段和降水云特征Table 1 Periods and characteristics of precipitation samples
Straka等[4-5]建立模糊规则的基础是用梯形隶属函数描述水凝物的特征,通过观测量与隶属函数比较获得分类的判据。隶属函数与相态特征越匹配,识别效果就越好。单一A的数值仅是识别相态的依据,而相态的A则可作为检验隶属函数的依据。统计高、中、低分段A的比例,中、低分段A比例越高说明隶属函数越不匹配。采用统计相态偏振参量特征的方法改进隶属函数,根据水凝物形成的物理规律,选出相态所在区域统计偏振参量的频率分布。同类相态偏振参量通常集中分布在一定范围内,取0.1和0.9频率对应的数值作为外边界,取0.2和0.8对作为内边界修改隶属函数。根据顾震潮[14]提出的层状云粒子分布,冰晶和干雪特征取7~11 km高度上统计冰雪混合区的偏振参量分布。根据冰雹的雷达回波特征[15-18],选取地面有降雹记录的对流过程,统计Zh>50 dBZ且出现三体散射回波的偏振参量分布。根据地物回波和晴空回波的特征[19-20],选取在距雷达不超过50 km、高度不超过5 km范围内无气象回波的数据,为晴空回波统计滤波前后Zh之差低于1 dBZ的参量分布,为地物回波统计滤波前后Zh之差大于20 dBZ的参量分布。
分析A可以检验方法准确性,也可以检验方法区分相似相态的能力。HCA方法的隶属函数之间重叠越多得到的A越接近,方法的区分能力越弱,轻微的扰动就可能改变识别的结果。为此检验方法的区分能力,统计相态A与次高A(A2)之差ΔA≤0.1的比例以及A2所属相态比例,差值偏小的识别结果比例越高说明方法的区分能力越弱。改进方法是分析相似相态区分能力弱的原因,提高区分能力较强的偏振参量的权重,验证修改权重后的识别效果并修改权重矩阵。
误差对相态识别的影响不可避免,而以往研究报道未见系统讨论误差的影响程度,因此,提出引入适量误差对比识别结果变化的误差敏感性检验法。雷达误差有系统偏差和随机噪声[21-22],向数据中加上或减去适当固定数值模拟系统偏差,用生成高斯白噪声的方式模拟随机噪声。高斯白噪声是按预设标准差生成满足高斯分布的随机值,白噪声的均值为0且68.27%的值在标准差范围内。因为不稳定识别结果通常在梯形隶属函数的斜边附近,所以用小于或等于斜边长度的数值即可检验方法敏感性。
相关研究引入一些附加步骤优化识别结果的空间分布,但未用统计方式检验订正效果,也未系统分析识别出的相态时空分布是否合理。为此,本文提出通过检验相态时空分布的连续性检验识别效果。统计各相态在不同高度的出现频率,讨论相态出现高度随时间变化,分析相态时空分布的合理性,针对不合理部分分析原因并找出修正方法。
对引入适当误差后的数据进行相态识别,统计各相态数量的变化比例,得到可保证识别结果稳定的误差范围:Zh为-0.5~+0.5 dBZ,ZDR为-0.1~+0.1 dB,ρhv为0~0.02,KDP为-0.3 ~+0.9 dB。各相态对Zh和ZDR误差敏感性的检验结果如图1所示,每种相态中减少的识别结果所占比例为负值,增加的识别结果所在比例为正值,ρhv和KDP误差的影响较小可忽略。
图1 相态识别对误差敏感性的检验(a)对Zh系统偏差的敏感性,(b)对ZDR系统偏差的敏感性,(c)对Zh噪声的敏感性,(d)对ZDR噪声的敏感性Fig.1 Algorithm sensitivity to data error(a)sensitivity to Zh system deviation,(b)sensitivity to ZDR system deviation,(c)sensitivity to Zh artificial noise,(d)sensitivity to ZDR artificial noise
对误差敏感的结果通常位于隶属函数斜边附近,占比不宜超过20%。对误差敏感的相态主要有地物、冰晶、干雪、大雨滴、湿雪、霰、大雨和冰雹。对比权重相近的Zh和ZDR误差的影响,发现同等程度下Zh误差影响比例远高于ZDR。Zh决定霰、大雨滴、大雨和冰雹相态ZDR和KDP的隶属函数,Zh负偏差也会使ZDR隶属函数范围缩小,部分识别结果因不满足更小范围的隶属函数而改变。这些特殊的ZDR隶属函数使Zh误差产生的影响远大于ZDR,且Zh和ZDR负偏差的影响程度远大于正偏差。因此,方法应用中可尽量提高Zh和ZDR的数据质量并注意负的系统误差对相态识别影响。另外,同等程度系统偏差的影响比噪声大,在实际中需注意雷达系统偏差对识别结果的影响。
用A检验隶属函数的合理性,不同分段A占比如表2所示。不合理的隶属函数使中低分段的比例偏高,如地物、晴空回波、干雪、冰晶、小到中雨和冰雹,这些相态中多数高分段比例偏低,说明隶属函数不适用,也有如干雪和小到中雨高分段A的比例很高,中间分段的比例较低,说明隶属函数适用但对相似相态的区分能力不足。
表2 A的分布和ΔA≤0.1的比例Table 2 Distribution of A and proportion of ΔA ≤ 0.1
用ΔA检验隶属函数的区分能力,差值越小说明方法区分能力越弱,识别结果越不理想。差值小于0.1的识别结果在各相态中的比例,以及主要易混淆的相态,如表2所示。可以看到,地物和晴空回波的区分效果不佳,这两种相态A均偏低,隶属函数不合理导致了方法对地物和晴空回波既不能有效识别也不能合理区分。区分能力不足的还有冰晶、干雪、大雨和冰雹,但这4种相态A中、高分段比例较高,可能是特征相似导致方法难以区分。
通常冰晶形成于云温低于0℃并且凝结核充足的区域,而后在一定条件下生长成雪或霰,霰和冻滴密集区有利于冰雹生成。因此,冰雹一般出现在零度层附近,周围有密集的霰和湿雪。而在识别结果中,低层常有异常出现冰雹相态的情况。以2017年4月21日08:00—17:00(北京时,下同)的强对流过程为例,统计冰晶、霰、干雪、湿雪、冰雹和大雨的出现频次,分别绘制6种相态时空分布的频率图(图2)。零度层附近冰雹与霰出现频率的变化基本一致,而低层识别出的冰雹不符合其生成规律。取13:54时次绘制水平和垂直结构剖面图(图3)。可见低空冰雹的环境温度远高于0℃且上方均为液相粒子,不符合生成冰雹的条件。大雨与冰雹的参量特征相似,方法对这两种相态的区分能力较弱,低层冰雹出现频率与大雨一致,可见方法存在将低层大雨误判为冰雹的问题。
图2 2017年4月21日08:00—17:00部分相态高度的时空分布(填色表示频率)Fig.2 Distribution of hydrometeor heights from 0800 BT to 1700 BT on 21 Apr 2017(the shaded denotes the frequency)
图3 2017年4月21日13:54水平和垂直方向的Zh和相态结构(a)Zh水平分布,(b)相态水平分布,(c)Zh垂直结构,(d)相态垂直结构Fig.3 Zh and hdyrometeors’ horizontal and vertical structure at 1354 BT 21 Apr 2017(a)Zh horizontal structure,(b)hdyrometeor horizontal structure,(c)Zh vertical structure,(d)hdyrometeor vertical structure
经过识别效果的检验分析,归纳识别效果不理想的相态得到表3。针对冰雹的A偏低改进隶属函数,为提高方法对冰雹和大雨的识别和区分能力改进权重,并利用冰雹的空间分布特征订正误差识别结果。针对冰晶和干雪的A偏低和不易区分改进高权重偏振参量的隶属函数和权重,针对地物和晴空回波的A偏低改进隶属函数。
表3 识别效果分析存在不足的相态Table 3 Insufficiency of hydrometeor classification results
修改非气象回波类相态的隶属函数,地物指雷达发射的电磁波遇到地物反射或散射回来的信号,晴空回波包括由大气湍流和空气折射率不均匀产生的非降水回波,以及由飞鸟、昆虫等生物产生的回波[19-20]。采用2017年3月18—20日和3月31日无气象回波的数据,统计地物和晴空回波的相态特征形成图4a~图4d,Zh-ρhv和σZh-σФDP二维频率分布图用不同颜色表示频率高低,边框的一维分布图由柱状体高低表示频率,实线和梯形上底是新隶属函数的内边界,虚线和梯形下底是新隶属函数的外边界。原隶属函数基于无地物抑制的样本确定,经过滤波地物的Zh,σZh和晴空回波的Zh小于原来的隶属函数值,导致地物和晴空回波的A过低。
图4 相态偏振参量特征和隶属函数改进(填色表示频率)(a)晴空回波的Zh-ρhv,(b)晴空回波的σZh-σФDP,(c)地物回波的Zh-ρhv, (d)地物回波的σZh-σФDP,(e)冰晶和干雪的Zh-ZDR(新隶属函数干雪为橙线,冰晶为黑线,下同), (f)冰晶和干雪的Zh-ρhv,(g)冰雹的Zh-ZDR,(h)冰雹的Zh-KDRFig.4 Hydrometeor variable feature and improvement of membership function(the shaded denotes the frequency) (a)Zh-ρhv feature for biological scatterers,(b)σZh-σФDP feature for biological scatterers,(c)Zh-ρhv feature for ground clutter,including that due to anomalous propagation,(d)σZh-σФDP feature for ground clutter,including that due to anomalous propagation,(e)Zh-ZDR feature for crystals of various orientations(the orange line) and dry aggregated snow(the black line),(f)Zh-ρhv feature for crystals of various orientations and dry aggregated snow, (g)Zh-ZDR feature for a mixture of rain and hail,(h)Zh-KDR feature for a mixture of rain and hail
续图4
干雪和冰晶的相态特征,采用两次层状云降水数据样本(2017年3月18日17:00—20:00和2016年6月4日13:00—14:00),统计相态特征绘制Zh-ZDR和Zh-ρhv二维频率分布图和Zh,ZDR和ρhv的频率柱状图,如图4e和图4f所示,黑线表示改进后的隶属函数冰晶,橙线表示干雪。利用有降雹记录的强对流过程(2016年5月9日16:00—18:00),统计冰雹回波的偏振参量Zh,ZDR和KDP的频率分布特征,绘制图4g和图4h。冰雹的ZDR原隶属函数(-0.3,0.0,f1,f1+0.5)的下边界过高,KDP原隶属函数(-10,-4,g1,g1+1)的上边界数过低,对冰雹特征描述不足使冰雹的A偏低(f1,g1[8]如式(2)、式(3)所示)。修改后的隶属函数曲线如图4所示,同时为保证冰雹和大雨隶属函数过渡不重叠,修改大雨隶属函数的下边界。以上改进的隶属函数如表4所示。
f1(Z)=-0.5+2.50×10-3Z+
7.50×10-4Z2,
(2)
g1(Z)=-44.0+0.8Z。
(3)
表4 改进前后隶属函数Table 4 Membership functions before and after modification
模糊规则的权重矩阵通常由综合理论值与研究者的经验确定[23],本文所用Park等[8]和Wu等[13]给出的权重矩阵是将权重划分成6个等级(1.0,0.8,0.6,0.4,0.2,0.0),依据偏振参量对相态识别的贡献程高低赋予权重。提高区分能力更强的偏振参量的权重能有效提高方法的区分能力,区分效果不好的冰晶与干雪主要在Zh隶属函数重叠的区域,该区域内ZDR更容易区分冰晶与干雪。因此,试验将冰晶的ZDR权重由0.6分别提高到0.65,0.70,0.75 和0.80,统计提高权重让冰晶各分段A和ΔA的变化(如图5所示)。在冰雹和大雨的Zh重叠区域主要依靠ZDR和KDP区分,考虑到KDP的权重已最高,仅修改ZDR权重,由0.8分别改为0.85,0.90,0.95 和1.0,统计大雨和冰雹A和ΔA的变化(图5)。提高ZDR权重后冰晶和大雨的识别效果有一定提高,对冰雹的识别效果没有提高,仅将冰晶和大雨的ZDR权重修改为0.8和1.0。
图5 冰晶、大雨和冰雹改变ZDR权重后A和ΔA的变化趋势Fig.5 Variability of A and ΔA for crystals,heavy rain and a mixture of rain and hail after changing ZDR weight
续图5
Zh为45 ~55 dBZ时,虽然很难根据偏振参量区分ZDR和KDP相近的大雨和冰雹,但两种相态的空间分布却有很大差异。因此,可利用冰雹分布的连续性检验并订正识别结果,附加步骤订正需分别检验冰雹与上方和下方粒子分布的连续性。首先检查识别出的冰雹上方是否为一致的冰相粒子,若是,则检查其下方冰雹的分布是否连续,存在一致的液相粒子间隔后再出现的冰雹无法通过检验。未通过检验的冰雹均订正为大雨相态,统计被订正冰雹的比例及平均出现高度,结果如表5所示。由表5可知,不满足冰雹连续性的识别结果出现高度均远低于零度层的高度,且远低于满足连续分布的冰雹高度。
表5 误判冰雹的比例和高度Table 5 Proportion and height of misclassified a mixture of rain and hail
本文利用2016—2017年暖季广州雷达降水数据进行了方法识别效果分析,提出了误差敏感性检验、识别效果检验和相态时空分布检验方法,并基于检验结果分析找出影响识别效果的关键因素加以改进。得到如下结论:
1) 通过A检验隶属函数的合理性以及方法对相似相态区分能力分析可知,之前的方法由于隶属函数不合理导致地物、晴空回波、干雪和冰晶的准确度较低,对3组相态地物和晴空回波、冰晶和干雪以及霰、大雨和冰雹区分能力弱。经过分析影响识别的关键因素,使用统计相态偏振参量特征的方法改进隶属函数,使用A统计检验法改进权重矩阵。
2) 通过引入的系统偏差和模拟噪声检验,找出保证全部相态识别稳定性数据质量要求:发现Zh误差范围在-0.5~+0.5 dBZ、ZDR误差在-0.1~+0.1 dB、ρhv误差在0~0.02、KDP误差在-0.3~+0.9 dB 范围内,识别结果稳定性较好。此外,检验发现部分相态(冰雹、大雨、大雨滴和湿雪)对误差影响的敏感性很高,主要是对偏振参量Zh和ZDR的误差表现敏感。多数相态对ρhv和KDP误差表现不敏感,因此,应用时应更注意Zh和ZDR参量的数据质量。
3) 提出相态时空分布合理性检验法,检验相态识别结果的空间分布是否符合物理规律。由低层冰雹面积异常增加的结果发现异常识别为冰雹的大雨相态,提出冰雹一致性订正方法对识别结果进行订正,经检验,该订正方法能有效筛除被异常识别的冰雹相态。