一種組合GAUSS-filter、SOBEL、NMS、OTSU 4種算法的圖像邊緣檢測的FPGA實現

張 灝，王素珍，鄭 宇，王 偉，李濰志，王 鵬，任貴珊，馬家麟

(青島大學 電子信息學院，山東 青島 266071)

1 引 言

邊緣檢測是分離圖像目標主體和背景的主要依據。通過檢測圖像上灰度變化顯著的像素點集合來確定圖像的邊緣。圖像的邊緣檢測在許多領域，比如醫學影像、地質勘探、衛星遙感、航天、交通、計算機視覺等方面都有重要的應用。目前對圖像清晰度的要求越來越高，圖像處理的信息量巨大。一般串行方式已經不能滿足對圖像大量信息的實時處理。大規模可編程片上系統(Field Programmable Gate Array，FPGA)芯片以其獨特的并行可編程處理信號的方式，成為處理圖像的重要器件。

邊緣檢測的常用方法是SOBEL算法。文獻[1-5]成功地在FPGA芯片上應用SOBEL算法實現了圖像的邊緣檢測。但是，其邊緣檢測的門限閾值需要反復的實驗確定，且圖像有噪聲的影響，同時邊緣輪廓也不夠細化。文獻[6]引入高斯濾波(GAUSS-filter)，對圖像先去噪聲處理，再送入SOBEL算法進行邊緣檢測，減少了圖像邊緣背景噪聲的干擾。文獻[7]把經過SOBEL算法檢測出的邊緣圖像送到非極大值抑制(Non-maximum Suppression，NMS)算法模塊處理，細化了圖像輪廓的邊緣。文獻[8-10]將大津法OTSU加入SOBEL算法中，使圖像邊緣的門限閾值能夠自適應圖像特點調整。

本文在FPGA器件上，綜合應用GAUSS-filter、SOBEL、NMS、OTSU等多種算法，進行圖像邊緣檢測。將原始采集到的彩色圖像轉換成灰度圖像，先對圖像進行高斯平滑處理，消除背景干擾性噪聲；應用SOBEL算法計算灰度梯度值，得到圖像邊緣；再應用NMS算法處理灰度，對比出最大的灰度像素集合；最后應用OTSU算法自適應圖像特點確定圖像門限閾值，提高圖像邊緣檢測的適應性。實驗結果證明本文方法對不同圖像都有很好的邊緣檢測效果。在運行速度上，比軟件實現時間少很多，提高了系統處理的速度。該實時圖像邊緣檢測方法能夠嵌入式地應用到許多領域中。

2 基于FPGA的實時邊緣檢測系統

本系統采用XILINX公司的ARTIX7系列的xc7a100tfgg484-2芯片作為邊緣檢測的主要器件，同時在FPGA主控芯片外圍配置雙倍速率同步動態隨機存儲器DDR3 SDRAM(Double Data Rate SDRAM)緩存圖像。基于FPGA的實時邊緣檢測系統如圖1所示。采用乒乓操作方法將來自CCD攝像機的分辨率為640×480的圖像先緩存到DDR3存儲器中，在圖像請求處理信號的作用下，圖像數據流將從DDR3存儲中被讀出，進行圖像的邊緣檢測處理。原始采集到的圖像為彩色圖，將圖像由彩色空間轉換到灰度空間；轉換后的灰度圖像混有噪聲，在進行圖像邊緣檢測之前，先將圖像送到高斯濾波器中處理，消除圖像中的噪聲。

圖1 整體結構框圖

高斯濾波處理后的圖像送到SOBEL算法模塊中，進行梯度的計算，包括X、Y方向的梯度計算，進而計算出像素的梯度值，這些梯度像素點對應的輪廓比較粗糙。為一步對圖像邊緣進行細化處理，要對梯度值進行非極大值抑制處理，即判斷當前像素點的梯度值是否是梯度方向上的極大值，如果是，則保留該梯度值；否則，則將其歸零。圖像的邊緣檢測需要確定閾值，采用OTSU算法，求取目標與背景之間的灰度類間方差來確定閾值。圖像的梯度值與閾值進行比較，得到圖像邊緣，將檢測到的邊緣輪廓送到顯示器進行顯示。

