侯晴宇 鞏晉南 樊志鵬 王一惠

(哈爾濱工業(yè)大學空間光學工程研究中心,哈爾濱 150001)

1 引 言

在天基光學空間目標態(tài)勢感知領域,目前的研究側重于遠距離點目標檢測、目標軌跡的預測與確認.而在遠距離成像無目標結構、紋理信息時,利用時序點信號對空間目標的在軌工作狀態(tài)及基本物理屬性參量進行反演具有重要的意義,可有效支撐空間目標狀態(tài)判別以及在軌維護的相關決策.

天基光學空間目標態(tài)勢感知的前提為目標的特性認知.目前,在軌目標光學成像特性的研究分為實測以及仿真兩類方法.由于在軌測量的成本較高,得到的測量數據相對有限,因此在該領域多采用地面材料特性測量與成像仿真結合的方式實現在軌目標成像特性研究.近些年來,對于目標材料表面可見光波段雙向反射分布函數(BRDF)測量以及光度信號、圖像仿真的研究較為深入,給出了多種空間目標材料BRDF的實驗室測量方法[1,2]以及建模方法[3],并建立了基于BRDF的空間目標圖像仿真方法[4?7].

空間目標的可見光波段光度信號由目標反射太陽光得到,它是目標外形結構、尺寸、姿態(tài)指向、表面材料等物理屬性參量的函數[8,9].由于耦合量多,基于光度信號的參量估計的難度較大,一般情況下需要給定先驗知識.在目標表面材料反射特性模型已知的條件下,國內外開展了一些基于地基光學探測系統(tǒng)對于空間目標姿態(tài)、外形的估計研究.Calef等[10]在假設目標姿態(tài)和指向已知的前提下,利用時序光度和熱輻射數據反演目標的三維外形.Hinks等[11]分析了姿態(tài)變化與光度信號變化之間的關系,以及利用時序光度信號推導姿態(tài)變化的可行性.Wetterer等[12,13]從濾波角度給出了姿態(tài)估計的流程.文獻[14]跟蹤了國內外基于光度數據反演目標特征信息的最新動態(tài),分析了目前基于地基平臺進行光度數據姿態(tài)反演的主要方法及特點.文獻[15]研究了國外幾種基于點目標測量信息的空間目標識別方法的基本原理和實現途徑.文獻[16]基于時序光度信號匹配的思想開展了衛(wèi)星形狀反演的研究.文獻[17]基于地面實驗測量獲得了目標光譜特性,反向提取目標的材料、大小和狀態(tài)等特征,具有一定的參考價值.

在上述文獻中,針對空間目標的相關狀態(tài)及特性反演研究多以目標材料的光學特性為已知輸入.但是,在實際應用中,非合作目標的本體材料組成未知,參數反演難度較大.因此,本文以含太陽翼的三軸穩(wěn)定衛(wèi)星為研究對象,基于雙面模型及BRDF的多級融合表征模型,建立了融合帆板在軌動態(tài)特性的光學特性宏觀表征模型,以在軌測量的時序光度信號為輸入,建立優(yōu)化方法估計模型參數,完成了對宏觀表征模型的重構;時序光度信號仿真以及參數反演仿真驗證了本文所提方法的有效性.

2 可見光波段時序光度建模

2.1 總體方案

為了體現空間目標在軌過程中帆板與本體運動特性上的差異,本文在文獻[4—7]的基礎上進行了相應完善,融入了帆板的對日指向計算模型.在軌空間目標的可見光波段時序光度特性建模流程如圖1所示.首先,根據空間目標的幾何結構模型進行面元劃分,確定本體面元的法向量及位置等相關屬性;然后,根據空間目標的三軸穩(wěn)定姿態(tài),確定該時刻的空間目標照明與觀測矢量.考慮帆板對日指向,確定帆板面元的法向量及位置等相關屬性;再根據照明與觀測矢量以及面元的相對位置信息,確定既能夠被光照又能夠被觀測的有效面元.結合有效面元的BRDF模型及太陽的輻射特性計算面元的反射特性,結合觀測星的位置確定光瞳處接收到的總的光度特性.

圖1 在軌空間目標時序光度特性建模流程Fig.1. Modeling process of sequence photometric characteristics of on-orbit space object.

2.2 時序光度信號模型

空間目標成像特性建模需要定義下列坐標系:衛(wèi)星本體坐標系Rb(OXY Z)、帆板坐標系Rp(OXpYpZp)、相機坐標系R′(OX′Y′Z′)及探測器坐標系R′′(OMN),如圖2所示.

