優(yōu)秀1數(shù)碼相機(jī)定位算法研究

1、2008 高教社杯大學(xué)生數(shù)學(xué)建模競賽編 號頁賽區(qū)評閱編號（由賽區(qū)評閱前進(jìn)行編號）：賽區(qū)評閱（可供賽區(qū)評閱時(shí)使用）：統(tǒng)一編號（由賽區(qū)送交前編號）：評閱編號（由評閱前進(jìn)行編號）：評閱人評分備注數(shù)碼相機(jī)算法研究摘要本文研究數(shù)碼相機(jī)中有標(biāo)定的相關(guān)問題。首先，本文建立了三個(gè)坐標(biāo)系：像素平面坐標(biāo)系、像物理平面坐標(biāo)系和相機(jī)坐標(biāo)系 。其中像素平面坐標(biāo)系和像物理平面坐標(biāo)系是同一個(gè)平面標(biāo)系是一個(gè)世界坐標(biāo)系，它以相機(jī)為參照物。不同需要而建立的；相機(jī)坐然后第一問確定圓心在像平面上的坐標(biāo)的問題，本文建立了兩個(gè)子模型：針孔相機(jī)模型和確定靶標(biāo)相對相機(jī)位置的模型，然后提出了運(yùn)用以上兩個(gè)子模型求解坐標(biāo)的方法。在第一個(gè)子模型針孔

2、相機(jī)模型中，本文對數(shù)碼相機(jī)進(jìn)行了適當(dāng)?shù)暮喕窗褦?shù)碼相機(jī)看成是一個(gè)針孔相機(jī)的結(jié)構(gòu)，利用射影幾何的有關(guān)知識建立了從相機(jī)坐標(biāo)到像物理坐標(biāo)的轉(zhuǎn)換關(guān)系模型。在第二個(gè)子模型確定靶標(biāo)相對相機(jī)位置的模型中，本文利用像平面上四個(gè)圖形公切線的交點(diǎn)建立了與靶標(biāo)平面的的模型。，并結(jié)合靶標(biāo)的、形狀，建立起了確定靶標(biāo)位置在建立了以上兩個(gè)子模型后，通過第二個(gè)子模型可以求出靶標(biāo)上圓心在相機(jī)坐標(biāo)系中的坐標(biāo)，再利用第一個(gè)子模型的轉(zhuǎn)換關(guān)系，就可以得到圓心在像平面上的坐標(biāo)。模型的求解，本文使用模擬退火算法計(jì)算出了像平面上四條公切線交點(diǎn)的坐標(biāo)，并使用基于最小二乘法的優(yōu)化工具箱的工具求解出靶標(biāo)的位置，進(jìn)一步求出了圓心在像平面上的坐標(biāo)，

3、 五個(gè)坐標(biāo)分別為： A0(-190.26,-196.77),B0 (-88.88,- 189.15),C0 (129.74,-172.72),D0 (72.85,119.30),E0 (-229.13,119.21)在模型的檢驗(yàn)?zāi)Ｐ椭?，本文分別討論了以上模型的精度和穩(wěn)定性。在精度檢驗(yàn)中，像平面上未被利用的圖形的輪廓上的點(diǎn)回靶標(biāo)平面上，并在靶標(biāo)平面上檢驗(yàn)這個(gè)輪廓是否與相應(yīng)的圓形重合。經(jīng)過檢驗(yàn)，輪廓上的點(diǎn)與相應(yīng)的圓形之間的平均偏差在 1 像素以內(nèi)，說明以上模型的精度很高。在隨后的穩(wěn)定性檢驗(yàn)中，我們通過計(jì)算機(jī)模擬的方式，隨機(jī)改變了像平面上圖形的輪廓，并對這些輪廓求解圓心，結(jié)果即使在輪廓損失了近 30

4、%的信息量時(shí)，圓心的平均偏移距離也只有 0.2682，不到一個(gè)像素，說明模型具有很好的穩(wěn)定性。最后，本文通過改變世界坐標(biāo)系，以靶標(biāo)作為參照物，給出了計(jì)算兩臺相機(jī)光學(xué)中心、像平面中心坐標(biāo)的方法，得出了兩臺相機(jī)相對位置的模型。： 系統(tǒng)標(biāo)定 射影幾何 針孔成像模型 模擬退火算法1一、 問題重述數(shù)碼相機(jī)在交通監(jiān)管（電子）等方面有廣泛的應(yīng)用。所謂數(shù)碼相機(jī)是指用數(shù)碼相機(jī)攝制物體的相片確定物體表面某些特征點(diǎn)的位置。最常用的方法是雙目，即用兩部相機(jī)來。對物體上一個(gè)特征點(diǎn)，用兩部固定于不同位置的相機(jī)攝得物體的像，分別獲得該點(diǎn)在兩部相機(jī)像平面上的坐標(biāo)。只要知道兩部相機(jī)精確的相對位置，就可用幾何的方法得到該特征點(diǎn)在

5、固定一部相機(jī)的坐標(biāo)系中的坐標(biāo)，即確定了特征點(diǎn)的位置。于是對雙目標(biāo)定。，精確地確定兩部相機(jī)的相對位置就是關(guān)鍵，這一過程稱標(biāo)定的一種做法是：在一塊平板上畫若干個(gè)點(diǎn)， 同時(shí)用這兩部相機(jī)照相，分別得到這些點(diǎn)在它們像平面上的像點(diǎn)，利用這兩組像點(diǎn)的幾何關(guān)系就可以得到這兩部相機(jī)的相對位置。然而，無論在物平面或像平面上我們都無法直接得到?jīng)]有幾何的“點(diǎn)”。實(shí)際的做法是在物平面上畫若干個(gè)圓（稱為靶標(biāo)），它們的圓心就是幾何的點(diǎn)了。而它們的像一般會變形，所以必須從靶標(biāo)上的這些圓的像中把圓心的像精確地找到，標(biāo)定就可實(shí)現(xiàn)。要求給出求解圓心的算法以及方法精度和穩(wěn)定性的檢驗(yàn)，最后討論兩部固定相機(jī)的相對位置關(guān)系。二、符號說明O

