福建省氣象信息中心 袁 偉 高 攀 鄭玉蘭

雷達拼圖射線狀雜波抑制算法及其應(yīng)用

福建省氣象信息中心 袁 偉 高 攀 鄭玉蘭

從雷達產(chǎn)品拼圖中出現(xiàn)的射線狀雜波出發(fā)，研究其產(chǎn)生的機理與特征，并在大量統(tǒng)計分析的基礎(chǔ)上，建立相應(yīng)的射線狀雜波篩選、去偽存真和逐點抑制的檢驗?zāi)Ｐ?，同時對模型中用到的檢驗參數(shù)進行了分析與論證。該文設(shè)計的射線狀雜波篩選游程算法，可以有效定位雷達產(chǎn)品數(shù)據(jù)中的射線狀雜波，應(yīng)用強相關(guān)性檢驗算法和T-值檢驗算法進一步得到更加可靠的射線狀雜波數(shù)據(jù)，最后設(shè)計格點平移窗口對射線狀雜波數(shù)據(jù)進行逐點抑制。三個步驟的抑制算法過程實施后，能夠改善雷達拼圖效果，在射線狀雜波對天氣過程判斷造成的干擾有一定的抑制作用。

雷達拼圖 射線狀雜波 檢驗?zāi)Ｐ?游程算法 強相關(guān)性檢驗

多普勒天氣雷達在短臨、臺風(fēng)、暴雨等天氣過程預(yù)報中起的作用已越來越重要。其觀測產(chǎn)品數(shù)據(jù)、小時累計降水產(chǎn)品、基本反射率、組合反射率、回波頂高、基本速度等十幾個直接產(chǎn)品和幾十種衍生產(chǎn)品已應(yīng)用到不同場合下天氣預(yù)報預(yù)測中[1]。與此同時，在雷達觀測數(shù)據(jù)的質(zhì)量控制研究方面，國內(nèi)外也做了大量的技術(shù)研究與應(yīng)用[2-5]，并取得了非常有價值的成效，由于雷達產(chǎn)品因種類繁多、應(yīng)用廣泛，仍存在部分沒有質(zhì)控或質(zhì)控算法未普及到的地方，給實際預(yù)報預(yù)警判斷上會帶來一定的困擾。

雷達基本反射率、組合反射率產(chǎn)品是雷達最基本的產(chǎn)品，常以一定區(qū)域范圍內(nèi)多部雷達進行拼圖的形式展示出來，在降水動態(tài)監(jiān)測預(yù)警中直觀、形象，對預(yù)報員而言，有很強的既視感，如圖1所示。

圖1 基本反射率多部雷達拼圖實況

雷達拼圖產(chǎn)品在實際應(yīng)用中，時有一些射線狀的雜波出現(xiàn)，對預(yù)報員準確判讀產(chǎn)生了一定干擾影響，同時也存在回波圖像遮擋等問題。本文將對雷達拼圖中出現(xiàn)的射線狀雜波進行研究，并提出相應(yīng)的抑制方案和實現(xiàn)過程，試圖最大可能的還原真實回波圖像。

1 雷達拼圖中射線狀雜波抑制概述

分析雷達拼圖產(chǎn)品中出現(xiàn)射線狀雷達回波的數(shù)據(jù)，可以發(fā)現(xiàn)它有幾個特點：（1） 常常以一根或多根孤立的線段存在，而周邊沒有回波。（2）線段出現(xiàn)在掃描線的末端居多。（3）線段上的回波強度一般不強，且大部分在中間值左右。（4）同一掃描線上，除該線段外，沒有更長的其他線段。

當多部雷達進行產(chǎn)品拼圖時，單站的射線狀雜波可能會遮蓋周邊雷達探測的有效數(shù)據(jù)，造成不必要的誤判，而且這類射線本身沒有實際上的天氣現(xiàn)象含義，但其往往與真實的降水回波疊加在一起，容易造成人工識別上的困擾。因此抑制這類非天氣過程引起的射線狀雜波，有一定的研究價值和意義。

根據(jù)上述特點，本文將追溯到拼圖產(chǎn)品來源的單站雷達產(chǎn)品，展開相應(yīng)的技術(shù)分析，試圖建立一個較為客觀的射線刻畫模型，遵循用模型算法定位射線，然后加以抑制，從而達成拼圖產(chǎn)品中剔除射線的研究目標。

2 射線狀雜波抑制算法

由于雷達拼圖的數(shù)據(jù)來自雷達單站產(chǎn)品，通過研究單站雷達產(chǎn)品的特征，采用射線狀雜波篩選、偽射線狀雜波剔除、相鄰射線狀雜波抑制、基準射線狀雜波抑制等措施，最終實現(xiàn)去偽存真、有效抑制射線狀雜波的算法。

下面的研究將以單站雷達基本反射率19號產(chǎn)品為例進行分析論證,其回波強度量化后從0～15分為16個等級，值越大表示回波越強，這里我們稱之為能量值。

2.1 射線狀雜波篩選

