CN109471053A - 一种基于双约束的介电特性迭代成像方法 - Google Patents
一种基于双约束的介电特性迭代成像方法 Download PDFInfo
- Publication number
- CN109471053A CN109471053A CN201811215000.8A CN201811215000A CN109471053A CN 109471053 A CN109471053 A CN 109471053A CN 201811215000 A CN201811215000 A CN 201811215000A CN 109471053 A CN109471053 A CN 109471053A
- Authority
- CN
- China
- Prior art keywords
- follows
- double
- dielectric property
- equation
- reaction equation
- 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
- G01—MEASURING; TESTING
- G01R—MEASURING ELECTRIC VARIABLES; MEASURING MAGNETIC VARIABLES
- G01R33/00—Arrangements or instruments for measuring magnetic variables
- G01R33/20—Arrangements or instruments for measuring magnetic variables involving magnetic resonance
- G01R33/44—Arrangements or instruments for measuring magnetic variables involving magnetic resonance using nuclear magnetic resonance [NMR]
- G01R33/48—NMR imaging systems
-
- G—PHYSICS
- G01—MEASURING; TESTING
- G01R—MEASURING ELECTRIC VARIABLES; MEASURING MAGNETIC VARIABLES
- G01R33/00—Arrangements or instruments for measuring magnetic variables
- G01R33/20—Arrangements or instruments for measuring magnetic variables involving magnetic resonance
- G01R33/44—Arrangements or instruments for measuring magnetic variables involving magnetic resonance using nuclear magnetic resonance [NMR]
- G01R33/48—NMR imaging systems
- G01R33/54—Signal processing systems, e.g. using pulse sequences ; Generation or control of pulse sequences; Operator console
Landscapes
- Physics & Mathematics (AREA)
- High Energy & Nuclear Physics (AREA)
- Condensed Matter Physics & Semiconductors (AREA)
- General Physics & Mathematics (AREA)
- Engineering & Computer Science (AREA)
- Signal Processing (AREA)
- Magnetic Resonance Imaging Apparatus (AREA)
Abstract
本发明公开了一种基于双约束的介电特性迭代成像方法,属于磁共振断层成像领域,首先通过磁共振成像设备获得成像所需要的MR收发相位数据;然后对MR收发相位数据进行高斯滤波处理,再将得到的收发相位代入对流反应方程中,利用有限差分方法将对流反应方程转换成线性系统方程;再将线性系统方程采用最小二乘方式构造全变分和小波的双约束优化问题,再通过Split Bregman方法对该优化问题进行迭代求解;最后将求解得到的结果取倒数从而得到介电特性中电导率分布的重建结果,本发明解决了传统介电特性成像方法存在误差较大,噪声影响严重以及临床适用性差的问题。
Description
技术领域
本发明属于磁共振断层成像领域,涉及基于磁共振相位信息获得组织电导率分布的成像方法,特别涉及一种基于双约束的介电特性迭代成像方法。
背景技术
介电特性(Electrical properties,Eps)主要包括电导率σ和电容率ε,它们是介质在电磁场环境下表达与电磁场相互作用的基本特性。在生物组织中,介电特性的值是与频率相关的,电导率的分布是计算特定吸收率(Specific absorption rate,SAR)的主要参数之一,SAR的分布是高场磁共振研究中的关键性问题之一,所以准确的估计出电导率的分布有助于高场或者超高场磁共振的发展对介电特性进行成像有着很深远且广阔的前景。
近年来,基于磁共振的介电特性断层成像(Magnetic resonance electricalproperties tomography,MR-EPT)方法受到研究者的注意。MR-EPT是一种无损的计算介电特性分布的成像方法,最早由Haccke等人提出,它是以亥姆霍兹方程为基础,通过MR系统采集得到的射频磁场的信息,进而求解出组织介电特性的分布。此外Wen等人注意到电导率的分布主要是受到射频磁场相位信息的影响,因此可以利用相位信息求解电导率的分布。
现有主要的MR-EPT的算法是基于组织介电特性分布均匀的假设前提下推导得出,相对应的方程为其中σ是电导率,为发射相位,μ0是真空磁导率,ω是角频率。这种算法存在很明显的缺点,首先均匀性的假设在很多的实际情况中并不满足,例如人体的大脑,这会使得重建的结果产生很大的误差;其次是收发相位假设,即根据互易定理可知发射线圈的相位和接收线圈的相位在使用正交线圈作为发射和接收线圈的情况下相等,因而采用扫描得到的收发相位的一半作为发射相位,但是这个假设只有在正交线圈同时作为发射和接收线圈的情况下才成立,因而对于医院临床中常用的阵列线圈来说并不适用;最后此算法还存在信噪比很低的问题,主要原因是公式中的Laplacian算子对噪声非常的敏感。
为了解决上述的问题,Balidemaj采用对比源反演算法进行介电常数分布重建,它主要是基于体积分方程实验,它不需要采用任何的假设,但到目前为止还未在实际数据中进行测试;Liu提出了基于梯度的MR-EPT,它将含有介电特性梯度的方程作为核心方程,通过多通道收发线圈采集得到的数据进行计算。这种方法虽然可以提高重建的质量,但是它却需要多通道线圈同时作为发射和接收线圈,这种线圈配置并不适用于现在的临床情况。
发明内容
本发明的目的在于:提供了一种基于双约束的介电特性迭代成像方法,解决了传统介电特性成像方法存在误差较大,噪声影响严重以及临床适用性差的问题。
本发明采用的技术方案如下:
一种基于双约束的介电特性迭代成像方法,包括以下步骤:
步骤1:获取MR收发相位数据,进行高斯滤波处理后得到收发相位
步骤2:将收发相位代入对流反应方程中,所述对流反应方程为:
其中,u为电导率的倒数,为u的梯度,为的梯度,为收发相位的拉普拉斯运算,
再通过有限差分方法对对流反应方程进行转换得到线性系统方程,所述线性系统方程为:
Aρ=b,
其中,将u写成列向量形式即为ρ,b是一个列向量,每一点的值都为2ωμ0;
步骤3:将线性系统方程采用最小二乘方式构造双约束优化问题,所述双约束优化问题的公式为:
其中,为全变分(total variation,TV),W是小波变换矩阵,λ1,λ2分别为TV项和小波变换项对应的参数,
再通过Split Bregman方法对双约束优化问题进行迭代求解,得到电导率的倒数ρ;
步骤4:对步骤3得到的ρ取倒数同时变为矩阵形式,从而得到介电特性中电导率分布的重建结果。
进一步地,所述步骤1的具体步骤为:
步骤1.1:通过磁共振成像(MRI)设备对成像物体进行扫描,利用体线圈作为发射线圈,相控阵列线圈作为接收线圈,选择自旋回波(spin echo,SE)序列或者梯度回波(gradient echo,GRE)序列,在参数设置中保存收发相位数据,扫描结束便可以获得重建所需要的MR收发相位数据;
步骤1.2:对MR收发相位数据进行高斯滤波处理,得到收发相位
进一步地,所述步骤2的具体步骤为:
步骤2.1:将收发相位代入对流反应方程中,所述对流反应方程为:
其中,u为电导率的倒数,u=1/σ,
为u的梯度,采用中心差分方式进行表达,
为收发相位的梯度,
为收发相位的拉普拉斯运算,采用中心差分方式计算得到,
因此,对流反应方程展开后为:
步骤2.2:设定成像物体的大小为矩阵N×M×L,对矩阵N×M×L中心的任意一点ui,j,k(i=2,3…N-1,j=2,3…M-1,k=2,3…L-1),采用中心差分方式对其梯度进行表示,得到该点的对流反应方程:
其中,Δx,Δy,Δz分别是x,y,z方向上的分辨率,将上述方程展开后为:
步骤2.3:对矩阵N×M×L处在边界上的点ui,j,k(i=1,N,j=1,M,k=1,L),采用前向差分或后向差分方式对其梯度进行表示,得到与步骤2.2相对应的对流反应方程,所述前向差分为:
所述后向差分为:
步骤2.4:将大小为N×M×L的矩阵u变为一个列向量后,根据步骤2.2和步骤2.3对每个点表达的结果,构建出与ρ对应的,大小为N*M*L的方阵A;
步骤2.5:根据步骤2.4得到的方阵A,通过有限差分方法对步骤2.1的对流反应方程进行展开,得到线性系统方程,所述线性系统方程为:
Aρ=b,
其中,b是一个列向量,每一点的值都为2ωμ0。
进一步地,所述步骤3的具体步骤为:
步骤3.1:将线性系统方程采用最小二乘方式构造双约束优化问题,所述双约束优化问题的公式为:
其中,为全变分(total variation,TV),W是小波变换矩阵,λ1,λ2分别为TV项和小波变换项对应的参数;
步骤3.2:通过Split Bregman方法对双约束优化问题进行求解,对TV项和小波变换项进行替换,将替换为d,将Wρ替换为w,所述双约束优化问题的公式变为:
其中,btv和bw来自于Bregman迭代;
步骤3.3:将步骤3.2的双约束优化问题分成三个子问题:
子问题1:
子问题2:
子问题3:
步骤3.4:对步骤3.3的三个子问题分别进行求解,对步骤3.2的btv和bw进行迭代更新,得到ρ。
更进一步地,所述步骤3.4的具体步骤为:
步骤3.4.1:选定参数λ1,λ2以及误差τ,对d,w,btv,bw进行初始化,初始值均可以设置为0;
步骤3.4.2:对步骤3.3的三个子问题分别进行求解,
因为子问题1是可导的,通过对ρ进行求导并令其等于0,得到ρ进行更新的公式:
其中,WTW=I,再通过梯度下降法求解,
通过软阈值的方式对子问题2和子问题3进行求解,得到各自迭代更新的公式:
其中,
步骤3.4.3:对btv,bw进行迭代更新:
步骤3.4.4:判断||ρk+1-ρk||2<τ是否满足,若满足,则停止迭代,输出ρk,若不满足,则返回步骤3.4.2继续进行迭代。
综上所述,由于采用了上述技术方案,本发明的有益效果是:
1.一种基于双约束的介电特性迭代成像方法,首先通过磁共振成像设备获得成像所需要的MR收发相位数据,在线圈的选择上更具灵活性,可以有效的适用与临床的应用;然后对MR收发相位数据进行高斯滤波处理,再将得到的收发相位代入对流反应方程中,利用有限差分方法将对流反应方程转换成线性系统方程,采用含有电导率梯度变化的方程进行计算,有效去除常规的均匀性假设带来的误差,提高重建结果的准确性;再将线性系统方程采用最小二乘方式构造全变分和小波的双约束优化问题,再通过Split Bregman方法对该优化问题进行迭代求解,TV和小波正则化可以解决求解过程中的振荡问题,有效抑制噪声的影响;最后将求解得到的结果取倒数从而得到介电特性中电导率分布的重建结果,本发明能够有效的提高电导率重建的准确性,抑制噪声对重建的影响,提高重建图像的质量。
附图说明
为了更清楚地说明本发明实施例的技术方案,下面将对实施例中所需要使用的附图作简单地介绍,应当理解,以下附图仅示出了本发明的某些实施例,因此不应被看作是对范围的限定,对于本领域普通技术人员来讲,在不付出创造性劳动的前提下,还可以根据这些附图获得其他相关的附图,其中:
图1是一种基于双约束的介电特性迭代成像方法的流程图;
图2是本发明实施例一中步骤1得到的真实电导率和收发相位数据图;
图3是本发明实施例一和传统方法进行介电特性成像的结果示意图;
图4是本发明和传统方法在不同的噪声情况下对介电特性成像得到的结果示意图;
图5是本发明和传统方法在不同的噪声情况下对介电特性成像结果的均值和方差统计图,统计区域为体模的中心区域;
图6是本发明和传统方法在不同的噪声情况下对介电特性成像结果的均值和方差统计图,统计区域为体模的外部区域;
图7是本发明和传统方法对人体大脑仿真数据进行介电特性成像的重建结果示意图;
图8是本发明和传统方法对添加了噪声的人体大脑仿真数据进行介电特性成像的重建结果示意图;
图9是本发明和传统方法对自制体模进行扫描得到的收发相位数据进行重建得到的电导率分布示意图。
具体实施方式
为了使本发明的目的、技术方案及优点更加清楚明白,以下结合附图及实施例,对本发明进行进一步详细说明。应当理解,此处所描述的具体实施例仅用以解释本发明,并不用于限定本发明,即所描述的实施例仅仅是本发明一部分实施例,而不是全部的实施例。通常在此处附图中描述和示出的本发明实施例的组件可以以各种不同的配置来布置和设计。
因此,以下对在附图中提供的本发明的实施例的详细描述并非旨在限制要求保护的本发明的范围,而是仅仅表示本发明的选定实施例。基于本发明的实施例,本领域技术人员在没有做出创造性劳动的前提下所获得的所有其他实施例,都属于本发明保护的范围。
需要说明的是,术语“第一”和“第二”等之类的关系术语仅仅用来将一个实体或者操作与另一个实体或操作区分开来,而不一定要求或者暗示这些实体或操作之间存在任何这种实际的关系或者顺序。而且,术语“包括”、“包含”或者其任何其他变体意在涵盖非排他性的包含,从而使得包括一系列要素的过程、方法、物品或者设备不仅包括那些要素,而且还包括没有明确列出的其他要素,或者是还包括为这种过程、方法、物品或者设备所固有的要素。在没有更多限制的情况下,由语句“包括一个……”限定的要素,并不排除在包括所述要素的过程、方法、物品或者设备中还存在另外的相同要素。
一种基于双约束的介电特性迭代成像方法,解决了传统介电特性成像方法存在误差较大,噪声影响严重以及临床适用性差的问题;
一种基于双约束的介电特性迭代成像方法,包括以下步骤:
步骤1:获取MR收发相位数据,进行高斯滤波处理后得到收发相位
步骤2:将收发相位代入对流反应方程中,所述对流反应方程为:
其中,u为电导率的倒数,为u的梯度,为的梯度,为收发相位的拉普拉斯运算,
再通过有限差分方法对对流反应方程进行转换得到线性系统方程,所述线性系统方程为:
Aρ=b,
其中,将u写成列向量形式即为ρ,b是一个列向量,每一点的值都是为2ωμ0;
步骤3:将线性系统方程采用最小二乘方式构造双约束优化问题,所述双约束优化问题的公式为:
其中,为全变分(total variation,TV),W是小波变换矩阵,λ1,λ2分别为TV项和小波变换项对应的参数,
再通过Split Bregman方法对双约束优化问题进行迭代求解,得到电导率的倒数ρ;
步骤4:对步骤3得到的ρ取倒数同时变为矩阵形式,从而得到介电特性中电导率分布的重建结果。
本发明首先通过磁共振成像设备获得成像所需要的MR收发相位数据,在线圈的选择上更具灵活性,可以有效的适用与临床的应用;然后对MR收发相位数据进行高斯滤波处理,再将得到的收发相位代入对流反应方程中,利用有限差分方法将对流反应方程转换成线性系统方程,采用含有电导率梯度变化的方程进行计算,有效去除常规的均匀性假设带来的误差,提高重建结果的准确性;再将线性系统方程采用最小二乘方式构造全变分和小波的双约束优化问题,再通过Split Bregman方法对该优化问题进行迭代求解,TV和小波正则化可以解决求解过程中的振荡问题,有效抑制噪声的影响;最后将求解得到的结果取倒数从而得到介电特性中电导率分布的重建结果,本发明能够有效的提高电导率重建的准确性,抑制噪声对重建的影响,提高重建图像的质量。
下面结合实施例对本发明的特征和性能作进一步的详细描述。
一种基于双约束的介电特性迭代成像方法,包括以下步骤:
步骤1:获取MR收发相位数据,进行高斯滤波处理后得到收发相位
步骤1.1:通过磁共振成像(MRI)设备对成像物体进行扫描,采用电磁仿真软件对鸟笼线圈的仿真数据,通过对鸟笼线圈上的两个激励端口添加激励源对体模或者其他对象进行仿真,发射相位和接收相位分别是在不同的激励下得到的,两次激励的不同在于激励源具有相反的相位信息,将两者相加即得到收发相位数据,利用体线圈作为发射线圈,相控阵列线圈作为接收线圈,选择自旋回波(spin echo,SE)序列或者梯度回波(gradientecho,GRE)序列,在参数设置中保存收发相位数据,扫描结束便可以获得重建所需要的MR收发相位数据,如图2所示;
步骤1.2:对MR收发相位数据进行高斯滤波处理,得到收发相位
步骤2:将收发相位代入对流反应方程中,再通过有限差分方法对对流反应方程进行转换得到线性系统方程;
步骤2.1:将收发相位代入对流反应方程中,所述对流反应方程为:
其中,u为电导率的倒数,u=1/σ,
为u的梯度,采用中心差分方式进行表达,
为收发相位的梯度,
为收发相位的拉普拉斯运算,采用中心差分方式计算得到,
因此,对流反应方程展开后为:
步骤2.2:设定成像物体的大小为矩阵N×M×L,对矩阵N×M×L中心的任意一点ui,j,k(i=2,3…N-1,j=2,3…M-1,k=2,3…L-1),采用中心差分方式对其梯度进行表示,得到该点的对流反应方程:
其中,Δx,Δy,Δz分别是x,y,z方向上的分辨率,将上述方程展开后为:
步骤2.3:对矩阵N×M×L处在边界上的点ui,j,k(i=1,N,j=1,M,k=1,L),采用前向差分或后向差分方式对其梯度进行表示,得到与步骤2.2相对应的对流反应方程,所述前向差分为:
所述后向差分为:
步骤2.4:将大小为N×M×L的矩阵u变为一个列向量后,根据步骤2.2和步骤2.3对每个点表达的结果,构建出与ρ对应的,大小为N*M*L的方阵A;
步骤2.5:根据步骤2.4得到的方阵A,通过有限差分方法对步骤2.1的对流反应方程进行展开,得到线性系统方程,所述线性系统方程为:
Aρ=b,
其中,b是一个列向量,每一点的值都是为2ωμ0。
步骤3:将线性系统方程采用最小二乘方式构造双约束优化问题,再通过SplitBregman方法对双约束优化问题进行迭代求解,得到电导率的倒数ρ;
步骤3.1:将线性系统方程采用最小二乘方式构造双约束优化问题,所述双约束优化问题的公式为:
其中,为全变分(total variation,TV),W是小波变换矩阵,λ1,λ2分别为TV项和小波变换项对应的参数;
步骤3.2:通过Split Bregman方法对双约束优化问题进行求解,对TV项和小波变换项进行替换,将替换为d,将Wρ替换为w,所述双约束优化问题的公式变为:
其中,btv和bw来自于Bregman迭代;
步骤3.3:将步骤3.2的双约束优化问题分成三个子问题:
子问题1:
子问题2:
子问题3:
步骤3.4:对步骤3.3的三个子问题分别进行求解,对步骤3.2的btv和bw进行迭代更新,得到ρ;
步骤3.4.1:选定参数λ1,λ2以及误差τ,对d,w,btv,bw进行初始化,初始值均可以设置为0;
步骤3.4.2:对步骤3.3的三个子问题分别进行求解,
因为子问题1是可导的,通过对ρ进行求导并令其等于0,得到ρ进行更新的公式:
其中,WTW=I,再通过梯度下降法求解,
通过软阈值的方式对子问题2和子问题3进行求解,得到各自迭代更新的公式:
其中,
步骤3.4.3:对btv,bw进行迭代更新:
步骤3.4.4:判断||ρk+1-ρk||2<τ是否满足,若满足,则停止迭代,输出ρk+1,若不满足,则返回步骤3.4.2继续进行迭代。
步骤4:对步骤3得到的ρ取倒数同时变为矩阵形式,从而得到介电特性中电导率分布的重建结果,如图3所示。
根据上述实施方式,对仿真和真实的收发相位数据进行介电特性成像的结果,本实施例选取其中部分实验结果进行分析比较,如图4至图9所示:
图4为本发明和传统方法在不同的噪声情况下对介电特性成像得到的结果示意图,由图可知,传统方法在组织的边界处产生很大的误差,同时随着噪声的增大,重建的结果也越来越差,而采用本发明重建的结果在不同组织的交界面没有明显的误差,同时对噪声的抑制也远远好于传统方法;
图5为本发明和传统方法在不同的噪声情况下对介电特性成像结果的均值和方差统计图,统计区域为体模的中心区域,真实值为1.5(S/m),传统方法对介电特性成像结果的均值分别为:1.2、1.21、1.2,而本发明对介电特性成像结果的均值分别为:1.46、1.47、1.5;图6为本发明和传统方法在不同的噪声情况下对介电特性成像结果的均值和方差统计图,统计区域为体模的外部区域,真实值1.5(S/m),传统方法对对介电特性成像结果的均值分别为:0.95、0.96、0.93,而本发明对介电特性成像结果的均值分别为:1、1.02、1.02,从图5和图6可以看出,本发明重建的结果更接近于真实值,同时误差也很小;
图7为本发明和传统方法对人体大脑仿真数据进行介电特性成像的重建结果示意图,由图可知,相较于真实的结果,传统方法只在均匀的部分能够重建出较好的结果,但在不同组织的交界处则会出现很大的错误,相反的是本发明则是很有效的重建出大脑内部的结构,并且在不同组织之间有明确的边界,且重建的结果更精确;
图8是本发明和传统方法对添加了噪声的人体大脑仿真数据进行介电特性成像的重建结果示意图,由图可知,传统方法已经无法重建出大脑的结构,并且重建结构的噪声也非常的大。相较而言,本发明则是有效了抑制了噪声对结果的影响,能够有效的分辨出大脑的结构,重建的结果更加的准确;
图9为本发明和传统方法对自制体模进行扫描得到的收发相位数据进行重建得到的电导率分布示意图,由图可知,本发明的结果更加的均匀,并且与真实测量得到的电导率的值相比较,本发明的结果更接近于真实值。
通过大量仿真和实际的实验证明本发明能够有效的解决传统方法中因为均匀性假设所造成的成像误差,提高成像的质量;此外该方法可以使用临床常见的线圈配置,因此该方法未来有助于推广到临床应用。
需要说明的是,由于说明书附图不得着色和涂改,所以本发明中部分区别明显的地方比较难以显示,若有必要,可提供彩色图片。
以上所述仅为本发明的较佳实施例而已,并不用以限制本发明的保护范围,任何熟悉本领域的技术人员在本发明的精神和原则之内所作的任何修改、等同替换和改进等,均应包含在本发明的保护范围之内。
Claims (5)
1.一种基于双约束的介电特性迭代成像方法,其特征在于,包括以下步骤:
步骤1:获取MR收发相位数据,进行高斯滤波处理后得到收发相位
步骤2:将收发相位代入对流反应方程中,所述对流反应方程为:
其中,u为电导率的倒数,▽u为u的梯度,为的梯度,为收发相位的拉普拉斯运算,
再通过有限差分方法对对流反应方程进行转换得到线性系统方程,所述线性系统方程为:
Aρ=b,
其中,将u写成列向量形式即为ρ,b是一个列向量,每一点的值都为2ωμ0;
步骤3:将线性系统方程采用最小二乘方式构造双约束优化问题,所述双约束优化问题的公式为:
其中,▽ρ为全变分(total variation,TV),W是小波变换矩阵,λ1,λ2分别为TV项和小波变换项对应的参数,
再通过Split Bregman方法对双约束优化问题进行迭代求解,得到电导率的倒数ρ;
步骤4:对步骤3得到的ρ取倒数同时变为矩阵形式,从而得到介电特性中电导率分布的重建结果。
2.根据权利要求1所述的一种基于双约束的介电特性迭代成像方法,其特征在于:所述步骤1的具体步骤为:
步骤1.1:通过磁共振成像(MRI)设备对成像物体进行扫描,利用体线圈作为发射线圈,相控阵列线圈作为接收线圈,选择自旋回波(spin echo,SE)序列或者梯度回波(gradientecho,GRE)序列,在参数设置中保存收发相位数据,扫描结束便可以获得重建所需要的MR收发相位数据;
步骤1.2:对MR收发相位数据进行高斯滤波处理,得到收发相位
3.根据权利要求1所述的一种基于双约束的介电特性迭代成像方法,其特征在于:所述步骤2的具体步骤为:
步骤2.1:将收发相位代入对流反应方程中,所述对流反应方程为:
其中,u为电导率的倒数,u=1/σ,
▽u为u的梯度,采用中心差分方式进行表达,
为收发相位的梯度,
为收发相位的拉普拉斯运算,采用中心差分方式计算得到,
因此,对流反应方程展开后为:
步骤2.2:设定成像物体的大小为矩阵N×M×L,对矩阵N×M×L中心的任意一点ui,j,k(i=2,3…N-1,j=2,3…M-1,k=2,3…L-1),采用中心差分方式对其梯度进行表示,得到该点的对流反应方程:
其中,Δx,Δy,Δz分别是x,y,z方向上的分辨率,将上述方程展开后为:
步骤2.3:对矩阵N×M×L处在边界上的点ui,j,k(i=1,N,j=1,M,k=1,L),采用前向差分或后向差分方式对其梯度进行表示,得到与步骤2.2相对应的对流反应方程,所述前向差分为:
所述后向差分为:
步骤2.4:将大小为N×M×L的矩阵u变为一个列向量后,根据步骤2.2和步骤2.3对每个点表达的结果,构建出与ρ对应的,大小为N*M*L的方阵A;
步骤2.5:根据步骤2.4得到的方阵A,通过有限差分方法对步骤2.1的对流反应方程进行展开,得到线性系统方程,所述线性系统方程为:
Aρ=b,
其中,b是一个列向量,每一点的值都为2ωμ0。
4.根据权利要求1所述的一种基于双约束的介电特性迭代成像方法,其特征在于:所述步骤3的具体步骤为:
步骤3.1:将线性系统方程采用最小二乘方式构造双约束优化问题,所述双约束优化问题的公式为:
其中,▽ρ为全变分(total variation,TV),W是小波变换矩阵,λ1,λ2分别为TV项和小波变换项对应的参数;
步骤3.2:通过Split Bregman方法对双约束优化问题进行求解,对TV项和小波变换项进行替换,将▽ρ替换为d,将Wρ替换为w,所述双约束优化问题的公式变为:
其中,btv和bw来自于Bregman迭代;
步骤3.3:将步骤3.2的双约束优化问题分成三个子问题:
子问题1:
子问题2:
子问题3:
步骤3.4:对步骤3.3的三个子问题分别进行求解,对步骤3.2的btv和bw进行迭代更新,得到ρ。
5.根据权利要求4所述的一种基于双约束的介电特性迭代成像方法,其特征在于:所述步骤3.4的具体步骤为:
步骤3.4.1:选定参数λ1,λ2以及误差τ,对d,w,btv,bw进行初始化,初始值均可以设置为0;
步骤3.4.2:对步骤3.3的三个子问题分别进行求解,
因为子问题1是可导的,通过对ρ进行求导并令其等于0,得到ρ进行更新的公式:
其中,▽T▽=Δ,WTW=I,再通过梯度下降法求解,
通过软阈值的方式对子问题2和子问题3进行求解,得到各自迭代更新的公式:
其中,
步骤3.4.3:对btv,bw进行迭代更新:
步骤3.4.4:判断||ρk+1-ρk||2<τ是否满足,若满足,则停止迭代,输出ρk+1,若不满足,则返回步骤3.4.2继续进行迭代。
Priority Applications (1)
Application Number | Priority Date | Filing Date | Title |
---|---|---|---|
CN201811215000.8A CN109471053B (zh) | 2018-10-18 | 2018-10-18 | 一种基于双约束的介电特性迭代成像方法 |
Applications Claiming Priority (1)
Application Number | Priority Date | Filing Date | Title |
---|---|---|---|
CN201811215000.8A CN109471053B (zh) | 2018-10-18 | 2018-10-18 | 一种基于双约束的介电特性迭代成像方法 |
Publications (2)
Publication Number | Publication Date |
---|---|
CN109471053A true CN109471053A (zh) | 2019-03-15 |
CN109471053B CN109471053B (zh) | 2020-01-31 |
Family
ID=65665721
Family Applications (1)
Application Number | Title | Priority Date | Filing Date |
---|---|---|---|
CN201811215000.8A Active CN109471053B (zh) | 2018-10-18 | 2018-10-18 | 一种基于双约束的介电特性迭代成像方法 |
Country Status (1)
Country | Link |
---|---|
CN (1) | CN109471053B (zh) |
Cited By (3)
Publication number | Priority date | Publication date | Assignee | Title |
---|---|---|---|---|
CN110426704A (zh) * | 2019-08-20 | 2019-11-08 | 中国科学院重庆绿色智能技术研究院 | 一种用于稀疏阵列的全变差快速成像算法 |
CN110687489A (zh) * | 2019-09-04 | 2020-01-14 | 华东师范大学 | 一种磁共振成像多回波磁场分布图拟合方法 |
CN112345989A (zh) * | 2020-11-18 | 2021-02-09 | 中国科学院电工研究所 | 一种肿瘤组织磁特性成像方法 |
Citations (14)
Publication number | Priority date | Publication date | Assignee | Title |
---|---|---|---|---|
US20070110290A1 (en) * | 2005-10-19 | 2007-05-17 | Siemens Corporate Research Inc. | Devices Systems and Methods for Processing Images |
JP2008008793A (ja) * | 2006-06-29 | 2008-01-17 | Tohoku Univ | 高温融体導電材料の熱物性測定方法及び測定装置 |
WO2010129711A1 (en) * | 2009-05-05 | 2010-11-11 | The Trustees Of Columbia University In The City Of New York | Devices, systems, and methods for evaluating vision and diagnosing and compensating impairment of vision |
US20110263478A1 (en) * | 2008-09-16 | 2011-10-27 | Simard J Marc | Sur1 inhibitors for therapy |
WO2012061475A3 (en) * | 2010-11-02 | 2012-07-19 | University Of Florida Research Foundation, Inc. | Systems and methods for fast magnetic resonance image reconstruction |
US20130228683A1 (en) * | 2011-08-10 | 2013-09-05 | Fei Company | Charged-Particle Microscope Providing Depth-Resolved Imagery |
CN105678822A (zh) * | 2016-01-13 | 2016-06-15 | 哈尔滨理工大学 | 一种基于Split Bregman迭代的三正则磁共振图像重构方法 |
CN105678715A (zh) * | 2016-02-15 | 2016-06-15 | 江苏大学 | 一种基于核回归全变分的去除图像噪声的方法 |
CN105814574A (zh) * | 2013-10-04 | 2016-07-27 | 塞昆纳姆股份有限公司 | 遗传变异的非侵入性评估的方法和过程 |
EP3171163A1 (en) * | 2015-11-18 | 2017-05-24 | FEI Company | Novel x-ray imaging technique |
CN107515424A (zh) * | 2017-07-26 | 2017-12-26 | 山东科技大学 | 一种基于vmd与小波包的微震信号降噪滤波方法 |
EP3346443A1 (en) * | 2015-09-04 | 2018-07-11 | Samsung Electronics Co., Ltd. | Method for restoring magnetic resonance image and magnetic resonance image processing apparatus |
CN108333628A (zh) * | 2018-01-17 | 2018-07-27 | 中国石油大学(华东) | 基于正则化约束的弹性波最小二乘逆时偏移方法 |
CN108416723A (zh) * | 2018-02-07 | 2018-08-17 | 南京理工大学 | 基于全变分正则化和变量分裂的无透镜成像快速重构方法 |
-
2018
- 2018-10-18 CN CN201811215000.8A patent/CN109471053B/zh active Active
Patent Citations (14)
Publication number | Priority date | Publication date | Assignee | Title |
---|---|---|---|---|
US20070110290A1 (en) * | 2005-10-19 | 2007-05-17 | Siemens Corporate Research Inc. | Devices Systems and Methods for Processing Images |
JP2008008793A (ja) * | 2006-06-29 | 2008-01-17 | Tohoku Univ | 高温融体導電材料の熱物性測定方法及び測定装置 |
US20110263478A1 (en) * | 2008-09-16 | 2011-10-27 | Simard J Marc | Sur1 inhibitors for therapy |
WO2010129711A1 (en) * | 2009-05-05 | 2010-11-11 | The Trustees Of Columbia University In The City Of New York | Devices, systems, and methods for evaluating vision and diagnosing and compensating impairment of vision |
WO2012061475A3 (en) * | 2010-11-02 | 2012-07-19 | University Of Florida Research Foundation, Inc. | Systems and methods for fast magnetic resonance image reconstruction |
US20130228683A1 (en) * | 2011-08-10 | 2013-09-05 | Fei Company | Charged-Particle Microscope Providing Depth-Resolved Imagery |
CN105814574A (zh) * | 2013-10-04 | 2016-07-27 | 塞昆纳姆股份有限公司 | 遗传变异的非侵入性评估的方法和过程 |
EP3346443A1 (en) * | 2015-09-04 | 2018-07-11 | Samsung Electronics Co., Ltd. | Method for restoring magnetic resonance image and magnetic resonance image processing apparatus |
EP3171163A1 (en) * | 2015-11-18 | 2017-05-24 | FEI Company | Novel x-ray imaging technique |
CN105678822A (zh) * | 2016-01-13 | 2016-06-15 | 哈尔滨理工大学 | 一种基于Split Bregman迭代的三正则磁共振图像重构方法 |
CN105678715A (zh) * | 2016-02-15 | 2016-06-15 | 江苏大学 | 一种基于核回归全变分的去除图像噪声的方法 |
CN107515424A (zh) * | 2017-07-26 | 2017-12-26 | 山东科技大学 | 一种基于vmd与小波包的微震信号降噪滤波方法 |
CN108333628A (zh) * | 2018-01-17 | 2018-07-27 | 中国石油大学(华东) | 基于正则化约束的弹性波最小二乘逆时偏移方法 |
CN108416723A (zh) * | 2018-02-07 | 2018-08-17 | 南京理工大学 | 基于全变分正则化和变量分裂的无透镜成像快速重构方法 |
Non-Patent Citations (6)
Title |
---|
LUSTIG M,ET AL: "M.Sparse MRI:The application of compressed sensing for rapid MR imaging", 《MAGNETIC RESONANCE IN MEDICINE》 * |
余丽红 等: "基于自适应正则化的全变分去噪算法", 《中国图象图形学报》 * |
刘晓晖 等: "基于Moreau_包络的近似平滑迭代磁共振图像重建算法", 《计算机应用》 * |
刘晓晖 等: "磁共振图像的原始-对偶近似迭代重建算法", 《计算机工程与应用》 * |
蒋明峰 等: "基于全变分扩展方法的压缩感知磁共振成像算法研究", 《电子与信息学报》 * |
邓红耀: "基于核范数和扩散的图像复原研究", 《中国博士学位论文全文数据库 信息科技辑》 * |
Cited By (6)
Publication number | Priority date | Publication date | Assignee | Title |
---|---|---|---|---|
CN110426704A (zh) * | 2019-08-20 | 2019-11-08 | 中国科学院重庆绿色智能技术研究院 | 一种用于稀疏阵列的全变差快速成像算法 |
CN110426704B (zh) * | 2019-08-20 | 2023-03-24 | 中国科学院重庆绿色智能技术研究院 | 一种用于稀疏阵列的全变差快速成像算法 |
CN110687489A (zh) * | 2019-09-04 | 2020-01-14 | 华东师范大学 | 一种磁共振成像多回波磁场分布图拟合方法 |
CN110687489B (zh) * | 2019-09-04 | 2021-08-20 | 华东师范大学 | 一种磁共振成像多回波磁场分布图拟合方法 |
CN112345989A (zh) * | 2020-11-18 | 2021-02-09 | 中国科学院电工研究所 | 一种肿瘤组织磁特性成像方法 |
CN112345989B (zh) * | 2020-11-18 | 2024-05-28 | 中国科学院电工研究所 | 一种肿瘤组织磁特性成像方法 |
Also Published As
Publication number | Publication date |
---|---|
CN109471053B (zh) | 2020-01-31 |
Similar Documents
Publication | Publication Date | Title |
---|---|---|
Buonincontri et al. | MR fingerprinting with simultaneous B1 estimation | |
Zhao et al. | Maximum likelihood reconstruction for magnetic resonance fingerprinting | |
US7800368B2 (en) | High field magnetic resonance | |
US8653818B2 (en) | Parallel transmission RF pulse design with local SAR constraints | |
EP1915089B1 (en) | Electric impedance imaging system, method, and computer program | |
US9983284B2 (en) | MRI with dixon-type water/fat separation and prior knowledge about inhomogeneity of the main magnetic field | |
Childs et al. | Impact of number of channels on RF shimming at 3T | |
JP6270709B2 (ja) | 磁気共鳴生データの再構成方法および装置 | |
Klauser et al. | Fast high‐resolution brain metabolite mapping on a clinical 3T MRI by accelerated H‐FID‐MRSI and low‐rank constrained reconstruction | |
RU2582474C2 (ru) | Эмуляция виртуальных катушек в мрт с параллельной передачей | |
US20180231626A1 (en) | Systems and methods for magnetic resonance fingerprinting for quantitative breast imaging | |
CN104781685A (zh) | 用于使用非相干采样和冗余haar子波的动态磁共振成像的图像重建 | |
US10018691B2 (en) | Correction of imaging methods in a magnetic resonance device | |
WO2013181703A1 (en) | Method of estimating specific absorption rate | |
CN109471053A (zh) | 一种基于双约束的介电特性迭代成像方法 | |
WO2011087847A2 (en) | Accelerated mri with nonlineear spatial encoding gradients | |
CN112370040B (zh) | 磁共振成像方法、装置、存储介质及电子设备 | |
Börnert et al. | A half-century of innovation in technology—preparing MRI for the 21st century | |
Liu et al. | Four‐dimensional diffusion‐weighted MR imaging (4D‐DWI): a feasibility study | |
van der Heide et al. | High‐resolution in vivo MR‐STAT using a matrix‐free and parallelized reconstruction algorithm | |
Liu et al. | MR-based electrical property tomography using a modified finite difference scheme | |
US20210166384A1 (en) | System, method and computer-accessible medium for tissue fingerprinting | |
Stijnman et al. | Transceive phase corrected 2D contrast source inversion‐electrical properties tomography | |
Shimron et al. | CORE‐PI: Non‐iterative convolution‐based reconstruction for parallel MRI in the wavelet domain | |
Rincón‐Domínguez et al. | Accelerated multi‐snapshot free‐breathing mapping based on the dual refocusing echo acquisition mode technique (DREAM): An alternative to measure RF nonuniformity for cardiac MRI |
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 |