附加相對(duì)約束關(guān)系的月表全景相機(jī)平差算法研究*

倪 濤，任 鑫，劉建軍，李春來，嚴(yán) 韋，陳王麗，張曉霞，高興燁

(1. 中國(guó)科學(xué)院大學(xué)，北京 100049；2. 中國(guó)科學(xué)院國(guó)家天文臺(tái)，北京 100101；3. 中國(guó)科學(xué)院月球與深空探測(cè)重點(diǎn)實(shí)驗(yàn)室，北京 100101)

隨著人類探索太空能力的不斷增強(qiáng)，深空探測(cè)領(lǐng)域已成為航天領(lǐng)域的熱點(diǎn)[1]，其中，以月球?yàn)槟繕?biāo)的探測(cè)任務(wù)達(dá)到了將近一半的比例。月表測(cè)繪制圖是月球探測(cè)的重要任務(wù)之一，主要通過月球探測(cè)器攜帶的科學(xué)載荷獲取探測(cè)數(shù)據(jù)，然后利用攝影測(cè)量技術(shù)進(jìn)行處理和計(jì)算，得到月表正射影像、月表三維地形圖等。我國(guó)近期發(fā)射的月球探測(cè)器，包括嫦娥三號(hào)、嫦娥四號(hào)和嫦娥五號(hào)都攜帶了兩臺(tái)全景相機(jī)，負(fù)責(zé)探測(cè)器著陸區(qū)影像數(shù)據(jù)采集的任務(wù)。全景相機(jī)獲取的影像數(shù)據(jù)用于制作月表著陸區(qū)厘米級(jí)精度的地形數(shù)據(jù)，該成果能夠?yàn)樯钊胙芯刻綔y(cè)點(diǎn)附近準(zhǔn)確的地質(zhì)構(gòu)造提供精確的形貌數(shù)據(jù)，為其它載荷(如礦物光譜儀)數(shù)據(jù)的綜合研究提供基礎(chǔ)數(shù)據(jù)，同時(shí)也是月球車探測(cè)點(diǎn)規(guī)劃的重要參考數(shù)據(jù)之一，因此，對(duì)于后續(xù)的研究意義重大。兩臺(tái)全景相機(jī)固定于巡視器的轉(zhuǎn)臺(tái)上，工作時(shí)采用雙目立體成像，相機(jī)間的相對(duì)位置和姿態(tài)關(guān)系(以下簡(jiǎn)稱相對(duì)外方位元素)保持不變，是后續(xù)攝影測(cè)量平差處理的重要約束條件之一。

目前主流的商業(yè)攝影測(cè)量軟件，包括smart 3D，photoscan等在稀疏匹配時(shí)都是以單幅影像為單位進(jìn)行光束法平差計(jì)算，成像模型基于攝影測(cè)量學(xué)的共線方程[2]：

(1)

(2)

其中，x0，y0，f為相機(jī)的內(nèi)方位元素；X0，Y0，Z0為相機(jī)的投影中心；R為旋轉(zhuǎn)矩陣。

該方法未加入相機(jī)間的相對(duì)約束關(guān)系，每幅影像都是自由地旋轉(zhuǎn)平移以達(dá)到光線束的最佳交會(huì)，因此解算得到的左右相機(jī)的相對(duì)外方位元素不是嚴(yán)格固定的，由左右相機(jī)構(gòu)成的立體像對(duì)解算的三維坐標(biāo)必然存在一定的誤差。

目前已經(jīng)有大量關(guān)于行星表面全景相機(jī)雙目立體影像的處理案例，文[3]總結(jié)了美國(guó) “勇氣號(hào)” 和 “機(jī)遇號(hào)” 的定位與制圖方法，提出平差方法主要來源于光束法平差及其改進(jìn)算法，文[4]利用嫦娥三號(hào)全景相機(jī)的雙目立體影像實(shí)現(xiàn)了月表地形的三維重建，算法部分主要利用基于攝影測(cè)量學(xué)共線方程的前方交會(huì)算法。經(jīng)過大量的文獻(xiàn)調(diào)研和總結(jié)發(fā)現(xiàn)，目前行星表面雙目立體成像數(shù)據(jù)處理大多采用基于傳統(tǒng)共線方程的平差算法，很少考慮左右相機(jī)相對(duì)外方位元素固定的約束條件。

