用动态胸部数字仿真模型检测图像重建算法性能的方法

摘要

一种用动态胸部数字仿真模型检测图像重建算法性能的方法，包括以下步骤：(1)用椭球和椭圆柱体排列组合构成三维动态胸部数字仿真模型；(2)获取仿真投影数据；(3)用需要检测的重建算法，对第(2)步中获得的仿真投影数据进行图像重建；(4)利用第(3)步的重建算法的按实际需要，计算各方面性能指标检测重建算法的重建性能；(5)改变模型中的参数，重复(1)－(4)步，进一步检测重建算法对动态物体在不同参数条件下性能变化情况。本发明对动态物体的重建性能，以确定重建算法是否适用于动态物体重建或对算法的部分参数进行调整始之适用于实际情况，或进一步确定不同算法间的优劣。

说明书

技术领域

本发明涉及的是一种生物医学成像技术领域的方法，具体地说，是一种用动态胸部数字仿真模型检测图像重建算法性能的方法。

背景技术

图像重建算法性能有：精度、密度分辨率、空间分辨率、时间分辨率、对噪声的敏感性、伪迹程度、鲁棒性等。在检测图像重建算法性能时，直接用人体或小动物测试，因人体或小动物内的密度分布的正确值无法事先知道，重建值无法与之比较，重建算法的许多性能指标无法定量表示。所以一般采用事先设定好密度分布的数字模型或实物模型，对其仿真投影、再用仿真投影重建，来检测重建算法性能。模型可分为二维和三维两种，最常见的如Shepp-Logan三维头部模型还有胸部、腹部等模型。

经对现有技术的文献检索发现，典型的采用事先设定好密度分布的胸部数字模型来检测重建算法性能的方法，如Jed D Pack和Frédéric Noo，Cone-beamreconstruction using 1D filtering along the projection of M-lines[使用延M线投影的一维滤波的锥形束重建算法]，Inverse Problem[逆问题]，Vol.21(2005)，pp.1105-1120。其中使用的模型为Katia Sourbelle提出的The FORBILD CT-simulation phantoms[FORBILD CT仿真模型]，Proc.of the1999 Int.Meeting on Fully 3D Image Reconstruction[1999年完全三维重建方法国际会议]。FORBILD仿真模型是由椭球，长方体，圆柱体及半环形圆管道等规则几何结构，通过空间的分布组合构成心，肺，骨骼及血管等人体各组织三维结构。采用FORBILD仿真模型来检测重建算法性能的方法的不足在于：FORBILD模型是静态的，采用此模型，不能检测重建算法对动态物体(如心脏)的重建性能。这与现今针对动态成像、小动物成像的算法研究日益重视的背景很不相称。

发明内容

本发明的目的在于针对现有技术的不足，提供一种用动态胸部数字仿真模型检测图像重建算法性能的方法，使其采用三维动态胸部(心脏模型是动态的)数字模型检测重建算法对动态物体的重建性能，以确定重建算法是否适用于动态物体重建或对算法的部分参数进行调整始之适用于实际情况，或进一步确定不同算法间的优劣。

本发明是通过以下技术方案实现的，本发明采用包含动态心脏的三维动态胸部数字模型，通过对该模型计算仿真投影，再用得到的仿真投影重建，来检测重建算法性能，特别对动态物体的重建性能。本发明包括以下步骤：

(1)用椭球和椭圆柱体排列组合构成三维动态胸部数字仿真模型。所述三维动态胸部数字仿真模型，包括模拟心脏的动态椭球，各轴长分别是一个随时间变化的周期函数，分别是A＝A(t)，B＝B(t)，C＝C(t)，以及由椭球组成的肺、手臂、胸骨和腹腔，由椭圆柱组成的肱骨，由椭球和椭圆柱合并而成的脊椎骨。

所述的用椭球和椭圆柱体排列组合构成三维动态胸部数字仿真模型，其具体过程如下：

设该数字模型整体为一长方体结构，定义其三维直角坐标系的原点位于该长方体的中心，X轴、Y轴、Z轴分别延长方体各棱，并以LX、LY、LZ分别代表长方体X轴、Y轴和Z轴的各棱长。以下所表示的各几何结构的参数均以该坐标系内归一化数值描述。

所述归一化数值a是指， $>>a>=>>P>L>>,>$>其中P指是指坐标系中某一分量的实际坐标值，L指该分量P所属的坐标轴的半棱长。

●定义空间椭球方程为： $>>>>x>2>>>A>2>>>+>>>y>2>>>B>2>>>+>>>z>2>>>C>2>>>\le >1>,>$>一般情况下，椭球中心坐标位于(EX，EY，EZ)，A、B、C分别代表椭球各轴长，椭球在X-Z平面内轴A与X轴夹角为α。