各個坐標系之間的矢量變換滿足:

式中RY(β)為繞Y(Yp)軸轉動β角時的旋轉矩陣;ξ,ψ,ζ分別為相機坐標系坐標方向在本體坐標系下的單位矢量;T0是相機坐標系原點在本體坐標系下的坐標;M,N分別為探測器的行列數;a為探測器像元尺寸;f為光學系統(tǒng)焦距.

圖2 坐標系之間的關系示意圖Fig.2.Relationship diagram among coordinate systems.

衛(wèi)星在軌運行過程中,反射能量絕大部分為太陽輻射,因此本文研究的衛(wèi)星反射特性主要是指衛(wèi)星反射直接入射太陽光的特性.設本體坐標系下空間目標第j個面中第k個面元在入瞳方向的輻射強度為反射光譜強度時波長和時間的函數,表示為

式中Esun(λ)為太陽輻照度;Aj,k為第j個面中第k個面元面積;面元矢量nj,k,oj,k為本體坐標系下第j個面中第k個面元-觀測星入瞳中心的單位矢量;s為本體坐標系下空間目標-太陽的單位矢量;fj(λ,nj,k,oj,k,s)表示第k個面元為第j種材料時的BRDF,〈·〉為點積運算.其中,oj,k,s,fj(λ,nj,k,oj,k,s)都為時變量.

探測系統(tǒng)入瞳接收到的面元輻照度為

光學系統(tǒng)入瞳接收到的目標輻通量為

式中Aaperture為光學系統(tǒng)的入瞳直徑.

光學系統(tǒng)入瞳接收到的目標輻通量為

考慮到遠距離點目標成像過程,面元-觀測星光學系統(tǒng)入瞳中心的單位矢量基本一致,記為o.(7)式簡化為

可以看出,空間目標光度信號與面元在本體坐標系下的法向量nj,k以及照明矢量s和觀測矢量o有關.

2.3 帆板的對日指向計算模型

對于三軸穩(wěn)定衛(wèi)星,本體面元的相關屬性如位置、法向量在本體坐標系下保持不變.但對于太陽帆板面元,由于保持對日指向運動,帆板與本體之間存在實時相對轉動,因此需要單獨計算太陽帆板面元在本體坐標系下的法向量.

由于太陽距地球距離遠,入射太陽光可近似為平行光.對于三軸穩(wěn)定衛(wèi)星,一般具有較大的太陽能帆板,從衛(wèi)星姿態(tài)控制的角度考慮,帆板通常為單自由度旋轉,即帆板繞衛(wèi)星本體OY軸旋轉,則帆板坐標系下照明矢量sp表示為

式中,s為本體坐標系中的太陽光入射方向矢量,s=[sx,sy,sz].

為確保太陽翼獲得最大能量,需保證帆板法線與太陽光入射矢量之間的夾角最小.為保證該條件,如圖2,在帆板坐標系內,sp的Xp分量為零,即

因此,本體坐標系下帆板的法向量表示為

3 光度模型的在軌重構方法

3.1 空間目標的等效雙面模型

對于正常工作狀態(tài)下的對地定向三軸穩(wěn)定衛(wèi)星,對其進行同軌道面或近軌道面觀測時,目標本體和觀測系統(tǒng)的相對姿態(tài)基本不發(fā)生變化,使得本體被觀測到的表面時序不變,假設該表面為平面(如圖3所示),并且本體和帆板不發(fā)生相互遮擋,可將觀測模型簡化為雙面模型,表示為

式中

其中Φp(t,λ),Φb(t,λ)分別為帆板和本體的光通量;Ap為帆板的面積;Ab為本體中被觀測表面的面積.

圖3 對三軸穩(wěn)定衛(wèi)星的在軌觀測過程Fig.3.On orbit observation process of three axis stabilized satellite.

在軌觀測過程中,光度信號Φ(t,λ),Esun(λ),Aaperture由光學成像系統(tǒng)獲得,為已知量;T0以及ζ,o,nb與觀測幾何相關,可基于空間定軌、定位方法獲得,可認為是已知量;np可以根據帆板對日指向的約束計算得到.Apfp(λ,np,o,s)和Abfb(λ,nb,o,s),即本體和帆板的面積-BRDF乘積,為未知量.

3.2 材料表面BRDF的多級融合表征模型

