李 海 孫 研 孟凡旺 白 錦 韓雁飛

(1.中國民航大學(xué)天津市智能信號(hào)與圖像處理重點(diǎn)實(shí)驗(yàn)室 天津 300300；2.中國航空工業(yè)集團(tuán)公司雷華電子技術(shù)研究所 江蘇無錫 214063；3.中國民用航空華北地區(qū)空中交通管理局天津分局 天津 300300)

0 引言

飛機(jī)在巡航階段(一般在22000ft以上)飛行時(shí)，高空冰晶是對(duì)飛機(jī)飛行安全威脅較大的天氣現(xiàn)象。在冰晶高聚集區(qū)域，發(fā)動(dòng)機(jī)葉片表面的冰晶堆積會(huì)造成發(fā)動(dòng)機(jī)失速與喘振，燃燒室熄火和損壞等。同時(shí)，冰晶還會(huì)進(jìn)入大氣數(shù)據(jù)傳感器，造成堵塞使得相關(guān)傳感器數(shù)據(jù)異常，造成飛機(jī)失控。近年來因遭遇高空冰晶而發(fā)生的航空飛行事故時(shí)有發(fā)生。美國聯(lián)邦航空管理局對(duì)近二十年來發(fā)生的發(fā)動(dòng)機(jī)空中事故進(jìn)行統(tǒng)計(jì)，統(tǒng)計(jì)結(jié)果表明這些事故大部分與高空冰晶相關(guān)。這些真實(shí)的事件表明高空冰晶對(duì)航空運(yùn)輸業(yè)和旅客飛行安全會(huì)產(chǎn)生嚴(yán)重影響，開展高空冰晶檢測(cè)的研究是十分重要的，而高空冰晶冰水含量估計(jì)是整個(gè)檢測(cè)流程的基礎(chǔ)，直接影響后續(xù)危險(xiǎn)級(jí)別判定的準(zhǔn)確性。因此針對(duì)高空冰晶冰水含量估計(jì)的研究具備重要的意義。

機(jī)載氣象雷達(dá)能探測(cè)飛機(jī)飛行前方空域的危險(xiǎn)氣象條件，可以用作飛行預(yù)警或規(guī)避危險(xiǎn)氣象的重要參考。一般情況下，機(jī)載氣象雷達(dá)探測(cè)降水目標(biāo)時(shí)是根據(jù)雷達(dá)回波的強(qiáng)度來區(qū)分的，并且認(rèn)為回波強(qiáng)度越大，含水成分越多，對(duì)飛行安全的影響越大。而高空冰晶檢測(cè)及其危險(xiǎn)級(jí)別量化可以通過所探測(cè)空間的冰水含量(Ice Water Content，IWC)，即單位立體空間的冰晶密度予以指示。IWC是指單位體積空氣中包含的固態(tài)水的質(zhì)量，單位為g/m。

國外學(xué)者提出可根據(jù)雷達(dá)反射率因子Z與IWC之間存在的統(tǒng)計(jì)關(guān)系，結(jié)合經(jīng)驗(yàn)公式通過雷達(dá)反射率因子反演估計(jì)IWC值。1954年空中劍橋研究中心的Atlas根據(jù)飛機(jī)實(shí)測(cè)的譜參數(shù)數(shù)據(jù)，獲得了雷達(dá)反射率因子-粒子分布-IWC之間的經(jīng)驗(yàn)關(guān)系，1992年美國的Liao等研究了Ka波段雷達(dá)和W波段雷達(dá)的Z-IWC關(guān)系，1995年Brown用毫米波機(jī)載雷達(dá)探測(cè)了冰晶云的Z-IWC經(jīng)驗(yàn)公式，2001年Austin等人通過采用數(shù)值模擬的方法，提出了一種基于毫米波雷達(dá)反射率和可見光層云微物理參數(shù)反演冰水含量的方法，2004年Matrosov等人研究了基于雷達(dá)反射率因子的估測(cè)海洋性層云含水量的方法，并比較了取不同的降水云與非降水云強(qiáng)度界限時(shí)云內(nèi)含水量的變化，2015年大氣介質(zhì)空間觀測(cè)實(shí)驗(yàn)室的Protat分別給出了35GHz 和95GHz雷達(dá)下的Z-IWC關(guān)系，并統(tǒng)計(jì)對(duì)比了實(shí)驗(yàn)誤差。但由于經(jīng)驗(yàn)公式形式過于簡(jiǎn)單具有較大局限性，難以準(zhǔn)確表征復(fù)雜的Z-IWC關(guān)系模型，特別是使用不同數(shù)據(jù)進(jìn)行冰水含量估計(jì)時(shí)會(huì)存在大誤差，魯棒性較差，因此如何獲得精確的Z-IWC關(guān)系模型是進(jìn)行高空冰晶冰水含量準(zhǔn)確估計(jì)亟待解決的問題。