在航空攝影領(lǐng)域，多視影像成像同樣存在相機(jī)間相對(duì)外方位元素固定的情況，與本文研究的月表全景相機(jī)成像方式有相似之處，且已經(jīng)有大量的研究成果。文[5]針對(duì)多鏡頭的MOMS-02相機(jī)提出了一種改進(jìn)的成像模型。文[6]提出了傾斜航空攝影中，考慮同一測(cè)站下視與斜視相機(jī)間的約束關(guān)系，采用6個(gè)偏心參數(shù)描述兩者之間的變換。文[7]提出多視影像的傾斜航空影像區(qū)域網(wǎng)可以采用3種平差方式，包括無約束的聯(lián)合定向方式、附加相對(duì)約束的定向方式以及直接定向方式。

本文在調(diào)研了大量航空攝影多視影像平差方法的文獻(xiàn)后，總結(jié)出一種附加相對(duì)約束關(guān)系的月表全景相機(jī)平差算法，主要是將同一測(cè)站左右相機(jī)的相對(duì)外方位元素保持不變的條件加入光束法平差模型中，將每個(gè)攝站的多個(gè)相機(jī)整體作為一個(gè)剛體單位進(jìn)行平差計(jì)算，該方法的平差模型更加嚴(yán)密，網(wǎng)型更加穩(wěn)定。經(jīng)過驗(yàn)證，該平差模型能夠保證同一測(cè)站兩臺(tái)全景相機(jī)的相對(duì)外方位元素保持不變，并且減少了平差模型中的未知數(shù)，能夠?yàn)槲覈?guó)探月任務(wù)中月表全景相機(jī)的攝影測(cè)量處理提供技術(shù)支持。

1 月表全景相機(jī)及實(shí)驗(yàn)數(shù)據(jù)簡(jiǎn)介

本次實(shí)驗(yàn)數(shù)據(jù)包括嫦娥三號(hào)巡視器全景相機(jī)環(huán)拍數(shù)據(jù)和嫦娥五號(hào)全景相機(jī)外場(chǎng)科學(xué)驗(yàn)證試驗(yàn)數(shù)據(jù)，全景相機(jī)參數(shù)見表1。兩組數(shù)據(jù)的成像方式均為全景相機(jī)雙目立體成像，影像數(shù)據(jù)獲取過程中保持左右相機(jī)相對(duì)外方位元素不變。

表1 全景相機(jī)參數(shù)Table 1 Parameters of panoramic cameras

嫦娥五號(hào)全景相機(jī)科學(xué)驗(yàn)證試驗(yàn)的目的在于評(píng)價(jià)相機(jī)影像質(zhì)量、驗(yàn)證全景相機(jī)幾何參數(shù)和安裝參數(shù)的測(cè)量方案。試驗(yàn)分為內(nèi)場(chǎng)和外場(chǎng)兩部分，其中外場(chǎng)試驗(yàn)于2015年11月29日開展，外場(chǎng)全景如圖1。全景相機(jī)在兩種不同姿態(tài)條件下完成了兩組圖像數(shù)據(jù)的獲取，圖像數(shù)據(jù)覆蓋了試驗(yàn)點(diǎn)附近180°水平范圍和-90°～0°俯仰范圍，總共獲取了195對(duì)影像數(shù)據(jù)，數(shù)據(jù)總量約3 GB。本文選取其中一組試驗(yàn)數(shù)據(jù)，共包括93對(duì)影像。

