基于卡爾曼預(yù)測的序列圖像橢圓檢測算法

趙江濤，張東亮，張鎖平，齊占輝

（國家海洋技術(shù)中心 天津 300112）

圖像平面上點(diǎn)、直線以及二次曲線的檢測是圖像處理、模式識別及機(jī)器視覺等領(lǐng)域的重要研究內(nèi)容。圓形目標(biāo)是現(xiàn)實(shí)生活中極為常見的目標(biāo)，由于透視投影變換圓在圖像坐標(biāo)系蛻變?yōu)闄E圓，快速魯棒的橢圓曲線提取有重要的意義。尋找一種快速、魯棒、精確的橢圓檢測算法一直是相關(guān)領(lǐng)域研究熱點(diǎn)。目前常用的橢圓檢測算法主要有最小二乘擬合法[1]、Hough變換法[2]以及基于橢圓幾何特征的算法[3]等。最小二乘擬合方法是橢圓擬合中最為常用的方法，該方法直觀、簡單、實(shí)時(shí)性好并且可以達(dá)到較高的擬合精度，能夠應(yīng)用于各種復(fù)雜的對象模型。目前大多數(shù)的橢圓檢測算法所針對的是靜態(tài)圖像，對序列圖像中動態(tài)橢圓目標(biāo)的實(shí)時(shí)檢測研究較少。在筆者所研究的計(jì)算機(jī)視覺測量系統(tǒng)中，需要對序列圖像中的橢圓目標(biāo)予以檢測及識別，通過分析各幀圖像中橢圓目標(biāo)的中心位置、長短軸等參數(shù)對目標(biāo)的運(yùn)動參數(shù)（位置、速度、姿態(tài)）予以估計(jì)。實(shí)驗(yàn)證明：本方法可以實(shí)現(xiàn)序列圖像中橢圓目標(biāo)的快速、魯棒檢測。

1 序列圖像中橢圓檢測算法

在二維平面上，橢圓和其他二次曲線方的一般形式可以描述為：

當(dāng)上式中的參數(shù)滿A2－4B<0時(shí)，該方程表示的是橢圓方程。另一種更直觀的橢圓表示方式是橢圓幾何參數(shù)法，假設(shè)橢圓的中心為（uc，vc），長軸和短軸分別為（2a 2b），橢圓長軸與坐標(biāo)系橫軸的夾角為θ，可以利用式（2）～（6）將橢圓的一般式轉(zhuǎn)換成幾何參數(shù)表示形式。

隨機(jī)抽樣一致性算法[4]RANSAC（RANdom SAmple Consensus）是一種魯棒性的模型估計(jì)算法，基于RANSAC的改進(jìn)的最小二乘橢圓擬合算法[5]能有效剔除錯誤樣本點(diǎn)的影響。對于序列圖像中的單幀圖像經(jīng)邊緣檢測后，邊緣二值圖像除包含待檢測橢圓的邊緣信息之外還包含大量的偽邊緣及其他干擾信息，分別對邊緣擬合之后會得到多個橢圓目標(biāo)，有必要利用先驗(yàn)知識對擬合得到的橢圓進(jìn)行篩選。用以下步驟對序列圖像中的橢圓目標(biāo)予以檢測，算法流程如圖1所示。

1）根據(jù)上一幀圖像中卡爾曼預(yù)測的結(jié)果設(shè)置本幀圖像中橢圓待檢測區(qū)域，使用Canny算子得到待檢測區(qū)域的邊緣點(diǎn)集；

2）利用邊緣特征點(diǎn)的連通性將圖像邊緣分為不同的輪廓，選取輪廓，若輪廓點(diǎn)的個數(shù)大于，則使用基于RANSAC的改進(jìn)最小二乘法擬合橢圓。否則，重新選取輪廓；

3）計(jì)算所擬合橢圓的橢圓率、擬合率以及面積，如果三者都滿足閾值條件，則將擬合的橢圓記為本幀圖像中的目標(biāo)橢圓。否則，返回步驟2）重新選擇輪廓；

4）返回步驟1）檢測下一幀圖像中的橢圓目標(biāo)。

其中，卡爾曼濾波預(yù)測的相關(guān)內(nèi)容在本文第二節(jié)予以介紹。橢圓率定義為λ=b/a，b為橢圓的短半軸，a為橢圓的長半軸，橢圓率描述了橢圓接近于圓的程度。擬合率定義為p=m/n，其中n為擬合橢圓所用輪廓邊界的總點(diǎn)數(shù)，m為該邊界上的點(diǎn)到擬合橢圓距離小于設(shè)定閾值的點(diǎn)集，擬合率描述了擬合橢圓所用輪廓與所得到橢圓之間的擬合程度。橢圓面積S=πab?？梢愿鶕?jù)先驗(yàn)知識及前期實(shí)驗(yàn)設(shè)定待檢測橢圓的面積、橢圓率以及擬合率的閾值范圍。

