蔣立輝，高 浩，莊子波，熊興隆

(1.民航氣象研究所中國民航大學，天津300300;2.智能信號與圖像處理天津市重點實驗室中國民航大學，天津300300)

1 引言

風切變表現(xiàn)為空間兩點之間風的矢量差，即在統(tǒng)一高度或不同高度短距離內(nèi)風向和風速的變化，飛機在風向風速不斷發(fā)生變化的環(huán)境中飛行，相應地就要發(fā)生突然性的空速變化，空速變化引起了升力變化，升力變化又引起了飛行高度的變化。如果遇到的是空速突然減小，而飛行員又未能立即采取措施，飛機就要掉高度，以至發(fā)生事故。通過對嚴重事故的分析確認后，低空風切變是引起這些飛機失事的主要原因。

目前對于低空風切變的識別這方面的課題研究還不是很成熟，如何有效地識別風切變?nèi)找娉蔀楸Ｕ厦窈桨踩\輸?shù)囊粋€重要課題，國內(nèi)外在此方面也取得了一定的進展。Potts認為下?lián)舯┝髯羁煽康奶卣魇堑孛孑椛ⅲ诙嗥绽账俣葓D上表現(xiàn)為在雷達徑向上的“牛眼”回波［1］。1987年Merritt根據(jù)下?lián)舯┝鞯竭_地面后向四周輻散的特性，提出了基于切變段提取的下?lián)舯┝魈綔y算法［2］。Hermes在研究自動識別、追蹤陣風鋒的方法中指出，當陣風鋒與雷達徑向不垂直時，輻合區(qū)在多普勒雷達上顯示為速度沿徑向減小的區(qū)域［3］。

近年來，由于激光雷達在探測低空風切變中具有精度高、分辨率高等特點，基于二維圖像的方法對風切變進行識別成為研究的新思路，高志光采用了Hu矩和灰度共生矩陣的串聯(lián)組合特征對風切變進行了識別［4］。然而Hu矩和灰度共生矩陣都只能刻畫圖像的全局特征，在同類風切變風向變化較大，尤其是風向與雷達徑向垂直時并不能得到理想的識別效果。

文中對激光雷達探測的二維風切變圖像進行了紋理分析，將風切變的紋理特征分為局部特征和全局特征，局部的紋理特征可以提取到風場局部風速變化較為微弱的信息，全局的紋理特征代表風場全局風速的變化，而對局部風速變化小的區(qū)域較為不敏感，二者都包含有代表風場類別的信息，文中采用典型相關分析對二者進行融合，最后通過最近鄰分類器進行匹配識別。識別的流程如圖1所示。

圖1 風切變識別的流程圖Fig.1 Flow chart of the proposed method

2 樣本的仿真

由于真實的風切變數(shù)據(jù)難以獲取，本節(jié)首先基于大氣流場理論，運用Fluent流場軟件對微下?lián)舯┝?、低空急流、海陸風進行了仿真模擬。進而通過多普勒激光雷達速度方位(VAD)掃描方式得到三種風切變的仿真激光雷達圖像。

2.1 風切變風場的建模

低空風切變是指600 m以下空氣層中風向和風速突然改變的現(xiàn)象。三種低空風切變的風場特性可概括如下:

下?lián)舯┝魇抢妆┰浦械膹娏蚁鲁翚饬?，在沖擊地面后向各個方向輻散。

低空急流是對流層下層的一個準管狀的狹長強風帶，在強風帶中心風速最大，上下和兩側(cè)逐漸減小。

海陸風切變是指海風從海面吹向大陸，與大陸上相反方向的背景風匯聚而形成的風切變，在海風鋒兩側(cè)風向相反。

根據(jù)以上三種低空風切變的特性，分別選用“壁面射流”模型和“平面自由射流”模型對風場進行數(shù)學描述，選用相應的湍流模型和邊界條件，具體參見文獻［5］，通過Fluent數(shù)值模擬仿真后便可得到三種低空風切變的三維風場，三種風場的尺寸均為8000 m×8000 m×1000m，如圖2所示。