對于太陽帆板,其表面材料為電池片,電池片的材料組分相差不大,其BRDF可以通過地面測量得到[3].對于空間目標本體,由于衛(wèi)星功能差異,可能存在多種包覆不同表面材料的構件,因此其BRDF未知,這里采用多級融合表征模型,該模型表示為

式中bm為與β′m(nb,o,s)對應的反射率,b={bm},bm需要滿足0≤bm≤1,bm≤1.表示為

m=0時,表示漫反射;m逐漸增大,表示鏡面反射的增強.在這里,提出應用Cook-Torrance BRDF作為參數化模型,即

式中α為n和b之間的夾角,b=(s+o)/2.

因為該模型考慮了面元微觀分布特性,并且對漫反射以及不同角分布的鏡面反射的表征能力較強,其形狀由μm調節(jié),μm表征了鏡反峰的寬度,可將其設定為0.04≤μm≤0.5,μm按照一定的間隔進行取值.此時,(14)式表示為

式中Q(t,λ)=〈np·s〉〈np·o〉·Apfpp(λ,np,o,s)+,則未知量為面積與反射率的乘積fpp及{Abbm}m=0,1,···,M.

3.3 光學特性模型的在軌重構方法

將以上Q(t,λ)式轉化為矩陣形式有

式中Q為T×1維向量,