圖1 橢圓檢測流程Fig.1 Flow of ellipse detection

2 卡爾曼濾波預(yù)測橢圓搜索區(qū)域

對整幅圖像進(jìn)行橢圓檢測需要較多的時(shí)間，如果能夠?qū)ο乱粠瑘D像中橢圓的位置有預(yù)先粗略的估計(jì)，將橢圓的待檢測區(qū)域限定在一個較小的范圍內(nèi)，可以減少算法所需要的存儲空間和運(yùn)行時(shí)間。另外，在全局范圍內(nèi)搜索橢圓目標(biāo)使得算法的抗干擾能力也不理想。利用卡爾曼濾波器的預(yù)測功能來預(yù)測下一幀圖像中待檢測橢圓目標(biāo)中心的位置，并設(shè)置下一幀圖像中橢圓目標(biāo)的大致搜索范圍。通過這種方法大大減小了算法的運(yùn)算量，保證了系統(tǒng)的實(shí)時(shí)性及魯棒性。

2.1 卡爾曼預(yù)測器

卡爾曼預(yù)測器[6]實(shí)質(zhì)上是一種遞歸預(yù)測算法，根據(jù)誤差無偏最小方差原則預(yù)測目標(biāo)下一時(shí)刻的運(yùn)動狀態(tài)?？柭鼮V波可以劃分為兩個過程：時(shí)間更新和測量更新。時(shí)間更新方程用于向前預(yù)測下一時(shí)刻狀態(tài)X （k）-及計(jì)算預(yù)測協(xié)方差P（k）-。測量更新方程用于估計(jì)當(dāng)前時(shí)刻狀態(tài)的最優(yōu)值X（k），計(jì)算卡爾曼增益Kg（k）并更新最優(yōu)值的協(xié)方差P（k）。在卡爾曼濾波器初始化后，通過圖2中卡爾曼濾波的5個基本方程就可以保證算法自回歸的運(yùn)行下去。在利用卡爾曼濾波的預(yù)測功能時(shí)系統(tǒng)所輸入的是k-1時(shí)刻系統(tǒng)狀態(tài)的測量值Z（k-1），輸出的是k時(shí)刻系統(tǒng)的預(yù)測值X（k）-。在沒有控制輸入的情況下Kalman濾波的遞歸運(yùn)算過程如圖2所示。

圖2 卡爾曼濾波預(yù)測流程Fig.2 Flow of Kalman filter prediction

2.2 橢圓檢測區(qū)域預(yù)測算法實(shí)現(xiàn)

考慮運(yùn)動速度、加速度對橢圓中心位置預(yù)測的影響，使用三階向量卡爾曼預(yù)測器對下一幀圖像中的橢圓中心位置進(jìn)行預(yù)測。為保證算法的實(shí)時(shí)性，設(shè)計(jì)了兩個卡爾曼預(yù)測器分別跟蹤橢圓中心在圖像坐標(biāo)系橫軸和縱軸方向，將一個二維預(yù)測降為兩個一維預(yù)測[7]，提高了算法的運(yùn)行速率。假設(shè)目標(biāo)物體運(yùn)動的軌跡投影在一定時(shí)間間隔T內(nèi)的加速度是常量，則：

上式中，yk-1、vk-1、ak-1分別表示 t=k-1時(shí)刻橢圓中心在圖像坐標(biāo)系縱軸方向上的位置、速度及加速度，T為兩幀圖像間的時(shí)間間隔。用示t=k時(shí)刻單個預(yù)測器的狀態(tài)向量，則橢圓中心運(yùn)動狀態(tài)的系統(tǒng)方程可以寫為：

系統(tǒng)量測方程為

上式中，C=[1 0 0]為量測矩陣，Z（k）為 t=k時(shí)刻所測量得到的橢圓中心在圖像坐標(biāo)系中縱軸方向的坐標(biāo)，V（k）為測量噪聲，其協(xié)方差矩陣為R。將基于卡爾曼濾波的橢圓中心位置預(yù)測算法描述為：

1）初始化卡爾曼濾波器，設(shè)置卡爾曼預(yù)測器的噪聲方差Q和R，并給定狀態(tài)初始值X（0）和對應(yīng)協(xié)方差矩陣初始值 P（0）；