6、p v xr OyR Ku Or yr x MtAB（ 粗體）像素平面坐標(biāo)系原點(diǎn) 像素平面坐標(biāo)系縱軸 像物理平面坐標(biāo)系橫軸相機(jī)坐標(biāo)系原點(diǎn)相機(jī)坐標(biāo)系縱軸旋轉(zhuǎn)矩陣相機(jī)標(biāo)定矩陣像素平面坐標(biāo)系橫軸 像物理平面坐標(biāo)系原點(diǎn)像物理平面坐標(biāo)系縱軸相機(jī)坐標(biāo)系橫軸投影矩陣平移矢量向量 AB，粗體表示向量三、問題分析題目第一問要求建立確定靶標(biāo)上圓的圓心在像平面的坐標(biāo)的模型。由于存在畸變， 所以像平面的圖形并不是規(guī)則的形狀，單純從像平面的圖形出發(fā)顯然不能找到這個(gè)圓心在像平面的坐標(biāo)。題目中給出了靶標(biāo)的等信息，但是靶標(biāo)像上的點(diǎn)卻很難與靶標(biāo)上的點(diǎn)建立起對應(yīng)關(guān)系。因此我們需要將模型簡化，簡化成一個(gè)小孔成像的模型，這樣就可以很

7、容易的將五個(gè)圓的外公切線的交點(diǎn)對應(yīng)起來。對應(yīng)起這四個(gè)點(diǎn)之后，配合靶標(biāo)的2、形狀等信息，我們就可以確定靶標(biāo)在三中與相機(jī)的位置關(guān)系。確定了靶標(biāo)與相機(jī)的位置關(guān)系后，就可以很容易的將靶標(biāo)上的點(diǎn)投射到靶標(biāo)像平面上，當(dāng)然也包括五個(gè)圓心。在求取出五個(gè)圓心的空間坐標(biāo)之后，將其投射到像平面上， 就得到了第一問需要求的坐標(biāo)。在求解時(shí)涉及到一個(gè)問題，就是像平面上怎樣確定切線。因?yàn)橄衿矫嫔系膱D形是不規(guī)則的，所以很難確定這些形狀的切線。因此我們考慮另外的方法，使用搜索的辦法， 利用模擬退火算法求解。在如何檢驗(yàn)?zāi)Ｐ偷膯栴}上，需要分兩方面進(jìn)行檢驗(yàn)，一是精度，而是穩(wěn)定性。按照以上的方法求圓心在像平面上的坐標(biāo)，并沒有充分利用像

8、平面上所有輪廓點(diǎn)的信息，因此可以利用這些點(diǎn)來檢驗(yàn)?zāi)Ｐ偷木取τ诜€(wěn)定性問題，可以采用計(jì)算機(jī)模擬的方法，隨機(jī)修改圖形的輪廓，并用以上的方法再次進(jìn)行求解，通過比較修改前后的結(jié)果來分析模型的穩(wěn)定性。最后，考慮另外一臺相機(jī)的相對位置問題。根據(jù)前面模型，我們應(yīng)能夠?qū)θ我庖慌_相機(jī)確定靶標(biāo)相對它的位置，因此可以以這個(gè)靶標(biāo)作為參照物，建立一個(gè)世界坐標(biāo)系，將這兩臺相機(jī)的位置在這個(gè)坐標(biāo)系里面表示出來，以此確定兩臺相機(jī)的相對位置。四、模型假設(shè)1、2、3、4、假設(shè)靶標(biāo)像的中心恰好在光軸上假設(shè)數(shù)碼相機(jī)中圖像平面與光軸垂直假設(shè)相機(jī)兩個(gè)方向上焦距相等假設(shè)透鏡的焦距很小，像距約等于焦距五、模型準(zhǔn)備(一) 靶圖像矩陣表示首先將

9、題目中的圖片保存出來，得到的圖像可以很方便的放到里面進(jìn)行處理。但在處理之前還要進(jìn)行進(jìn)一步：(1)(2)，使用 imread()函數(shù)將文件讀入將矩陣變?yōu)?0-1 矩陣對于用以上方式得到的矩陣，有兩個(gè)值：0、15。其中 0 代表像素為白色的點(diǎn)，15 代表像素為黑色的點(diǎn)。為了方便下面處理，對需要把以上像素為 15 的點(diǎn)值全部變?yōu)?1。以上兩步的源代碼見附錄一。3(二) 圖像輪廓的提取在提取圖像輪廓時(shí)，首先要引入計(jì)算機(jī)圖像處理技術(shù)中四鄰域的概念。四鄰域：某個(gè)像素的上、下、左、右四個(gè)像素成為該像素的四鄰域。如下圖所示：圖 1 四鄰域示意圖則 A 的四鄰域?yàn)?B、C、D、E 四個(gè)像素。在尋找邊界時(shí)，對任意

10、一個(gè)值為 1 的像素，只要其四鄰域有一個(gè)值為 0，則認(rèn)為這個(gè)像素為邊界上的像素。此外，經(jīng)過這樣處理后還要剔除孤立的點(diǎn)。方法是若一個(gè)像素的周圍（包括像素四個(gè)角上的像素）值都為 0，則此點(diǎn)為孤立點(diǎn)，予以剔除。提取圖像輪廓的程序見附錄二。提取出的輪廓見下圖：圖 2提取出的輪廓示意圖4FBGCADHEI(三) 坐標(biāo)定義1.像素平面坐標(biāo)系像素平面坐標(biāo)即像素的行、列所標(biāo)識的坐標(biāo)，如下圖所示：123456nOp 1u23456mv圖 3像素平面坐標(biāo)示意圖其中每個(gè)像素表示為矩陣的行、列數(shù)，如(u,v)=(1,1)代表第一行第一列的像素。2.像物理平面坐標(biāo)系像物理平面坐標(biāo)系是以圖像的中心為坐標(biāo)原點(diǎn)的坐標(biāo)系，如下

11、圖所示：5OpuOrxrvyr圖 4 物理平面坐標(biāo)示意圖圖中 Op-uv 坐標(biāo)系為像素平面坐標(biāo)系，Or-xryr 為物理平面坐標(biāo)系。兩者的轉(zhuǎn)換關(guān)系為1：é 1u ùs1êd0 úéuùúéxr ùêxêvú = êv úêy ú0（1）ê úêë1úû0 úê r úêêêëdúê

