李瀟 胡維多

(中國空間技術研究院,北京100094)(北京航空航天大學宇航學院,北京100191)

1 引言

小行星探測將是人類未來深空探測最復雜的任務之一。光學導航就是利用航天器自身的光學傳感器設備探測目標小行星,然后規(guī)劃和處理探測到的目標天體傳感器數(shù)據(jù),結合其他星載設備獲得的數(shù)據(jù),確定航天器相對目標小行星的相對位置和姿態(tài)信息。國內外學者很早就開展了關于軟著陸小行星的導航算法研究。目前結合圖像處理的導航技術發(fā)展較為成熟的一種方案是:利用機載導航相機對小行星表面拍照,從圖像中提取出可視自然特征點,精確地跟蹤這些特征點,用激光測距儀測量探測器到特征點的距離,并利用EKF濾波處理特征點圖像和距離信息,從而確定探測器與目標天體之間的相對位置與姿態(tài)。文獻[1]基于特征點檢測的GNC方法,給出了探測器軟著陸小行星表面的導航算法。文獻[2-3]以日本的MUSES_C計劃為背景,詳細介紹了探測器逼近并著陸小行星的指令控制方案。文獻[4]則根據(jù)圖像序列匹配得到的特征點對,給出了由本質矩陣求取空間探測器的運動比例參數(shù),實現(xiàn)軟著陸過程的導航。

上述文獻都是在假定特征點跟蹤不受小行星自轉影響下展開研究,但是為了確保安全性,盡管探測器在實際著陸小行星表面的過程中下降速度不會很大(如“隼鳥”號探測器的最終下降速度僅為0.1m/s)[2],而考慮到小行星一般自轉周期較小[5],即小行星自轉很快(如Eros小行星的自轉周期僅有5.27h),因此CCD相機有時很難跟蹤拍攝得到小行星表面地形特征點,從而也就不能得到特征點精確的位置信息。本文為了解決特征點跟蹤受小行星自轉因素的影響,特征點溢出相機視野的問題,提出了一種新的方案。首先建立全景著陸地區(qū)地形特征圖像數(shù)據(jù)庫,圖像數(shù)據(jù)庫中包含了地表特征點的經(jīng)緯度信息、圖像紋理信息以及各特征點之間的位置相關性;然后針對CCD實拍得到的特征點在數(shù)據(jù)庫中進行匹配;最后結合匹配點的位置信息,通過計算機視覺原理測量獲取探測器和小行星之間的相對位置和相對姿態(tài)信息。

2 圖像特征庫的建立

2.1 特征點提取

從機載導航相機拍攝得到的目標小行星表面地形特征圖像中獲取可視自然特征點是后續(xù)推演導航算法的關鍵前提。本文采用了圖像處理技術領域一種較為成熟的特征點采集算法——Forstner算子。在一幅灰度圖像中,類似于邊緣、角點、紋理等都可以被認為是該幅圖像的特征信息,而對于一幅小行星地表特征圖像來說,最容易被提取的特征莫過于隕石坑、大石塊之類的自然特征點,這是由于它們一般實際尺寸較大,形狀較為規(guī)則且亮度變化大更易于從地形背景中分割出來。因此,基于上述考慮本文選取的圖像特征信息就是小行星表面的自然特征點。具體的數(shù)據(jù)庫建立步驟如下:

1)利用Forstner算子提取特征區(qū)域(這里不詳細闡述提取步驟,可參考文獻[6])。

2)計算該特征區(qū)域的7個不變中心距[6]。選擇中心距作為數(shù)據(jù)庫特征之一的理由是:不變中心距具有平移、比例和旋轉不變性,針對后續(xù)圖像匹配時可以很好的克服由于CCD相機的平移、旋轉以及圖像比例變化帶來的誤差影響,有效提高匹配精度。

