大角速度條件下星像運動軌跡建模及誤差評估

何貽洋，王宏力，馮磊，*，由四海，陳志侃

（1.火箭軍工程大學 導彈工程學院，西安710025； 2.中國人民解放軍66133部隊，北京100144）

星敏感器作為一種姿態測量設備，具有精度高、自主性好以及無誤差漂移等優點，因而在航空航天領域具有重要的應用價值［1-2］。星敏感器通常捷聯安裝于彈道導彈、衛星等載體上，若存在大姿態角機動或者受到外界環境影響時，由于目前星敏感器的動態性能仍然較差，星敏感器拍攝的星圖會模糊，這可能嚴重影響星點提取成功率和精度，進而降低星圖識別成功率，甚至可能導致星敏感器無法正常定姿［3-5］。

經典的模糊圖像復原方法主要分為以下兩大類：盲復原方法和非盲復原方法。由于星圖的先驗信息可以基于慣導信息等獲取，因而在模糊星圖復原的研究中，非盲復原方法應用更為廣泛［6-10］。模糊核函數是非盲復原方法的基礎，而分析星圖模糊的機理是獲得模糊核函數的第一步。文獻［11-13］將星圖模糊機理近似為成像平面的勻速直線運動，利用Radon變換、倒頻譜分析等方法得到星圖模糊的尺度和方向，從而估計模糊星圖的點擴散函數。文獻［14-15］研究了轉動條件下成像平面上的星像運動軌跡，給出星像運動軌跡近似為一段直線的結論，但對于星像運動軌跡形式分析并不全面。李楠和路小波［16］研究了散焦和運動模糊的復合模型，提出了一種基于自相關的參數估計方法，在一定程度上減輕了混疊現象帶來的干擾。文獻［17］中采用功率譜密度函數的方法分析了載體振動對星敏感器成像的影響，可以較好地模擬星圖的振動模糊效果，這也為星圖模糊的研究提供了新思路。文獻［18］提出了一種基于線段擴散函數的星像能量分布模型，為星圖模糊機理的研究開辟了新思路，但仍是基于星像運動軌跡為直線這一前提的?？梢姴簧傥墨I中都將星像運動軌跡當做勻速直線運動來處理，但實際上星像運動軌跡隨著載體運動形式的不同而各異。故而吳小娟和王新龍［3］對于繞光軸和非光軸轉動時成像平面上的星像運動軌跡進行了建模分析（光軸指星敏感器測量坐標系Osxsyszs的Oszs軸，非光軸指垂直光軸方向的Osxs軸和Osys軸），但對于繞任意方向轉動時的星像運動軌跡形式未做進一步分析。目前對星圖模糊的研究主要聚焦于模糊核函數估計和模糊復原算法，而星像運動對星點提取精度影響是怎樣的？將產生多大影響？這些問題鮮有提及，但實際上這是研究星圖模糊對星敏感器性能影響時應該首先回答的問題。

本文針對大角速度條件下星敏感器星像拖尾的問題，對星像運動軌跡的數學模型進行了研究，并定量評估了星像運動對星點提取精度影響的大小，為提高星敏感器在復雜工作環境下的動態性能提供參考。

1 星像運動軌跡建模

根據文獻［3］的研究結果，在曝光時間內，當繞非光軸轉動時，星像運動軌跡近似為勻速直線運動的一段；當繞光軸轉動時，星像運動軌跡為一段圓弧，而實際中星敏感器轉動可繞任意方向。故本文對繞星敏感器測量坐標系內任意方向轉動時的星像運動軌跡進行分析建模。

