基于卷積曲面的動態實時星圖模擬

閆勁云， 劉慧， 趙偉強， 江潔

(1. 中國計量科學研究院, 北京 100029； 2. 北京航空航天大學儀器科學與光電工程學院, 北京 100083)

星敏感器是航天器中最重要的姿態敏感設備，在設計、調試和測試等環節中，需要提供輸入激勵以完善和評估其各項功能和性能，星模擬器是實現這一任務的最佳途徑[1-2]。星模擬器包括單星模擬器、靜態多星模擬器、動態多星模擬器[3-4]。動態多星模擬器可模擬全天區任意指向的星圖，然后由光學系統將星圖輸出為適合星敏感器接收的星光信號[5]。

隨著星敏感器的動態性能和數據更新率等各種指標的不斷提高，對多星模擬器的要求也隨之提高[6-7]。星敏感器的動態性能是指當其隨著載體運動時，星點成像發生拖尾的情況下，保持正常姿態輸出的能力。為了測試動態性能，要求多星模擬器可以模擬拖尾的星點，目前多星模擬器對于動態拖尾星點模擬的算法研究還較少[5,8]。當前星敏感器的數據更新率普遍可達10 Hz，一些星敏感器可達30 Hz[6]，多星模擬器要求動態刷新顯示變化的星空，其刷新率不低于星敏感器的數據更新率，通常大于30 Hz，因此對算法的實時性也提出了較高要求[9]。

多星模擬器算法的核心任務是高精度實時星圖生成，主要包括導航星搜索、星點坐標變換、星點光斑彌散及像素離散化等步驟。導航星搜索是指根據星敏感器的視軸指向在全天區范圍內搜索當前視場內的星點信息。星點坐標變換將搜索到的星點信息變換為模擬星圖中的圖像坐標和灰度信息。星點光斑彌散模擬了鏡頭的彌散作用，將星點變換為星點光斑[10]。像素離散化過程將連續的星點光斑變為離散的數字圖像。通過星表分區[3,5]和重新組織構建星表[11]等方法，可以解決搜索導航星耗時較高的問題，本文對此不再贅述。靜態下的星點坐標變換模型較為成熟，但在動態條件下，星點在模擬星圖中的坐標構成一條軌跡。劉海波[8]和劉朝山[12]等推導了星點在像面的速度與星敏感器角速度的關系，但未給出星點軌跡曲線的表達式。Sun等[7]在推導過程中使用了一階近似，因此僅適用于曝光時間較短的情況。本文在此基礎上，給出了完備的星點軌跡模型。對于模擬星點光斑彌散過程，Sun等[7]利用運動模糊圖像中的概念，將動態拖尾星點在頻域中表示為鏡頭彌散函數和運動退化函數的乘積，但運動退化函數在時域中的物理含義仍有待進一步明確。Yan等[6]利用線擴散函數的概念，得出了拖尾星點光斑的解析表達式模型，但僅適用于星點在像面做勻速直線運動的情況。本文針對這些問題，建立了拖尾星點的卷積曲面模型，更準確地描述了動態拖尾星點的成像過程。對于像素離散化過程，由于涉及對光斑曲面的二重積分，普遍存在計算復雜、實時性較差的問題[13]。Yan等[6]的方法僅適用于勻速直線運動的情況，劉朝山等[12]的方法假定星點掃過每個像素的時間一致，因此不適合非勻速運動的情況。本文在卷積曲面模型的基礎上提出了一種快速像素離散算法，適用于星點做任意運動的情況，并且大大加快仿真速度，以滿足動態多星模擬器的實時星圖仿真要求。

1 動態拖尾星點模擬

1.1 星點運動軌跡

當星敏感器隨載體發生角位移時，在曝光時間內，星點在星敏感器像面移動成像，星點中心的運動軌跡為一條曲線段。

設星表中恒星i在天球坐標系下的映像為Vi。在t時刻，此恒星在星敏感器坐標系下的映像為Wi(t)，則Vi和Wi(t)可分別表示為

(1)

(2)

式中：αi和δi分別為赤經和赤緯；(xi(t),yi(t))為星點在星敏感器像面的投影坐標；F為星敏感器的鏡頭焦距。

Wi(t)和Vi是同一個矢量在不同坐標系下的映像，二者關系為

Wi(t)=R(t)Vi

(3)

式中：R(t)為星敏感器的姿態矩陣。

設星敏感器相對于天球坐標系的旋轉角速度為ω=[ωx,ωy,ωz]T，式(3)對時間t求導得[8]

-[ω×]Wi(t)

(4)

(5)

若曝光起始時刻t0星點在像面坐標為(xi0,yi0)，則動態下星點在像面的運動軌跡C為

(6)

1.2 基于卷積曲面的光斑彌散模型

