一种迭代演化模型的CT图像稀疏重建方法与流程

本发明涉及ct图像重建领域，特别涉及一种迭代演化模型的ct图像稀疏重建方法。

背景技术：

计算机断层成像(computedtomography，ct)是根据x射线穿透待检测物体横截面时被介质吸收衰减的现象，利用接收到的投影数据重建得到物体的断层图像。传统的重建方法解析法要求投影数据满足数据完整性要求(tuy-smith条件)，即在180°角度范围内进行均匀密集的投影，否则重建结果会存在严重伪影。然而在众多应用场景下，往往只能得到采样角度受限、采样间隔较大的稀疏投影数据。例如医疗应用的脑部ct中，过大的x射线照射剂量会对大脑造成二次伤害，因此常通过采用加大投影角度间隔的方式减小ct射线的照射剂量；在工业设备的无损探伤中，大型的不规则印制电路板或在役的管道限制了ct成像系统的几何位置，因此只能在有限角度下对目标进行扫描，得到投影角度范围受限的数据。因此，提高稀疏投影数据的重建质量是ct重建领域十分重要且亟需解决的问题。

针对ct稀疏投影数据重建问题，现有技术如凸集投影-全变分最小化(projectiononconvexsets-totalvariation，pocs-tv)算法等[1][2]根据ct图像的分片平滑特性，通过在正则化框架下引入tv约束，建立目标函数(其中a为投影系数矩阵，x为待重建图像的列向量形式，b为投影数据)的方式，将ct重建问题转化为凸优化求解问题，通过相应的凸优化迭代法求解得到重建结果。

在投影角度范围受限、投影角度间隔较大等采样条件下获得的ct稀疏投影数据重建问题上，现有技术存在重建结果伪影多、图像质量差的缺点，本发明提出一种新型迭代演化模型的ct图像稀疏重建方法，提高重建图像质量，减少重建伪影。

由于投影角度范围受限、投影角度间隔较大等采样条件下获得的投影数据对称性、完备性极低，重建算法的搜索空间巨大。现有技术如凸集投影-全变分最小化(projectiononconvexsets-totalvariation，pocs-tv)算法基于凸优化求解思路，从单一的初始解出发，按照目标函数在当前解的负梯度方向路径进行搜索。一方面重建质量受初始解质量的严重影响，在单一解上的迭代优化使得重建过程无法在有限时间内遍历巨大的搜索空间，往往只能得到局部最优解，重建图像质量差。

另外，特别是投影角度范围受限的情况下，一方面，投影数据不具有对称性，因此不同方向的梯度对图像的信息恢复作用不同，全变分约束并不能保证正确的到达全局最优的迭代路径；另一方面，pocs-tv算法采用的梯度下降法将迭代路径固定在了目标函数的负梯度方向上，在全变分约束并不能保证到达全局最优的迭代路径的前提下，也无法跳出当前方向去寻找更优的解，无法得到高质量的解。

[1]emilys,kaocm,panx.accurateimagereconstructionfromfew-viewsandlimited-angledataindivergent-beamct[j].j.x-rayscitechnol,2006,14(2):119-139.

[2]sidkyey,panx.imagereconstructionincircularcone-beamcomputedtomographybyconstrained,total-variationminimization[j].physicsinmedicineandbiology,2008,53(17):4777-4807.

技术实现要素：

本发明的主要目的在于克服现有技术的缺点与不足，提供一种迭代演化模型的ct图像稀疏重建方法。

本发明的目的通过以下的技术方案实现：

一种迭代演化模型的ct图像稀疏重建方法，包括以下步骤：

s1、采集固定投影角度范围或投影角度间隔达到设定值的稀疏投影数据；

s2、生成四个初始解：伪逆、随机值、全0值、全1值；

s3、通过随机漫步获得四组解集；

s4、对每一个解集中的每一个解进行有效性约束；

s5、依据适应度评价函数评价每组随机解集中的每一个解；

s6、按新型迭代演化模型对每一个解进行转移；

s7、根据基于凸优化的迭代算法对每一个解进行更新；

s8、依据适应度评价挑选每一组的当前最优解更新；

s9、判断是否达到最大迭代次数：若达到，则依据适应度评价选取最优解作为重建结果；若未达到，则返回步骤s3。

步骤s1中，所述固定投影角度范围为(0°,180°)，所述投影角度间隔大于等于3°小于等于6°。

所述步骤s2，具体为：

伪逆：根据离散化的投影方程ax＝b，求其伪逆作为初始解：其中t表示当前迭代次数，其中a为投影系数矩阵，x为待重建图像的列向量形式，b为投影数据；

随机值：随机生成初始解：

全0值：元素全为0的初始解

全1值：元素全为1的初始解

所述步骤s3，具体为：

对每一个初始解执行随机漫步操作，获得四组随机解集

其中，t表示当前迭代次数；i＝1,…,4，i表示解集组别；k表示第i组解集中的第k个解；rand[length(x)]表示维度与x相同的随机量，服从(0,∑)的高斯分布；通过生成随机解集，对解空间进行充分的搜索，并推动迭代过程跳出局部最优点。

所述步骤s4，具体为：

对每组随机解集中的每一个解进行有效性约束：

其中，j表示解的第j个元素；有效性约束使得解中每个元素的值都是有意义的图像像素值。

所述步骤s5，具体为：

依据适应度评价函数评价每组随机解集中的每一个解：

在适应度评价函数中，第一项为数据保真项，评价重建图像与投影数据的一致性；第二项为tv正则化约束项，用于评价重建图像的质量；适应度评价函数值越低，说明解的质量越好。