嫦娥三號(hào)巡視器搭載了全景相機(jī)的科學(xué)載荷，主要用于獲取巡視區(qū)月表的三維光學(xué)影像，在著陸器與巡視器實(shí)現(xiàn)兩器月面分離后，巡視器圍繞著陸器一圈獲取了月表影像。全景相機(jī)在第1個(gè)月晝期間，在兩器月面互拍點(diǎn)A、D上以彩色模式對(duì)著陸器成像51對(duì)，在N0106點(diǎn)上分別以俯仰角-7°和-19°，偏航角175°～176°的探測(cè)模式對(duì)巡視器周圍月面進(jìn)行360°環(huán)拍，分別獲取了112對(duì)全色圖像和56對(duì)彩色圖像。第2個(gè)月晝期間，全景相機(jī)以彩色成像模式在N0203和N0205兩個(gè)探測(cè)點(diǎn)對(duì)月面進(jìn)行環(huán)拍，共獲取環(huán)拍數(shù)據(jù)112對(duì)。本文選取第1月晝?cè)贓點(diǎn)和S3點(diǎn)的全色圖像環(huán)拍數(shù)據(jù)，包括56對(duì)影像。圖2為嫦娥三號(hào)在N0106點(diǎn)獲取的部分影像數(shù)據(jù)。

2 雙目立體相機(jī)的平差方法

為方便推導(dǎo)，首先說明攝影測(cè)量學(xué)中幾個(gè)重要的坐標(biāo)系。

“工作坊”英文為“workshop”，最早源于德國(guó)包豪斯學(xué)院現(xiàn)代建筑設(shè)計(jì)奠基人之一格拉皮烏斯的“工廠學(xué)徒制”教育理念[6]?！肮ぷ鞣弧笔且环N開放的交流方式，其構(gòu)成要素包括：主題內(nèi)容；學(xué)習(xí)方法；活動(dòng)序列和社會(huì)環(huán)境[7]?！肮ぷ鞣弧钡倪\(yùn)行包括四個(gè)環(huán)節(jié)：工作坊創(chuàng)設(shè)-分組實(shí)訓(xùn)-成果交流-實(shí)訓(xùn)評(píng)價(jià)[8]，具體形式包括案例分析、角色扮演、集體分享、團(tuán)體討論、頭腦風(fēng)暴、教師點(diǎn)評(píng)等[6]。

圖1 嫦娥五號(hào)全景相機(jī)科學(xué)驗(yàn)證試驗(yàn)外場(chǎng)區(qū)域
Fig.1 Outfield of scientific verification test for panoramic cameras in Chang′E-5

圖2 嫦娥三號(hào)全景相機(jī)在N0106點(diǎn)獲取的部分影像
Fig.2 Partial images taken by panoramic cameras in Chang′E-3 in point N0106

(1)像平面坐標(biāo)系：攝影方向與影像平面的交點(diǎn)稱為像主點(diǎn)，以像主點(diǎn)為原點(diǎn)，兩對(duì)邊機(jī)械框標(biāo)的連線為x軸和y軸，其中與航線方向一致的連線為x軸，航線方向?yàn)檎颉?/p>

(2)像空間坐標(biāo)系：用于表示像點(diǎn)在像方空間的過渡坐標(biāo)系，以攝影中心為原點(diǎn)，主光軸為z軸，正向?yàn)閿z影的反方向，x，y軸與像平面坐標(biāo)系的x，y軸平行。

(3)像空間輔助坐標(biāo)系：以攝影中心為原點(diǎn)，以鉛垂方向?yàn)閆軸，航線方向?yàn)閄軸，可以通過外方位元素構(gòu)成的旋轉(zhuǎn)矩陣對(duì)像空間坐標(biāo)系進(jìn)行旋轉(zhuǎn)變換得到。

