国际生物医学成像杂志》上

PDF
国际生物医学成像杂志》上/2011年/文章
特殊的问题

并行计算在医疗成像应用

把这个特殊的问题

研究文章|开放获取

体积 2011年 |文章的ID 952819年 | https://doi.org/10.1155/2011/952819

安德斯·埃克伦,Mats Andersson汉斯Knutsson, 真正的四维图像去噪在GPU上”,国际生物医学成像杂志》上, 卷。2011年, 文章的ID952819年, 16 页面, 2011年 https://doi.org/10.1155/2011/952819

真正的四维图像去噪在GPU上

学术编辑器:哈立德z Abd-Elmoniem
收到了 2011年3月31日
修改后的 2011年6月23日
接受 2011年6月24日
发表 2011年10月01

文摘

图像去噪技术的使用是许多医疗成像应用的一个重要组成部分。一个常见的应用是提高图像质量低剂量(噪声)计算机断层扫描(CT)的数据。而三维图像去噪之前已经应用于独立几卷,没有很多工作真正的4 d图像去噪,在算法同时考虑几卷。问题4 d图像去噪与2 d和3 d去噪相比,是计算复杂度成倍增加。在本文中,我们描述一个小说真正的4 d图像去噪的算法,基于局部自适应滤波,以及如何实现它的图形处理单元(GPU)。该算法应用于心脏4 d CT数据集的分辨率512×512×445×20。结果是GPU可以完成在大约25分钟如果使用空间滤波去噪和大约8分钟如果FFT-based过滤。CPU实现需要几天的处理时间为空间滤波和大约50分钟FFT-based过滤。短的处理时间增加真实的临床价值显著4 d图像去噪。

1。介绍

图像去噪是常用的医学成像,以帮助医生看到异常图像。图像去噪是第一个应用于二维图像(1- - - - - -3),然后扩展到三维数据(4- - - - - -6),3 d数据可以收集一些2 d图像随着时间的推移,或作为一个3 d体积。大量的医学成像模式(例如,计算机断层扫描(CT)、超声(美国)和磁共振成像(MRI))现在收集4 d数据提供的可能性,也就是说,时间分辨的体积数据。这使得有可能,例如,检查哪些部分的大脑是活跃在特定任务(功能性磁共振成像(fMRI))。虽然4 d CT数据可以看到3 d的心跳,缺点是较低的x射线曝光用于4 d CT数据收集,而3 d CT数据采集,为了不伤害病人。当接触的数量减少,噪声的数据量大幅增加。

三维图像去噪曾被应用于独立几个时间点,但没有工作真正的4 d图像去噪的算法同时考虑几个卷(而不是单个卷一次)。Montagnat et al。7)四维各向异性扩散滤波应用于超声波卷和Jahanian et al。8)应用4 d小波去噪扩散张量磁共振成像数据。对于CT数据,它可以额外的有益的使用时间维度去噪,作为重建的一些文物随时间。从而可以充分利用移除这些文物的4 d数据。虽然真正的4 d图像去噪是非常强大的,缺点是处理时间对维数呈指数增长。

图形处理单元(GPU)的快速发展导致了很多算法在医学成像领域一直在GPU上实现,为了节省时间,能够运用更高级的分析。给GPU快速发展的一个例子,比较三个消费者从表中给出了Nvidia图形卡1。每一代GPU之间的时间是2 - 3年。在医学成像领域的一些例子,利用GPU的计算能力是图像配准9- - - - - -13,图像分割14- - - - - -16和功能磁共振成像分析17- - - - - -20.]。


属性/ GPU 9800 GT GTX 285 GTX 580

数量的处理器核心 112年 240年 512年
正常的全局内存的大小 512 MB 1024 MB 1536 MB
全局内存带宽 57.6 GB / s 159.0 GB / s 192.4 GB / s
不变的记忆 64 KB 64 KB 64 KB
共享内存/议员 16 KB 16 KB 48/16 KB
每个议员浮点寄存器 8192年 16384年 32768年
每个议员L1缓存 没有一个 没有一个 16/48 KB
L2高速缓存 没有一个 没有一个 768 KB

在图像去噪领域,一些算法也已在GPU上实现。已经在2001年Rumpf和Strzodka [21)描述了如何应用各向异性扩散3在GPU)。杰夫里(22)做了一个比较不同的GPU实现三维各向异性扩散和双边滤波的数据。苏和徐23)在2010年提出了如何利用GPU加速小波图像去噪。Zhang et al。24]描述了基于gpu的图像处理和增强技术动态体积医学图像可视化,但增强的在这种情况下是指增强可视化,而不是4 d数据。最近,GPU用于实时图像去噪。2007年,陈等人。25)使用双边滤波(26在GPU上]edge-aware图像的实时处理。丰特斯et al。27)2011年GPU用于实时超声数据的去噪和古森斯et al。28]在2010年成功地运行常用的外地意味着算法(29日在真正的时间。

据我们所知,没有任何工作做些真正的4 d图像去噪在GPU上。在这个工作我们因此提出一个新颖的算法,基于局部自适应过滤、四维去噪和描述如何在GPU上实现它,为了减少处理时间,从而大大增加的临床价值。

2。方法

在本节中,该算法用于真正的4 d图像去噪将被描述。

2.1。维度的力量

展示更多的维度,维度的力量,可以改善去噪的结果,一个小测试是首先进行合成数据。4 d数据的大小 ,但是没有信号变化在过去的两个维度。数据包含一个大的一步,细线,阴影从左上角到右下角。大量的4 d加性噪声终于添加到数据中。图像去噪不同维度的应用。对于二维的情况,进行了去噪 形象,3 d的情况下,去噪是一个完成的 体积和4 d的情况下使用的所有数据。一个各向异性使用低通滤波器去噪,过滤器和相同的维度数据和面向结构。原始测试数据,测试数据与噪声去噪结果图1。很明显,明显改善去噪结果为每个新维度。

2.2。自适应滤波在4 d

去噪的方法,我们的工作是基于自适应滤波。它介绍了2 d Knutsson等人于1983年(21992年],然后扩展到3 d (4]。在这部作品中,相同的基本原则是用于自适应滤波在4 d。主要的思想是先估计局部结构张量(30.](通过使用第一套过滤器)在每一个地区的数据,然后让张量控制重建过滤器(第二组过滤器)。重建这一术语在本文中不应混淆CT数据的重建。当地的结构张量 在4 d 对称矩阵在每次体素, 和包含数据的局部结构信息,可以用来控制重建滤波器的权重。的结果自适应滤波,平滑完成结构(如线条和边缘在2 d),但不垂直。

2.3。自适应滤波比其他图像去噪的方法

最近开发的图像去噪方法相比(例如,外地意味着[29日),各向异性扩散3)和双边滤波(26]),自适应滤波在我们的示例中用于图像去噪的4 d原因有三。首先,自适应滤波计算比其他方法更有效。外地意味着可以得到很好的结果,但是该算法可以非常耗时(即使使用gpu)。各向异性扩散是一个迭代算法,因此可以相当缓慢。自适应滤波是一种直接方法不需要迭代。双边滤波不仅需要每个滤波器系数的乘法和每个数据值,还有一个强度范围的评价函数(例如,一个指数)这是昂贵得多比乘法来执行。第二,我们去噪算法的调优参数的,而容易理解和探索。当第一次去噪结果已经获得,通常明显如何改变参数来改善结果。这对其他方法并非总是如此。第三,自适应滤波方法已经被证明是非常健壮的(非常很少,一个奇怪的结果是获得)。自适应滤波是用于二维图像去噪在商业软件20多年临床和3 d最近的一项研究[31日]证明了其潜力、健壮性和临床验收。外地意味着算法只适用如果数据包含几个街区具有相似属性。