3)統(tǒng)計各個特征區(qū)域中心點的經(jīng)緯度信息。由于小行星的立體圖是可以得到的[5],根據(jù)立體圖和經(jīng)度、緯度可推斷在任意點的半徑,從而得到特征點間的距離。因此,針對每個特征區(qū)域計算距離其最近的兩個特征區(qū)域的歐氏距離,將這三個特征區(qū)域構成一封閉三角形,一同存入數(shù)據(jù)庫。

圖1為著陸地區(qū)全景地形圖,圖2為所有特征區(qū)域(用白框表示)。圖3為利用Forstner算子提取得到的特征區(qū)域中心點,所有原圖均取自NASA的NEAR計劃。

圖1 著陸地區(qū)表面地形原圖

圖2 特征區(qū)域圖像(白框)

圖3 特征區(qū)域中心點(圓圈顯示)

2.2 建立特征庫

建立特征庫所要保存的信息如下:1)特征區(qū)域中心點的經(jīng)緯度;2)特征區(qū)域圖像的不變中心距;3)特征區(qū)域之間的位置相關性;4)特征區(qū)域在圖像中的行、列位置;5)特征區(qū)域的編號(編號排列順序依照圖3中的標記)。表1為所建立的特征區(qū)域位置信息,取著陸區(qū)域的經(jīng)緯度范圍:經(jīng)度120°～210°和緯度-20°～-40°。表2為特征區(qū)域的7個不變中心距信息(其中 φ7只有比例和平移不變性)。表3為特征區(qū)域位置相關性信息,該表包含了距離每個特征區(qū)域最近和次近的兩點信息,每個特征區(qū)域的跨度為1km。

表1 特征區(qū)域位置信息

表2 特征區(qū)域不變中心距信息

表3 特征區(qū)域位置相關性信息

3 特征區(qū)域匹配

探測器下降過程中,CCD相機拍攝著陸區(qū)域圖像,由于小行星自轉速率快,而探測器自身下降速度較小,因此,在前一時刻拍攝的圖像中提取得到的特征點在下一時刻不能完全保證還出現(xiàn)在相機視野中,即發(fā)生了特征點溢出相機視野的現(xiàn)象,這樣就不能精確跟蹤得到特征點。為解決這個問題,本文在上文中建立了著陸區(qū)域全景特征數(shù)據(jù)庫,在探測器下降過程中,用相同的特征點提取方法提取特征點,計算這些提取得到的特征區(qū)域的不變中心距、特征區(qū)域之間位置相關性等信息,然后在數(shù)據(jù)庫里查找進行匹配得到這些特征區(qū)域的精確位置信息。具體的步驟如下:

1)同樣采用Forstner算子提取下降過程中圖像中的特征區(qū)域。

2)建立一個特征信息向量v=(v1,v2,v3,v4,v5,v6),這六個向量分量由特征區(qū)域的前6個不變中心距組成(去掉φ7是因為考慮到φ7只有比例和平移不變性,不具有旋轉不變性,若代入計算后會有較大誤差)。同樣將數(shù)據(jù)庫中每個特征區(qū)域的前6個不變中心距也構成一特征向量v′=(v1,v2,v3,v4,v5,v6)

3)計算這兩個向量差的2-范數(shù),即‖v′-v‖2,若值小于預先設定的閾值,則可判斷為可能匹配正確,將這特征區(qū)域標定為候選特征區(qū)域。

4)基于最大互信息配準方法計算數(shù)據(jù)庫中特征區(qū)域和對應的候選特征區(qū)域之間的互信息。用A表示數(shù)據(jù)庫中特征區(qū)域圖像,B表示候選特征區(qū)域圖像,則圖像A,B之間的互信息可表示為

Ⅰ(A,B)=H(A)+H(B)-H(A,B)

其中H(A)為 A的信息熵,H(B)為B的信息熵,H(A,B)為 A,B的聯(lián)合信息熵。當Ⅰ(A,B)達到最大時,可判定匹配正確。圖4～6展示了匹配過程和結果。

圖4 下降過程拍攝原圖

圖5 特征區(qū)域(白框)

圖6 正確匹配得到的結果(白框)

4 光學測距與定姿算法