如圖3，o-xy為像平面坐標(biāo)系，S1-xyz為像空間坐標(biāo)系，S1-XYZ為像空間輔助坐標(biāo)系。設(shè)左相機(jī)的投影中心為S1，S1的物方坐標(biāo)為(XS1,YS1,ZS1)，外方位元素構(gòu)成的旋轉(zhuǎn)矩陣為R1，右相機(jī)的投影中心為S2，S2的物方坐標(biāo)為(XS2,YS2,ZS2)，外方位元素構(gòu)成的旋轉(zhuǎn)矩陣為R2。設(shè)右相機(jī)的像空間坐標(biāo)系旋轉(zhuǎn)至與左相機(jī)的像空間坐標(biāo)系平行時(shí)，旋轉(zhuǎn)矩陣為M。

(3)

圖3 3個(gè)坐標(biāo)系的位置關(guān)系

Fig.3 Position relationships of three coordinate systems

圖4 左相機(jī)旋轉(zhuǎn)過程
Fig.4 Rotation process of left camera

圖5 右相機(jī)旋轉(zhuǎn)過程
Fig.5 Rotation process of right camera

設(shè)S1到S2的平移量為(Δx, Δy, Δz)，該平移量以左相機(jī)的像空間坐標(biāo)系為參照。于是，有以下關(guān)系：

(4)

(5)

其中，x0，y0，f為右相機(jī)的內(nèi)方位元素；X，Y，Z和x，y為觀測(cè)點(diǎn)的物方坐標(biāo)和像平面坐標(biāo)；k為比例系數(shù)。

將(3)式、 (4)式代入(5)式，得到：

(6)

(7)

對(duì)于每一對(duì)相片，旋轉(zhuǎn)矩陣M和平移量(Δx, Δy, Δz)都保持不變。通過上述式子，每個(gè)右相機(jī)影像共線方程中的外方位元素可以用對(duì)應(yīng)的左相機(jī)外方位元素以及旋轉(zhuǎn)矩陣M和平移量(Δx, Δy, Δz)轉(zhuǎn)換得到。實(shí)際的平差計(jì)算中，對(duì)于左相機(jī)的誤差模型采用原始的共線方程，對(duì)于右相機(jī)的誤差模型采用本文推導(dǎo)的方程。

3 實(shí)驗(yàn)驗(yàn)證及分析

首先對(duì)影像進(jìn)行特征點(diǎn)提取和同名像點(diǎn)匹配。試驗(yàn)采用的尺度不變特征變換(Scale-Invariant Feature Transform, SIFT)算子[8]是用于圖像處理領(lǐng)域的一種描述，對(duì)旋轉(zhuǎn)、尺度縮放、亮度變化保持不變，對(duì)視角變化、仿射變換、噪聲也具有一定的穩(wěn)定性。一般來說，尺度不變特征變換算子的計(jì)算結(jié)果存在一定數(shù)量的誤匹配點(diǎn)，本文對(duì)匹配結(jié)果使用隨機(jī)抽樣一致性算法進(jìn)一步篩選，最終嫦娥三號(hào)巡視器環(huán)拍數(shù)據(jù)選取其中54對(duì)影像，共提取42 085個(gè)匹配點(diǎn)，嫦娥五號(hào)科學(xué)驗(yàn)證試驗(yàn)數(shù)據(jù)選取其中75對(duì)影像，共提取141 498個(gè)匹配點(diǎn)。

由于光束法平差屬于非線性優(yōu)化，需要提供未知數(shù)的初始值。在地球平臺(tái)，通?？梢酝ㄟ^全球定位系統(tǒng)測(cè)量和計(jì)算接收機(jī)載體的運(yùn)動(dòng)方位信息[9]。但是在深空探測(cè)領(lǐng)域，由于客觀因素的制約，直接獲取的相機(jī)外方位元素精度較低，不能作為平差計(jì)算的初始值。本文利用計(jì)算機(jī)視覺領(lǐng)域的運(yùn)動(dòng)恢復(fù)結(jié)構(gòu)問題(Structure From Motion, SFM)的相關(guān)理論，根據(jù)匹配結(jié)果使用五點(diǎn)法計(jì)算本質(zhì)矩陣，再用矩陣分解法分解該矩陣得到相機(jī)外方位元素的估算值，最后利用基于共線方程的前方交會(huì)方法估算像點(diǎn)對(duì)應(yīng)的物方三維坐標(biāo)。