2.4。使用正交滤波器估计局部结构张量

当地的结构张量可以,例如,使用正交滤波器估计(5,30.]。正交滤波器 一分之零一半的频域(过滤器)的定义的方向,可以表示为两个极性可分离函数,径向函数 和一个定向功能 , 在哪里 是频率变量。径向函数是一个对数正态函数 在哪里 滤波器的中心频率和吗 在八度是带宽()。定向功能取决于角度 方向向量之间的过滤器 和归一化频率坐标向量 作为 , 正交过滤器是笛卡尔不可分的和复杂的价值在空间域,甚至真正的一部分是在2 d作为探测器,而虚部是奇数和2 d作为边缘检测器。在3 d,偶数和奇数过滤器对应于一个平面探测器和一个3 d边缘检测器。4 d的平面和三维边缘可能除了时变。的复数滤波器响应 是一个带通滤波估计的版本的分析信号大小 和相位 , 张量乘以计算正交滤波器的幅度响应 过滤方向向量的外积 然后求和结果对所有过滤器 , 在哪里 标量常数,取决于数据的维数(5,30.), 正交滤波器的数量吗 是单位矩阵。由此产生的张量是不变的阶段正交滤波器的幅度响应,不变的是当地社区的类型(例如,在2 d明亮的线条,黑暗,黑暗,光明的边缘,等等)。这与在当地的结构张量估计利用梯度算子,如Sobel过滤器。

过滤器的数量需要估计张量取决于数据的维数和由独立的组件的数量给出对称局部结构张量。所需数量的过滤器是3对2 d, 6 10 3 d和4 d。给定的张量公式,然而,假设过滤器是均匀分布。可以6 3 d过滤器均匀扩散,但它是不可能传播10过滤器均匀4 d。因为这个原因,12正交滤波器用于4 d(即。,a total of 24 filters in the spatial domain, 12 real valued and 12 complex valued). To apply 24 nonseparable filters to a 4D dataset requires a huge number of multiplications. In this paper a new type of filters,单项过滤器(32),因此代替。

2.5。估计使用单项过滤器局部结构张量

单项过滤器也有一个径向函数 和一个定向功能 。的定向部分单项过滤器产品的正整数频率变量的组件 。o(1)的单项滤波器矩阵, ,两个 在频域定义为 的单项过滤器是第一次描述了2 d,然后推广到四维。

2.5.1。在2 d单项过滤器

在2 d,频率变量定义为在此工作 。方向一阶单项过滤器的一部分 在空间域 在频域。二维图给出了一阶的单项过滤器2。方向二阶单项过滤器的一部分 在空间域 在频域。二维二阶的单项过滤器在图3

单项滤波器响应矩阵 要么是卷积计算的空间域或频域乘法。一个简单信号阶段 (例如, );订单1和2的单项滤波器响应矩阵可以写成 一阶产品有奇怪的功能,从而与奇怪的正弦函数,二阶相关产品是偶函数,从而即使余弦函数(注意正交滤波器的相似之处,甚至有一个实部和一个奇怪的虚部)。通过使用这一事实 外产品的滤波器响应矩阵 当地的结构张量 然后计算 从这个表达式,很明显,估计张量,正如前面,相位不变的广场是一个奇怪的部分,甚至是一个部分的平方的总和。如何计算高阶的张量信息单项,看到我们最近的工作(32]。

2.5.2。单项4 d的过滤器

总共14不可分的4 d单项过滤器(4一阶的奇怪 甚至十二级 )空间的支持 体素应用于CT卷。过滤器有一个径向函数对数正态和中心频率 和2.5个八度的带宽。过滤器内核优化对理想的频率响应,空间位置,和预期的信噪比(5,33]。

利用方程(10)4 d的情况下,更换频率与单项滤波器响应变量,10个组件的结构张量计算 在哪里 表示滤波器响应单项过滤器 。第一项有关 ,与其余的条款 ,在总

如果单项过滤器是用来代替正交过滤器,所需数量的4 d过滤器从而减少从24 - 14所示。另一个优势是,单项过滤器需要一个较小的空间支持,这使得它更容易保留细节和对比的处理。一个较小的空间支持也导致更低的数字滤波器系数,减少处理时间。

2.6。控制张量

当局部结构张量 据估计,然后映射到一个控制张量 级,通过映射(能量)的各向同性张量。这个映射的目的是进一步提高去噪。2 d和3 d图像去噪,该映射可以通过先计算结构张量的特征值和特征向量的每个元素数据。描述的映射是第一个2 d,然后4 d。

2.6.1。映射在二维张量的大小