3 算法原理

3.1 SOBEL算法

SOBEL算法是通過SOBEL算子矩陣與圖像鄰域矩陣進行卷積，計算X方向和Y方向的梯度分量，來獲取圖像的邊緣信息。圖2為邊緣檢測模板。

(a)3×3 neighborhood window (b)Edge detection openator

圖2 邊緣檢測模板

Fig.2 Edge detection template

圖2(a)為中心像素點的3×3鄰域窗口，圖2(b)為SOBEL算法水平方向、垂直方向的梯度算子。將像素鄰域與X、Y方向梯度算子進行卷積運算，分別可以得到水平方向的梯度值Gx、垂直方向的梯度值Gy為

(1)

3×3鄰域的梯度卷積值為

(2)

其中G為當前鄰域內的梯度值。因為FPGA不容易執行開平方的操作，為了方便FPGA內部電路的計算，此處將鄰域的梯度卷積值簡化表示為

G=|Gx|+|Gy|.

(3)

3.2 NMS算法

用GAUSS filter和SOBEL算法進行邊緣檢測可以得到粗線條的圖像輪廓。對圖像輪廓細化處理的算法是非極大值抑制NMS。NMS是將鄰域中心像素點的梯度值與它所在梯度方向上的梯度值的線性插值相比較，判斷其是否為最大值，如果是，則保留該像素點梯度值；否則，直接將鄰域中心像素點的梯度值歸為0。這種方法能夠去除圖像中存在的偽邊緣像素點[11]。

非極大值抑制必須要進行梯度方向的判斷，梯度方向計算公式一般為θ(i,j)=arctan(Gx/Gy)，使用FPGA實現反正切函數相對較難。在論文中是通過判斷梯度值在X方向Gx的分量和在Y方向的Gy分量之間大小，以及兩者的正負來判斷當前梯度方向在哪一個“方向象限”，這樣去處理當前的梯度方向會減少一定的計算量。

(a)θ(i,j)∈[0°,45°) (b)θ(i,j)∈[45°,90°)

(c)θ(i,j)∈[90°,135°) (d)θ(i,j)∈[135°,180°)

如圖3所示，其中Gi為3×3鄰域中的像素點的梯度值，mi為中心像素點G4的梯度方向與鄰域的交點。圖中只畫出[0°,180°]上4個“方向象限”內的梯度方向圖，[180°,360°]上梯度方向由此類推。梯度方向上的值可以由鄰域內相應的梯度值進行線性插值求得。例如：圖3(b)的m4點的值可以由G6和G7進行線性插值[12-13]求出。G4所有梯度方向上的線性插值公式如表1所示，其中x和y分別代表|Gx|、|Gy|。

表1 線性插值公式表

從表1中可以看出，通過線性插值公式可得到[0°,360°]梯度方向上的mi值，NMS方法是G4與mi的值進行比較，如果G4比mi的值大，則保留G4;否則將G4置為0。這樣經過NMS處理后的任意像素點的梯度值可表示為公式(4)，其中KnmsG4為非極大值抑制結果。

(4)

3.3 OTSU算法

設置一幅圖像的灰度閾值就能進行邊緣檢測。如果設置閾值過大，會丟失一些邊緣信息；如果設置的閾值較小，那么會出現一些假的邊緣，使得輪廓不準確。為了使邊緣最大程度地體現目標信息，就得恰當地確定灰度閾值。本文使用OTSU算法自適應圖像灰度來確定當前的閾值。OTSU算法通過計算背景和目標灰度之間的類間方差來區分圖像邊緣，類間方差越大，證明背景和目標之間的差別越大[14-16]。

假設w(D0)為背景D0出現的概率，w(D1)為目標D1出現的概率，u(D0)和u(D1)分別為D0和D1的灰度均值。w(D0)、w(D1)、u(D0)以及u(D1)計算公式如公式(5)所示：

(5)

其中：P(i)為灰度級i出現的概率，T為背景和目標的分界灰度級，L為最大灰度級。則整幅圖像的背景與目標之間的類間方差g為：

g=w(D0)*w(D1)*(u(D0)-u(D1))2，

