基于CUDA平臺的DR圖像增強(qiáng)處理加速算法

【作 者】何祥彬，周荷琴，李方勇

中國科學(xué)技術(shù)大學(xué)自動化系，安徽，合肥，230027

直接數(shù)字化X線攝影（Digital Radiography, DR）是現(xiàn)代醫(yī)療診斷中的一種先進(jìn)醫(yī)學(xué)成像技術(shù)，具有曝光劑量小、成像速度快和動態(tài)范圍大等優(yōu)點。但DR成像過程易受噪聲等的影響，使原始的DR圖像對比度達(dá)不到醫(yī)學(xué)診斷要求，需要對其進(jìn)行一定的增強(qiáng)處理。多尺度增強(qiáng)算法較適用于DR圖像的增強(qiáng)處理，其方法主要包括小波變換和在小波分析基礎(chǔ)上提出的塔型多尺度增強(qiáng)算法[1]。后者雖然簡化了圖像分解和重構(gòu)的計算過程，但由于DR圖像的分辨率極高，其算法運(yùn)行仍需要較長的時間，因此如何提高其運(yùn)行速度是能否將它應(yīng)用到實際DR系統(tǒng)上的關(guān)鍵問題。有作者在不影響圖像整體質(zhì)量前提下，通過去掉一些冗余的計算來提高塔型算法的速度[2]，但加速效果仍不甚明顯。隨著計算機(jī)技術(shù)的快速發(fā)展，不斷有更高性能的部件問世，如何深入挖掘它們的潛力，將它們用到新型醫(yī)學(xué)成像設(shè)備中去，也是生物醫(yī)學(xué)工程的課題。例如，CPU的多線程和并行計算功能已被成功應(yīng)用于圖形繪制和圖像處理，但由于操作系統(tǒng)的流水線作業(yè)調(diào)度方式和CPU本身在浮點運(yùn)算上的限制[3]，使基于CPU的多線程操作對加速效果仍不夠明顯。Intel奔騰處理器的單指令多數(shù)據(jù)(Single-Instruction Multiple Data, SIMD)技術(shù)，能讓一條指令同時對4個數(shù)據(jù)進(jìn)行相同的操作，我們小組曾將其用于基于CT圖像的體繪制中，明顯提高了繪制的速度[4]。然而，對于分辨率很高的DR圖像，這兩種技術(shù)均難以滿足速度的要求。為此，我們嘗試用GPU（Graphic Processing Unit，圖形處理器）來提高多尺度塔型算法的運(yùn)行速度，采用CUDA（Compute Unified Device Architecture, 計算統(tǒng)一設(shè)備體系結(jié)構(gòu)）架構(gòu)，將算法中包含的大量卷積運(yùn)算交給具有超強(qiáng)并行運(yùn)算能力的GPU來執(zhí)行，顯著提高了算法的運(yùn)行速度。

1 加速方法設(shè)計

1.1 CUDA架構(gòu)

GPU是專為計算密集型、高度并行化的算法而設(shè)計的處理器。它將更多的晶體管用于數(shù)據(jù)處理，而非數(shù)據(jù)緩存和流控制，具有極高的計算密度，對浮點數(shù)的計算能力遠(yuǎn)勝于CPU和其它處理器，已進(jìn)入計算主流[5]。CUDA是顯卡廠商N(yùn)VIDIA公司推出的一種基于GPU的軟硬結(jié)合的通用計算構(gòu)架，并提供硬件的直接訪問接口，不必靠圖像API接口來實現(xiàn)對GPU的訪問。它用CPU做總控制器，與GPU結(jié)合，以C為編程語言，并進(jìn)行了適度擴(kuò)展，能大量開發(fā)高性能計算指令。因此，能在GPU強(qiáng)大計算能力的基礎(chǔ)上，實現(xiàn)高效的密集數(shù)據(jù)計算，且非常適合并行計算[6]。

CUDA的運(yùn)行模型如圖1所示，在Host（CPU）上分配了不同的任務(wù)，一個個線程Kernel按照線程網(wǎng)格(Grid)的概念在Device（GPU）上執(zhí)行。這樣，每個Grid都得到資源，但不同Grid間相互獨立。每個Grid又可包含很多個線程塊(Block)，各Block之間不相互通信，也不并行工作，但共享同一個任務(wù)分配的資源。每個Block又把任務(wù)分配給多個線程(Thread)，每個Thread都有自己的編號，我們可以很方便地給它分配具體的任務(wù)，Thread之間共享資源，可相互通信，并行工作，并按既定的規(guī)則同步。每個Block下的線程均可讀寫自己的寄存器和本地內(nèi)存(local memory)和Block的共享(share)內(nèi)存，也可以讀寫Grid的全局(global )內(nèi)存、常量(constant)內(nèi)存和紋理(texture)內(nèi)存。

