一种稀疏采样角度下基于变分不等式的CT图像重建方法技术

本发明专利技术公开了一种稀疏采样角度下基于变分不等式的CT图像重建方法，包括以下步骤：(1)针对扇束CT，在0到180度角度范围内进行等角度间隔投影扫描，投影方向角度间隔的大小在2到6度之间，获取稀疏投影数据y；(2)通过X射线源，探测器以及待重建物体的位置信息计算投影矩阵A；(3)根据步骤(1)中得到的投影数据y以及步骤(2)中得到的投影矩阵A，同时引入图像梯度的稀疏性以及非负性作为先验知识得到稀疏采样角度下图像重建问题的重建模型；(4)将步骤(3)中的重建模型转化为变分不等式形式；(5)求解步骤(4)中的变分不等式，得到重建后的图像。本发明专利技术在保证图像重建质量的前提下，可以提升重建的收敛速度并降低单次迭代的时间。

【技术实现步骤摘要】
一种稀疏采样角度下基于变分不等式的CT图像重建方法
本专利技术涉及精密电子封装过程中元器件内部图像重建以及医学CT图像重建领域，特别是涉及稀疏采样角度下基于变分不等式的CT图像重建方法。
技术介绍
计算机断层成像技术(ComputedTomography,CT)是X射线照相技术与复杂的计算机信号处理方法相结合的产物，在不同的方向对物体进行X射线投影，通过测量的投影数据获得物体的横截面信息，可以准确而又直观地重构出物体内部的结构。目前CT技术已被广泛应用于安全检查，工业无损探伤以及医学诊断等领域。在CT重建的过程中，我们往往只能获得不完全的，稀疏采样角度下的投影数据。如考虑到集成电路大规模生产的需求，要求短时间内完成封装元器件的内部图像的重建并进行质量分析，因此，投影数据采集以及重建的速度对于提高集成电路检测的效率至关重要，相比于每隔1度进行一次投影，增大投影角度的间隔，采集稀疏角度的投影数据无疑会大幅度提升数据采集的速度，有助于缺陷检测过程的快速进行；在医学领域，考虑到X光辐射对人体的危害性，高剂量的X光辐射会增加人体患疾病的风险，在满足CT重建质量要求的情况下，希望能够尽可能减少投影的角度，降低X光射线穿过的次数来降低X光的辐射剂量。对于这类稀疏采样角度下图像重建问题，通常使用迭代算法进行求解，但是使用传统的迭代算法，如最大期望(ExpectationMaximization，EM)算法，有序子集最大期望(OrderedSubsetExpectationMaximization，OSEM)算法，代数重建算法(AlgebraicReconstructionTechnique，ART)、联合迭代重建算法(SimultaneousIterativeReconstructionTechnique，SIRT)并不能得到高质量的重建图像。为了提升重建的质量，在重建的过程中考虑正则化的概念，以约束、罚项等形式引入图像的先验知识，如图像梯度的稀疏性以及非负性质。然而这些考虑了正则化的迭代重建方法如凸集投影-全变分(ProjectionOnConvexSets-TotalVariation,POCS-TV)算法由于投影数据的不完全性，对重建的速度产生了重大的影响，往往需要经过很多次的迭代才能获得高质量的重建图像，收敛速度慢，单次迭代时间长。因此，对于稀疏采样角度下的图像重建问题，必须专利技术新的算法在保证重建质量的前提下，提高重建的收敛速度，降低单次迭代的时间。
技术实现思路
本专利技术的目的在于克服现有方法的缺点与不足，提供了稀疏采样角度下基于变分不等式的CT图像重建方法，能够在保证重建质量的前提下，有效提高收敛的速度，并降低单次迭代的时间。为了达到上述目的，本专利技术采用以下技术方案：一种稀疏采样角度下基于变分不等式的CT图像重建方法，包括下述步骤：(1)针对扇束CT，在0到180度角度范围内进行等角度间隔投影扫描，投影方向角度间隔的大小在2到6度之间，获取稀疏投影数据y；(2)通过X射线源，探测器以及待重建物体的位置信息计算投影矩阵A；(3)根据步骤(1)中得到的投影数据y以及步骤(2)中得到的投影矩阵A，同时引入图像梯度的稀疏性以及非负性作为先验知识进行约束得到稀疏采样角度下图像重建问题的重建模型：其中u为待重建图像，被离散化为一个n×n的图像矩阵；表示图像梯度；矩阵A中的每一个元素Ai,j表示为第i条射线穿过第j个像素的长度；N为待重建图像像素的个数；M为投影数据的个数；(4)将步骤(3)中的稀疏采样角度下图像重建问题的重建模型转化为变分不等式形式：<f(u*),u-u*>≥0(2)其中u*,u∈C，u*为变分不等式的解，F(u)＝||u||TV，f(u)＝F'(u)，C1＝{u|||Au-y||2≤ε}，C2＝{u|u≥0}，C＝C1∩C2；(5)求解步骤(4)中的变分不等式，得到重建图像，即找到一个待重建的图像u*，使得对于任意的u∈C，<f(u*),u-u*>≥0成立。优选的，步骤(3)中，||u||TV的梯度是一个图像，该图像的每个像素值都是||u||TV相对应像素的偏微分，并用如下表示：优选的，在重建的过程中，待重建图像u以及||u||TV的梯度都被转化为列向量，图像中第s行，第t列的像素点对应于列向量中第(s-1)·n+t个像素点。优选的，步骤(4)中稀疏采样角度下图像重建模型与变分不等式之间的转化，转化过程包括以下步骤：变分不等式的定义：设G:Rn→Rn为连续映射，为非空闭凸集，存在x*∈S，使得对于<G(x*),x-x*>≥0成立；(4-1)令C1＝{u|||Au-y||2≤ε}，C2＝{u|u≥0}，C＝C1∩C2，F(u)＝||u||TV，f(u)＝F'(u)公式(1)变为argminF(u)s.t.u∈C(3)其中F(u)是一个凸函数，因为C1，C2皆为凸集，所以C＝C1∩C2也为凸集(4-2)令u*为公式(3)的最优解，那么对于F(u)≥F(u*)；(4-3)对于任何的γ∈[0,1]，当u，u*属于凸集C时，u*+γ(u-u*)∈C，于是根据步骤(4-2)可得:F(u*+γ(u-u*))-F(u*)≥0(4)因此对于公式(2)成立；这表明，u*也是变分不等式，即公式(2)的解，公式(1)等价于公式(2)。优选的，步骤(5)所述的变分不等式求解方法，求解过程具体包含以下步骤：(5-1)初始化：u0＝0，d0＝1，ζ＝0.2，设置总的迭代次数NTotal，当前迭代次数k＝0；(5-2)对当前的迭代解uk应用投影算法，得到中间解(5-3)计算下一次迭代中的梯度下降法步长dk+1；(5-4)将列向量转化为二维图像矩阵二维图像矩阵中第s行，第t列的像素点与列向量中第(s-1)·n+t个像素点对应；(5-5)对进行平滑操作，得到平滑过后的结果(5-6)将转化为列向量得到uk+1，图像中第s行，第t列的像素点对应于列向量中第(s-1)·n+t个像素点；(5-7)k＝k+1，判断k是否小于设置的总迭代次数NTotal；若k小于NTotal，跳转至步骤(5-2)继续迭代；否则重建结束，将重建结果转化为二维图像矩阵。优选的，步骤(5-2)所述的投影算法，具体为：(5-2-1)pk＝uk-τf(uk)；(5-2-2)zk＝PC(pk)，其中为投影算子，计算pk在可行集C上的投影。优选的，步骤(5-2-1)所述pk的计算，pk＝uk-τf(uk)是梯度下降法，使用如下步骤来实现梯度下降法的过程，其中t为当前迭代的次数，设置固定迭代次数为NGD，dk为下降的步长，具体为：(5-2-1-1)初始化：s0＝uk,t＝0；(5-2-1-2)计算梯度：(5-2-1-3)梯度下降法：st+1＝st-dk·v/|v|；(5-2-1-4)t＝t+1，判断t是否小于设置的迭代次数NGD；若t小于NGD，跳转至步骤(5-2-1-2)继续迭代；否则迭代停止，优选的，步骤(5-2-2)所述对pk在C上投影的计算，C＝C1∩C2，由两部分组成，因此求pk在可行集C上的投影也分为两部分：先求pk在可行集C1上的投影zk，λ为正则化参数，并设置迭代次数t为当前迭代的次数；再限制zk为非负，得到具体过程为：(本文档来自技高网...

【技术保护点】
一种稀疏采样角度下基于变分不等式的CT图像重建方法，其特征在于，包括下述步骤：(1)针对扇束CT，在0到180度角度范围内进行等角度间隔投影扫描，投影方向角度间隔的大小在2到6度之间，获取稀疏投影数据y；(2)通过X射线源，探测器以及待重建物体的位置信息计算投影矩阵A；(3)根据步骤(1)中得到的投影数据y以及步骤(2)中得到的投影矩阵A，同时引入图像梯度的稀疏性以及非负性作为先验知识进行约束得到稀疏采样角度下图像重建问题的重建模型：argminu||u||TV,s.t.||Au-y||2≤ϵ,u≥0---(1)]]>其中u为待重建图像，被离散化为一个n×n的图像矩阵；||u||TV=Σ|▿us,t|=Σ(us,t-us-1,t)2+(us,t-us,t-1)2+ξ]]>表示图像梯度；矩阵A中的每一个元素Ai,j表示为第i条射线穿过第j个像素的长度；N为待重建图像像素的个数；M为投影数据的个数；(4)将步骤(3)中的稀疏采样角度下图像重建问题的重建模型转化为变分不等式形式：<f(u*),u‑u*&gt;≥0   (2)其中u*,u∈C，u*为变分不等式的解，F(u)＝||u||TV，f(u)＝F'(u)，C1＝{u|||Au‑y||2≤ε}，C2＝{u|u≥0}，C＝C1∩C2；(5)求解步骤(4)中的变分不等式，得到重建图像，即找到一个待重建的图像u*，使得对于任意的u∈C，<f(u*),u‑u*>≥0成立。...

