基于KNN的DSA图像去噪及GPU的快速实现

王光磊，裴晨辉，苑　昊，王　斌，刘秀玲

(河北大学 电子信息工程学院 河北省数字医疗工程重点实验室，河北 保定 071002)

基于KNN的DSA图像去噪及GPU的快速实现

王光磊，裴晨辉，苑昊，王斌，刘秀玲

(河北大学 电子信息工程学院 河北省数字医疗工程重点实验室，河北 保定 071002)

摘要:为快速地去除或减少DSA(Digital Subtraction Angiography)图像的噪声，对比评价KNN(K Nearest Neighbors)算法对高斯噪声、泊松噪声、斑点噪声、椒盐噪声4种噪声去除或减少的效果，帮助医生快速准确地为病人诊断疾病。提出的算法主要贡献在于构建了基于GPU(Graphics Processing Unit)的加速方法，使传统图像去噪的运算速度得到大幅提升。基于图像降质、图像还原过程建模，使用KNN算法对4种噪声去除或减少，并对算法做并行化处理，利用GPU加速实现去噪的过程。通过实验得出，KNN算法能较好地去除或减少高斯噪声、泊松噪声来还原DSA图像，使用CUDA(Compute Unified Device Architecture)编写可在GPU上运行的程序，利用GPU对1 024×1 024像素的24位深度的DSA图像去噪，平均渲染帧率能达到190.53 f/s(帧/秒)，较传统CPU(Central Processing Unit)串行，平均处理速度提高70.86倍。使用GPU加速能够快速地处理数据量较大、计算密集的DSA噪声图像，实现有效并且快速的高斯噪声去除，帮助医生精、准、快地诊断疾病。

关键词:DSA；KNN；GPU；CUDA；图像去噪

近几年，心脑血管疾病发病率呈逐年上升趋势，其危害已超过传染性疾病，成为我国威胁居民健康最大类疾病。最近的一次统计数据显示，中国心脑血管疾病死亡已占总死亡构成的41%。在医学界，对于心脑血管疾病能够准确诊断也是卫生工作人员和医疗技术人员一直孜孜追求的。然而在各类疾病诊断中，脑血管疾病的诊断比较困难。

1980年，美国威斯康星大学Mistretta小组、亚利桑纳大学Nadelman小组研制成功数字减影血管造影机，在血管疾病诊断方面得到了很好的应用，尤其是心脑血管疾病的诊断。自此，医学领域有了伤害性小、方便简洁、安全实用，并且影像清楚分辨率高的血管造影系统。然而，DSA发展到今天仍存在一些问题，比如它在采集图像数据的时候会受到一些干扰，从而会使图像产生噪声。改善图像质量可以通过添加滤光片、对近红外光敏感的探测器、改良CCD感光元件等改进前端设备的方法进行预处理的调节。由于前端预处理的局限性，为了提高系统的鲁棒性与可靠性，降低成本，需要在后端利用软件进行图像去噪处理。然而，对于医学成像又有其特点，数据量巨大，计算密集，还有对实时性的要求也很高，传统的串行计算很难满足这些要求。

近些年来，国内外研究人员在图像去噪方面做了许多探索与研究。现如今已有的去噪方法大类几十种，小类上百种，如直方图均衡算法、经典逆滤波法、基于Retinex的相关算法、小波及曲波变换方法、维纳滤波、空间域的核回归方法、非局部均值、基于参数估计的算法、双边滤波、基于概率框架的盲图像复原算法等常见方法[1-4]。BUADES A，DELEDALLE C，COUPE P，GUO Y等人首先提出和应用了非局部均值(Non-local Means， NLM)去噪算法[5-9]。时至今天，在NLM提出十年间里也得到了许多改进，如DOWSON N，DARBON J，ORCHARD J，LIU Y L，WANG J，WANG Z M等人将算法改进后提高了算法实现的速度[10-14]。CHEN Y，DONOHO D L，CHANG S G等人运用小波阈值算法实现了图像去噪[15-17]。ZHANG J等人将小波与双边滤波算法结合实现了医学超声图像的去噪[18]。自2007年NVIDIA推出CUDA以来，研究人员才逐渐重视到并行计算在医学领域存在着巨大潜力。美国通用电气医疗集团利用GPU实现了锥形束反投影与透视CT重建。西门子公司BATENBURG，SIJBERS等人利用GPU加速了CT重建速度，并将硬件成本压缩了上百倍；STONE等人使用GPU加快高级MRI重建速度14.35倍。VERONICA GIDEN等人利用GPU依照场景图对用于XIP的CUDA内核管线进行成功构建，并提高了构建速度[19]。WOLFGANG WEIN等人使用GPU成功实现了针对诊断成像和成像引导介入的自动CT超声配准。2014年德克萨斯大学由STEVE博士带领的小组，将GPU应用在癌症治疗领域，原本需要计算70个多小时才能够完成的复杂质子放射疗法缩短到10 s，缩短时间的同时还提高了计算精度。GPU在医学领域的应用近十年间得到了空前的发展，为人类医疗事业做出了很大的贡献，GPU在医学方面的重要性也不言而喻。