在二维情况下,大小 张量的计算 在哪里 是两个特征值。大小 规范化是在0和1之间,然后映射到不同吗 称所谓的辅助功能 在哪里 , , 是用于控制参数映射。的 变量成正比数据的信噪比(信噪比)和作为软噪声阈值, 主要控制过度(可用于动态范围压缩或放大区域级略高于噪声阈值),和 主要控制曲线的斜率/柔软。这个映射的目的是控制高通滤波的一般用法信息。高通的信息应该只用于哪里有定义良好的数据结构。如果结构张量的大小很低,可以假定一个地区只包含噪音。给出了函数m图的一些示例4

2.6.2。在二维映射的各向同性张量

各向同性的 在2 d计算吗 和映射到 所谓的mu-function根据 在哪里 是用于控制参数映射, 主要是控制曲线的过渡和 主要控制斜率/柔软。这个映射的目的是控制高反差保留信息的使用非惯用的方向,也就是说,方向是由最小的特征值对应的特征向量。这是通过提高张量的各向同性如果是轻微的各向同性,甚至使其更多的各向异性,如果是各向异性。给出了mu-function图的一些示例5。给出了各向同性映射图的一些示例6。辅助功能和mu-function进一步解释(5]。

2.6.3。二维张量的映射

控制张量 最终计算 在哪里 是最大的特征值对应的特征向量 是最小的特征值对应的特征向量 。因此保存的映射eigensystem,但改变了特征值,从而的形状张量。

2.6.4。完整的四维张量的映射

矩阵的大小 ,有直接的公式如何计算特征值和特征向量,但是 矩阵,没有公式和这样的复杂映射。当然是有可能通过其他方法来计算特征值和特征向量,如迭代算法,但这将是非常耗时的映射控制张量在每次体素。映射的局部结构张量在这个工作因此控制张量的方式不明确需要计算特征值和特征向量。张量的大小是第一次计算 在哪里 弗罗贝尼乌斯范数表示。指数将决定如何接近 估计张量大小将;指数会给更好的精度高,但一个指数8在实践中被证明是足够的。为了减少计算负荷, 是计算 在哪里 表示矩阵乘法。 然后计算 最大的运营商是整个数据集,这样的最大价值 1, 然后映射到 通过使用辅助功能。

首先绘制出各向同性,结构张量归一化 这样张只携带信息的各向异性(形状)。这一事实 使用eigensystem相同,这样控制张量可以计算吗 在哪里 是单位矩阵。下面的公式是一个临时修改这个基本思想,没有明确需要计算的各向同性 为我们的CT数据,给好的结果。基本的想法是,张量的特征值的比值变化本身乘以一个张量时的次数,从而的形状张量也变化。这种方法并没有给出相同的结果作为初始各向同性的映射,但是它避免了显式计算特征值和特征向量。一个帮助变量 首先是计算的 然后张量的控制 是计算

由此产生的传递函数映射每个特征值在图给出7。小成为更小的特征值,特征值大的变得更大。这个特征值映射的结果是类似于给出各向同性映射示例图6

2.7。计算去噪数据

11不可分的重建过滤器,一个低通滤波器 的零级和10个高通滤波器 二阶,空间的支持 体素应用于CT卷。去噪4 d数据 计算lowpass-filtered数据之和, 为每一个高通滤波器和高通滤波过滤数据 , 加权的组件 张量的控制 , 结果是4 d数据低通滤波器过滤在四面八方,然后高通的信息放回哪里有一个明确的结构。高通的信息放回在当地社区的主要方向(给出的最大特征向量与特征值),如果张量大小是高。高通的信息放回在非惯用的方向(由特征向量与特征值越小),如果张量大小是高各向异性低。

2.8。完整的算法

所有加工步骤的去噪算法给出了表2。在本例中CT数据不包含任何重要的结构信息的频率 因此低通空间维度,卷过滤然后downsampled 2的一个因素 。当局部结构张量据估计,它是低通滤波器过滤、分离低通滤波器的尺寸 ,改善估计每次体素和确保结果重建滤波器不同顺利。注意,这种平滑不降低分辨率图像数据,但只有的分辨率张量场。张量的映射后,控制张量的原始分辨率CT数据的插值。


处理步骤 分辨率,科幻小说 内存消耗,科幻小说 分辨率,FFTBF 内存消耗,FFTBF

低通滤波,将采样的CT卷 406 MB 294 MB

与14个单项过滤器和过滤计算局部结构张量(10),(11)) 1376 MB 1791 MB

低通滤波的局部结构张量组件(规范化卷积(25)) 1276 MB 720 MB

计算张量大小和映射与函数m ((17),(18),(19),(13)) 1376 MB 770 MB

映射的局部结构张量控制张量(20.),(22),(23)) 1376 MB 770 MB

低通滤波控制张量的组件(规范化卷积(25)) 1476 MB 820 MB

与11重建过滤器,过滤插值控制
张量,并计算
去噪数据(24)
2771 MB 2110 MB
(三轮×6片= 18总共去噪片)

虽然提出了算法易于实现,空间滤波有11个重建过滤器的大小 (14 641滤波器系数)应用到数据集的决议 需要大约375 0000亿乘法。这就是为什么GPU是必要的为了做四维去噪在合理的时间内。

2.9。规范化卷积

算法的一个主要的缺点是,使用标准的卷积,有效元素的数量 方向(即。,片)decreases rapidly. If the algorithm is applied to a dataset of the resolution 掉,两片第一由于低通滤波器卷积的大小 。将采样后,有16片数据。单项过滤器的大小 滤波器响应的,因此只有10片是有效的。在每个结构张量的低通滤波组件,另一个四片丢失,然后另外四个低通滤波过程中丢失的张量的控制。结果是这样只剩下两片有效毕竟曲线玲珑。时间维度可能存在同样的问题,但由于心脏循环周期是很自然的在时间方向,使用循环卷积,从而所有的时间点是有效的。

损失的有效切片可以避免使用规范化卷积(34),低通滤波的数据之前,将采样和张量的低通滤波组件。在规范化卷积,肯定是附加到每个信号值。certainty-weighted滤波器响应计算轨 在哪里 是肯定的, 是信号, 过滤器, 表示点态乘法 表示卷积。数据内的确定性设置为1和0以外的数据。注意,这个简单的版本的规范化卷积(归一化平均)不能申请重建的单项过滤器和过滤器,这些过滤器有正面和负面的系数。可以应用这些过滤器的完整规范化卷积的方法,但它会显著增加计算负载。

