CN102068271A - 一种基于呼吸相位的胸部或腹部ct图像的回顾性分类方法 - Google Patents

Abstract

本发明涉及一种基于呼吸相位的胸部或腹部CT图像的回顾性分类方法，该方法的步骤为：测量预设呼吸周期，设置CT扫描的参数；进行电影模式CT扫描获得图像，并给每张图片赋予床位号、层位号和相位号；对每一图片进行阈值分割，获得人体体表轮廓；利用人体体表轮廓的像素点的高度值建立每个相位的轮廓矩阵；累加每个相位的轮廓矩阵中元素的高度值，在每个床位中挑出累加和最大值所对应的一组图片作为一个完整周期的图片阵列；以任一床位为基准床位，对两相位间轮廓矩阵的每一列向量进行三次样条曲线光滑拟合，再计算两相位相邻层位间轮廓高度差分的距离和，挑出距离和最小的图片序号构建分类矩阵；选择呼吸周期差异程度最小的分类矩阵作为回顾性分类结果。

Description

一种基于呼吸相位的胸部或腹部CT图像的回顾性分类方法

技术领域

本发明涉及一种医学CT图像的处理方法，具体涉及根据呼吸相位对人体CT图像进行分类的方法，该方法所得到的CT图像可用于制订肿瘤放射治疗计划。

背景技术

计算机断层成像技术，即CT(Computed Tomography)是肿瘤放射治疗必不可少的工具和技术。4D-CT是一种评价肿瘤和器官运动的新技术，即在CT扫描的过程中，同步记录代表病人呼吸运动的时间信息，所以在扫描结果中除三维空间信息外还包括了病人呼吸运动的时间信息。它可以实现对运动肿瘤形状的动态观察，为研究肿瘤的运动提供了条件。对比与3D-CT，4D-CT不仅可以减少器官的运动伪影而且可以体现肿瘤和正常器官的运动信息，根据不同病人的呼吸运动的个体差异性进行个体化的放疗计划设计。该技术应用于放疗定位、计划制定及实施中，从而进入图象引导的四维放射治疗时代。

现代医学研究结果显示，影响肿瘤放射治疗时精确定位的众多因素中，呼吸运动最为显著。因此，CT图像的回顾性分类是将所得的各床位的CT图像按一个呼吸周期的时相进行归类。然而，要进行上述分类就需要在带PRM系统的计算机断层成像设备上进行。由于带PRM系统的计算机断层成像设备具有一套获得两个荧光标记点的光学系统(PRM系统)，因此我们且不说它的复杂性和设备造价高的不足，就肿瘤放射治疗时精确定位来说还存在下述缺陷：1、由于制订肿瘤放射治疗计划的CT扫描与计划的执行往往不是在同一张床上进行，因此难免存在换床的误差；2、即使是在同一张床上进行，由于所述的光学系统中的具有荧光标记点的塑料盒放在人体的不同位置也会影响实际放疗效果；3、CT扫描的过程是将患者置于热塑固定体模内进行的，患者呼吸过程中的起伏运动受到所述体模的限制，使得PRM系统不能准确地测出呼吸幅度，进而不能正常重建4D-CT。

对于带PRM(Realtime Position Management，RPM)系统的计算机断层成像设备重建4D-CT方法(如Keall P J，Starkschall G，Shukla H，Forster K M，Ortiz V，Stevens C W，Vedam S S，George R，Guerrero T andMohan R 2004Acquiring 4D thoracic CT scans using a multislice helical method Phys.Med.Biol.492053-67)，在CT扫描时，RPM系统与CT机相互通信，通过CT机的出射线信号实现CT扫描与跟踪信号的联动，获得每张CT图像成像时两个荧光标记点的高度信号，这个信号作为呼吸特征信号并绘制出每个床位的呼吸特征曲线，对呼吸特征曲线中两个波峰和中间的波谷分别设为0°，180°，360°相位角度，其余相位通过线性插值获得相位角度实现4D CT回顾性分类。带PRM(Realtime Position Management，RPM)系统的计算机断层成像设备除了上述已经讨论的缺陷外，在该方法中涉及的呼吸周期被认为是在扫描前测得的固定不变的，没有考虑扫描状态时呼吸频率的变化的事实，无法挑选出扫描状态时周期变化差异最小的4D-CT分类，该方法实现的回顾性分类也就不能精确反映扫描状态时呼吸周期的规律，其不准确性是显然的。

此外，张书旭等提出了一种在多螺旋CT上重建4D-CT的方法(中国组织工程研究与临床康复，第14卷，第26期，2010-06-25，作者：张书旭等)，该方法如下所述：分割出读入的Cine模式扫描CT图像并进行目标组织的体积计算，然后根据目标组织随呼吸过程变化，确定每一张图片在一个完整呼吸周期中的相位，最后将各床位的CT图像按相位归类。由以上描述可见，上述方法没考虑呼吸频率的变化，即将所得到的各床位的图像按一个固定周期的相位进行分类，因此难免出现某一图像的分类相位与实际相位不符。众所周知，受自身因素和外部条件的影响，人的每一个呼吸周期是不等的，因此用不变的呼吸周期对每一床位的图像进行呼吸相位分类，所造成的不准确性是可预见的。故，将采用上述方法所得到的回顾分类的图像用于制订肿瘤放射治疗计划，同样存在不精确性。

发明内容

鉴于现有技术存在上述不足，本发明所要解决的技术问题是提高基于呼吸相位的CT图像的回顾性分类的准确性。

