數(shù)碼相機(jī)定位問(wèn)題研究

1、PAGE PAGE 37數(shù)碼相機(jī)定位問(wèn)題研究趙志剛 薛江堂 冷振鵬摘要基于雙目CCD立體測(cè)量系統(tǒng)標(biāo)定技術(shù)被廣泛用于交通監(jiān)管中，該技術(shù)的核心是攝影測(cè)量。本文利用MATLAB對(duì)圓的邊界提取，再用最小二乘法擬合橢圓曲線，并借助攝影測(cè)量基本公式建立二維模型求得靶標(biāo)圓心像坐標(biāo)。并對(duì)該模型作檢驗(yàn)，再利用針孔模型分別求得兩部相機(jī)像平面對(duì)應(yīng)于物平面的旋轉(zhuǎn)向量和平移向量，進(jìn)而確定兩相機(jī)的相對(duì)位置。首先，我們假設(shè)：已知四點(diǎn)的物平面坐標(biāo)和像平面坐標(biāo)，借助攝影測(cè)量基本公式建立二維模型求得物平面與像平面的對(duì)應(yīng)關(guān)系，但有.八個(gè)未知數(shù)待定（將在第二問(wèn)中給予求解）。其次，我們根據(jù)A、B、C、D、E五個(gè)圓形靶標(biāo)的成像情況利用M

2、ATLAB軟件提取其邊緣閾值。在此基礎(chǔ)上利用最小二乘法擬合橢圓求得其中心坐標(biāo)。在不考慮畸變影響的前提下，該中心坐標(biāo)即為圓心的像坐標(biāo)。任取A、B、C、D、E中四點(diǎn)代入模型一中，即可求得物平面和像平面的對(duì)應(yīng)關(guān)系。再次，我們?cè)诘谌龁?wèn)中分別以A、B、C為研究對(duì)象求出其在模型一的條件下的圓心像坐標(biāo)。利用檢驗(yàn)，比較擬合橢圓中心坐標(biāo)與模型一求得結(jié)果的差異，在置信度情況下，這三組結(jié)果無(wú)顯著差異，從而檢驗(yàn)了模型的精度及穩(wěn)定性。最后利用線性相機(jī)模型（針孔模型）確定世界坐標(biāo)系和計(jì)算機(jī)數(shù)字圖像坐標(biāo)的對(duì)應(yīng)關(guān)系，從而分別求得兩部攝像機(jī)的旋轉(zhuǎn)矩陣和平移向量，從而我們可以求得兩相機(jī)相機(jī)坐標(biāo)系間的關(guān)系：?jiǎn)栴}二的求解結(jié)果如下表：

3、A、B、C、D、E五點(diǎn)在像平面上的像坐標(biāo)坐標(biāo)點(diǎn)像坐標(biāo)（mm）像坐標(biāo)（像素）A（-49.83，51.37）（-188，194）B（-23.51，49.43）（-89，187）C（34.16，45.13）（129，171）D（18.89，-31.47）（71，-119）E（-60.09，-31.21）（-227，-118）注：該坐標(biāo)系原點(diǎn)為像平面的中點(diǎn)，向右為x軸正方向，向上為y軸正方向。問(wèn)題重述（1）問(wèn)題的背景攝影技術(shù)自 20 世紀(jì) 40 年代開始應(yīng)用于交通事故分析，已經(jīng)得到廣泛的應(yīng)用，但僅用做簡(jiǎn)單定性分析，隨著計(jì)算機(jī)視覺(jué)和圖像處理技術(shù)的發(fā)展,攝影測(cè)量技術(shù)在交通事故現(xiàn)場(chǎng)測(cè)量中的應(yīng)用研究已經(jīng)成為熱

4、點(diǎn)，國(guó)內(nèi)外許多學(xué)者已經(jīng)做了大量研究，使定量分析成為可能。如在 80 年代，Kerkoff 對(duì)透視投影發(fā)展的歷史、透視繪圖原理和透視成像原理進(jìn)行了詳細(xì)闡述，根據(jù)透視原理研究了利用攝影圖像確定拖痕長(zhǎng)度等的方法，逐步形成了二維攝影測(cè)量方法，開發(fā)出 Pc-rect 等現(xiàn)場(chǎng)測(cè)量軟件。在 1994 年，Nicholas 等提出了交通事故現(xiàn)場(chǎng)的反投影照片三維重建法 。近年來(lái)，隨著計(jì)算機(jī)視覺(jué)原理的突破和飛速發(fā)展，國(guó)內(nèi)學(xué)者李江教授、許洪國(guó)教授等相繼提出了利用立體視覺(jué)原理的多照片重建交通事故現(xiàn)場(chǎng)的方法，并初步開發(fā)了軟件。為改變標(biāo)定參考物的限制，進(jìn)行照相機(jī)自標(biāo)定研究，魯光泉等提出了基于基礎(chǔ)矩陣的交通事故現(xiàn)場(chǎng)三維重建

5、方法。當(dāng)前主要的交通事故現(xiàn)場(chǎng)重建方法，有二維方法、三維方法等。（2）問(wèn)題的提出用數(shù)碼相機(jī)攝制物體的相片確定物體表面某些特征點(diǎn)的位置，目前最常用的定位方法是雙目定位，即用兩部相機(jī)來(lái)定位。對(duì)物體上一個(gè)特征點(diǎn)，用兩部固定于不同位置的相機(jī)攝得物體的像，分別獲得該點(diǎn)在兩部相機(jī)像平面上的坐標(biāo)。只要知道兩部相機(jī)精確的相對(duì)位置，就可用幾何的方法得到該特征點(diǎn)在固定一部相機(jī)的坐標(biāo)系中的坐標(biāo)，即確定了特征點(diǎn)的位置。于是對(duì)雙目定位，精確地確定兩部相機(jī)的相對(duì)位置是關(guān)鍵，這一過(guò)程稱為系統(tǒng)標(biāo)定。本問(wèn)題的第四小問(wèn)就是解決該問(wèn)題。系統(tǒng)標(biāo)定最常用的一種做法是：在一塊平板上畫若干個(gè)點(diǎn)，同時(shí)用這兩部相機(jī)照相，分別得到這些點(diǎn)在它們像平