單站雷達產(chǎn)品數(shù)據(jù)上射線狀雜波的特征，主要體現(xiàn)在：孤立的一條射線，射線上回波較強，周邊回波非常弱（量化后一般為零），射線長短不一，粗細不均等。根據(jù)這些特點，設(shè)計回波游程統(tǒng)計，相對集中的能量判決算法，初步過濾出射線狀雜波。

針對單站雷達產(chǎn)品，每條掃描線上，統(tǒng)計游程：

其中N為掃描距離數(shù)；ri為掃描線上第i個掃描距離的回波值；n0、n1分別為掃描線上回波值為0和非0時的連續(xù)統(tǒng)計計數(shù)；Y[n0]、Y[n1]分別為統(tǒng)計連續(xù)回波為0和非0出現(xiàn)的距離數(shù)長度計數(shù)。p為連續(xù)能量計數(shù)器，MP =8為能量中值，P為連續(xù)非0距離數(shù)的回波能量累加和。

記游程Y[n0]最大值所對應(yīng)的距離長度為Md0。記游程Y[n1]最大值所對應(yīng)的距離長度為Md1，此時相應(yīng)的能量為P[Md1]。記游程Y[n1]第二大值所對應(yīng)的距離長度為d1。

則滿足如下條件的掃描線可定為疑似射線狀雜波：

(1)Md1/N >0.5，有回波的最大距離長度超過總長度的一半以上。

(2)P[Md1] < 1，能量值控制在一般回波強度范圍內(nèi)，即要求不出現(xiàn)特別強的回波。

或

(1)（Md0+ Md1）/N > 0.8，有回波值和無回波值的距離長度之和占據(jù)大部分掃描線。

(2)P[Md1] < 1能量值控制在一般回波強度范圍內(nèi)，即要求不出現(xiàn)特別強的回波。

(3)d1/ Md1 < 0.5，有回波的第二長度的距離數(shù)不超過有回波的最大距離長度一半。

上述判決條件是基于2017年2月份中15000個文件（總計5788562條掃描線）的數(shù)據(jù)統(tǒng)計分析情況而定。設(shè)定幾個中值進行比較，統(tǒng)計結(jié)果見表1～表3。

表1 有回波的最大距離長度與總長度占比分析

表2 有回波值和無回波值的距離長度之和與總長度占比分析

表3 有回波值的第二距離數(shù)與第一距離數(shù)比率分析

統(tǒng)計表中的偏好系數(shù)為：估算射線數(shù)（樣本比率乘以疑似文件數(shù)）與疑似射線的比值，其值越大表明篩選效果越好。綜合表1～3結(jié)果分析，在篩選射線狀雜波時，為盡可能的選入存在射線狀雜波的數(shù)據(jù)（偏好系數(shù)相差不大的情況下，選稍小比率或占比），可選取如下門限參數(shù)：①有回波的最大距離長度與總長度占比=0.5；②有回波值和無回波值的距離長度之和與總長度占比=0.8；③有回波值的第二距離數(shù)與第一距離數(shù)比率=0.5。

2.2 偽射線狀雜波剔除

經(jīng)過射線狀雜波篩選預(yù)處理后，大部分射線狀雜波文件能夠選入，但仍然存在部分偽射像雜波，需要進一步進行剔除過濾，盡管會錯過不少真實射線狀雜波的選入，但可以確保極大概率剩下真實的射線狀雜波。

根據(jù)單站產(chǎn)品數(shù)據(jù)文件分析，射線狀雜波的周邊能量值一般不超過1。設(shè)計基于當前射線狀雜波掃描線為基準，統(tǒng)計與周邊掃描線的相關(guān)性,如強相關(guān)性、T-值分布等[6]，進行進一步濾除。

以當前射線狀雜波所在掃描線為基準(稱為：基準線)，按掃描方位角前后設(shè)定相關(guān)統(tǒng)計窗口(窗口內(nèi)的掃描線稱為：窗口線)，窗口數(shù)記為W。

（1）強相關(guān)性統(tǒng)計。

用來做兩條掃描線的相關(guān)性檢測。強相關(guān)系數(shù)越大，表明兩者之間的相似性越強。

計算基準線與窗口線的強相關(guān)系數(shù)：

其中xi為基準線上的射線狀雜波段的數(shù)據(jù)，m為該雜波段數(shù)據(jù)的長度（即點數(shù)），為該段數(shù)據(jù)的均值。yki為第k個窗口線上的對應(yīng)基準線掃描距的數(shù)據(jù)，為相應(yīng)的均值。

這樣得到W個強相關(guān)系數(shù)Sk。

（2）T-值檢測。

用來做兩條掃描線的顯著性差異檢測，T-值越大，表明兩者之間的關(guān)系越小。

這樣得到W個T-值系數(shù)kT。

根據(jù)上述統(tǒng)計分析得到的兩種檢測標準，我們用于檢測基準線是否為符合一定檢測水平的射線狀雜波，另一方面可以識別窗口線是否可以歸類到射線狀雜波。

檢驗水平門限值的選定統(tǒng)計如圖2、圖3所示。圖2中，橫軸為強相關(guān)系數(shù)門限范圍[0-1]，按0.1刻度放大10倍，縱軸為符合檢測標準的射線數(shù)，總共統(tǒng)計6757條疑似射線數(shù)，實際射線數(shù)7條?？梢妔0取值越大，檢測門檻就越低，當s0=1時，有970條通過檢測。s0取值越小，檢測門檻越高，當s0=0時，疑似射線全部排除。