針對(duì)數(shù)學(xué)模型相對(duì)復(fù)雜，所需確定參數(shù)較多的Z-IWC關(guān)系模型，利用一些智能算法進(jìn)行自適應(yīng)的參數(shù)估計(jì)從而得到準(zhǔn)確的模型是解決問題的新方案。遺傳算法是一種在處理各種復(fù)雜優(yōu)化問題時(shí)效果顯著，且具備較強(qiáng)魯棒性的算法。從1995年開始，國內(nèi)對(duì)遺傳算法的研究也日益增加，作為啟發(fā)式搜索算法的一員，遺傳算法在最近幾十年內(nèi)在信號(hào)處理、最優(yōu)控制、生產(chǎn)規(guī)劃等各個(gè)領(lǐng)域都有著長(zhǎng)足的進(jìn)步與發(fā)展。遺傳算法在求解函數(shù)優(yōu)化問題時(shí)，相比其它優(yōu)化算法，最大的特點(diǎn)是處理對(duì)象不是參數(shù)本身，而是對(duì)參數(shù)集編碼獲得基因個(gè)體，因此該算法不會(huì)受到函數(shù)連續(xù)性的限制。另外該算法不需要問題的梯度信息等輔助信息，僅根據(jù)適應(yīng)度函數(shù)的數(shù)值大小就可以選出最優(yōu)個(gè)體，且適應(yīng)度函數(shù)也不會(huì)受到連續(xù)可微的條件限制，所以適用范圍更廣。針對(duì)高空冰晶冰水含量估計(jì)問題，如能設(shè)計(jì)遺傳算法進(jìn)行參數(shù)估計(jì)，發(fā)揮其魯棒性較好的優(yōu)點(diǎn)，確定較為精確的Z-IWC關(guān)系模型將具有很好的應(yīng)用價(jià)值。

基于此，本文提出了一種基于遺傳算法的機(jī)載氣象雷達(dá)高空冰晶冰水含量估計(jì)方法。該方法對(duì)雷達(dá)反射率因子與高空冰晶冰水含量關(guān)系模型中待定參數(shù)進(jìn)行精細(xì)尋優(yōu)得到最優(yōu)解，然后使用估計(jì)的參數(shù)得到確定的模型后利用測(cè)試數(shù)據(jù)的雷達(dá)反射率因子進(jìn)行反演，估計(jì)高空冰晶的冰水含量。針對(duì)2016年3月17日美國俄克拉荷馬州地區(qū)一次真實(shí)氣象場(chǎng)景的數(shù)值天氣預(yù)報(bào)模式輸出數(shù)據(jù)進(jìn)行實(shí)驗(yàn)，利用遺傳算法估計(jì)冰水含量，對(duì)比遺傳算法與經(jīng)驗(yàn)公式法冰水含量估計(jì)結(jié)果，并給出不同時(shí)間點(diǎn)、不同高度層的冰水含量估計(jì)結(jié)果，與真實(shí)值的對(duì)比結(jié)果表明，該方法可有效地估計(jì)高空冰晶的冰水含量，估計(jì)結(jié)果較為準(zhǔn)確，統(tǒng)計(jì)分析的指標(biāo)參量(包括相關(guān)性系數(shù)、RMSE、誤差統(tǒng)計(jì)占比等)較為理想。

1 雷達(dá)反射率因子與冰水含量關(guān)系模型

雷達(dá)反射率因子與粒子直徑的6次方成正比，一般用來反映單位體積內(nèi)云、冰晶、雨滴等粒子群的濃度情況。冰水含量(IWC)反映了單位立體空間的冰晶密度。對(duì)雷達(dá)反射率因子與冰水含量關(guān)系進(jìn)行建模是利用雷達(dá)反射率因子反演冰水含量的關(guān)鍵。

根據(jù)Heymsfield等人的研究，冰晶粒子群的等效雷達(dá)反射率因子和冰水含量表達(dá)式分別如公式(1)和公式(2)所示。

(1)

(2)

