倪 濤,任 鑫,劉建軍,李春來,嚴 韋,陳王麗,張曉霞,高興燁
(1. 中國科學院大學,北京 100049;2. 中國科學院國家天文臺,北京 100101;3. 中國科學院月球與深空探測重點實驗室,北京 100101)
隨著人類探索太空能力的不斷增強,深空探測領域已成為航天領域的熱點[1],其中,以月球為目標的探測任務達到了將近一半的比例。月表測繪制圖是月球探測的重要任務之一,主要通過月球探測器攜帶的科學載荷獲取探測數據,然后利用攝影測量技術進行處理和計算,得到月表正射影像、月表三維地形圖等。我國近期發射的月球探測器,包括嫦娥三號、嫦娥四號和嫦娥五號都攜帶了兩臺全景相機,負責探測器著陸區影像數據采集的任務。全景相機獲取的影像數據用于制作月表著陸區厘米級精度的地形數據,該成果能夠為深入研究探測點附近準確的地質構造提供精確的形貌數據,為其它載荷(如礦物光譜儀)數據的綜合研究提供基礎數據,同時也是月球車探測點規劃的重要參考數據之一,因此,對于后續的研究意義重大。兩臺全景相機固定于巡視器的轉臺上,工作時采用雙目立體成像,相機間的相對位置和姿態關系(以下簡稱相對外方位元素)保持不變,是后續攝影測量平差處理的重要約束條件之一。
目前主流的商業攝影測量軟件,包括smart 3D,photoscan等在稀疏匹配時都是以單幅影像為單位進行光束法平差計算,成像模型基于攝影測量學的共線方程[2]:
(1)
(2)
其中,x0,y0,f為相機的內方位元素;X0,Y0,Z0為相機的投影中心;R為旋轉矩陣。
該方法未加入相機間的相對約束關系,每幅影像都是自由地旋轉平移以達到光線束的最佳交會,因此解算得到的左右相機的相對外方位元素不是嚴格固定的,由左右相機構成的立體像對解算的三維坐標必然存在一定的誤差。
目前已經有大量關于行星表面全景相機雙目立體影像的處理案例,文[3]總結了美國 “勇氣號” 和 “機遇號” 的定位與制圖方法,提出平差方法主要來源于光束法平差及其改進算法,文[4]利用嫦娥三號全景相機的雙目立體影像實現了月表地形的三維重建,算法部分主要利用基于攝影測量學共線方程的前方交會算法。經過大量的文獻調研和總結發現,目前行星表面雙目立體成像數據處理大多采用基于傳統共線方程的平差算法,很少考慮左右相機相對外方位元素固定的約束條件。
在航空攝影領域,多視影像成像同樣存在相機間相對外方位元素固定的情況,與本文研究的月表全景相機成像方式有相似之處,且已經有大量的研究成果。文[5]針對多鏡頭的MOMS-02相機提出了一種改進的成像模型。文[6]提出了傾斜航空攝影中,考慮同一測站下視與斜視相機間的約束關系,采用6個偏心參數描述兩者之間的變換。文[7]提出多視影像的傾斜航空影像區域網可以采用3種平差方式,包括無約束的聯合定向方式、附加相對約束的定向方式以及直接定向方式。
本文在調研了大量航空攝影多視影像平差方法的文獻后,總結出一種附加相對約束關系的月表全景相機平差算法,主要是將同一測站左右相機的相對外方位元素保持不變的條件加入光束法平差模型中,將每個攝站的多個相機整體作為一個剛體單位進行平差計算,該方法的平差模型更加嚴密,網型更加穩定。經過驗證,該平差模型能夠保證同一測站兩臺全景相機的相對外方位元素保持不變,并且減少了平差模型中的未知數,能夠為我國探月任務中月表全景相機的攝影測量處理提供技術支持。
本次實驗數據包括嫦娥三號巡視器全景相機環拍數據和嫦娥五號全景相機外場科學驗證試驗數據,全景相機參數見表1。兩組數據的成像方式均為全景相機雙目立體成像,影像數據獲取過程中保持左右相機相對外方位元素不變。