3所示。GPU实现

在本节中,GPU实现去噪算法的描述。CUDA(统一计算设备架构)编程语言,英伟达(35),解释了柯克和胡36),已被用于实现。开放计算语言(OpenCL) [37)可能是一个更好的选择,因为它可以在任何硬件运行相同的代码。

3.1。创建4 d指数

CUDA编程语言可以很容易地为每个线程生成二维指标,例如,通过使用算法1。生成3 d指数比较难,因为每个线程阻塞三个维,但网格只能两个维。一种方法来生成3 d指数给出了算法2。生成4 d指数是更加困难。在4 d导航数据,使用上述3 d索引方法,和内核就曾呼吁每个时间点。

/ /代码执行之前内核推出了
intthreadsInX=32;
intthreadsInY=16;
intblocksInX=DATA_W / threadsInX;
intblocksInX=DATA_H / threadsInY;
dimGrid=dim3 (blocksInXblocksInY);
dimBlock=dim3 (threadsInXthreadsInY,1);
/ /代码执行内部内核
intx=blockIdx.x*blockDim.x+threadIdx.x;
inty=blockIdx.y*blockDim.y+threadIx.yd;

/ /代码执行之前内核推出了
intthreadsInX=32;
intthreadsInY=16;
intthreadsIn = 1;
intblocksInX = (DATA_W + threadsInX-1) / threadsInX;
intblocksInY = (DATA_H + threadsInY-1) / threadsInY;
intblocksInZ = (DATA_D + threadsInZ-1) / threadsInZ;
dim3dimGrid=dim3 (blocksInX blocksInY * blocksInZ);
dim3dimBlock=dim3 (threadsInX threadsInY threadsInZ);
/ /代码执行内部内核
intblockIdxz=__float2uint_rd (blockIdx.y*invBlocksInY);
intblockIdxy=blockIdx.yblockIdxz*blocksInY;
intx=blockIdx.x*blockDim.x+threadIdx.x;
inty=blockIdxy*blockDim.y+threadIdx.y;
intz=blockIdxz*blockDim.z+threadIdx.z;

3.2。空间与基于FFT的过滤

基于快速傅里叶变换(FFT)时可以非常有效的过滤大不可分的高维度的过滤器适用于大数据集,但空间过滤通常更快如果过滤器是小或笛卡尔分离。FFT-based过滤的主要优势在于处理时间是一样的,不管空间大小的过滤器。一个小免费奖金是循环过滤实现。然而FFT-based过滤的主要缺点是内存需求,随着过滤器需要存储在同一分辨率的数据,并为每个元素复值数量。

看哪一种最适合GPU的过滤,因此实现空间和FFT-based过滤。过滤与小的可分离的低通滤波器(应用数据之前downsampled和光滑的张量组件),只可分空间过滤实现。

3.3。空间滤波

空间滤波可以实现,而许多方面,尤其是在四维空间。一个简单的方法来实现2 d和3 d的GPU是利用过滤纹理内存的缓存,把过滤器内核在内存。然而,这种方法的缺点是实现将内存带宽非常有限,而不是计算性能。另一个问题是,不可能在CUDA编程语言中使用4 d纹理。人会来存储4 d数据作为一个大的1 d纹理或几个2 d或3 d纹理。更好的方法是利用共享内存,增加一个因素3大小之间的Nvidia GTX 285和Nvidia GTX 580。第一次读到的数据共享内存,然后并行滤波器响应计算。通过使用共享内存,线程可以共享数据在一个非常有效的方式,这是有益的过滤结果计算了两个相邻元素主要使用相同的数据。

多维过滤器可以分离或不可分的(单项过滤器和重建过滤器是不可分的,而不同的低通滤波器可分)两种不同的空间滤波函数实现的。

3.3.1。分离过滤

我们可分4 d卷积器实现,首先做过滤的所有行,然后所有的列,然后所有的棒,最后所有的时间点。数据加载到共享内存和有效的并行滤波器响应计算。过滤器内核内存存储在常数。四个内核,使用16 KB的共享内存,这样3线程块可以并行运行在每个多处理器的Nvidia GTX 580。

3.3.2。不可分的过滤

共享内存的方法很好地适合不可分的2 d过滤但不亦不可分的3 d和4 d过滤。共享内存的大小在Nvidia GTX 580 48 KB为每个多处理器,从而只可能,例如,健康 浮点值。如果大小的4 d过滤器 ,只有 可以为每个多处理器生成有效的滤波器响应。更好的方法在4 d是不可分的过滤而不是使用一个优化的二维滤波内核,然后积累的滤波器响应对其他维度求和通过调用每个片二维滤波函数,每个时间点的过滤器。中给出的方法是用伪代码描述算法3

/ /过滤所有时间数据
(intt = 0;t < DATA_T;t + +)
{
/ /过滤所有数据
(intz = 0;z < DATA_D; z + +)
{
/ /过滤器反应GPU0
重置< < <dimGrid, dimBlock> > > (d_Filter_Responses);
/ /过滤所有时间过滤器
(inttt = 0;tt < FILTER_T;tt + +)
{
/ /过滤所有过滤器
(intzz = 0;zz < FILTER_D;zz + +)
{
/ /复制当前的过滤器系数
/ /常数内存GPU
zz CopyFilterCoefficients (tt);
/ /二维过滤GPU
/ /增量过滤器反应
/ /内部过滤函数
Conv2D< < <dimGrid,dimBlock> > >(d_Filter_Responses);
}
}
}
}

我们不可分的二维卷积器首先读取 像素到共享内存,然后计算有效的滤波器响应所有的14个单项过滤器或11重建过滤器同时,最后将结果写入全局内存。两个版本的卷积器实施,最大限度地支持 过滤器和最大限度地支持 过滤器。第一个计算 有效的滤波器响应,和第二个计算 有效的滤波器响应。作为 浮点值只需要16 KB的内存,三个线程块可以每个多处理器上同时运行。这将导致 每个多处理器有效过滤反应。为获得最佳性能,2 d过滤循环完全展开通过Matlab生成的代码脚本。