6、面上的像點(diǎn)，利用這兩組像點(diǎn)的幾何關(guān)系就可以得到這兩部相機(jī)的相對(duì)位置。然而，無(wú)論在物平面或像平面上我們都無(wú)法直接得到?jīng)]有幾何尺寸的“點(diǎn)”。所以我們實(shí)際的做法是在物平面上畫若干個(gè)圓（稱為靶標(biāo)），它們的圓心就是幾何的點(diǎn)了。而它們的像一般會(huì)變形，所以必須從靶標(biāo)上的這些圓的像中把圓心的像精確地找到，系統(tǒng)標(biāo)定就可實(shí)現(xiàn)。題目中第一問(wèn)和第二問(wèn)就是解決該問(wèn)題。該題目不但要求我們建立模型和算法，而且要求我們對(duì)所見(jiàn)的模型進(jìn)行精度和穩(wěn)定性分析。在前三個(gè)問(wèn)題的基礎(chǔ)上，第四問(wèn)要求我們通過(guò)物平面與像平面的對(duì)應(yīng)關(guān)系找到像平面相對(duì)物平面的旋轉(zhuǎn)矩陣和平移矩陣，從而建立兩部相機(jī)像平面的夾角和距離關(guān)系，以便確定兩部相機(jī)的相對(duì)位置。條

7、件假設(shè)（1）不考慮相機(jī)自身因素導(dǎo)致的誤差（2）對(duì)第二問(wèn)中抽樣選取部分邊界點(diǎn)擬合的橢圓誤差不予考慮（3）兩相機(jī)透鏡光心處于同一水平高度（4）兩部相機(jī)的內(nèi)部參數(shù)是相同（5）本題中所給圖像信息和數(shù)據(jù)真實(shí)準(zhǔn)確（6）不考慮人為因素造成的誤差符號(hào)說(shuō)明：物點(diǎn)在二維圖像平面坐標(biāo)系上的坐標(biāo)：物點(diǎn)在二維物平面坐標(biāo)系上的坐標(biāo)：表示圖像的灰度值：處理前圖像的信息矩陣：處理后圖像的信息矩陣：相機(jī)的像距（這里假定兩部相機(jī)的焦距一樣）：像平面上x軸方向每一像素的物理尺寸：像平面上y軸方向每一像素的物理尺寸：若圖像中心在u、v坐標(biāo)系中的坐標(biāo)為（）（以像素為單位）：被攝像點(diǎn)在世界坐標(biāo)系中坐標(biāo)：被攝像點(diǎn)在攝像機(jī)1中的u、v坐標(biāo)系

8、中的坐標(biāo)：被攝像點(diǎn)在攝像機(jī)1中的u、v坐標(biāo)系中的坐標(biāo)：被攝像點(diǎn)在攝像機(jī)1的相機(jī)坐標(biāo)系中的坐標(biāo)：被攝像點(diǎn)在攝像機(jī)2的相機(jī)坐標(biāo)系中的坐標(biāo)：系統(tǒng)內(nèi)部參數(shù)矩陣：系統(tǒng)外部參數(shù)矩陣：投影矩陣：旋轉(zhuǎn)矩陣：平移向量模型的建立及分析問(wèn)題一的模型：對(duì)于給定的像平面上的靶標(biāo)圓的圖像，我們可以先提取其邊界坐標(biāo)，再根據(jù)最小二乘法進(jìn)行橢圓擬合（在問(wèn)題二中將進(jìn)行詳細(xì)討論），用擬合橢圓的中心近似地代替靶標(biāo)上圓心像的坐標(biāo)。但是，當(dāng)我們只有靶標(biāo)上的圓而沒(méi)有它的圖像時(shí)，是否也可以通過(guò)某種替代關(guān)系，找出它在相機(jī)平面上的坐標(biāo)呢？本文對(duì)其展開如下討論：圖 SEQ 圖表 * ARABIC 1如上圖所示：在靶標(biāo)上建立實(shí)際坐標(biāo)系，在相機(jī)平面內(nèi)

9、建立平面坐標(biāo)系，透視中心表示相機(jī)光心。若在實(shí)際坐標(biāo)系中各待求像坐標(biāo)點(diǎn)在同一平面上，則可以根據(jù)直接線性變換，推出從實(shí)際坐標(biāo)平面到圖像平面的對(duì)應(yīng)關(guān)系，若已知實(shí)際坐標(biāo)系中4個(gè)或4個(gè)以上共面點(diǎn)的空間坐標(biāo)和像坐標(biāo)，這種對(duì)應(yīng)關(guān)系就可以求得。根據(jù)對(duì)應(yīng)關(guān)系，可以求平面上任一靶標(biāo)圓心的像坐標(biāo)了。其原理如下：設(shè)空間坐標(biāo)系中任意點(diǎn)和它在二維圖像平面坐標(biāo)系上的坐標(biāo)為，根據(jù)攝影測(cè)量基本公式空間點(diǎn)和圖像點(diǎn)滿足假設(shè)四個(gè)已知點(diǎn)的空間坐標(biāo)和圖像坐標(biāo)分別是 和 可得如下方程組：通過(guò)對(duì)上述方程組的求解，我們可以得到到，這樣就建立起空間坐標(biāo)系中的點(diǎn)和二維圖像平面對(duì)應(yīng)點(diǎn)之間的關(guān)系，也就是說(shuō)如果知道空間坐標(biāo)系中一點(diǎn)，通過(guò)以下公式就能得到