12、5; 1 úûy001 úû其中 dx 和 dy 為像素點(diǎn)在 x 軸與 y 軸方向上的物理(在這里使用像素為)，u0和 v0 為點(diǎn) Or 在像素坐標(biāo)系中的坐標(biāo)，s 為圖像傾斜度，對于數(shù)碼相機(jī)，可以認(rèn)為圖像傾斜度為一個(gè)極小量，一般為 0。1再次，我們選擇 dx 和 dy 分別為 1，s=0，由于行、列像素?cái)?shù)均為偶數(shù)，因此選擇u0=512.5，v0=384.5，此時(shí)，要計(jì)算 u、v 時(shí)，需要代入像素的中心的坐標(biāo)來計(jì)算。3.相機(jī)坐標(biāo)系界的坐標(biāo)，現(xiàn)根據(jù)假設(shè) 1、假設(shè) 2 做如下坐標(biāo)：以相機(jī)光學(xué)相機(jī)坐標(biāo)即整個(gè)中心為原點(diǎn) O，O 與像的中心 Or 連成的直線為 z

13、軸，以過 Or 平行于像平面 Or-xryr 且平行于 Orxr 的直線為 x 軸，以過 Or 平行于像平面 Or-xryr 且平行于 Oryr 的直線為 y 軸，建立直角坐標(biāo)系，如下圖所示：6zOrOxryrxy圖 5 相機(jī)坐標(biāo)系示意圖(四) 分別提取五個(gè)圓輪廓的像素坐標(biāo)為了便于下面求解，還要將五個(gè)圓輪廓的坐標(biāo)分別提取出來。由于在上一個(gè)處理中 ， 已經(jīng)將圖像中可能出現(xiàn)的孤立的點(diǎn)除去，因此現(xiàn)在圖像中只存在輪廓上的點(diǎn)，這意味著圖像中每一個(gè)像素值為 1 的點(diǎn)其四周至少有一個(gè)像素的值為 1。因此我們可以使用如下方法提取五個(gè)圓輪廓的像素坐標(biāo)：i.令 k=1ii. iii. iv.v.令 i=1，j=1

14、(i,j)像素值是否為 1，若不為 1，進(jìn)入第 v 步，否則，進(jìn)入第 iv 步令 i=i+1,j=j+1，跳至第 ii 步。若 i=1024 且 j=768，終止程序。將(i,j)至 Bk，并(i-1,j)、(i+1,j)、(i,j+1)、(i,j-1)、(i-1,j-1)、(i-1,j+1)、(i+1,j-1)、(i+1,j+1)是否等于 1，并令(i,j)為以上檢驗(yàn)出等于 1 的數(shù)值，跳至vi.第 v 步(i,j)像素值是否為 1，若不為 1，將(i,j)至 Bk，進(jìn)入，否則，進(jìn)入vii. 令 k=k+1，對整幅圖像剔除孤立點(diǎn)，進(jìn)入第 ii 步。通過以上步驟即可得到五個(gè)圓的輪廓的坐標(biāo)（源程

15、序見附錄三）。7六、模型建立(一) 針孔相機(jī)模型的建立相機(jī)即一個(gè)從三維世界點(diǎn)的集合到二維世界點(diǎn)的集合的一個(gè)。針孔相機(jī)模型是一個(gè)沒有引入任何非線性畸變的，空間點(diǎn)P 通過一個(gè) 3×4 的矩陣到影像點(diǎn) p。1.射影幾何射影幾何是本文討論的基礎(chǔ)。射影幾何是遠(yuǎn)處元素的問題。幾何的擴(kuò)展，解決了幾何中無窮對于三維坐標(biāo)(X,Y,Z)，其齊次坐標(biāo)表示為(X1,X2,X3,X4)，其中 X=X1/X4，Y=X2/X4， Z=X3/X4。對于非零數(shù)，(X1,X2,X3,X4)與(X1,X2,X3,X4)表示同一個(gè)點(diǎn)。若X4=0，則該點(diǎn)位無窮遠(yuǎn)點(diǎn)。22.針孔相機(jī)模型空間點(diǎn) P 通過一個(gè) 3×4 的

16、投影矩陣到像平面上的點(diǎn) p：3p = MP投影矩陣 M 可以分為三部分：é- 10ù0- 1 000úéR tùM = KTG = Kê 00（2）úê 01úêêë 0- 10úûëû其中 R 為旋轉(zhuǎn)矩陣，t 為平移矢量。中間矩陣兩個(gè)負(fù)號代表通過小孔物體呈倒像。K 為相機(jī)標(biāo)定矩陣，它包含相機(jī)內(nèi)部參數(shù)數(shù)據(jù)：- fu cotqé fuu0 ùK = ê 0v úf / sinq（3）ê&#

17、234;ë 00 ú1 úûv0其中 5 個(gè)參數(shù)為：fu,fv：相機(jī)焦距，分別以像素為這里取 fu=fv=f=1577。沿水平方向和垂直方向測量。根據(jù)假設(shè) 3、4，u0,v0：相機(jī)主點(diǎn)，即 Or 點(diǎn)，根據(jù)假設(shè) 1，這里均為 0。：相機(jī)水平與垂直軸向間的夾角，根據(jù)假設(shè) 2，取為/2。則標(biāo)定矩陣可以改寫為：8é- f0 ù0- fK = ê 00 úêêë 0ú0- 1úû因此，相機(jī)坐標(biāo)(x,y,z)到像物理平面坐標(biāo)(xr,yr)就可以寫成如下形式：é

18、;xùé- fé x ù0- f000- 10ùúéR tùêyú1úêêúê úêx2 ú = ê000úê（4）1ûêzú0úë 0êx úêê1úë 3 ûëûë û其中：x= x2=1x, yrrxx33(二) 確定靶標(biāo)平