由于公式(1)與公式(2)均含有斜率參量成分，因此可利用共有成分將兩公式進(jìn)行聯(lián)立從而得到Z-IWC關(guān)系模型。對(duì)公式(1)進(jìn)行變換，提出斜率參量成分可得

(3)

對(duì)公式(2)進(jìn)行變換，提出斜率參量成分可得

(4)

將公式(3)與公式(4)進(jìn)行聯(lián)立可得到

(5)

對(duì)公式(5)進(jìn)行變換可得

(6)

將確定的、、等參數(shù)帶入公式(6)可得

(7)

公式(7)即為求得的Z-IWC關(guān)系模型，模型中的待定參數(shù)為、、。在Z-IWC模型中，只要知道了、、三個(gè)參數(shù)，就可以得到準(zhǔn)確的Z-IWC模型，即可以通過雷達(dá)反射率因子準(zhǔn)確反演冰水含量。其中，是冰晶質(zhì)量?jī)缏煞植己瘮?shù)(D)中表征冰晶粒子成分的參數(shù)；是冰晶粒子數(shù)濃度分布函數(shù)()中表征冰晶粒子尺寸的參數(shù)；是冰晶粒子數(shù)濃度分布函數(shù)()中決定分布形狀的參數(shù)。

2 基于遺傳算法的Z-IWC參數(shù)估計(jì)

本文使用遺傳算法對(duì)Z-IWC模型中的三個(gè)未知參數(shù)、、進(jìn)行參數(shù)估計(jì)，以確定具體的Z-IWC關(guān)系模型，從而完成冰水含量的估計(jì)。遺傳算法的核心部分包括：

1)種群初始化與編碼；

2)適應(yīng)度函數(shù)的選?。?/p>

3)遺傳選擇；

4)交叉和變異；

5)迭代終止。

利用遺傳算法進(jìn)行參數(shù)估計(jì)的具體實(shí)現(xiàn)流程如圖1所示。下面對(duì)利用遺傳算法進(jìn)行Z-IWC模型參數(shù)估計(jì)的主要操作過程進(jìn)行具體論述。

圖1 基于遺傳算法的Z-IWC模型參數(shù)估計(jì)流程圖

1)編碼與種群初始化。遺傳算法進(jìn)行參數(shù)估計(jì)，首先對(duì)要優(yōu)化問題的解進(jìn)行編碼，編碼的目的主要是為了改變優(yōu)化問題解的表現(xiàn)形式，使其方便遺傳運(yùn)算。本文采用二進(jìn)制編碼的方式對(duì)Z-IWC關(guān)系模型中的三個(gè)參數(shù)、、分別進(jìn)行編碼，將每一個(gè)染色體(解)編碼為基因型。完成編碼操作后，接著需要設(shè)置初始種群。對(duì)于模型中每一個(gè)待定參數(shù)，本文通過隨機(jī)生成的形式，來創(chuàng)造一個(gè)初始種群，文中是隨機(jī)生成具有均勻分布的200個(gè)個(gè)體，作為初始種群。

2)適應(yīng)度函數(shù)的選取。遺傳算法的收斂速度以及能否找到最優(yōu)解，直接依賴于適應(yīng)度函數(shù)的選取。為了使估計(jì)的冰水含量值更加接近冰水含量真值，本文構(gòu)造滿足下列優(yōu)化準(zhǔn)則的目標(biāo)函數(shù)(,,)作為適應(yīng)度函數(shù)如式(8)所示。

(8)

式(8)中，IWC(,,)為利用氣象場(chǎng)景數(shù)據(jù)空氣密度、水汽混合比等參數(shù)計(jì)算出來的冰水含量值，作為IWC的真實(shí)值，IWC′(,,)為利用式(7)的Z-IWC關(guān)系模型求得的冰水含量估計(jì)值，為使用的樣本個(gè)數(shù)。

4)遺傳交叉與變異。遺傳算法的收斂性，主要取決于變異和交叉兩大核心操作。變異的實(shí)質(zhì)是在群體中發(fā)掘個(gè)體的多樣性，提高算法的局部隨機(jī)搜索能力，而交叉操作決定產(chǎn)生新個(gè)體的速度，對(duì)實(shí)現(xiàn)算法的全局搜索性能具有關(guān)鍵性的作用。為了更好地完成交叉和變異操作，以盡可能得到最優(yōu)解，本文設(shè)置交叉概率為06，設(shè)置變異概率為008。