10、其對(duì)應(yīng)點(diǎn)在二維圖像平面的坐標(biāo)。問(wèn)題二的求解：在問(wèn)題一中，我們建立了求解沒(méi)有圖像的空間點(diǎn)的像坐標(biāo)，它是基于4個(gè)以上含有像坐標(biāo)的空間點(diǎn)的。因此，只有當(dāng)我們把四個(gè)（或以上）像坐標(biāo)求出，模型一才是有效的?，F(xiàn)在就已知空間點(diǎn)和它的圖像求其像坐標(biāo)的問(wèn)題進(jìn)行如下討論： 我們認(rèn)為空間中一個(gè)圓孔經(jīng)過(guò)相機(jī)成像后，在圖像平面上呈現(xiàn)為橢圓形狀。對(duì)于圖像中的圓中心像坐標(biāo)，可以用如下方法進(jìn)行提取：自動(dòng)閾值分割圖像中的橢圓，提取橢圓輪廓并用輪廓數(shù)據(jù)進(jìn)行最小二乘橢圓擬合，求出橢圓的中心作為空間圓孔所對(duì)應(yīng)的圖像坐標(biāo)，下面是具體算法。1.取ROI圖像灰度范圍內(nèi)的中值作為初始閾值： 其中，I表示圖像的灰度值。2.根據(jù)將圖像分割成目標(biāo)

11、和背景兩個(gè)部分，求出兩部分的平均灰度值，。其中為權(quán)重系數(shù)，一般取=1(1)求出新的閾值：如果，結(jié)束迭代。去迭代結(jié)束時(shí)的為最佳分割閾值對(duì)ROI圖像進(jìn)行二值化。但是，在本題中，圖像的坐標(biāo)已進(jìn)行二值化，因此上述步驟本文并未進(jìn)行。(2)提取橢圓輪廓并擬合橢圓包含源口特征的ROI圖像進(jìn)行輪廓提取，其提取方法如下： 首先將本題中所給的分別率為的圖像信息讀入到MATLAB的一數(shù)組A中，然后利用MATLAB中的邊界函數(shù)edge對(duì)數(shù)組A進(jìn)行處理，通過(guò)這一過(guò)程我們可將圖像中邊界的相關(guān)信息讀入另一數(shù)組B中，再通過(guò)對(duì)數(shù)組B的簡(jiǎn)單處理，我們便可得到圖像的輪廓信息。得到如下圖所示的圖像：圖 SEQ 圖表 * ARABIC

12、 2 處理前的圖像（即本文中數(shù)組A對(duì)應(yīng)的圖像）圖 SEQ 圖表 * ARABIC 3 處理后的圖像（即本文中數(shù)組B對(duì)應(yīng)的圖像）注：數(shù)組A、B都是的二維數(shù)組，其對(duì)應(yīng)元素值為0和1 其中圖像上白色的區(qū)域值為1，黑色的區(qū)域值為0從而得到橢圓輪廓數(shù)據(jù)，下表為其中一個(gè)橢圓的輪廓上的坐標(biāo)值對(duì)以E點(diǎn)為圓心的圓在像中部分邊界點(diǎn)坐標(biāo)提取編號(hào)坐標(biāo)編號(hào)坐標(biāo)1（502，246）12（493，324）2（509，246）13（502，324）3（483，254）14（517，318）4（475，262）15（527，310）5（471，270）16（532，302）6（468，278）17（535，294）7（467，

13、286）18（537，286）8（467，294）19（537，278）9（469，302）20（536，270）10（472，310）21（533，262）11（480，318）22（526，254）注：該坐標(biāo)系為以像左上角為原點(diǎn)，向下為x軸，向右為y軸，以像素為坐標(biāo)單位這樣通過(guò)這些邊界點(diǎn)坐標(biāo)就可以利用最小二乘法實(shí)現(xiàn)橢圓擬合。(3)求解橢圓的中心作為空間圓孔中心所對(duì)應(yīng)的圖像坐標(biāo)橢圓的中心坐標(biāo)可以通過(guò)如下公式求得：對(duì)于曲線的中心坐標(biāo)可以由下式得出：因此，擬合橢圓后，可以將橢圓方程的系數(shù)代入上式求出圓孔中心圖像坐標(biāo)。我們依然以E點(diǎn)為研究對(duì)象，我們把上文表格中的E值橫坐標(biāo)和縱坐標(biāo)分存在兩個(gè)數(shù)組中:

14、X=502,509,483,475,526;Y=246,246,254,262,254。 在本文中，我們直接采用美國(guó)人Andrew Fitzgibbon, Maurizio Pilu, Bob Fisher編的最小二乘法擬合橢圓，采用上述原理直接返回橢圓的中心E（-227，-118）。依照此方法我們求得A、B、C、D、E的像素坐標(biāo)和實(shí)際坐標(biāo)如下表：A、B、C、D、E的像素坐標(biāo)A點(diǎn)坐標(biāo)（-188，194）B點(diǎn)坐標(biāo)（-89，187）C點(diǎn)坐標(biāo)（129，171）D點(diǎn)坐標(biāo)（71，-119）E點(diǎn)坐標(biāo)（-227，-118）注：該坐標(biāo)系原點(diǎn)為像平面的中點(diǎn)，向右為x軸，向上為y軸，以像素為坐標(biāo)單位A、B、C、D

15、、E五點(diǎn)在像平面上的像坐標(biāo)A點(diǎn)坐標(biāo)（-49.83，51.37）B點(diǎn)坐標(biāo)（-23.51，49.43）C點(diǎn)坐標(biāo)（34.16，45.13）D點(diǎn)坐標(biāo)（18.89，-31.47）E點(diǎn)坐標(biāo)（-60.09，-31.21）注：該坐標(biāo)系原點(diǎn)為像平面的中點(diǎn)，向右為x軸，向上為y軸，以毫米為坐標(biāo)單位結(jié)果的檢驗(yàn)：因?yàn)槲锲矫鎯?nèi)A、B、C三點(diǎn)共線，這三點(diǎn)在像平面內(nèi)也應(yīng)該共線即可以近似地認(rèn)為。從而驗(yàn)證了我們第二問(wèn)結(jié)果的科學(xué)性。問(wèn)題三：基于成對(duì)數(shù)據(jù)的檢驗(yàn)：為了比較像平面中的擬合橢圓中心點(diǎn)坐標(biāo)和建立多元線性回歸擬合出來(lái)的像平面上的靶標(biāo)中心坐標(biāo)的差異，我們?nèi)∪龑?duì)坐標(biāo) ， ， 作對(duì)比試驗(yàn)，得到一批成對(duì)的觀察值，然后分析觀察數(shù)據(jù)作出