19、面相對相機(jī)位置的模型通過以上分析，界的點(diǎn)可以通過一次線性變換轉(zhuǎn)換成像物理平面坐標(biāo)上的點(diǎn)，然而這個(gè)過程并不是可逆的。三應(yīng)的關(guān)系，所以單純從像平面上的點(diǎn)上的點(diǎn)不能與二上的點(diǎn)建立對回三上的點(diǎn)是不可能的。盡管如此，如果我們一些參數(shù)（比如靶標(biāo)大?。屵@些參數(shù)變?yōu)橐阎?，或者從像平面上提取一些有用的信息，在這樣嚴(yán)格的條件下，仍然可以將像平面上的點(diǎn)近似還原到三中。本問題，我們可以從像平面中提取一些點(diǎn)，這些點(diǎn)可以很容易找到三中與之對應(yīng)的點(diǎn)，然后根據(jù)已知的靶標(biāo)的長度等數(shù)據(jù)，就可以近似求出靶標(biāo)相對于相機(jī)的位置。9i.在靶像上找到四個(gè)定點(diǎn)圓的外切線，如下圖所示，確定四個(gè)交點(diǎn)分別為 A、B、C、D。圖 6 切線示意圖

20、到 z=-1 平面上：ii.將四個(gè)點(diǎn)按照如下公式-1éxùé- f0 ùé x ù0- f1êyú = ê0 úêx ú0（5）ê úêúê2 úêëzúû0êë0-1úûêë- 1úû得到四個(gè)新的點(diǎn) A、B、C、D。10iii.得到的四個(gè)點(diǎn) A、B、C、D既為靶標(biāo)在 z=-1 平面上的射影，如下圖

21、所示：DDBBCC'OAA圖 7 靶標(biāo)在 z=-1 上的投影示意圖在靶標(biāo)上與這四個(gè)點(diǎn)相對應(yīng)的是A、B、C、D四個(gè)點(diǎn)，它們是 A、C、D、E四個(gè)圓外公切線的交點(diǎn)，并與 A、B、C、D一一對應(yīng)。由此我們可以設(shè) OA=1OA, OOB=2OB, OOC=3OC, OOD=4OD，進(jìn)而得到 AB、AC、BD、CD的表參數(shù)1,2,3,4：。這樣我們就可以建立以下方程求解ìïïïïl OB- lOA = l21l OC- lOA = l31l OD- l OB = lï42ïl OC- l OD = l（6）í3OA)

22、 × (OA) × (4ï (OB-OB-OC-OD-OA) = 0OB) = 0ïï (ïï(l3 OC-OA) × (OC-OD) = 0ï(l OC- l OD) × (l OD- l OB) = 0î3442其中 l 為靶標(biāo)上線段 AB的長度。需要注意的是雖然只有四個(gè)未知數(shù)，但是方程仍然使用了 8 個(gè)方程，主要是這個(gè)模型沒有考慮畸變的因素，實(shí)際的數(shù)據(jù)帶入并不是完全吻合的，所以要在這里求解一個(gè)偏差最小的解，所以這 8 個(gè)方程都不能舍棄。利用以上方程組求出 4 個(gè)參數(shù)之后，就可以確

23、定靶標(biāo)平面上四個(gè)點(diǎn) A、B、C、11D的坐標(biāo)，靶標(biāo)平面在相機(jī)坐標(biāo)中的位置也隨之確定。(三) 尋找圓心的方法經(jīng)過如上處理后，我們就得到了靶標(biāo)在界中相對于相機(jī)的位置。這時(shí)只要將靶標(biāo)平面的圓心到相平面上即可找到靶標(biāo)平面上圓心的位置。按照針孔相機(jī)透視模型的法則進(jìn)行，即：éxùé- fé x ù0- f000ùúéR tùêyú1úêêúê úêx2 ú = ê0000úê（7）1

24、51;êzú- 10úë 0êx úêê1úë 3 ûëûë û由于已經(jīng)在以上的模型中求出了靶標(biāo)所在的相對位置，因此其中的R、t 分別變?yōu)榫仃嚭土阆蛄浚琭 已知，利用這些條件就可以求出像物理平面上圓心的坐標(biāo)。假設(shè)靶標(biāo)上的圓心坐標(biāo)為 Oo(xo,yo,zo)，則圓心的坐標(biāo)計(jì)算公式為：0ùéxo ùé x ùé- f0- f000êy ú1êx ú =

25、 ê0úê o ú00（8）ê2 úêúêzúêx úê- 10úêo úë 3 ûëûë 1 ûìX =x1ïxí3ïY = x2ïîx3(X,Y)既為要求的圓心坐標(biāo)。七、模型的求解(一) 像物理坐標(biāo)系中切線的求解給出的像物理坐標(biāo)系的圖像并不是完美的圓形甚至不是橢圓形，給求解公切線帶來了一定的難度。并且計(jì)算機(jī)中的圖形

26、表示是離散形式的，這就更加難以處理。為此，我們引入了模擬退火算法（SA）進(jìn)行求解。模擬退火算法是一種搜索算法，相比于遺傳算法、蟻群算法等算法，這種算法比較容易陷入局部最優(yōu)。但高效而又準(zhǔn)確的算法。本問題，不存在局部最優(yōu)問題，因此模擬退火算法仍是一個(gè)121.求解思路求解切線方程并不是最優(yōu)問題，但是我們可以通過一定方法將其轉(zhuǎn)化成一個(gè)最優(yōu)問題。首先是約束條件。我們所研究的圖形的切線跟圖形理論上應(yīng)該只有一個(gè)交點(diǎn)，但是在計(jì)算機(jī)中，這個(gè)形狀的邊界點(diǎn)是離散的而非連續(xù)的，所以很難保證切線與圖形只有一個(gè)交點(diǎn)，即使在離散的坐標(biāo)中的確只有一個(gè)交點(diǎn)。因此我們這里并不希望切線與圖形之間有交點(diǎn)，而是圖形與切線的最近的點(diǎn)是相