圖2 三種風切變的三維流場仿真圖Fig.2 The three － dimensional field simulation figure

2.2 仿真激光雷達圖像的生成

多普勒激光雷達速度方位掃描(VAD)方式［6］如圖3所示。

圖3 雷達圓錐掃描示意圖Fig.3 The scanning way of lidar

雷達波束在一個固定仰角α上，沿方向角θ進行掃描。設去向的徑向速度為正，來向的徑向速度為負，雷達位于坐標原點。其中，Vk為水平風，Vf為垂直風，θ0為水平風和X軸的夾角。

由圖4可知，徑向速度v:

式(1)描述了雷達徑向風速與三維實際風場各個方向上風速分量的關系。這樣就可以由仿真得到的三維風場數(shù)據(jù)得到基于激光雷達的變化風場探測數(shù)據(jù)。待測風場的激光雷達PPI圖像如圖4所示，周圍的藍色區(qū)域代表無有效地風場數(shù)據(jù)，圖像的高度代表雷達徑向距離庫，水平寬度代表方位角，單位庫長為100 m，方位角分辨率為1°。

圖4 三種風切變的激光雷達PPI圖像Fig.4 The lidar PPI image of three kinds of wind shear

3 組合紋理特征提取

3.1 局部紋理特征

局部二進制模式(Local Binary Pattern，LBP)［7］是一種基于灰度描述圖像紋理特征的不相關算子，它通過對圖像任意一點與其周圍點的灰度值大小關系來表征圖像的局部紋理特征，具有對灰度值的線性變化不敏感的特性。

LBP算子通常由參數(shù)(P，R)來表示，其中P表示鄰域內(nèi)包含的像素個數(shù)，R表示鄰域半徑，對于不同的(P，R)值，LBP算子也不相同，圖5為兩種不同LBP算子。

圖5 兩種不同的LBP算子Fig.5 Two different LBP operator

對于任意的LBP算子，它的編碼公式為:

式(2)中，gc為中心像素;gi(i=0，…，P －1)為以中心像素為圓心;R為半徑的圓周上的P個像素值。LBP算子的計算過程如圖6所示。

圖6 LBP算子的計算過程Fig.6 The calculation process of LBP operator

為了得到旋轉(zhuǎn)不變的LBP算子，文獻［8］對原始的LBP算子做了改進，即不斷旋轉(zhuǎn)圓形鄰域得到一系列的LBP值，取其最小值作為該鄰域的LBP值，計算公式如下:

式(3)中，LBPriP，R表示旋轉(zhuǎn)不變得 LBP 算子，ROR(x，i)為旋轉(zhuǎn)函數(shù)，表示將x循環(huán)右移i位。風切變圖像經(jīng)過LBP算子運算后，可得到對應的LBP編碼圖像，如圖7所示。對得到的LBP編碼做直方圖統(tǒng)計，將歸一化后統(tǒng)計直方圖做為最終的局部特征向量。

圖7 三種風切變的LBP編碼圖像Fig.7 The LBP coded image of three kinds of wind shear

3.2 全局紋理特征

灰度－梯度共生矩陣(Gray－Gradinent Co－occurrence Matrix，GGCM)是利用整幅圖像的灰度和梯度綜合信息來提取紋理特征?；叶龋荻裙采仃嚨脑豀(i，j)定義為在歸一化的灰度圖像F(m，n)和歸一化的梯度圖像G(m，n)中共同具有灰度為i和梯度為j的總像素點數(shù)。以圖像的總像素點數(shù)歸一化后得到概率為p(i，j)。從灰度－梯度共生矩陣中，可以提取15個二次特征參量［9］。

在雷達風切變圖像中，灰度代表徑向風速值，梯度代表徑向風速的切變值，為了克服傳統(tǒng)的梯度算子對噪聲敏感的不足，文中采用最小二乘擬合法分別計算雷達的徑向速度的徑向切變和徑向速度的切向切變(除以r是為了要歸一化)，并將徑向方位組合切變G做為最終的梯度值，計算公式如下:

式(4)中，vi為第i個點的徑向速度;ri為相應點距雷達的距離;n為參與計算的同一徑向上的像素點數(shù)。

式(5)中，vj為第j個點的徑向速度;θj為相應地方位角;r為距離雷達的距離;m為參與計算的相同徑向距離上像素點數(shù)。

由最小二乘求得的三種風切變圖像的切變值圖像如圖8所示。

圖8 三種風切變的切變值圖像Fig.8 The shear value images obtained by LS

3.3 特征融合

典型相關分析(Canonical Correlation Analysis，CCA)是一種提取兩個特征集間線性相關的強有力的工具。通過局部特征和全局特征的提取，我們得到兩組特征集:X=［x1，…，xn］，Y=［y1，…，yn］，其中 xi∈Rp，yi∈Rq，i=1，2，…，n。n 為樣本的個數(shù)。CCA的目標可看為分別對特征集X和特征集Y尋找兩組基向量wx∈Rp和wy∈Rq。使得隨機變量u=wTxx和v=wTyy之間的相關性達到最大，即相關系數(shù)ρ(u，v)最大。CCA可以通過下面的式子描述:

式(7)是一種最優(yōu)化問題，可以通過Lagrange乘子法求得最優(yōu)解(wxi，wyi)，i=1，…，d(d≤min(p，q))，具體的求解方法參見文獻［10］，這樣便求出了兩組投影矩陣 Wx=［wx1，…，wxd］，Wy=［wy1，…，wyd］，融合后的特征可以定義為:

其中，ti即為全局紋理特征和局部紋理特征融合后的組合紋理特征(Global－Local Texture Features，GLTF)。

4 應用最近鄰分類器進行分類識別

從激光雷達風切變圖像中提取到LBP特征，灰度－梯度共生矩陣特征并構(gòu)成融合特征后，可以通過最近鄰分類器進行分類識別。

最近鄰分類器是按距離來分類的，這種分類器較簡單且花費時間很小，其分類的思想為:如果未知樣本f的特征和已知類別樣本fk的特征之間的歐氏距離最小，則判斷f與fk同類，假設有c個類別w1，w2，…，wc，每類有表明類別的樣本 Ni個，則可以規(guī)定wi類的判別函數(shù)為:

決策規(guī)則為:若gj(f)=(f)，i=1，2，…，c，則 f∈wj。

5 結(jié)果分析

文中采用自建立的雷達掃描圖像樣本庫對常見的三種低空風切變進行分類識別，每種激光雷達風切變圖像80幅，分別為待測風場位于激光雷達不同方位下的掃描圖像，總計240張圖像，每幅圖像的尺寸分別代表徑向距離和方位角。識別過程采用交叉驗證的策略，以避免人為選擇樣本所導致的識別結(jié)果的偶然性，即隨機選取20幅作為訓練樣本，剩下的作為測試樣本，重復進行50次。

為了突出文中算法的優(yōu)勢，將其與文獻［4］中采用的Hu+GLCM特征、LBP特征、GGCM特征所得的識別率進行對比。Hu+GLCM包括Hu矩的7個形狀特征和GLCM的4個紋理特征;LBP算法采用p=8，R=1的旋轉(zhuǎn)不變LBP算子，包含36個特征;GGCM包括小梯度優(yōu)勢、大梯度優(yōu)勢、灰度分布的不均勻性、梯度平均、灰度熵、混合熵、慣性、逆差矩等15個全局統(tǒng)計特征;GLTF由 LBP特征和GGCM特征經(jīng)CCA融合而成，特征維數(shù)為15。實驗結(jié)果如圖9和表1所示。

圖9 四種算法的識別率曲線Fig.9 The recognition rate curve of four algorithms

