基于改進Delaunay三角剖分的空間目標紅外輻射成像方法

王彩云,趙煥玥,李曉飛,劉純勝,王立國,王佳寧

(1.南京航空航天大學航天學院,江蘇 南京 210001;2.北京電子工程總體研究所,北京 100854)

1 引 言

由于航天技術(shù)的發(fā)展,各軍事強國發(fā)射進入太空的空間目標數(shù)量劇增。空間目標的體積愈加小型化,結(jié)構(gòu)更加復雜化,形狀越來越多樣化,空間目標探測被賦予新的挑戰(zhàn)。紅外探測作為一種重要的光學探測手段,因其對目標背景溫差的高敏感性,探測的高隱蔽性以及可以進行多光譜探測的優(yōu)點成為了理想的空間探測方式[1]。與外場試驗相比較,紅外成像仿真可以獲得特定位置、時間和天氣下的紅外圖像,不僅能減少試驗成本,還能提高試驗效率,具有廣泛的應(yīng)用前景與研究價值。

國內(nèi)在空間目標的紅外成像起步較晚,但是也取得了很大進步。例如文獻[2]～[4]中闡述了空間目標的紅外輻射特性分析建模方法,對其紅外輻射強度進行了數(shù)值求解,并將計算結(jié)果轉(zhuǎn)化成灰度值,最后采用計算機圖形學算法生成紅外圖像,為空間目標紅外探測提供了有力的工具。但是在紅外成像方面,已有的研究大都基于理論數(shù)值計算階段,尚未開發(fā)一套從目標表面溫度分布求解到目標紅外輻射特性分析求解,再轉(zhuǎn)入目標紅外成像的系統(tǒng)仿真軟件,而國外成熟的紅外成像商業(yè)軟件由于版權(quán)限制,價格昂貴,也不是最佳選擇[5]。

所以基于經(jīng)濟因素和可行性的綜合考慮,本文提出了一種基于三角剖分的空間目標紅外成像仿真方法,建立了一套完整的空間目標紅外成像仿真方法。此方法首先根據(jù)空間目標熱平衡方程[6],利用節(jié)點網(wǎng)絡(luò)法求解出目標的各節(jié)點溫度,然后利用STK得出探測器、目標、地球與太陽的位置關(guān)系,隨后由目標節(jié)點溫度求得目標紅外輻射強度,最后利用改進Delaunay三角剖分的方法重建目標得到紅外圖像。以某衛(wèi)星[7]為例進行了成像仿真實驗,成像結(jié)果驗證了本文方法的有效性,為進一步的空間目標的探測與識別提供了基礎(chǔ)。

2 空間目標表面溫度求解模型

2.1 理論計算模型

空間目標的紅外輻射特性主要取決于目標的表面溫度分布。要對空間目標進行紅外成像,首先應(yīng)該建立目標表面溫度求解模型,然后根據(jù)溫度分布求解其紅外輻射特性。

為了簡化計算,忽略空間目標自身工作的發(fā)熱,目標表面溫度的分布主要由太陽直接輻射[8]、地球反射太陽輻射[9]、地球直接輻射和目標表面向宇宙空間的輻射四個因素決定。

(1)太陽直接輻射

太陽的半徑為6.960×105km,可以被看作是溫度為5900 K的黑體。把太陽與地球的距離當作其平均距離(1.4961×108km),地球大氣層外法線方向的太陽總輻照度Es為:

(1)

式中,Eλ是太陽的光譜輻射照度。

另外,因為地球半徑遠遠小于太陽與地球距離,可以將太陽光視為輻射度均勻的平行光,則在目標表面的第i個單元吸收的太陽直接輻射qi1為:

qi1=εsEsSiFs,i

(2)

式中,εs表示材料對太陽輻射的吸收率;Si為表面面積;Fs,i是外表面對于太陽直接輻射的角系數(shù)。

(2)地球反射太陽輻射

除了太陽直接輻射,目標表面還會接收地球反射太陽輻射,被視為地球表面的漫反射,則該表面單元所吸收的地球反射太陽輻射qi2為:

qi2=εrEsρFse,i

(3)

式中,εr為吸收率;ρ為地球反射率;Fse,i為角系數(shù)。

(3)地球直接輻射

地球可以被視為一個溫度280 K的黑體,其空間某處的總輻射照度Ee為:

(4)

式中,Re為 6450 km,是地球及其大氣系統(tǒng)總半徑;h為目標和大氣層的距離,則該表面單元吸收的地球直接輻射qi3為:

qi3=εeEeSiFe,i

(5)

式中,εe為吸收率；Fe,i為角系數(shù)。

(4)目標向宇宙空間的輻射

若目標的溫度不低于熱力學零度,則目標表面一定會向外輻射能量。第i個單元表面向宇宙空間的輻射能qi4為:

qi4=εσTi4(t)Si

(6)

式中,ε為發(fā)射率;σ為斯蒂芬-玻爾茲曼常數(shù);Ti(t) 為該部分在時間t的溫度。

(5)內(nèi)能的增加

空間目標的內(nèi)能會根據(jù)溫度的變化而發(fā)生變化,在任一瞬間,第i個單元內(nèi)能的增加量為:

(7)

式中,c表示表面的比熱容;mi為質(zhì)量。

2.2 熱平衡方程的求解

根據(jù)上文的描述,熱平衡方程如下:

qi1+qi2+qi3-qi4=qi5

(8)

(9)

(2)將高次項變?yōu)橐淮雾?

Ti4(t)=4Ti3(t-Δt)Ti(t)-3Ti4(t-Δt)

(10)

如果給定太陽、目標、地球的位置關(guān)系,就能計算出式(8)的角系數(shù)Fs,i、Fse,i和Fe,i的值,然后由式(8)～(10),可知對目標任意表面單元,只要知道其在t-Δt時刻的溫度值Ti(t-Δt),其在t時刻的溫度值Ti(t)就可以推算得到。所以,若知道了目標表面的初始溫度分布以及目標與太陽、地球的方位,目標的表面溫度分布就能求出。

本文采用節(jié)點網(wǎng)絡(luò)法對空間目標的表面溫度分布進行計算。節(jié)點網(wǎng)絡(luò)法是把空間表面劃分為許多熱力學特性基本一致的單元,把單元作為節(jié)點,對節(jié)點和背景的熱交換進行分析,建立各節(jié)點的瞬時熱平衡方程。在給出的初始條件下,遞推計算這些方程組,便可以計算出目標表面節(jié)點的溫度值和其變化規(guī)律[9]。

3 空間目標的紅外輻射特性計算

空間目標對空間某處的紅外輻射包含兩部分:目標自身紅外輻射和目標反射背景輻射。

(1)目標自身紅外輻射

在空間目標表面溫度分布被計算出以后,可以采用普朗克公式,在其紅外波段范圍的積分得到:

(11)

式中,λ1,λ2為紅外波段范圍的上下限;C1為第一輻射常數(shù),其值為3.742×108W·μm4/m2;T為單元表面溫度;C2為第二輻射常數(shù),其值為1.439×104μm·K。

目標自身紅外輻射強度采用漫輻射源的輻射特性得出:

(12)

式中,ΔA為目標表面單元面積;θi為目標表面單元法線與輻射方向的夾角。

(2)目標反射輻射

目標對背景的反射輻射亮度為:

(13)

式中,ρe為等效反射系數(shù);E為外界對目標輻射量的總和。目標反射強度可采用漫反射源的輻射特性得到:

Ir=∑LrΔAcosθi

(14)

式中,ΔA為目標表面單元面積;θi為目標表面單元面法線與反射方向的夾角。

因此,空間目標的總輻射強度為:

I=Ifa+Ir

(15)

4 空間目標的紅外輻射特性成像

4.1 圖像灰度值計算

前一節(jié)計算出的輻射強度并不是紅外成像仿真需要的最終結(jié)果,需要把經(jīng)過成像系統(tǒng)效應(yīng)作用后的輻射強度轉(zhuǎn)化成灰度值。

現(xiàn)在通常使用的是均勻量化,輻射強度和灰度值之間是一個線性關(guān)系。一幅圖像中像面微元對應(yīng)的最高和最低輻射強度分別為Imin和Imax,設(shè)圖像灰度等級為0～255級,則對應(yīng)于輻射強度為I(j)的像面微元的灰度值為:

(16)

4.2 基于改進Delaunay三角剖分的紅外輻射特性成像方法

前面已經(jīng)計算出紅外圖像中節(jié)點的灰度值,下一步是將這些離散的節(jié)點導出,重建為空間目標,并且選取探測器所在角度進行成像,這一步要采用Delaunay三角剖分[11]方法。

Delaunay三角剖分算法常用的分為兩種:基于投影的三角剖分算法[12],此算法容易出現(xiàn)拼接錯誤并且不利于復雜曲面的重構(gòu)；另一種為基于Delaunay的生長算法,即在三維離散點中構(gòu)造一個種子三角形,隨后根據(jù)Delaunay準則尋找最優(yōu)點擴展形成一個新的三角形,以此三角形為基礎(chǔ)再次向外搜索最優(yōu)點,直到所有的散射點都加入到三角網(wǎng)中構(gòu)成三角形[13]。此算法的優(yōu)點是思想簡單,能對復雜三維曲面進行重構(gòu),但是此方法在搜索滿足條件的第三點時計算復雜,耗費大量時間,且易出現(xiàn)曲面重疊、空洞等現(xiàn)象[14]。所以本文將采用三維柵格法[15]和自適應(yīng)空間外接球方法[16]對其進行改進,提出一種改進Delaunay三角剖分算法。

三維柵格法首先讀取三維點云坐標,求得點云在三個坐標軸X、Y和Z軸上的最大值和最小值,確定能將所有數(shù)據(jù)包含在內(nèi)的最小長方體,然后選區(qū)合適的邊長值,將長方體分割成若干個相同大小的立方體。最后依據(jù)各點的三位坐標值將點存儲在相對應(yīng)的小立方體編號數(shù)組中,完成了點云數(shù)據(jù)的空間劃分。該方法可以根據(jù)點云的疏密程度自由調(diào)整子區(qū)域的大小,在領(lǐng)域點搜索的過程中,可以減小搜索范圍,提高效率。