16、推斷，由于A、B、C三點(diǎn)坐標(biāo)與其位置有關(guān)，不能看作三個(gè)獨(dú)立的隨機(jī)變量觀察值，因此我們采用逐對(duì)比較法。如果不考慮橫縱坐標(biāo)的差異，我們可以將試驗(yàn)對(duì)象看作6對(duì)相互獨(dú)立的觀察結(jié)果 。其中 其他類似。再令 則相互獨(dú)立，.相互獨(dú)立，又由于，.是由同一因素引起的，可認(rèn)為它們服從同一分布。123456橢圓中心坐標(biāo)值M-49.8351.37-23.5149.4334.1645.13模型求取值N-50.1651.42-23.2549.4933.5145.18差值Q0.26-0.057-0.33-0.060.65-0.05今假設(shè)，=1,2,3,4,5,6 這就是說(shuō)，.構(gòu)成正態(tài)分布的一個(gè)樣本，其中未知，我們需要基于這

17、一樣本檢驗(yàn)假設(shè)。原假設(shè)，備擇假設(shè)分別假設(shè)，.的樣本均值和樣本方差的觀察值為，按基于單個(gè)正態(tài)分布總體均值的檢驗(yàn)，通過(guò)概率知識(shí)知檢驗(yàn)的拒絕域?yàn)?，其中，為置信度，這里取0.05即可知拒絕域?yàn)橛捎^察值，現(xiàn)的值不落在拒絕域內(nèi)，故接受。認(rèn)為用多元線性回歸擬合出來(lái)像平面上靶標(biāo)中心坐標(biāo)與像平面的擬合橢圓中心坐標(biāo)無(wú)明顯差異。模型的精度和穩(wěn)定性討論：對(duì)于圓形靶標(biāo)，經(jīng)過(guò)透視投影變換，在攝像機(jī)平面所成的像為橢圓，在圓心像點(diǎn)位置測(cè)定過(guò)程中，我們模型的作法是通過(guò)圖像處理技術(shù)得到投影橢圓的幾何中心，并近似認(rèn)為它是空間圓的圓心在攝像機(jī)像平面上的準(zhǔn)確二維像點(diǎn)，然而由于透視變換的固有特性，只有當(dāng)圓所在平面與攝像機(jī)圖像平面平行時(shí)，

18、這種近似才是正確的。在二者不平行時(shí)，用像平面上投影橢圓的幾何中心近似空間圓心的實(shí)際像點(diǎn)，必然在測(cè)量中引入誤差，這種誤差就是所謂的畸變誤差。隨著圓所在平面與攝像機(jī)圖像平面之間夾角的增加，畸變誤差將增大。從某種意義上說(shuō)該誤差影響了模型的精度和穩(wěn)定性。但考慮實(shí)際情況，用于交通監(jiān)管的數(shù)碼相機(jī)一般距拍攝目標(biāo)距離不大，而且畸變誤差的范圍不過(guò)幾微米（圓平面與像平面夾角為0時(shí)除外）到幾十微米，因此該模型可近似看作穩(wěn)定且準(zhǔn)確的，畸變誤差不影響該模型的應(yīng)用。問(wèn)題四的模型 在前三問(wèn)的基礎(chǔ)上，我們對(duì)空間點(diǎn)像坐標(biāo)的求解已不再陌生。在圖像特征點(diǎn)坐標(biāo)及其在像平面的坐標(biāo)已知的前提下，如何確定兩部相機(jī)的相對(duì)位置呢？以下對(duì)該問(wèn)題

19、展開討論：建模準(zhǔn)備：1.圖像坐標(biāo)系、照相機(jī)坐標(biāo)系與世界坐標(biāo)系的定義圖 4 相機(jī)成像原理圖4為一個(gè)位于三維空間的照相機(jī)成像幾何模型，其中各坐標(biāo)系的定義如下：稱為世界坐標(biāo)系，該坐標(biāo)系的設(shè)置一般是根據(jù)實(shí)際情況加以考慮。稱為照相機(jī)坐標(biāo)系，其原點(diǎn)位于鏡頭的物方主點(diǎn)上，軸與鏡頭的光軸重合，方向由鏡頭指向物方。稱為成像平面二維坐標(biāo)系。在理想的光學(xué)系統(tǒng)中，成像平面的原點(diǎn)位于光軸上，但實(shí)際中存在一定的偏移量。稱為計(jì)算機(jī)數(shù)字圖像平面二維坐標(biāo)系。該平面是對(duì)鏡頭成像平面數(shù)字量化的結(jié)果，其原點(diǎn)一般位于圖像的左上角。照相機(jī)采集的圖像以標(biāo)準(zhǔn)電視信號(hào)的形式輸入計(jì)算機(jī)，經(jīng)計(jì)算機(jī)中的專用數(shù)模轉(zhuǎn)換板變換成數(shù)字圖像。每幅數(shù)字圖像在計(jì)

20、算機(jī)內(nèi)為MN數(shù)組M行N列的圖像中的每一個(gè)元素(稱為象素)的數(shù)值即是圖像點(diǎn)的亮度(或稱灰度。若為彩色圖像，則圖像上象素的亮度將由紅、綠、藍(lán)三種顏色的亮度表示。)，如圖5所示，在圖像上定義直角坐標(biāo)系u、v，每一象素的坐標(biāo)(u，v)分別是該象素在數(shù)組中的列數(shù)與行數(shù)，所以(u，v)是以象素為單位的圖像坐標(biāo)系的坐標(biāo)。由于(u，v)只表示象素位于數(shù)組中的列數(shù)與行數(shù)，并沒(méi)有用物理單位表示出該象素在圖像中的位置，因而，需要再建立以物理單位(例如毫米)表示的圖像坐標(biāo)系。該坐標(biāo)系以圖像內(nèi)某一點(diǎn)為原點(diǎn)，x軸與y軸分別與u，v軸平行，(x，y)表示以毫米為單位的圖像坐標(biāo)系的坐標(biāo)。在x、y坐標(biāo)系中，原點(diǎn)定義在相機(jī)光軸與

