基于非局部均值的低剂量CT图像处理并行化优化

摘要

低剂量CT (LDCT)图像常因斑驳噪声/伪影的严重增加而严重退化，导致临床诊断精度降低。非局部均值(non - local means, NLM)滤波利用LDCT图像中的大尺度斑块相似性信息，可以有效地去除斑点噪声/伪影。但NLM滤波在LDCT成像中的应用也需要较高的计算成本，因为通常需要在大搜索窗口内进行密集的斑块相似度计算，以包含足够的结构相似度信息来抑制噪声/伪影。为了提高其临床可行性，本研究进一步优化了NLM滤波的并行化，避免了行强度计算和对称权计算的重复计算。共享内存与快速速度也用于提出的方法的行明智强度计算。定量实验表明，可以相对于传统的直像素明智的平行化实现显着的加速度。

1.介绍

x射线计算机断层扫描(CT)可以在毫米级反映人体衰减图，为临床诊断提供丰富的组织、器官或病变的三维信息。虽然CT在临床上应用广泛，但CT检查时对患者的辐射一直是人们广泛关注的问题。据[1CT辐射可能会增加甚至癌症发育代谢异常的风险。降低辐射剂量的最实用的手段是减少管电流（Milliampere（MA））或管电流时间产品（毫安第二（MAS））。但是，降低MA / MAS设置通常会导致具有增加的斑点噪声和条纹伪像的CT图像[2那3.，会影响诊断的准确性[4.-7.］．研究人员可以通过开发新的重建或后处理算法来抑制低剂量CT (LDCT)图像中的噪声和伪影。目前提高LDCT图像质量的方法大致可以分为三类:预处理方法、迭代重建方法和后处理方法。

第一种是指在常规FBP重建之前，通过抑制投影原始数据中的噪声来改善CT成像的技术。这些技术的关键是找到投影数据的准确统计分布，并设计有效的恢复算法[5.那6.］．第二种是迭代重建方法，将LDCT成像视为不适定逆问题，通过一些迭代优化解将问题作为先验正则化代价函数来求解[7.那8.］．虽然迭代重建能有效地获得良好的重建图像质量，但迭代重建最广为人知的限制是迭代优化需要大量的计算。另外，出于专利保护的考虑，目前主流CT设备供应商通常没有提供格式良好的投影数据，这严重限制了这两个研究方向的研究和可能的临床应用。

第三是指后处理方法，可以直接应用于改善LDCT图像。CT图像中的噪声，伪像和正常组织的分布和规模特征需要在设计有效的后处理算法中共同考虑[9.那10］．它被指向[9.-14指出利用信息冗余特性的非局部均值(non - local means, NLM)滤波可以有效地抑制噪声和伪影，而不使图像细节明显模糊。我们还将注意到NLM中的patch相似性度量也被用于构建层析重建的正则化项[15那16］．

然而，由于在LDCT图像中噪声和伪影往往以显著的幅度分布，抑制噪声和伪影实际上需要一个较大的搜索窗口来包含更多的结构信息，这意味着计算代价很大。鉴于目前放射科工作量大，这将极大地限制其临床应用。为了克服这一问题，本文提出了一种改进的基于gpu的并行化方法来加速NLM滤波。该方法通过逐行强度计算和权值计算避免重复计算，优化了NLM滤波中的计算。快速GPU共享内存的数据访问速度也被很好地利用来减少数据操作成本。在二维LDCT图像上的实验结果表明，改进的并行化方法可以显著缩短计算时间，使其成为一种潜在的适用于LDCT成像的处理方法。

2.基于非局部均值的低剂量CT图像处理

与基于强度梯度信息的恢复算法相比，NLM滤波可以在不模糊图像结构的情况下抑制边缘保持的噪声/伪影。在NLM滤波中，一个图像小块与一组邻近大区域的相似小块进行匹配，这样可以利用邻近大尺度下更多的结构相似信息来抑制LDCT图像中的噪声和伪影。NLM算法用搜索窗口内强度的加权平均值来代替像素强度．每个权重表示搜索窗口中中心像素与相邻像素之间的相似度，由这两个像素周围的patch之间的欧氏距离计算。让 表示待处理像素，令表示搜索邻域窗口中的像素，令表示处理后的图像，令表示待处理的图像;二维NLM滤波算法可以表示为[17]： 在哪里为中心的搜索窗口;表示以为中心的两个patch的相似度和，分别为半径;为距离相关的高斯核函数;patch中的像素数为．我们通常使用这个参数在(2)来控制平滑效果。

在图中1-4.，图中给出了4幅二维LDCT图像的NLM滤波处理结果1（a）那2（a）那3(一个),4(一)图中给出了两张相应的标准剂量CT (SDCT)图像1（b）那2 (b)那3 (b),4 (b)作为引用。分别用降低的管电流80 mA和常规管电流240 mA采集LDCT和SDCT图像。我们可以看到，CT图像主要由强度范围有限的像素组成，代表不同组织的强度分布在整个图像域。其他扫描参数均为默认设置。与参考SDCT图像相比，我们可以看到管电流降低会导致LDCT图像噪声和伪影严重增加。数据1（c）那2 (c)那3 (c),4 (c)说明了搜索窗口大小为的NLM滤波处理结果，补丁大小为（),和参数设置为10。数字5.通过对一组二维胸部LDCT图像进行处理，显示三维体积的结果。图中的插图1-4.在合适的窗口显示。所有参数均在放射科资深医生指导下设置。可以看出，NLM滤波可以有效地抑制LDCT图像中的斑纹噪声和伪影，而不会导致明显的结构模糊。

要突出大型搜索窗口的重要性，我们还列出了图中的1（d）NLM处理的结果搜索窗口(其他参数设置为与图中的结果相同1（c）)．我们可以看到这个较小的处理搜索窗口没有给出满意的工件抑制(见箭头)。因此，一个大的搜索窗口(到)，以包含足够大的相似度信息来抑制LDCT图像中的噪声和伪影，如图所示1（c）[10］．但是，在较大的搜索窗口内进行patch相似度计算，往往会带来较大的计算量。对于一个图像大小，搜索窗口的像素数为patch半径为，我们得到了计算复杂度对于原始的基于CPU的串行处理，其总计算复杂度达为大小的图像。这种计算成本过高，无法为放射科常规提供实时CT成像;因此，我们需要加快NLM的滤波速度，以便在临床中得到快速的应用。

3.基于cuda的NLM算法GPU加速

3.1。基于CUDA的GPU加速简介

利用基于GPU的技术实现并行算法已经成为并行计算领域的一个显著趋势。基于GPU的并行化是通过对可并行化的原始网格计算任务中的粗尺度补丁和细尺度线程联合并行来实现的[18-20.］．CUDA(计算统一设备架构)技术为开发者提供了一个软件平台，可以使用c风格代码设计并行任务，直接访问虚拟指令集和GPU内存。每个运行在CUDA上的并行化函数被称为内核，我们使用表示基于[的连接的并行级联功能20.］．核函数的输出的输入．我们使用表示核函数的输入数据在处理和表示功能的输出数据．我们可以用(5.)表示内核函数，如下所示: 在哪里表示图像中的像素位置。需要指出的是，在某些情况下，并非所有的输入数据都需要更新(例如，)．

３．２．基于GPU并行化的NLM滤波算法

传统的GPU并行化通过直接逐像素并行化加速NLM滤波算法。基于上述(1） - （4.)，我们通常将算法分为以下四个部分(6.） - （9.），在循环中计算。循环的数量被设置为搜索窗口大小遍历窗口中所有相邻的点．第一个核函数(6.)通过GPU并行计算强度差异，具有计算复杂度．这里我们用并行运算时间来量化计算复杂度。在(6.)，表示相邻像素在以为中心的搜索窗口中的空间位置，也可以用空间位置表示．我们设置初始化。考虑

对数据，第二个核函数使用(7.)，根据(6.） 在．我们可以看到，第二个核函数的计算复杂度大约是因为有对两个比较块中的每个像素对进行加权求和操作。考虑

第三个核函数(8.）计算权重和强度的总和和，计算复杂度是对于这个操作。考虑

在最后一个循环中，中的最后一个内核函数9.)用于计算最终输出图像考虑 这里表示最后一个循环编号。运算的计算复杂度(9.)也．最终图像输出为．从(6.） - （9.)，可以看出，传统并行化算法的整个计算复杂度为．

3．3．改进的GPU加速NLM滤波算法

在上述的常规并行化方法中，(7.）串行应用以计算贴片相似度，这导致使用大搜索窗口时的大计算成本。因此，我们的第一个改进致力于降低本部分的计算复杂性。数字6.说明了补丁的大小（)，红点表示中心点。方程(1） - （4.）显示，NLM滤波中的贴片相似度可以通过两个补丁中的对应像素的相应像素的加权和的加权之和来量化。在图中6.，可以看出，对于两个斑块中位于绿色点的中心点，同行蓝色点的强度差之和与中心点向下移动到红色点时的强度差之和是相同的。这意味着当中心点移动时，行之间的强度差被反复计算像素的距离。因此，我们可以通过以下逐行计算来有效地计算patch差异: 两个独立像素之间的强度差在哪里和表示搜索窗口中相邻的点位置。的数据集垂直方向上的不同点。因此，对于大小的补丁，我们从(10),取值为不同的行对。在存储到共享内存中之前，只需要计算一次行差，然后再计算另一次 业务(10)可以很容易地通过从共享内存加载数据，然后执行高斯加权来获得．对于具有快速单精度浮点处理的GPU，主要计算成本为10）在于全局存储器的数据访问操作，因为与全局存储器访问相比，共享内存访问的时间成本是微不足道的。（计算复杂性）（11）可以大致估计[21］．

与传统的GPU并行化相似，我们也将算法分为以下四个部分(11） - （14)和在循环中计算。假设输入图像的大小，我们设置的大小是．的数据那那那是尺寸的,表示居中的搜索窗口中的相邻点．初始化也是用data设置的那．

第一个核函数计算每个行对的强度差之和，乘以基于行到中心点的垂直距离计算的高斯权值。这个核函数的计算复杂度是．考虑

第二个核函数基于(12)．在(12)，我们通过将第一个核函数(11)．核函数的计算复杂度(12） 是．考虑

第二个改进是利用()中计算权重的对称性，节省了一半的计算成本。2)．显然,我们有（表示像素的位置偏移量在搜索窗口中)。基于这种对称性我们也在积累累积加权强度时的位置．这样，我们只需要遍历搜索窗口中一半的像素。第三个核函数为(13)，其计算复杂度为．考虑

然后，最后一个内核函数(14)，可得到最终处理后的图像:

现在，关于搜索窗口所需的操作号的最终循环编号变为（表示搜索窗口的半径），大约是．最终输出图像为．综上所述，改进算法的总计算复杂度约为，这近似等于．我们可以看到计算复杂度已经降低到相对于传统的并行化。

4.实验结果与分析

在本节中，我们比较了不同方法的计算成本。为了验证所提出的加速方法对NLM滤波带来的改进，我们进行了同样的处理LDCT图像见图1（a）使用串行算法（基于CPU），传统的并行化算法（基于GPU）和改进的并行化算法（基于GPU）。在本节中，我们没有说明处理的图像，因为与图中相同的图像1（a）被获得。我们将补丁大小设置为并记录相对于搜索窗口大小的计算时间。数字7.图示了串行算法和传统并行化算法的计算时间。我们可以观察到传统的并行化通过直的像素明智的并行化显着降低计算成本，并实现了原始串行算法的100多次的加速度。我们的实验的系统配置如下给出。

4．1.硬件环境

CPU: Inter(R) Core(TM) i7-3770 @ 3.40 GHz;记忆:8 GB;显卡:NVIDIA GeForce GTX 680, 1536 CUDA核;有效存储时钟:6008 MHz;内存带宽:192gb /s;内存:2gb;内存总线类型:256位。

4．2．软件环境

操作系统:Win7 64位;Matlab: R2011a;CUDA: 4.0。

然后比较了传统并行算法和改进并行算法在搜索窗口大小上的计算时间。补丁大小固定为．如图所示8.，当搜索窗口大小变得大于，加速度比近似等于，这一观察结果与上述推导的加速度比是一致的．此外，我们比较了传统并行算法和改进并行算法在patch大小上的计算时间。搜索窗口大小固定为．由于较大的patch大小往往会导致图像模糊，因此我们设最大patch大小为．我们可以从图中看到9.当补丁尺寸增加时，加速度的明显增量，并且再次验证上述推导率．

5.讨论与结论

本文进一步优化了CT图像处理中NLM滤波的并行化。该方法通过逐行强度计算和权值计算避免重复计算，优化了NLM滤波中的并行计算。快速GPU共享内存的数据访问速度也得到了很好的利用。我们将改进算法应用于LDCT图像处理，发现改进算法比传统的并行化算法有显著的加速比。现在，处理一个需要大约0.8秒CT图像与搜索窗口,patch，并发现该参数设置在NLM滤波中能够提供有效的处理。本文只提供2 d NLM过滤结果,我们会强调,同样的并行策略可以很容易地扩展到更多的计算量的3 d加速NLM过滤,和相同的加速度比2 d情况下可以预期的,因为他们有相同的计算结构。具体来说，这个扩展可以通过替换(11)通过平面优化。然而，我们也需要指出，对于大层厚(>2 mm)的CT切片，由于层间连续性较差，不建议采用三维NLM滤波处理。

目前，NLM中的结构相似思想在图像处理的其他领域(如图像分割和图像重建)得到了广泛的应用[15那16那22那23］．本文提出的并行优化方法可直接用于加速这些应用中的斑块相似度计算。在目前的并行化方法中，反映patch相似度的权值是通过(7.)，可以通过内核间的操作进一步并行化，以实现进一步的加速。将多核CPU策略与GPU并行化技术相结合，提高计算速度。这种优化策略可以很容易地使用补丁相似类型度量来加速其他重建或恢复任务[24-26］．我们也可以考虑通过将分数阶度量加入到patch相似度的计算中来提高NLM滤波的性能[27］．还需要评估拟议的处理可能带来的分割/配准(与CT图像相关)的潜在精度提高[28-30.］．这些问题都将在今后的工作中加以解决。

利益冲突

作者声明本文的发表不存在利益冲突。

致谢

辽宁省科技攻关项目(no . 2013225089);国家自然科学基金资助项目(no . 81370040, no . 31100713);江苏省青兰项目。国家重点基础研究发展计划(no . 2010CB732503)资助。

参考文献

1. 张刚，孙东，闫平，赵慧，李志强，“基于单尺度retinex理论的LDCT图像对比度增强算法”国际计算智能建模控制与自动化会议论文集，第1282-1287页，奥地利维也纳，2008年12月。视图:出版商的网站|谷歌学术搜索
2. M. K. Kalra, M. M. Maher, T. L. Toth等人，“CT辐射剂量优化策略”，放射学号，第230卷3，页619-628,2004。视图:出版商的网站|谷歌学术搜索
3. E. Angel, N. Yaghmai, C. M. Jude等，“蒙特卡罗模拟评估多探测器CT中管电流调制对乳腺剂量的影响，”医学与生物学的物理学第54卷第5期3，第497-511页，2009。视图:出版商的网站|谷歌学术搜索
4. R. NELSON，“肺癌的低剂量CT筛选产生高误报率”第45届美国临床肿瘤学会年会(ASCO '09)论文集, 2009年。视图:谷歌学术搜索
5. 王建军，吕海华，温建文，“基于多尺度加权最小二乘正弦图的低剂量x射线计算机断层成像恢复”，生物医学工程上的IEEE交易，第55卷，第55期3，第1022-1031页，2008。视图:出版商的网站|谷歌学术搜索
6. T. Kubo, Y. Ohno, S. Gautam等，“在肺部CT中使用3D自适应原始数据过滤器:对减少辐射剂量的影响”，美国x射线杂志号，第191卷。4, pp. W167-W174, 2008。视图:出版商的网站|谷歌学术搜索
7. “基于自适应加权非局部先验的低剂量x射线计算机断层成像贝叶斯统计重构，”计算机化的医学影像和图形第33卷第3期7，第495-500页，2009。视图:出版商的网站|谷歌学术搜索
8. I. A. Elbakri和J.A.Fessier，“X射线计算机断层扫描中的迭代图像重建的高效和准确的可能性”医学影像2003:图像处理，第5032卷SPIE会议记录，页1839-1850，圣地亚哥，加利福尼亚州，美国，2003年2月。视图:出版商的网站|谷歌学术搜索
9. 陈勇，杨振宇，胡玉华等，“基于伪影抑制的胸椎低剂量CT图像处理，”医学与生物学的物理学(第57卷)9, pp. 2667-2688, 2012。视图:出版商的网站|谷歌学术搜索
10. 陈勇，陈伟，尹旭东等，“基于加权强度平均的低剂量腹部CT图像改进，”欧洲放射学杂志，第80卷，第2期。2，页42 - e49, 2011。视图:出版商的网站|谷歌学术搜索
11. F. P. X. de Fontes, G. A. Barroso, P. Coupé，和P. Hellier，“实时超声图像去噪”，实时图像处理杂志，第6卷，第2期1，pp。15-22，2011。视图:出版商的网站|谷歌学术搜索
12. Xu和K. Mueller，“一种用于NLM滤波正则化CT重建的参考图像数据库方法”全部3D的诉讼程序，第116-119页，2011。视图:谷歌学术搜索
13. 贾旭东，田振中，娄烨。Sonke, S. B. Jiang，“基于时间非局部方法的四维锥形束CT重建和增强”，医学物理学，卷。39，没有。9，pp。5592-5602,2012。视图:出版商的网站|谷歌学术搜索
14. “基于局部噪声水平自适应非局部均值滤波的CT去噪”，“基于局部噪声水平自适应非局部均值滤波的CT去噪”，医学物理学号，第41卷。1、文章编号011908,2014。视图:出版商的网站|谷歌学术搜索
15. “基于自适应加权非局部先验的低剂量x射线计算机断层成像贝叶斯统计重构，”计算机化的医学影像和图形第33卷第3期7，第495-500页，2009。视图:出版商的网站|谷歌学术搜索
16. 陈勇，马军，冯青，罗磊，石鹏，陈伟，“非局部先验贝叶斯层析重建，”数学成像与愿景，第30卷，第2期2，页133-146,2008。视图:出版商的网站|谷歌学术搜索
17. 布ades, B. Coll和j - m。莫尔，“图像去噪的非局部算法”IEEE计算机愿景和模式识别关于IEEE计算机协会会议的诉讼程序（CVPR '05），第60-65页，2005年6月。视图:出版商的网站|谷歌学术搜索
18. 陈勇，罗莉，陈伟等，“基于混合先验模型的联合map层析重建，”多尺度建模与仿真，第9卷，第5期。4，第1399-1419页，2011。视图:出版商的网站|谷歌学术搜索
19. D. Gembris, M. Neeb, M. Gipp, A. Kugel，和R. Männer，“使用NVIDIA CUDA的GPU系统的相关性分析”，实时图像处理杂志，第6卷，第2期4, pp. 275 - 280,2011。视图:出版商的网站|谷歌学术搜索
20. B. Goossens, H. Luong, J. Aelterman, A. Pižurica，和W. Philips，“用于彩色视频序列去噪的gpu加速实时NLMeans算法”，在智能视觉系统的先进概念，第6475卷计算机科学课堂讲稿，第46-57页，施普林格，柏林，德国，2010。视图:谷歌学术搜索
21. Z. Zhuang，Y. Chen，H. Shu，L. Luo，C.Toumoulin和J.-L。Coatrieux，“使用改进的并行化非局部的快速低剂量CT图像处理意味着过滤”国际医学生物特征学会议论文集， pp. 147-150，中国深圳，2014年6月。视图:出版商的网站|谷歌学术搜索
22. B. Caldairou, N. Passat, P. A. Habas, C. Studholme, F. Rousseau，“非局部模糊分割方法在脑MRI中的应用”，模式识别，第44卷，第5期。第9页，1916-1927,2011。视图:出版商的网站|谷歌学术搜索
23. M. Protter，M. Elad，H. Takeda和P. Milanfar，“概括为超级分辨率的重建，”IEEE图像处理汇刊第18卷第2期1，pp。36-51,2009。视图:出版商的网站|谷歌学术搜索|MathSciNet
24. “基于非局部算子的欠采样核磁共振图像重建”，“基于非局部算子的欠采样核磁共振图像重建”，医学图像分析第18卷第2期6, pp. 843-856, 2014。视图:出版商的网站|谷歌学术搜索
25. J.Ma，J. Huang，Q. Feng等人，“使用以前的正常剂量扫描的低剂量计算断层扫描图像恢复”医学物理学第38卷第2期10, pp. 5713-5731, 2011。视图:出版商的网站|谷歌学术搜索
26. 陈颖，石磊，冯秋琴等，“低剂量CT图像处理的伪迹抑制字典学习”，IEEE医学影像汇刊第33卷第3期12, pp. 2271-2292, 2014。视图:出版商的网站|谷歌学术搜索
27. 基于截断左手Grünwald-Letnikov公式的分数阶导数的新定义α< 1和中位数校正，”摘要与应用分析， 2014年，第914386号，9页，2014。视图:出版商的网站|谷歌学术搜索
28. 杨军，王勇，陈伟，“基于冠状动脉图像的三维重建方法”，IEEE图像处理汇刊第18卷第2期7，pp。1563-1572，2009。视图:出版商的网站|谷歌学术搜索|MathSciNet
29. Q.冯，M. Foskey，W. Chen和D. Shen，使用人口和患者特异性统计分割CT前列腺图像进行放射治疗，“医学物理学，第37卷，第2期8，页4121 - 4132,2010。视图:出版商的网站|谷歌学术搜索
30. 杨军，王勇，唐胜，周胜，刘勇，陈伟，“基于薄板样条的x射线血管造影图像的多分辨率弹性配准，”IEEE核科学学报第54卷第5期1，页152-166,2007。视图:出版商的网站|谷歌学术搜索

版权

版权所有©2015 Libo Zhang等人。这是分布下的开放式访问文章创意公共归因许可证，允许在任何媒介上不受限制地使用、传播和复制，但必须正确引用原作。