CN108280859A - 一种采样角度受限下的ct稀疏投影图像重建方法及装置 - Google Patents
一种采样角度受限下的ct稀疏投影图像重建方法及装置 Download PDFInfo
- Publication number
- CN108280859A CN108280859A CN201711420601.8A CN201711420601A CN108280859A CN 108280859 A CN108280859 A CN 108280859A CN 201711420601 A CN201711420601 A CN 201711420601A CN 108280859 A CN108280859 A CN 108280859A
- Authority
- CN
- China
- Prior art keywords
- solution
- random
- formula
- optimal solution
- projection
- Prior art date
- Legal status (The legal status is an assumption and is not a legal conclusion. Google has not performed a legal analysis and makes no representation as to the accuracy of the status listed.)
- Granted
Links
Classifications
-
- G—PHYSICS
- G06—COMPUTING; CALCULATING OR COUNTING
- G06T—IMAGE DATA PROCESSING OR GENERATION, IN GENERAL
- G06T11/00—2D [Two Dimensional] image generation
- G06T11/003—Reconstruction from projections, e.g. tomography
- G06T11/006—Inverse problem, transformation from projection-space into object-space, e.g. transform methods, back-projection, algebraic methods
-
- G—PHYSICS
- G06—COMPUTING; CALCULATING OR COUNTING
- G06F—ELECTRIC DIGITAL DATA PROCESSING
- G06F17/00—Digital computing or data processing equipment or methods, specially adapted for specific functions
- G06F17/10—Complex mathematical operations
- G06F17/11—Complex mathematical operations for solving equations, e.g. nonlinear equations, general mathematical optimization problems
-
- G—PHYSICS
- G06—COMPUTING; CALCULATING OR COUNTING
- G06T—IMAGE DATA PROCESSING OR GENERATION, IN GENERAL
- G06T2211/00—Image generation
- G06T2211/40—Computed tomography
- G06T2211/424—Iterative
Landscapes
- Engineering & Computer Science (AREA)
- Physics & Mathematics (AREA)
- General Physics & Mathematics (AREA)
- Mathematical Physics (AREA)
- Theoretical Computer Science (AREA)
- Mathematical Analysis (AREA)
- Mathematical Optimization (AREA)
- Pure & Applied Mathematics (AREA)
- Algebra (AREA)
- Computational Mathematics (AREA)
- Data Mining & Analysis (AREA)
- Operations Research (AREA)
- Databases & Information Systems (AREA)
- Software Systems (AREA)
- General Engineering & Computer Science (AREA)
- Image Processing (AREA)
- Apparatus For Radiation Diagnosis (AREA)
Abstract
本发明公开了一种采样角度受限下的CT稀疏投影图像重建方法及装置。所述方法包括根据投影数据,获得投影方程的解的伪逆矩阵,根据伪逆矩阵,生成随机解集,对当前随机解集中的各个解进行相应的保留或替换,迭代次数达预设最大值时从当前随机解集中选择出最优解作为所需获得的重建结果等步骤。所述装置包括用于存储至少一个程序的存储器和用于加载至少一个程序以执行本发明方法的处理器。本发明通过求离散化投影重建方程的伪逆作为算法的初始值,保证了初始解的质量,通过随机游走生成一组解并分别进行迭代优化,保证了优化路径的多样性,从而克服了传统重建方法由初始解和和迭代路径的缺陷带来的困难。本发明应用于图像处理技术领域。
Description
技术领域
本发明涉及图像处理技术领域,尤其是一种采样角度受限下的CT稀疏投影图像重建方法及装置。
背景技术
术语解释:
CT:Computed Tomography(计算机断层成像)
TV:Total Variation(全变分)
POCS-TV:Projection On Convex Sets-Total Variation(凸集投影-全变分最小化算法)
CT成像是通过利用X射线穿透待检测物体时发生的强度衰减信息,重建得到物体的断层图像。而提高角度受限下稀疏投影数据的重建质量是CT重建领域十分重要且亟需解决的难题。如医学中,在对人体大脑进行扫描时,过大的射线照射剂量可能会对大脑造成二次伤害,因此需要减小CT射线的照射剂量,只能得到角度受限的稀疏投影数据;在工业中,为了对大型的不规则印制电路板或在役的管道进行无损探伤,受CT成像系统几何位置限制,只能在有限角度下对待重建物体进行扫描,无法得到完整角度投影数据。在采样角度受限情况下,因缺失连续角度范围的投影数据,不满足解析法对获得精确重建结果的数据完整性要求(Tuy-Smith条件),直接采用解析法重建存在严重伪影。针对角度受限投影数据的重建,现有技术通常在TV约束的正则化框架下利用基于凸优化的迭代法进行求解。现有技术在迭代重建法的原目标函数(其中A为投影系数矩阵,x为待重建图像的列向量形式,b为投影数据)基础上,根据CT图像的分片光滑特性,添加待重建图像的TV项作为正则化约束,得到目标函数然后通过基于凸优化的迭代法求解得到重建结果。
但是,由于采样角度受限下的稀疏投影数据对称性、完备性极低,重建算法的搜索空间巨大,传统的基于凸优化的迭代重建方法从单一的初始解出发,按照目标函数在当前解的负梯度方向路径进行搜索,一方面重建质量受初始解质量的严重影响,另一方面,单一的固定迭代路径使得优化过程无法在有限时间内遍历巨大的搜索空间,易陷入局部最优且无法跳离,重建图像质量差。由此可见,现有的重建方法面临的困难,主要是由初始解和迭代路径引起的。
发明内容
为了解决上述技术问题,本发明的第一目的在于提供一种采样角度受限下的CT稀疏投影图像重建方法,第二目的在于提供一种采样角度受限下的CT稀疏投影图像重建装置。
本发明所采取的第一技术方案是:
一种采样角度受限下的CT稀疏投影图像重建方法,包括:
S1.获取投影方程的解的伪逆矩阵,所述投影方程根据投影数据而建立;
S2.根据伪逆矩阵,生成与本轮迭代相对应的随机解集;
S3.若当前为第一轮迭代过程,则返回步骤S2开始执行,反之,将当前随机解集与上一轮迭代过程的随机解集进行对应解的比较,根据比较结果,对当前随机解集中的各个解进行相应的保留或替换;
S4.判断当前迭代次数是否已达预设最大值,若是,则根据适应度评价,从当前随机解集中选择出最优解,然后执行步骤S5,反之,返回从步骤S2开始执行;
S5.根据最优解重建得到CT图像。
进一步地,所述步骤S2具体为:根据伪逆矩阵,通过以下公式获得多个随机解,将获得的所有解构成与该轮迭代相对应的随机解集:
式中,t为迭代轮次,i表示序号,N表示自然数集,为第t轮迭代所得的随机解集中第i个解,xini为投影方程的解的伪逆矩阵,ε是一维度与xini相同的服从(0,σ)分布的随机变量。
进一步地,步骤S3具体包括:
S31.对当前随机解集中的每一个解进行有效性约束;
S32.计算当前随机解集中每一个解的后验概率密度;
S33.根据当前随机解集中每一个解的后验概率密度与上一轮迭代过程的随机解集中每一个解的后验概率密度,计算梅特罗波利斯比率;
S34.判断梅特罗波利斯比率是否大于第一阈值,若是,则保留当前随机解集中的解,反之,将当前随机解集中的解替换成上一轮迭代过程的随机解集中的对应解。
进一步地,步骤S31具体为:利用下式,对当前随机解集中的每一个解进行有效性约束:
式中,为的第j个元素;
步骤S32具体为:利用下式,计算当前随机解集中每一个解的后验概率密度:
式中,为的后验概率密度,C为的协方差矩阵;
步骤S33具体为:利用下式,计算梅特罗波利斯比率:
式中,α为梅特罗波利斯比率;
步骤S34具体为:利用下式,对梅特罗波利斯比率是否大于第一阈值进行判断,从而保留当前随机解集中的解,或者将当前随机解集中的解替换成上一轮迭代过程的随机解集中的对应解:
式中,rand()为一值域为[0,1]的随机函数,rand()的值为第一阈值。
进一步地,步骤S4中,选择最优解所用的公式为:
式中,为所要获取的最优解,F(·)为适应度评价函数,A为投影系数矩阵,b为投影数据,为的TV正则化约束项。
进一步地,步骤S4与步骤S5之间还包括以下步骤:
S4A.利用基于凸优化的迭代重建算法,对最优解进行更新。
进一步地,步骤S4A具体包括:
S4A1.采用ART迭代算法对最优解进行处理,从而得到一致最优解;
S4A2.对一致最优解进行图像非负性约束,从而得到非负最优解;
S4A3.根据非负最优解的TV项的梯度,以非负最优解作为迭代变量进行迭代,将得到的迭代结果作为更新后的最优解。
进一步地,所述步骤S4A1具体为:使用以下公式对最优解进行处理,从而得到一致最优解:
式中,M为投影数据总数,bm为第m个投影数据,Am为投影系数矩阵A的第m行,根据公式求得的fM为一致最优解;
所述步骤S4A2具体为:使用以下公式对一致最优解进行图像非负性约束,从而得到非负最优解:
式中,g0(j)为g0的第j个元素,根据公式求得的g0为非负最优解。
进一步地,所述步骤S4A3具体为:使用以下公式进行迭代,将得到的迭代结果作为更新后的最优解:
式中,Ngrad为梯度下降迭代总次数,
表示gn中坐标为(p,q)的点的像素值,
根据公式计算求得的为更新后的最优解。
本发明所采取的第二技术方案是:
一种采样角度受限下的CT稀疏投影图像重建装置,包括:
存储器,用于存储至少一个程序;
处理器,用于加载所述至少一个程序以执行第一技术方案所述一种采样角度受限下的CT稀疏投影图像重建方法。
本发明的有益效果是:本发明通过求离散化投影重建方程的伪逆作为算法的初始值,保证了初始解的质量,通过随机游走生成一组解并分别进行迭代优化,保证了优化路径的多样性,从而克服了传统重建方法由初始解和迭代路径的缺陷带来的困难。
进一步地,本发明还利用高斯马尔科夫模型对解的进化能力进行预测,基于解的进化潜能对其进行随机更新,使得算法能够在提高全局搜索能力的同时保证效率;对基于适应度评价下的最优解使用基于凸优化的迭代更新,并将更新后的最优解作为最终的重建结果,进一步提高了包含随机因素下算法的效率。本发明全局搜索能力强,重建图像不存在严重伪影,显著提高了重建图像的质量,降低重建误差,能够在数据缺失严重的情况下实现图像重建,特别在处理含噪数据时表现出很强的鲁棒性和适应性。
附图说明
图1为本发明方法的流程图;
图2为Shepp-Logan模型原图;
图3为投影角度间隔为3°时应用POCS-TV算法的重建效果图;
图4为投影角度间隔为3°时应用本发明方法的重建效果图;
图5为投影角度间隔为3°时重建结果第64行的像素值曲线图;
图6为投影角度间隔为3°时重建结果第64列的像素值曲线图;
图7为投影角度间隔为6°时应用POCS-TV算法的重建效果图;
图8为投影角度间隔为6°时应用本发明方法的重建效果图;
图9为投影角度间隔为6°时重建结果第64行的像素值曲线图;
图10为投影角度间隔为6°时重建结果第64列的像素值曲线图;
图11为本发明装置的结构框图。
具体实施方式
实施例1
一种采样角度受限下的CT稀疏投影图像重建方法,如图1所示,包括:
S1.获取投影方程的解的伪逆矩阵,所述投影方程根据投影数据而建立;
S2.根据伪逆矩阵,生成与本轮迭代相对应的随机解集;
S3.若当前为第一轮迭代过程,则返回步骤S2开始执行,反之,将当前随机解集与上一轮迭代过程的随机解集进行对应解的比较,根据比较结果,对当前随机解集中的各个解进行相应的保留或替换;
S4.判断当前迭代次数是否已达预设最大值,若是,则根据适应度评价,从当前随机解集中选择出最优解,然后执行步骤S5,反之,返回从步骤S2开始执行;
S5.根据最优解重建得到CT图像。
在应用本发明方法之前,可以使用各种方法采集角度受限下的稀疏投影数据。所述角度受限下的稀疏投影数据,可以是投影角度范围为(0°,90°),投影间隔大于等于3°、小于等于6°所得的投影数据。根据投影数据重建图像,实质是要求解离散化的投影方程Ax=b,其中A为投影系数矩阵,x为待重建图像的列向量形式,b为投影数据。求得方程Ax=b的解x,便能根据x重建图像。
x的伪逆矩阵可以使用以下公式xini=(ATA)-1ATb求得。由于在本发明方法中,x的伪逆矩阵用于作为后续一系列迭代运算的初值,因此下标ini;本发明方法中只使用到了x的伪逆矩阵,因此本发明方法中“伪逆矩阵”专指x的伪逆矩阵。
根据伪逆矩阵,可以使用随机游走等算法产生一系列矩阵,由于这些矩阵都有可能是方程Ax=b的解或者解的近似,对这些矩阵使用迭代算法,将收敛得到方程Ax=b的解,因此将一系列矩阵的集合作为随机解集。
在以下方法中,需要使用到迭代算法或者循环的处理步骤,在一轮迭代过程当中,可能会对随机解集或者随机解集当中特定的解做处理,但在此轮迭代过程当中,处理前和处理后的随机解集都称作当前随机解集。
在执行第一轮迭代的情况下,由于不存在上一轮迭代,因此步骤S3将无法直接执行。可以在执行S1和S2从而得到一个随机解集之后,再次执行S2,在执行后一次S2时便进入第二轮迭代,此时便可以执行步骤S3。预设的迭代次数最大值可以为1,这样,可以在执行S1和S2后,跳过S3,然后执行S4中根据适应度评价,从当前随机解集中选择出最优解的步骤,然后跳到步骤S5执行。
步骤S3中,将本轮迭代的随机解集,也就是当前随机解集,与上一轮迭代过程的随机解集一起针对各对对应解分别进行比较。对应解指在两个解集中索引号、序号、编号相同或者地位相同的两个解。经过各个对应解的比较,按照进化能力的标准,在当前随机解集中的解保留下来,或者将当前随机解集中的部分或者全部解替换成上一轮迭代过程的随机解集中的对应解,所得的随机解集仍称作当前随机解集。上述操作使得当前随机解集中的各个解具有更佳的进化能力。
步骤S4中判断当前迭代次数是否已达预设最大值,如果没有,那么就返回步骤S2重新执行,使得当前随机解集不断完善。如果已达预设最大值,那么就使用适应度评价的方法,根据对当前随机解集的适应度评价结果,在当前随机解集中选择出一个最优解作为方程Ax=b的解,也就是所需获得的重建结果。
进一步作为优选的实施方式,所述步骤S2具体为:根据伪逆矩阵,通过以下公式获得多个随机解,将获得的所有解构成与该轮迭代相对应的随机解集:
式中,t为迭代轮次,i为解的编号,为第t轮迭代所得的随机解集中第i个解,xini为投影方程的解的伪逆矩阵,ε是一维度与xini相同的服从(0,σ)分布的随机变量。
由于ε是随机变量,每一次执行S2,ε的取值都是随机的,其与xini叠加后所得的结果也是随机的,由此每一次执行S2都可以得到一个新的随机解集。可以根据实际计算的需要,确定随机解集中解的总个数。
下面使用高斯马尔科夫模型对解的进化能力进行定量评价,然后基于解的进化能力对当前随机解集进行更新。
进一步作为优选的实施方式,步骤S3具体为:
S31.对当前随机解集中的每一个解进行有效性约束;
S32.计算当前随机解集中每一个解的后验概率密度;
S33.根据当前随机解集中每一个解的后验概率密度与上一轮迭代过程的随机解集中每一个解的后验概率密度,计算梅特罗波利斯比率;
S34.判断梅特罗波利斯比率是否大于第一阈值,若是,则保留当前随机解集中的解,反之,将当前随机解集中的解替换成上一轮迭代过程的随机解集中的对应解。
步骤S31-S33可以计算出梅特罗波利斯(Metropolis)比率,梅特罗波利斯比率可以反映出当前随机解集与上一轮迭代过程的随机解集的进化能力孰大孰小。根据实际需要,确定第一阈值,即梅特罗波利斯比率的对比值,如果梅特罗波利斯比率小于或等于第一阈值,表明当前随机解集的进化能力不如上一轮迭代过程的随机解集,那么就将当前随机解集中的部分解或者全部解替换成上一轮迭代过程的随机解集的对应解;如果梅特罗波利斯比率大于第一阈值,表明当前随机解集的进化能力比上一轮迭代过程的随机解集好,可以将当前随机解集中的部分解或全部解保留下来。经过这样处理后的当前随机解集仍称作当前随机解集。反复执行步骤S34,直至对当前随机解集中的所有解进行保留或替换。
进一步作为优选的实施方式,步骤S31具体为:利用下式,对当前随机解集中的每一个解进行有效性约束:
式中,为的第j个元素。
由公式可知,有效性约束可以将当前随机解集中的各个解的各元素都约束在(0,1)的范围内。即使当前随机解集中各个解是经过随机过程得到的,通过有效性约束,也可以保证各个解中每个元素的值对应于图像像素值都是有意义的。
步骤S32具体为:利用下式,计算当前随机解集中每一个解的后验概率密度:
式中,为的后验概率密度,C为的协方差矩阵。
步骤S33具体为:利用下式,计算梅特罗波利斯比率:
式中,α为梅特罗波利斯比率。
步骤S34具体为:利用下式,对梅特罗波利斯比率是否大于第一阈值进行判断,从而保留当前随机解集中的解,或者将当前随机解集中的解替换成上一轮迭代过程的随机解集中的对应解:
式中,rand()为一值域为[0,1]的随机函数,rand()的值为第一阈值。
在步骤S34中,可以分别对每对对应解进行判断,而且第一阈值不是一个固定不变的值,由于其是随机函数,其对于不同的对应解可以取不同的值。由于梅特罗波利斯比率进行对比的阈值不同,也就是对每对对应解的进化能力的评价标准是不同的,于是可能出现当前随机解集中部分解被保留、部分解被替换的情况。而且,由于阈值的随机性,还可能出现当小于也就是的进化能力不及时,仍有可能被保留下来的情况,这样,本发明方法可以对更广阔的解空间进行探索,提高全局搜索能力,可以缓解陷入局部最优的情况。
进一步作为优选的实施方式,步骤S4中,选择最优解所用的公式为:
式中,为所要获取的最优解,F(·)为适应度评价函数,A为投影系数矩阵,b为投影数据,表示的TV正则化约束项。
适应度评价函数中,为数据保真项,评价重建图像与投影数据的一致性,为TV正则化约束项,评价重建图像的分片光滑性,适应度评价函数值越小,表明解的质量越好。
在经过预定轮次的迭代后,便可以从当前随机解集中选择出最优解作为所需获得的重建结果。优选地,选择最优解的标准是,使得适应度评价函数取得最小值时的xi t的值。
进一步作为优选的实施方式,步骤S4与步骤S5之间还包括以下步骤:
S4A.利用基于凸优化的迭代重建算法,对最优解进行更新。
步骤S4中得到的最优解,可以直接执行步骤S5,也就是直接用于作为所需获得的重建结果。为了取得更好的效果,可以在执行步骤S4得到最优解后,接着执行步骤S4A,对最优解进行更新,然后执行步骤S5,根据经过步骤S4A处理的最优解,也就是更新后的最优解,来重建得到CT图像。
进一步作为优选的实施方式,步骤S4A具体包括:
S4A1.采用ART迭代算法对最优解进行处理,从而得到一致最优解;
S4A2.对一致最优解进行图像非负性约束,从而得到非负最优解;
S4A3.根据非负最优解的TV项的梯度,以非负最优解作为迭代变量进行迭代,将得到的迭代结果作为更新后的最优解。
进一步作为优选的实施方式,所述步骤S4A1具体为:使用以下公式对最优解进行处理,从而得到一致最优解:
式中,M为投影数据总数,bm为第m个投影数据,Am为投影系数矩阵A的第m行,根据公式求得的fM为一致最优解。
上述公式是ART迭代算法的一种实现方法,最终得到的结果fM为一致最优解。
进一步作为优选的实施方式,所述步骤S4A2具体为:使用以下公式对一致最优解进行图像非负性约束,从而得到非负最优解:
式中,g0(j)为g0的第j个元素,根据公式求得的g0为非负最优解。
进一步作为优选的实施方式,所述步骤S4A3具体为:使用以下公式进行迭代,将得到的迭代结果作为更新后的最优解:
式中,Ngrad为梯度下降迭代总次数,
表示gn中坐标为(p,q)的点的像素值,根据公式计算求得的为更新后的最优解。经过步骤S4A对最优解使用基于凸优化的更新处理,可以保证包含随机因素下算法的效率。
实施例2
在本实施例中,使用实施例1的方法,对二维Shepp-Logan模型(128×128)的模拟投影数据进行重建,并在相同的迭代次数下与POCS-TV算法进行比较。假设模拟投影数据获取及重建过程的参数如下:
(1)X射线源与待重建物体中心的距离为256mm;
(2)X射线检测器与待重建物体中心的距离为256mm;
(3)X射线检测器单元的个数为256,每单元的宽度为0.5mm;
(4)投影角度范围为(0°,90°),投影角度间隔为3°和6°;
(5)为说明本发明在含噪数据的重建上具有更强的鲁棒性和适应性,在投影数据上添加0均值,0.15方差的高斯噪声;
(6)设置总的迭代次数为5000次。
Shepp-Logan模型原图如图2所示。
投影角度间隔为3°时,用POCS-TV算法的重建结果如图3所示,用本发明方法的重建结果如图4所示,重建结果第64行的像素值曲线如图5所示,第64列的像素值曲线如图6所示。
投影角度间隔为6°时,用POCS-TV算法的重建结果如图7所示,用本发明方法的重建结果如图8所示,重建结果第64行的像素值曲线如图9所示,第64列的像素值曲线如图10所示。
从图3至图10可以看出,本发明方法的重建结果轮廓完整,细节更加清晰,更接近于模型原图,重建质量远胜于POCS-TV算法。
实施例3
一种采样角度受限下的CT稀疏投影图像重建装置,如图11所示,包括:
存储器,用于存储至少一个程序;
处理器,用于加载所述至少一个程序以执行实施例1和实施例2所述一种采样角度受限下的CT稀疏投影图像重建方法。
以上是对本发明的较佳实施进行了具体说明,但对本发明创造并不限于所述实施例,熟悉本领域的技术人员在不违背本发明精神的前提下还可做作出种种的等同变形或替换,这些等同的变形或替换均包含在本申请权利要求所限定的范围内。
Claims (10)
1.一种采样角度受限下的CT稀疏投影图像重建方法,其特征在于,包括:
S1.获取投影方程的解的伪逆矩阵,所述投影方程根据投影数据而建立;
S2.根据伪逆矩阵,生成与本轮迭代相对应的随机解集;
S3.若当前为第一轮迭代过程,则返回步骤S2开始执行,反之,将当前随机解集与上一轮迭代过程的随机解集进行对应解的比较,根据比较结果,对当前随机解集中的各个解进行相应的保留或替换;
S4.判断当前迭代次数是否已达预设最大值,若是,则根据适应度评价,从当前随机解集中选择出最优解,然后执行步骤S5,反之,返回从步骤S2开始执行;
S5.根据最优解重建得到CT图像。
2.根据权利要求1所述的一种采样角度受限下的CT稀疏投影图像重建方法,其特征在于,所述步骤S2具体为:根据伪逆矩阵,通过以下公式获得多个随机解,将获得的所有解构成与该轮迭代相对应的随机解集:
式中,t为迭代轮次,i表示序号,N表示自然数集,为第t轮迭代所得的随机解集中第i个解,xini为投影方程的解的伪逆矩阵,ε是一维度与xini相同的服从(0,σ)分布的随机变量。
3.根据权利要求2所述的一种采样角度受限下的CT稀疏投影图像重建方法,其特征在于,步骤S3具体包括:
S31.对当前随机解集中的每一个解进行有效性约束;
S32.计算当前随机解集中每一个解的后验概率密度;
S33.根据当前随机解集中每一个解的后验概率密度与上一轮迭代过程的随机解集中每一个解的后验概率密度,计算梅特罗波利斯比率;
S34.判断梅特罗波利斯比率是否大于第一阈值,若是,则保留当前随机解集中的解,反之,将当前随机解集中的解替换成上一轮迭代过程的随机解集中的对应解。
4.根据权利要求3所述的一种采样角度受限下的CT稀疏投影图像重建方法,其特征在于:
步骤S31具体为:利用下式,对当前随机解集中的每一个解进行有效性约束:
式中,为的第j个元素;
步骤S32具体为:利用下式,计算当前随机解集中每一个解的后验概率密度:
式中,为的后验概率密度,C为的协方差矩阵;
步骤S33具体为:利用下式,计算梅特罗波利斯比率:
式中,α为梅特罗波利斯比率;
步骤S34具体为:利用下式,对梅特罗波利斯比率是否大于第一阈值进行判断,从而保留当前随机解集中的解,或者将当前随机解集中的解替换成上一轮迭代过程的随机解集中的对应解:
式中,rand()为一值域为[0,1]的随机函数,rand()的值为第一阈值。
5.根据权利要求1所述的一种采样角度受限下的CT稀疏投影图像重建方法,其特征在于,步骤S4中,选择最优解所用的公式为:
式中,为所要获取的最优解,F(·)为适应度评价函数,A为投影系数矩阵,b为投影数据,为的TV正则化约束项。
6.根据权利要求1-5任一项所述的一种采样角度受限下的CT稀疏投影图像重建方法,其特征在于,步骤S4与步骤S5之间还包括以下步骤:
S4A.利用基于凸优化的迭代重建算法,对最优解进行更新。
7.根据权利要求6所述的一种采样角度受限下的CT稀疏投影图像重建方法,其特征在于,步骤S4A具体包括:
S4A1.采用ART迭代算法对最优解进行处理,从而得到一致最优解;
S4A2.对一致最优解进行图像非负性约束,从而得到非负最优解;
S4A3.根据非负最优解的TV项的梯度,以非负最优解作为迭代变量进行迭代,将得到的迭代结果作为更新后的最优解。
8.根据权利要求7所述的一种采样角度受限下的CT稀疏投影图像重建方法,其特征在于:
所述步骤S4A1具体为:使用以下公式对最优解进行处理,从而得到一致最优解:
式中,M为投影数据总数,bm为第m个投影数据,Am为投影系数矩阵A的第m行,根据公式求得的fM为一致最优解;
所述步骤S4A2具体为:使用以下公式对一致最优解进行图像非负性约束,从而得到非负最优解:
式中,g0(j)为g0的第j个元素,根据公式求得的g0为非负最优解。
9.根据权利要求8所述的一种采样角度受限下的CT稀疏投影图像重建方法,其特征在于,所述步骤S4A3具体为:使用以下公式进行迭代,将得到的迭代结果作为更新后的最优解:
式中,Ngrad为梯度下降迭代总次数,
表示gn中坐标为(p,q)的点的像素值,
根据公式计算求得的为更新后的最优解。
10.一种采样角度受限下的CT稀疏投影图像重建装置,其特征在于,包括:
存储器,用于存储至少一个程序;
处理器,用于加载所述至少一个程序以执行权利要求1-9任一项所述一种采样角度受限下的CT稀疏投影图像重建方法。
Priority Applications (1)
Application Number | Priority Date | Filing Date | Title |
---|---|---|---|
CN201711420601.8A CN108280859B (zh) | 2017-12-25 | 2017-12-25 | 一种采样角度受限下的ct稀疏投影图像重建方法及装置 |
Applications Claiming Priority (1)
Application Number | Priority Date | Filing Date | Title |
---|---|---|---|
CN201711420601.8A CN108280859B (zh) | 2017-12-25 | 2017-12-25 | 一种采样角度受限下的ct稀疏投影图像重建方法及装置 |
Publications (2)
Publication Number | Publication Date |
---|---|
CN108280859A true CN108280859A (zh) | 2018-07-13 |
CN108280859B CN108280859B (zh) | 2021-03-30 |
Family
ID=62802211
Family Applications (1)
Application Number | Title | Priority Date | Filing Date |
---|---|---|---|
CN201711420601.8A Active CN108280859B (zh) | 2017-12-25 | 2017-12-25 | 一种采样角度受限下的ct稀疏投影图像重建方法及装置 |
Country Status (1)
Country | Link |
---|---|
CN (1) | CN108280859B (zh) |
Cited By (6)
Publication number | Priority date | Publication date | Assignee | Title |
---|---|---|---|---|
CN109903353A (zh) * | 2019-01-28 | 2019-06-18 | 华南理工大学 | 一种迭代演化模型的ct图像稀疏重建方法 |
CN110210443A (zh) * | 2019-06-11 | 2019-09-06 | 西北工业大学 | 一种优化投影对称性近似稀疏分类的手势识别方法 |
WO2020151424A1 (zh) * | 2019-01-24 | 2020-07-30 | 浙江大学 | 一种基于各向异性全变分的有限角度ct重建算法 |
CN111956180A (zh) * | 2019-05-20 | 2020-11-20 | 华北电力大学(保定) | 一种重建光声内窥层析图像的方法 |
CN112288832A (zh) * | 2020-12-24 | 2021-01-29 | 中国人民解放军国防科技大学 | 一种角度受限及稀疏采样的层析成像图像重构方法 |
CN115619890A (zh) * | 2022-12-05 | 2023-01-17 | 山东省计算中心(国家超级计算济南中心) | 基于并行随机迭代求解线性方程组的断层成像方法及系统 |
Citations (7)
Publication number | Priority date | Publication date | Assignee | Title |
---|---|---|---|---|
CN103824295A (zh) * | 2014-03-03 | 2014-05-28 | 天津医科大学 | 一种肺部ct图像中粘连血管型肺结节的分割方法 |
CN104091355A (zh) * | 2014-06-06 | 2014-10-08 | 华南理工大学 | 一种采样角度受限下的ct图像重建方法 |
CN104103086A (zh) * | 2014-06-06 | 2014-10-15 | 华南理工大学 | 一种稀疏采样角度下基于变分不等式的ct图像重建方法 |
CN105118078A (zh) * | 2015-09-24 | 2015-12-02 | 中国科学院苏州生物医学工程技术研究所 | 欠采样的ct图像重建方法 |
US20160242721A1 (en) * | 2015-02-20 | 2016-08-25 | Kabushiki Kaisha Toshiba | Apparatus and method for fast iterative reconstruction in computed tomography |
CN106163403A (zh) * | 2013-12-18 | 2016-11-23 | 伊利克塔股份有限公司 | Ct图像中的目标特定剂量和散射估计 |
CN107122725A (zh) * | 2017-04-18 | 2017-09-01 | 深圳大学 | 一种基于联合稀疏判别分析的人脸识别方法及其系统 |
-
2017
- 2017-12-25 CN CN201711420601.8A patent/CN108280859B/zh active Active
Patent Citations (7)
Publication number | Priority date | Publication date | Assignee | Title |
---|---|---|---|---|
CN106163403A (zh) * | 2013-12-18 | 2016-11-23 | 伊利克塔股份有限公司 | Ct图像中的目标特定剂量和散射估计 |
CN103824295A (zh) * | 2014-03-03 | 2014-05-28 | 天津医科大学 | 一种肺部ct图像中粘连血管型肺结节的分割方法 |
CN104091355A (zh) * | 2014-06-06 | 2014-10-08 | 华南理工大学 | 一种采样角度受限下的ct图像重建方法 |
CN104103086A (zh) * | 2014-06-06 | 2014-10-15 | 华南理工大学 | 一种稀疏采样角度下基于变分不等式的ct图像重建方法 |
US20160242721A1 (en) * | 2015-02-20 | 2016-08-25 | Kabushiki Kaisha Toshiba | Apparatus and method for fast iterative reconstruction in computed tomography |
CN105118078A (zh) * | 2015-09-24 | 2015-12-02 | 中国科学院苏州生物医学工程技术研究所 | 欠采样的ct图像重建方法 |
CN107122725A (zh) * | 2017-04-18 | 2017-09-01 | 深圳大学 | 一种基于联合稀疏判别分析的人脸识别方法及其系统 |
Non-Patent Citations (2)
Title |
---|
王林元,等: "基于稀疏优化的计算机断层成像图像不完全角度重建综述", 《物理学报》 * |
陈海,等: "基于缺陷投影的稀疏角度数据CT图像迭代重建仿真研究", 《医疗卫生装备》 * |
Cited By (10)
Publication number | Priority date | Publication date | Assignee | Title |
---|---|---|---|---|
WO2020151424A1 (zh) * | 2019-01-24 | 2020-07-30 | 浙江大学 | 一种基于各向异性全变分的有限角度ct重建算法 |
US11727609B2 (en) | 2019-01-24 | 2023-08-15 | Zhejiang University | Limited-angle CT reconstruction method based on anisotropic total variation |
CN109903353A (zh) * | 2019-01-28 | 2019-06-18 | 华南理工大学 | 一种迭代演化模型的ct图像稀疏重建方法 |
CN109903353B (zh) * | 2019-01-28 | 2023-02-14 | 华南理工大学 | 一种迭代演化模型的ct图像稀疏重建方法 |
CN111956180A (zh) * | 2019-05-20 | 2020-11-20 | 华北电力大学(保定) | 一种重建光声内窥层析图像的方法 |
CN110210443A (zh) * | 2019-06-11 | 2019-09-06 | 西北工业大学 | 一种优化投影对称性近似稀疏分类的手势识别方法 |
CN110210443B (zh) * | 2019-06-11 | 2022-03-15 | 西北工业大学 | 一种优化投影对称性近似稀疏分类的手势识别方法 |
CN112288832A (zh) * | 2020-12-24 | 2021-01-29 | 中国人民解放军国防科技大学 | 一种角度受限及稀疏采样的层析成像图像重构方法 |
CN115619890A (zh) * | 2022-12-05 | 2023-01-17 | 山东省计算中心(国家超级计算济南中心) | 基于并行随机迭代求解线性方程组的断层成像方法及系统 |
CN115619890B (zh) * | 2022-12-05 | 2023-04-07 | 山东省计算中心(国家超级计算济南中心) | 基于并行随机迭代求解线性方程组的断层成像方法及系统 |
Also Published As
Publication number | Publication date |
---|---|
CN108280859B (zh) | 2021-03-30 |
Similar Documents
Publication | Publication Date | Title |
---|---|---|
CN108280859A (zh) | 一种采样角度受限下的ct稀疏投影图像重建方法及装置 | |
US11610346B2 (en) | Image reconstruction using machine learning regularizers | |
Gregor et al. | Computational analysis and improvement of SIRT | |
Liu et al. | Adaptive-weighted total variation minimization for sparse data toward low-dose x-ray computed tomography image reconstruction | |
Tian et al. | Low-dose CT reconstruction via edge-preserving total variation regularization | |
EP2881039A1 (en) | X-ray computer tomography image pick-up device and image reconstruction method | |
Li et al. | DECT-MULTRA: Dual-energy CT image decomposition with learned mixed material models and efficient clustering | |
CN110651302B (zh) | 用于图像重建的方法和设备 | |
Hamelin et al. | Design of iterative ROI transmission tomography reconstruction procedures and image quality analysis | |
CN110335325A (zh) | 一种ct图像重建方法及其系统 | |
Hu et al. | An extended simultaneous algebraic reconstruction technique (E‐SART) for X‐ray dual spectral computed tomography | |
Majee et al. | Multi-slice fusion for sparse-view and limited-angle 4D CT reconstruction | |
Yu et al. | Weighted adaptive non-local dictionary for low-dose CT reconstruction | |
Abascal et al. | A novel prior-and motion-based compressed sensing method for small-animal respiratory gated CT | |
Wang et al. | Sparse-view cone-beam CT reconstruction by bar-by-bar neural FDK algorithm | |
CN110298900A (zh) | 一种基于各向异性基函数建立spect重构的方法 | |
Yan et al. | EM+ TV based reconstruction for cone-beam CT with reduced radiation | |
Qiu et al. | Evaluating iterative algebraic algorithms in terms of convergence and image quality for cone beam CT | |
Iborra et al. | Noise analysis in computed tomography (CT) image reconstruction using QR-Decomposition algorithm | |
CN116188615A (zh) | 一种基于正弦域和图像域的稀疏角度ct重建方法 | |
Gong et al. | Structure-guided computed tomography reconstruction from limited-angle projections | |
Saha et al. | Multi-axial CT reconstruction from few view projections | |
Gopi et al. | Iterative computed tomography reconstruction from sparse-view data | |
Flores et al. | Application of a modified LSQR method for CT imaging reconstruction with low doses to patient | |
US20220172412A1 (en) | Medical image reconstruction method, computer device and storage medium |
Legal Events
Date | Code | Title | Description |
---|---|---|---|
PB01 | Publication | ||
PB01 | Publication | ||
SE01 | Entry into force of request for substantive examination | ||
SE01 | Entry into force of request for substantive examination | ||
GR01 | Patent grant | ||
GR01 | Patent grant |