21、圖像平面的交點(diǎn)，該點(diǎn)一般位于圖像中心處，但由于相機(jī)制作的原因，也會(huì)有些偏離，若在u、v坐標(biāo)系中的坐標(biāo)為（），每一個(gè)象素在x軸與y軸方向上的物理尺寸為dx,dy，則圖像中任意一個(gè)象素在兩個(gè)坐標(biāo)系下的坐標(biāo)有如下關(guān)系： 圖5 圖像坐標(biāo)系為以后使用方便，我們用齊次坐標(biāo)與矩陣形式將上式表示 （*1）逆關(guān)系可寫為： （*2）照相機(jī)成像幾何關(guān)系可由圖5表示。其中O點(diǎn)稱為照相機(jī)光心，照相機(jī)成像幾何關(guān)系可由圖5表示。其中O點(diǎn)稱為照相機(jī)光心，軸和軸與圖像的x軸與y軸平行，軸為相機(jī)的光軸，它與圖像平面垂直。光軸與圖像平面的交點(diǎn)，即為圖像坐標(biāo)系的原點(diǎn)，由點(diǎn)O與軸組成的直角坐標(biāo)系稱為相機(jī)坐標(biāo)系，為相機(jī)焦距。 圖 6 相

22、機(jī)坐標(biāo)系與世界坐標(biāo)系由于照相機(jī)可位于環(huán)境中的任何位置，可以在環(huán)境中再選擇一個(gè)基準(zhǔn)坐標(biāo)系來(lái)描述照相機(jī)的位置，并用它描述環(huán)境中任何物體的位置，該坐標(biāo)系稱為世界坐標(biāo)系。它由，軸組成。照相機(jī)坐標(biāo)系與世界坐標(biāo)系之間的關(guān)系可以用旋轉(zhuǎn)矩陣R與平移向量t來(lái)描述。因此，空間中某一點(diǎn)P在世界坐標(biāo)系與照相機(jī)坐標(biāo)系下的齊次坐標(biāo)如果分別是和，于是存在如下關(guān)系：（*3）其中，R為33正交單位矩陣；t為三維平移向量；0，為44矩陣。2.線性相機(jī)模型空間任何一點(diǎn)P在圖像上的成像位置可以用針孔模型近似描述，即任何點(diǎn)P在圖像上的投影位置p為光心O與P點(diǎn)的連線OP與圖像平面的交點(diǎn)。這種關(guān)系也稱為中心射影或透視投影。由比例關(guān)系有如下

23、關(guān)系式 （*4）其中，(x，y)為p點(diǎn)的圖像坐標(biāo)；為空間點(diǎn)P在相機(jī)坐標(biāo)系下的坐標(biāo)。用齊次坐標(biāo)與矩陣表示上述透視投影關(guān)系 （*5）將式(*2)與式(*3)代入式(*5)，得到以世界坐標(biāo)系表示的P點(diǎn)坐標(biāo)與其投影點(diǎn)p的坐標(biāo)(u，v)的關(guān)系： (*6)其中 ，為34矩陣，稱為投影矩陣； 完全由決定，由于 只與照相機(jī)內(nèi)部結(jié)構(gòu)有關(guān)，稱這些參數(shù)為照相機(jī)內(nèi)部參數(shù)； 完全由照相機(jī)相對(duì)于世界坐標(biāo)系的方位決定，稱為照相機(jī)外部參數(shù)，確定某一照相機(jī)的內(nèi)外參數(shù)，稱為照相機(jī)定標(biāo)。在，（u，v），和（）已知的情況向下，運(yùn)用matlab強(qiáng)大的矩陣處理功能可以求得：= (*7)由矩陣的對(duì)應(yīng)關(guān)系可的和模型的建立: 對(duì)于該問(wèn)題，如果

24、給出如下數(shù)據(jù)，即可以確定出兩部相機(jī)的相對(duì)位置：相機(jī)的像距（這里假定兩部相機(jī)的焦距一樣）：像平面上x軸方向每一像素的物理尺寸：像平面上y軸方向每一像素的物理尺寸：若圖像中心在u、v坐標(biāo)系中的坐標(biāo)為（）（以像素為單位）：被攝像點(diǎn)在世界坐標(biāo)系中坐標(biāo)：被攝像點(diǎn)在攝像機(jī)1中的u、v坐標(biāo)系中的坐標(biāo)：被攝像點(diǎn)在攝像機(jī)1中的u、v坐標(biāo)系中的坐標(biāo)：被攝像點(diǎn)在攝像機(jī)1的相機(jī)坐標(biāo)系中的坐標(biāo)：被攝像點(diǎn)在攝像機(jī)2的相機(jī)坐標(biāo)系中的坐標(biāo)通過(guò)下面的公式可分別求出和 (*8) (*9)圖 7接下來(lái)我們建立兩相機(jī)相機(jī)坐標(biāo)系間的關(guān)系：通過(guò)以上兩公式便可得到兩相機(jī)的相機(jī)坐標(biāo)系間的旋轉(zhuǎn)矩陣和平移向量，也就是這里需要的兩相機(jī)的相對(duì)位置。

25、模型進(jìn)一步討論基于透視投影變換,空間圓在攝像機(jī)像平面上所成的像是橢圓曲線5 . 為了準(zhǔn)確獲取空間圓中心,目前視覺(jué)檢測(cè)中本文采用的方法的方法是:先獲得空間圓中心在攝像機(jī)像平面上的二維像點(diǎn),然后基于視覺(jué)檢測(cè)模型求出其空間三維點(diǎn)坐標(biāo),即為空間圓中心. 并且認(rèn)為通過(guò)對(duì)相機(jī)像平面上的橢圓曲線進(jìn)行擬合得到的橢圓中心就是空間圓中心在相機(jī)像平面上的準(zhǔn)確二維像點(diǎn). 然而,這種處理方法只是一種近似.由于透視投影變換所固有的特性,只有在空間圓截面與相機(jī)像平面平行時(shí),這種處理方法才是準(zhǔn)確的. 在大多數(shù)情形下,存在一定的偏差,稱其為畸變誤差. 如果在實(shí)際應(yīng)用中直接用像平面上擬合得到的橢圓中心像點(diǎn)代替空間橢圓中心的實(shí)際像