繞星敏感器測量坐標系Osxsyszs內任意方向轉動時，可將角速度分解到Osxsyszs3個坐標軸上，結合文獻［3］的相關成果，可知在曝光時間內，星像在成像平面上的運動為勻速直線運動和定軸勻角速度運動的疊加。為了方便描述星像運動軌跡，將星像在成像平面上勻速直線運動速度分量 記 為v0x＞0，v0y＞0，v0x≈fωy／dCCD［3］，v0y≈fωx／dCCD。f為星敏感器焦距；ωx、ωy和ωz分別為繞Osxs、Osys和Oszs轉動角速度；dCCD為星敏感器成像平面像元物理尺寸。不妨設星像在成像平面上的定軸勻角速度運動沿順時針方向，角速度記為ωz，星像初始位置位于成像平面的Oy軸上，星像運動示意圖如圖1所示。

圖1 繞任意軸方向轉動時的星像運動示意圖Fig.1 Schematic diagram of star spotmotion when rotating around any axis

首先給出圓的微分方程，不妨設圓心為（a，b）、半徑為R的圓方程為

對x求導，整理可得圓的微分方程為

由于Δt非常小，則可近似認為在Δt時間內，星點在成像平面上的運動速度不變，則從第i→i＋1時刻，星點的位移量為

由式（4）可得

式（5）實質上為導數的定義，因此有

將式（6）中參數t消去，得到運動軌跡的微分方程為

若繞光軸的角速度接近于零，根據式（7），可知星像運動軌跡的圓心位置位于無窮遠處，那么此時形成的圓弧即為一段直線，這與文獻［3］中繞非光軸轉動時的星像運動軌跡為一段直線的結論相一致。所以，無論轉動軸指向如何，星像在成像平面上的運動軌跡都可認為是一段圓弧。

2 星像運動對星點提取精度影響評估

2.1 星像運動時的能量分布建模

對于靜態條件下的星圖，星像在像平面上的能量分布符合二維高斯分布，數學表示為

式中：E0為曝光時間內星點的總能量；（x0，y0）為星點位置；（x，y）為星圖上的某個像素點；σ為高斯彌散斑半徑，代表了能量彌散程度，σ越大，星像能量越分散，根據二維高斯分布的特點以及3σ準則，星像99%以上的能量集中于距離星點的3σ范圍以內。

當星敏感器處于動態環境下，星像在像平面上會產生拖尾現象，星像的能量不再符合二維高斯分布。值得注意的是，在星敏感器鏡頭參數確定的情況下，在一次曝光時間內星像的總能量可認為不變，但在星像運動軌跡上能量將按照某種規律進行分布。下面對星敏感器處于動態情況下，星像的能量分布進行建模分析。

將星敏感器曝光時間Te等分為N個區間［ti，ti＋1］，i＝0，1，…，N－1，Δt＝ti＋1－ti＝Te／N，t0＝0，tN＝Te，Te內星像的總能量為E0，由于每個時間區間［ti，ti＋1］非常小，因而可以認為在［ti，ti＋1］內星像的能量仍服從二維高斯分布，不妨設在［t0，t1］的星像質心坐標為（x0，y0），那么在［ti，ti＋1］內的星像質心點坐標為（xi，yi），i＝0，1，…，N－1。在［ti，ti＋1］時間內，星像能量分布為

式中：Ei，i＋1＝E0／N。

假設上述高斯分布的彌散斑半徑相等，將各個微小曝光時間內的星像能量分布函數進行疊加，得到星像總的能量分布函數為

為了更直觀展示星圖模糊情況下星像的能量分布情況，在繞星敏感器測量坐標系Osxsyszs的坐標軸轉動的條件下，假設E0＝1，σ＝0.5，得到繞Osxs、Osys和Oszs軸轉動時，模糊星圖上星像能量分布的仿真示意圖如圖2所示。由圖2可知，在（a）、（b）條件下，星像運動軌跡為一段直線；在（c）、（d）條件下，星像運動軌跡為一段圓弧。模糊星像不再是符合二維高斯分布的高斯彌散斑，而是呈帶狀分布。

圖2 模糊星圖上星像能量分布仿真示意圖Fig.2 Schematic diagram of star spot energy distribution simulation on blurred star image

2.2 星像運動對星點提取精度影響的評估分析