βp為T×1向量,其中向量元素=p();βb為T×M矩陣,其中矩陣元素[為,

將(19)式簡化為

(20)式的誤差向量可以用矩陣向量形式表示:r(B′)=Q?βB′,由于B′中元素的物理意義為反射率面積的乘積,所以對B′中的元素做如下約束:b′≥ 0.

特征矩陣β為帆板和本體的BRDF信息矩陣,當T＞M時,即特征矩陣β行數大于列數時,只需適當選取B′使誤差r(B′)在2范數意義下最小,即‖r(B′)‖=‖Q?βB′‖最小,即可對B′進行求解,即完成了空間目標復雜材料表面的物性參數反演以及宏觀表征模型的重構.之所以稱之為宏觀表征模型,是因為該模型利用BRDF的基函數重構得到,是微觀BRDF模型的宏觀光度體現.

4 仿真驗證

4.1 時序光度信號的仿真

4.1.1 仿真輸入

衛(wèi)星本體尺寸為1.2 m×1.2 m×1.2 m,帆板尺寸為6.8 m×1.2 m×0.1 m;上圓柱的直徑為0.1 m,高0.3 m;下圓柱的直徑為1 m,高0.5 m;與下圓柱相連接的半球高0.3 m.衛(wèi)星的結構如圖4所示.

圖4 衛(wèi)星結構圖 (a)衛(wèi)星結構簡圖;(b)面元劃分后的結果Fig.4.Structure of satellites:(a)Satellite structure;(b)result of surface division.

仿真計算過程中,衛(wèi)星本體表面包覆黃色熱控材料,帆板正面貼滿太陽能電池片,帆板背面涂有機黑漆,本體上下兩個組件涂有機白漆.對于有機黑漆(反射率0.04)和有機白漆(反射率0.9),其反射特性遵從朗伯漫反射定律.黃色熱控材料和電池片具有較強的鏡反射特性,利用雙向反射分布函數進行描述.對于黃色熱控材料及太陽電池片,其BRDF的測量及建模方法見文獻[5],采用改進的Sun模型對太陽電池片的BRDF進行描述,能夠較好地擬合測量數據,如圖5所示.圖5給出了入射角為30°,45°,60°時測量得到的雙向反射分布函數與反射角的關系.可以看出,兩種材料表現出了很強的鏡反射特性,反射能量主要集中在鏡面反射方向±10°的范圍.

圖5 (網刊彩色)黃色熱控材料與太陽電池片的BRDF測量與建模結果 (a)入射角30°;(b)入射角45°;(c)入射角60°Fig.5.(color online)BRDF measurement and modeling results of yellow thermal control materials and solar cells:(a)Incident angle of 30°;(b)incident angle of 45°;(c)incident angle of 60°.

觀測衛(wèi)星探測系統(tǒng)參數如表1所列.

選擇兩種軌道參數作為仿真參數,表2表示近軌,表3表示觀測星與目標星同軌.

表1 觀測星探測系統(tǒng)參數Table 1.Parameters of observing satellite detection system.

表2 目標衛(wèi)星與觀測衛(wèi)星的軌道參數1Table 2.Orbital elements 1 of the target satellite and observation satellite.

表3 目標衛(wèi)星與觀測衛(wèi)星的軌道參數2Table 3.Orbital elements 2 of the target satellite and observation satellite.

4.1.2 仿真結果

對于表2的軌道參數進行仿真得到的空間目標時序光度信號如圖6所示,本體和帆板的反射峰值均在38幀,帆板能量占較大比例.

針對表3中的軌道參數仿真得到的目標整體時序光度信號如圖7(a)所示,圖中信號有兩個峰值,分別為衛(wèi)星帆板和本體的反射峰值.帆板的鏡面反射峰值在30幀,本體反射峰值點在42幀;此時,光度信號呈現了多峰特性.目標本體和帆板的時序光度信號如圖7(b)和圖7(c)所示.

4.2 光學特性宏觀模型的重構

4.2.1 軌道參數1對應的重構結果

針對上述衛(wèi)星模型、衛(wèi)星軌道參數1、探測器參數仿真所得到的目標時序光度曲線,進行基于雙面法的多級BRDF反演,結果列于表4.根據反演的結果,分別對本體、帆板、目標整體的光度信號進行重構,結果如圖8所示.

由圖8可見,觀測星與目標星在同一軌道時,時序光度信號呈現單峰特性,重構與實測數據整體擬合效果很好,但帆板、本體重構誤差較大,主要原因為帆板和本體信號發(fā)生耦合,本文方法對混合信號曲線中各組分的特征辨識度低,無法精確分離和提取各個信號,故在同軌觀測情況下,算法誤差較大.

表4 衛(wèi)星軌道1對應的參數反演結果Table 4.Parameters inversion results on orbit 1.

4.2.2 軌道參數2對應的重構結果

針對上述衛(wèi)星模型、衛(wèi)星軌道參數2、探測器參數仿真所得到的目標時序光度曲線,進行基于雙面法的多級BRDF反演,結果列于表5.根據反演的結果,分別對本體、帆板、目標整體的光度信號進行重構,結果如圖9所示.

圖6 軌道1目標時序光度信號 (a)整體時序光度信號;(b)本體時序光度信號;(c)帆板時序光度信號Fig.6.Space object photometric sequence signal on orbit 1:(a)Photometric signal of space object;(b)photometric signal of body;(c)photometric signal of panel.

表5 衛(wèi)星軌道參數2對應的參數反演結果Table 5.Parameters inversion result on orbit 2.

圖7 軌道2目標時序光度信號 (a)整體時序光度信號;(b)本體時序光度信號;(c)帆板時序光度信號Fig.7.Space object photometric sequence signal on orbit 2:(a)Photometric signal of space object;(b)photometric signal of body;(c)photometric signal of panel.

由圖9可見,觀測星和目標星處于相近軌道時,時序光度信號為多峰特性,分別體現了本體和帆板分別占優(yōu)時的特性,對于二者的信號解混起到了決定性的作用,此時針對二者的參數反演及信號重構精度較高.

4.2.3 重構結果的對比分析

重構精度用確定系數(R-square)來表征,確定系數是通過數據的變化來表示兩組數據的相關性,取值范圍在0—1之間,越接近1,表明數據相關性越好,對原始數據的解釋能力越強;確定系數表達式如下:

圖8 (網刊彩色)軌道1對應的時序光度信號重構結果(a)目標整體重構結果;(b)本體重構結果;(c)帆板重構結果Fig.8.(color online)Reconstructing photometric sequence signal for orbit 1:(a)Reconstructing result of space object;(b)reconstructing result of body;(c)reconstructing result of panel.

其中SSE為擬合數據與原始數據對應點的誤差的平方和,SST為原始數據和均值之差的平方和.

經過計算軌道參數1對應的目標衛(wèi)星整體、帆板、本體的R-square分別為0.99913,0.97883,0.86196.軌道參數2對應的目標衛(wèi)星整體、帆板、本體的R-square分別為0.99753,0.99978,0.97304.

對于軌道1,帆板和本體的確定系數均低于軌道2的反演結果,但是目標整體光度信號的確定系數卻高于軌道2的反演結果.說明該反演方法對于目標星和觀測星在同一軌道的情況下,帆板本體信號的分離有一定誤差,但對目標整體光學模型具有很強的描述能力.

圖9 (網刊彩色)軌道2對應的時序光度信號重構結果(a)目標整體重構結果;(b)本體重構結果;(c)帆板重構結果Fig.9.(color online)Reconstructing photometric sequence signal for orbit 2:(a)Reconstructing result of space object;(b)reconstructing result of body;(c)reconstructing result of panel.

對于軌道2,帆板的重構精度最高,幾乎與原始數據重合.重構精度高是因為帆板的主要材料為太陽能電池片,且沒有其他構件,因此光學特性顯著,可以進行精確的反演與重構.本體的重構精度較帆板略低,由圖4衛(wèi)星結構可知,衛(wèi)星本體除了一個立方體構件,還有其他構件.本文反演采用雙面法理論,即將本體、帆板各視為一個面,所以本體的其他構件對本體重構的精度有一定影響,但重構誤差在可允許的范圍內.

綜上,采用基于雙面法的多級BRDF融合模型對空間目標光學特性宏觀模型具有很強的描述能力,當目標信號有多峰特性的條件下,能夠進行精準的衛(wèi)星帆板、本體信號分離;在無多峰特性時,能夠較為精確地重現目標整體的光度信號.

5 結 論

針對在軌非合作空間目標的光學特性反演問題,基于在軌觀測過程的雙面假設以及BRDF多級表征,提出了基于時序光度信號分析的光學特性模型參數反演及模型重構方法.利用在軌空間目標光度模型進行了時序光度信號的仿真,結果表明提出的宏觀表征模型的在軌重構方法能夠實現近軌條件下針對本體、帆板97%以上的信號重構精度.該方法可為天基平臺空間目標光學態(tài)勢感知提供一種解決途徑,并可為空間目標姿態(tài)、形狀反演提供技術支撐.

[1]Yuan Y,Sun C M,Zhang X B 2010Acta Phys.Sin.59 2097(in Chinese)[袁艷,孫成明,張修寶 2010物理學報59 2097]

[2]Wang H Y,Zhang W,Dong A T 2013Measurement46 3654

[3]Hou X Y,Zhi X Y,Zhang H L,Zhang W 2014Proc.SPIE9299 929914

[4]Yuan Y,Sun C M,Zhang X B,Zhao H J,Wang Q 2010Acta Opt.Sin.30 2748(in Chinese)[袁艷,孫成明,張修寶,趙慧潔,王潛2010光學學報30 2748]

[5]Wang H Y,Zhang W 2012J.Mod.Opt.59 547

[6]Wang H Y,Zhang W,Dong A T 2012Appl.Opt.51 7810

[7]Wang H Y,Zhang W,Wang F G 2012Sci.China:Technol.Sci.55 982

[8]Linares R,Shoemaker M,Walker A,Mehta P M,Palmer D M,Thompson D C,Koller J,Crassidis J L 2013Proceedings of the Advanced Maui Optical and Space Surveillance Technologies ConferenceMaui,Hawaii,United States,September 10–13,2013 p1889

[9]Coughlin J 2014Proceedings of the Advanced Maui Optical and Space Surveillance Technologies ConferenceMaui,Hawaii,United States,September 9–12,2014

[10]Calef B,Africano J,Birge B,Hall D,Kervin P 2006Proc.SPIE6307 6307E

[11]Hinks J C,Linares R,Crassidis J L 2013AIAA Guidance,Navigation,and Control(GNC)ConferenceBoston,MA,United States,August 19–22,2013 p1

[12]Wetterer C J,Jah M 2009J.Guid.Control.Dynam.32 1648

[13]Holzinger M J,Alfriend K T,Wetterer C J,Luu K K,Sabol C,Hamada K 2014J.Guid.Control.Dynam.37 921

[14]Gou X R,Du X P,Liu H 2016Laser Opt.Progress.53 1(in Chinese)[茍瑞新,杜小平,劉浩 2016激光與光電子學進展53 1]

[15]Zhuang X X,Ruan N J,Zhao S S 2016Infrared Laser Eng.45 201(in Chinese)[莊緒霞,阮寧娟,趙思思 2016紅外與激光工程45 201]

[16]Du X P,Liu H,Chen H,Gou X R,Du M J 2016Acta Opt.Sin.36 251(in Chinese)[杜小平,劉浩,陳杭,茍瑞新,杜明江2016光學學報36 251]

[17]Sun C M,Zhao F,Yuan Y 2015Acta Phys.Sin.64 034202(in Chinese)[孫成明,趙飛,袁艷 2015物理學報64 034202]