【技术特征摘要】
1.一种稀疏采样角度下基于变分不等式的CT图像重建方法，其特征在于，包括下述步骤：(1)针对扇束CT，在0到180度角度范围内进行等角度间隔投影扫描，投影方向角度间隔的大小在2到6度之间，获取稀疏投影数据y；(2)通过X射线源，探测器以及待重建物体的位置信息计算投影矩阵A；(3)根据步骤(1)中得到的投影数据y以及步骤(2)中得到的投影矩阵A，同时引入图像梯度的稀疏性以及非负性作为先验知识进行约束得到稀疏采样角度下图像重建问题的重建模型：其中u为待重建图像，被离散化为一个n×n的图像矩阵；表示图像梯度；矩阵A中的每一个元素Ai,j表示为第i条射线穿过第j个像素的长度；N为待重建图像像素的个数；M为投影数据的个数；(4)将步骤(3)中的稀疏采样角度下图像重建问题的重建模型转化为变分不等式形式：<f(u*),u-u*>≥0(2)其中u*,u∈C，u*为变分不等式的解，F(u)＝||u||TV，f(u)＝F'(u)，C1＝{u|||Au-y||2≤ε}，C2＝{u|u≥0}，C＝C1∩C2；(5)求解步骤(4)中的变分不等式，得到重建图像，即找到一个待重建的图像u*，使得对于任意的u∈C，<f(u*),u-u*>≥0成立。2.根据权利要求1所述的稀疏采样角度下基于变分不等式的CT图像重建方法，其特征在于，步骤(3)中，||u||TV的梯度是一个图像，该图像的每个像素值都是||u||TV相对应像素的偏微分，并用如下表示：3.根据权利要求1所述的稀疏采样角度下基于变分不等式的CT图像重建方法，其特征在于，在重建的过程中，待重建图像u以及||u||TV的梯度都被转化为列向量，图像中第s行，第t列的像素点对应于列向量中第(s-1)·n+t个像素点。4.根据权利要求1所述的稀疏采样角度下基于变分不等式的CT图像重建方法，其特征在于，步骤(4)中稀疏采样角度下图像重建模型与变分不等式之间的转化，转化过程包括以下步骤：变分不等式的定义：设G:Rn→Rn为连续映射，为非空闭凸集，存在x*∈S，使得对于<G(x*),x-x*>≥0成立；(4-1)令C1＝{u|||Au-y||2≤ε}，C2＝{u|u≥0}，C＝C1∩C2，F(u)＝||u||TV，f(u)＝F'(u)公式(1)变为argminF(u)s.t.u∈C(3)其中F(u)是一个凸函数，因为C1，C2皆为凸集，所以C＝C1∩C2也为凸集(4-2)令u*为公式(3)的最优解，那么对于F(u)≥F(u*)；(4-3)对于任何的γ∈[0,1]，当u，u*属于凸集C时，u*+γ(u-u*)∈C，于是根据步骤(4-2)可得:F(u*+γ(u-u*))-F(u*)≥0(4)因此对于公式(2)成立；这表明，u*也是变分不等式，即公式(2)的解，公式(1)等价于公式(2)。5.根据权利要求1所述的稀疏采样角度下基于变分不等式的CT图像重建方法，其特征在于，步骤(5)所述的变分不等式求解方法，求解过程具体包含以下步骤：(5-1)初始化：u0＝0，d0＝1，ζ＝0.2，设置总的迭代次数NTotal，当前迭代次数k＝0；(5-2)对当前的迭代解uk应用投影算法，得到中间解(5-3)计算下一次迭代中的梯度下降法步长dk+1；(5-4)将列向量转化为二维图像矩阵二维图像矩阵中第s行，第t列的像素点与列向量中第(s-1)·n+t个像素点对应；(...

【专利技术属性】
技术研发人员：高红霞，陈科伟，吴丽璇，胡跃明，
申请(专利权)人：华南理工大学，
类型：发明
国别省市：广东;44