14单项过滤器的大小 ,这需要135 KB的内存存储为浮动,但是不断的记忆只有64 KB。由于这个原因, 滤波器系数存储一次,然后为每个时间点和每一块更新。它可能商店 滤波器系数,而是只存储 系数,过滤器的大小(2.75 KB)足够小,总是在不断记忆的缓存(8 KB)。同样的方法用于11重建过滤器的大小

3.4。FFT-Based过滤

虽然由Nvidia CUFFT图书馆支持1 d, 2 d和3 d fft算法,没有直接支持4 d fft算法。FFT笛卡儿可分的,但是可以做一个四维FFT应用连续四1 d FFT算法。CUFFT库支持推出一批1 d FFT算法,这样,许多一维FFT的并行运行。批1 d应用fft算法在第一维度的数据存储(例如, 如果数据存储 )。之间1 d FFT,从而需要更改的数据(例如,来自 )。这种方法的缺点是所花费的时间变化的数据可以比实际执行1 d FFT。CUFFT库的最新版本支持推出一批二维FFT。通过两个连续的二维FFT,它足以改变数据的顺序一次,而是三次。

转发4 d FFT首次应用于卷。过滤器与零填充相同的分辨率数据,然后转换到频域。过滤,数据和过滤器之间的复数乘法,然后应用逆4 d FFT应用于滤波器的响应。逆变换后,FFT的转变是必要的;但是没有这样的功能在CUFFT图书馆。当张量分量和去噪数据计算,每个四坐标转移通过使用帮助功能,看到算法4

__device__intShift_FFT_Coordinate (int坐标,intDATA_SIZE)
{
如果(协调>(ceilf (DATA_SIZE / 2)1))
{
返回协调ceilf (DATA_SIZE / 2);
}
其他的
{
返回协调+floorf (DATA_SIZE / 2);
}
}

作为单项过滤器只有一个实部或虚部在空间域,节省一些额外的时间,把一个单项过滤器的实部和虚部的另一个单项过滤前的4 d FFT应用于在过滤器。当执行复杂的乘法在频域中,因此两个过滤器应用在同一时间。倒数4 d FFT后,第一个过滤器响应提取提取实部和第二滤波器响应的虚部。同样的技巧用于10高通滤波重建过滤器。

3.5。内存方面的考虑

实现四维去噪算法的主要问题在GPU上的有限大小全局内存(3 GB在我们的例子中)。这是更加困难的GPU驱动程序可以使用多达100 - 200 MB的全局内存。在GPU上存储所有CT数据在同一时间是不可能的,一个CT体积的决议 大约需要467 MB的内存如果使用32位浮点数。存储滤波器响应的问题更多。举一个例子,来存储所有的11重建滤波器响应作为数据集的浮动大小 需要大约103 GB的内存。因此去噪的完成片(例如,16或32)。

空间滤波的算法从数据的分辨率 并downsampled 。控制张量计算 体素,和去噪数据计算 压的时间。处理所有的445片需要12分。

FFT-based过滤的算法从数据的分辨率 并downsampled 。然后计算出控制张量 体素,和去噪数据计算 压的时间。处理所有的445片需要26分。

存储的10组件控制张量在同一分辨率作为一个运行的原始CT数据空间滤波( )将需要约12.2 GB的内存。张量的控制需要插值因子2在每一个空间维度,因为据估计downsampled数据,另一种方法是使用。插值GPU的张量是一个完美的任务,由于线性插值的硬件支持。10张量分量,计算,因此存储在10纹理,然后进行插值计算去噪数据时。通过使用这种方法,只有10变量的决议 需要存储在同一时间。

2国家在解决数据,使用的方程,和内存消耗每一步的去噪算法,用于空间滤波和FFT-based过滤。出分辨率指的是有效解决数据在每个处理步骤之后,一些数据过滤后被认为是无效的操作。的内存消耗大的原因FFT-based滤波是空间滤波可以做一片或一次一个卷,虽然FFT-based过滤应用到数量足够大的切片和时间点在同一时间。我们不能使用超过2 GB的内存FFT-based过滤;原因之一可能是CUFFT函数内部使用临时变量使用的一些记忆。自从CUFFT库的源代码不可用,很难进一步调查这一假设。

4所示。数据

的4 d CT数据集用于测试我们的GPU实现收集与西门子SOMATOM定义Flash双能CT扫描仪中心的医学影像科学可视化(CMIV)。数据集包含近9000 DICOM文件和数据的分辨率为512×512×445×20像素点的时间。各体素的空间大小 毫米。在图像采集管电流调制在心动周期和减少辐射暴露在心脏收缩阶段。因此,噪声随时间的数量。

5。结果

5.1。处理时间

比较我们的GPU实现处理时间和CPU实现。使用GPU是Nvidia GTX 580,配备512处理器核心和3 GB的内存(Nvidia GTX 580通常配备1.5 GB的内存)。使用的CPU是一个英特尔至强2.4 GHz 4处理器核心和12 MB的L3缓存,使用了12 GB的内存。所有的实现使用32位浮点数。所使用的操作系统是Linux Fedora 14 64位。

对于CPU实现,OpenMP(开放多处理)库(38,39),这样所有的4个处理器核心并行工作。没有其他类型的CPU的优化,如SSE2,。我们十分清楚,可以使一个更好的CPU实现。这个比较的目的是给一个指示CPU和GPU的性能。如果CPU的代码将矢量化,CPU处理时间可以除以一个因素3或4 (FFT已经非常优化的除外)。

表中给出了处理时间3,4,5,6。空间滤波的CPU处理时间估计,因为它需要几天在整个数据集上运行算法。multi-GPU实现的处理时间尺度,而线性的GPU的,因为每个GPU并行可以工作在不同子集的片。我们的电脑包含三个GPU,所有GPU的处理时间从而可以除以3的一个因素。


数据大小 空间滤波CPU 空间滤波GPU GPU加速 FFT滤波CPU FFT滤波GPU GPU加速

17.3分钟 5.7秒 182年 25岁的年代 1.8秒 13.9
2.3 h 36.0秒 230年 3.3分钟 14.3秒 13.9


数据大小 CPU GPU GPU加速

42岁的年代 1.0秒 42
292年代 7.3秒 40


数据大小 空间滤波CPU 空间滤波GPU GPU加速 FFT滤波CPU FFT滤波GPU GPU加速