26、點(diǎn),必然在測(cè)量中引入誤差.針對(duì)此問(wèn)題,應(yīng)該考慮畸變對(duì)的影響。1， 空間圓的攝像機(jī)針孔成像原理空間橢圓的攝像機(jī)針孔成像原理即透視投影變換如圖1 所示： 圖 8 空間圓的像機(jī)針孔成像原理 在圖8 中,O 為投影中心,即光學(xué)中心. O1 為空間圓的中心. Oxyz 為攝像機(jī)坐標(biāo)系, Oxyz為新建坐標(biāo)系,稱之為空間圓坐標(biāo)系, 其坐標(biāo)原點(diǎn)O與O 重合, z軸的正向指向正交于空間圓O1 A1 B1 C1 D1 所在平面,即為其法線方向. 為攝像機(jī)像平面,與z 軸正交. 為O1 A1 B1 C1 D1 在上所成的曲線像. O2 A2 B2 C2 D2 為O1 A1 B1 C1 D1在xOy坐標(biāo)面上的正交投

27、影. 設(shè)O1 在Oxyz中的坐標(biāo)為( x , y , z) ,矢量OO1在Oxyz中的方向角分別為,建立空間圓中心O1 在攝像機(jī)像平面上的成像畸變誤差模型： 空間橢圓在上的像橢圓中心坐標(biāo)模型，為了建模方便起見(jiàn), 不失一般性, 考慮O1 A1 B1 C1 D1 的如下位置情形:O1 A1 B1 C1 D1 的某一軸A1 C1 與矢量OO1 正交,另一軸為B1 D1。通過(guò)透視中心O 形成一個(gè)斜橢圓錐面, 在Oxyz下記為OA1 B1 C1 D1 , 在Oxyz 下記為OA1 B1 C1 D1 ,則OA1 B1 C1 D1 的外形線A1 OC1 為一等腰三角形. 設(shè)= A1 OO1 =C1 OO1

28、,又假設(shè),當(dāng)OO1 D1 90時(shí),令 90時(shí), 令 =B1 OO1 ;當(dāng)OO1 B1 D1時(shí),令 = D1 OO1 =B1 OO1 . 令矢量OB2 與x軸正向沿逆時(shí)針?lè)较虻膴A角為.O1 A1 B1 C1 D1 的像平面曲線 的方程建立過(guò)程如下:在Oxyz 下,為與OA1 B1 C1 D1 的截交線. 設(shè)的方程為F ( x , y , z) = 0 , OA1 B1 C1 D1的方程為G( x , y , z) = 0 ,則的方程為F( x , y , z) = 0G( x , y , z) = 0(1)在平面z= z內(nèi),O1 A1 B1 C1 D1 在Oxyz下的方程可以表達(dá)為（2）對(duì)于矢量

29、OO2 ,有如下關(guān)系成立: = = = (3)所以有 =，,則= = = =- (4)由式(3) 和式(4) 、式(2) 進(jìn)一步整理,并令 為變量, 即用z代替 , 則得到在Oxyz下OA1 B1 C1 D1的方程為+ （5）由圖1 可見(jiàn),Oxyz 與Oxyz只存在旋轉(zhuǎn)關(guān)系,設(shè)由下式給出 (6)令,利用式(6) 、式(5)轉(zhuǎn)換為Oxyz 下OA1 B1 C1 D1 的方程為依式(1) ,將: z = f 帶入式(7) ,并令k = u sin( r1 cos - r7cos) -cos( r4 cos - r7 cos) l = u sin( r2 cos - r8 cos) -cos( r5

30、 cos - r8 cos) m = u sin( r3 cos - r9 cos) -cos( r6 cos - r9 cos) fn = tan cos( r1 cos - r7 cos) +sin( r4cos - r7 cos) p = tan cos( r2 cos - r8 cos) +sin( r5cos - r8 cos) q = tan cos( r3 cos - r9 cos) +sin( r6 cos - r9 cos) fr = ur7 tans = ur8 tant = ur9 f tan （7）則得到在Oxyz 下,O1 A1 B1 C1 D1 在上的透視投影的方程

31、為( k2 + n2 - r2 ) x2 + 2 ( kl + np - rs) xy +( l2 + p2 - s2 ) y2 + 2 ( mk + nq - rt) x +2 ( ml + pq - st) y + m2 + q2 - t2 = 0 (8)式(8) 為一二次曲線表達(dá)式,它可以表示橢圓、圓、雙曲線、拋物線或者直線. 但在實(shí)際的視覺(jué)檢測(cè)系統(tǒng)中,被測(cè)物體總是處于像平面的前方. 因此,O1 A1 B1 C1 D1 在上所成的像 為一閉合的曲線,且為橢圓或者圓. 因此, 的中心坐標(biāo)( a , b)由下式給出:a = ( l2 + p2 - s2 ) ( mk + np - rt) -

32、( kl + np - rs) ( ml + pq - st) / ( kl + np - rs) 2 - ( k2 + n2 -r2 ) ( l2 + p2 - s2 ) b = ( k2 + n2 - r2 ) ( ml + pq - st) -( kl + np - rs) ( mk + nq - rt) / ( kl + np - rs) 2 - ( k2 + n2 -r2 ) ( l2 + p2 - s2 ) (9)2. 2 空間橢圓中心在上的實(shí)際像點(diǎn)坐標(biāo)模型O1 A1 B1 C1 D1 的中心O1 在上的實(shí)際像點(diǎn)的坐標(biāo)可以利用直線OO1 與的交點(diǎn)來(lái)求取. 首先建立OO1 的方程.O