圖1 CUDA運(yùn)行模型Fig.1 The running model of CUDA

1.2 多尺度塔型增強(qiáng)算法的加速

多尺度塔型增強(qiáng)算法包括圖像分解、子帶增強(qiáng)和圖像重建等三個主要步驟，其流程如圖2所示。圖像分解是多尺度塔型增強(qiáng)算法不可缺少的部分，算法由上而下的分解過程，即為拉普拉斯金字塔的建立過程[7]。圖中LP代表低通濾波器，一般采用5×5的高斯低通濾波器。Fi-1是i級輸入圖像（i=1,2,……L+1），經(jīng)過濾波并以2步長進(jìn)行抽樣，得到第i級圖像的近似Fi。接著再對Fi進(jìn)行內(nèi)插和濾波，把得到的圖像與Fi-1相減，便可得到第i級的拉普拉斯金字塔Gi。在圖像分解過程中得到的每一級圖像Gi，代表著圖像中不同尺度的高頻信息，即圖像的細(xì)節(jié)信息。圖2中虛線框部分是子帶增強(qiáng)，它是塔型算法的關(guān)鍵步驟。采用非線性變換對分解得到的每一級圖像即不同尺度的細(xì)節(jié)信息進(jìn)行子帶增強(qiáng)。所選用的非線性變換函數(shù)直接影響著圖像的增強(qiáng)效果。子帶增強(qiáng)后，緊接著進(jìn)行圖像重建，它是圖像分解的逆過程，由下而上，將經(jīng)過增強(qiáng)后的圖像高頻信息相加，便可得到處理后的圖像。

DR圖像的分辨率很高，每幅圖片包含的像素一般在2000×3000左右，因此在多尺度塔型增強(qiáng)的分解和重構(gòu)的每一層上進(jìn)行高斯低通濾波時，都需要進(jìn)行大量的卷積運(yùn)算，這也是塔型算法耗時多的主要原因。若能提高這部分卷積運(yùn)算的速度，就會顯著縮短整個算法的運(yùn)行時間。

很多文獻(xiàn)都對多尺度塔形算法進(jìn)行了介紹，這里不再贅述，僅考察算法中涉及的包含大量運(yùn)算的空間域內(nèi)的二維卷積公式：

圖2 塔型多尺度增強(qiáng)算法流程Fig.2 Flow diagram of Laplacian pyramid multi-scale contrast enhancement

式中s(i,j)為像素點的值，k(i,j)為濾波器的值。由式（1）可知，對于n*n的二維濾波器，每個輸出點的計算需要作n*n次乘法及n*(n-1)次加法運(yùn)算。在塔型增強(qiáng)算法中常用的5×5的高斯低通濾波器為:

該濾波器可以由兩個濾波器相乘得到，即

式中h=[0.05 0.25 0.40 0.25 0.05]。因此，我們可以對該二維濾波進(jìn)行分解，先使用h濾波器進(jìn)行一次行濾波，再使用hT進(jìn)行一次列濾波。容易證明，采取這樣的濾波形式可將算法復(fù)雜度從原來的O(n2)降為O(n)，從而顯著減少計算量，縮短算法的運(yùn)行時間，而濾波效果又沒有改變。此外，先做行卷積，再作列卷積，每次讀取顯存都可以利用其尋址模式順序讀取，進(jìn)一步提高了GPU的運(yùn)算速度。

1.3 基于CUDA的加速

根據(jù)CUDA的運(yùn)行模型，將要處理的卷積任務(wù)劃分為多個 Block，Block中的每個線程執(zhí)行一個輸出點的卷積運(yùn)算。為了不執(zhí)行過多的內(nèi)存拷貝和邊界控制，我們在實現(xiàn)加速方法時選用紋理內(nèi)存來進(jìn)行線程的數(shù)據(jù)存取。下面介紹主要的實現(xiàn)步驟：

(1) 首先聲明紋理變量及紋理模式texture＜usigned short, 2, cudaReadModeElementType＞ texData；

(2) 開辟CUDA數(shù)組空間cudaArray*a_Data，并將圖像數(shù)據(jù)從顯存拷入CUDA數(shù)組空間 cudaMemcpy

ToArray (filData[i],0, 0, pCrt, crtS*sizeof(unsigned sh ort),cudaMemcpyDeviceToDevice)；