5)迭代終止。本文設(shè)置進(jìn)化代數(shù)計(jì)數(shù)器=0，設(shè)置最大進(jìn)化代數(shù)=100，運(yùn)行終止條件設(shè)置為若優(yōu)化準(zhǔn)則目標(biāo)函數(shù)的目標(biāo)值滿足()≤=001，或者運(yùn)行迭代次數(shù)=100次，則停止遺傳進(jìn)化過程，選取最后一代的最優(yōu)個(gè)體作為最優(yōu)解，以得到確定的Z-IWC模型。利用確定的模型對(duì)雷達(dá)數(shù)據(jù)進(jìn)行處理，得到高空冰晶冰水含量估計(jì)結(jié)果。

3 冰水含量估計(jì)方法流程與步驟

基于遺傳算法的機(jī)載氣象雷達(dá)高空冰晶冰水含量估計(jì)方法具體實(shí)現(xiàn)步驟如下：

1)步驟1：對(duì)雷達(dá)數(shù)據(jù)進(jìn)行質(zhì)量控制；

2)步驟2：建立高空冰晶冰水含量與雷達(dá)反射率關(guān)系模型(Z-IWC關(guān)系模型)；

3)步驟3：利用遺傳算法對(duì)Z-IWC關(guān)系模型中的未知參數(shù)進(jìn)行精細(xì)尋優(yōu)；

4)步驟4：根據(jù)確定后的參數(shù)得到最終確定的Z-IWC關(guān)系模型；

5)步驟5：利用確定的Z-IWC關(guān)系模型對(duì)雷達(dá)數(shù)據(jù)進(jìn)行計(jì)算，反演IWC得到高空冰晶冰水含量估計(jì)結(jié)果。

總體流程框圖如圖2所示。

圖2 總體流程圖

4 實(shí)驗(yàn)及結(jié)果分析

為了驗(yàn)證所提方法性能，本文將2016年3月17日美國俄克拉荷馬州地區(qū)西經(jīng)95°至100°，北緯30°至34°范圍內(nèi)的一次真實(shí)氣象場(chǎng)景的數(shù)值天氣預(yù)報(bào)模式輸出數(shù)據(jù)作為樣本進(jìn)行檢測(cè)分析，分別選取俄克拉荷馬州地區(qū)2016年3月17日14時(shí)與2016年3月17日16時(shí)兩個(gè)時(shí)間點(diǎn)的數(shù)據(jù)進(jìn)行實(shí)驗(yàn)。

為了驗(yàn)證所提方法的有效性和準(zhǔn)確性，針對(duì)兩個(gè)不同時(shí)刻點(diǎn)的不同高度層的數(shù)據(jù)分別給出了本文方法進(jìn)行冰水含量估計(jì)與經(jīng)典經(jīng)驗(yàn)公式法進(jìn)行冰水含量估計(jì)的實(shí)驗(yàn)結(jié)果對(duì)比圖，如圖3和圖4所示。本文方法與1954年Atlas等人提出的經(jīng)驗(yàn)公式IWC=0.064Z，1997年Aydin等人提出的經(jīng)驗(yàn)公式IWC=0.104Z以及2002年Mace等人提出的經(jīng)驗(yàn)公式IWC=0.104Z進(jìn)行對(duì)比。從圖3和圖4的實(shí)驗(yàn)結(jié)果中可以看出，經(jīng)驗(yàn)公式法因受經(jīng)驗(yàn)取值的影響，冰水含量估計(jì)結(jié)果均產(chǎn)生了不同程度的偏差，誤差較大，而采用遺傳算法完成冰水含量估計(jì)，誤差較小，具有更優(yōu)越的魯棒性。

圖3 2016年3月17日14時(shí)不同高度層Z-IWC參數(shù)估計(jì)結(jié)果對(duì)比圖

圖4 2016年3月17日16時(shí)不同高度層Z-IWC參數(shù)估計(jì)結(jié)果對(duì)比圖

從真值與估值相關(guān)性的角度，下面給出了俄克拉荷馬州地區(qū)2016年3月17日14時(shí)與2016年3月17日16時(shí)真值IWC與估值IWC的相關(guān)性散點(diǎn)圖，分別如圖5和圖6所示。圖中橫坐標(biāo)是利用遺傳算法估計(jì)得到的IWC值，縱坐標(biāo)是IWC真值，當(dāng)估計(jì)值完全等于真實(shí)值時(shí)數(shù)據(jù)點(diǎn)應(yīng)該全部分布在對(duì)角線上，那么越靠近對(duì)角線則說明估計(jì)的越準(zhǔn)確，從圖5和圖6可以看出數(shù)據(jù)點(diǎn)都在對(duì)角線附近，說明誤差較小，估值與真值之間具有很強(qiáng)的相關(guān)性，證明了本文方法進(jìn)行高空冰晶冰水含量估計(jì)的準(zhǔn)確性。

