CN102376097A - 非对称检测器ct迭代重建方法 - Google Patents

Abstract

本发明涉及一种非对称检测器CT迭代重建方法，其特征在于包括以下步骤：从CT扫描机获得的生数据经过空气校正、非线性校正、串扰校正、散焦校正以及定标校正获得投影数据；判断上述投影数据为断层扫描数据还是螺旋扫描数据；如为断层扫描数据，则对断层扫描数据进行分组；计算投影权重，得到权值矩阵；根据上述权值矩阵进行OSEM迭代正、反投影；判断是否满足迭代终止条件，如满足，则终止迭代正、反投影，得到初始CT图像；对上述初始CT图像进行后处理，得到重建的CT图像。本发明方法直接使用不完整投影数据进行迭代反投影重建，可以减轻伪影，抑制由插值造成的数据不一致性带来的噪声水平差异，主要适用于四排以下CT检测器。

Description

非对称检测器CT迭代重建方法

技术领域

本发明涉及医学图像处理技术，具体的说是一种非对称检测器CT迭代重建方法。

背景技术

在CT扫描中，为了增大扫描视野或减少辐射剂量，通常在检测器中心通道的一侧增加或者减少部分检测器模块，从而形成非对称检测器阵列，如图1所示，减少了部分检测器；相对于对称检测器结构，图中的非对称结构相当于有部分投影数据缺失，导致重建视野缩小。如图2所示，扇束几何重组为平行束几何，其中FOV_In为[-γ_m+α，γ_m-α]所能够覆盖的视野范围，FOV_Out为完整扫描视野与FOV_In之间的视野部分。该图仅为示意图，一个通道和它的对面通道不一定能完全对应，要补偿缺失部分的数据，目前的处理方法都是通过插值，利用其对面通道数据(Complementarydata)反投影，但这种方法在数据不一致性较大时会引入伪影，需要寻找适当的平滑过渡函数来抑制。

在统计计算中，最大期望(EM，Expectation-Maximization)算法是在概率模型中寻找参数最大似然估计的算法，其中概率模型依赖于无法观测的隐藏变量。最大期望算法经过两个步骤交替进行计算，第一步是计算期望(E)，也就是将隐藏变量像能够观测到的一样包含在内从而计算最大似然的期望值；另外一步是最大化(M)，也就是最大化在E步上找到的最大似然的期望值从而计算参数的最大似然估计。M步上找到的参数然后用于另外一个E步计算，这个过程不断交替进行。

传统的MLEM(Maximum Likelihood Expectation Maximization)方法计算式为：

${\mathrm{\μ}}^{\left(k\right)}\left(j\right)=\frac{{\mathrm{\μ}}^{(k-1)}\left(j\right)}{\underset{i}{\mathrm{\Σ}}w(i,j)}\underset{i}{\mathrm{\Σ}}\frac{w(i,j)P\left(i\right)}{\underset{j}{\mathrm{\Σ}}w(i,j^{\′}){\mathrm{\μ}}^{(k-1)}\left(j^{\′}\right)}---\left(1\right)$

其中：

μ^(k)(j)：第j个像素在第k次迭代中的衰减系数

P(i)：第i个检测器单元的测量投影值

w(i，j)：第j个像素对第i个投影的贡献权重

MLEM方法旨在寻找与测量的投影数据具有最大似然性(ML)的估计解，其迭代过程是由EM算法来实现的。由于是以统计规律为基础，MLEM重建法具有很好的抗噪声能力，尤其是在处理统计性差的数据时，更能显示出它相对于解析法的优越性，但是这种方法仍然存在迭代法的运算量大，运算时间长等缺点。

发明内容

针对现有技术中在CT扫描中心通道两侧检测器不对称的情况下(去掉一侧的部分通道)，无法获得完整的一圈反投影数据的不足之处，本发明要解决的技术问题是提供一种可以减轻伪影和噪声水平差异的非对称检测器CT迭代重建方法。

为解决上述技术问题，本发明采用的技术方案是：

本发明一种非对称检测器CT迭代重建方法包括以下步骤：

从CT扫描机获得的生数据经过空气校正、非线性校正、串扰校正、散焦校正以及定标校正获得投影数据；

判断上述投影数据为断层扫描数据还是螺旋扫描数据；

如为断层扫描数据，则对断层扫描数据进行分组；

计算投影权重，得到权值矩阵；

根据上述权值矩阵进行OSEM迭代正、反投影；

判断是否满足迭代终止条件，如满足，则终止迭代正、反投影，得到初始CT图像；

对上述初始CT图像进行后处理，得到重建的CT图像。

如果不满足迭代终止条件，则返回根据上述权值矩阵进行迭代正、反投影步骤。

如为螺旋扫描数据，则对螺旋扫描数据经z方向滤波后重组到建像平面后再进行分组；