基于灰度的算法將星點等同于星像灰度分布區域的灰度極值點，利用灰度分布信息求取其質心的精確坐標，是目前常用的星點提取算法，如質心法、高斯曲面擬合法等［5］。

質心法將星像區域內像元的灰度值作為其坐標的權重，然后計算星像區域的一階矩，得到星點坐標。高斯曲面擬合法要求星像灰度分布區域近似服從二維高斯分布，才能獲得較高星點提取精度，高斯曲面擬合方法顯然不再適用于星像模糊嚴重的情形。因此，本文著重研究由轉動引起的星圖運動模糊現象對質心法提取精度的影響。

質心法的計算窗口取為一矩形區域，這一區域的左上角像元坐標設為（1，1），右下角的像元坐標為（m，n）。不妨設G（x，y）為像元（x，y）處的灰度，與該像元的能量成正比，有G0～E0，Gi～Ei，i＋1；根據2.1節可知，每個微小曝光時間區間對應的星像灰度近似服從二維高斯分布，且星像光斑絕大部分能量集中于距離星點坐標3σ范圍以內。那么質心法的計算窗口應能夠包含所有的（xi－3σ，xi＋3σ）＆（yi－3σ，yi＋3σ）區域。

結合式（9）、式（10）則有式（11）～式（14）成立：

（xi，yi）為［ti，ti＋1］內的星像質心坐標，利用質心法可以得到（xi，yi）的計算公式如下：

將式（17）、式（18）分別代入式（15）、式（16），則有

下面 分 析（x′0，y′0）與 靜 態 條 件 下 星 點 坐 標（x0，y0）之間的誤差與星像運動軌跡長度之間的關系。

星點坐標誤差的模定義如下：

式（20）表示實際提取的星點與理想的星點之間的距離。

下面對繞非光軸和光軸轉動的情況下星圖模糊引起的星點坐標誤差進行評估分析。

1）繞非光軸轉動時星圖模糊引起的星點坐標誤差評估

繞非光軸轉動時，星像在成像平面上的運動近似為勻速直線運動，假設星像在成像平面的Ox及Oy軸上的運動速度分別為vx、vy，則在曝光時間 Te內，星 像 模 糊 的 長 度 為 L1＝

式中：i＝0，1，…，N－1；Δt＝ti＋1－ti。

可得

得到星像運動引起的星點坐標誤差為

2）繞光軸轉動時星圖模糊引起的星點坐標誤差評估

繞光軸轉動時，星像在成像平面上作圓周運動，設繞光軸運動的角速度為ωz；星像逆時針做圓周運動，星像圓周運動半徑記為r，曝光時間Te分割為N個等距的微小時間區間Δt；在Δt內，星像模糊的長度l＝rωzΔt，那么Te內星像運動軌跡長度L2＝（N－1）l（也可根據式（4）迭代計算L2），不妨設初始質心（x0，y0）與像平面Ox軸的夾角為β0，對于第i個Δt對應的質心坐標（xi，yi），i＝0，1，…，N－1可表示為

由于星像在Te內所轉動的角度ωzTe為小量，可做如下近似：cos（iωzΔt）≈1，sin（iωzΔt）≈iωzΔt，將式（24）三角函數展開后簡化為

將式（25）代入式（19）中，則有

得到星圖運動模糊引起的星點坐標誤差為

根據式（27）、式（28）可知，當繞光軸轉動時，利用質心法獲得的星點坐標（x′0，y′0）與（x0，y0）之間的誤差（Δx，Δy）與β0有關，星點坐標誤差的模L′2近似為星像運動軌跡長度L2的一半，不隨β0而改變，但與（x0，y0）與主點之間的距離有關，距離主 點 越 遠，L2越 大，星 點 誤 差 的 模也 就越大。

當繞任意方向轉動時，星像在成像平面上的軌跡仍然為一段圓弧，根據上述分析推導可知，星點坐標誤差的模仍可用星像運動軌跡長度L3的二分之一來衡量。在此不再贅述。