2）由時(shí)間更新方程并根據(jù)k-1時(shí)刻的最優(yōu)估計(jì)值X（k-1）預(yù)測k時(shí)刻系統(tǒng)狀態(tài)X（k）-，同時(shí)更新k時(shí)刻狀態(tài)預(yù)測值的協(xié)方差 P（k）-，并將預(yù)測值 X（k）-作為系統(tǒng)的輸出；

3）根據(jù)測量更新方程可以計(jì)算得到k時(shí)刻的卡爾曼增益Kg（k），輸入k時(shí)刻的狀態(tài)測量值Z（k）并結(jié)合k時(shí)刻狀態(tài)預(yù)測值X（k）-可以得到 k時(shí)刻的最優(yōu)估計(jì)值X（k）。計(jì)算k時(shí)刻最優(yōu)估計(jì)值的協(xié)方差P（k），并返回步驟2）使算法自回歸的運(yùn)行下去。

假設(shè)由第k-1幀圖像預(yù)測的第k幀圖像中橢圓的中心坐標(biāo)為（uc，vc），在 k-1幀計(jì)算的橢圓的長軸為 r，則將第 k幀圖像中橢圓的待檢測區(qū)域設(shè)置為以（uc，vc）為中心，以2（r-Tr）為邊長的矩形區(qū)域[8]。圖3中實(shí)線框?yàn)橛煽柭鼮V波預(yù)測的下一幀圖像中橢圓的待檢測區(qū)域。其中，閾值可以根據(jù)目標(biāo)橢圓運(yùn)動情況進(jìn)行調(diào)整，對于運(yùn)動變化較快的目標(biāo)可以適當(dāng)增大，此時(shí)橢圓檢測范圍增大Tr，穩(wěn)定性提高。由卡爾曼濾波的性質(zhì)可知，濾波器初始化后的前幾次預(yù)測并不可靠，此時(shí)應(yīng)適當(dāng)增大Tr保證算法的穩(wěn)定性，經(jīng)過若干幀跟蹤之后，以較小的Tr即可保證待檢測橢圓位于預(yù)測的檢測區(qū)域內(nèi)。

圖3 卡爾曼預(yù)測橢圓區(qū)域Fig.3 Ellipse location prediction with kalman predictor

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

3.1 仿真實(shí)驗(yàn)

利用MATLAB生成含有噪聲的標(biāo)準(zhǔn)橢圓邊界點(diǎn)并添加部分錯誤點(diǎn)，分別利用最小二乘法和改進(jìn)的方法擬合橢圓，兩者檢測結(jié)果如圖4所示。由圖可以看出由于樣本點(diǎn)中含有一定比例的不屬于待擬合橢圓邊界的點(diǎn)集，最小二乘法所擬合橢圓并不準(zhǔn)確，而基于RANSAC方法能排除這些錯誤點(diǎn)的干擾，準(zhǔn)確的擬合出目標(biāo)橢圓。

3.2 實(shí)際圖像檢測結(jié)果

圖4 橢圓擬合結(jié)果Fig.4 Result of ellipse fitting

在序列圖像動態(tài)橢圓目標(biāo)檢測實(shí)驗(yàn)中，對自由下落并反彈的球體的投影圖像進(jìn)行橢圓檢測。所錄制視頻圖像分辨率為320*240，幀頻為20 frame/s，并設(shè)置待檢測區(qū)域閾值Tr=30。圖5為從該段序列圖像中所截取的幾幀圖像，圖像中的方形區(qū)域?yàn)橥ㄟ^卡爾曼濾波器所預(yù)測的橢圓待檢測區(qū)域。各幀圖像中所預(yù)測的方形區(qū)域邊長最大不超過50像素，待檢測區(qū)域的面積僅僅是全局圖像面積的1/30，大大減小了下一幀圖像中橢圓的搜索區(qū)域。在不同的實(shí)驗(yàn)狀況下，可以根據(jù)目標(biāo)運(yùn)動的劇烈情況設(shè)置不同的Tr來保證橢圓目標(biāo)位于待檢測區(qū)域內(nèi)。圖5中的橢圓為按本文方法所擬合的橢圓曲線與原二值邊緣圖像疊加的結(jié)果。