星敏感器為了獲得亞像素級的定位精度，無窮遠處的星光經過星敏感器鏡頭在圖像傳感器像面并非匯聚于單個像素，而是成像為一個彌散光斑。靜態下星點光斑能量分布可表示為[6-7]

f(x,y)=ΦTh(x-xc,y-yc)

(7)

式中：Ф為星光被鏡頭收集并入射到像面的輻射通量；T為曝光時間；(xc,yc)表示星點光斑中心位置；h(x,y)為鏡頭的點擴散函數，通常以二維高斯函數來描述[6]。

(8)

(9)

式中：E0為0等星在地球大氣外表面處的輻照度[1]，E0= 2.96×10-14W/mm2；m為入射星光的星等；D為星敏感器的光學鏡頭孔徑；τ為鏡頭的透射率；ρ為鏡頭的彌散半徑。

動態情況下，星點在像面發生移動，星點光斑的能量分布構成一個二維卷積曲面[14-15]，如圖1所示。圖1中，星點的軌跡曲線為C，r為軌跡曲線 上的點，r=(xc(t),yc(t))，其參數方程見式(6)。

圖1 拖尾星點光斑成像過程Fig.1 Imaging process of smeared star spot

p為像平面內的任意一點，p=(x,y)。h(p)為星敏感器光學成像系統的點擴散函數，h(|p-r|)表示此函數的中心位于點r。記權重函數g(r)表示星點經過軌跡曲線上的弧微元|dr|所需的時間，則有

(10)

在點r附近取弧微元|dr|，星點掃過此微元期間的能量為Φg(r)|dr|。顯然，當星點的速度越快，星光累積的能量越少；反之，當星點的速度越慢，能量累積越多。由于光學鏡頭的彌散作用，此星點能量按照點擴散函數描述的方式彌散到了像面，能量分布可表示為Φg(r)h(|p-r|)|dr|，那么沿曲線積分后可得星點在像面的總能量分布為

Φg(p)?h(p)

(11)

式中：?為卷積運算符。星點能量分布為一卷積曲面，即f(p)為點擴散函數h(p)與權重函數g(r)的卷積。拖尾星點光斑的卷積曲面模型具有清晰的物理含義，但不便于運算，為此進一步導出動態下星點能量分布積分模型為

(12)

式(11)和式(12)所表示的模型精確描述了動態條件下的拖尾星點成像過程。

1.3 卷積曲面的快速像素離散算法

入射到像面的光斑能量，在像元內轉化為電荷信號，再經過讀出和放大，轉換為數字圖像。那么位于第i行和第j列的像素的圖像灰度值為[1,5-8,12-13]

(13)

式中：K為轉換系數，是圖像傳感器的量子效率和轉換增益的乘積。圖像傳感器對光斑能量的作用不僅僅體現為離散采樣，還包括積分作用，此過程稱為像素離散。顯然，像素灰度值為一個三重積分表達式，因此仿真一幅星圖需要對每個星點光斑范圍內的每一個像素計算三重積分，運算的時間開銷較大，不能滿足動態多星模擬器實時性的要求。

由拖尾星點的卷積曲面模型的物理含義出發，提出一種快速算法，以解決運算量較大的問題。具體思路為：將曝光時間T等間隔地分為N段，每段時間間隔ΔT=T/N。相應的，軌跡曲線也被分為了N段。當ΔT足夠短時，地軌跡上的每一段可近似為一點，此時的星點光斑成像近似為靜態下的情況。靜態下對二維高斯函數的二重積分可推導得其原函數，將耗時的數值積分轉化為函數求值。最后將這N點沿軌跡累加，即可得到

(14)

(15)

式中：(xck,yck)為第k時間段時星點光斑中心在軌跡上的坐標；erf(x)為誤差函數。式(15)將復雜的三重積分化簡為函數求值、乘法和加法運算，大大提高了拖尾星點的仿真速度。

2 仿真及分析

仿真計算機的配置如下。CPU：Inter Core i7，內存8 GB，仿真環境：MATLAB 2015b。仿真星敏感器的主要參數為：視場20°×20°，鏡頭直徑27 mm，焦距37.92 mm，像素尺寸5.5 μm × 5.5 μm，圖像分辨率2 048像素×2 048像素。采用SAO星表，選取亮于6.0等星的星作為仿真星表。

2.1 精度對比

本節首先對比3種星點仿真算法的精度，再研究對軌跡劃分到何種精度，本文算法的精度可以達到滿意結果。式(13)表示的三重積分法是未優化的原始算法，速度慢，但精度高，得到的星點光斑可認為是真值。本文算法的精度依賴于對積分劃分的細度，隨著N增大，精度逐漸逼近真值。第3種算法為運動模糊圖像卷積法，拖尾星點由原始星點圖像與運動模糊核函數圖像卷積而得[7,10]。定義星點圖像灰度誤差為