在實際應用中，將星敏感器相關參數代入式（4）計算星像運動軌跡長度，根據本文結論，可迅速估計出由于星圖模糊引起的星點坐標誤差，若星圖模糊引起的星點坐標誤差較小，不影響后續的星圖識別和定姿，則不需進行模糊星圖復原，否則需對模糊星圖進行復原，以保證星敏感器正常工作。所以本文為評估星敏感器是否需要啟用模糊星圖復原算法提供了參考。

3 仿真實驗驗證

3.1 基本條件設置

星敏感器視場大小設置為9°×9°，焦距為65.76 mm，像元大小為20μm，面陣大小為512 pixel×512 pixel，曝光時間為100ms，星敏感器鏡頭孔徑為0.04m；從SKY2000主星表篩選出星等小于6的導航星，并剔除其中的雙星與變星，剩余4 908 顆導航星組成本文的星表。利用Monte Carlo方法隨機產生星敏感器的一個光軸指向，其在天球坐標系下的赤經、赤緯為（36°，50°）；星圖模擬仿真實驗中加入的高斯白噪聲均值為0，均方差為6（灰度范圍為0～255）。

3.2 實驗分析

在3.1節的仿真條件下，靜態、無噪聲情況下獲得一幅理想的星圖如圖3所示；靜態、有噪聲情況下獲得一幅星圖如圖4所示。選取星圖中3個具有代表性的星像進行比較，記為star1、star2和star3。star1位于中心區域，靠近主點；star2位于邊緣區域；star3介于star1和star2之間的區域，這3個星像可反映成像平面上的不同區域受星圖模糊影響的情況。

圖3 靜態無噪聲條件的理想星圖Fig.3 Ideal star image under static noise-free conditions

圖4 靜態有噪聲條件的星圖Fig.4 Star image under static and noisy conditions

設置以下4組角速度大小和方向，每組條件下仿真5次：①繞星敏感器測量坐標系的Osxs軸轉動，角速度為1，2，3，4，5（°）／s；②繞星敏感器測量坐標系的Osys軸轉動，角速度為－1，－2，－3，－4，－5（°）／s；③繞星敏感器測量坐標系的Oszs軸 轉 動，角 速 度 為8，10，12，14，16（°）／s；④繞星敏感器測量坐標系的三軸轉動，繞Osxs軸轉動角速度為1，2，3，4，5（°）／s，繞Osys軸轉動角速度為－1，－2，－3，－4，－5（°）／s，繞Oszs軸轉動角速度為8，10，12，14，16（°）／s。其余仿真條件不變，得到star1、star2和star3的運動模糊情況，分別如圖5～圖8所示。

圖5 繞O s x s軸轉動時star1、star2和star3的運動軌跡Fig.5 Motion trajectory of star1，star2 and star3 when rotating around O s x s axis

圖6 繞O s y s軸轉動時star1、star2和star3的運動軌跡Fig.6 Motion trajectory of star1，star2 and star3 when rotating around O s y s axis

圖7 繞O s z s軸轉動時star1、star2和star3的運動軌跡Fig.7 Motion trajectory of star1，star2 and star3 when rotating around O s z s axis

圖8 繞O s x s、O s y s和O s z s軸轉動時star1、star2和star3的運動軌跡Fig.8 Motion trajectory of star1，star2 and star3 when rotating around O s x s，O s y s and O s z s axis

為了驗證本文關于星像運動對星點提取影響評估研究所得結論的正確性，進行了仿真實驗驗證。設置了以下3種轉動的仿真條件：

1）繞Osxs軸角速度ωx＝5sin（πt／2）；

2）繞Oszs軸角速度ωz＝16sin（πt）；

3）繞Osxs軸角速度ωx＝5sin（πt／2），繞Osys軸角速度ωy＝－6sin（πt／3），繞Oszs軸轉動角速度ωz＝16sin（πt）。