现代医学形态学研究结果表明，受呼吸运动影响的人体器官，虽然每一个呼吸周期是不等的，但是其胸部扩张和收缩的变化曲线近似于正弦波，且每一个周期的变化曲线基本相似。基于上述理论，本发明解决上述问题的技术方案如下所述：

一种基于呼吸相位的胸部或腹部CT图像的回顾性分类方法，该方法由以下步骤组成：

(1)首先对需要进行CT扫描的对象进行呼吸周期的时长标定，以对象的平均呼吸周期为预设呼吸周期，然后根据CT扫描器的转速，设定一个床位CT扫描器转动的转数，使一个床位中探测器所得到的CT图片张数大于探测器在一个预设呼吸周期内获得的CT图片张数，小于探测器在两个预设呼吸周期内获得的CT图片张数；

(2)启动CT扫描器对人体的胸部或腹部进行电影模式CT扫描，并将所获得的每一张CT图片按人体进入扫描区的先后顺序给CT图像的每一图片赋予床位号和层位号，按扫描成像的时间顺序给CT图像的每一图片赋予相位号；

(3)将CT图像的每一图片进行阈值分割并获取人体体表轮廓曲线；

(4)以通过人体轮廓曲线在垂直方向的最低点的水平线为基准，对每一图片自左向右依次竖向逐点扫描，当图片的一列像素点的像素值全为零时记为零，当图片的一列由上而下有一个像素点的像素值不为零时，便记下第一个像素值不为零的像素点的高度值；然后将所记录的数据按自左向右的顺序排列，再以所得到的排列为行，按床位层位号的顺序由上向下排列，得到每一个床位的轮廓矩阵；

(5)将每一矩阵中所有元素的高度值累加求和，从一个床位中挑出所求和的最大值所对应的一个相位的一组图片，然后以该组图片为起点将相位号大于该组图片相位号的每一组图片按相位号由小到大排列，得到每一个床位内一个完整呼吸周期的图片阵列；

(6)选择任一床位为基准床位，以其完整呼吸周期作为回顾性分类的呼吸周期构建分类矩阵，该矩阵的构建方法是，自被选中的那床位起向上或向下依次以前一床位为基准，逐列从本床位对应的半个完整呼吸周期内的图片阵列中挑出一列相邻层位间轮廓高度的差分的距离和最小的一列图片并按相位号的顺序排，将所得到图片阵列按床的顺序上下排列，分别得到一个与基准床位内一个完整呼吸周期的图片阵列的列数相等的矩阵；其中，所述的相邻层位间轮廓高度的差分的距离和的计算方法如下所述：

(6.1)从所述的轮廓矩阵T中挑出每个床的一个完整呼吸周期轮廓矩阵，以所挑出的每一个轮廓矩阵的每列为一个列向量，并对每一列向量中相邻两点进行拟合，建立所挑出的每一个轮廓矩阵的每个列向量中相邻两点间的三次样条函数；

(6.2)构建每个列向量中的光滑样条曲线的目标函数，然后对每一目标函数求导，并令其为零，再结合自然边界条件求解出每个三次样条函数的多项式系数，得到每一拟合点的拟合值，然后把每一拟合点的拟合值分别代入所对应的每个床的一个完整呼吸周期轮廓矩阵，得到相应的拟合矩阵；

(6.3)累加相邻床位的两拟合矩阵间两相邻层位所有列的二阶前向差分的绝对值，求出两床位间每两个相位相邻层位间轮廓高度差分的距离和；

(7)以所构建的分类矩阵的每一行图片所对应的相位号为一行依次排列，然后删除相位号完全相同的行，得到一个以相位号排列的阵列，再用该阵列的每一列所对应的相位号分别减去第一列所对应的相位号，每一分类矩阵便得到一个第一列相位号为零的相位矩阵；然后用基准床位所对应的那一行的相位号减去相位矩阵中其它床位的每一行的对应相位号，再把所有差值的绝对值累加，所得和除以本相位矩阵的列数，比较所得商的大小，选择所得商最小(商最小代表相位分类的周期变化差异程度最小，相位分类也就最能反映扫描状态时呼吸周期的规律)的相位矩阵对应的分类矩阵为回顾性分类结果。

本发明所述的方法根据CT扫描过程中CT扫描对象呼吸周期的实际变化，计算以每个床位为基准床位获得的分类矩阵的周期变化差异程度，选择周期变化差异程度最小的分类矩阵作为回顾性分类结果，因此所得到的分类矩阵结果最能代表对象在进行扫描时的实际呼吸周期，准确性显著提高。此外，本发明所述的方法不是使用肺体积、空气量等高度概括的呼吸特征信号，而是对两床位中两个相位轮廓矩阵的每一列轮廓曲线向量进行分析处理，因此就能清楚地知道所得到的4D-CT图像中伪影的分布情况，为重新修订肿瘤放射治疗计划提供必要的依据和信息。

附图说明

图1为本发明所述回顾性分类的流程图。

图2为进行分割并提取后的人体轮廓曲线图。

图3为对人体轮廓曲线进行竖向逐点扫描示意图。

图4为CT图片的相位与呼吸曲线的对应关系示意图，图中，A、B、C表示床位号，水平箭头表示时间坐标，正弦波曲线表示呼吸曲线，垂直于的竖线表示CT扫描器旋转一周所得到的一组图片(下述实施例中为8张)，图(a)、(b)和(c)上方的阿拉伯数字分别表示本床位中该组图片所对应的相位号，图(a)、(b)和(c)下方的角度分别表示本床位中该组图片所对应的相位角。