因此，综合考虑后，本文提出了一种针对DSA图像去噪的GPU加速方案，并很好地应用于KNN去噪算法。实验表明，使用GPU加速能够快速处理数据量较大、计算密集的噪声图像，实现有效的噪声去除。

1KNN算法及其实现

1.1DSA噪声

DSA系统最主要的组成部分是控制处理机和快速图像处理机，以及模拟数字转换器、存储器等电子器件，整个图像的摄制、存储、处理和传递都是用数字形式进行的。所以，在这其中不可避免地会引入图像的降质、模糊以及噪声的干扰，采集到的图像信息及特征会衰减、被掩盖，甚至出现模糊，严重影响了医生对病人病情的诊断。其中最值得引起关注的是，在所有噪声中高斯噪声的比重最大，也就是说高斯噪声对影像的影响最大。因此，若对高斯噪声能够较好地滤除或者减弱，对最终的成像将会有较大的意义，对医生准确诊断病人的病情帮助也是最大的。

1.2基于KNN算法的DSA图像去噪

1.2.1图像降质

图像增强的主要目的是要迎合人们的视觉感官，使人尽可能地从图像中得到更多的信息或令视觉效果更好，具有艺术性。如图像增强和图像去噪主要是以预先确定的目标来针对性地还原图像。尽管两者有相重叠的领域，但图像增强主要是主观处理，而图像去噪原则大部分是客观处理。图像去噪是利用退化现象的某种先验知识来复原一幅带噪声的图像，而不能够使图像有较大失真，最大程度地从客观上进行最真实的还原。因此，图像去噪是面向一个退化的模型，而采用降质过程的逆过程做出针对性的并且目的很明确的处理，来恢复图像。

通过建模[20-21]，降质过程可以看作一个退化函数和一个加性噪声项，对一幅输入图像f(x,y)进行处理，产生一幅退化后的图像g(x,y)

g(x,y)=H[f(x,y)]+η(x,y)

(1)

给定函数g(x,y)和退化函数H以及加性噪声项函数η(x,y)，从而图像去噪的目的便成为了获得一个预知的原图像f^(x,y)。原则上，这种预知应尽可能地接近原图像，而且退化函数H、加性噪声项函数η(x,y)所含信息越丰富，f^(x,y)，f(x,y)越相近。

g(x,y)=h(x,y)*f(x,y)+η(x,y)

(2)

式中：h(x,y)为H的空间表示。空间域的卷积和频率的乘法组成了一个傅里叶变换对，所以可以用等价的频率域来写出前面的模型，如图1所示。

图1　降质—去噪处理模型

假设H为线性、空间不变过程，便可以得到在空间域中表示的降质图像

G(u,v)=H(u,v)F(u,v)+N(u,v)

(3)

式中：用大写字母表示空降质间域中相应项的傅里叶变换。由于线性、空间不变性，所以退化函数H可以由卷积来建模。在此次建模中，假设H是恒等算子，并只处理由噪声造成的退化。

为了更好地观察KNN算法对DSA图像的去噪效果，分析算法的针对性、有效性、稳定性以及评价本算法对去噪的成效，故实验将人为针对性地引入高斯噪声、泊松噪声、斑点噪声、椒盐噪声4种噪声来实验。这些噪声都是空间噪声，其值是随机数，由概率密度函数p(z)或等价的积累分布函数F(z)表征，它们均由均值m和方差σ2决定。a，b均为常数。这4种噪声的各参数如下：

1)高斯噪声

概率密度函数

(4)

积累分布函数

(5)

均值和方差

m=a，σ2=b2

(6)

2)泊松噪声

概率密度函数

(7)

积累分布函数

(8)

均值和方差

m=a，σ2=a

(9)

3)斑点噪声

概率密度函数

(10)

积累分布函数

(11)

均值和方差

(12)

4)椒盐噪声

概率密度函数

(13)

积累分布函数

(14)

均值和方差

m=(0)Pp+k(1-Pp-Ps)+(2n-1)Ps

σ2=(0-m)2Pp+(k-m)2(1-Pp-Ps)+

