基于GPU的數(shù)學(xué)形態(tài)學(xué)運算并行加速研究

張 聰，邢同舉，羅 穎，張 靜，孫 強(qiáng)

（電子科技大學(xué) 光電信息學(xué)院，四川 成都 610054）

數(shù)學(xué)形態(tài)學(xué)的應(yīng)用幾乎遍及計算機(jī)圖形圖像處理的所有方面，包括圖像濾波、圖像分割分類、圖像測量、模式識別以及紋理分析與合成，其應(yīng)用還涉及遙感監(jiān)測、工業(yè)自動化檢測測量、生物醫(yī)學(xué)影像、圖像壓縮、軍事、航空航天等眾多領(lǐng)域。在運用數(shù)學(xué)形態(tài)學(xué)處理氣象預(yù)報、工程計算和模擬等問題時，需要對巨大數(shù)據(jù)量進(jìn)行多次重復(fù)計算，這些計算又必須在有限的時間內(nèi)完成，比如在幾秒或幾分鐘內(nèi)完成。這就對數(shù)學(xué)形態(tài)學(xué)在處理這些問題時的運算速度提出了較高的要求。在數(shù)學(xué)形態(tài)學(xué)中，腐蝕和膨脹是最基本的兩個運算，這兩種對偶運算好比形態(tài)學(xué)字母表中的字母，其他運算都可以根據(jù)這兩種運算表述[1]，因此形態(tài)學(xué)處理速度決定于膨脹、腐蝕處理速度。

人們通過在腐蝕、膨脹處理中采用改進(jìn)算法來提高其運算速度[2]，但其速度也只能提高幾倍，這在某些大數(shù)據(jù)量、實時性要求較高的應(yīng)用中仍不能滿足需要。針對這種情況，文章中提出了通過圖形處理器（GPU）并行處理的方式對數(shù)學(xué)形態(tài)學(xué)運算進(jìn)行加速。實驗結(jié)果表明，通過GPU對數(shù)學(xué)形態(tài)學(xué)運算進(jìn)行加速，其運算速度可達(dá)到1～2個數(shù)量級的提高。

1 膨脹腐蝕的基本理論

正文內(nèi)容腐蝕運算：構(gòu)造一個結(jié)構(gòu)元素，結(jié)構(gòu)元素的原點定位于待處理的目標(biāo)元素上，結(jié)構(gòu)元素完全包含在目標(biāo)區(qū)域中時，則結(jié)構(gòu)元素的原點所對應(yīng)的目標(biāo)像素點被保留，負(fù)責(zé)被腐蝕掉，置為背景點。膨脹運算：構(gòu)造一個結(jié)構(gòu)元素，結(jié)構(gòu)元素的原點定位在待處理的目標(biāo)元素上，當(dāng)有目標(biāo)點被結(jié)構(gòu)元素覆蓋時該點被膨脹為目標(biāo)點。

膨脹腐蝕運算中存在一些重要的運算性質(zhì)，在設(shè)計算法時使用這些性質(zhì)往往能簡化處理過程和提高運算速度。對偶性便是其中的一個重要性質(zhì)，對偶性可用下面公式來表示：（XΘB）C=XC⊕B；（X⊕B）C=XCΘB。 這兩個公式表明，圖像膨脹等價于該圖像的補(bǔ)被相同結(jié)構(gòu)元素腐蝕所得結(jié)果的補(bǔ)，反之亦然[3]。在本文中討論主要集中在腐蝕運算的并行加速，而膨脹運算可以通過待處理圖像補(bǔ)的腐蝕運算求得。

2 GPU簡介

2.1 GPU硬件架構(gòu)

與CPU相比，GPU有著很高的計算吞吐量，這主要在于GPU和CPU不同的設(shè)計理念。GPU最初設(shè)計就是為了圖形處理，由于圖像渲染的高度并行性，GPU設(shè)計者通過增加數(shù)學(xué)邏輯單元（ALU）和控制單元（Control）的方式提高處理能力和存儲器帶寬。這樣就在GPU硬件架構(gòu)中形成了眾多的流處理器，Nvidia GTX 280的架構(gòu)[4]如圖1所示。