图5为一个完整呼吸周期内CT图片的相位与呼吸曲线的对应关系示意图，图中，A、B、C表示床位号，水平箭头表示时间坐标，正弦波曲线表示呼吸曲线，垂直于的竖线表示CT扫描器旋转一周所得到的一组图片(下述实施例中为8张)，图(a′)、(b′)和(c′)上方的阿拉伯数字分别表示本床位中该组图片所对应的相位号，图(a′)、(b′)和(c′)下方的角度分别表示本床位中该组图片所对应的相位角。

具体实施方式

本例的实施对象为一患肺癌的病人，所使用的CT扫描器为GE lightspeed 16排CT，该CT机每转一周的时长为0.5秒，每一张图片的大小为512×512。如图1所示，所述患肺癌的病人的CT图像的回顾性分类过程如下：

(1)首先让患者躺在CT床上，用呼吸机测出其平均呼吸周期为5秒，并以此为预设呼吸周期，然后拟定扫描计划如下：扫描3个床位，每床8个层位，每层的层厚为1.5mm，每一个预设呼吸周期CT扫描器转10圈，每床扫描器转16圈，采用电影模式扫描。

(2)启动CT扫描器对人体的胸部进行电影模式CT扫描，并将所获得的每一张CT图片，按人体进入扫描区的先后顺序给每一图片赋予床位号和层位号，按扫描成像的时间顺序给每一图片赋予相位号；然后，以A、B和C表示三个床位的床位号，以每一层图片按时序排列为行，并按人体进入扫描区的先后顺序由上向下排列为列，便得到以下三个图片矩阵：

$\begin{array}{cccccc}{A}_{\mathrm{1,1}}& A_{\mathrm{1,2}}& K& A_{1,j}& K& A_{\mathrm{1,16}}\\ A_{\mathrm{2,1}}& A_{\mathrm{2,2}}& K& A_{2,j}& K& A_{\mathrm{2,16}}\\ M& M& M& A_{i,j}& M& M\\ A_{\mathrm{8,1}}& A_{8,2}& K& A_{8,j}& K& A_{\mathrm{8,16}}\end{array}---\left(A\right)$

$\begin{array}{cccccc}{B}_{\mathrm{1,1}}& B_{\mathrm{1,2}}& K& B_{1,j}& K& B_{\mathrm{1,16}}\\ B_{\mathrm{2,1}}& B_{\mathrm{2,2}}& K& B_{2,j}& K& B_{\mathrm{2,16}}\\ M& M& M& B_{i,j}& M& M\\ B_{\mathrm{8,1}}& B_{8,2}& K& B_{8,j}& K& B_{\mathrm{8,16}}\end{array}---\left(B\right)$

$\begin{array}{cccccc}{C}_{\mathrm{1,1}}& C_{\mathrm{1,2}}& K& C_{1,j}& K& C_{\mathrm{1,16}}\\ C_{\mathrm{2,1}}& C_{\mathrm{2,2}}& K& C_{2,j}& K& C_{\mathrm{2,16}}\\ M& M& M& C_{i,j}& M& M\\ C_{\mathrm{8,1}}& C_{8,2}& K& C_{8,j}& K& C_{\mathrm{8,16}}\end{array}---\left(C\right)$

三个图片矩阵中，床位号的第一个下标表示层位号，即8层图像按人体进入扫描区的先后顺序排列的序号，床位号的第二个下标表示相位号，即CT扫描器转动圈数的时序编号，如A_i，j表示第一个床位第i层位第j相位的图片。三个床位共得到384张大小为512×512的图片。

(3)将步骤(2)所得到的384张图片分别进行阈值分割并提取人体体表轮廓曲线，其中一张图片的人体体表轮廓曲线如图2所示。

(4)参见图3，以通过人体轮廓曲线在垂直方向的最低点的水平线为基准，对每一图片自左向右依次竖向逐点扫描，当图片的一列像素点的像素值全为零时记为零，当图片的一列由上而下有一个像素点的像素值不为零时，便记下第一个像素值不为零的像素点的高度值，即该像素点至所述水平线的垂直距离；然后将所记录的数据t按自左向右的顺序排列，再将所得到的排列以床位为单位按层位号的顺序由上向下排列，一个床位便得到16个大小为8×512的轮廓矩阵T，三个床位共得到48个大小为8×512的轮廓矩阵T。以某一床位第j个相位的图片为例，该相位的图片对应的轮廓矩阵可表达为：

该矩阵中t_i，ch为第i层位图片中第ch列的数据。

(5)将每一矩阵T中所有元素的高度值累加求和，从一个床位中挑出一个最大值所对应的一组图片，然后以该组图片为起点将相位号大于该组图片相位号的每一组图片按相位号由小到大排列，得到每一个床位内一个完整呼吸周期的图片阵列。

如图4所示，假设A床位得到的16个相位的图片与呼吸曲线的对应关系如图4中(a)图所示，相位1-16分别有8张图片，分别对应于步骤(2)中矩阵(A)中的一列；假设B床得到的16个相位的图片与呼吸曲线的对应关系如图4中(b)图所示，相位1-16分别有8张图片，分别对应于步骤(2)中矩阵(B)中的一列；假设C床得到的16个相位的图片与呼吸曲线的对应关系如图4中(c)图所示，相位1-16分别有8张图片，分别对应于步骤(2)中矩阵(C)中的一列。那么经过本步骤后的结果如图5所示，其中，A床位参与分类的图片的相位号为3～11的9组图片，具体见下面的矩阵(A′)；B床位参与分类的图片的相位号为1～11的11组图片，具体见下面的矩阵(B′)；C床位参与分类的图片的相位号为4～8的5组图片，具体见下面的矩阵(C′)。