(2n-1-m)2Ps

(15)

1.2.2KNN算法实现

针对一幅退化图像的去噪，最简单方法是在上文介绍的模型中忽略噪声项，并形成如下形式的估计

(16)

式中：G(u,v)是降质图像的傅里叶变换，然后取F^(u,v)的傅里叶反变换得到图像的相应估计，来直接对噪声图像进行逆滤波，考虑噪声时，可以将估计表示为

(17)

KNN算法用来减少图像噪声也是基于高斯噪声滤波器，但比其更复杂。假设让u(x)作为最原始的噪声图像，KNNh,ru(x)是经过KNN滤波器和参数h和r之后的结果。让Ω(P)在空间上以一定的距离环绕像素点P，本文将它视为N×N大小的一系列像素，此时N=2M+1(自然数M=0，1，2，3…)，这样P就是Ω(P)的中心像素，于是

(18)

式中：C(x)是正常化系数。由公式得知，周围像素权重决定于h、r和Ω。r是一个传统的高斯噪声系数，由带噪图像决定，在用KNN算法处理时，r=N。h的选择具有主观性，它是用户自己来选择的，选择的不同最后的结果也会不同，用户可以根实验结果来调节h。对难以量化的h值，会根据图像默认一个值，最后的结果也是较为可观的。Ω的大小取决于图像的大小，为了达到一个良好的视觉效果，通常是5×5或7×7块的像素块儿。在本实验中，选择了7×7块的像素块进行实验，如图2所示。

图2　7×7像素块

1.3GPU加速实现

1.3.1耗时分析

在KNN算法中，K近邻搜索方法是暴力搜索法。本实验要处理的图像是1 024×1 024像素点的RGB图像，也就是说一幅图像是一个1 024×1 024×3的彩色像素组，其中每个彩色像素是一个三值组，这三个值分别对应一个特定空间位置处该RGB图像的红、绿、蓝分量，如图3所示。

图3　RGB分量示意图

RGB图像可以看作是三色图像的合成，把它们输入显示器对应的三端时，会在屏幕上生成一幅人类可视的多彩图像。其合成的三幅图像分别称作红色、绿色、蓝色分量图像，每个分量的数据类型又是判断三色分量值大小唯一的依据。本实验所采用的图像是unit8类的，其灰度值为0～255。同时用来表示这些分量图像像素值的比特数决定了一幅图RGB图像的比特深度，所以本实验的RGB图像的比特深度为24。因此，本实验的计算量非常大，传统的串行计算模式难以满足对速度的要求，为了提高计算效率，采用了GPU加速。

1.3.2GPU加速实现过程

在本实验中，为针对每个参考点的计算分配了一个线程，以发挥GPU优势。一个线程块分配了1 024个线程，大于192个线程，这样就可以在寄存器和共享内存内存取数据，1个时钟周期的延迟也可以忽略掉了。无论CPU还是GPU在计算的时候主要是存取数据耗时，对于GPU而言，若一个线程块有192个线程，正好可以在32个线程从启动到结束的时候，其他5组32个线程会有一组准备好数据，这样就可以无延时的进行计算，从而最大程度上利用了硬件资源。根据实验选取DSA图像的大小正好分配1 024个线程块，这样一个SM(Streaming Multiprocessor)上能够运行多个线程块，又线程数目大于了192个，这样延迟被忽略掉，从而最大程度上利用共享内存。参考点与线程对应图如图4所示。

实验时，首先将数据传输到GPU中。然后再将一个SM里需要计算的线程的数据放进共享内存中，共享内存再将运行线程直接使用的数据放进寄存器里面，来提高计算速度，最大程度挖掘GPU硬件性能。GPU端计算完毕后，再将计算结果传回到CPU端，这样就实现了GPU处理带噪声的DSA图像。其实现过程图如图5所示。

2实验结果

2.1实验环境

硬件方面，本文选用英特尔公司出品的专门面向服务器和工作站的高端处理器Xeon E5-2680，同时，为了更好地使实验效果对比明显，在原有的显卡NVIDIA Quadro K5200的基础上再加一块NVIDIA Tesla K40c。

图4　参考点与线程的对应图

图5　GPU实现过程

软件方面，系统选择64位的Windows 7专业版，CUDA平台选择Visual Studio 2013 + CUDA 7.0这样较新的组合搭配。环境的主要参数如表1。

表1实验环境参数