接续计算投影权重，得到权值矩阵步骤。

所述分组是对一圈投影数据中0～90度范围内的投影数据进行间隔方式的子集分组，对子集中的每一个角度做正、反投影的同时进行与该角度相隔90度、180度以及270度的其他三个角度的正、反投影。

所述OSEM迭代正、反投影为：在每一次迭代过程中将投影数据分成L个子集，每一个子集对重建图像各像素点值校正后，重建图像被更新一次，所有的子集运算一遍。

所述OSEM迭代正、反投影为直接使用不完整投影数据进行迭代正、反投影重建。

本发明具有以下有益效果及优点：

1.本发明方法直接使用不完整投影数据进行迭代反投影重建，可以减轻伪影，抑制由插值造成的数据不一致性带来的噪声水平差异，主要适用于四排以下CT检测器。

附图说明

图1为现有技术中减少部分检测器的非对称检测器阵列；

图2为现有技术中非对称平行束几何示意图；

图3为本发明方法流程图；

图4为本发明方法涉及的加权函数示意图；

图5为本发明方法涉及的平行束几何示意图

图6为本发明方法中FBP缺失数据一圈反投影；

图7为本发明方法中FBP缺失数据插值成半圈反投影；

图8为本发明方法中OSEM缺失数据一圈反投影。

具体实施方式

如图3所示，本发明非对称检测器CT迭代重建方法包括以下步骤：

从CT扫描机获得的生数据(Raw Data)经过空气校正、非线性校正、串扰校正、散焦校正以及定标校正获得投影数据；

判断上述投影数据为断层扫描数据还是螺旋扫描数据；

如为断层扫描数据，则对断层扫描数据进行分组；

计算投影权重，得到权值矩阵；

根据上述权值矩阵进行OSEM迭代正、反投影；

判断是否满足迭代终止条件，如满足，则终止迭代正、反投影，得到初始CT图像；

对上述初始CT图像进行后处理，得到重建的CT图像。

如果不满足迭代终止条件，则返回根据上述权值矩阵进行迭代正、反投影步骤；

如为螺旋扫描数据，则对螺旋扫描数据经z方向(即三维坐标系中的z轴方向)滤波后重组到建像平面后再进行分组；

接续计算投影权重，得到权值矩阵步骤。

1.投影数据的重组及分组

对于断层扫描数据，不需要进行重组，只需进行分组；而对于螺旋扫描数据，通过传统的z方向滤波的方法将不同z位置的投影数据重组到建像平面。

本发明采用正交的子集排序方法，将0～90度范围内的投影数据进行间隔方式的子集分组(将投影划分为N个子集，每隔M/N个投影选取一个投影，其中M为0～90度范围内的投影数)，每次进行四个相互正交角度(即与该角度相隔90度、180度以及270度)的迭代正、反投影。

2.OSEM迭代正、反投影

为了加快收敛速度，减少运算时间，提高图像质量，人们提出了很多快速算法，其中OSEM(Ordered Subsets Expectation Maximization有序子集期望最大化)是很有应用前景的一种快速迭代重建算法，它是在最大似然期望法MLEM的基础上发展起来的。

所述OSEM迭代正、反投影为直接使用不完整投影数据进行迭代正、反投影重建，替代通过插值补充数据后进行滤波反投影重建的方法，可以抑制由数据不一致性带来的噪声水平差异，减轻伪影。

OSEM方法属于“块”迭代法中的一种，在每一次迭代过程中将投影数据分成L个子集，每一个子集对重建图像各像素点值校正以后，重建图像便被更新一次，所有的子集运算一遍，称为一次迭代过程。在MLEM方法一次迭代过程中，使用所有的投影数据对重建图像每一个像素点的值进行校正，重建图像被更新一次，而在OSEM方法一次迭代过程中重建图像被更新L次，所以OSEM方法具有加快收敛的作用。

令

为初始图像，

为k次迭代后的图像，P₁，P₂......P_L为全部投影数据划分而成的有序子集，具体算法过程如下：

k＝0，初始化

重复下面步骤直到迭代到指定迭代次数为止。

A)

k＝k+1；

B)对子集l＝1，2，.....L计算投影：计算重建图像在每个子集对应的投影角度下探测器的期望正投影值

${\hat{P}}^l=\underset{j}{\mathrm{\Σ}}w(i,j){\mathrm{\μ}}^l\left(j\right),i\∈P_l;---\left(2\right)$

w(i，j)为投影权重，μ^l(j)为像素衰减值；

对每一个像素计算：

${\mathrm{\μ}}^{(l+1)}\left(j\right)={\mathrm{\μ}}^{\left(l\right)}\left(j\right)\underset{i\∈P_l}{\mathrm{\Σ}}\frac{P_l\·w(i,j)}{{\hat{P}}^l\left(i\right)}/\underset{i\∈P_l}{\mathrm{\Σ}}w(i,j)---\left(3\right)$