自適應(yīng)空間外接球方法的思想是找到點集中距離生長邊的中心點最近的點組成新三角形,確定此三角形外接球并找到球內(nèi)部的點云數(shù)據(jù),最后通過迭代找到最優(yōu)擴展點。

綜上,本文的改進Delaunay三角剖分算法的步驟為:

(1)利用上述三維柵格法完成某個點的k鄰域?qū)ふ遥?/p>

(2)構(gòu)造種子三角形,利用種子三角形來生長三角網(wǎng)；

(3)采用自適應(yīng)空間外接球方法尋找三角形生長的最優(yōu)點；

(4)對未使用的點集進行迭代生長。

本文采用這種方法先對空間目標進行重建,將紅外輻射強度轉(zhuǎn)化的灰度值分布于空間目標的表面,隨后采用STK計算出探測器方位,根據(jù)其方位矢量得出二維圖像。

5 空間目標的紅外輻射成像仿真實驗

本文以某衛(wèi)星為例,進行紅外輻射成像仿真實驗,首先進行目標三維建模并計算衛(wèi)星表面溫度分布,隨后利用第二節(jié)中的公式求出目標紅外輻射特性大小,最后采用改進的Delaunay三角剖分成像方法得到目標的紅外圖像。

5.1 實驗參數(shù)

首先表1中列出了探測目標的物性參數(shù)。表2列出了目標軌道與探測軌道參數(shù)。

表1 目標物性參數(shù)表

表2 目標與探測軌道參數(shù)

目標和探測器軌道均為太陽同步軌道,其中ad為軌道半長軸;ed為軌道偏心率;id為軌道傾角;Ωd為升交點赤經(jīng);ωd為近地點輻角;τd為升交點地方時。

5.2 結(jié)果及分析

為了驗證本文計算空間目標紅外輻射強度的模型與成像方法是否合理有效。圖1、圖2首先比較了在3～5μm與8～12μm兩個波段下目標本體面元與太陽帆板正面面元總紅外輻射強度變化,從2009 年 6 月 5 日 15∶30開始,橫坐標為距開始時間的秒數(shù),縱坐標為紅外輻射強度大小。由于紅外輻射強度隨目標周期性運動也呈現(xiàn)周期性變化,圖1展現(xiàn)了紅外輻射強度的兩個周期變化曲線。

圖1 兩種波段下目標本體面元總輻射強度隨時間變化曲線

由圖1(a)與圖1(b)對比可知,8～12 μm波段探測的目標本體面元輻射強度比3～5 μm更強,采用8～12 μm波段更易于探測目標。在某一時間段,目標本體面元輻射強度達到最小保持不變,這時面元處于陰影區(qū)。

圖2 兩種波段下目標太陽帆板面元總輻射強度隨時間變化曲線

對比圖1、圖2,可以看出在光照區(qū),目標太陽帆板正面面元總輻射強度在8～12 μm波段大于目標本體面元的總輻射強度,在3～5 μm波段并不明顯,這是因為太陽輻射能量主要集中在中短波,目標本體反射率較大,所以在3～5 μm波段反射了更多的太陽能量。

圖3為目標紅外輻射成像仿真結(jié)果,仿真從5400 s(17∶00)開始,直到114000 s(18∶40)結(jié)束,每隔1200 s成一次像。

圖3 不同時間目標雙波段成像結(jié)果

由圖3可以看出,同一時間8～12 μm波段輻射強度大于3～5 μm波段。目標的成像角度隨時間不斷變化,在光照區(qū),目標的太陽帆板正面的亮度大于目標本體,滿足理論實際。綜上所述,本文提出的基于改進Delaunay三角剖分的空間目標紅外輻射成像方法有效。

為了比較傳統(tǒng)的Delaunay三角剖分算法與本文改進后的Delaunay三角剖分算法的重建效率,進行了兩種算法對同一目標的重建時間對比實驗,實驗平臺均為Matlab R2017a,且均在安裝 64 位Windows 10 操作系統(tǒng)的同一臺計算機上完成,表3是運行時間比較結(jié)果。

表3 算法運行時間比較

由表3可以看出,本文改進后的算法重建時間比傳統(tǒng)Delaunay三角剖分算法重建時間短,效率更高。

6 結(jié) 論

本文提出了一種基于改進Delaunay三角剖分的空間目標紅外輻射成像方法,建立了一套完整的空間目標紅外輻射成像方法。此方法首先根據(jù)熱平衡方程,利用節(jié)點網(wǎng)絡(luò)法求解出空間目標的各節(jié)點溫度,然后利用STK得出探測器、目標、地球與太陽的位置關(guān)系,隨后由目標節(jié)點溫度求得目標紅外輻射強度,最后利用改進的Delaunay三角剖分的方法重建目標得到紅外圖像。本文以某衛(wèi)星為例進行了成像仿真實驗,并對實驗結(jié)果進行了分析,驗證了本文方法的有效性,此方法對空間目標的紅外探測以及識別具有重要意義。