(3) 將CUDA數(shù)組綁定到紋理內(nèi)存cudaBind ToArray(texData, filData[i])。這樣就可以利用紋理內(nèi)存，使GPU快速的完成所有計算，最后將計算結(jié)果輸出到CPU。

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

我們用2560×3072的DR胸部圖片為實驗數(shù)據(jù)。實驗環(huán)境為CPU AMD 4200+，主頻4.4GHz，2GB內(nèi)存，Windows XP操作系統(tǒng)，編譯工具為Visual C++2005。GPU選用NVIDIA公司的GeForce 9500GT，驅(qū)動的版本為NVIDAI Driver 2.2。

實驗中采用文獻(xiàn)[8]提出的非線性變換公式作為子帶增強(qiáng)函數(shù)，分別在GPU和CPU下實現(xiàn)塔型算法。圖3是結(jié)果圖像，其中（a）為原始圖像，（b）、(c)分別為利用GPU和CPU分解得到的圖像細(xì)節(jié)信息，(d)、(e)分別為GPU和CPU下增強(qiáng)后的圖像。從圖3可以直觀地看出，無論是增強(qiáng)過程中的分解結(jié)果還是增強(qiáng)后的最終結(jié)果，在GPU和CPU下均是一致的，且處理后的圖像清晰度相對原圖均有較大的提高。

為了進(jìn)一步驗證在GPU下做塔型增強(qiáng)算法的正確性，對GPU和CPU下的運(yùn)算結(jié)果的每一像素點計算差，然后求平方和，結(jié)果為0，表明算法在GPU和CPU下執(zhí)行后得到的結(jié)果是一致的。

實驗統(tǒng)計了20次塔型增強(qiáng)算法的運(yùn)行時間，然后取平均值，得到表1所示的結(jié)果。可見，在CPU下執(zhí)行單線程對DR圖像實施增強(qiáng)，計算非常耗時；采用多線程方案可以稍微提高計算速度，但其加速效果不是非常明顯；而在基于GPU的CUDA平臺下，DR圖像的增強(qiáng)處理耗時僅為CPU單線程方式的1/30，加速效果非常明顯。

圖3 GPU和CPU下DR胸部圖像增強(qiáng)結(jié)果比較Fig.3 The comparison of processing results based on GPU and CPU

表1 CPU和GPU下塔型算法的運(yùn)行耗時Tab.1 The time-consuming of CPU and GPU

3 結(jié)論

圖像增強(qiáng)是DR系統(tǒng)圖像處理過程中的重要環(huán)節(jié)，由于DR圖像包含的數(shù)據(jù)量很大，因此提高增強(qiáng)算法的運(yùn)行速度直接關(guān)系到DR系統(tǒng)的性能。本文基于CUDA平臺，利用GPU來完成多尺度塔型增強(qiáng)算法中圖像的分解和合成過程中的大量卷積運(yùn)算，能夠在保持原有增強(qiáng)效果的基礎(chǔ)上，將算法運(yùn)行速度提高近30倍，為在微型機(jī)系統(tǒng)上快速處理DR圖像提供了一條新途徑，具有較大的實用價值。而且，所介紹的方法也可用于其他大容量數(shù)據(jù)處理的場合。

[1] P. Vuylsteke, E. Schoeters. Multiscale Image Contrast Amplification (MUSICATM)[J]. Proc SPIE Image Processing,1994, 2167:551 - 560.

[2] 李名慶, 高新波, 許晶. 多尺度塔型醫(yī)學(xué)圖像增強(qiáng)算法[J]. 中國生物醫(yī)學(xué)工程學(xué)報, 2006, 25(2):178 – 181.

[3] 宋曉麗, 王慶. 基于GPGPU的數(shù)字圖像并行化預(yù)處理[J]. 計算機(jī)測量與控制, 2009, 17(6):1169 – 1171.

[4] 袁非牛，諸葛斌，周荷琴，等. 基于奔騰SIMD和分割技術(shù)的快速體繪制[J]. 中國圖形圖像學(xué)報，2003，8（12）：1438 – 1443.

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

[6] CUDA-走向GPGPU新時代. 程序員, 2008, 10(3):36 – 39.

[7] (美)岡薩雷斯等. 數(shù)學(xué)圖像處理[M]. 北京：電子科學(xué)出版社,2004.

[8] Xuanqin Mou, Min Zhang. Nonlinear Multi-scale Contrast Enhancement for Chest Radiograph[J]. ICIP, 2008, 3184 –3187.