所述步骤s6，具体为：

以为当前估计值，为预测值，为观测值，k＝2,3,4,5，对预测值建立新型迭代演化模型：

由于预测值由随机漫步操作得到，无法保证其质量，但同时为了保证本方法在局部最优点的跳离能力，不能仅靠适应度评价的好坏来决定完全保留或舍弃预测值。因此本方法设计的新型迭代演化模型考虑了的引导作用以及对预测值的评价结果，将转移到更有可能收敛到全局最优解的路径上，同时保持对解空间的探索性。

所述步骤s7，具体为：

根据基于凸优化的迭代算法对所有解进行更新：

(1)将作为基于凸优化的迭代重建算法初值，按目标函数进行更新：

(2)采用art迭代算法保证数据一致性：

其中，m为投影数据总数；bm表示第m个投影数据；am表示投影稀疏矩阵a的第m行；

(3)图像非负性约束：

(4)基于梯度下降法求解tv最小化：

其中，表示中坐标为(p,q)的像素值；

(5)得到更新后的解：

所述步骤s8，具体为：

依据适应度评价挑选每一组的当前最优解：

其中，i＝1,…,4，k＝1,2,…,5。

所述步骤s9，具体为：

判断是否达到最大迭代次数：若已达到最大迭代次数，则对4个解进行适应度评价，选取最优解作为重建结果：若未达到，则令返回步骤s3。

本发明与现有技术相比，具有如下优点和有益效果：

1、本发明通过多点初始化的方法获取多个初始解，使得初始解能够有效分布在解空间的不同位置。

2、通过随机漫步生成解集，一方面使得本方法能够对巨大的解空间进行最大程度的探索，另一方面能够在解的迭代陷入局部最优时通过随机漫步的推动跳出局部最优点。

3、本发明提出新型迭代演化模型，建立估计、预测、观测模型，对解的演化过程进行转移控制，使得算法能够在保持对解空间的探索性的同时提高迭代效率。

4、本发明在投影角度范围受限、投影角度间隔较大等采样条件下获得的ct稀疏投影数据重建问题上，现有技术无法对巨大的解空间进行全面的搜索，存在重建结果伪影多、图像质量差的缺点，本发明通过设计多点初始化、随机漫步、新型迭代演化模型等策略，对解空间进行有效的探索，提高重建图像质量，减少重建伪影。

5、本发明设计不同的迭代演化模型也可以实现对解的转移控制。

附图说明

图1是本发明所述一种迭代演化模型的ct图像稀疏重建方法的流程图。

具体实施方式

下面结合实施例及附图对本发明作进一步详细的描述，但本发明的实施方式不限于此。

如图1所示，一种迭代演化模型的ct图像稀疏重建方法，包括以下步骤：

a、采集投影角度范围受限或投影角度间隔较大的稀疏投影数据：投影数据角度范围小于(0°,180°)，投影间隔大于等于3°小于等于6°。

b、生成四个初始解：

①伪逆：根据离散化的投影方程ax＝b(其中a为投影系数矩阵，x为待重建图像的列向量形式，b为投影数据)，求其伪逆作为初始解：其中t表示当前迭代次数。

②随机值：随机生成初始解：

③全0值：元素全为0的初始解

④全1值：元素全为1的初始解

c、对每一个初始解执行随机漫步操作，获得四组随机解集

其中，t表示当前迭代次数；i＝1,…,4，表示解集组别；k表示第i组解集中的第k个解。rand[length(x)]表示维度与x相同的随机量，服从(0,∑)的高斯分布。通过生成随机解集，保证优化路径的多样性，对解空间进行充分的搜索，并推动迭代过程跳出局部最优点。

d、对每组随机解集中的每一个解进行有效性约束：

其中，j表示解的第j个元素。有效性约束使得解中每个元素的值都是有意义的图像像素值。

e、依据适应度评价函数评价每组随机解集中的每一个解：

在适应度评价函数中，第一项为数据保真项，评价重建图像与投影数据的一致性；第二项为tv正则化约束项，用于评价重建图像的质量。适应度评价函数值越低，说明解的质量越好。

f、以为当前估计值，为预测值，为观测值，对预测值建立新型迭代演化模型：

由于预测值由随机漫步操作得到，无法保证其质量，但同时为了保证本方法在局部最优点的跳离能力，不能仅靠适应度评价的好坏来决定完全保留或舍弃预测值。因此本方法设计的新型迭代演化模型考虑了的引导作用以及对预测值的评价结果，将转移到更有可能收敛到全局最优解的路径上，同时保持对解空间的探索性。

g、根据基于凸优化的迭代算法对所有解进行更新：

g1、将作为基于凸优化的迭代重建算法初值，按目标函数进行更新：

g2、采用art迭代算法保证数据一致性：

其中，m为投影数据总数；bm表示第m个投影数据；am表示投影稀疏矩阵a的第m行。

g3、图像非负性约束：

g4、基于梯度下降法求解tv最小化：

其中，

表示中坐标为(p,q)的像素值。

g5、得到更新后的解：

h、依据适应度评价挑选每一组的当前最优解：

其中，i＝1,…,4，k＝1,2,…,5。

i、若已达到最大迭代次数，则进入步骤j；否则，令跳回步骤c；

j、对4个解进行适应度评价，选取最优解作为重建结果：

上述实施例为本发明较佳的实施方式，但本发明的实施方式并不受上述实施例的限制，其他的任何未背离本发明的精神实质与原理下所作的改变、修饰、替代、组合、简化，均应为等效的置换方式，都包含在本发明的保护范围之内。