27、鄰的。這樣最大的誤差就圍以內(nèi)。比如在下圖中：在一個(gè)像素范圖 8 切線示意圖灰色代表切線，黑色代表圖形輪廓。如果切線與形狀相交，只會造成更大的誤差。然后是目標(biāo)函數(shù)。由于約束條件選擇了讓切線不與圖形相交，因此目標(biāo)函數(shù)要讓切線盡量向圖形靠近。我們在這里選取四條切線圍成的面積最小為目標(biāo)函數(shù)。在計(jì)算四邊形面積時(shí)，將四邊形分為兩個(gè)三角形，用以下公式計(jì)算三角形面積：ab df111S = 1 c2 e其中(a,b)、(c,d)、(e,f)為三角形三個(gè)頂點(diǎn)的坐標(biāo)。此外，用于求解的參數(shù)也不直接用四條切線的方程系數(shù)，而是用四個(gè)交點(diǎn)作為參數(shù)進(jìn)行求解。四個(gè)交點(diǎn)坐標(biāo)可以表示四條切線，并且減少了參數(shù)數(shù)目，省去了計(jì)算交點(diǎn)坐

28、標(biāo)的步驟。132.初始解的給出在這里給出初始解只是為了加快算法收斂的速度，應(yīng)用中也可以隨機(jī)取四個(gè)可行的點(diǎn)，此做法不失一般性。為了給出初始解，我們現(xiàn)在原圖中畫四條靠近圖形但是不與圖形相交的直線，然后找到四個(gè)交點(diǎn)的坐標(biāo)，如下圖所示：圖 9 初始解尋找由此得到的四個(gè)初始解（像素坐標(biāo)）為：(256,139),(235,546),(690,172),(617,544)。143.模擬退火算法過程i.ii. iii.設(shè)定初始溫度 T0，最終溫度 T1，擾動最大最小范圍 Lmax,Lmin，初始化迭代次數(shù) cs令 T=T0 計(jì)算初始解 r0 的目標(biāo)函數(shù)值 E0產(chǎn)生一個(gè)0,L內(nèi)的隨機(jī)擾動加至 r0 上，加入擾動

29、后是否在可行域上，若否，則繼續(xù)產(chǎn)生可行域，直到加入擾動后在可行域上，加入擾動后的計(jì)算 r1 的目標(biāo)函數(shù)值 E1計(jì)算E=E1-E0若E<0，則令 r0=r1,E0=E1 跳至第 viii 步，否則，跳至第 iii 步產(chǎn)生一個(gè)隨機(jī)數(shù) u，若 u<exp(E/T)，則令 r0=r1,E0=E1r1iv.v.vi. vii. viii.迭代次數(shù) k 是否為一個(gè)整數(shù)常數(shù) K 的整數(shù)倍，若是，則令 T=aT，a 為一個(gè)常數(shù)若 T<T1，則進(jìn)入 x，否則令 k=k+1，跳至第 iii 步結(jié)束ix.x.為了防止程序在最后收斂可行域邊緣產(chǎn)生隨機(jī)擾動時(shí)經(jīng)常超出可行域，我們在一開始給定了一個(gè)最大最

30、小范圍 Lmax,Lmin，在每次迭代時(shí)，L 都做如下計(jì)算：- LLmax - LminL= LmaxminT + L- TT - TminminT - T1010我們選取 T0=100,T1=0.01,Lmax=10,Lmin=1,k=10,a=0.9 進(jìn)行計(jì)算（源程序見附錄四），得到的結(jié)果為：A(285,143),B(240,539),C(685,177),D(613,537)，如下圖所示：圖 10 求解切線結(jié)果圖（放大后）15(二) 圓心坐標(biāo)的計(jì)算1.靶標(biāo)平面位置的確定前面的模型建立中已經(jīng)給出了確定靶標(biāo)平面的方程組的模型。由于像存在畸變，因此方程組很有可能無解。但是可以根據(jù)這個(gè)方程組求一

31、個(gè)近似解，這個(gè)近似解也就是我們要取的近似的靶標(biāo)平面。具體求法我們使用里面的 fsolve()函數(shù)。fsolve 函數(shù)是軟件優(yōu)化工具箱中用于求解非線性方程組的函數(shù)，其算法基于最小二乘法，可用于求解非線性方程組的零點(diǎn)。如果一個(gè)系統(tǒng)能表示成一個(gè)非線性方程組，便可用 fsolve 函數(shù)求零解。4在使用 fsolve 函數(shù)時(shí)，初始解的選取非常重要。若初始解不合理，則fsolve 函數(shù)很容易陷入局部最優(yōu)或者提前結(jié)束迭代。我們可以根據(jù)像的大小和像距來大體估計(jì)靶標(biāo)平面距離相機(jī)的距離：在像中，四邊形的長大約為 400 像素，即大約有 105mm。而靶標(biāo)的真實(shí)長度有100mm，因此可以得出結(jié)論，靶標(biāo)距離相機(jī)坐標(biāo)