(6)

其中：w(D1)=1-w(D0)，將公式(5)中變量代入公式(6)可得：

(7)

其中N(i)為灰度級i的統計數。不同的圖像對應不同的輪廓邊緣，其對應的最大類間方差不同，當每幅圖像的g達到最大時，此時的灰度級i即為我們所需要的最佳閾值。

4 FPGA實現

4.1 基于FPGA的SOBEL算法實現

在FPGA芯片上實現SOBEL算法的過程中，產生3×3鄰域領域滑動窗口是關鍵技術。本文采用XILINX公司的SHIFT_RAM模塊進行3×3鄰域窗口設計，其深度設置為一幅圖像的一行像素點數，寬度為8位。

圖4 水平方向和垂直方向梯度計算

如圖4所示，構建3×3鄰域需要2個SHIFT_RAM緩存模塊，每行3個D觸發器，共3行。3×3鄰域窗口滑窗式地遍歷整幅圖像實現SOBEL算法，從3×3鄰域窗輸出的數據對應圖2(a)中的像素灰度值。圖4右是實現公式(1)的結構框圖，計算輸出水平方向Gx和垂直方向Gy。圖5是實現公式(3)的電路結構圖，判斷X、Y方向梯度值正負，再分別相加計算SOBEL算法的梯度值。

圖5 SOBEL算子梯度值計算

4.2 基于FPGA的NMS算法實現

4.2.1 NMS公式改進

由表1可看出，在進行線性插值計算時，會涉及除法，因為FPGA計算除法需要消耗較多的資源，本文將除法轉化為乘法，實現插值計算。將“方向象限”線性插值公式兩邊同時乘以分母，變除法為乘法運算。例如線性插值m1的公式可表示為xm1=yG2+(x-y)G5。在m1的計算中，x>y，取Mbig=max(x,y)，Msmall=min(x,y)，則有：

Mbig*m1=Msmall*G2+(Mbig-Msmall)*G5，

(8)

同理，在m2～m8的計算中，通過比較x、y之間的大小關系，將線性插值公式總結為式(9)：

Mbig*mi=Msmall*GA+(Mbig-Msmall)*GB，

(9)

其中GA和GB分別為相應的3×3鄰域內的兩個梯度值。相應地在非極大值抑制判斷表達式中進行改變，可得改進后的判斷表達式：

(10)

4.2.2 方向象限編碼

在FPGA實現公式(9)的計算中，首先通過判斷Gx和Gy的大小以及正負來確定鄰域內的中心像素點的梯度方向，為此用b2b1b0對8個“方向象限”方向進行編碼，用info_data表示，便于程序運行時進行查找表計算。可得到表2所示內容。

表2 方向象限編碼

如表2所示，“方向象限”編碼信號info_data的b2位為Gy的符號位，b1位為Gx的符號位，b0位是兩者之間的大小關系，如果|Gx|>|Gy|，則為0，反之，則為1。例如，第5、6“方向象限”，其Gx、Gy均為負值，則info_data的b2、b1位都為1，但確定梯度方向具體在二者中哪個，則要根據Gx和Gy的大小進行判斷，第5“方向象限”|Gx|>|Gy|，“方向象限”編碼為110，第6“方向象限”|Gy|>|Gx|，“方向象限”編碼為111。根據表2中的“方向象限”信息看出，每個“方向象限”的Gx和Gy大小情況不同，為統一計算，將其中較大的值賦值給Mbig，較小的值賦值給Msmall。圖6為“方向象限”信息仿真圖，out_x_data=Mbig為梯度方向上梯度值的較大值，out_x_data=Msmall為梯度方向上梯度值的較小值，out_info=info_data為“方向象限”編碼信息。

圖6 “方向象限”信息仿真圖

4.2.3 NMS內插公式的FPGA實現