●定义空间椭圆柱方程为： $>\left(\begin{array}{}>>>>x>2>>>A>2>>>+>>>y>2>>>B>2>>>\le >1>,>>>>|>z>|>\le >>Z>Lim>>>\end{array}\right)$>一般情况下，椭圆中心坐标位于(CX，CY，CZ)，A、B分别代表椭圆各轴长，椭圆柱在X-Z平面内轴A与X轴夹角为β。

●以椭圆柱建立身体部分，椭圆柱中心坐标(0，0，0)；轴长A＝0.6，B＝0.35；Z_Lim＝1.0；β＝0°；密度d＝1.00。

所述的密度d代表各组织的相对平均密度结构，不同的密度将反映不同的灰度级，这里设身体的密度为1.00，其它各组织密度参考此密度给出。

●以椭球建立左手臂部分，椭球中心坐标(0.8，0，0)；轴长A＝0.17，B＝0.17，C＝2.2；α＝4°；d＝1.00。

●以椭球建立右手臂部分，椭球中心坐标(-0.8，0，0)；轴长A＝0.17，B＝0.17，C＝2.2；α＝-4°；d＝1.00。

●以内外相套的两椭圆柱建立左臂肱骨部分，外椭圆柱：中心坐标(0.79，0，0)；轴长A＝0.06，B＝0.06；Z_Lim＝2.2；β＝4°；d＝1.50。内椭圆柱：中心坐标(0.79，0，0)；轴长A＝0.04，B＝0.04；Z_Lim＝2.2；β＝4°；d＝0.98。

●以内外相套的两椭圆柱建立右臂肱骨部分，外椭圆柱：中心坐标(-0.79，0，0)；轴长A＝0.06，B＝0.06；Z_Lim＝2.2；β＝-4°；d＝1.50。内椭圆柱：中心坐标(-0.79，0，0)；轴长A＝0.04，B＝0.04；Z_Lim＝2.2；β＝-4°；d＝0.98。

●以椭球建立左肺部分，椭球中心坐标(0.3，0.01，0.08)；轴长A＝0.25，B＝0.2，C＝0.8；α＝3°；d＝0.26。

●以椭球建立右肺部分，椭球中心坐标(-0.28，0.01，0.08)；轴长A＝0.25，B＝0.2，C＝0.8；α＝-5°；d＝0.26。

●以内外相套的两椭球建立胸骨部分，外椭球：中心坐标(0，0.25，0)；轴长A＝0.05，B＝0.02，C＝0.6；α＝0°；d＝1.25。内椭球：中心坐标(0，0.25，0)；轴长A＝0.04，B＝0.01，C＝0.55；α＝0°；d＝0.98。

●以内外相套的两椭圆柱及两相切的椭球共同构成一节脊椎骨，外椭圆柱：中心坐标(0，-0.2，Z)；轴长A＝0.05，B＝0.05；Z_Lim＝0.06；β＝0°；d＝1.92。内椭圆柱：中心坐标(0，-0.2，Z)；轴长A＝0.03，B＝0.03；Z_Lim＝0.04；β＝0°；d＝1.18。椭球1：中心坐标(0，0.26，Z)；轴长A＝0.08，B＝0.01，C＝0.06；α＝0°；d＝1.92。椭球2：中心坐标(0，0.29，Z)；轴长A＝0.01，B＝0.02，C＝0.06；α＝0°；d＝1.92。

●按上步中的方法在不同的坐标位置Z处，建立互不重叠的若干脊椎骨。

●以椭球建立腹腔部分，椭球中心坐标(0，0.07，-1)；轴长A＝0.52，B＝0.2，C＝0.7；α＝0°；d＝0.95。

●以椭圆球建立心脏部分，椭球中心坐标(0.05，0，0.07)；α＝30°；d＝1.10；为满足心脏的动态变化，这里椭球各半径将是关于时间变化的一个周期函数，A＝A(t)，B＝B(t)，C＝C(t)。

(2)根据需要检测的CT(计算机断层成像术)，SPECT(单光子发射计算机断层成像术)，PET(正电子发射计算机断层成像术)等重建算法，按照算法要求的扫描模式(如在CT重建算法中，使用扇形束圆轨迹的滤波反投影重建算法和由Feidkamp，Dewis和Kress提出的基于锥束圆轨迹的重建算法，简称FDK算法)，根据实际中数据采集的原理，获取仿真投影数据。如在CT重建算法中，数据是通过X射线管采得的，其原理如公式p＝∫ρ·dl(其中，dl为空间内某一X射线经过的路径，ρ为该路径上所经过的物体密度，p则为采得的数据)，根据这一物理采集过程，用数字方法计算X射线路经，并用模型密度替换实际物体密度，对公式离散化处理后按扫描模式获得三维动态胸部数字仿真模型在不同时刻下的仿真投影。在SPECT和PET中，其仿真投影过程类似于CT，不同处在于数据发生原理和采集方式，根据具体情况将上述的原理公式替换成相对应的原理公式，并作数字化处理即可。