圖1中：Nvidia GTX 280中有30個完整前端流多處理器（SM），每個流多處理器包含8個流處理器（SP）。每個SP又可看成一個SIMD處理器。每個SM擁有獨立的完整前端，包括取指、譯碼、發(fā)射和執(zhí)行單元等，但每兩個SM共享一條存儲器流水線。SP不具有獨立的處理器核，它們有獨立的寄存器和指令指針，但沒有取指和調(diào)度單元構(gòu)成的完整前端。

圖1 Nvidia GTX 280架構(gòu)Fig.1 Nvidia GTX 280 architecture

2.2 CUDA軟件開發(fā)環(huán)境

由于GPU具有很高的計算吞吐量，利用GPU進(jìn)行通用計算已經(jīng)引起了人們的關(guān)注[5－8]。CUDA架構(gòu)正是專門為了GPU通用計算而設(shè)計的一種全新模塊，它使得開發(fā)人員無需再過多地考慮嚴(yán)格的資源限制和編程限制。CUDA C是在CUDA GPU上的編程語言，CUDA C語言本質(zhì)上是對C的簡單擴(kuò)展，它的出現(xiàn)使得開發(fā)人員不需要學(xué)習(xí)計算機(jī)圖形學(xué)和著色語言，就可以編寫出高效的GPU并行程序。

CUDA 編程模型[9]將 CPU 作為主機(jī)（Host），GPU 作為設(shè)備（Device）[10]。在這個模型中GPU和CPU協(xié)同工作，CPU上負(fù)責(zé)邏輯性較強(qiáng)的事務(wù)處理和串行計算，GPU則專注于執(zhí)行高度并行化處理任務(wù)。運行在GPU上的CUDA函數(shù)稱為內(nèi)核函數(shù)（kernel），一個完整的CUDA程序是由Host串行代碼和Device上的kernel函數(shù)共同組成的，如圖2所示。

運行在設(shè)備上的kernel函數(shù)定義和普通C語言函數(shù)基本一致，當(dāng)需要被主機(jī)調(diào)用時，只需在函數(shù)定義前面加上__global__修飾符__global__void Copy（char*src， char*dst）。在主機(jī)調(diào)用kernel函數(shù)時，需要傳遞額外的參數(shù)，如：Copy＜＜＜blocksPerGrid， threadsPerBlock＞＞＞（ptr_src， ptr_dst）。kernel函數(shù)是以線程網(wǎng)格（Grid）的形式組織的，每個Grid由若干個線程塊（block）組成，每個block由若干線程（thread）組成。調(diào)用kernel函數(shù)所需傳遞的額外參數(shù)blocksPerGrid和threadsPerBlock分別是Grid和block的維度設(shè)計。

3 基于GPU的腐蝕運算加速實現(xiàn)

圖像腐蝕采用（2n+1）×（2n+1）且系數(shù)全為 1的模板，即進(jìn)行n尺寸腐蝕。要確定某一像素點是否被腐蝕掉，只需將模板中心對準(zhǔn)當(dāng)前待處理像素點，查看被模板覆蓋的所有像素點，如果全為目標(biāo)點則當(dāng)前待處理像素點保留為目標(biāo)點，否則腐蝕掉。腐蝕可進(jìn)一步轉(zhuǎn)化為鄰域像素點灰度值求和的問題，每一個像素點都由被模板覆蓋所有像素點灰度值的和決定。由此，確定每一點是否被腐蝕掉，所要執(zhí)行的操作是相同的，只是所處理的數(shù)據(jù)不同而已。這是典型的單指令多數(shù)據(jù)問題，具有很好的并行性。

圖2 CUDA編程模型Fig.2 CUDA programming structure

3.1 實現(xiàn)時需要注意問題

