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

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

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

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

1 膨脹腐蝕的基本理論

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

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

2 GPU簡(jiǎn)介

2.1 GPU硬件架構(gòu)

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

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

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

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

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

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

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

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

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

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

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

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

3.2 實(shí)現(xiàn)過程

上面已經(jīng)提到CUDA編程模型是GPU和CPU協(xié)同工作的，所以GPU腐蝕運(yùn)算的實(shí)現(xiàn)過程是CPU和GPU共同完成的。CPU負(fù)責(zé)內(nèi)存分配、Grid和block的維度設(shè)定、啟動(dòng)在GPU上運(yùn)行的kernel函數(shù)并向其傳遞參數(shù)、最后從GPU獲取處理結(jié)果。GPU負(fù)責(zé)大量并行運(yùn)算，并將處理數(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è)計(jì)grid和block的維度。為了提高性能將 block 維度固定為 dim3 blockDim（16，16），則 grid的維度則根據(jù)圖像寬度width、高度height和block的維度設(shè)定，可設(shè)定為 dim3 gridDim（ （width+15）/16，（height+15）/16）；

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

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

4 實(shí)驗(yàn)結(jié)果及分析

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

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

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

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

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

5 結(jié) 論

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

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

[2]陸宗騏，朱煜.數(shù)學(xué)形態(tài)學(xué)腐蝕膨脹運(yùn)算的快速算法[C]//第十三屆全國圖象圖形學(xué)學(xué)術(shù)會(huì)議，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].計(jì)算機(jī)輔助設(shè)計(jì)與圖形學(xué)學(xué)報(bào)，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邊緣檢測(cè)算法[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.