33、O1 即直線OO1 在Oxyz下的方程為 (10)由式(3) 和式(6) ,得到OO1 在Oxyz 下的方程為( r1 x + r2 y + r3 z) (11)OO1 與: z = f 的交點(diǎn)即為O1 點(diǎn)的實(shí)際像點(diǎn). 把z = f 帶入式(11) ,得( r1 cos - r7 cos) x + ( r2 cos - r8 cos) y =r9 f cos - r3 fcos( r4 cos - r7 cos) x + ( r5 cos - r8 cos) y =r9 f cos - r6 f cos (12)由式(12) 得到O1 點(diǎn)的實(shí)際像點(diǎn)的坐標(biāo)( a, b) 的表達(dá)式為(13)其中k

34、1 = r1cos - r7 cosl1 = r2 cos - r8 cosm1 = ( r3 cos - r9 cos) fn1 = r4 cos - r7 cosp1 = r5cos - r8 cosq1 = ( r6 cos - r9 cos) f 成像畸變誤差模型式(9) 和式(13) 只是在O1 A1 B1 C1 D1 所在的平面與平行時(shí)才有a = a, b = b. 當(dāng)二者之間存在一定夾角的時(shí)候,它們并不相等. 用這2 個(gè)點(diǎn)之間的距離d 來(lái)描述它們之間偏差的大小,即空間橢圓中心的攝像機(jī)針孔成像(即透視投影變換)畸變誤差數(shù)學(xué)模型為 (14)式(14) 中, a , b 和a, b已分

35、別由式(9) 和式(13)具體給出.模型的改進(jìn)：我們?cè)诘谒膯?wèn)建立的模型是線性模型不能準(zhǔn)確地描述成像幾何關(guān)系，尤其在使用廣角鏡頭時(shí)，在遠(yuǎn)離圖像中心處會(huì)有較大的畸變?？紤]到各種畸變影響，我們采用下列公式描述非線性畸變： （15）其中，為由小孔線性模型計(jì)算出來(lái)的圖像坐標(biāo)的理想值；(x,y)是實(shí)際的圖像點(diǎn)的坐標(biāo)，與是非線性畸變值，它與圖像點(diǎn)在圖像中的位置有關(guān)，可用以下公式表達(dá)： （16）是非線性畸變系數(shù)。一般情況下，非線性模型的徑向畸變已經(jīng)能夠滿足描述非線性畸變，國(guó)內(nèi)外研究表明引入過(guò)多的非線性參數(shù)不僅不能提高精度，反而引起解的不穩(wěn)定，根據(jù)需要選擇合適的非線性模型。 參考文獻(xiàn)1 阮秋琦編著數(shù)字圖像處理學(xué)

36、M北京：電子工業(yè)出版社，2001.2 邱茂林，馬頌德，李毅.計(jì)算機(jī)視覺(jué)中攝像機(jī)定標(biāo)綜述J.自動(dòng)化學(xué)報(bào)，2006， 26（I）：43：55.3 葉建軍.直接線性變換解法的適用性J.測(cè)繪技術(shù)，1990,(2):16-19.4 李欣菊. 攝像機(jī)標(biāo)定算法研究D, 中科院自動(dòng)化研究所, 2003.5 孟曉橋, 胡占義.一種新的基于圓環(huán)點(diǎn)的攝像機(jī)自標(biāo)定方法J.軟件學(xué)報(bào), 2002, 13(5):957964.6 李德仁, 鄭肇葆. 解析攝影測(cè)量. 北京: 測(cè)繪出版社,2002.7 張劍清, 張祖勛. 數(shù)字?jǐn)z影測(cè)量. 武漢: 武漢測(cè)繪科技大學(xué)出版社, 2006.8 方德植, 陳奕培. 射影幾何. 北京: 高

37、等教育出版社,2003.9 李鵬, 王軍寧.攝像機(jī)標(biāo)定方法簡(jiǎn)述J.技術(shù)與應(yīng)用, 2007: 21- 24.10毛劍飛, 鄒細(xì)勇, 諸靜.改進(jìn)的平面模板兩步法標(biāo)定攝像機(jī)J.中國(guó)圖象圖形學(xué)報(bào), 2004, 9( 7) : 846- 852.11秦大輝. 結(jié)構(gòu)光二維測(cè)量系統(tǒng)關(guān)鍵技術(shù)的研究D . 華中科技大學(xué), 2005: 2 - 12.12馬頌德,張正友. 計(jì)算機(jī)視覺(jué)計(jì)算理論與算法基礎(chǔ)M . 北京:科學(xué)出版社, 1998: 52 - 93.13吳曉良.數(shù)字影像明顯目標(biāo)的精確定位D.武漢:武漢測(cè)繪科技大學(xué),1999.14楊啟帆，何勇，談之奕.數(shù)學(xué)建模競(jìng)賽。浙江大學(xué)出版社，2005.15盛驟，謝式千，潘

38、承毅.概率論與數(shù)理統(tǒng)計(jì)（第三版）.北京：高等教育出版社，2000.16David R. Anderson,Dennis J. Sweeney,Thomas A. Williams等著，于淼等譯.數(shù)據(jù)模型與決策.機(jī)械工業(yè)出版社，2002.17葉其孝.大學(xué)生數(shù)學(xué)建模競(jìng)賽輔導(dǎo)教材（一）.長(zhǎng)沙：湖南教育出版社，1998.19張培強(qiáng).MATLAB語(yǔ)言演算紙式的科學(xué)工程計(jì)算語(yǔ)言.合肥:中國(guó)科技大學(xué)出版社，1995.20成禮智等.數(shù)學(xué)建模的理論與實(shí)踐.北京：中國(guó)科技大學(xué)出版社，1999.21數(shù)碼相機(jī)定位.www. /view/1870113.html附錄最小二乘法擬合橢圓的MATLAB程序function

39、a = fitellipse(X,Y)% FITELLIPSE Least-squares fit of ellipse to 2D points.% A = FITELLIPSE(X,Y) returns the parameters of the best-fit% ellipse to 2D points (X,Y).% The returned vector A contains the center, radii, and orientation% of the ellipse, stored as (Cx, Cy, Rx, Ry, theta_radians)% Example:

40、Run fitellipse without any arguments to get a demoif nargin = 0 % Create an ellipse t = linspace(0,2);Rx = 300 Ry = 200 Cx = 250 Cy = 150 Rotation = .4 % Radians x = Rx * cos(t); y = Ry * sin(t); nx = x*cos(Rotation)-y*sin(Rotation) + Cx; ny = x*sin(Rotation)+y*cos(Rotation) + Cy; % Draw it plot(nx,

41、ny,o); % Fit it fitellipse(nx,ny) % Note it returns (Rotation - pi/2) and swapped radii, this is fine. returnend% normalize datamx = mean(X);my = mean(Y);sx = (max(X)-min(X)/2;sy = (max(Y)-min(Y)/2;x = (X-mx)/sx;y = (Y-my)/sy;% Force to column vectorsx = x(:);y = y(:);% Build design matrixD = x.*x x

42、.*y y.*y x y ones(size(x) ;% Build scatter matrixS = D*D;% Build 6x6 constraint matrixC(6,6) = 0; C(1,3) = -2; C(2,2) = 1; C(3,1) = -2;% Solve eigensystemgevec, geval = eig(S,C);% Find the negative eigenvalueI = find(real(diag(geval) 0.2)=0;Y=fft2(double(IA);Y=fftshift(Y);Ya=Y.*Hd;Ya=ifftshift(Ya);I

43、a=ifft2(Ya);figuresubplot(2,2,1),imshow(uint8(IA);subplot(2,2,2),imshow(uint8(Ia);figuresurf(Hd,Facecolor,interp,Edgecolor,none,Facelighting,phong); 二、理想高通濾波器IA=imread(lena.bmp);f1,f2=freqspace(size(IA),meshgrid);Hd=ones(size(IA);r=sqrt(f1.2+f2.2);Hd(r0.2)=0;Y=fft2(double(IA);Y=fftshift(Y);Ya=Y.*Hd;

44、Ya=ifftshift(Ya);Ia=real(ifft2(Ya);figuresubplot(2,2,1),imshow(uint8(IA);subplot(2,2,2),imshow(uint8(Ia);figuresurf(Hd,Facecolor,interp,Edgecolor,none,Facelighting,phong); Butterworth低通濾波器IA=imread(lena.bmp);f1,f2=freqspace(size(IA),meshgrid);D=0.3;r=f1.2+f2.2;n=4;for i=1:size(IA,1) for j=1:size(IA,

45、2) t=r(i,j)/(D*D); Hd(i,j)=1/(tn+1); endendY=fft2(double(IA);Y=fftshift(Y);Ya=Y.*Hd;Ya=ifftshift(Ya);Ia=real(ifft2(Ya);figuresubplot(2,2,1),imshow(uint8(IA);subplot(2,2,2),imshow(uint8(Ia);figuresurf(Hd,Facecolor,interp,Edgecolor,none,Facelighting,phong); Butterworth高通濾波器IA=imread(lena.bmp);f1,f2=fr

46、eqspace(size(IA),meshgrid);D=0.3;r=f1.2+f2.2;n=4;for i=1:size(IA,1) for j=1:size(IA,2) t=(D*D)/r(i,j); Hd(i,j)=1/(tn+1); endendY=fft2(double(IA);Y=fftshift(Y);Ya=Y.*Hd;Ya=ifftshift(Ya);Ia=real(ifft2(Ya);figuresubplot(2,2,1),imshow(uint8(IA);subplot(2,2,2),imshow(uint8(Ia);figuresurf(Hd,Facecolor,int

47、erp,Edgecolor,none,Facelighting,phong); 高斯低通濾波器IA=imread(lena.bmp);IB=imread(babarra.bmp);f1,f2=freqspace(size(IA),meshgrid);D=100/size(IA,1);r=f1.2+f2.2;Hd=ones(size(IA);for i=1:size(IA,1) for j=1:size(IA,2) t=r(i,j)/(D*D); Hd(i,j)=exp(-t); endendY=fft2(double(IA);Y=fftshift(Y);Ya=Y.*Hd;Ya=ifftshif

48、t(Ya);Ia=real(ifft2(Ya);figuresubplot(2,2,1),imshow(uint8(IA);subplot(2,2,2),imshow(uint8(Ia);figuresurf(Hd,Facecolor,interp,Edgecolor,none,Facelighting,phong); 高斯高通濾波器IA=imread(lena.bmp);IB=imread(babarra.bmp);f1,f2=freqspace(size(IA),meshgrid);%D=100/size(IA,1);D=0.3;r=f1.2+f2.2;for i=1:size(IA,1)

49、 for j=1:size(IA,2) t=r(i,j)/(D*D); Hd(i,j)=1-exp(-t); endendY=fft2(double(IA);Y=fftshift(Y);Ya=Y.*Hd;Ya=ifftshift(Ya);Ia=real(ifft2(Ya);figuresubplot(2,2,1),imshow(uint8(IA);subplot(2,2,2),imshow(uint8(Ia);figuresurf(Hd,Facecolor,interp,Edgecolor,none,Facelighting,phong); 梯形低通濾波器IA=imread(lena.bmp)

50、;IB=imread(babarra.bmp);f1,f2=freqspace(size(IA),meshgrid);%D=100/size(IA,1);D0=0.1;D1=0.4;r=sqrt(f1.2+f2.2);Hd=zeros(size(IA);Hd(r=D0 & r(i,j)=D1 Hd(i,j)=(D1-r(i,j)/(D1-D0); end endendY=fft2(double(IA);Y=fftshift(Y);Ya=Y.*Hd;Ya=ifftshift(Ya);Ia=real(ifft2(Ya);figuresubplot(2,2,1),imshow(uint8(IA);s

51、ubplot(2,2,2),imshow(uint8(Ia);figuresurf(Hd,Facecolor,interp,Edgecolor,none,Facelighting,phong); 梯形高通濾波器IA=imread(lena.bmp);IB=imread(babarra.bmp);f1,f2=freqspace(size(IA),meshgrid);%D=100/size(IA,1);D0=0.1;D1=0.4;r=sqrt(f1.2+f2.2);Hd=ones(size(IA);Hd(r=D0 & r(i,j)=D1 Hd(i,j)=(D0-r(i,j)/(D0-D1); end endendY=fft2(double(