所述的仿真投影是根据实际投影数据采集原理，用解析式或数值积分的方法获取仿真投影数据的过程。

(3)用需要检测的CT，SPECT，PET重建算法，对第(2)步中获得的仿真投影数据进行图像重建。

(4)利用第(3)步中获得的重建结果对需要检测的CT，SPECT，PET重建算法的按实际需要，计算各方面性能指标检测重建算法的重建性能，如SSP(层灵敏度)用于检测体积CT算法的Z轴分辨率；重建图像与模型之间的方差以检测重建图像的空间分辨率，方差越小说明图像质量越高；通过在不同极端条件下重建图像与模型间方差的比较，检测重建算法的鲁棒性；通过不同时刻的重建图像与模型的误差分析，检测重建算法的时间分辨率。具体的性能数值随重建算法不同和实际需要而变化。

(5)改变模型中的用于描述心脏的动态方程的参数和仿真投影中的参数，重复(1)-(4)步，进一步检测重建算法对动态物体在不同参数条件下性能变化情况。

本发明采用包含动态心脏的三维动态胸部数字模型可检测重建算法对动态物体的重建性能的原理在于：用一个类比说明，用照相机拍摄静物，照相机用三脚架固定，不论快门开启时间是多少，拍摄的照片总是清晰的(对比度可能有差异)，就不能体现出曝光时间短的优点；如用照相机拍摄动物，在曝光量足够的条件下，曝光时间短的，拍摄的照片就清晰，否则照片会因为物体与相机的相对运动而模糊。总之，用动物代替静物可测试照相机某种设置下对动态物体成像的性能。同样，采用三维静态胸部数字模型，不论投影数据采集多快，不能反映重建算法对动态物体的优良重建性能，采用包含动态心脏的三维动态胸部数字模型可检测重建算法对动态物体的重建性能。

本发明的有益效果是：(1)由于采用包含动态心脏的三维动态胸部数字模型，可以检验重建算法对动态物体的时间分辨率，以判断被检测算法对于动态物体是否实际有效或在特定参数下有效。(2)本发明中采用的数字模型仅有两种几何机构组成，椭球和椭圆柱，在实际应用中便于计算机仿真模拟的实现，提高算法检测效率，节约时间。(3)本发明中的模型更接近人体生理结构，对算法实用性的检测更具参考价值。

附图说明

图1空间直角坐标系定义示意图

图2模型截面三视图，其中：图2-a为模型冠状截面图；图2-b为模型失状截面图；图2-c为模型横截面图

其中，1：身体；2：右臂；3：左臂；4：右肱骨；5：左肱骨；6：心脏(t＝0s时刻)；7：右肺；8：左肺；9：胸骨；10：腹腔；11：脊椎骨。

图3单源锥束圆轨迹CT扫描方式示意图

具体实施方式

为更好地理解本发明的技术方案，以下结合附图及具体的实施例作进一步描述，实施例中的被检测算法为FDK算法，该算法适用于单源锥束圆轨迹扫描方式的CT，实施例的具体步骤如下：

1.建立三维动态胸部数字模型。

(1)如图1，建立空间直角坐标系，设该数字模型整体为一长方体结构，各棱长分别为：X轴方向棱长256(像素)，Y轴方向棱长256(像素)，Z轴方向棱长140(像素)。

按下列各表参数依次建立模型各主要部位：