$\begin{array}{cccccc}{A}_{\mathrm{1,3}}& A_{\mathrm{1,4}}& K& A_{1,j}& K& A_{\mathrm{1,11}}\\ A_{\mathrm{2,3}}& A_{\mathrm{2,4}}& K& A_{2,j}& K& A_{\mathrm{2,11}}\\ M& M& M& A_{i,j}& M& M\\ A_{8,3}& A_{8,4}& K& A_{8,j}& K& A_{\mathrm{8,11}}\end{array}---\left(A^{\′}\right)$

$\begin{array}{cccccc}{B}_{\mathrm{1,1}}& B_{\mathrm{1,2}}& K& B_{1,j}& K& B_{\mathrm{1,11}}\\ B_{2,1}& B_{2,2}& K& B_{2,j}& K& B_{\mathrm{2,11}}\\ M& M& M& B_{i,j}& M& M\\ B_{8,1}& B_{8,2}& K& B_{8,j}& K& B_{\mathrm{8,11}}\end{array}---\left(B^{\′}\right)$

$\begin{array}{cccccc}{C}_{\mathrm{1,4}}& C_{\mathrm{1,5}}& K& C_{1,j}& K& C_{\mathrm{1,8}}\\ C_{2,4}& C_{2,5}& K& C_{2,j}& K& C_{2,8}\\ M& M& M& C_{i,j}& M& M\\ C_{8,4}& C_{8,5}& K& C_{8,j}& K& C_{8,8}\end{array}---\left(C^{\′}\right)$

(6)选择任一床位为基准床位，以其完整呼吸周期作为回顾性分类的呼吸周期构建分类矩阵，该矩阵的构建方法是，自基准床位起向上或向下依次以前一床位为基准，逐列从本床位对应的半个完整呼吸周期内的图片阵列中挑出一列相邻层位间轮廓高度的差分的距离和最小的一列图片并按相位号的顺序排列，将所得到图片阵列按床的顺序上下排列，分别得到一个与基准床位内一个完整呼吸周期的图片阵列的列数相等的矩阵。例如，

以A床位为基准床位所构建的分类矩阵为：

$\begin{array}{ccccccccc}{A}_{\mathrm{1,3}}& A_{\mathrm{1,4}}& A_{\mathrm{1,5}}& A_{\mathrm{1,6}}& A_{\mathrm{1,7}}& A_{\mathrm{1,8}}& A_{\mathrm{1,9}}& A_{\mathrm{1,10}}& A_{\mathrm{1,11}}\\ A_{\mathrm{2,3}}& A_{\mathrm{2,4}}& A_{\mathrm{2,5}}& A_{\mathrm{2,6}}& A_{\mathrm{2,7}}& A_{\mathrm{2,8}}& A_{\mathrm{2,9}}& A_{\mathrm{2,10}}& A_{\mathrm{2,11}}\\ M& M& M& M& A_{i,7}& M& M& M& M\\ A_{\mathrm{8,3}}& A_{\mathrm{8,4}}& A_{8,5}& A_{\mathrm{8,6}}& A_{\mathrm{8,7}}& A_{\mathrm{8,8}}& A_{\mathrm{8,9}}& A_{\mathrm{8,10}}& A_{\mathrm{8,11}}\\ B_{\mathrm{1,1}}& B_{\mathrm{1,2}}& B_{\mathrm{1,3}}& B_{\mathrm{1,5}}& B_{\mathrm{1,6}}& B_{\mathrm{1,7}}& B_{\mathrm{1,8}}& B_{\mathrm{1,10}}& B_{\mathrm{1,11}}\\ B_{\mathrm{2,1}}& B_{\mathrm{2,2}}& B_{\mathrm{2,3}}& B_{\mathrm{2,5}}& B_{2,6}& B_{\mathrm{2,7}}& B_{2,8}& B_{\mathrm{2,10}}& B_{\mathrm{2,11}}\\ M& M& M& M& B_{i,6}& M& M& M& M\\ B_{\mathrm{8,1}}& B_{\mathrm{8,2}}& B_{\mathrm{8,3}}& B_{\mathrm{8,5}}& B_{\mathrm{8,6}}& B_{\mathrm{8,7}}& B_{\mathrm{8,8}}& B_{\mathrm{8,10}}& B_{\mathrm{8,11}}\\ C_{\mathrm{1,4}}& C_{\mathrm{1,5}}& C_{\mathrm{1,5}}& C_{\mathrm{1,6}}& C_{\mathrm{1,6}}& C_{\mathrm{1,6}}& C_{\mathrm{1,7}}& C_{\mathrm{1,7}}& C_{\mathrm{1,8}}\\ C_{\mathrm{2,4}}& C_{\mathrm{2,5}}& C_{\mathrm{2,5}}& C_{\mathrm{2,6}}& C_{\mathrm{2,6}}& C_{\mathrm{2,6}}& C_{\mathrm{2,7}}& C_{\mathrm{2,7}}& C_{\mathrm{2,8}}\\ M& M& M& M& C_{i,6}& M& M& M& M\\ C_{\mathrm{8,4}}& C_{\mathrm{8,5}}& C_{\mathrm{8,5}}& C_{\mathrm{8,6}}& C_{\mathrm{8,6}}& C_{\mathrm{8,6}}& C_{\mathrm{8,7}}& C_{\mathrm{8,7}}& C_{\mathrm{8,8}}\end{array}---\left(I\right)$