2 7.5 h 3.3米 136年 5.6分钟 1.1分钟 5.1
2.5天 23.9米 150年 45分钟 8.6分钟 5.2


数据大小 空间滤波CPU 空间滤波GPU GPU加速 FFT滤波CPU FFT滤波GPU GPU加速

7.8 h 3.5米 133年 6.7米 1.2米 5.6
2.6天 26.3米 144年 52.8米 8.9米 5.9

5.2。去噪结果

显示四维去噪的结果,原始的CT数据去噪相比,数据采用体绘制。免费MeVisLab软件开发项目(http://www.mevislab.de/)使用。两个体积渲染器,一个用于原始数据,一个用于数据去噪,同时运行,并同步视角和传递函数。图8显示了原始的体积渲染和去噪数据对不同时间点和视图天使。很明显,很多噪声去噪被移除,但由于数据的去噪算法改变了直方图,很难做一个客观的比较,即使相同的传递函数。

电影原和去噪数据探索和两个体积渲染器也。视频,数据downsampled因子2在空间维度,以减少内存使用量。体积渲染器自动循环遍历所有的时间点。视频可以发现http://www.youtube.com/watch?v=wflbt2sV34M

通过看视频,很容易看到原始数据中的噪声随时间。

6。讨论

我们提出了如何实现真正的4 d图像去噪在GPU上。结果是4 d图像去噪成为实际上可能如果使用GPU,从而显著临床价值增加。

6.1。处理时间

让一个完全公平的比较在CPU和GPU之间相当困难。它一直在讨论40)如果GPU加速已报告在文献中是合理的或如果他们与未优化CPU实现的结果。在我们看来,理论和实践的处理性能,可以实现不同的硬件并不是唯一一个有趣的话题。在研究环境中,之间可实现的加工性能和所花费的时间的实现也很重要。从这个角度看,我们认为CPU-GPU比较是相当公平的,因为几乎在同一时间花在做CPU和GPU实现。CUDA编程语言设计和开发了并行计算从一开始,虽然不同的插件添加到C编程语言能够做并行计算。虽然它很容易使CPU多线程实现,例如,通过使用OpenMP图书馆,更先进的CPU优化通常是包括和经常需要汇编编程更加困难。

而空间滤波可以显著低于FFT-based不可分的过滤器过滤,有一些优势(较低的内存使用除外)。一是可以选择感兴趣的区域(ROI)去噪,而在整个数据集进行去噪。另一个优势是,过滤网络(41,42)可以应用,这样的滤波器响应许多小过滤器结合相同的滤波器响应从一个大的过滤器。过滤网络可以减少乘法因子5在2 d, 3 d和300年25在4 d (43]。然而网络设计和优化一个过滤器需要更多的工作比优化一个过滤器(33]。另一个问题是,当过滤内存使用量会显著提高网络使用,因为许多滤波器响应需要存储在内存中。过滤网络GPU是一个有前途的未来的研究领域。

从我们的结果,很明显,FFT-based滤波速度比空间滤波对大型不可分的过滤器。数据大小,不是一个两个在每个维度的力量,然而基于FFT的方法可能不那么有效。因为医生通常不会看3 d和4 d数据作为体积渲染,而是2 d切片,然而空间滤波方法的优点是,可以做去噪对感兴趣的区域(例如,一个特定的切片或体积)。这是一个浪费时间来提高数据的部分不是由医生使用。空间滤波方法比FFT-based还可以处理更大的数据集的方法,因为它足以存储滤波器响应为一片或一次一个卷。最近,我们获得的CT数据集与100的时间点,相比20时间点。不可能使用FFT-based方法对于这个数据集。

有几个原因的GPU加速FFT-based过滤远小于GPU加速的空间滤波。首先,CUFFT图书馆不包括任何直接支持4 d FFT的,而且我们必须实现自己的4 d FFT后两个二维FFT的应用。二维FFT的存储顺序之间的数据发生了变化。它可以花更长的时间变化的顺序数据比实际执行FFT。如果英伟达包括直接支持4 d FFT CUFFT图书馆,我们确信他们的实现会比我们的更有效。第二,FFT的CPU非常优化,因为它在很多应用程序中,使用CPU和卷积器是没有完全优化。CUDA编程语言只有几岁,和GPU标准库不优化CPU标准库。gpu的硬件设计也变化很快。一些工作已经完成,以进一步优化CUFFT图书馆。Nukada et al。44,45)已经创建了自己的GPU FFT库已被证明比CUFFT图书馆提供更好的性能。他们绕过的问题改变数据的顺序,从而实现一个更有效的实现。2008年,他们的3 d FFT的5 - 6倍的3 d FFT CUFFT图书馆。第三,由于更大的内存需求FFT-based过滤不可能实现大加速对GPU实现CPU实现。如果一个GPU具有更高的全局内存会被使用,FFT-based实现更高效。

6.2。4 d和CUDA图像处理

正如前面讨论的,4 d图像处理在CUDA更难实现2 d和3 d图像处理。例如,没有4 d纹理,没有4 d fft算法,没有直接支持4 d(或三维)指数。然而,由于功能磁共振成像数据也是4 d,我们曾获得一些经验如何与CUDA (4 d图像处理18- - - - - -20.]。我们得出的结论实施后的4 d图像去噪算法CUDA编程语言是这样CUDA不是非常适合4 d图像处理,但由于其灵活性,它还可能实现该算法很容易。

6.3。真正的5 d图像去噪

似乎不可能的医学图像数据有超过4个维度,但有些工作已经完成关于如何收集5 d数据(46]。三维空间和五个维度两个时间维度,一个用于呼吸的节奏,一个心脏节律。5 d数据的一个主要优势是,病人可以在数据采集正常呼吸,当病人呼吸在收集4 d数据。5 d数据,它是可能的,例如,注视心脏和肺部只看到移动,或注视肺部只看到心跳。如果提出算法扩展到5 d,它将需要使用共有20个单项过滤器和16重建过滤器。5 d数据集的大小 所需数量的乘法,空间滤波重建过滤器将比375年增加0000亿的4 d到11.9亿亿( 5 d)。重建滤波器响应的大小会增加103 GB的4 d为5 d 2986 GB。这仍然是只有一个数据集的一个耐心,我们预计,空间和时间分辨率的医学成像模式在未来将进一步增加。除了5维的,还可以收集数据与不止一个内心的维度。这是,例如,如果心脏的血液流动进行了研究。对于流数据,需要存储在一个三维向量每次体素,而不是一个强度值。

