锥束CT扫描方式对ART算法图像重建的影响研究

锥束CT扫描方式对ART算法图像重建的影响研究 第 34 卷 第 10 期 核技术Vol. 34, No.10 2011 年 10 月NUCLEAR TECHNIQUES October 2011 锥束 CT 扫描方式对 ART 算法图像重建的影响研究 伍建辉 李公平 周 龙 潘小东 林俏露 (兰州大学 核科学与技术学院 兰州 730000) 摘要 平板探测器锥束 CT 的扫描方式对 ART(Algebraic Reconstruction Techniques)算法图像重建会产生影响。 本文用 128 ， 128 ， 128 三维 Shepp-Logan 头部模型作为仿真对象，探讨了不同源-扫描中心距离、投影数、缺 角度下扫描方式对图像重建的影响。结果明，源-扫描中心距离与扫描中心-探测器距离成一固定比例时，其 重建图像较好;投影数固定情况下，均匀扫描 180?重建图像质量较好;缺失角度情况下，投影数尽可能以1? –2? 角度间隔均匀填满第一角度间隔区域，剩余投影数均匀填满第二角度间隔区域时，其重建图像质量较好。 关 键词 锥束 CT，扫描方式，ART 算法，三维 Shepp-Logan 头部模型，图像重建质量 中图分类号 R814.42 计算机断层成像技术(Computed Tomography, CT 图像重建过程，探讨锥束 CT 扫描方式对 ART [1] CT)是公认的最佳无损检测手段，已广泛应用于工 算法图像重建的影响。 业、医学等领域。射线穿过被测物会衰减，CT 用 1 三维 ART 算法图像重建原理 不同角度的射线衰减重构物体内部密度，从而检测 [2]物体内部的缺陷、材料成分等。CT 可分为传统 三维 ART 算法图像重建原理与二维 ART 算法 CT、螺旋 CT、锥束 CT。锥束 CT 的扫描速度快、 一致，但其被测物体素长度计算及投影方式不同。 图像分辨率高、辐射利用率高，作一次圆周扫描即 1.1 锥束 CT 成像原理 可完成传统 CT 及螺旋 CT 的数百个断层扫描，且 其 X、Y、Z 方向的分辨率一致，能明显减少传统图 1 为锥束 CT 扫描示意图，其中 i–o 为穿过 CT 及螺旋 CT 检测中常见的容积效应。 图像重建被测物的射线，D为射线源-旋转轴距离， source-axis 算法是锥束 CT 检测的重要环节，直 D为旋转轴-平板探测器距离，被测物绕中 axis-detection 接关系到图像重建质量。其算法大体可分为分析法 心轴转动。穿越被测物的 X 射线锥束被探测器接收， 和迭代法 。现常用 的分析法 是 FDK 算法 所得数据(该角度的投影)传输给计算机。将被测物 (Feldkamp-Davis-Kress algorithm)，经深入研究已提 转过一个设定角度，进行相同的操作，直至完成预 出一些减小误差的方法。但解析算法要求有均匀、 定角度的扫描，由计算机进行图像重建及图像处理。 全角度的投影，而有些应用对象并无可能作全角度 扫描。ART 算法(Algebra Reconstruction Technique) 能很好解决这一问题，ART 算法可在缺角度、不均 匀投影情况下很好重建图像，这在工业及医学中有 非常重要的价值。但 ART 算法不像 FDK 算法那样 能用硬件实现，这在一定程度上限制了其应用。而 随着计算机技术的发展，重建速度不再是大问题。 一些学者对 ART 算法的硬件化及并行化已取得了 [3–5]一定成果，这将加速 ART 算法的实际应用。对 加速锥束 ART 算法权因子计算、锥束 ART 算法图 图 1 锥束 CT 扫描示意图 [6–9]像重建的投影排序等，已有深入研究。本文基于 Fig.1 Diagram of cone beam CT scanning. Visual C ，， 6.0 及 MATLAB 7.0 实验平台，仿真 —————————————— 第一作者:伍建辉，男，1985 年生，2011 年 7 月于兰州大学核科学与技术学院获硕士学位，主要从事工业 CT 研究 通讯作者:李公平，E-mail:ligp@lzu.edu.cn 收稿日期:2011-02-20，修回日期:2011-07-03 j , n 1.2 ART 算法图像重建 (k ) p , l ， u , i ij jj ,1 ( k ，1) ( k ) 强度为 I的 X 射线穿过衰减系数为 μ 长度为 l 0 (3) ， l U , U ， ,ij j , n 2 的体素后强度变为 I，则其对该体素的投影值为l , ijj ,1 (1) /I) = μlln(I 0其中，,为松弛因子，k 为迭代次数，p为第 i 条射 i 投影值一般写为 p。若有一系列强度已知的射 线的投影值。 线 1, 2, 3, …, m 穿过被测物，则有: 式(2)中，n 个体素涉及 n 个方程，若其中某些 ,l， ,l， ,l， ， ,l, p 111 221 331 n n1 1方程线性相关，则会影响着方程组的解。同样，ART ,l， ,l， ,l， ， ,l, p 112 222 332 n n 2 2算法中方程组建立的方式对图像重建的精度影响也 (2) 不可忽视。因此，CT 的投影方式(即方程建立方 ,l， ,l， ,l， ， ,l, p 式)会影响 ART 算法图像重建。 11m 22 m 33m n nm m 其中，l表示第 m 条射线经过第n 个体素时的长度。 nm 2 不同扫描方式对图像重建的影响 图像重建，就是通过上述方程组求解各体素的衰减 采用 C ，， 6.0 编写程序，仿真锥束 CT 扫描过 系数 μ 值。但因射线方程个数巨大(如 128 ， 128 ， 128 [11]3 程。扫描参数为:扫描轨迹为圆形;X 射线源为 个体素，就要建立 128 个方程)，不易直接求解， 理想点光源，单能光谱;射线源最大锥角 9.6?; 而是采用不同算法(如本文的 ART 算法)求解。 mm;探测器矩阵 128 ， 128;探测器 D=780 ART 算法图像重建的思路为:任取一初始值 source-axis(0)(0)像元尺寸 0.512 mm ; 三 维图像 重 建尺寸 U(U为 n ， 1 维向量)，通过式(3)迭代求一次近似 (1)(1)(2)128 ， 128 ， 128 ;重建体 素尺寸 0.256 mm 。用 图像 U，再基于 U用式(3)求一次近似图像 U， [10]MATLAB 7.0 构造 128 ， 128 ， 128 的 Shepp-Logan 如此继续，直至所求的近似图像达到要求为止。 [12,13]头部模型，由十个椭球组成，其参数见表 1。 表 1 头部模型参数 Parameters of Shepp-Logan model Table 1 , , , , No xyza b c 000 1 0.000 0.000 0.000 0.690 0.920 0.900 0 0 0 1 2 0.000 0.000 0.000 0.662 0.874 0.880 0 0 0 –0.8 3 –0.220 0.000 –0.250 0.410 0.160 0.210 0 0 –108 –0.2 4 0.220 0.000 –0.250 0.310 0.110 0.220 0 0 –72 –0.2 5 0.000 0.350 –0.250 0.210 0.250 0.500 0 0 0 0.1 6 0.000 0.100 –0.250 0.046 0.046 0.046 0 0 0 0.1 7 –0.080 –0.650 –0.250 0.046 0.023 0.020 0 0 0 0.1 8 0.060 –0.650 –0.250 0.046 0.023 0.020 0 0 90 0.1 9 0.060 –0.105 0.625 0.056 0.040 0.100 0 0 90 0.1 10 0.000 0.100 0.625 0.056 0.056 0.100 0 0 0 0.1 注: x, y, z为椭球的中心坐标; a、b、c 分别为X、Y、Z 方向的半轴长; ,、,、,分别为椭球绕X、Y 、Z 轴逆时针旋转 000 的角度; ,为椭球密度, 用来模拟人体断层组织对X 射线的衰减系数 Note: x, yand zis the coordinates to ellipsoid center; a, b and c are the axle length of X, Y and Z direction, respectively; ,, , 00 0 and , are angles of X, Y and Z axis against clockwise rotation around ellipsoid; , is the density of ellipsoid to simulate human tissue on the X-ray tomography attenuation coefficient. 采用以下两个判据来评价重建后的图像: 第 10 期 伍建辉等:锥束 CT 扫描方式对 ART 算法图像重建的影响研究 761 其中，t，r分别表示测试模型和重建模型中第 i,j,ki,j,k i 层、j 行、k 列的体素密度; t 表示测试模型密度 的平均值，r 表示重建模型的密度平均值;,为图像 相似系数，表示重建所得图像与原始图像模型之间的相似程度，其值越大，表示其重建图像质量越好; d 为归一化均方根距离测量值，它表示重建后图像 真实再现测试模型图像的情况，其值越小表示其重 建图像质量越好。 2.1 不同射线源-扫描中心距对图像重建的影响 在本仿真初始程序中，D=D=source-axisaxis-detection 780 mm，扫描 1–10 次的 d 值分别为 0.47828、 0.38996、0.31866、0.261、0.21575、0.18073、0.1536、0.13235、0.1155 和 0.101919。为探讨 D对图 source-axis 像重建的影响，先取 D=D=460、540、source-axisaxis-detection 620、700、780、860、940、1020 和 1100 mm,不同 距离对图像重建影响的结果见图 2(a);再将射线源- 探测器距离固定为 1560 mm，取 D= 660、 source-axis 680、700、720、740、760、800、820、840、860、 880 和 900 mm，对图像重建影响的结果见图2(b )。以 D=D= 780 mm 时迭代 8 次时的 d source-axisaxis-detection 值(d)为基准，得到上述两类条件下迭代 8 次所 780-8 得 d 值(d)的差值，即 Δd= d– d，见图 2(a)、 i-8i-8 780-8 (b)的插图。下文差值图形均如此法进行处理。 由图 2 插图，D=D时，D source-axisaxis-detection source-axis 对图像重建的影响不是很大，迭代 8 次，d=0.13–0.14;但当 D= 1560 mm–D时迭代 8 axis-detectionsource-axis 次，d=0.13–0.2。而 D=D=780 mm source-axisaxis-detection 的重建质量最好。 图 2 D对图像重建的影响 source-axis Fig.2 Effect of Don image reconstruction. source-axis (a) D=D= 460–1100 mm, (b) D= source-axisaxis-detection source-axsi2.2 投影数对图像重建的影响 660–900 mm and D= 1560 mm –D axis-detection source-axis 选取松弛参数为 0.25，迭代 8 次，扫描时，每 转动 1?扫描一次。投影数从6 0–360 每 20 个投影数 为一间隔时的重建结果分析图如图 3 所示。可见随 着投影角度的增加，其重建图像误差逐渐减小，相 似性不断增大。投影数大于 180，d 值减小幅度变缓， 几乎成一水平线。故实际应用中，投影数在 180 左 右即可。 2.3 投影数一定时，扫描方式对图像重建的影响 取投影数为 60、120、180，投影数为 60 时， 扫描间隔分别为 1?、1.5?、2?、2.5?、3?、3.5?、4?、 4.5?、5?、5.5?、6?;投影数为 120 时，扫描间隔分 别为 1?、1.5?、2?、2.5?、3?;投影数为 180 时，扫 图 3 不同投影数图像重建分析图 Fig.3 Analysis of image reconstructing as different projections. 描间隔分别为 1?、1.5?、2?。 投影数为 60，扫描间隔为 3?迭代 1–10 次所得 0.25373、0.21583、0.18585、0.16201、0.14286 和 d 值分别为 0.58835、0.51332、0.43952、0.37278、0.12727，图 4(b)给出了各扫描间隔所得结果的 d 值 0.31687、0.27176、0.23574、0.20688、0.18356 和和差值分析(插图)。 ,值见图投影数为 180，迭代 1–10 次所得 d 和0.16447，图 4(a)给出了各扫描间隔所得结果的d 值 和差值分析(插图)。 4(c)。从各图可知，投影数为60 、120、180 扫描间 投影数为 120，扫描间隔为 2?迭代 1–10 次所得隔分别为 3?、1.5?–2.0?、1?时其重建质量较好。 d 值分别为 0.50035、0.42618、0.35938、0.30133、 图 4 投影数为 60(a)和 120 (b)时，不同扫描间隔的图像重建结果，以及投影数为 180(c)时不同扫描间隔的d 与,值 Fig.4 Analysis of image reconstruction with (a) 60 and (b) 120 projections at under different scanning intervals, and d and , values of reconstructed image with 180(c) projections under scanning intervals of 1º, 1.5º and 2.5º. 0.44196、0.38044、0.32956、0.28819、0.25456、0.22705 2.4 缺角度下扫描方式对图像重建的影响 0.20432，, 值分别为 0.76141、0.81127、0.85867、 和 设 330?到 30?之间不能扫描。CT 检测时其探测 0.89659、0.92419、0.94353、0.95704、0.96661、0.97350 器与射线源保持相对位置一致，则对应的 150?到 和 0.97858，投影数为 60、120、180 的扫描采用表210?也不能扫描。取松弛参数为0 .25。投影数为 240 2 所示的多种方式。图 5 给出了这些扫描方式与投 时，扫描间隔 1?，扫描 120 次，迭代 1–10 次所得 影数为 240 时所得结果的 d 值及 , 值的差值。 结果的 d 值分别为 0.65339、0.58678、0.51281、 表 2 扫描方式参数(每 n 度扫描一次，共 m 次) Table 2 Parameters of scanning mode (m scans at an interval of n degree per scan) 扫描角60 projections 120 projections 180 projections Scan Angles n m n m n m n m n m n m n m n m n m n m 30?–150? 1 60 2 60 2 30 3 30 4 30 1 120 2 60 1 60 1 120 1 120 210?–330? – – – – 2 30 3 30 4 30 1 120 2 60 1 60 1 60 2 60 标号 Symbols 601 602 603 604 605 1201 1202 1203 1801 1802 由图 5，在缺失一定角度不能扫描的情况下， 重建质量较好;投影数为 180 时，在 30?–150?扫描 1?，在 210?–330?扫描间隔为 2?时为最好。 间隔为 投影数越多，重建图像的质量一般越好，如投影数 180 的重建图像质量比投影数为 60 与 120 的要好。 由此可以看出，在投影数一定情况下，投影数以 同时，投影数为 60 时，不同的扫描方式重建的图像 1?–2?角度间隔均匀填满第一角度间隔后，剩余投影 也有很大差别，在 30?–150?间扫描间隔为 2?时，图 数均匀填满第二角度间隔时其重建质量较好。图 6 像重建效果最好;投影数为 120 时，在 30?–150?扫 为 602、1202、1802 扫描方式重建图像，可见重建 描间隔为 2?，在 210?–330?扫描间隔为 2?时，图像 图像随投影数的增加其质量越好。 第 10 期 伍建辉等:锥束 CT 扫描方式对 ART 算法图像重建的影响研究 763 图 5 缺角度情况下不同扫描方式重建图像分析图(重建图像的 d 差值(a)，,差值(b)) The value of d (a) and , (b) of image reconstruction at different modes of missing angle scan. Fig.5 Model 602 1202 1802 (a) (b) (c) 图 6 不同投影数下不同扫描方式重建图 (a)Z 轴切片 48, (b)Y 轴切片 64, (c)X 轴切片 50 Fig.6 Reconstructed image as different scanning mode in different projections. th th th Reconstructed image of (a) 48slice in the Z-axis, (b)48slice in the Z-axis, (c) 50slice in the X-axis 于 180?范围之内。缺角度时，在条件许可的情况下， 3 小结 应尽量增加投影数，同时使投影数尽可能地以 1?–2? 角度间隔均匀填满第一角度间隔，然后尽可能地使 综上所述，不同的扫描方式，其图像重建质量剩余投影数均匀分布于第二角度间隔之内。由此可 不同。在本仿真过程中，当射线源到扫描中心的距 知，其 CT 的扫描方式对 ART 图像重建具有重要的 离与扫描中心到探测器的距离相等时，其重建图像 影响。 较好。所以在实际中，应尽量找到射线源到扫描中 参考文献 心与扫描中心到探测器距离在重建图像质量最好时 的比例关系。同时也应尽量考虑射线源强度对图像 先武, 李时光, 王珏. 最佳无损检测手段——工业 CT 1 重建质量的影响，兼顾这两方面的因素选择最好的 技术的发展[J]. 光电工程, 1995, 22(4):52–57 距离。在投影数一定的情况下，使投影数均匀分布 XIAN Wu, LI Shiguang, WANG Jue. Development of the 刘勇奎, 沈红, 石教英. 一个有效的沿三维直线的体素 optimal nondestructive testing approach-industrial CT[J]. 8 遍历整数算法[J]. 计算机学报. 2002, 25(11):1257–1262 Opto-Electronic Engineering, 1995, 22(4):52–57 叶云长. 计算机层析成像检测[M]. 北京: 机械工业出 Liu Yongkui, Shen Hong, Shi Jiaoying. An efficient 2 版社, 2006 integer algorithm for traversing voxels along a 3D line[J]. YE Yunchang. Computed Tomography Test[M]. Beijing: Chinese Journal of Computers. 2002, 25(11):1257–1262 Mechanical Industry Press, 2006 锥束 ART 算法快速图 张顺利, 张定华, 黄魁东, 等. 9 张雷. 代数类三维 CT 图像重建算法及其硬件加速技术 像重建[J]. 仪器仪表学报. 2009, 30(4):887-8923 研究:硕士学位[D]，清华大学, 2007 ZHANG Shunli, ZHANG Dinghua, HUANG Kuidong, et ZHANG Lei. The study of three-dimensional CT image al. Fast image reconstruction with cone-beam ART reconstruction using algebraic methods and its hardware algorithm[J]. Chinese Journal of Scientific Instrument. acceleration: M. S. Dissertation[D], Tsinghua University, 2009, 30(4):887–892 2007 庄天戈. CT 原理与算法[M]. 上海:上海交通大学出版 10 刘远. CT 代数重建的并行加速和图像质量优化技术研 社, 1992 4 究:硕士学位论文[D], 西北工业大学, 2006 ZHANG Tiange. Principle and algorithm of computed LIU Yuan. Research on techniques of CT ART's parallel tomography[M]. Shanghai:Shanghai Jiaotong University acceleration and optimization of image quality: M. S. Press, 1922 Dissertation[D], Northwestern Polytechnical University, 刘远朋, 张丰收. 锥束计算机断层成像系统[M]. 北京: 11 2006 知识产权出版社, 2008, 199. 刘春华, 潘晋孝. 锥束 ART 算法的并行运算实现[J]. LIU Yuanpeng, ZHANG Fengshou. Cone beam computed CT 理论与应用研究. 2007, 16(3):13 , 19 5 tomography system[M]. Beijing:Intellectual Property LIU Chunhua, PAN Jinxiao. The study of cone-beam Right Press, 2008, 199 ART algorithm using parallel computing[J], CT Theory 刘泽, 孙丰荣, 李艳玲, 等. 基于 3D Shepp-Logan 头部 and Applications. 2007, 16(3):13 –19 12 模型的三维医学图像重建仿真[J].生物医学工程学杂志. Siddon R L. Fast calculation of the exact radiological path 2006, 23(5): 938–943 for a three dimensional CT array[J]. Med Phys. 1986, 6 LIU Ze, SUN Fengrong, LI Yanling, et al. The 12:252–255 three-dimension medical image reconstruction simulation 王宏钧, 路宏年, 傅健. 代数重建技术中投影序列选择 on 3D Shepp-Logan head phantom[J]. Journal of 次序的研究[J]. 光学技术. 2006, 32(3):389–395 7 Biomedical Engineering. 2006, 23(5): 938–943 WANG Hongjun, LU Hongnian, FU Jian. Research of Hengyong Yu, Shiying Zhao, Ge Wang. A differentiable projection access order in algebraic reconstruction 13 Shepp - Logan phantom and its applications in exact cone- technique[J]. Optical Technique. 2006, 32(3):389–395 beam CT[J]. Phys Med Biol. 2005, 50:5583–5595 Effect of scanning mode of cone beam CT on the image reconstruction of ART algorithm WU Jianhui LI Gongping ZHOU Long PAN Xiaodong LIN Qiaolu (School of nuclear science and technology, Lanzhou University, Lanzhou 730000, China) Abstract Cone beam computed tomography (CT) based on flat detector has interested many CT researchers. Different CT scanning modes have different effect on ART image reconstruction. In this paper, using a 128 ， 128 ， 128 3D Shepp-Logan head model as the object of simulation, impacts on image reconstruction are studied with different geometries of the radiation source, scanning center and detectors; different number of projections; and different scanning modes under fixed projections. An object to which scans cannot be done at certain angles is investigated, too. The following results are obtained. The image reconstruction is better when the source-center and the center-detector distances are in a certain proportion. Under fixed projections, the reconstructed image quality is better with a uniform 180º scan. To obtain a better image reconstruction in occasions where some angles are unable to scan, the projections should fill the first interval angle range uniformly with the angle interval 1?–2?, and the remaining projections should fill the second interval angle range uniformly. Key words Cone beam CT, Scanning mode, ART algorithm, 3D Shepp-Logan head model, Quality of reconstructed image CLC R814.42