以B床位为基准床位所构建的分类矩阵为：

$\begin{array}{ccccccccccc}{A}_{\mathrm{1,3}}& A_{\mathrm{1,4}}& A_{\mathrm{1,4}}& A_{\mathrm{1,5}}& A_{\mathrm{1,6}}& A_{\mathrm{1,7}}& A_{\mathrm{1,8}}& A_{\mathrm{1,8}}& A_{\mathrm{1,9}}& A_{\mathrm{1,10}}& A_{\mathrm{1,11}}\\ A_{\mathrm{2,3}}& A_{\mathrm{2,4}}& A_{\mathrm{2,4}}& A_{\mathrm{2,5}}& A_{\mathrm{2,6}}& A_{\mathrm{2,7}}& A_{\mathrm{2,8}}& A_{\mathrm{2,8}}& A_{\mathrm{2,9}}& A_{\mathrm{2,10}}& A_{\mathrm{2,11}}\\ M& M& M& M& A_{i,6}& M& M& M& M& M& M\\ A_{\mathrm{8,3}}& A_{\mathrm{8,4}}& A_{\mathrm{8,4}}& A_{\mathrm{8,5}}& A_{\mathrm{8,6}}& A_{\mathrm{8,7}}& A_{\mathrm{8,8}}& A_{\mathrm{8,8}}& A_{\mathrm{8,9}}& A_{\mathrm{8,10}}& A_{\mathrm{2,11}}\\ B_{\mathrm{1,1}}& B_{\mathrm{1,2}}& B_{\mathrm{1,3}}& B_{\mathrm{1,4}}& B_{\mathrm{1,5}}& B_{\mathrm{1,6}}& B_{\mathrm{1,7}}& B_{\mathrm{1,8}}& B_{\mathrm{1,9}}& B_{\mathrm{1,10}}& B_{\mathrm{1,11}}\\ B_{\mathrm{2,1}}& B_{\mathrm{2,2}}& B_{\mathrm{2,3}}& B_{\mathrm{2,4}}& B_{\mathrm{2,5}}& B_{\mathrm{2,6}}& B_{\mathrm{2,7}}& B_{\mathrm{2,8}}& B_{\mathrm{2,9}}& B_{\mathrm{2,10}}& B_{\mathrm{2,11}}\\ M& M& M& M& B_{i,5}& M& M& M& M& M& M\\ B_{\mathrm{8,1}}& B_{\mathrm{8,2}}& B_{\mathrm{8,3}}& B_{\mathrm{8,4}}& B_{\mathrm{8,5}}& B_{\mathrm{8,6}}& B_{\mathrm{8,7}}& B_{\mathrm{8,8}}& B_{\mathrm{8,9}}& B_{\mathrm{8,10}}& B_{\mathrm{8,11}}\\ C_{\mathrm{1,4}}& C_{\mathrm{1,4}}& C_{\mathrm{1,5}}& C_{\mathrm{1,5}}& C_{\mathrm{1,6}}& C_{\mathrm{1,6}}& C_{\mathrm{1,6}}& C_{\mathrm{1,7}}& C_{\mathrm{1,7}}& C_{\mathrm{1,8}}& C_{\mathrm{1,8}}\\ C_{\mathrm{2,4}}& C_{\mathrm{2,4}}& C_{\mathrm{2,5}}& C_{\mathrm{2,5}}& C_{\mathrm{2,6}}& C_{\mathrm{2,6}}& C_{\mathrm{2,6}}& C_{\mathrm{2,7}}& C_{\mathrm{2,7}}& C_{\mathrm{2,8}}& C_{\mathrm{2,8}}\\ M& M& M& M& C_{i,6}& M& M& M& M& M& M\\ C_{\mathrm{8,4}}& C_{\mathrm{8,5}}& C_{\mathrm{8,5}}& C_{\mathrm{8,5}}& C_{\mathrm{8,6}}& C_{\mathrm{8,6}}& C_{\mathrm{8,7}}& C_{\mathrm{8,7}}& C_{\mathrm{8,8}}& C_{\mathrm{8,8}}& C_{\mathrm{8,8}}\end{array}---\left(\mathrm{II}\right)$

以C床位为基准床位所构建的分类矩阵为：

$\begin{array}{ccccc}{A}_{\mathrm{1,3}}& A_{\mathrm{1,5}}& A_{\mathrm{1,7}}& A_{\mathrm{1,9}}& A_{\mathrm{1,11}}\\ A_{\mathrm{2,3}}& A_{\mathrm{2,5}}& A_{\mathrm{2,7}}& A_{\mathrm{2,9}}& A_{\mathrm{2,11}}\\ M& M& A_{i,7}& M& M\\ A_{\mathrm{8,3}}& A_{\mathrm{8,5}}& A_{\mathrm{8,7}}& A_{\mathrm{8,9}}& A_{\mathrm{8,11}}\\ B_{\mathrm{1,1}}& B_{\mathrm{1,4}}& B_{\mathrm{1,6}}& B_{\mathrm{1,9}}& B_{\mathrm{1,11}}\\ B_{\mathrm{2,1}}& B_{\mathrm{2,4}}& B_{\mathrm{2,6}}& B_{\mathrm{2,9}}& B_{\mathrm{2,11}}\\ M& M& B_{i,6}& M& M\\ B_{\mathrm{8,1}}& B_{\mathrm{8,4}}& B_{\mathrm{8,6}}& B_{\mathrm{8,9}}& B_{\mathrm{8,11}}\\ C_{\mathrm{1,4}}& C_{\mathrm{1,5}}& C_{\mathrm{1,6}}& C_{\mathrm{1,7}}& C_{\mathrm{1,8}}\\ C_{\mathrm{2,4}}& C_{\mathrm{2,5}}& C_{\mathrm{2,6}}& C_{\mathrm{2,7}}& C_{\mathrm{2,8}}\\ M& M& C_{i,6}& M& M\\ C_{\mathrm{8,4}}& C_{\mathrm{8,5}}& C_{\mathrm{8,6}}& C_{\mathrm{8,7}}& C_{\mathrm{8,8}}\end{array}---\left(\mathrm{III}\right)$