7所示。结论

最后,利用GPU,真正的4 d图像去噪变得实际可行。我们的实现当然可以应用到其他模式,如超声和MRI,不仅CT数据。短的处理时间也几乎可以进一步提高去噪算法和调整所使用的参数。

实际可行的发展之间的运行时间2 d (2)和三维(4)图像去噪技术是10年,从3 d、4 d花费的时间大约是20年。由于gpu的快速发展,它是希望需要等待10 - 20年5 d图像去噪。

确认

这项工作得到了林奈CADICS和研究中心批准号2008 - 3813年由瑞典研究理事会。CT数据收集中心的医学影像科学可视化(CMIV)。作者要感谢NovaMedTech林雪平大学项目的财政支持GPU硬件、功能Johan但对于支持CUDA的安装,和王Chunliang设置传输函数体绘制。

引用

  1. js。李,“数字图像增强和噪声过滤利用当地的统计,“IEEE模式分析与机器智能,卷2,不。2、165 - 168年,1980页。视图:谷歌学术搜索
  2. h·e·Knutsson r·威尔逊和g . h . Granlund“各向异性非平稳图像估计及其applications-part我:恢复噪声图像,”IEEE通信没有,卷。31日。3、388 - 397年,1983页。视图:谷歌学术搜索
  3. p . Perona一起和j·马利克尺度空间,边缘检测使用各向异性扩散。”IEEE模式分析与机器智能,12卷,不。7,629 - 639年,1990页。视图:出版商的网站|谷歌学术搜索
  4. h . Knutsson l . Haglund h .酒吧间招待员,g . Granlund”框架的各向异性自适应过滤和分析图像序列和交易量,”《IEEE国际会议音响、语音信号处理,(ICASSP),第472 - 469页,1992年。视图:谷歌学术搜索
  5. g . Granlund和h . Knutsson对计算机视觉信号处理,波士顿Kluwer学术质量,美国,1995年。
  6. 张炳扬。威斯汀,l . Wigstrom t . Loock l . Sjoqvist r . Kikinis h . Knutsson,“三维自适应滤波在磁共振血管造影术,”磁共振成像杂志》上,14卷,第71 - 63页,2001年。视图:谷歌学术搜索
  7. j . Montagnat m . Sermesant h . Delingette g . Malandain和n . Ayache“各向异性过滤基于模型的分割的4 d圆柱形超声心动图图像,”模式识别的字母,24卷,不。4 - 5,815 - 825年,2003页。视图:出版商的网站|谷歌学术搜索
  8. h . Jahanian a Yazdan-Shahmorad, h . Soltanian-Zadeh”先生的4 d小波噪声抑制扩散张量数据,”《IEEE国际会议音响、语音信号处理,(ICASSP)2008年4月,页509 - 512。视图:出版商的网站|谷歌学术搜索
  9. k . Pauwels和m . m . Van Hulle“实时分阶段光学流GPU,”《IEEE计算机学会学报计算机视觉与模式识别研讨会、会议(CVPR)2008年6月,页1 - 8,。视图:出版商的网站|谷歌学术搜索
  10. p . Muyan-Ozcelik j·d·欧文斯j .夏和s . s . Samant”快速变形登记在GPU上:恶魔的CUDA实现,”《计算机科学及其应用国际会议上,(ICCSA)2008年7月,页223 - 233。视图:出版商的网站|谷歌学术搜索
  11. Bui和j·p·布洛克曼”,使用GPGPU加速图像配准,性能分析”第二届研讨会通用图形处理单元处理,(GPGPU-2)2009年3月,页38-45,。视图:出版商的网站|谷歌学术搜索
  12. a·埃克伦·m·安德森和h . Knutsson“阶段基于卷使用CUDA登记,”IEEE国际会议上声学学报》,演讲,和信号处理,(ICASSP)2010年3月,页658 - 661。视图:出版商的网站|谷歌学术搜索
  13. r·夏姆斯·Sadeghi r·肯尼迪,r·哈特利”的调查在多核医学图像配准和GPU,”IEEE信号处理杂志,27卷,不。2、文章ID 5438962, 50 - 60, 2010页。视图:出版商的网站|谷歌学术搜索
  14. a . e . Lefohn j·e·盖茨,r·t·惠特克”3 d互动,基于gpu的水平集分割,“在计算机科学的课堂讲稿卷,2878年,第572 - 564页,2003年。视图:谷歌学术搜索
  15. 诉Vineet和p . j . Narayanan CUDA削减:快图GPU,削减”《IEEE计算机学会学报计算机视觉与模式识别研讨会、会议(CVPR)2008年6月,页1 - 8,。视图:出版商的网站|谷歌学术搜索
  16. a .阿布拉莫夫t . Kulvicius f . Worgotter, b . Dellen“GPU的实时图像分割,”面对Multicore-Challenge学报》上卷,6310在计算机科学的课堂讲稿施普林格,页131 - 142年,2011年。视图:谷歌学术搜索
  17. d . Gembris m . Neeb m . Gipp a .面食和r的方式,“相关分析使用NVIDIA的CUDA GPU系统上,“实时图像处理》杂志上,2010年,页1 - 6。视图:出版商的网站|谷歌学术搜索
  18. 雷欧,o . Friman、m·安德森和h . Knutsson”GPU加速交互式界面功能磁共振成像数据的探索性功能连通性分析,”诉讼IEEE国际会议上的图像处理,(ICIP),第1624 - 1621页,2011年。视图:谷歌学术搜索
  19. a·埃克伦·m·安德森和h . Knutsson“GPU-possibilities功能磁共振成像分析和挑战,”计算机在生物医学方法和项目。在出版社。视图:出版商的网站|谷歌学术搜索
  20. a·埃克伦·m·安德森和h . Knutsson”快速随机排列测试使客观评价方法单一主题功能磁共振成像分析,“国际生物医学成像杂志》上文章ID 627947卷,2011年,2011年。视图:出版商的网站|谷歌学术搜索
  21. m . Rumpf和r . Strzodka”,在图形硬件非线性扩散,”程序如/ IEEE TCVG研讨会上可视化,第84 - 75页,2001年。视图:谷歌学术搜索
  22. m .杰夫,“比较的GPU实现双边和各向异性扩散过滤器对3 d生物医学数据集“技术。众议员劳伦斯- 3425 e,劳伦斯伯克利国家实验室,伯克利,加州,美国。视图:谷歌学术搜索
  23. 许苏y和z,”并行实现小波图像去噪的可编程PC-grade图形硬件,“信号处理,卷90,不。8,2396 - 2411年,2010页。视图:出版商的网站|谷歌学术搜索
  24. 问:张、r . Eagleson和t·m·彼得斯”基于gpu的图像处理和增强技术动态体积医学图像可视化,”《第四届IEEE国际研讨会上生物医学成像:从纳米到宏观,(位ISBI)2007年4月,页1168 - 1171。视图:出版商的网站|谷歌学术搜索
  25. j·陈,美国巴黎,f .勾勒出“edge-aware图像的实时处理双边网格,ACM交易图形”特殊利益集团的程序在计算机图形和交互技术会议,没有。103年,p . 2007。视图:谷歌学术搜索
  26. c预和r . Manduchi灰色和彩色图像双边滤波,”《IEEE 6计算机视觉国际会议1998年1月,页839 - 846。视图:谷歌学术搜索
  27. f .丰特斯·g·巴罗佐,p .轿跑车和p的命运,“实时超声图像去噪,”实时图像处理》杂志上》第六卷,第15 - 22,2010页。视图:谷歌学术搜索
  28. b·古森斯h .陈德良j . Aelterman A . Pizurica和w·飞利浦“GPU-accelerated实时NLMeans算法去噪彩色视频序列,”学报》第12届国际会议上先进的智能视觉系统,概念(ACIVS)卷,6475在计算机科学的课堂讲稿页,46-57施普林格,2010年。视图:谷歌学术搜索
  29. A . Buades b·科尔和j·m·莫雷尔,“非本地图像去噪算法,”《IEEE计算机学会学报计算机视觉与模式识别、会议(CVPR)2005年6月,页60 - 65,。视图:谷歌学术搜索
  30. h . Knutsson“使用张量,代表当地结构”学报北欧图像分析会议(SCIA),第251 - 244页,1989年。视图:谷歌学术搜索
  31. f . Forsberg诉Berghella d·a·默顿k . Rychlak j . meier和比比戈德堡,“比较改进的三维超声成像,图像处理技术”超声波在医学杂志》上卷,29号4、615 - 619年,2010页。视图:谷歌学术搜索
  32. h . Knutsson张炳扬。威斯汀,m·安德森”代表局部结构使用张量二世”学报北欧图像分析会议(SCIA)卷,6688在计算机科学的课堂讲稿施普林格,页545 - 556年,2011年。视图:谷歌学术搜索
  33. h . Knutsson m·安德森和功能j .但,“先进的过滤设计,”学报北欧图像分析会议(SCIA),第193 - 185页,1999年。视图:谷歌学术搜索
  34. h . Knutsson和c·f·威斯汀,”规范化和卷积差:不完整和不确定数据的插值和滤波的方法,”《IEEE计算机学会学报计算机视觉与模式识别、会议(CVPR)1993年6月,页515 - 523。视图:谷歌学术搜索
  35. 英伟达,CUDA编程指南,4.0版。,2011年。
  36. d·柯克和w .研究会会长,编程大规模并行处理器,Handson方法沃尔瑟姆,摩根考夫曼,质量,美国,2010年。
  37. Khronos集团& OpenCL, 2010,http://www.khronos.org/opencl/
  38. OpenMP并行编程的API规范,2011年,http://www.openmp.org/
  39. b·查普曼g . Jost, r·范德不是使用OpenMP,便携式共享内存并行编程美国剑桥,麻省理工学院出版社,质量,2007年。
  40. 诉w·李,c . Kim j . Chhugani et al .,“揭穿100 x GPU与CPU神话:吞吐量计算CPU和GPU的评价,”学报》第37届国际研讨会上计算机体系结构(ISCA)2010年6月,页451 - 460。视图:出版商的网站|谷歌学术搜索
  41. m·安德森、功能j .但和h . Knutsson“过滤网络,”信号和图像处理的程序,(SIP),第217 - 213页,1999年。视图:谷歌学术搜索
  42. b . Svensson m·安德森和h . Knutsson”过滤网络有效评估当地的三维结构,”诉讼IEEE国际会议上的图像处理,(ICIP)2005年9月,页573 - 576。视图:出版商的网站|谷歌学术搜索
  43. m·安德森、功能j .但和h . Knutsson”顺序滤波器树高效2 d、3 d和4 d方向估计,“技术。众议员石- isy - r - 2070,电气工程系,林雪平,瑞典林雪平大学,1998。视图:谷歌学术搜索
  44. a . Nukada理事长绪方y、t . Endo和s .松岗”带宽密集三维FFT内核使用CUDA gpu,”学报》国际会议为高性能计算、网络、存储和分析,(SC)2008年11月,页1 - 11,。视图:出版商的网站|谷歌学术搜索
  45. a . Nukada s松岗,“自动三维FFT库CUDA gpu,”高性能计算国际会议网络、存储和分析,(SC)2009年11月,页1 - 10,。视图:出版商的网站|谷歌学术搜索
  46. a . Sigfridsson j . p . e . Kvitting h . Knutsson和l . Wigstrom“五维MRI结合同时解决心脏和呼吸阶段容积成像,”磁共振成像杂志》上,25卷,不。1,第121 - 113页,2007。视图:出版商的网站|谷歌学术搜索

版权©2011安德斯·埃克伦等。这是一个开放的分布式下文章知识共享归属许可,它允许无限制的使用、分配和复制在任何媒介,提供最初的工作是正确引用。


更多相关文章

PDF 下载引用 引用
下载其他格式更多的
订单打印副本订单
的观点3281年
下载1169年
引用

相关文章

我们致力于分享发现相关COVID-19尽快。我们将提供无限的出版费用豁免接受研究文章以及案例报告和案例系列COVID-19有关。评论文章被排除在这个豁免政策。注册在这里作为一个评论家,帮助快速新提交。