在GPU上編碼實現(xiàn)腐蝕運算時，為了正確地執(zhí)行和盡可能地提高執(zhí)行效率需要注意幾個問題：1）避免當(dāng)前處理的像素點被已處理像素點所干擾。在具體實現(xiàn)時，生成一幅和原圖同樣大小的新圖，每個線程將處理結(jié)果存放在這幅新圖中，而運算數(shù)據(jù)則從原圖中獲取；2）充分利用GPU硬件資源，盡可能為每個像素點都開辟一個線程，目前GPU每個block允許開辟的線程數(shù)量已達(dá)到1 024個，每個Grid允許每一維開辟的block數(shù)量為65 535。所以對幾乎所有的二值圖像，都能實現(xiàn)為每一個像素點開辟一個線程；3）腐蝕運算訪問內(nèi)存模式具有很強(qiáng)的空間局部性，GPU紋理緩存是專門為了加速這種訪問模式而設(shè)計的，在這種內(nèi)存訪問模式下，使用紋理內(nèi)存能減少內(nèi)存流量是性能得到提高。4）block維度設(shè)計時，每個block的線程數(shù)量應(yīng)該為64的倍數(shù)[11]，這樣GPU的性能才能最大限度地提高。

3.2 實現(xiàn)過程

上面已經(jīng)提到CUDA編程模型是GPU和CPU協(xié)同工作的，所以GPU腐蝕運算的實現(xiàn)過程是CPU和GPU共同完成的。CPU負(fù)責(zé)內(nèi)存分配、Grid和block的維度設(shè)定、啟動在GPU上運行的kernel函數(shù)并向其傳遞參數(shù)、最后從GPU獲取處理結(jié)果。GPU負(fù)責(zé)大量并行運算，并將處理數(shù)據(jù)保存下來。其具體顯現(xiàn)步驟如下：

1）讀取待處理圖像，根據(jù)圖像大小在device上開辟顯存devSrc，并將圖像數(shù)據(jù)復(fù)制到這塊顯存上，并將其綁定二維紋理內(nèi)存；

2）在device開辟和devSrc同樣大小的顯存devDst，用于GPU保存處理結(jié)果。在host上開辟和devDst同樣大小的內(nèi)存hostDst，用于接收GPU處理的結(jié)果；

3）根據(jù)待處理圖像大小，設(shè)計grid和block的維度。為了提高性能將 block 維度固定為 dim3 blockDim（16，16），則 grid的維度則根據(jù)圖像寬度width、高度height和block的維度設(shè)定，可設(shè)定為 dim3 gridDim（ （width+15）/16，（height+15）/16）；

4）啟動kernel函數(shù)，并將設(shè)定好的參數(shù)傳遞到kernel函數(shù)中；

5）kernel函數(shù)執(zhí)行完畢，將處理結(jié)果傳回host。host保存并顯示結(jié)果。

4 實驗結(jié)果及分析

文章所采用實驗平臺為：Interl（R）Xeon（R）CPU 5160，主頻 2.99 GHz，Windows7操作系統(tǒng)，顯卡為 GForce GTX 570。文章中所處理的為 1 024×1 024、7 350×5 700 的圖像。表格1為在對大小為1 024×1 024的圖像進(jìn)行不同尺寸腐蝕運算時，CPU串行腐蝕和GPU并行腐蝕所耗用時間對比。

表1 在不同腐蝕尺寸下CPU和GPU腐蝕時間對比Tab.1 The comparison of the erosion time between CPU and GPU in different erosion sizes

從表格1中可看出，隨著腐蝕尺寸的增大，即隨著運算時間復(fù)雜度的增加，CPU串行腐蝕運算所耗用時間大幅度增加，GPU并行腐蝕運算所耗用時間沒有大幅增加，GPU并行運算加速比越來越大。