32、xoy 平面與像距離相機(jī)坐標(biāo) xoy 平面的距離大體相等。而 O 點(diǎn)到Or 點(diǎn)有 1577 像素，因此不妨令初始值為1500,1500,1500,1500。如此，使用 fsolve 函數(shù)求解（見附錄五），求得：1,2,3,4=1657.63,1815.08,1829.54,1984.68同時(shí)計(jì)算所得的四條邊平均長度誤差為 3.19 像素，相對誤差 0.844%，誤差比較小， 可以接受。因此靶標(biāo)上四條外公切線交點(diǎn) A、B、C、D的坐標(biāo)分別為：(239.13,253.85,-1657.63),(313.64,-177.82,-1815.08),(-200.12,240.73,-1829.54),(

33、-126.48,-191.92,-1984.68)。用以上四個(gè)點(diǎn)中的任意三個(gè)就可以表示靶標(biāo)所在的平面。2.圓心的確定對于下圖所示的靶標(biāo)平面：16圖 11 靶標(biāo)平面如果以 AB和 AC作為平面的基向量，則五個(gè)圓的圓心分別可以表示為：1212A' '' A=A'''C''' +A''' B'' '1241244212A' '' B=A'''C''' +A''' B''

34、'124124A'' 'C = 112 A'''C''' +A''' B'''12（9）124124A'''D = 112 A'''C''' + 112 A' ''B ' '124124A'' 'E =A'''C''' + 112 A'''B''&

35、#39;12124124中，上述表依然成立。我們已經(jīng)知道了 A,B,C三點(diǎn)在相同樣，在三機(jī)坐標(biāo)系中的坐標(biāo)，因此很容易就可以將 A、B、C、D、E 的坐標(biāo)表示出來。然后通過真空相機(jī)模型的轉(zhuǎn)換矩陣就可以很容易的找到像物理平面上對應(yīng)的圓心的坐標(biāo)。經(jīng)過計(jì)算機(jī)程序計(jì)算（見附錄五），得到的圓心坐標(biāo)為：（像物理坐標(biāo)系）A0(-190.26,-196.77),B0 (-88.88,-189.15),C0 (129.74,-172.72)D0 (72.85,119.30),E0 (-229.13,119.21)（像素坐標(biāo)系）A0 (322,188),B0 (424,195),C0 (642,212),D0(58

36、5,504),E0 (283,504)圓心位置如下圖所示：17圖 12 圓心位置八、模型的精度分析與檢驗(yàn)(一) 誤差產(chǎn)生的分析1.畸變產(chǎn)生的誤差以上本文討論的是線性變換的模型，并沒有考慮圖形產(chǎn)生的畸變。實(shí)際的相機(jī)拍的照都是有畸變的，這些畸變的產(chǎn)生一部分是由于透鏡成像的，另一部分是由于空氣等的折射因素，此外還有極其復(fù)雜的其他免會對精度產(chǎn)生影響。以上模型建立在沒有畸變的基礎(chǔ)上，難2.連續(xù)圖形離散化產(chǎn)生的誤差被拍攝的物體一般是連續(xù)的，然而物體拍攝到數(shù)碼相機(jī)后，原來連續(xù)的點(diǎn)被離散化 ， 會丟失很多圖像邊緣的信息。在以上模型中，求解切線就嚴(yán)重依賴于圖像的邊緣。而這些圖像邊緣丟失的信息會影響切線的形狀，進(jìn)

37、而影響求得的靶標(biāo)的相隨位置以及最后求得的圓心。183.邊緣圖像顏色不確定帶來的誤差這種誤差不同于以上離散化的誤差。以上離散化產(chǎn)生的誤差是指續(xù)的線離散化時(shí)需要丟棄一些點(diǎn)，以至于形狀改變。而邊緣化圖像顏色不確定產(chǎn)生的主要是拍攝時(shí)對焦不精確等產(chǎn)生的。實(shí)際拍攝中，對焦不可能完全精確，因此在圖像中顏色變化之間的點(diǎn)的顏色經(jīng)常會顯示出一些過渡顏色，比如在黑色與白色之間，并不是黑白分明，而是會有一些像灰色。在以上我們的處理中，有一在 Photoshop 中用“閾值”工具將圖像處理為黑白顏色。然而這往往是不精確的。比如我們選擇閾值為 205，當(dāng)選擇 206 時(shí)，又有一些像素將被處理成白色，依次類推。閾值無法準(zhǔn)確

38、確定，并且也不存在精確的閾值，因此這種操作會給圖像邊緣帶來一定誤差。4.求解誤差以上所講所有導(dǎo)致誤差的都是客觀條件上產(chǎn)生的，而求解誤差則存在于模型本身當(dāng)中。比如我們以上的模型是建立在“光線垂直攝入相機(jī)的光學(xué)中心，出剛好射在圖像的中心上”等一系列假設(shè)條件之上的，這些假設(shè)不一定符合事實(shí)。在大多數(shù)情況下，這些假設(shè)總是近似符合事實(shí)的（否則模型就會出現(xiàn)不符合實(shí)際的錯(cuò)誤，而不是誤差 ）， 這里的“近似”符合事實(shí)也會產(chǎn)生誤差。(二) 精度檢驗(yàn)的基本思路在以上模型的求解過程中，我們并沒有利用圖形中給出的全部信息就得到了結(jié)果， 圖像的邊緣并沒有被完全利用，這些沒有被完全利用的信息就是我們建立精度檢驗(yàn)?zāi)Ｐ偷幕A(chǔ)。

39、我們以上求解圓心坐標(biāo)時(shí)，是直接以三內(nèi)靶標(biāo)上的圓心投射到像上，并以此作為圓心。這里面暗含了一個(gè)假設(shè)，即像上的點(diǎn)都是由我們求得的靶標(biāo)位置上的圓的點(diǎn)形成的。但是求得的靶標(biāo)的位置是不精確的，因此三內(nèi)圓上的點(diǎn)與像上圓上的點(diǎn)可能會有偏差。在理想情況下，若不存在這種偏差，則我們求得的圓心就是要求的精確的圓心；這種偏差越大，圓心的位置就越不可靠。因此精度檢驗(yàn)可以轉(zhuǎn)換為像上的圖形與靶標(biāo)上的圓的關(guān)系的檢驗(yàn)。具體方法是，首先將像上的輪廓的點(diǎn)投射到以上模型求得的靶標(biāo)上，然后計(jì)算靶標(biāo)上相應(yīng)的圓心距離這些點(diǎn)的平均距離和標(biāo)準(zhǔn)差，平均距離越貼近于靶標(biāo)上圓的半徑，則模型越準(zhǔn)確；標(biāo)準(zhǔn)差越小，說明模型越準(zhǔn)確。(三) 檢驗(yàn)算法及求解

40、1.從像物理平面坐標(biāo)到相機(jī)坐標(biāo)下面我們考慮一個(gè)在像上的點(diǎn)怎樣如下圖所示：到求出的靶標(biāo)平面上：19AAOA圖 13示意圖現(xiàn)假設(shè)像平面上任意一點(diǎn) A 坐標(biāo)為(a,b,c)，則 AO=(-a,-b,-c)，若假設(shè) OA=AO， 那么，OA=-(a,b,c)。現(xiàn)在我們?nèi)“袠?biāo)平面上的點(diǎn) A(m,n,p)，同時(shí) AB=(ux,uy,uz),AC=(vx,vy,vz)。則 OA 還可以表示為 OA=OA+ 1AB+ 2AC=(m,m,p)+ 1(ux,uy,uz)+ 2(vx,vy,vz)，用矩陣可以表示為：é abc ùl ll êuu ú = -m npê

41、;uxêëvxyz ú12vz úûvy要求解，只要計(jì)算如下式：c ù-1é ap êuxêëvxbêu úl ll = -m nu（10）yz ú12vz úûvy即可。這樣，就可以得到 A的坐標(biāo)為(-a,-b,-c)。2.計(jì)算點(diǎn)與圓心距離的均值與方差到靶標(biāo)平面上的點(diǎn)位A，預(yù)期對應(yīng)的圓心為 O，則假設(shè)計(jì)算兩個(gè)點(diǎn)之間的距離|OA|。對于原邊緣的整個(gè)序列上的點(diǎn)，計(jì)算所有的距離，并計(jì)算出所有這些 點(diǎn)距離相應(yīng)圓心的距離的均值與標(biāo)準(zhǔn)差。經(jīng)過計(jì)算（源程序見

42、附錄六），得到如下的結(jié)果（以像素為）：20表一 精度檢驗(yàn)結(jié)果從以上表可以看出，以|OA|的平均偏差都不大于一個(gè)像素，標(biāo)準(zhǔn)差也在 2 像素左右標(biāo)準(zhǔn)差系數(shù)非常小，說明像平面上的點(diǎn)到靶標(biāo)平面上，得到的點(diǎn)比較好的擬合了原靶標(biāo)平面上的五個(gè)圓，因此有理由相信，經(jīng)過以上算法計(jì)算得到的像平面上的圓心精度很高。九、模型的穩(wěn)定性模擬與檢驗(yàn)由于種種，像的輪廓上的點(diǎn)很多情況下是難以確定的，而這些難以確定的點(diǎn)常常會影響模型求解的結(jié)果。一個(gè)好的模型不僅要精確，而且在輪廓出現(xiàn)變化時(shí)，模型的結(jié)果不能變化太大，即要求模型的穩(wěn)定性要好。為了檢驗(yàn)本模型的穩(wěn)定性，我們采用計(jì)算機(jī)模擬的方法來破壞像的輪廓，再采用以上建立的模型求解，比較

43、兩個(gè)結(jié)果就可以模型穩(wěn)定性的好壞。為了簡便起見，我們只考慮輪廓縮小的情況，不考慮輪廓放大的情況。對于輪廓縮小的情況，我們可以隨機(jī)減少輪廓上的點(diǎn)來模擬圖像輪廓的變化。在輪廓上減少點(diǎn)，利用模擬退火算法求解外公切線交點(diǎn)坐標(biāo)時(shí)，其最優(yōu)解目標(biāo)函數(shù)值只會比減少錢更小，因此不妨令模擬退火算法的初始點(diǎn)為以上我們模型確定的4 個(gè)點(diǎn) 。令 p 為輪廓上點(diǎn)的保留概率，則輪廓上的點(diǎn)被刪除的概率為 1-p，我們在 3 個(gè)不同的 p 值下利用以上模型計(jì)算圓心的坐標(biāo)。每個(gè) p 值做了 10 次模擬（模擬程序見附錄七）， 模擬的結(jié)果如下：表二 穩(wěn)定性檢驗(yàn)表其中每個(gè)點(diǎn)的偏移距離即新的圓心與前面模型求解出的圓心之間的距離；總偏移距

44、21p總偏移距離平均偏移距離單個(gè)點(diǎn)最大偏移距離0.910.0210.850.110.713.4140.26821.414圓心偏差偏差/半徑標(biāo)準(zhǔn)差系數(shù)(/)A45.980.620.01371.900.0413B46.030.670.01491.630.0354C46.350.990.02192.260.0487D46.280.920.02032.490.0538E45.850.490.01081.900.0417離即 10 在一個(gè) p 下 10 次模擬每個(gè)點(diǎn)偏移距離的綜合；每個(gè)p 下做 10 次模擬，每次模擬有 5 個(gè)坐標(biāo)，因此每個(gè) p 下有 50 個(gè)點(diǎn)，平均偏移距離既為總偏移距離除以 50。表

45、中各數(shù)據(jù)為像素。從以上數(shù)據(jù)可以看出，在輪廓發(fā)生改變后，即使在輪廓損失了近 30%的信息量時(shí)， 圓心的平均偏移距離也只有 0.2682，不到一個(gè)像素，因此該求解圓心的模型具有很好的穩(wěn)定性。十、兩部固定相機(jī)相對位置的確定為了確定兩步固定相機(jī)的相對位置，需要先建立一個(gè)相機(jī)之外的世界坐標(biāo)系。在這里不妨以靶標(biāo)上的點(diǎn) A 為原點(diǎn)，AC 的x 軸正方向，AE 的y 軸正方向，以 AC 與 AE 的叉積的z 軸正方向建立直角坐標(biāo)系。建立了世界坐標(biāo)系后，就可以確定兩個(gè)相機(jī)的光學(xué)中心 O1、O2 以及像平面中心 Or1、Or2 的坐標(biāo)。利用這四個(gè)點(diǎn)的坐標(biāo)，就可以確定兩臺相機(jī)的相對位置。利用問題 1 中的模型，兩個(gè)

46、相機(jī)都可以根據(jù)各自在像平面上呈的像來確定靶標(biāo)平面與的位置關(guān)系。利用這個(gè)關(guān)系就可以求出四個(gè)點(diǎn)在世界坐標(biāo)系中的坐標(biāo)。下面將給出 O1 的世界坐標(biāo)系的求法，其他三個(gè)點(diǎn)方法類似。如下圖所示，回到針孔相機(jī)模型的世界坐標(biāo)系中，O1 仍為坐標(biāo)原點(diǎn)，A,B,C,D,E 各點(diǎn)的坐標(biāo)均為已知：DEO1CO1A圖 14 確定 O1 點(diǎn)坐標(biāo)示意圖過 O1 做靶標(biāo)平面的垂線 O1O1，并與靶標(biāo)平面相交于 O1。設(shè) O1O1的長度為 h，根據(jù)體積守恒，可以根據(jù)下式計(jì)算 h：221/6×|(O1A×O1C)·O1E|=1/2×h×|AE|×|AC|根據(jù) O1O1的

47、長度 h 以及 O1O1與 AE,AC 的垂直關(guān)系，求 O1的坐標(biāo)。得到坐標(biāo)后， 計(jì)算 a 和 b，使得：AO1=aAC+bAE。最后，返回到相機(jī)之外的以 A 為原點(diǎn)的世界坐標(biāo)系中，O1 的坐標(biāo)就可以表示為(a,b,h)。利用上述方法，得到四個(gè)點(diǎn) O1、O2、Or1、Or2 在世界坐標(biāo)系中的坐標(biāo)，就可以得到兩部相機(jī)的位置關(guān)系。十一、模型擴(kuò)展以上我們建立的是線性的模型，然而實(shí)際上，由于鏡頭畸變、大氣折射等因素的影響，實(shí)際成像過程并不滿足共線條件，而是一種非線性幾何關(guān)系。因此可以考慮建立一個(gè)多項(xiàng)式模型來綜合考慮以上因素的影響。以下我們建立一個(gè)有理函數(shù)模型：在有理函數(shù)模型中，像點(diǎn)坐標(biāo)(r,c)表達(dá)為

48、以相應(yīng)地面點(diǎn)坐標(biāo)(X,Y,Z)為自變量的多項(xiàng)式的比值，基本方程為：p(Xn,Yn,Zn)1ìïrn =p (Xn,Yn,Zn)2ïíï（11）p(,)X Y Znnnïcn=3p (,)X Y Zïînnn4式中rn和cn 分別是標(biāo)準(zhǔn)化的像素坐標(biāo)， Xn ，Yn 和Zn 分別是標(biāo)準(zhǔn)化的地面坐標(biāo)。所謂標(biāo)準(zhǔn)化是指 2 個(gè)像點(diǎn)坐標(biāo)和 3 個(gè)地面點(diǎn)坐標(biāo)分別平移和縮放使它們的值落在- 1.0,1.0中。一般地，多項(xiàng)式次數(shù)最高取 3，地面點(diǎn)用三維坐標(biāo)表示，則每個(gè)多項(xiàng)式有以下形式：m1 m2 m3p = åå

49、åa* Xi *Yj * Zk（12）ijki=0 j=0 k=0其中aijk 為有理函數(shù)多項(xiàng)式系數(shù)， 0 £ m1 ,m2 ,m3 £ 3 ， (m1 + m2 + m3) £ 3 。當(dāng)函數(shù) p= 1 時(shí)，上式變?yōu)橐话愣囗?xiàng)式。光學(xué)投影引起的畸變表示為一階多項(xiàng)式，而鏡p=24頭畸變、大氣折射和地面曲率等變形可由多項(xiàng)式來描述，高階部分和其它未知畸變可用三階多項(xiàng)式來趨近。(1)式是有理函數(shù)模型的正解形式，同樣通過逆變換可得到Xn ，Yn 的反解形式。答解有理函數(shù)模型的實(shí)質(zhì)就是確定多項(xiàng)式的系數(shù)，而系數(shù)的數(shù)目隨多項(xiàng)式階數(shù)、分母多項(xiàng)式的不同組合而變化，共有 9 種

50、情況，如下表：23表三 系數(shù)數(shù)目決定由于(1)式為非線性方程，可對其進(jìn)行線性化處理建立其誤差方程：ìp(Xn,Yn,Zn)1p (Xn,Yn,Zn)2ïvr =- rïíï（13）p(,)X Y Znnnïvc=- c3p (,)X Y Zïînnn4然后用最小二乘法求得相應(yīng)多項(xiàng)式系數(shù)。解算過程如下所示：像方格網(wǎng)模型變換物方格網(wǎng)圖 15 解算過程24創(chuàng)建像格模型變換生成 物 格解算參數(shù)多項(xiàng)式階數(shù)分母組合多項(xiàng)式系數(shù)數(shù)目至少需點(diǎn)個(gè)數(shù)3p2 ¹ p478393p2=p459303p2=p4=140203p2 &

51、#185; p438192p2=p429152p2=p4=120102p2 ¹ p41471p2=p41161p2=p4=184十二、模型評價(jià)優(yōu)點(diǎn)：1. 用靶標(biāo)像的四條公切線的交點(diǎn)對應(yīng)于靶標(biāo)上的圓的公切線的交點(diǎn)，由此確定靶標(biāo)平面在相機(jī)坐標(biāo)中的位置，用模擬退火算法和最小二乘法求解，誤差較??；2. 精度檢驗(yàn)的模型中充分考慮了靶標(biāo)像邊界的信息，從整體上把握了圓的作用；3. 穩(wěn)定性檢驗(yàn)中，將靶標(biāo)像上的離散點(diǎn)分別以三個(gè)不同的概率隨機(jī)去掉部分點(diǎn)，其異概率性和隨機(jī)性增強(qiáng)了檢驗(yàn)的可靠性和說服力。缺點(diǎn)：1. 本文模型沒有考慮相機(jī)拍攝時(shí)的畸變情況，使得結(jié)果存在一定的誤差。十三、參1月 1 日。234，基

52、于編隊(duì)遙感圖像的對地算法研究，國防科技大學(xué)，2004 年 10，高等幾何，：大學(xué)，2005 年。，基于多像滅點(diǎn)進(jìn)行相機(jī)標(biāo)定的方法研究，大學(xué)，2004 年 11 月 1 日。戰(zhàn)，基于的非線性方程組求解的方法，科技資訊，2008年 5 月 13 日。十四、附錄附錄一原圖像數(shù)據(jù)在中的處理%data.m a=imread('image3.bmp','bmp'); dat=zeros(length(a(:,1),length(a(1,:); dat=dat+1;for i=1:length(a(:,1)for j=1:length(a(1,:) if a(i,j)=15d

53、at(i,j)=0;endendendclear a;clear i;clear j;25附錄二 提取輪廓源程序%boundary.m%提取邊界data;bound=zeros(length(dat(:,1),length(dat(1,:); for i=1:length(dat(:,1)for j=1:length(dat(1,:) if dat(i,j)=1if dat(i-1,j)=1 & dat(i+1,j)=1 & dat(i,j-1)=1 & dat(i,j+1)=1 bound(i,j)=0;elsebound(i,j)=1;endelsebound(i,j)=0;endendend%檢測孤立點(diǎn)for i=1:length(dat(:,1)for j=1:length(dat(1,:) if dat(i,j)=1if dat(i-1,j)=0 & dat(i+1,j)=0 & dat(i,j-1)=0 1,j+1