圖5 2016年3月17日14時(shí)不同高度層真值與估值冰水含量相關(guān)性散點(diǎn)圖

圖6 2016年3月17日16時(shí)不同高度層真值與估值冰水含量相關(guān)性散點(diǎn)圖

在給出真值IWC與估值IWC相關(guān)性散點(diǎn)圖的同時(shí)，對(duì)利用本文方法估計(jì)得到的IWC值的誤差大小及出現(xiàn)的概率情況進(jìn)行統(tǒng)計(jì)，分別給出2016年3月17日14時(shí)和2016年3月17日16時(shí)不同高度層數(shù)據(jù)的誤差統(tǒng)計(jì)直方圖如下圖7和圖8所示。圖7和圖8中的誤差為估計(jì)的IWC值減去IWC真值的絕對(duì)誤差，由圖可以看出兩個(gè)時(shí)刻點(diǎn)不同高度層數(shù)據(jù)的真值與估計(jì)值的誤差出現(xiàn)在0附近的概率較大，均占到80%左右，故驗(yàn)證了算法的準(zhǔn)確性。

圖7 2016年3月17日14時(shí)不同高度層誤差分布直方圖

圖8 2016年3月17日16時(shí)不同高度層誤差分布直方圖

另外，分別統(tǒng)計(jì)本文方法與經(jīng)驗(yàn)公式法進(jìn)行冰水含量估計(jì)時(shí)估值與真值之間的相關(guān)性系數(shù)與均方根誤差(Root Mean Square Error，RMSE)，2016年3月17日14時(shí)數(shù)據(jù)的統(tǒng)計(jì)結(jié)果如表1和表2所示，2016年3月17日16時(shí)數(shù)據(jù)的統(tǒng)計(jì)結(jié)果如表3和表4所示。由此可見本文方法與經(jīng)驗(yàn)公式法相比，相關(guān)性系數(shù)更接近1，表明具有很強(qiáng)的相關(guān)性。同時(shí)本文方法的RMSE都小于0.1，說明估計(jì)值與真值之間偏差較小，而經(jīng)驗(yàn)公式法的RMSE偏差較大，對(duì)比表明本文算法準(zhǔn)確性更高。

表1 IWC估值與真值之間的相關(guān)性系數(shù)(2016年3月17日14時(shí)數(shù)據(jù))

表2 IWC估值與真值之間的RMSE(2016年3月17日14時(shí)數(shù)據(jù))

表3 IWC估值與真值之間的相關(guān)性系數(shù)(2016年3月17日16時(shí)數(shù)據(jù))

表4 IWC估值與真值之間的RMSE(2016年3月17日16時(shí)數(shù)據(jù))

5 結(jié)束語

本文提出了一種基于遺傳算法的機(jī)載氣象雷達(dá)高空冰晶冰水含量估計(jì)方法。首先對(duì)高空冰晶冰水含量與雷達(dá)反射率因子關(guān)系進(jìn)行建模，得到理論推導(dǎo)的基礎(chǔ)模型之后，利用遺傳算法對(duì)模型中待定參數(shù)進(jìn)行精細(xì)尋優(yōu)獲得最優(yōu)解，然后利用估計(jì)的參數(shù)得到確定的模型后對(duì)測(cè)試數(shù)據(jù)的雷達(dá)反射率因子進(jìn)行反演，估計(jì)高空冰晶的冰水含量。通過對(duì)仿真數(shù)據(jù)進(jìn)行試驗(yàn)驗(yàn)證，利用遺傳算法估計(jì)冰水含量，對(duì)比遺傳算法與經(jīng)驗(yàn)公式法冰水含量估計(jì)結(jié)果，并給出不同時(shí)間點(diǎn)、不同高度層的冰水含量估計(jì)結(jié)果，與真實(shí)值的對(duì)比結(jié)果表明，該方法可有效地估計(jì)高空冰晶的冰水含量，估計(jì)結(jié)果較為準(zhǔn)確，統(tǒng)計(jì)分析的指標(biāo)參量(包括相關(guān)性系數(shù)、RMSE、誤差統(tǒng)計(jì)占比等)較為理想。