3種仿真條件中的角速度變化如圖9所示。

以曝光時間0.1 s作為采樣周期，仿真時間為6 s，得到如下仿真結果。

1）第1組：仿真條件1）下star1、star2和star3的星點坐標誤差結果和星點坐標誤差的模與星像運動軌跡長度二分之一的對比，如圖10、圖11所示。

2）第2組：仿真條件2）下star1、star2和star3的星點坐標誤差結果和星點坐標誤差的模與星像運動軌跡長度二分之一的對比，如圖12、圖13所示。

3）第3組：仿真條件3）下star1、star2和star3的星點坐標誤差結果和星點坐標誤差的模與星像運動軌跡長度二分之一的對比，如圖14、圖15所示。

圖9 仿真條件下的角速度變化曲線Fig.9 Angular velocity curves under simulation conditions

圖10 仿真條件1）的星點坐標誤差Fig.10 Star spot coordinate error under simulation condition 1）

圖11 仿真條件1）的星點坐標誤差的模與星像運動軌跡長度二分之一的對比Fig.11 Comparison ofmodulus of star spot coordinate error and one half of star spotmotion trajectory length under simulation condition 1）

圖12 仿真條件2）的星點坐標誤差Fig.12 Star spot coordinate error under simulation condition 2）

圖13 仿真條件2）的星點坐標誤差的模與星像運動軌跡長度二分之一的對比Fig.13 Comparison ofmodulus of star spot coordinate error and one half of star spotmotion trajectory length under simulation condition 2）

圖14 仿真條件3）的星點坐標誤差Fig.14 Star spot coordinate error under simulation condition 3）

圖15 仿真條件3）的星點坐標誤差的模與星像運動軌跡長度二分之一的對比Fig.15 Comparison ofmodulus of star spot coordinate error and one half of star spotmotion trajectory length under simulation condition 3）

繞星敏感器測量坐標系Osxs軸轉動時的星點坐標誤差結果如圖10所示，此時star1、star2和star3的x軸坐標誤差幾乎都為零，3個星像的y軸坐標誤差曲線的變化規律一致，可見繞非光軸轉動對星圖不同區域星像的影響大小是一樣的，與星像在星圖上的位置無關；繞星敏感器測量坐標系Oszs軸轉動時的星點提取誤差結果如圖12所示，star1、star2和star3坐標誤差曲線的幅值和相位各不相同，這取決于星點與主點的相對位置，star1距離主點最近，其坐標誤差曲線的幅值最小，而star2距離主點最遠，其坐標誤差曲線的幅值最大；繞星敏感器測量坐標系三軸轉動時的星點提取誤差結果如圖14所示，3個星點的x軸坐標誤差變化周期與繞Osys軸角速度ωy的變化周期接近，3個星點的y軸坐標誤差變化周期與繞Osxs軸角速度ωx的變化周期接近，可見在本文仿真條件下繞Oszs軸角速度ωz帶來的影響遠遠小于繞Osxs、Osys軸角速度的影響。圖11、圖13和圖15是3種仿真條件下，star1、star2和star3坐標誤差的模與其運動軌跡長度二分之一的對比，可見無論角速度大小和方向如何，星像運動對星點提取精度的影響都可近似用星點運動軌跡長度二分之一來衡量評估。這也驗證了本文結論的正確性和有效性。

4 結 論

本文建立了星敏感器繞任意方向轉動時在成像平面上星像運動軌跡的數學模型，得到了此時的星像運動軌跡仍是一段圓弧的結論；在此基礎上，基于微元思想建立了星像運動時的能量分布模型，推導出質心法的星點提取誤差的模近似為星點運動軌跡長度二分之一，這可為定量評估星像運動對星點提取影響的大小提供重要的理論參考。本文結論也說明了星圖模糊確實會影響星敏感器的工作性能，開展模糊星圖復原的研究是十分必要的。