趙棟

(沈陽(yáng)市勘察測(cè)繪研究院有限公司，遼寧 沈陽(yáng) 110000)

1 引 言

20世紀(jì)80年代，光束法平差以及附加參數(shù)的引入，促使很多研究者進(jìn)行了相機(jī)標(biāo)定方法的各種改進(jìn)。比如Ziemann提出的將鉛垂線和光束法平差相結(jié)合的兩步標(biāo)定法，即在標(biāo)定過(guò)程中將畸變差和主距、像主點(diǎn)分開求解，通過(guò)反復(fù)迭代完成標(biāo)定。目前室內(nèi)三維控制場(chǎng)標(biāo)定在攝影測(cè)量和計(jì)算機(jī)視覺(jué)界中使用比較多、精度較高，大量實(shí)驗(yàn)證明，基于高精度室內(nèi)三維控制場(chǎng)的標(biāo)定精度是最高的，而基于平面格網(wǎng)的標(biāo)定方法是最方便的。

隨著普通數(shù)碼相機(jī)的普及，研究各種方便實(shí)用、靈活度高，精度好的相機(jī)標(biāo)定系統(tǒng)已經(jīng)成為當(dāng)今攝影測(cè)量界和計(jì)算機(jī)視覺(jué)界的研究熱點(diǎn)之一。并且問(wèn)題主要集中在使用何種標(biāo)定參考對(duì)象和采用何種算法上。自標(biāo)定算法和自檢校光束法平差是兩種非常有效的標(biāo)定方法，前者利用同名點(diǎn)和約束條件進(jìn)行相機(jī)標(biāo)定靈活方便，后者利用高精度控制點(diǎn)進(jìn)行標(biāo)定，適用于特殊環(huán)境下的高精度測(cè)量。

2 研究趨勢(shì)

目前，找到能適應(yīng)不同投影類型(包括中心投影)，標(biāo)定簡(jiǎn)單、可靠性高、精度高、算法快速的軟件是主要的研究趨勢(shì)。傳統(tǒng)標(biāo)定模型的畸變項(xiàng)已不適合用來(lái)矯正超廣角鏡頭(以下簡(jiǎn)稱魚眼鏡頭)拍攝影像的畸變，因此需要用新的模型和畸變項(xiàng)取代傳統(tǒng)的標(biāo)定模型。投影模型的確定是核心，很多文章提到了有關(guān)投影模型的確定，例如：英向華等提出的一種基于球面透視投影約束的魚眼鏡頭校正方法，此方法通過(guò)擬合一個(gè)與選取的離散點(diǎn)距離最近的圓曲面，輔以多項(xiàng)式的切向和徑向畸變項(xiàng)來(lái)確定相機(jī)內(nèi)參數(shù)。缺陷是不同相機(jī)的投影方式各異，投影面可能并不是圓曲面，而是一個(gè)不規(guī)則曲面[1]。基于幾何特征的標(biāo)定方法通過(guò)提取特征物體并矯正到實(shí)際的幾何模型上來(lái)標(biāo)定內(nèi)參數(shù)，如邱志強(qiáng)等提出的用射影不變性糾正魚眼鏡頭畸變，利用直線物體在糾正影像中也應(yīng)該是直線這樣一種約束來(lái)標(biāo)定相機(jī)內(nèi)參數(shù)，這對(duì)特征點(diǎn)的選取和提取精度要求較高[2]。還有通過(guò)立體標(biāo)定來(lái)確定內(nèi)參，這種方法對(duì)像點(diǎn)匹配精度要求很高。目前用非專業(yè)測(cè)量相機(jī)進(jìn)行低精度測(cè)量已大量普及，用其部分替代傳統(tǒng)測(cè)量可以大大降低測(cè)量成本，研究適用于非專業(yè)測(cè)量相機(jī)的高精度標(biāo)定算法，可使普通相機(jī)在測(cè)量領(lǐng)域的應(yīng)用更加廣泛。

3 標(biāo)定模型

3.1 經(jīng)典標(biāo)定模型