本文引入計算機視覺成像原理,通過三維景物和二維圖像的轉換關系來獲得空間三維自然特征點在探測器本體坐標系下的位置矢量。首先建立坐標系:1)探測器本體坐標系Ot-xyz:該坐標系原點位于探測器質心,z軸指向目標小行星。2)導航相機坐標系Oc-xyz:該坐標系原點位于相機光心,z軸沿著光軸方向指向目標小行星。3)著陸點坐標系Of-xyz:該坐標系原點位于著陸點,x軸指向正東,y軸指向正北,z軸指向當?shù)厮矫娣ㄏ蚍较颉?)圖像物理坐標系Oi-xyz:該坐標系原點位于攝像機光軸和圖像平面的交點。根據(jù)計算機視覺成像測量原理:在Oc-xyz系下像點p(x,y)與空間物點p(x,y,z)的幾何關系如下(相機焦距為f)

將式(1)用矩陣形式表示為

本文采用以下的測距與定姿算法[7]:

設定Pi,Pj是Of-xyz系中的兩個特征點,其在相機成像平面上的像點為 pi,pj,其在Oc-xyz系下的對應點為P′i,P′j。很顯然,由于Pi,Pj是三維空間中的兩定點,因此Pi,Pj兩點在Of-xyz系下的距離應等于在Oc-xyz系下的距離,即:

將式(3)經(jīng)過變形可構造特征點間距關系方程,由式(4)表示為

設Pi,Pj在Oc-xyz系下的向量為Ri,Rj,則

式中

在得到測距后,采用QUEST姿態(tài)求解算法[8]得到探測器本體坐標系相對著陸點坐標系的姿態(tài)轉移矩陣和姿態(tài)歐拉角。

5 數(shù)值仿真實驗

為了驗證測距與定姿算法的可行性,建立如下的仿真系統(tǒng),設定圖1中的點(222,444)為探測器下降著陸點(圖中箭頭所指),CCD導航相機的分辨率為300×300像素,相機焦距為10mm,一個像素的實際大小為0.1mm×0.1mm,光心位于(150,150)處,經(jīng)過特征點提取,在圖像中可以獲得5個特征區(qū)域,經(jīng)過與圖像數(shù)據(jù)庫(表1～3)的匹配得到序號為4、6、9、12、13的5個特征區(qū)域。建立著陸點坐標系Of-xyz,著陸點即位于坐標原點,各坐標軸如前所述。查找數(shù)據(jù)庫得到這5個特征區(qū)域中心點的位置坐標分別如表4所示。

表4 Of-xyz系下特征區(qū)域中心點坐標值

將原圖像在圖像平面內分別旋轉0°、5°、8°、10°、15°,用來模擬探測器在下降過程中相對于著陸表面,滾動角處于上述角度時所拍攝得到的5幅不同的圖像,得到像點坐標如表5所示。

表5 Oc-xyz系下像點坐標值mm

經(jīng)過第4節(jié)的測距與定姿算法計算得到測距值和姿態(tài)角結果如表6所示。

表6 測距值和姿態(tài)角值

結合二維運動估計相關理論,本文提出的基于地形特征庫的測距與定姿算法可以很好的解決因小行星自轉速度快導致特征點溢出相機視野的問題。表3建立了特征點位置相關性信息,針對每個特征點都可以在表中查找得到與之距離最近的兩個特征點的位置坐標及位置相關向量。探測器在著陸過程中,部分特征點可能會溢出導航相機視野,此時只需要跟蹤與這些溢出特征點位置靠近且仍停留在相機視野里的另一組特征點,通過提取這些特征點,在數(shù)據(jù)庫中查找得到與之距離靠近的溢出特征點的空間三維位置信息,隨后結合二維運動估計相關理論,依靠這組未溢出視野特征點的像素坐標變化情況,建立像點二維運動模型,從而估計得到溢出特征點的像點坐標,這樣就得到了溢出特征點的空間三維坐標和像點坐標,最后運用上述的測距與定姿算法進行求解。