相機(jī)鏡頭畸變模型采用布朗(brown)模型，只考慮徑向畸變的部分，誤差模型如下：

Δx=-x(k0+k1r2+k2r4)，

(8)

Δy=-y(k0+k1r2+k2r4)，

(9)

誤差模型中未知數(shù)與觀測(cè)值為非線性關(guān)系，故先對(duì)方程進(jìn)行線性化，采用泰勒級(jí)數(shù)展開并保留一次項(xiàng)，然后采用經(jīng)典的高斯-牛頓法(Gauss-Newton)進(jìn)行迭代求解[10]。

考慮到迭代計(jì)算可能存在粗差的影響，本文加入了Huber loss函數(shù)，該函數(shù)是一種常用的魯棒的回歸損失函數(shù)，具體形式如下：

(10)

其中，a為誤差方程算出的殘差值；δ為函數(shù)的參數(shù)，可以根據(jù)具體殘差的大小給定。

利用嫦娥五號(hào)數(shù)據(jù)對(duì)δ的取值進(jìn)行實(shí)驗(yàn)，最終通過對(duì)比結(jié)果選取δ=1.0，即觀測(cè)點(diǎn)反投影殘差大于1個(gè)像素時(shí)，減弱其對(duì)整體誤差的影響。

對(duì)試驗(yàn)數(shù)據(jù)分別進(jìn)行了附加相對(duì)約束關(guān)系和未附加相對(duì)約束關(guān)系的平差計(jì)算，并統(tǒng)計(jì)了平差計(jì)算的相關(guān)數(shù)據(jù)，結(jié)果見表2。另外，通過共線方程計(jì)算像點(diǎn)反投影的殘差，并將像點(diǎn)中誤差以影像為單位進(jìn)行統(tǒng)計(jì)，兩種方法得到結(jié)果的像點(diǎn)殘差分布規(guī)律如圖6。

表2 平差實(shí)驗(yàn)結(jié)果Table 2 Results of adjustment tests

從以上的統(tǒng)計(jì)結(jié)果可以看出，兩種方法得到結(jié)果的像點(diǎn)反投影誤差規(guī)律非常接近，而在未知數(shù)個(gè)數(shù)方面，附加相對(duì)約束關(guān)系的平差方法能夠減少部分外方位元素，對(duì)于節(jié)約計(jì)算機(jī)內(nèi)存具有一定意義。

通過平差計(jì)算得到的影像外方位元素反求左右相機(jī)的相對(duì)姿態(tài)，轉(zhuǎn)角系統(tǒng)按照x-y-z軸的順序，3個(gè)旋轉(zhuǎn)角分別為ω,φ,κ，統(tǒng)計(jì)了3個(gè)轉(zhuǎn)角的平均值、中誤差、最大誤差等，結(jié)果見表3。通過相機(jī)的外方位線元素求取左右相機(jī)間的距離，并統(tǒng)計(jì)了相關(guān)的誤差指標(biāo)，結(jié)果見表4。

表3和表4顯示，附加相對(duì)約束關(guān)系的平差方法能夠確保模型內(nèi)部的剛性關(guān)系，而傳統(tǒng)方法由于是完全自由平差，在這方面產(chǎn)生一定誤差，為了估算該誤差產(chǎn)生的影響，利用瞬時(shí)視場(chǎng)角對(duì)其進(jìn)行估算，兩組實(shí)驗(yàn)數(shù)據(jù)的全景相機(jī)焦距約為50 mm，像元大小為7.4 μm，取兩者比值得到瞬時(shí)視場(chǎng)角約為0.148 mrad，再通過與旋轉(zhuǎn)角的誤差值取比值，可以估算得出立體像對(duì)相對(duì)旋轉(zhuǎn)角誤差造成的像點(diǎn)誤差大小，最后以像對(duì)為單位進(jìn)行統(tǒng)計(jì)，結(jié)果如圖7。