經(jīng)典的模型以中心投影為基礎(chǔ)，加上線性和非線性的畸變項(xiàng)對(duì)圖像邊緣的徑向和切向畸變進(jìn)行糾正。經(jīng)典模型的求解首先通過(guò)二維直接線性變換(DLT)或者求同形矩陣給出外方位元素的初值，然后用光束平差進(jìn)行參數(shù)的優(yōu)化，最后得到可靠的相機(jī)內(nèi)參數(shù)[3，4，5]，式(1)描述了經(jīng)典模型所依賴的共線方程。

(1)

通過(guò)平面格網(wǎng)控制點(diǎn)將其轉(zhuǎn)化為8個(gè)參數(shù)相關(guān)的二維DLT模型，如下：

(2)

式(2)中：L1…L8為二維DLT模型的8個(gè)參數(shù)，由于8個(gè)參數(shù)相關(guān)，故得到的參數(shù)精度不夠高，故只能作為迭代解算的初始值。由于切向畸變和徑向畸變的影響對(duì)高精度測(cè)量不可忽略，如公式(3)所示△x，△y為像點(diǎn)坐標(biāo)的畸變改正數(shù)，其中：K1，K2，K3…為徑向畸變系數(shù)，p1，p2為切向畸變系數(shù)。

(3)

3.2 魚眼相機(jī)成像模型

魚眼相機(jī)由于其超大視場(chǎng)角的原因，不遵守透視投影，故物方點(diǎn)、投影中心和像方點(diǎn)不在一條直線上(除了與主光軸平行的那條光線)，基于共線方程的經(jīng)典模型不適用于魚眼相機(jī)的成像模型，故有必要討論魚眼相機(jī)的投影模型的性質(zhì)。

(1)投影模型

人工選取的區(qū)域，并由計(jì)算機(jī)自動(dòng)識(shí)別的圓心

(ⅰ)中心投影：r=f×tan(θ)

(4)

(ⅱ)立體投影：r=2×f×tan(θ/2)

(5)

(ⅲ)等距投影：r=f×θ

(6)

(ⅳ)等角投影：r=2×f×sin(θ/2)

(7)

(ⅴ)正交投影：r=f×sin(θ)

(8)

式中：f為焦距，θ為入射角，r為投影半徑

5種投影模型中投影半徑r與入射角θ的函數(shù)關(guān)系如圖1所示。

圖1 以上5種投影模型中投影半徑r與入射角θ的函數(shù)關(guān)系

圖2 魚眼相機(jī)投影模型中投影半徑與入射角的關(guān)系圖示

由圖2可以看出，P為物方點(diǎn)，經(jīng)過(guò)投影中心，彎曲投射到了p，而不是沿直線到達(dá)p′。

目前魚眼鏡頭的成像模型大都可以用以上其中一種公式近似表達(dá)，等距投影是理想的球面投影，但是魚眼鏡頭的等效投影面是一個(gè)近似球面的不規(guī)則曲面而非理想曲面，考慮到擬合度，可以選擇更一般的多項(xiàng)式作為投影模型即：

r=k1×θ+k2×θ2+k3×θ3+k4×θ4+k5×θ5+…+kn×θn

(9)

將以上5個(gè)典型的投影模型展開為泰勒多項(xiàng)式：

(ⅰ)中心投影：

(10)

(ⅱ)立體投影：

(11)

(ⅲ)等距投影：r=f×θ

(12)

(ⅳ)等角投影：

(13)

(ⅴ)正交投影：

(14)

由以上式子可以看出理想的魚眼成像模型的泰勒展開式都是奇次項(xiàng)，在θ接近0時(shí)，三次項(xiàng)以上的項(xiàng)可以忽略不計(jì)。隨著θ增大，高次項(xiàng)的作用越來(lái)越明顯。以上投影模型的泰勒展開式呈現(xiàn)出高度的一致性，故可以選擇一個(gè)一般的奇次多項(xiàng)式開擬合各種成像模型的參數(shù)，定為r(θ)，k1、k2、k3等代表成像模型系數(shù)。選擇多少次的模型取決于控制點(diǎn)數(shù)量，精度和特征點(diǎn)提取的精度。對(duì)于多項(xiàng)式擬合，次數(shù)越高，理論上也就越精確，但是需要一定數(shù)量和足夠高精度的觀測(cè)值。否則，高次項(xiàng)沒(méi)有意義。