其中，所述的相邻层位间轮廓高度的差分的距离和的计算方法如下所述：

(6.1)从所述的轮廓矩阵T中挑出每个床的一个完整呼吸周期轮廓矩阵，以所挑出的每一个轮廓矩阵的每列为一个列向量，并对每一列向量中相邻两点进行拟合，建立所挑出的每一个轮廓矩阵的每个列向量中相邻两点间的三次样条函数。本例中，所述三次样条函数的具体建立过程如下：

将步骤(4)中所述的48个轮廓矩阵中的每一个轮廓矩阵的每列作为一个列向量，共得到48×512个列向量，其中，某一床位第j个相位所对应的轮廓矩阵T_j中的第ch列的列向量为

t_ch＝(t_1，ch，t_2ch，.....，t_8，ch)′；那么对每一列向量中相邻两点进行拟合，即可得到48×512×7条样条曲线，其中，任一列向量中任两点间的光滑样条曲线可表达为

s_i，ch(x)＝d_i，ch(x-p_i，ch)³+c_i，ch(x-p_i，ch)²+b_i，ch(x-p_i，ch)¹+a_i，ch所示的三次样条函数，式中，x表示i与i+1层两图片中人体轮廓曲线上边缘上任一对应点间样条线上任一点至相位号相同且层位号为1的图片间的距离，P表示任一层图片至相位号相同且层位号为1的图片间的距离，a、b、c和d分别表示待求的三次样条函数的多项式系数。

(6.2)构建每个列向量中的光滑样条曲线的目标函数，然后对每一目标函数求导，并令其为零，再结合自然边界条件求解出每个三次样条函数的多项式系数，得到每一拟合点的拟合值，然后把每一拟合点的拟合值分别代入所对应的每个床的一个完整呼吸周期轮廓矩阵，得到相应的拟合矩阵。本例中，所述拟合矩阵的具体获得过程如下：

为求解上述的a、b、c和d的值，首先构建出每一光滑样条曲线的目标函数，其中，某一列向量中任两点间的光滑样条曲线s_i，ch(x)＝d_i，ch(x-p_i，ch)³+c_i，ch(x-p_i，ch)²+b_i，ch(x-p_i，ch)¹+a_i，ch的目标函数为