ferr=|E(j,i)-I(j,i)|

(16)

進一步定義最大灰度誤差為星點區域內誤差絕對值最大的灰度值；定義灰度RMS誤差為星點區域內每一個像素的誤差灰度值的均方差。

設[ωx,ωy,ωz]=[5(°)/s, 5(°)/s,0]，曝光時間T=20 ms，鏡頭的彌散半徑為0.7像素，星點初始位置為像面中心，則由式(6)可得，星點拖尾長度為14.33像素。采用三重積分法得到的拖尾星點如圖2(a)所示；采用本文算法，將星點軌跡劃分為143段，仿真得到的拖尾星點如圖2(b)所示；采用圖像卷積法得到拖尾星點如圖2(c)所示；后兩者的星點圖像灰度誤差分別如圖2(d)、(e)所示。為便于比較，圖中星點灰度值在最大值處進行了歸一化。

對比圖2(a)、(b)可知，本文仿真結果與真值較為一致。由圖2(d)可知，僅起始點和終點附近存在一定誤差，星點拖尾中間部位的誤差可忽略，其灰度RMS誤差為0.005 4。由圖2(a)、(c)可知，采用圖像卷積法得到的拖尾星點，其星點能量彌散程度較大，圖2(e)所示的星點誤差也較大，其灰度RMS誤差為0.102 0。

由式(14)可知，當軌跡劃分越細，則仿真精度越高。本文算法在不同劃分間隔下，星點灰度誤差如表1所示。

由表1可知，隨著軌跡劃分間隔增加，星點灰度誤差逐漸增大。由于成像過程中存在的背景噪聲及其他噪聲，星點圖像灰度存在起伏，若仿真星點最大誤差小于噪聲引起的灰度起伏，則仿真誤差對最終結果可忽略。若以最大灰度誤差為限，軌跡劃分間隔為0.05像素時，仿真星點精度滿足要求；若以RMS誤差為限，軌跡劃分為0.10像素即可滿足要求。

圖2 拖尾星點光斑仿真精度對比Fig.2 Simulation accuracy comparison of smeared star spots

劃分間隔/像素0.020.050.100.20最大灰度誤差0.00530.01290.02600.0525灰度RMS誤差0.00110.00270.00540.0109

2.2 仿真速度對比

令多星模擬器產生1 000個隨機視軸指向，分別利用三重積分法、運動模糊圖像卷積法和本文算法仿真生成拖尾星圖，統計每幅仿真星圖的仿真時間，3種算法仿真耗時的直方圖統計如圖3所示。

由圖3可知，三重積分法仿真一幅星圖所需時間為幾秒，運動模糊圖像卷積法所需時間為幾百毫秒，本文算法耗時為幾十毫秒，3種算法仿真一幅星圖的平均耗時分別為3 030、309和22 ms。顯然，本文算法的仿真速度大大優于前兩種算法。

圖3(c)所示的個別星圖的仿真耗時大于33 ms，不滿足多星模擬器刷新率的要求。為此，繪制模擬星圖仿真耗時與星點數量的關系圖(見圖4)，以進一步分析超時原因。

由圖4可知，星圖仿真耗時與星點數量是線性關系。星點數量較多時，仿真耗時也越多。對于單個拖尾星點，平均仿真耗時為0.45 ms，那么為保證單幅星圖仿真時間小于33 ms，只要星圖中 星點個數小于70即可。為此，對于模擬星圖中星點數量大于70顆星的情況，根據星點的亮度信息，篩選其中較亮的70顆星，其余亮度較暗的星予以剔除。這樣的處理對實際使用影響較小，原因在于星敏感器為獲得較高的姿態精度，僅采用較亮的若干顆星進行姿態解算。較暗的星形成的星點光斑圖像的信噪比低，因此對其定位的精度也較低，應予以剔除。增加對多余暗星剔除的步驟，得到的仿真耗時直方圖如圖5所示。

圖3 不同算法仿真耗時直方圖Fig.3 Histogram of simulation time consumption for different algorithm

圖4 星圖仿真耗時與星圖中星點數量的關系Fig.4 Star map simulation time consumption vs. star number of star image

由圖5可知，剔除多余暗星后，本文的快速星點仿真算法可滿足動態多星模擬器的實時星圖仿真的要求。

圖5 剔除暗星后仿真耗時直方圖Fig. 5 Histogram of simulation time consumption after dark star is excluded

3 結 論

針對動態多星模擬器對拖尾星點仿真難以滿足實時性的問題，本文首先建立了動態拖尾星點的卷積曲面成像模型，然后提出了一種基于卷積曲面的快速像素離散算法。得到的仿真星點誤差低于圖像噪聲水平，仿真精度高；同時仿真一幅星圖耗時小于33 ms，比現有算法速度提高14倍，滿足了多星模擬器刷新率30 Hz的要求。