魚眼鏡頭也存在非線性的對(duì)稱和非對(duì)稱畸變，對(duì)稱畸變可以由成像模型來(lái)修復(fù)，非對(duì)稱畸變可以通過(guò)增加適當(dāng)?shù)膮?shù)來(lái)加以改正。

這里自然而然得從形式上給出不考慮非線性畸變的不同于經(jīng)典模型的魚眼相機(jī)成像模型：

u=r(θ)×cos(φ)×mu+u0

(15)

v=r(θ)+sin(φ)×mv+v0

(16)

以上兩式中：(u，r)為像點(diǎn)坐標(biāo)，θ為像點(diǎn)對(duì)應(yīng)的光線入射角，φ為像點(diǎn)與主點(diǎn)的連線與像坐標(biāo)系中u方向軸的夾角，(u0，v0)為主點(diǎn)坐標(biāo)，mu為影像水平方向上單位長(zhǎng)度所包含的像素個(gè)數(shù)，mv為影像垂直方向上單位長(zhǎng)度所包含的像素個(gè)數(shù)[6～9]。

4 實(shí)驗(yàn)過(guò)程

4.1 特征點(diǎn)提取

特征點(diǎn)的提取在模式識(shí)別中有很多方法例如：二值分割，邊緣算子提取等，這里我們選用計(jì)算機(jī)視覺(jué)開源庫(kù)里OPENCV里的Canny算子進(jìn)行圓形物邊緣的提取以便計(jì)算圓的中心[10]。Canny算子是John F. Canny1986年開發(fā)出來(lái)的一個(gè)多級(jí)邊緣檢測(cè)算法，實(shí)踐證明，這種邊緣檢測(cè)算法是眾多邊緣檢測(cè)算子中表現(xiàn)較突出的一種方法，下面我們?cè)敿?xì)闡述特征點(diǎn)提取的步驟：

(1)對(duì)原圖進(jìn)行中值濾波，然后用Canny算子提取邊緣，如圖3所示。

圖3 使用Canny算子得到的邊緣二值圖像

(2)人工選擇區(qū)域并且由計(jì)算機(jī)識(shí)別圓心，如圖4所示。

圖4 人工選取的區(qū)域，并由計(jì)算機(jī)自動(dòng)識(shí)別的圓心(綠色環(huán)狀物)

4.2 內(nèi)參數(shù)初始化

對(duì)于主點(diǎn)(u0，v0)的初始值，可以由像片的長(zhǎng)度和高度來(lái)決定，即：

(17)

式中：widht為影像以像素為單位的寬度，height為影像以像素為單位的高度，對(duì)于mu，mv，可以由像片大小，視場(chǎng)角和焦距算出，式中符號(hào)含義與前文相同：

(18)

對(duì)于5個(gè)投影參數(shù)k1、k2、k3、k4、k5，我們采用r=k1×θ+k2×θ3進(jìn)行擬合。θ為0°～75°的等間隔采樣點(diǎn)(單位為弧度)，投影半徑r=f×θ，并且暫時(shí)忽略k3、k4、k5。得到至少2個(gè)方程，解出k1、k2作為整體迭代求解的初始值。

4.3 外參數(shù)初始化

由于魚眼影像的變形越靠近邊緣越大，我們只將離主點(diǎn)較近的9個(gè)圓心點(diǎn)作為像方點(diǎn)用二維DLT進(jìn)行外參數(shù)的初始化。

(19)

(X，Y)為控制點(diǎn)坐標(biāo)，(x，y)為相應(yīng)的像點(diǎn)坐標(biāo)，當(dāng)控制點(diǎn)大于等于4個(gè)時(shí)，將式(18)線性化，通過(guò)線性方程求解二維DLT的8個(gè)參數(shù)，由于8個(gè)參數(shù)是內(nèi)外方位元素的函數(shù)值，由這8個(gè)參數(shù)解出外方位元素。當(dāng)平面控制點(diǎn)的z坐標(biāo)為0時(shí)，經(jīng)典的共線方程可表示為式(19)：

(20)

將式(20)轉(zhuǎn)化成式(21)：

(21)

式中：λ=(a3×XS+b3×YS+c3×ZS)

比較式(19)與式(21)，可得出以下8式：

(22)

由上式可得：

(23)

由上式最終可得到：

(24)