式中，为光滑参数(与层厚和像素点面积的大小相关，根据经验取值；本例层厚为1.5mm，像素点面积0.5830mm²，根据经验取

然后对所述的目标函数求导，并令其为零，结合三次样条函数的自然边界条件(即列向量中头尾两端点的样条函数的二次和三次多项式系数为0，也就是d_1，ch＝c_1，ch＝d_7，ch＝c_7，ch＝0)求解出a、b、c和d，再把所得到的a、b、c和d的值代入所构建的每一光滑样条曲线，得到每一拟合点的拟合值，然后把每一拟合点的拟合值代入相应的列向量t_ch＝(t_1，ch，t_2ch，.....，t_8，ch)′，得到相应的轮廓拟合向量s_ch＝(s_1，ch(p_1，ch)，s_1，ch(p_2，ch)，s_2，ch(p_3，ch)，，s_3，ch(p_4，ch)....，s_i，ch(p_i+1，ch)，...，s_7，ch(p_8，ch))′，再用所得到的轮廓拟合向量替换矩阵T中对应的列，那么，A床位便得到9个大小为8×512的拟合矩阵L，B床位便得到11个大小为8×512的拟合矩阵L，C床位便得到5个大小为8×512的拟合矩阵L，三个床位共得到25个大小为8×512的拟合矩阵L，其中，某一床位相位j所对应的拟合矩阵可表达为：

$L_j=\left(\begin{array}{cccccc}{s}_{1,1}\left(p_{\mathrm{1,1}}\right)& s_{\mathrm{1,2}}\left(p_{\mathrm{1,2}}\right)& K& s_{1,\mathrm{ch}}\left(p_{1,\mathrm{ch}}\right)& K& s_{\mathrm{1,512}}\left(p_{\mathrm{1,512}}\right)\\ s_{\mathrm{1,1}}\left(p_{\mathrm{1,2}}\right)& s_{\mathrm{1,2}}\left(p_{\mathrm{2,2}}\right)& K& s_{1,\mathrm{ch}}\left(p_{2,\mathrm{ch}}\right)& K& s_{\mathrm{1,512}}\left(p_{\mathrm{2,512}}\right)\\ s_{\mathrm{2,1}}\left(p_{\mathrm{3,1}}\right)& s_{\mathrm{2,2}}\left(p_{\mathrm{3,2}}\right)& K& s_{\mathrm{2,2}}\left(p_{\mathrm{3,2}}\right)& K& s_{\mathrm{2,512}}\left(p_{\mathrm{2,512}}\right)\\ s_{\mathrm{3,1}}\left(p_{\mathrm{4,1}}\right)& s_{\mathrm{3,2}}\left(p_{\mathrm{4,2}}\right)& K& s_{3,\mathrm{ch}}\left(p_{5,\mathrm{ch}}\right)& K& s_{7,512}\left(p_{\mathrm{8,512}}\right)\\ M& M& M& s_{i,\mathrm{ch}}\left(p_{i+1,\mathrm{ch}}\right)& M& M\\ s_{\mathrm{7,1}}\left(p_{\mathrm{8,1}}\right)& s_{\mathrm{7,2}}\left(p_{\mathrm{8,2}}\right)& K& s_{7,\mathrm{ch}}\left(p_{8,\mathrm{ch}}\right)& K& s_{\mathrm{7,512}}\left(p_{\mathrm{8,512}}\right)\end{array}\right),$

其中s_i，ch(p_i+1，ch)的值对应样条曲线s_i，ch(x)在该相位图片第i+1层第ch列的拟合值。

(6.3)累加相邻床位的两拟合矩阵间两相邻层位所有列的二阶前向差分的绝对值，求出两床位间每两个相位相邻层位间轮廓高度差分的距离和。本例中，所述两床位间每两个相位相邻层位间轮廓高度差分的距离和的具体计算过程如下：

首先分别令，前一床位中完整呼吸周期内相位号相同的任一组图片所对应的拟合矩阵为(h₁，h₂.，h_i..，h₈)′(其h_i表示前一床位中拟合矩阵L_j的第i行向量)，后一床位中完整呼吸周期内相位号相同的任一组图片所对应的拟合矩阵为(k₁，k₂.，k_i..，k₈)′(其k_i表示后一床位中拟合矩阵L_j的第i行向量)，那么累加前一床位和后一床位间两相邻层位的二阶前向差分，求出前一床位和后一床位间每两个相位相邻层位间轮廓高度的差分的距离和，(式中k_i，ch表示后一床位中拟合矩阵L_j的第i行第ch列的元素值，h_i，ch表示前一床位中拟合矩阵L_j的第i行第ch列的元素值)。

(7)以任一床位为基准床位所构建的分类矩阵的每一行图片所对应的相位号为一行依次排列，然后删除相位号完全相同的行，得到一个以相位号排列的阵列，其中，以A床位为基准床位所构建的分类矩阵得到的阵列见下述的阵列(A”)，以B床位为基准床位所构建的分类矩阵得到的阵列见下述的阵列(B”)，以C床位为基准床位所构建的分类矩阵得到的阵列见下述的阵列(C”)；再用所得到的每个阵列的每一列所对应的相位号分别减去第一列所对应的相位号，每一分类矩阵便得到一个第一列相位号为零的相位矩阵，然后用基准床位所对应的那一行的相位号减去相位矩阵中其它床位的每一行的对应相位号，这样以A床位为基准床位所构建的分类矩阵得到的阵列(A”)、(B”)和(C”)便相应变为(A^*)、(B^*)，和(C^*)；再把每个阵列中所有差值的绝对值累加，所得和除以相位矩阵的列数，比较所得商的大小(商最小代表相位分类的周期变化差异程度最小，相位分类也就最反映呼吸周期规律)，选择所得商最小的相位矩阵对应的分类矩阵为回顾性分类结果。本例中，以A床位为基准床位所构建的分类矩阵对应所得的商为2.89，以B床位为基准床位所构建的分类矩阵对应所得的商为3.91，以C床位为基准床位所构建的分类矩阵对应所得的商为5.20，可见以A床位为基准床位所构建的分类矩阵对应所得的商最小，故选择上述分类矩阵(I)，即以A床位为基准床位所构建的分类矩阵，作为回顾性分类结果。

$\begin{array}{ccccccccc}3& 4& 5& 6& 7& 8& 9& 10& 11\\ 1& 2& 3& 5& 6& 7& 8& 10& 11\\ 4& 5& 5& 6& 6& 6& 7& 7& 8\end{array}---\left(A^{,,}\right)$

$\begin{array}{ccccccccccc}3& 4& 4& 5& 6& 7& 8& 8& 9& 10& 11\\ 1& 2& 3& 4& 5& 6& 7& 8& 9& 10& 11\\ 4& 5& 5& 5& 6& 6& 7& 7& 8& 8& 8\end{array}---\left(B^{,,}\right)$

$\begin{array}{ccccc}3& 5& 7& 9& 11\\ 1& 4& 6& 9& 11\\ 4& 5& 6& 7& 8\end{array}---\left(C^{,,}\right)$

$\begin{array}{ccccccccc}0& 1& 2& 3& 4& 5& 6& 7& 8\\ 0& 0& 0& -1& -1& -1& -1& -2& -2\\ 0& 0& 1& 1& 2& 3& 3& 4& 4\end{array}---\left(A^{*}\right)$

$\begin{array}{ccccccccccc}0& 0& 1& 1& 1& 1& 1& 2& 2& 2& 2\\ 0& 1& 2& 3& 4& 5& 6& 7& 8& 9& 10\\ 0& 0& 1& 2& 2& 3& 3& 4& 4& 5& 6\end{array}---\left(B^{*}\right)$

$\begin{array}{ccccc}0& -1& -2& -3& -4\\ 0& -2& -3& -5& -6\\ 0& 1& 2& 3& 4\end{array}---\left(C^{*}\right).$

在本例中，由于CT扫描对象是患肺癌的病人，因此其胸部表面轮廓的起伏变化与患者吸气和呼气的过程完全吻合。基于上述事实，本例中的步骤(4)和(5)即找到了患者在扫描中每一个床位的一个完整呼吸周期中吸气末端的相位(见图5)，那么通过步骤(5)得到的每一个床位内一个完整呼吸周期的图片阵列的第一列的起点相同，即同处于一个完整呼吸周期的吸气末端[见矩阵(A′)、(B′)和(C′)]，从为步骤(6)所得到的分类矩阵(I)、(Ⅱ)和(Ⅲ)之间具有可比性奠定了基础；步骤(6)在步骤(5)所得到的矩阵(A′)、(B′)和(C′)的基础上，分别自基准床位起向上或向下依次以前一床位为基准，逐列从本床位对应的半个完整呼吸周期内的图片阵列中挑出一列相邻层位间轮廓高度的差分的距离和最小的一列图片并按相位号的顺序排列，从而保证了相邻层位间轮廓的平滑性，可显著减少重建图象中的伪影；步骤(7)又对步骤(6)所得到的分类矩阵(Ⅰ)、(Ⅱ)和(Ⅲ)中找出与基准床位呼吸周期差异程度最小的分类矩阵(Ⅰ)作为最终分类结果，从进一步缩小了A、B、C三个床之间呼吸周期差异程度(结合图4比较分类矩阵(Ⅰ)、(Ⅱ)和(Ⅲ)也可直观地看出)。

综上所述，本发明所得到回顾性分类结果，不仅每一个床位内一个完整呼吸周期的图片阵列的第一列的起点同处于一个完整呼吸周期中吸气末端，而且相邻层位间轮廓平滑，各床之间呼吸周期差异程度小，因此其精确性明显高于现有技术。

Claims (1)

1.一种基于呼吸相位的胸部或腹部CT图像的回顾性分类方法，该方法由以下步骤组成：

(1)首先对需要进行CT扫描的对象进行呼吸周期的时长标定，以对象的平均呼吸周期为预设呼吸周期，然后根据CT扫描器的转速，设定一个床位CT扫描器转动的转数，使一个床位中探测器所得到的CT图片张数大于探测器在一个预设呼吸周期内获得的CT图片张数，小于探测器在两个预设呼吸周期内获得的CT图片张数；

(2)启动CT扫描器对人体的胸部或腹部进行电影模式CT扫描，并将所获得的每一张CT图片按人体进入扫描区的先后顺序给CT图像的每一图片赋予床位号和层位号，按扫描成像的时间顺序给CT图像的每一图片赋予相位号；

(3)将CT图像的每一图片进行阈值分割并获取人体体表轮廓曲线；

(4)以通过人体轮廓曲线在垂直方向的最低点的水平线为基准，对每一图片自左向右依次竖向逐点扫描，当图片的一列像素点的像素值全为零时记为零，当图片的一列由上而下有一个像素点的像素值不为零时，便记下第一个像素值不为零的像素点的高度值；然后将所记录的数据按自左向右的顺序排列，再以所得到的排列为行，按床位层位号的顺序由上向下排列，得到每一个床位的轮廓矩阵；

(5)将每一矩阵中所有元素的高度值累加求和，从一个床位中挑出所求和的最大值所对应的一个相位的一组图片，然后以该组图片为起点将相位号大于该组图片相位号的每一组图片按相位号由小到大排列，得到每一个床位内一个完整呼吸周期的图片阵列；

(6)选择任一床位为基准床位，以其完整呼吸周期作为回顾性分类的呼吸周期构建分类矩阵，该矩阵的构建方法是，自被选中的那床位起向上或向下依次以前一床位为基准，逐列从本床位对应的半个完整呼吸周期内的图片阵列中挑出一列相邻层位间轮廓高度的差分的距离和最小的一列图片并按相位号的顺序排列，将所得到图片阵列按床位号的顺序上下排列，分别得到一个与基准床位内一个完整呼吸周期的图片阵列的列数相等的矩阵；其中，所述的相邻层位间轮廓高度的差分的距离和的计算方法如下所述：

(6.1)从所述的轮廓矩阵T中挑出每个床的一个完整呼吸周期轮廓矩阵，以所挑出的每一个轮廓矩阵的每列为一个列向量，并对每一列向量中相邻两点进行拟合，建立所挑出的每一个轮廓矩阵的每个列向量中相邻两点间的三次样条函数；

(6.2)构建每个列向量中的光滑样条曲线的目标函数，然后对每一目标函数求导，并令其为零，再结合自然边界条件求解出每个三次样条函数的多项式系数，得到每一拟合点的拟合值，然后把每一拟合点的拟合值分别代入所对应的每个床的一个完整呼吸周期轮廓矩阵，得到相应的拟合矩阵；

(6.3)累加相邻床位的两拟合矩阵间两相邻层位所有列的二阶前向差分，求出两床位间每两个相位相邻层位间轮廓高度差分的距离和；

(7)以所构建的分类矩阵的每一行图片所对应的相位号为一行依次排列，然后删除相位号完全相同的行，得到一个以相位号排列的阵列，再用该阵列的每一列所对应的相位号分别减去第一列所对应的相位号，每一分类矩阵便得到一个第一列相位号为零的相位矩阵；然后用基准床位所对应的那一行的相位号减去相位矩阵中其它床位的每一行的对应相位号，再把所有差值的绝对值累加，所得和除以本相位矩阵的列数，比较所得商的大小，选择所得商最小的相位矩阵对应的分类矩阵为回顾性分类结果。

Images