环境参数电脑型号戴尔PrecisionTower7910Tower操作系统Windows7专业版64位SP1(DirectX11)软件平台VisualStudio2013+CUDA7.0MATLABR2014a处理器英特尔Xeon(至强)E5-2680v3@2.50GHz(X2)主板戴尔0215PR(英特尔Haswell-EDMI2)显卡NVIDIAQuadroK5200(4Gbyte)加速卡NVIDIATeslaK40c(12Gbyte)内存192Gbyte(三星DDR42133MHz)主硬盘LSI(512Gbyte)(X4)

2.2去噪效果

由于本实验在数据选取时选择的DSA图像较大，通常图像的细节是医生诊断疾病的关键，故本文将放大图像的细节来做对比。截取部分如图6所示。

图6　DSA图像细节

实验分别引进了两组高斯噪声(m=0，σ2=0.05和σ2=0.03)、一组泊松噪声、一组斑点噪声(σ2=0.2)、一组椒盐噪声(d=0.03，d代表噪声密度)。其5组实验结果如图7所示。

从实验结果来看：KNN算法对选取的两组带高斯噪声的DSA图像具有较好去除效果，对泊松噪声在细节方面去噪较好，对斑点噪声几乎无任何功效，对椒盐噪声也有一定的噪声消除成效。

2.3加速效果

实验在对KNN去噪算法用CUDA实现后，做了一个对比实验，以展示算法用GPU运算的优势。实验分别选取了5位病人的DSA图像，每位病人选取2张，共10张，以准确测算去噪速度。GPU和CPU对比实验数据如表2所示。

从5个病人的10组对比实验结果来看，GPU的处理速度明显高于CPU的处理速度。其中，用GPU计算平均渲染帧率为190.53 f/s(帧/秒)，用CPU计算平均渲染帧率为2.69 f/s(帧/秒)，平均加速比为70.83，即应用本文的GPU加速策略可以将计算速度提高70倍以上。因此，应用本文的方法通过GPU的处理速度可以满足医学成像对实时性的要求，实现医生对病人的实时监测诊断。

e　椒盐噪声去噪前后图7　5组DSA图像去噪结果

比较项图像1图像2图像3图像4图像5图像编号se003se004se028se040se019se021se011se016se006se009图像规格1024×1024GPU渲染帧率/(f·s-1)190.9190.8190.6189.7190.6190.3190.4190.7190.8190.5CPU渲染帧率/(f·s-1)2.72.62.72.72.82.72.72.72.62.7加速比70.7073.3870.5970.2668.0770.4870.5270.6373.3870.56

3结论

DSA图像在采集过程中不可以避免地会引入一些噪声，为解决此问题，本实验应用了KNN算法来进行去噪。但由于DSA图像数据量大、计算密集，又考虑到实际情况，去噪的实时性也要好，即使对KNN算法进行改进仍很难满足这些要求。于是可以将解决的重点放在处理速度上，考虑到KNN算法的特点，可并行程度较高，使用CUDA编写可在GPU上运行的程序，用GPU处理后，效果较为理想。

参考文献：

[1]AI D， YANG J， FAN J， et al. Denoising filters evaluation for magnetic resonance images [J]. Optik-international journal for light and electron optics，2015(23)：3844-3850.

[2]SINGH C， AGGARWAL A. An efficient approach for image sequence denoising using Zernike moments-based nonlocal means approach[J]. Computers & electrical engineering，2015(11)：56-60.

[3]IFTIKHAR M A，JALIL A，RATHORE S，et al. An extended non-local means algorithm： application to brain MRI[J]. International journal of imaging systems and technology，2014(24)：293-305.

[4]KUANG Y，ZHANG L，ZHANG Y. An adaptive rank-sparsity K-SVD algorithm for image sequence denoising[J]. Pattern recognition letters，2014(45)：16-54.

[5]BUADES A，COLL B，MOREL J M. A review of image denoising algorithms，with a new one [J]. Multiscale model simul，2005，4(2)：490-530.

[6]BUADES A，COLL B，MOREL J M. Nonlocal image and movie denoising [J]. International journal of computer vision， 2008，76(2)：123-139.

[7]DELEDALLE C，DENIS L，TUPIN F. Iterative weighted maximum likelihood denoising with probabilistic patch based weights[J]. IEEE transactions on image processing，2009，18(12)：2661-2672.

[8]COUPE P，HELLIER P，KERVRANN C，et al. Nonlocal means-based speckle filtering for ultrasound Images [J].IEEE transactions on image processing，2009，18(10)：2221-2229.

[9]GUO Y，WANG Y，HOU T. Speckle filtering of ultrasonic images using a modified non local-based algorithm [J]. Biomedical signal processing and control，2011(6)：129-138.