圖7為非極大值實現框圖，上半部分為非極大值抑制部分，下半部分為3×3鄰域生成部分。要進行非極大值抑制計算，需要構建3×3鄰域窗口，窗口構建仿照SOBEL算法模塊實現。進入窗口的數據是經過“方向象限”信息處理后的Mbig、Msmall、out_info與Gi拼接的數據。經3×3鄰域窗口緩存后，分別提取Mbig、Msmall、out_info與Gi對應的鄰域中心點的值data_max，data_min,out_info_4以及G4，將這些數據送入非極大值抑制模塊進行非極大值抑制計算。將“方向象限”編碼值out_info_4作為LUT查找表的地址，進而尋找到線性插值的3×3鄰域梯度值數據GA和GB。例如，當3×3滑動窗口鄰域的中心梯度值的out_info為001時，插值公式對應梯度值數據分別為G2、G1、G6、G7，則梯度插值索引值{c3,c2,c1,c0}={2,1,6,7}，根據索引值也可找到梯度值數據。圖中減法器用來計算Mbig與Msmall的差值Mbig-Msmall，乘法器主要負責公式(9)中的乘法操作，加法器實現式公式(9)左端兩數據的相加，比較器COM、與門&、二選一數據選擇器MUX則是用來實現公式(10)非極大值的判斷，使用D觸發器實現數據延時使數據之間進行對齊。DOUT為非極大值抑制結果，且nmsG4=DOUT。

圖7 非極大值實現原理圖

4.3 OTSU全局閾值分割

要實現OTSU全局閾值分割，首先是要統計當前灰度級所對應的像素點數；對于灰度級的統計，采用直方圖統計方法來進行。送入到直方圖統計的數據是經過非極大值抑制后的數據DOUT。直方圖統計主要由雙口RAM及檢測控制器組成，其中雙口RAM位寬設置為8位，深度為32位，雙口RAM用于存儲每個灰度級所對應的統計值。直方圖統計的思想如下：利用每個像素點所對應的灰度值作為雙口RAM的存儲及讀取地址，當有當前圖像的梯度值輸入時，雙口RAM中相應地址內的統計值會被讀出與遞增計數器數據進行相加，直到輸入的梯度值不是該地址時，會將此梯度值對應的統計值緩存到RAM，并進行下一個梯度值的統計。

圖8 直方圖時序設計圖

直方圖統計的時序設計圖如圖8所示。對進入直方圖統計模塊的灰度值數據DIN延遲得到DIN_1、DIN_2。在論文中，控制器中會有一個遞增計數器CNT，遞增計數器使能信號INC_EN拉高的條件是DIN_1與DIN_2相等，使得CNT遞增，用來計算連續灰度值的個數；否則該計數器對應的復位使能信號RST_CNT拉高，CNT置為1，不再計數。對于涉及到的雙口RAM，B端口的地址address_b為DIN_1，A端口的地址address_a為DIN_2，兩者之間存在時鐘延遲是為了讓B口讀出的統計值q_b剛好與A口地址對應的遞增計數器CNT對齊。A端口的寫使能WR_A與B端口的讀使能RD_B拉高的條件是DIN_1與DIN_2數據不相等，此時需要將前一灰度值已統計好的數據CNT_VLE寫入RAM，CNT_VLE為q_b與CNT的和，同時將當前灰度值的統計值q_b從RAM中讀出，進行當前灰度值的統計計算。例如，當DIN為222時，延遲后在DIN_1、DIN_2上分別出現灰度值222，因為灰度值222與前一灰度值76不相等，需要將前一灰度值76的統計數據CNT_VLE寫入RAM，將上一次灰度值222出現時已統計好的數據q_b讀出，同時在INC_EN拉高的控制下q_b與CNT相加，將統計值CNT_VLE從38變到42。

圖9 最大類間方差實現圖

直方圖統計完成后，每一個灰度級所對應的統計值N(i)就存儲在雙口RAM中對應的地址中，背景和目標之間的最大類間方差需要使用這些統計值進行計算。在類間方差計算時，按照地址從0～255的順序將梯度的統計值讀出，將梯度值與梯度值的統計值進行對齊，方便兩者的相乘運算。圖9為OTSU算法的公式(7)的計算框圖。為了在FPGA中實現OTSU算法過程中省去兩個除法器，此處去掉公式(7)右半部分的除數N×N，這樣處理不會影響類間方差之間的比較，節約了資源。圖10為OTSU計算閾值仿真圖。

圖10 OTSU閾值仿真圖