在圖3中，可以看到在對圖像進(jìn)行同一腐蝕尺寸的腐蝕處理時，GPU腐蝕處理大小為7 350×5 700的圖像比處理1 024×1 024的圖像加速比要高。當(dāng)腐蝕尺寸增加到6時，對于1 024×1 024的圖像加速比可達(dá)到 80左右，對于7 350×5 700的圖像加速比可達(dá)150左右。GPU并行運算之所以有這么高的加速效果，主要因為它有如此多的數(shù)學(xué)邏輯單元，所以能同時運行較多的線程。它通過大量線程的切換來隱藏訪存延遲，當(dāng)某個線程被阻塞了，馬上切換到其他線程進(jìn)行處理，所以當(dāng)待處理數(shù)據(jù)規(guī)模越大時，GPU資源就越能被充分利用，加速效果就更好。

圖3 對不同大小圖片，GPU加速比隨著腐蝕尺寸的變化而變化Fig.3 In different size of images，the speedup of GPU changes when the erosion size changes

5 結(jié) 論

通過GPU并行處理的方式，使形態(tài)學(xué)運算速度大幅提高。實驗結(jié)果表明，隨著待處理數(shù)據(jù)規(guī)模的增加和處理復(fù)雜度的增加，GPU并行腐蝕運算的加速效果顯著增加。在本實驗平臺上GPU并行腐蝕運算可達(dá)到150倍的性能提升，目前，還未見任何優(yōu)化算法達(dá)到如此高的加速效果的報道。本文還在支持CUDA的較低檔的型號為GT220的顯卡上進(jìn)行了測試，在其上進(jìn)行并行腐蝕運算也可達(dá)到10倍以上的性能提升?？梢?，在數(shù)學(xué)形態(tài)學(xué)運算時，同時考慮性能和成本兩方面，GPU相對于CPU都有絕對的優(yōu)勢[12]。所以，在涉及數(shù)學(xué)形態(tài)學(xué)運算的領(lǐng)域，GPU并行運算是一個必然趨勢。

[1]Solille P.形態(tài)學(xué)圖像分析原理與應(yīng)用[M].北京：清華大學(xué)出版社，2008：47-48.

[2]陸宗騏，朱煜.數(shù)學(xué)形態(tài)學(xué)腐蝕膨脹運算的快速算法[C]//第十三屆全國圖象圖形學(xué)學(xué)術(shù)會議，2006-11-06，2006：306-311.

[3]孫即祥.圖像分析[M].北京：科學(xué)出版社，2005：183-184.

[4]Garland M，Grand S L，Nickolls J，et al.Parallel computing experiences with CUDA[J].Micro，IEEE，2008，28（4）:13-27.

[5]Nickolls J，Dally W J.The GPU computing era[J].IEEE Micro，2010，30（2）:56-69.

[6]Hawick K A，Leist A，Playne D P.Parallel graph component labelling with GPUs and CUDA[J].Parallel Computing，2010，36（12）:655-678.

[7]Ufimtsev J S，Schulten K.GPU-accelerated molecular modeling coming of age[J].Journal of Molecular and Modelling，2010，29（2）:116-125.

[8]吳恩華，柳有權(quán).基于圖形處理器（GPU）的通用計算[J].計算機(jī)輔助設(shè)計與圖形學(xué)學(xué)報，2004， 16（5）：601-612.

WU En-hua，LIU You-quan.General purpose computation on GPU[J].Journal of Computer-aided Design& Computer Graphics， 2004， 16（5）:601-612.

[9]Harish P，Narayanan P J.Accelerating large graph algrithms on the GPU using CUDA[J].Computer Science，2007（4873）:197-208.

[10]NVIDIA Corporation，NVIDIA C 3.1，NVIDIA CUDAC Programming Guid Version 3.1[S].2010.

[11]左顥睿，張啟衡，徐勇，等.基于GPU的快速Sobel邊緣檢測算法[J].光電工程，2009， 36（1）：8-12.

ZUO Hao-rui， ZHANG Qi-heng， XU Yong， et al.Fast sobel edge detection algorithm based on GPU[J].Opto-electronic Engineering， 2009， 36（1）:8-12.

[12]Seung I P，Ponce S P，HUANG Jing，et al.Low-cost，highspeed computer vision using NVIDIA's CUDA architecture[C]//2008 37th IEEE Applied Imagery Pattern Recognition Workshop，2008:1-7.