C) ${\widehat{\mathrm{\μ}}}^k={\mathrm{\μ}}^{L+1}$

式(3)对i的求和是在第l个子集内进行的。

3.计算投影权重

本发明采用传统的函数加权方式：

射线γ_i进入物体并在经过的所有点处采样该物体(如图4所示)，由于每一点并非都与图像网格点重合，所以就要以每一采样点为中心的加权函数h(x，y)来加权计算采样点值，落在加权函数h范围中的网格点x_j都通过加权函数适当加权而结合进来。

w(i，j)＝∫h(X(s_i)-X(x_j)，Y(s_i)-Y(x_j))ds_i                   (4)

采用贝塞尔-凯塞函数作为加权函数，其方程为：

$h^m\left(d\right)=\frac{{\left(\sqrt{1-{(d/r)}^2}\right)}^m{I}_m\left(\mathrm{\α}\sqrt{1-{(d/r)}^2}\right)}{I_m\left(\mathrm{\α}\right)},0\≤d\≤r---\left(5\right)$

式中：I_m为第一类修正贝塞尔函数，m为贝塞尔函数的阶数，r为函数半径，d为距函数中心点距离(从像素x_j到射线P_i的最短距离)，α为控制函数形状的参数。

权值的闭合形式为：

$w(i,j)=r\sqrt{2\mathrm{\π}/\mathrm{\α}}\·{(1-{(d/r)}^2)}^{m+1/4}{I}_{m+1/2}\left(\mathrm{\α}\sqrt{1-{(r/a)}^2}\right)I_m\left(\mathrm{\α}\right)---\left(6\right)$

非对称反投影几何示意图如图5所示。被扫描物体的阴影部分在投影P1下未被照射到，而在P2下被照射到，如果采取一圈滤波反投影的方式重建，这部分CT值由于滤波采样的不均衡偏低。图6是一个典型的水模非对称反投影的例子，去掉中心通道左侧部分通道，同时水模偏心放置，此时模体中的一部分不能在全部投影下都被完全照射到，因此造成CT值明显偏低。

要解决此问题，可以将数据重组到完整半圈，需进行插值，用对面通道数据补偿缺失数据，同时要在截断处需做平滑处理，否则由于数据的不一致性会引入明显的伪影，如图7所示：

如果用OSEM迭代重建方法进行缺失数据的一圈反投影，由于它是使用EM算法求参数的最大似然估计，而不是投影采样贡献累加，投影不均衡而造成部分CT值偏低的现象消失，同时也不需要进行插值，可以一定程度上抑制数据不一致带来的伪影及噪声差异，非对称OSEM迭代反投影图像结果如图8所示。

图像后处理包括图像定标、去环、去金属伪影、骨硬化校正等。

Claims (6)

1.一种非对称检测器CT迭代重建方法，其特征在于包括以下步骤：

从CT扫描机获得的生数据经过空气校正、非线性校正、串扰校正、散焦校正以及定标校正获得投影数据；

判断上述投影数据为断层扫描数据还是螺旋扫描数据；

如为断层扫描数据，则对断层扫描数据进行分组；

计算投影权重，得到权值矩阵；

根据上述权值矩阵进行OSEM迭代正、反投影；

判断是否满足迭代终止条件，如满足，则终止迭代正、反投影，得到初始CT图像；

对上述初始CT图像进行后处理，得到重建的CT图像。

2.按权利要求1所述的非对称检测器CT迭代重建方法，其特征在于：

如果不满足迭代终止条件，则返回根据上述权值矩阵进行迭代正、反投影步骤。

3.按权利要求1所述的非对称检测器CT迭代重建方法，其特征在于：

如为螺旋扫描数据，则对螺旋扫描数据经z方向滤波后重组到建像平面后再进行分组；

接续计算投影权重，得到权值矩阵步骤。

4.按权利要求1或3所述的非对称检测器CT迭代重建方法，其特征在于：所述分组是对一圈投影数据中0～90度范围内的投影数据进行间隔方式的子集分组，对子集中的每一个角度做正、反投影的同时进行与该角度相隔90度、180度以及270度的其他三个角度的正、反投影。

5.按权利要求1所述的非对称检测器CT迭代重建方法，其特征在于：所述OSEM迭代正、反投影为：在每一次迭代过程中将投影数据分成L个子集，每一个子集对重建图像各像素点值校正后，重建图像被更新一次，所有的子集运算一遍。

6.按权利要求1所述的非对称检测器CT迭代重建方法，其特征在于：所述OSEM迭代正、反投影为直接使用不完整投影数据进行迭代正、反投影重建。

Images