飛機的紅外圖像仿真

王 彪，童中翔，王超哲，馬 榜

(空軍工程大學(xué)航空航天工程學(xué)院，陜西 西安710038)

1 引言

隨著紅外制導(dǎo)技術(shù)的不斷發(fā)展，紅外制導(dǎo)導(dǎo)彈的抗干擾能力不斷提高，目前最新一代的紅外制導(dǎo)導(dǎo)彈普遍采用了紅外成像制導(dǎo)技術(shù)［1］。紅外成像制導(dǎo)探測器能夠利用目標(biāo)的形狀特征、灰度分布規(guī)律、運動特點等諸多物理特征形成對目標(biāo)的識別，進而形成制導(dǎo)信號，其抗干擾能力得到了極大提高，因此，對目標(biāo)的紅外圖像特征進行研究已變得十分重要。

目前紅外圖像的獲取方式有三種，一是在真實條件下進行實際拍攝，二是在地面進行模擬測試，三是理論建模與仿真。理論建模與仿真法與前兩種方法相比，其所需經(jīng)費較少，且不易受場地、設(shè)備、天氣等因素的限制，可以在一定的精度下模擬出飛機在任意飛行條件下的紅外圖像，在研究目標(biāo)紅外圖像的過程中具有重要意義。

針對飛行器的紅外圖像仿真方法，廣大專家學(xué)者進行了廣泛的研究。在機體紅外圖像方面，文獻［2］借助Fluent軟件生成了蒙皮的紅外圖像，文獻［3］、［4］通過建立熱平衡方程生成了機體紅外圖像，但建立的模型不能反映飛機外形對氣動加熱的影響。在排氣系統(tǒng)紅外圖像方面，文獻［5］對發(fā)動機熱噴流進行了紅外成像仿真，文獻［6］采用粒子系統(tǒng)對尾焰圖像進行了仿真研究。然而由于飛機在空中受到多種因素的影響，其紅外輻射模型的建立比較復(fù)雜，模型的精確性與實時性往往不能同時滿足，因此在建立仿真模型的過程中，如何快速地生成滿足精度要求的圖像是模型的關(guān)鍵。本文將整機輻射分解為機體輻射與排氣系統(tǒng)輻射兩部分，建立了一種可用于任意視線方向的飛機紅外圖像的生成模型，并利用OpenGL圖形庫生成了可視化的紅外圖像，通過驗證模型計算速度較快，且滿足一定的精確性。

2 機體紅外輻射計算

2．1 蒙皮溫度場建模

飛機在空中飛行的過程中要受到高速氣流的氣動加熱、內(nèi)熱源加熱、環(huán)境輻射熱量與蒙皮向外輻射熱量的綜合作用，如圖1所示，使得機體的各部位呈現(xiàn)出不均勻的溫度分布，因此，求解蒙皮的溫度場是生成機體紅外圖像的關(guān)鍵。

圖1 飛機熱交換示意圖

本文中首先將機體表面分成n個面元，認(rèn)為每個面元具有均勻的溫度與輻射量。綜合考慮各因素對面元溫度變化的影響，建立式(1)的熱平衡方程。

式中，mi為面元的質(zhì)量;ci為面元的比熱容;Ti為面元溫度;qcv，qin，qenv，qcon，qrad分別為氣動加熱、內(nèi)熱源加熱、環(huán)境輻射熱、相鄰面元之間的熱傳遞與蒙皮向外輻射熱量的熱流密度。

由于機體蒙皮很薄，面元之間的熱傳導(dǎo)面積很小，且蒙皮溫度場是連續(xù)的，相鄰面元之間的溫差不大，因此忽略面元之間的熱傳遞，即認(rèn)為:

2．1．1 氣動加熱

高速氣流流過飛機表面時，會因氣流壓縮與摩擦產(chǎn)生大量的熱量，熱量大小與面元所在位置的機體形狀相關(guān)。根據(jù)機體各部位形狀不同，將機體分為機身、機頭、翼前緣、翼面四個部分分別求解其氣動加熱大小。

a)機身及翼面氣動加熱計算

機身可近似等價為圓柱形，氣流的流動方向沿圓柱軸向，其熱流密度可按照平板氣動加熱計算。翼面的熱流密度也可近似按照平板氣動加熱公式計算，其熱流密度公式為:

式中，Ti為面元i的溫度;Tr與αox分別為恢復(fù)溫度和熱交換系數(shù)，其計算式為:

式中，Te，ue，ρe分別為附面層邊緣氣體的溫度、速度和密度;pr為普朗特數(shù);cp為空氣的定壓比熱;Rex為特征長度x處的雷諾數(shù);ρ*是參考溫度T*對應(yīng)的密度，可由氣體狀態(tài)方程求解:

其中，R為氣體常數(shù);T*可取Eckert參考溫度:

b)機頭氣動加熱計算

機頭可近似等價為圓錐形，其熱流密度公式和平板熱流密度公式的區(qū)別在于換熱系數(shù)不同，圓錐的換熱系數(shù)公式為［7］:

c)翼前緣氣動加熱計算

翼前緣可近似等價為一個無限后掠圓柱，其駐點線上的熱流密度qor與半球的駐點熱流密度qob有如下近似關(guān)系:

半球的駐點熱流密度qob可按照Kemp－Riddle經(jīng)驗公式計算:

式中，R0為半球半徑;ρsl為海平面大氣密度;Ts為駐點溫度。

對于后掠機翼，當(dāng)后掠角在一定范圍內(nèi)時(0°≤Λ≤60°)，對后掠角的修正因子如下［8］:

2．1．2 內(nèi)熱源加熱

內(nèi)熱源對蒙皮面元的加熱以兩種方式進行:熱輻射與通過金屬結(jié)構(gòu)的熱傳導(dǎo)。將內(nèi)熱源看作是溫度為T0、發(fā)射率為ε0、有一定形狀的輻射灰體，通過第一種方式輻射到面元i的熱流密度為:

式中，αi為蒙皮面元的吸收率;F0－i為內(nèi)熱源到面元i的輻射角系數(shù)。

內(nèi)熱源通過第二種方式傳導(dǎo)到面元i的熱流密度為:

式中，λc為綜合導(dǎo)熱系數(shù);l0為內(nèi)熱源到面元的熱傳導(dǎo)距離。

則內(nèi)熱源對蒙皮加熱的總熱流密度為:

2．1．3 環(huán)境輻射熱

對于高空飛行的飛機，天空散射太陽輻射、地面輻射及地面反射太陽輻射影響不大，本文中予以忽略，環(huán)境輻射對蒙皮溫度的影響只考慮太陽直射輻射與大氣輻射，即:式中，Esun為太陽常數(shù)(Esun=1353 W/m2);τh為從大氣層外邊界到飛行高度h的路徑上大氣對太陽輻射的透過率，計算方法可參考文獻［9］;θs為太陽直射光線與面元法線方向的夾角，當(dāng)0≤θs＜π/2時，qsun按式(17)計算，π/2＜θs≤π時，qsun=0。

2．1．4 蒙皮自身輻射熱量Ti

將蒙皮看作發(fā)射率為εi，溫度為Ti的輻射灰體，則蒙皮面元自身輻射熱量為:

2．1．5 熱平衡方程求解

根據(jù)以上描述，可完成式(1)中熱平衡方程的建立，方程的右側(cè)實際上是Ti的函數(shù)，其中含有T4

i、Ti、T－1i項及與Ti無關(guān)的常數(shù)項，記作:

其中，C1，C2，C3，C4為與Ti無關(guān)的系數(shù)。

因為本文所求的是穩(wěn)態(tài)時的溫度場，因此令dTi/dt=0，得到以下方程:

方程(20)雖然是一個5次方程，但根據(jù)其物理意義，方程有唯一解，且Ti＞T∞。因此在求解時，設(shè)定合適的門限值σ與步長B，以T∞為初始值，當(dāng)滿足式(21)時，認(rèn)為T∞+kB即為方程(20)的解。

2．2 機體紅外輻射亮度計算

得到蒙皮溫度場后，將蒙皮看作灰體，應(yīng)用普朗克定律和基爾霍夫定律可得蒙皮自身輻射的紅外輻射亮度，其計算公式為:

式中，ε為蒙皮發(fā)射率，與蒙皮材料有關(guān)，本文取ε=0．8;c1與c2為輻射常數(shù)。

蒙皮的輻射亮度不僅與其自身輻射有關(guān)，還與反射的環(huán)境輻射有關(guān)，本文主要考慮蒙皮面元對太陽輻射的反射，其反射輻射亮度為:

式中，θD為面元法向與視線方向的夾角;其余各參數(shù)與式(16)一致;F(λ1，λ2，Tsun)為黑體輻射函數(shù)，其計算式為:

式中，M(λ，Tsun)為波長為λ，溫度為Tsun的黑體光譜輻射出射度，可由普朗克公式計算得出，Tsun=5900 K。

蒙皮面元的紅外輻射亮度等于自身輻射亮度與反射太陽輻射亮度之和，即:

3 排氣系統(tǒng)紅外輻射計算

3．1 尾焰的溫度和組分分布計算

本文采用半經(jīng)驗法計算尾焰流場，將尾焰分成初始段與主體段兩部分，其中初始段可分為核心區(qū)與混合區(qū)，如圖2所示。設(shè)核心區(qū)長度為LC，外邊界錐角為α，尾噴管半徑為R0，初始段內(nèi)邊界半徑為rc，尾焰外邊界半徑為rm，飛機飛行速度為Vplane，尾焰流速為Vnozzle，當(dāng)?shù)芈曀贋閂sonic。

圖2 尾焰流場模型

各幾何參數(shù)的計算公式如下:核心區(qū)長度:

邊界錐角(α0=8．3°):

本文將大氣簡化成由CO2、H2O、N2、O2四種氣體構(gòu)成，根據(jù)渦流的傳播規(guī)律，對于坐標(biāo)為(x，r)的點，其溫度計算式為:

其中，T0為發(fā)動機噴口的氣體溫度;Ta為飛機所在環(huán)境的大氣溫度;Tw為x軸上溫度，其計算式為:

氣體的重量分?jǐn)?shù)計算式為:

式中，gi為尾焰氣體第i種組分的重量成分;g0i為發(fā)動機噴口處第i種氣體的重量成分;gai為飛機所在環(huán)境的大氣中第i種氣體的重量成分。

3．2 排氣系統(tǒng)紅外輻射亮度計算

由于紅外成像設(shè)備在某一時刻拍攝到的尾焰紅外圖像實際上是一幅二維圖像，圖像中每個像素點的灰度值取決于尾焰經(jīng)由包絡(luò)面發(fā)射向成像設(shè)備的紅外輻射亮度大小，因此只需對尾焰包絡(luò)面劃分網(wǎng)格，如圖3所示，網(wǎng)格的數(shù)量取決于計算的精度。

圖3 尾焰包絡(luò)面網(wǎng)格

計算前首先需要對每個結(jié)點進行可見性判斷，判斷方法如下:設(shè)結(jié)點P的坐標(biāo)為(xp，yp，zp)，外法向為 珗np=(npx，npy，npz);探測器所在位置的坐標(biāo)為(xD，yD，zD)，法向為 珗nD=(nDx，nDy，nDz)，探測器接收角為φ，點P到探測器的光束方向為珗ndir，如圖4所示。若 珗np與 珗ndir的內(nèi)積(珗np，珗ndir)＜0且 珗ndir與 珗nD的夾角滿足:

則點P可見，反之點P不可見。

圖4 可見性判斷

若點P可見，求取珗ndir所在直線與包絡(luò)面的另一個交點Q，交點Q有兩種情況，一是Q點在尾焰擴散邊界面上，二是Q點在尾噴口截面上。第一種情況下，Q點的輻射亮度為周圍環(huán)境大氣的光譜輻射亮度。第二種情況下，將尾噴管看作灰色輻射體，Q點的輻射亮度按照下式計算:

式中，εQ為尾噴管熱空腔的發(fā)射率;TQ為熱空腔溫度。計算得到LQ后，將QP路徑劃分為n段，如圖5所示。近似認(rèn)為每一小段路徑Rk，k+1中氣體的溫度和組分相同，對QP路徑應(yīng)用Malkmus統(tǒng)計窄譜帶模型與C－G近似［10］，得到尾焰通過點P的紅外輻射亮度大小。

圖5 QP路徑分段

4 探測器成像仿真

4．1 灰度值量化

分別計算出機體和尾焰的紅外輻射亮度后，對所得的結(jié)果進行比較，得出輻射亮度的最大值Lmax與最小值Lmin，設(shè)任意一點i的輻射亮度為Li，i點的灰度可設(shè)置為:

Gmax與Gmin分別為灰度值上限與下限，取值范圍在［0，255］。

4．2 噪聲模型

噪聲的灰度值變化近似符合高斯分布，噪聲仿真方法如下:

假設(shè)成像系統(tǒng)的感光元件由m×n個像素組成，需要首先生成m×n個均值為μ、方差為σ的隨機高斯數(shù)，作為噪聲灰度值矩陣，然后將圖像分為m×n個單元，在每個單元原有灰度值的基礎(chǔ)上疊加噪聲灰度值矩陣，即可得到疊加噪聲的紅外圖像。若噪聲類型是不隨幀改變的空間噪聲，在生成下一幀圖像時，不重新生成噪聲灰度值矩陣。若噪聲類型是隨機三維噪聲，則在生成下一幀圖像時，重新生成噪聲灰度值矩陣。

4．3 三維圖像投影

真實成像系統(tǒng)所成的圖像是二維圖像，因此需要對三維圖像進行投影，本文采取透視投影法對三維圖像進行投影，投影時定義圖所示的平截頭體ABCD－A'B'C'D'，將三維圖像置于平截頭體中，然后將平截頭體投影到探測器所在位置。圖中w/h為圖像寬高比，θ1與θ2為探測器的視場角。

圖6 透視投影示意圖

5 仿真算例及分析

本文模擬了F16飛機在不同條件下的紅外圖像，飛行條件設(shè)置如下:飛行高度3 km，飛行速度1．2 Ma，太陽位于目標(biāo)體軸系的y軸正方向，發(fā)動機尾噴口氣體流速600 m/s，探測器與目標(biāo)之間的距離為900 m，像素數(shù)為640×480，探測器視場角為3°×2．25°。圖7(a)與圖7(b)分別為探測器位于目標(biāo)正上方與目標(biāo)側(cè)后方時所成的8～12μm波段的紅外圖像;圖7(c)與圖7(d)分別為探測器位于目標(biāo)正上方與目標(biāo)側(cè)后方時所成的3～5μm波段的紅外圖像;圖7(e)與圖7(f)分別為探測器位于目標(biāo)正上方，但不考慮太陽輻射時所成的3～5μm與8～12μm波段的紅外圖像。

保持探測器的視場角為3°×2．25°不變，改變探測器與目標(biāo)之間的距離，圖8(a)為探測器與目標(biāo)距離1800 m時所成的8～12μm波段的紅外圖像;圖8(b)為探測器與目標(biāo)距離3200 m時所成的8～12μm波段的紅外圖像，

按照4．2節(jié)中的噪聲模型生成640×480像素的隨機噪聲，圖9(a)與圖9(b)分別為疊加噪聲前后飛機的側(cè)向紅外圖像。

分別對比圖7(a)與圖7(c)、圖7(b)與圖7(d)可知，在8～12μm波段，目標(biāo)輪廓較明顯，機體輻射占主要部分，而在3～5μm波段，機體輻射幾乎不可見，尾焰輻射占主要部分。

分別對比圖7(a)與圖7(b)、圖7(c)與圖7(d)可知，探測器與目標(biāo)之間的距離不變時，探測器角度不同，目標(biāo)紅外圖像所占像素點的數(shù)量及目標(biāo)輪廓存在很大差別。

分別對比圖7(c)與圖7(e)、圖7(a)與圖7(f)可知，太陽輻射在3～5μm波段對于機體輻射的影響較大，而在8～12μm波段對于機體輻射的幾乎沒有影響。而對比圖7(c)與圖7(d)，圖7(d)中的機體輻射亮度相對于圖7(c)中較低，說明機體反射太陽輻射的能量與機體角度相關(guān)。

對比圖8(a)與圖8(b)可知，視場角不變的情況下，探測器距目標(biāo)越遠(yuǎn)，目標(biāo)圖像所占的像素點越少，目標(biāo)的輪廓特征及灰度變化越不明顯。

對比圖9(a)與圖9(b)可知，疊加噪聲后，目標(biāo)與背景的對比變?nèi)?，不利于探測器對目標(biāo)的識別。

圖7 探測器與目標(biāo)距離400 m時的紅外圖像

圖8 改變探測器與飛機距離時的紅外圖像

圖9 疊加噪聲前后圖像對比

6 結(jié)語

本文研究了飛機紅外圖像的仿真方法，然后以F16飛機為例，生成了飛機的紅外圖像，并對不同條件下的紅外圖像進行了對比分析。在計算蒙皮表面溫度場時，采用求解熱平衡方程法，并基于傳熱學(xué)原理提出了蒙皮氣動加熱的計算方法;計算排氣系統(tǒng)紅外輻射時，采用基于包絡(luò)面的紅外輻射計算方法，結(jié)合非均勻熱氣體的計算模型，有效減小了計算量。與傳統(tǒng)方法相比，該方法能在準(zhǔn)確性的基礎(chǔ)上以較快的速度生成飛機的紅外圖像，對于分析飛機的紅外圖像及進行紅外成像制導(dǎo)導(dǎo)彈的仿真特征具有較強的應(yīng)用價值。