注：表中各参数均为归一化数值，其中表二中t的单位为秒。

                    表1  模型各主要部分(脊椎、心脏除外)参数表

  序号  名称  形状                               参数  密度  1  身体  CX  CY  CZ  A  B  Z_Lim  β  椭圆柱  0  0  0  0.6  0.35  1.0  0°  1.00  2  左臂  EX  EY  EZ  A  B  C  α  椭球  0.8  0  0  0.17  0.17  2.2  4°  1.00  3  右臂  椭球  -0.8  0  0  0.17  0.17  2.2  -4°  1.00  4  左肱骨  CX  CY  CZ  A  B  Z_Lim  β  外圆柱  0.79  0  0  0.06  0.06  2.2  4°  1.50  内圆柱  0.79  0  0  0.06  0.06  2.2  4°  0.98  5  右肱骨  外圆柱  -0.79  0  0  0.06  0.06  2.2  -4°  1.50  内圆柱  -0.79  0  0  0.06  0.06  2.2  -4°  0.98  6  左肺  EX  EY  EZ  A  B  C  α  椭球  0.3  0.01  0.08  0.25  0.2  0.8  3°  0.26  7  右肺  椭球  -0.28  0.01  0.08  0.25  0.2  0.8  -5°  0.26  8  胸骨  外椭球  0  0.25  0  0.05  0.02  0.6  0°  1.25  内椭球  0  0.25  0  0.04  0.01  0.55  0°  0.98  9  腹腔  椭球  0  0.07  -1  0.52  0.2  0.7  0°  0.95

                     表2  心脏模型参数表

  序号  名  称  形  状  β  CX  CY  CZ  10  心  脏  椭  球  30°  0.05  0  0.07  密度  A(t)  B(t)  C(t)  1.10  0.13+  0.02sin(2πt)  0.11+  0.02sin(2πt)  0.3+  0.05sin(2πt)

                                        表3  脊椎骨模型参数表

  序号  名称  形状                                      参数  密度  11  脊椎  CX  CY  CZ  A  B  Z_Lim  β  外圆柱  0  -0.2  Z  0.05  0.05  0.06  0°  1.92  内圆柱  0  -0.2  Z  0.03  0.03  0.04  0°  1.18  EX  EY  EZ  A  B  C  α  椭球1  0  0.26  Z  0.08  0.01  0.06  0°  1.92  椭球2  0  0.29  Z  0.01  0.02  0.06  0°  1.92

                                            表4  各脊椎骨Z轴坐标位置表

  脊椎骨序号  1  2  3  4  5  6  7  8  9  Z轴坐标  0.8  0.6  0.4  0.2  0.0  -0.2  -0.4  -0.6  -0.8

模型实际效果如图2中所示，按标号，依次为：1为身体；2、3分别为右、左臂；4、5分别为右、左肱骨，其中浅色部分为骨质，深色部分为髓质；6为t＝0s时刻的心脏，心脏大小将随时间周期性变化；7、8为右、左肺；9为胸骨，其中浅色部分为骨质，深色部分为髓质；10为腹腔；11为脊椎骨，其中浅色部分为骨质，深色部分为髓质。

2.按单源锥束圆轨迹CT的扫描方式(如图3所示)和原理，用公式：p＝∫ρ·dl(其中，dl为空间内某一X射线经过的路径，ρ为该路径上所经过的模型密度，p则为采得的数据)对每条X射线计算仿真投影p(x，y，θ)，(x，y)代表投影平面内某一点的位置，θ为投影平面内射线的角度参数。设锥形束旋转一圈的时间为45秒。

3.采用FDK锥形束虑波反投影重建算法，用仿真投影数据p(x，y，θ)进行直接三维体数据重建。

4.将重建后的三维体数据与t＝0s，t＝0.5s时的数字模型进行比较，发现重建结果中，心脏部分由于模型的运动，比较模糊，伪迹严重，进一步计算重建体数据与这两个时刻下模型的方差，得 $>sup>>\delta >>t>=>0>s>>2sup>>=>0.0235>,sup>>\delta >>t>=>0>.>5>s>>2sup>>=>0.0314>,>$>而一般静态物体的方差平均 $>sup>>ver>>\delta >‾>>>>s>2sup>>=>0.0001>$>左右，为说明该算法的时间分辨率低，在当前仿真参数条件下不适于动态物体的重建。

5.增大模型中描述心脏的周期函数的周期至4秒，锥形束旋转一圈的时间仍为45s，重复1-4步，发现心脏部分的模糊减小，重新计算方差，得 $>sup>>\delta >>t>=>0>s>>2sup>>=>0.009>8>,sup>>\delta >>t>=>0>.>5>s>>2sup>>=>0.0157>,>$>可见两时刻的方差都有所减小，说明重建图像更接近原模型。

6.设定描述心脏的周期函数的周期仍为1秒，锥形束旋转一圈的时间改为20秒，重复1-4步，发现心脏伪迹减小，重新计算方差，得 $>sup>>\delta >>t>=>0>s>>2sup>>=>0.0083>,sup>>\delta >>t>=>0>.>5>s>>2sup>>=>0.0132>,>$>两时刻的方差相比步骤4中的值也有所减小。从而可以确定，该算法在扫描速度低，物体运动快时，重建质量低；当扫描速度快，物体运动慢时，重建质量有所提高；当物体静止不动时，重建质量较高。

通过此例可发现，本发明可以检测算法的重建图像在不同时刻下与模型间的误差和在不同参数设置下对动态模型重建的性能变化，从而判定重建算法对动态物体的适用性。对模型的仿真投影过程中，以现有的软件和硬件，其耗时分别为，硬件：57s，软件：308s，同等条件下，对Katia Sourbelle提出的模型进行仿真投影的耗时分别为，硬件：104s，软件：632s。说明本发明在实现上效率更高。

对于SPECT和PET重建算法，也可用类似的方法进行比较分析，但在仿真投影时，须采用与数据发生和采集方式相对应的原理公式进行模拟仿真。