圖6為橫向和縱向的卡爾曼預(yù)測結(jié)果及誤差，從圖6（a）和圖6（c）可以看出經(jīng)過幾次迭代之后，卡爾曼濾波器在兩個方向上都可以準(zhǔn)確預(yù)測出下一幀圖像中橢圓中心的實(shí)際位置。圖6（b）和圖6（d）也表明濾波器可以較為準(zhǔn)確的預(yù)測出下一幀圖像中橢圓目標(biāo)的中心位置，最大預(yù)測誤差在20個像素之內(nèi)。通過卡爾曼濾波預(yù)測可以較為準(zhǔn)確的預(yù)測下一幀圖像中橢圓目標(biāo)的待檢測區(qū)域。

4 結(jié) 論

為實(shí)現(xiàn)對序列圖像中的動態(tài)橢圓實(shí)時(shí)檢測，結(jié)合RANSAC算法和最小二乘法設(shè)計(jì)了魯棒性的橢圓檢測流程，能有效剔除橢圓擬合過程中干擾邊緣和噪聲點(diǎn)的影響。序列圖像中橢圓目標(biāo)的檢測對系統(tǒng)的運(yùn)算速度要求較高，采用卡爾曼濾波預(yù)測能夠充分的利用歷史信息，縮小下一幀圖像中對待檢測橢圓目標(biāo)的搜索區(qū)域，提高了系統(tǒng)的運(yùn)算速度和穩(wěn)定性。此外，當(dāng)在某一幀圖像中未能檢測出橢圓目標(biāo)時(shí)可以根據(jù)預(yù)測結(jié)果繼續(xù)對下一幀圖像中的橢圓進(jìn)行檢測，提高了系統(tǒng)的可靠性。目前已成功將該算法應(yīng)用于所開發(fā)的水下目標(biāo)運(yùn)動參數(shù)視覺測量系統(tǒng)中，取得了令人滿意的效果。

[1]鄒益民，汪渤.一種基于最小二乘的不完整橢圓擬合算法[J].儀器儀表學(xué)報(bào)，2006，27（7）:809-810.ZOU Yi-min，WAN Bo.Fragmental ellipse fitting based on least square algorithm[J].Chinese Journal of Scientific Instrument，2006，27（7）:809-810.

[2]袁立，葉露，賈建祿.基于Hough變換的橢圓檢測算法[J].中國光學(xué)與應(yīng)用光學(xué)，2010，3（4）:380-382.YUAN Li，YE Lu，JIA Jian-lu.Ellipse detection algorithm based on hough transform[J].Chinese Journal of Optics and Applied Optics，2010，3（4）:380-382.

[3]于海濱，劉濟(jì)林.基于中心提取的RHT在橢圓檢測中的應(yīng)用[J].計(jì)算機(jī)輔助設(shè)計(jì)與圖形學(xué)學(xué)報(bào)，2007，19（9）:1108-1112.YU Hai-bin，LIU Ji-lin.Ellipse detection by the RHT based on center extraction[J].Hournal of Computer-Aided Desigh&Computer Graphics，2007，19（9）:1108-1112.

[4]Fischler M A，Bolles R C.Random sample consensus:a paradigm for model fitting with applications to image analysis and automated cartography[J].Graphics and Image Processing，1981，24（6）:381-383.

[5]閆蓓，王斌，李媛.基于最小二乘法的橢圓擬合改進(jìn)算法[J].北京航空航天大學(xué)學(xué)報(bào)，2008，34（3）:295-298.YAN Pei，WANG Bin，LI Yuan.Optimal ellipse fitting method based on least-square principle[J].Journal of Beijing University of Aeronautics and Astronautics，2008，34 （3）:295-298.

[6]郝利云，吳玲達(dá).基于卡爾曼預(yù)測器的數(shù)據(jù)預(yù)取方法[J].計(jì)算機(jī)仿真，2012，29（10）:10-12.HAO Li-yun，WU Ling-da.Kalman predictor based data perfecting method[J].Computer Simulation，2012，29（10）:10-12.

[7] 姚紅閣，耿軍雪. 基于卡爾曼預(yù)測的視頻目標(biāo)實(shí)時(shí)跟蹤[J].西安工業(yè)大學(xué)學(xué)報(bào)，2007，27（2）:174.YAO Hong-ge，GENG Jun-xue. A real-time tracking for videotarget based on kalman prediction[J]. Journal of Xi’an TechnoloficalUniversity，2007，27（2）:174.

[8] 王忠立，高文. 基于最小二乘預(yù)測的動態(tài)快速橢圓檢測[J].信息與控制，2003，32（7）:732.WANG Zhong-li，GAO Wen. Dynamic fast ellipse detectionbased on location prediction with least square[J]. Informationand Control，2003，32（7）:732.