如圖7、圖8是探測器在平移運動過程中拍攝得到的前后兩幀圖像,對圖7經(jīng)過特征提取,經(jīng)過與數(shù)據(jù)庫中匹配得到序號為4、6、9、12、13的特征點(圖中白框表示),用同樣的操作在圖8中得到序號為4、6、13的特征點(圖中白框表示)。經(jīng)過在數(shù)據(jù)庫中查找得知溢出特征點序號為9、12。

圖7 前一幀圖像

圖8 后一幀圖像

根據(jù)前后兩幀圖像中序號為4、6、13特征點的圖像坐標,計算這三個點橫、縱坐標平移量的均值作為圖像行、列方向的平移量,從而可估計得到序號為9、12的溢出特征點在圖8中的圖像坐標,如表7所示,測距值如表8所示。

表7 前后兩幀特征點圖像坐標

表8 前后兩幀測距值

經(jīng)過大量計算結果表明,姿態(tài)角在0°～30°范圍之內,誤差可控制在正負0.5°之內,超過30°,誤差將增大,對測距和定姿結果有一定影響?？紤]到探測器在下降過程中,為了保證著陸安全,發(fā)動機噴管要保證垂直于著陸平面,即探測器只有垂直速度,沒有水平速度,因此探測器相對著陸平面的姿態(tài)角變化不會太大,在較小的姿態(tài)角變化范圍之內,本文的測距與定姿算法可以保證較高的精度。

6 結束語

本文針對小行星自轉速度快,表面地形特征點較難跟蹤,容易溢出導航相機視野的情況,設計了基于建立圖像數(shù)據(jù)庫的測距、定姿算法。本算法首先采用Forstner算子提取著陸點全景圖像的特征區(qū)域,將特征區(qū)域圖像的位置信息、位置相關性、圖像紋理信息、熵信息存入圖像數(shù)據(jù)庫,然后利用計算機視覺原理針對匹配得到的特征區(qū)域在CCD圖像中的像素位置,運用2D到3D的轉化關系得到探測器和各個特征區(qū)域的直線距離以及探測器相對著陸平面的姿態(tài)角信息,為后續(xù)的自主實時導航提供了較高精度的位置和姿態(tài)信息。最后通過仿真實驗驗證了這一方案的可行性。

本文僅針對著陸地區(qū)建立了地形特征庫,適用于探測器最終逼近階段,為了將本文的測距與定姿方法應用到探測器軟著陸小行星的全過程,需要建立目標小行星全景范圍的地形特征庫,同時考慮到小行星形狀大小并不規(guī)則,且受光照條件影響較大,尤其是在光照角度較小時圖像較暗,很難精確跟蹤特征點,因此要建立全景的地形數(shù)據(jù)庫將是一項巨大且艱巨的工程。為此,本文作者正在做這方面研究的一些嘗試和探索。

[1]MISU T,HASHIMOTO T.Optical Guidance for Autonomous Landing of Spacecraf t[R].IEEE Transactions on Aerospace and Electronic System,1999,35(2):459-472.

[2]KAWAGUCHI J,HASHIMOTO T,KUBOTA T.Autonomous Optical Guidance and Navigation Strategy around a Small Body[J].Journal of Guidance,Control and Dynamics.1997,20(5):1010-1220.

[3]KUBOTA T,HASHIMOTO T.An Autonomous Navigation and Guidance System for MUSES-C Asteroid Landing[J].Acta Astronautica.2003,52:125-131.

[4]翟冬麗.基于圖像序列的小行星軟著陸導航方法研究[D].哈爾濱:哈爾濱工業(yè)大學,2006.

[5]胡維多,等.飛行器近小行星軌道動力學的特點及研究意義[J].天文學進展,2009,27(2):152-166.

[6]楊帆,等.數(shù)字圖像處理與分析[M].北京:北京航空航天大學出版社,2007.

[7]朱仁璋.航天器交會對接技術[M].長沙:國防出版社,2007.

[8]SHUSTER M D.Three-Axis Attitude Determination from Vector Observations[J].J.Guidance and Control,1981,4(1):70-77.