表1 全景相機參數Table 1 Parameters of panoramic cameras
嫦娥五號全景相機科學驗證試驗的目的在于評價相機影像質量、驗證全景相機幾何參數和安裝參數的測量方案。試驗分為內場和外場兩部分,其中外場試驗于2015年11月29日開展,外場全景如圖1。全景相機在兩種不同姿態條件下完成了兩組圖像數據的獲取,圖像數據覆蓋了試驗點附近180°水平范圍和-90°~0°俯仰范圍,總共獲取了195對影像數據,數據總量約3 GB。本文選取其中一組試驗數據,共包括93對影像。
嫦娥三號巡視器搭載了全景相機的科學載荷,主要用于獲取巡視區月表的三維光學影像,在著陸器與巡視器實現兩器月面分離后,巡視器圍繞著陸器一圈獲取了月表影像。全景相機在第1個月晝期間,在兩器月面互拍點A、D上以彩色模式對著陸器成像51對,在N0106點上分別以俯仰角-7°和-19°,偏航角175°~176°的探測模式對巡視器周圍月面進行360°環拍,分別獲取了112對全色圖像和56對彩色圖像。第2個月晝期間,全景相機以彩色成像模式在N0203和N0205兩個探測點對月面進行環拍,共獲取環拍數據112對。本文選取第1月晝在E點和S3點的全色圖像環拍數據,包括56對影像。圖2為嫦娥三號在N0106點獲取的部分影像數據。
為方便推導,首先說明攝影測量學中幾個重要的坐標系。
“工作坊”英文為“workshop”,最早源于德國包豪斯學院現代建筑設計奠基人之一格拉皮烏斯的“工廠學徒制”教育理念[6]。“工作坊”是一種開放的交流方式,其構成要素包括:主題內容;學習方法;活動序列和社會環境[7]。“工作坊”的運行包括四個環節:工作坊創設-分組實訓-成果交流-實訓評價[8],具體形式包括案例分析、角色扮演、集體分享、團體討論、頭腦風暴、教師點評等[6]。

圖1 嫦娥五號全景相機科學驗證試驗外場區域
Fig.1 Outfield of scientific verification test for panoramic cameras in Chang′E-5

圖2 嫦娥三號全景相機在N0106點獲取的部分影像
Fig.2 Partial images taken by panoramic cameras in Chang′E-3 in point N0106
(1)像平面坐標系:攝影方向與影像平面的交點稱為像主點,以像主點為原點,兩對邊機械框標的連線為x軸和y軸,其中與航線方向一致的連線為x軸,航線方向為正向。
(2)像空間坐標系:用于表示像點在像方空間的過渡坐標系,以攝影中心為原點,主光軸為z軸,正向為攝影的反方向,x,y軸與像平面坐標系的x,y軸平行。
(3)像空間輔助坐標系:以攝影中心為原點,以鉛垂方向為Z軸,航線方向為X軸,可以通過外方位元素構成的旋轉矩陣對像空間坐標系進行旋轉變換得到。
如圖3,o-xy為像平面坐標系,S1-xyz為像空間坐標系,S1-XYZ為像空間輔助坐標系。設左相機的投影中心為S1,S1的物方坐標為(XS1,YS1,ZS1),外方位元素構成的旋轉矩陣為R1,右相機的投影中心為S2,S2的物方坐標為(XS2,YS2,ZS2),外方位元素構成的旋轉矩陣為R2。設右相機的像空間坐標系旋轉至與左相機的像空間坐標系平行時,旋轉矩陣為M。

(3)

圖3 3個坐標系的位置關系
Fig.3 Position relationships of three coordinate systems

圖4 左相機旋轉過程
Fig.4 Rotation process of left camera

圖5 右相機旋轉過程
Fig.5 Rotation process of right camera
設S1到S2的平移量為(Δx, Δy, Δz),該平移量以左相機的像空間坐標系為參照。于是,有以下關系:
(4)
(5)
其中,x0,y0,f為右相機的內方位元素;X,Y,Z和x,y為觀測點的物方坐標和像平面坐標;k為比例系數。
將(3)式、 (4)式代入(5)式,得到:
(6)

(7)
對于每一對相片,旋轉矩陣M和平移量(Δx, Δy, Δz)都保持不變。通過上述式子,每個右相機影像共線方程中的外方位元素可以用對應的左相機外方位元素以及旋轉矩陣M和平移量(Δx, Δy, Δz)轉換得到。實際的平差計算中,對于左相機的誤差模型采用原始的共線方程,對于右相機的誤差模型采用本文推導的方程。
首先對影像進行特征點提取和同名像點匹配。試驗采用的尺度不變特征變換(Scale-Invariant Feature Transform, SIFT)算子[8]是用于圖像處理領域的一種描述,對旋轉、尺度縮放、亮度變化保持不變,對視角變化、仿射變換、噪聲也具有一定的穩定性。一般來說,尺度不變特征變換算子的計算結果存在一定數量的誤匹配點,本文對匹配結果使用隨機抽樣一致性算法進一步篩選,最終嫦娥三號巡視器環拍數據選取其中54對影像,共提取42 085個匹配點,嫦娥五號科學驗證試驗數據選取其中75對影像,共提取141 498個匹配點。
由于光束法平差屬于非線性優化,需要提供未知數的初始值。在地球平臺,通常可以通過全球定位系統測量和計算接收機載體的運動方位信息[9]。但是在深空探測領域,由于客觀因素的制約,直接獲取的相機外方位元素精度較低,不能作為平差計算的初始值。本文利用計算機視覺領域的運動恢復結構問題(Structure From Motion, SFM)的相關理論,根據匹配結果使用五點法計算本質矩陣,再用矩陣分解法分解該矩陣得到相機外方位元素的估算值,最后利用基于共線方程的前方交會方法估算像點對應的物方三維坐標。
相機鏡頭畸變模型采用布朗(brown)模型,只考慮徑向畸變的部分,誤差模型如下:
Δx=-x(k0+k1r2+k2r4),
(8)
Δy=-y(k0+k1r2+k2r4),
(9)