表1 四種算法識別率的對比Tab.1 Recognition rate comparison of the four algorithms %

圖9給出了50次分別應用Hu+GLCM特征，LBP特征，GGCM特征和GLTF特征進行風切變識別時的識別曲線，表1給給出了50次識別的平均識別率，S1、S2、S3分別代表微下?lián)舯┝鳌⒌涂占绷?、海陸風三種風切變，“±”后面的數(shù)字是識別率的標準差。

從表1中可以看出，對三種低空風切變的識別中，Hu+GLCM算法、LBP算法和GGCM算法的平均識別率只有83.31%，93.52%和 92.53%，運用文中所述的組合紋理特征進行識別，平均識別率達到99.02%，識別效果得到了很大改善。

6 結(jié)論

本文提出了基于組合紋理的低空風切變算法，從激光雷達探測到的低空風切變圖像中提取LBP特征和GGCM特征，通過CCA融合得到組合紋理特征，LBP特征反應圖像局部的紋理，可以提取到風場局部變化較為微弱的信息，GGCM特征反應圖像的全局紋理，代表風場全局風速的變化。通過對三種低空風切變的識別，該算法取得了較高的識別率，而且特征維數(shù)較少，達到了較好的分類識別效果。

［1］ Potts R.Microburst precursors observed with Doppler radar［C］//Conference on Radar Meteorology，24th，Tallahassee，F(xiàn)L，1989:158 －162.

［2］ Merritt M W.Automated detection of microburst windshear for terminal doppler weather radar［C］//1987 Cambridge Symposium.International Society for Optics and Photonics，1988:61 －69.

［3］ Hermes L G，Witt A，Smith SD，et al.The gust－ front detection and wind－shift algorithms for the terminal doppler weather radar system［J］.Journal of Atmospheric and O-ceanic Technology，1993，10(5):693 －709.

［4］ Jiang Lihui，Gao Zhiguang，Xiong Xinglong，et al.The study on type recognition of low altitude wind shear based on lidar image processing.［J］.Infrared and Laser Engineering，2012，41(12):3410 －3415.(in Chinese)蔣立輝，高志光，熊興隆.基于激光雷達圖像處理的低空風切變類型識別研究［J］.紅外與激光工程，2012，41(12):3410－3415.

［5］ Liu Shuo.The research on numerical simulation of high rise building outdoor flow field［D］.Harbin:Harbin Institute of Technology，2007.(in Chinese)劉朔.高層建筑室外氣流場的數(shù)值模擬研究［D］.哈爾濱:哈爾濱工業(yè)大學，2007.

［6］ Wang Bangxin，Shen Fahua，Sun Dongsong，et al.Beam scanning and wind field measurement of direct－detection doppler lidar［J］.Infrared and Laser Engineering，2007，36(1):69 －72.(in Chinese)王邦新，沈法華，孫東松，等.直接探測多普勒激光雷達的光束掃描和風場測量［J］.紅外與激光工程，2007，36(1):69 －72.

［7］ Ojala T Pietikainen M，Harwood D.A comparative study of texture measures with classification based on feature distribution［J］.Pattern Recognition，1996，29(1):51－59.

［8］ Ojala T Pietikainen M，Maenpaa T.Multiresolution gray scale and rotation invariant texture analysis with local binary patterns［J］.IEEE Transactions on Pattern Analysis and Machine Intelligence，2002，24(7):971 －987.

［9］ Xia Deshen，Jin Sheng，Wang Jian.Fractal dimension and GGCM meteorology cloud pictures recognition［J］.Journal of Nanjing University of Science and Technology，1999，23(4):289 －292.(in Chinese)夏德深，金盛，王建.基于分數(shù)維與灰度梯度共生矩陣的氣象云圖識別［J］.南京理工大學學報，1999，23(4):289－292.

［10］ DR Hardoon，SSzedmak，JShawe － Taylor.Canonical correlation analysis:an overview with application to learning methods［J］.Neural Computation，2004，16(12):2639－2664.