主點(diǎn)(x0，y0)與4.2節(jié)提到的主點(diǎn)(u0，v0)等同，f在前文4.2節(jié)已有論述。

(25)

通過(guò)sin(ω)=-b3可求得ω，再由旋轉(zhuǎn)矩陣的正交性：

(26)

由式(22)，可得到：

(27)

XS，YS，ZS的初值可通過(guò)解(27)方程組得到。

4.4 整體求解

觀測(cè)方程如下：

u=r(θ)×cos(φ)×mu+u0

(28)

v=r(θ)×sin(φ)×mv+v0

(29)

由4.2節(jié)求出的3個(gè)旋轉(zhuǎn)角可得到a1，a2，a3，b1，b2，b3，加上3個(gè)平移參數(shù)XS，YS，ZS，根據(jù)坐標(biāo)齊次變換：

(30)

式(30)(X，Y)為物方坐標(biāo)點(diǎn)，(x，y，z)為近似像空間坐標(biāo)點(diǎn)：

(31)

式(19)，式(20)建立了物方點(diǎn)和像方點(diǎn)之間的關(guān)系，由上述關(guān)系，對(duì)于每一對(duì)物方點(diǎn)和像方點(diǎn)可得間接平差的誤差方程：

v=Ax-L

(32)

式(32)中：A為觀測(cè)值對(duì)每個(gè)待求參數(shù)的偏導(dǎo)，V為殘差，L為觀測(cè)值減去近似值。對(duì)每一點(diǎn)對(duì)依次求出像點(diǎn)坐標(biāo)對(duì)內(nèi)參、外參、畸變系數(shù)的偏導(dǎo)數(shù)。但是這里對(duì)內(nèi)外參求偏導(dǎo)數(shù)的解析式較為復(fù)雜。我們應(yīng)用數(shù)值微分，通過(guò)8點(diǎn)的函數(shù)值求得近似偏導(dǎo)數(shù)，由數(shù)值微分理論，求偏導(dǎo)數(shù)的精度可達(dá)o(x^4)。

最后用非線性鄰域里較為可靠的優(yōu)化方法L_M進(jìn)行非線性最小二乘求解，此方法亦將“殘差平方和最小”作為準(zhǔn)則進(jìn)行優(yōu)化。

4.5 影像糾正

4.5.1 直接糾正

通過(guò)原影像坐標(biāo)的反解，得到入射角θ和極角φ，再代入中心投影方程得到糾正后影像的坐標(biāo)，由于重投影后，糾正影像的高度和寬度比原圖像大，會(huì)有空白點(diǎn)，可以對(duì)其進(jìn)行格網(wǎng)插值，選擇雙三次卷積進(jìn)行插值。

具體步驟如下：

求解多項(xiàng)式r=k1×θ+k2×θ3+…+k5×θ9的根，選擇合理的值φ=acos(x/r)，注意r=0的情形；

4.5.2 間接糾正

直接糾正會(huì)導(dǎo)致糾正圖像上存有空白點(diǎn)，間接糾正是從糾正影像上的坐標(biāo)點(diǎn)通過(guò)中心投影反算得到(θ，φ)，再應(yīng)用模型正算映射到原影像，通過(guò)雙三次卷積插值恢復(fù)其灰度。具體介紹其步驟：

5 結(jié)果與討論

基于MFC框架設(shè)計(jì)了一個(gè)從數(shù)據(jù)輸入到結(jié)果輸出操作過(guò)程流程化的工具，對(duì)算法進(jìn)行了實(shí)現(xiàn)，采用了University of Oulu的Juho Kannala and Sami S. Brandt提供的影像數(shù)據(jù)，其相機(jī)為Watec 221S CCD color camera，魚眼鏡頭的焦距為 1.178 mm，糾正結(jié)果如圖5、圖6所示。

圖5 魚眼鏡頭拍攝的影像

圖6 其中一個(gè)通道進(jìn)行糾正后的影像

對(duì)于糾正影像的精度，放棄了用觀測(cè)值殘差來(lái)評(píng)價(jià)，而直接用特殊物體幾何特征。選取糾正影像上共線的一組圓心坐標(biāo)，進(jìn)行直線擬合，得到的標(biāo)準(zhǔn)差為 0.415 7個(gè)像素，糾正精度達(dá)到了攝影測(cè)量的要求。