誤差模型中未知數與觀測值為非線性關系,故先對方程進行線性化,采用泰勒級數展開并保留一次項,然后采用經典的高斯-牛頓法(Gauss-Newton)進行迭代求解[10]。
考慮到迭代計算可能存在粗差的影響,本文加入了Huber loss函數,該函數是一種常用的魯棒的回歸損失函數,具體形式如下:
(10)
其中,a為誤差方程算出的殘差值;δ為函數的參數,可以根據具體殘差的大小給定。
利用嫦娥五號數據對δ的取值進行實驗,最終通過對比結果選取δ=1.0,即觀測點反投影殘差大于1個像素時,減弱其對整體誤差的影響。
對試驗數據分別進行了附加相對約束關系和未附加相對約束關系的平差計算,并統計了平差計算的相關數據,結果見表2。另外,通過共線方程計算像點反投影的殘差,并將像點中誤差以影像為單位進行統計,兩種方法得到結果的像點殘差分布規律如圖6。

表2 平差實驗結果Table 2 Results of adjustment tests
從以上的統計結果可以看出,兩種方法得到結果的像點反投影誤差規律非常接近,而在未知數個數方面,附加相對約束關系的平差方法能夠減少部分外方位元素,對于節約計算機內存具有一定意義。
通過平差計算得到的影像外方位元素反求左右相機的相對姿態,轉角系統按照x-y-z軸的順序,3個旋轉角分別為ω,φ,κ,統計了3個轉角的平均值、中誤差、最大誤差等,結果見表3。通過相機的外方位線元素求取左右相機間的距離,并統計了相關的誤差指標,結果見表4。
表3和表4顯示,附加相對約束關系的平差方法能夠確保模型內部的剛性關系,而傳統方法由于是完全自由平差,在這方面產生一定誤差,為了估算該誤差產生的影響,利用瞬時視場角對其進行估算,兩組實驗數據的全景相機焦距約為50 mm,像元大小為7.4 μm,取兩者比值得到瞬時視場角約為0.148 mrad,再通過與旋轉角的誤差值取比值,可以估算得出立體像對相對旋轉角誤差造成的像點誤差大小,最后以像對為單位進行統計,結果如圖7。

圖6 平差結果像點中誤差統計。(a) 嫦娥五號數據處理結果; (b) 嫦娥三號數據處理結果
Fig.6 RMS statistics of points of adjustment results

表3 立體像對相對姿態誤差統計Table 3 Error statistics of relative pose between stereo images

表4 立體像對距離的誤差統計Table 4 Error statistics of distance between stereo images
由圖7可以看出,嫦娥五號科學驗證試驗數據的處理結果中,像點誤差約為1個像素,最大超過3個像素,嫦娥三號巡視器全景相機環拍數據的處理結果中,像點誤差約為1個像素,最大超過2個像素。相對姿態誤差導致的像點誤差達到了像素級別,對于最終解算得到的三維坐標必然會產生一定影響。
最后利用附加相對約束關系的平差結果對影像數據進行了密集匹配,通過前方交會的計算方法得到像點對應的物方三維坐標,并生成影像區域的三維模型。這部分內容借助軟件Smart 3D完成,結果如圖8、圖9。

圖7 相對姿態誤差導致的像點誤差分布情況。(a) 嫦娥五號數據處理結果; (b) 嫦娥三號數據處理結果
Fig.7 Distribution of image point error caused by inaccurate relative pose

圖8 嫦娥五號全景相機科學驗證試驗外場區域三維模型
Fig.8 3D model of outfield for scientific verification test of panoramic cameras in Chang′E-5

圖9 嫦娥三號巡視器月表環拍區域三維模型
Fig.9 3D model of lunar surface around Chang′E-3 rover
本文分析了當前月表雙目立體相機平差處理存在的問題,參考了計算機視覺和航空攝影測量領域的處理經驗后使用一種改進的平差模型,利用嫦娥五號科學驗證試驗數據和嫦娥三號巡視器環拍數據進行了驗證。通過結果比較證明,本文使用的方法能夠確保不同立體像對之間的相對外方位元素保持一致,提高左右相機間的內部吻合精度,得到的平差結果更加可靠。另外,該方法還能減少部分外方位元素未知數,在節約計算機內存方面有一定意義。因此,本文提出的平差方法可以應用于相關的月表影像處理工作,為后續的研究提供支持。