圖6 平差結(jié)果像點(diǎn)中誤差統(tǒng)計(jì)。(a) 嫦娥五號(hào)數(shù)據(jù)處理結(jié)果； (b) 嫦娥三號(hào)數(shù)據(jù)處理結(jié)果
Fig.6 RMS statistics of points of adjustment results

表3 立體像對(duì)相對(duì)姿態(tài)誤差統(tǒng)計(jì)Table 3 Error statistics of relative pose between stereo images

表4 立體像對(duì)距離的誤差統(tǒng)計(jì)Table 4 Error statistics of distance between stereo images

由圖7可以看出，嫦娥五號(hào)科學(xué)驗(yàn)證試驗(yàn)數(shù)據(jù)的處理結(jié)果中，像點(diǎn)誤差約為1個(gè)像素，最大超過3個(gè)像素，嫦娥三號(hào)巡視器全景相機(jī)環(huán)拍數(shù)據(jù)的處理結(jié)果中，像點(diǎn)誤差約為1個(gè)像素，最大超過2個(gè)像素。相對(duì)姿態(tài)誤差導(dǎo)致的像點(diǎn)誤差達(dá)到了像素級(jí)別，對(duì)于最終解算得到的三維坐標(biāo)必然會(huì)產(chǎn)生一定影響。

最后利用附加相對(duì)約束關(guān)系的平差結(jié)果對(duì)影像數(shù)據(jù)進(jìn)行了密集匹配，通過前方交會(huì)的計(jì)算方法得到像點(diǎn)對(duì)應(yīng)的物方三維坐標(biāo)，并生成影像區(qū)域的三維模型。這部分內(nèi)容借助軟件Smart 3D完成，結(jié)果如圖8、圖9。

圖7 相對(duì)姿態(tài)誤差導(dǎo)致的像點(diǎn)誤差分布情況。(a) 嫦娥五號(hào)數(shù)據(jù)處理結(jié)果； (b) 嫦娥三號(hào)數(shù)據(jù)處理結(jié)果
Fig.7 Distribution of image point error caused by inaccurate relative pose

圖8 嫦娥五號(hào)全景相機(jī)科學(xué)驗(yàn)證試驗(yàn)外場(chǎng)區(qū)域三維模型

Fig.8 3D model of outfield for scientific verification test of panoramic cameras in Chang′E-5

圖9 嫦娥三號(hào)巡視器月表環(huán)拍區(qū)域三維模型

Fig.9 3D model of lunar surface around Chang′E-3 rover

4 總結(jié)與展望

本文分析了當(dāng)前月表雙目立體相機(jī)平差處理存在的問題，參考了計(jì)算機(jī)視覺和航空攝影測(cè)量領(lǐng)域的處理經(jīng)驗(yàn)后使用一種改進(jìn)的平差模型，利用嫦娥五號(hào)科學(xué)驗(yàn)證試驗(yàn)數(shù)據(jù)和嫦娥三號(hào)巡視器環(huán)拍數(shù)據(jù)進(jìn)行了驗(yàn)證。通過結(jié)果比較證明，本文使用的方法能夠確保不同立體像對(duì)之間的相對(duì)外方位元素保持一致，提高左右相機(jī)間的內(nèi)部吻合精度，得到的平差結(jié)果更加可靠。另外，該方法還能減少部分外方位元素未知數(shù)，在節(jié)約計(jì)算機(jī)內(nèi)存方面有一定意義。因此，本文提出的平差方法可以應(yīng)用于相關(guān)的月表影像處理工作，為后續(xù)的研究提供支持。