圖2 S-檢驗水平門限值統(tǒng)計圖

圖3中，橫軸為T-值門限范圍[0-10]，按0.5刻度放大2倍，縱軸為符合檢測標準的射線數(shù)，總共統(tǒng)計6757條疑似射線數(shù)，實際射線數(shù)7條。可見t0取值越大，檢測門檻就越高，t0≥4，疑似射線全部排除；t0取值越小，檢測門檻越低，當t0=0時，有2326條通過檢測。

圖3 T-檢驗水平門限值統(tǒng)計圖

綜合實際射線數(shù)，以及上述兩種檢測門限的趨勢，取s0=0.2，t0=3.0，可以大概率剔除偽射線狀雜波。

同理統(tǒng)計分析后，我們?nèi)1=0.6，t0=7.0，可以顯著判定相鄰的窗口線是否屬于射線狀雜波。

2.3 雜波抑制的逐點判決方法

對于一個單站雷達產(chǎn)品文件，通過上文所提供的算法計算后，可以大概率找到射線狀雜波回波掃描線。如前所述，我們可以定位出射線狀雜波的位置和長度，為盡量避免可能有用的回波數(shù)據(jù)，采用逐點判決。

假設(shè)xi（i=1,m）為基準線上的雜波點能量值，以xi為中心點，設(shè)定一個7×7的網(wǎng)格窗口，如果窗口內(nèi)的（除中心點）點yijk的平均值低于中心點到一定閾值時，則將該中心點能量值置為0。移動窗口遍歷基準線上的m個雜波點，從而實現(xiàn)雜波的抑制：

Pi越小，表明剔除中心點的可能性越大。

實驗表明，一般取0.5作為Pi的判別閾值。當Pi<0.5時，中心點可剔除（能量值置為0），否則我們保留該點數(shù)據(jù)。

3 實際應(yīng)用效果

上述算法應(yīng)用到實際雷達拼圖處理中，有效地對射線狀雜波進行了抑制，在減少了一些射線狀回波對短臨預(yù)報的干擾方面有一定作用（圖4，圖5）。此外，該方法還能抑制2016年9月出現(xiàn)的一次“餅圖”現(xiàn)象（圖6）。

圖4 射線狀回波抑制前后對比實況1

圖5 射線狀回波抑制前后對比實況2

圖6 一次抑制“餅圖”現(xiàn)象的前后對比

4 結(jié)論

射線狀雜波的抑制想法主要來源于工作中發(fā)現(xiàn)諸類現(xiàn)象不少，且或多或少在短臨預(yù)報、氣象服務(wù)等方面存在一定的影響，目前也確實沒有完全能夠抑制的方法，基于盡量發(fā)現(xiàn)，準確抑制的想法，經(jīng)過大量數(shù)據(jù)分析及實驗驗證，試圖尋找一些技術(shù)方法來達成該目標。

目前設(shè)計的抑制算法能有效識別并剔除射線狀雜波。但算法在實際應(yīng)用過程中，還存在一些難以定位的射線狀雜波，現(xiàn)有的做法是不對其進行抑制處理。此外對這些射線雜波，經(jīng)過分析，如果一些檢驗標準放寬些，可以很容易尋找到，但同時也會帶來負面影響，導(dǎo)致一些非射線狀的雜波被誤判，因此算法還有很大的改進空間，如設(shè)計更多的檢驗算法，設(shè)計更靈活、能自適應(yīng)的檢驗標準，引入其他氣象觀測數(shù)據(jù)作為參照等等。

此外，本算法適用于射線狀雜波的抑制，也為對于地物回波、鳥群干擾等非天氣現(xiàn)象回波方面的研究提供了一種參考思路，即從數(shù)據(jù)自身出發(fā)，研究其特征與非天氣現(xiàn)象回波的關(guān)聯(lián)性，從而提出一些假設(shè)檢驗標準。

[1] 俞小鼎,姚秀萍, 熊廷南, 等. 多普勒天氣雷達原理與業(yè)務(wù)應(yīng)用[M]. 北京：氣象出版社, 2006.

[2] 陳媛, 陳江民, 王紫陽, 等. 天氣雷達反射率基數(shù)據(jù)質(zhì)量控制的幾種算法[J].天氣與減災(zāi)研究, 2007,30(3): 48-51.

[3] 馬中元.CINRAD雷達數(shù)據(jù)質(zhì)量控制方法初探[J]. 氣象. 2010,36(8): 134-141.

[4] Waldteufel P，Corbin H. On the analysis of single-doppler radar data[J]. Journal of Applied Meteorology, 1979, 18(4):532-542.

[5] Meng Z. Methods for improving data quality for WSR-98D weather radar[J]. Meteorological Science & Technology, 2006, 34(21): 85-89.

[6] 梁之舜,鄧集賢, 楊維權(quán). 概率論與數(shù)理統(tǒng)計[M].北京: 高等教育出版社,1992.