5 測試結果分析

如圖11所示為在FPGA芯片上實現的未進行高斯濾波的SOEBL算法、進行高斯濾波后的不同閾值下的SOBEL算法、濾波后SOBEL和OTSU算法結合、濾波后的SOBEL和NMS算法結合、綜合GAUSS-filter、SOBEL、NMS、OTSU算法后的圖像邊緣檢測比較結果。

(a)實際拍攝的人物圖片處理效果

(b)實際拍攝的物體結構圖片處理效果

(c)麗達圖片處理效果

圖11(a)是對實際拍攝的人物圖片處理的效果；圖11(b)是對實際拍攝的物體結構圖片處理的效果；圖11(c)是對畫面層次豐富的麗達圖片處理的效果。其中A圖為原始灰度圖;B圖為未進行高斯濾波的SOEBL算法邊緣檢測效果圖，因為未進行高斯濾波，所以有很多噪聲存在圖像中；C圖、D圖、E圖為只用SOBEL算法實現的不同閾值下的邊緣檢測效果，其中C圖是因為閾值過小，導致過多非圖像邊緣信息顯示出來的效果；D圖為閾值過大，圖像邊緣不能全部顯示的效果，輪廓不連續；E圖是經過多次實驗找出的恰當閾值下所顯示圖像邊緣。由此可見，單一的SOBEL算法需要反復測試來確定一個恰當的閾值，而且不同的圖像需要不同的閾值，不便于對圖像邊緣的快速檢測。F圖為SOBEL和OTSU算法結合的圖像邊緣檢測效果，克服了多次尋找閾值的缺點；G圖為SOEBL算子與NMS算法的結合，將E圖像邊緣進行細化,但對圖像缺少自適應性；H圖為綜合GAUSS-filter、SOBEL、NMS、OTSU算法后的圖像邊緣檢測結果。比較圖11中的3組圖片的邊緣檢測效果可得，綜合應用4種算法的H圖不僅能夠畫出圖像邊緣，而且還能對邊緣進行細化，對不同圖像都有好的邊緣檢測效果。

表3 運行時間對比

經過多次實驗，得出論文中使用FPGA實現邊緣檢測的效果與軟件實現效果差別不大，但是在實現時間上兩者有較大差別。表3是圖像處理軟、硬件運行的平均時間比較，從表3看出，使用FPGA實現邊緣檢測比軟件實現時間少很多。

文中提到的硬件實現，是FPGA設計軟件上的功能驗證，將圖像處理過程在硬件仿真工具modelsim上實現；提取txt文件，然后在軟件工作環境QT上進行圖像處理效果比較。

FPGA硬件實現時間僅與工作時鐘頻率有關。文中軟件實現與硬件實現速度的比較是指圖像處理過程分別在硬件仿真工具modelsim上運行和在軟件工作環境QT運行時間的比較。

在實驗中我們使用的PC機配置參數為：64位 win7操作系統，運行內存為6 G，主頻為2.7 GHz。FPGA的晶振為200 MHz，其中圖像工作頻率為40 MHz。

6 結 論

在FPGA上綜合應用GAUSS-filter、SOBEL、NMS、OTSU等多種算法，進行圖像邊緣檢測。將原始采集到的彩色圖像轉換成灰度圖像，先對圖像進行高斯平滑處理，消除干擾性噪聲；應用SOBEL算法計算圖像梯度值，再應用NMS算法處理灰度，挑出輪廓邊緣梯度值最大的灰度像素集合，最后應用OTSU算法自適應地確定圖像門限閾值進行邊緣檢測。單一SOBEL算法需要人為設置閾值，閾值容易出現過大或者過小的現象，導致圖像出現假邊緣或者不連續邊緣；SOBEL和OTSU算法結合又不能對圖像進行細化處理。4種算法綜合的檢測不僅準確找到圖像邊緣，而且還能對邊緣進行細化，彌補了兩種邊緣檢測算法的不足。實驗結果顯示對不同的圖像都有很好的檢測效果。文中所提出的在FPGA上實現的圖像邊緣檢測方法具有準切、實時、高速、嵌入式等特點，能夠應用到許多領域中。