[10]DOWSON N， SALVADO O. Hashed nonlocal means for rapid image filtering [J].IEEE transactions on pattern analysis and machine intelligence，2011，33(3)：485-499.

[11]DARBON J，CUNHA A，CHAN T F，et al. Fast nonlocal filtering applied to electron cryomicroscopy[C]//Proc. IEEE International Symposium on Biomedical Imaging： From Nano to Macro， 2008. [S.l.]：IEEE， 2008： 1331-1334.

[12]ORCHARD J，EBRAHIMI M，WONG A. Faster nonlocal-means image denoising using the FFT[J]. IEEE transactions on image processing，2008(6)：45-49.

[13]LIU Y L，WANG J， CHEN X，et al. A robust and fast non-local means algorithm for image denoising [J]. Journal of computer science and technology， 2008，23(2)：270-279.

[14]WANG Z M，ZHANG L. An adaptive fast nonlocal image

denoising algorithm[J]. Journal of image and graphics，2009，14(4)：669-675.

[15]CHEN Y，HAN C. Adaptive wave let threshold for image denoising[J]. Electronics letters，2005，41(10)： 586-587.

[16]DONOHO D L. Denoising by soft-thresholding [J]. IEEE transactions on information theory，1995，41(3)： 613-627.

[17]CHANG S G，YU B，VETTERLI M. Adaptive wave let thresholding for image denoising and compression[J]. IEEE transactions on image processing，2000，9(9)： 1532-1546.

[18]ZHANG J，WANG C，CHENG Y. Despeckling for medical ultrasound images based on wavelet and bilateral filter [J]. Journal of image and graphics，2014，19(1)：126-132.

[19]WOLFGANG W，SHELBY B，ALI K，et al. Automatic CT-ultrasound registration for diagnostic imaging and image-guided intervention[J]. Medical image analysis，2008(12)：577-585.

[20]FU B，LI W W，FU Y P. An image topic model for image denoising. [J]. Neurocomputing，2015(169)：119-123.

[21]HE N，WANGA J B，ZHANG L L，et al. An improved fractional-order differentiation model for image denoising. [J]. Signal processing，2015(112)：180-188.

责任编辑：时雯

DSA image denoising based on KNN and rapid implementation of GPU

WANG Guanglei，PEI Chenhui，YUAN Hao，WANG Bin，LIU Xiuling

(KeyLaboratoryofDigitalMedicalEngineeringofHebeiProvince，CollegeofElectronicandInformationEngineering，HebeiUniversity，HebeiBaoding071002，China)

Abstract:The purpose of this paper is to efficiently reduce or eliminate the noise of DSA (Digital Subtraction Angiography) image and comparatively evaluate the efficiency of noise reduction of including Gaussian noise， Poisson noise， Speckle noise and Salt and Pepper noise based on KNN (K Nearest Neighbors) algorithm. The results can be used to provide high quality image to doctors for quick and accurate disease diagnosis for the patient. The main contribution of the proposed paper constructs an acceleration structure using GPU， which can dramatically reduce the time cost of the traditional image denoising algorithm. The process modeling of image quality and image reconstruction is supplemented in KNN algorithm to reduce or eliminate the four kinds of noises. Additionally， the parallel processing based on GPU (Graphics Processing Unit) is executed to accelerate the denoising procedure. The experiments prove that Gaussian noise， Poisson noise can be removed or reduced effectively to reconstruct DSA image by KNN algorithm. The developed programming codes using CUDA (Compute Unified Device Architecture) works well in GPU to realize the average rendering frame rate 190.53 f/s for denoising of 1 024×1 024 DSA image with 24 bit depth. Compare with the traditional CPU (Central Processing Unit) serial， the average processing speed rate is increased by 70.86 times. Using GPU acceleration can quickly process DSA noise image with large amount of data and intensive computational mission. The results can provide accurate and fast diagnosis to doctors with the efficient Gaussian noise reduction.

Key words:DSA；KNN；GPU；CUDA；image denoising

中图分类号：TN911.73

文献标志码:A

DOI：10.16280/j.videoe.2016.06.003

基金项目：国家自然科学基金项目(61473112;61203160);河北省教育厅项目(QN2014166)；河北省自然科学项目(F2015201196)

收稿日期：2015-12-10

文献引用格式：王光磊，裴晨辉，苑昊，等.基于KNN的DSA图像去噪及GPU的快速实现[J].电视技术，2016，40(6)：10-16.

WANG G L，PEI C H，YUAN H, et al.DSA image denoising based on KNN and rapid implementation of GPU[J].Video engineering，2016,40(6)：10-16.