本文只研究了飛機目標(biāo)本身的紅外圖像，并未對背景的紅外圖像以及紅外輻射的大氣傳輸過程進行研究，在下一步的研究中，應(yīng)將飛機目標(biāo)本身的紅外圖像、背景的紅外圖像與紅外輻射的大氣傳輸過程相結(jié)合，從而模擬出更加真實的紅外場景。

［1］ GE Wei，CAO Dongjie，HAO Hongxu．Application of IR control and guidance technology in precise attack weapons［J］．Acta Armamentarll，2010，36(S2):35－38．(in Chinese)葛煒，曹東杰，郝宏旭．紅外制導(dǎo)技術(shù)在精確打擊武器中的應(yīng)用［J］．兵工學(xué)報，2010，36(S2):35－38．

［2］ ZHANG Ke，LI Ning，LIU Fumei，et al．Method of aircraft skin infrared radiation image generating［J］．Laser＆Infrared，2011，41(3):272－277．(in Chinese)張可，黎寧，劉福美，等．一種飛行器蒙皮紅外輻射圖像生成方法［J］．激光與紅外，2011，41(3):272－277．

［3］ WANG Fei，HE Jing，WANG Xin－sai，et al．Simulation model of IR imaging for the aeroplane［J］．Infrared and Laser Engineering，2007，36(3):352－356．(in Chinese)王飛，賀菁，王新賽，等．空中飛機紅外成像仿真模型研究［J］．紅外與激光工程，2007，36(3):352－356．

［4］ ZHANG Ke．Research on technology of creating aircraft’s skin infrared images［D］．Nanjing:Nanjing University of Aeronautics and Astronautics，2009．(in Chinese)張可．飛行器蒙皮紅外輻射圖像生成技術(shù)研究［D］．南京:南京航空航天大學(xué)，2010．

［5］ MEI Fei，JIANG Yong，CHEN Shiguo，et al．Infrared imaging prediction model for aero－engine exhaust plume［J］．Laser＆Infrared，2012，42(8):909－913．(in Chinese)梅飛，江勇，陳世國，等．一種新的三維實時紅外尾焰仿真方法［J］．紅外技術(shù)，2012，42(8):909－913．

［6］ YU Yang，TANG Xinyi，LIU Peng，et al．A new way of real－time 3D simulation of infrared plume［J］．Infrared Technology，2009，31(10):577－580．(in Chinese)于洋，湯心溢，劉鵬，等．一種航空發(fā)動機噴流紅外成像仿 真 模 型［J］．激 光 與 紅 外，2009，31(10):577－580．

［7］ 姜貴慶，劉連元．高速氣流傳熱與燒蝕熱防護［M］．北京:國防工業(yè)出版社，2003．

［8］ REN Qingmei，YANG Zhibin，CHENG Zhu，et al．Develop－ment of the platform for analysis coupling aero heating and structural temperature field［J］．Structure＆Environment Engineering，2009，36(5):33－38．(in Chinese)任青梅，楊志斌，成竹，等．氣動加熱與結(jié)構(gòu)溫度場耦合分析平臺研發(fā)技術(shù)［J］．強度與環(huán)境，2009，36(5):33－38．

［9］ MAO Xia，HU Haiyong，HUANG Kang．Calculation method for airplane IR radiation and atmospheric transmit－ttance［J］．Journal of Beijing University of Aeronautics and Astronautics，2009，35(10):1228－1231．(in Chinese)毛峽，胡海勇，黃康．飛機紅外輻射及大氣透過率計算方法［J］．北京航空航天大學(xué)學(xué)報，2009，35(10):1228－1231．

［10］LI Jianxun，TONG Zhongxiang，WANG Chaozhe，et al．Calculation and simulation on infrared radiation of hot jet from engine［J］．Spectroscopy，2013，33(1):7－13．(in Chinese)李建勛，童中翔，王超哲，等．發(fā)動機熱噴流紅外輻射計算與仿真［J］．光譜學(xué)與光譜分析，2013，33(1):7－13．