CN110618463B - 一种核磁共振数据反演方法、装置、存储介质及设备 - Google Patents

一种核磁共振数据反演方法、装置、存储介质及设备 Download PDF

Info

Publication number
CN110618463B
CN110618463B CN201910822315.7A CN201910822315A CN110618463B CN 110618463 B CN110618463 B CN 110618463B CN 201910822315 A CN201910822315 A CN 201910822315A CN 110618463 B CN110618463 B CN 110618463B
Authority
CN
China
Prior art keywords
vector
magnetic resonance
iteration
nuclear magnetic
norm
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.)
Active
Application number
CN201910822315.7A
Other languages
English (en)
Other versions
CN110618463A (zh
Inventor
郭江峰
谢然红
肖立志
徐陈昱
谷明宣
吴勃翰
Current Assignee (The listed assignees may be inaccurate. Google has not performed a legal analysis and makes no representation or warranty as to the accuracy of the list.)
China University of Petroleum Beijing
Original Assignee
China University of Petroleum Beijing
Priority date (The priority date 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 date listed.)
Filing date
Publication date
Application filed by China University of Petroleum Beijing filed Critical China University of Petroleum Beijing
Priority to CN201910822315.7A priority Critical patent/CN110618463B/zh
Publication of CN110618463A publication Critical patent/CN110618463A/zh
Application granted granted Critical
Publication of CN110618463B publication Critical patent/CN110618463B/zh
Active legal-status Critical Current
Anticipated expiration legal-status Critical

Links

Images

Classifications

    • GPHYSICS
    • G01MEASURING; TESTING
    • G01VGEOPHYSICS; GRAVITATIONAL MEASUREMENTS; DETECTING MASSES OR OBJECTS; TAGS
    • G01V3/00Electric or magnetic prospecting or detecting; Measuring magnetic field characteristics of the earth, e.g. declination, deviation
    • G01V3/14Electric or magnetic prospecting or detecting; Measuring magnetic field characteristics of the earth, e.g. declination, deviation operating with electron or nuclear magnetic resonance
    • GPHYSICS
    • G01MEASURING; TESTING
    • G01VGEOPHYSICS; GRAVITATIONAL MEASUREMENTS; DETECTING MASSES OR OBJECTS; TAGS
    • G01V3/00Electric or magnetic prospecting or detecting; Measuring magnetic field characteristics of the earth, e.g. declination, deviation
    • G01V3/38Processing data, e.g. for analysis, for interpretation, for correction

Landscapes

  • Life Sciences & Earth Sciences (AREA)
  • Physics & Mathematics (AREA)
  • Engineering & Computer Science (AREA)
  • Environmental & Geological Engineering (AREA)
  • Geology (AREA)
  • Remote Sensing (AREA)
  • General Life Sciences & Earth Sciences (AREA)
  • General Physics & Mathematics (AREA)
  • Geophysics (AREA)
  • High Energy & Nuclear Physics (AREA)
  • Magnetic Resonance Imaging Apparatus (AREA)

Abstract

本申请实施方式公开了一种核磁共振数据反演方法、装置、存储介质及设备。所述方法包括:获取核磁共振回波数据;基于向量的Lp范数,建立目标函数;其中,所述向量满足预设的约束方程,所述约束方程通过所述核磁共振回波数据确定,所述Lp范数中的参数p的取值范围为大于等于0且小于等于1;根据所述目标函数的解,生成所述核磁共振回波数据的反演谱。本申请可以提高核磁共振反演谱的分辨率。

Description

一种核磁共振数据反演方法、装置、存储介质及设备
技术领域
本申请涉及油气勘探中的核磁共振测井数据处理技术领域,特别涉及一种核磁共振数据反演方法、装置、存储介质及设备。
背景技术
核磁共振测井利用氢原子核与磁场相互作用发生的共振现象来实现井下油气储层信息的观测,是一种有效的探测岩石物理性质的方法。核磁共振测井测量地层的原始数据是成千上万个回波组成的回波串,通常将测量得到的回波串反演成核磁共振谱用于储层岩石物理参数的计算,从而反演核磁共振谱的精度直接关系着计算的岩石物理参数的准确性。所以,研究高精度的核磁共振数据反演方法非常重要。
目前,已有的核磁共振数据反演方法包括正则化方法、迭代方法以及智能优化方法。这些方法在反演核磁共振数据时,目标函数都包括残差L2范数这一项,导致反演得到的解的稀疏性差,从而得到的核磁共振反演谱的分辨率较低。因此,亟需一种核磁共振数据反演方法,以提高核磁共振反演谱的分辨率。
发明内容
本申请实施例的目的是提供一种核磁共振数据反演方法、装置、存储介质及设备,以提高核磁共振反演谱的分辨率。
为达到上述目的,本申请实施例提供一种核磁共振数据反演方法,包括:
获取核磁共振回波数据;
基于向量的Lp范数,建立目标函数;其中,所述向量满足预设的约束方程,所述约束方程通过所述核磁共振回波数据确定,所述Lp范数中的参数p的取值范围为大于等于0且小于等于1;
根据所述目标函数的解,生成所述核磁共振回波数据的反演谱。
在一个实施例中,所述目标函数为:
Figure BDA0002187902520000011
其中,
Figure BDA0002187902520000012
为向量的Lp范数;向量s为所述目标函数的解,向量s中的所有元素都大于等于0;p为Lp范数的参数;
所述约束方程为:
Ks=y
其中,向量y为不同时刻的回波幅度组成的向量,所述不同时刻的回波幅度通过所述核磁共振回波数据得到;K为核矩阵,所述核矩阵K中的元素为:
Figure BDA0002187902520000021
其中,t为回波幅度衰减时间,T2为横向弛豫时间。
在一个实施例中,利用以下步骤得到所述目标函数的解:
设置参数ε、阈值tol、Lp范数中的参数p、最大迭代次数itmax,并将迭代次数it的初始值设置为1;所述参数ε和所述阈值tol的取值范围为大于0;
计算向量s(0)=K\y,并确定向量s(0)中的元素个数n;
计算初始残差
Figure BDA0002187902520000022
进行第it次迭代计算,所述第it次迭代计算包括:计算第it次迭代后的解s(it),以及第it次迭代后的残差r(it);所述第it次迭代后的解s(it),以及第it次迭代后的残差r(it)通过如下方法得到:
s(it)=WKT/(KWKT)y
Figure BDA0002187902520000023
其中,对角矩阵W=diag(1./w),向量w=[w1,w2,…wi,…wn],
Figure BDA0002187902520000024
在第1次迭代时,
Figure BDA0002187902520000025
为所述向量s(0)中的第i个元素;在迭代次数it的值大于等于2时,
Figure BDA0002187902520000026
为前一次迭代后的解s(it-1)中的第i个元素;
在|r(it)-r|<tol;或,|r(it)-r|≥tol且it=itmax的情况下,令s=s(it),得到所述向量s。
在一个实施例中,还包括:
在|r(it)-r|≥tol且it<itmax的情况下,令s=s(it)、r=r(it)
将迭代次数it增加1后,进行第it次迭代计算。
在一个实施例中,所述令s=s(it),得到向量s,包括:
令s=s(it)
若向量s中含有小于0的元素,将所述小于0的元素设置为0,得到向量s。
本申请实施例还提供一种核磁共振数据反演装置,所述装置包括:
回波数据获取模块,用于获取核磁共振回波数据;
目标函数建立模块,用于基于向量的Lp范数,建立目标函数;其中,所述向量满足预设的约束方程,所述约束方程通过所述核磁共振回波数据确定,所述Lp范数中的参数p的取值范围为大于等于0且小于等于1;
反演谱生成模块,用于根据所述目标函数的解,生成所述核磁共振回波数据的反演谱。
在一个实施例中,所述反演谱生成模块,包括:
参数设置单元,用于设置参数ε、阈值tol、Lp范数中的参数p、最大迭代次数itmax,并将迭代次数it的初始值设置为1;所述参数ε和所述阈值tol的取值范围为大于0;
初始解计算单元,用于计算向量s(0)=K\y,并确定向量s(0)中的元素个数n;
初始残差计算单元,用于计算初始残差
Figure BDA0002187902520000031
迭代计算单元,用于进行第it次迭代计算,所述第it次迭代计算包括:计算第it次迭代后的解s(it),以及第it次迭代后的残差r(it);所述第it次迭代后的解s(it),以及第it次迭代后的残差r(it)通过如下方法得到:
s(it)=WKT/(KWKT)y
Figure BDA0002187902520000032
其中,对角矩阵W=diag(1./w),向量w=[w1,w2,…wi,…wn],
Figure BDA0002187902520000033
在第1次迭代时,
Figure BDA0002187902520000034
为所述向量s(0)中的第i个元素;在迭代次数it的值大于等于2时,
Figure BDA0002187902520000035
为前一次迭代后的解s(it-1)中的第i个元素;
第一确定单元,用于在|r(it)-r|<tol;或,|r(it)-r|≥tol且it=itmax的情况下,令s=s(it),得到所述向量s。
在一个实施例中,所述反演谱生成模块还包括:
第二确定单元,用于在|r(it)-r|≥tol且it<itmax的情况下,令s=s(it)、r=r(it);并将迭代次数it增加1后,进行所述第it次迭代计算。
本申请实施例还提供一种计算机设备,包括处理器以及用于存储处理器可执行指令的存储器,所述处理器执行所述指令时实现上述任意实施例中所述核磁共振数据反演方法的步骤。
本申请实施例还提供一种计算机可读存储介质,其上存储有计算机指令,所述指令被执行时实现上述任意实施例中所述核磁共振数据反演方法的步骤。
由以上本申请实施例提供的技术方案可见,本申请中的目标函数仅含有Lp范数(Lp范数中的参数p的取值范围为大于等于0且小于等于1),不含残差L2范数,因此,利用本申请提供的核磁共振数据反演方法,反演的解稀疏性好,得到的反演谱分辨率较高。
附图说明
为了更清楚地说明本申请实施例或现有技术中的技术方案,下面将对实施例或现有技术描述中所需要使用的附图作简单地介绍,显而易见地,下面描述中的附图仅仅是本申请中记载的一些实施例,对于本领域普通技术人员来讲,在不付出创造性劳动性的前提下,还可以根据这些附图获得其他的附图。
图1是本申请实施例提供的一种核磁共振数据反演方法的流程图;
图2是本申请实施例提供的一个含有双峰的T2谱模型图;
图3是本申请实施例提供的正演的回波间隔为0.15ms,回波个数为3000,信噪比为150的回波串;
图4是本申请实施例提供的正演的回波间隔为0.15ms,回波个数为3000,信噪比为50的回波串;
图5是本申请实施例提供的正演的回波间隔为0.15ms,回波个数为3000,信噪比为10的回波串;
图6是本申请实施例提供的三种不同反演方法分别处理压缩后的图3中回波数据得到的核磁共振T2谱与模型的对比结果图;
图7是本申请实施例提供的三种不同反演方法分别处理压缩后的图4中回波数据得到的核磁共振T2谱与模型的对比结果图;
图8是本申请实施例提供的三种不同反演方法分别处理压缩后的图5中回波数据得到的核磁共振T2谱与模型的对比结果图;
图9是本申请实施例提供的一种核磁共振数据反演装置的模块结构图;
图10是本申请实施例提供的计算机设备的示意图。
具体实施方式
本申请实施例提供一种核磁共振数据反演方法、装置、存储介质及设备。
为了使本技术领域的人员更好地理解本申请中的技术方案,下面将结合本申请实施方式中的附图,对本申请实施方式中的技术方案进行清楚、完整地描述,显然,所描述的实施方式仅仅是本申请一部分实施方式,而不是全部的实施方式。基于本申请中的实施方式,本领域普通技术人员在没有作出创造性劳动前提下所获得的所有其他实施方式,都应当属于本申请保护的范围。
参考图1所示,为本申请实施例提供的一种核磁共振数据反演方法的流程图,可以包括如下步骤:
S101:获取核磁共振回波数据。
S102:基于向量的Lp范数,建立目标函数;其中,所述向量满足预设的约束方程,所述约束方程通过所述核磁共振回波数据确定,所述Lp范数中的参数p的取值范围为大于等于0且小于等于1。
在一维核磁共振数据反演问题可以理解为求解如下公式所示的积分方程中的s(T2),公式为:
Figure BDA0002187902520000051
其中,y(t)为t时刻的回波幅度,t为回波幅度衰减时间,s(T2)为待求的核磁共振T2谱,即核磁共振反演谱。
上面的公式可以用矩阵方程的形式表示成如下形式:
y=Ks
其中,向量y为不同时刻的回波幅度组成的向量,所述不同时刻的回波幅度通过所述核磁共振回波数据得到;K为核矩阵,所述核矩阵K中的元素为:
Figure BDA0002187902520000052
其中,t为回波幅度衰减时间,T2为横向弛豫时间。
上面的矩阵方程可以作为约束条件,即约束方程,对核磁共振反演的解进行约束,进一步的,基于Lp范数建立用于核磁共振反演的目标函数,并求解出使得该目标函数取到最小值的解,目标函数为如下形式:
Figure BDA0002187902520000053
其中,
Figure BDA0002187902520000054
为向量的Lp范数;向量s为所述目标函数的解,向量s中的所有元素都大于等于0;p为Lp范数的参数。
S103:基于所述目标函数的解,得到所述核磁共振回波数据的反演谱。
在得到目标函数后,还需要对其进行迭代求解,具体的,可以采用IRLS(Iteratively Reweighted Least Squares,迭代重加权最小二乘)算法求解,包括如下步骤:
设置参数ε(ε>0)、阈值tol、Lp范数的参数p(0≤p≤1)、最大迭代次数itmax,并将迭代次数it的初始值设置为1;
计算向量s(0)=K\y,并确定向量s(0)中的元素个数n;
计算初始残差
Figure BDA0002187902520000061
进行第it次迭代计算,所述第it次迭代计算包括:计算第it次迭代后的解s(it),以及第it次迭代后的残差r(it);所述第it次迭代后的解s(it),以及第it次迭代后的残差r(it)通过如下方法得到:
s(it)=WKT/(KWKT)y
Figure BDA0002187902520000062
其中,对角矩阵W=diag(1./w),向量w=[w1,w2,…wi,…wn],
Figure BDA0002187902520000063
在第1次迭代时,
Figure BDA0002187902520000064
为所述向量s(0)中的第i个元素;在迭代次数it的值大于等于2时,
Figure BDA0002187902520000065
为前一次迭代后的解s(it-1)中的第i个元素;
在|r(it)-r|<tol;或,|r(it)-r|≥tol且it=itmax的情况下,令s=s(it),若向量s中含有小于0的元素,将所述小于0的元素设置为0,从而得到所述向量s。
在|r(it)-r|≥tol且it<itmax的情况下,令s=s(it)、r=r(it);将迭代次数it增加1后,进行第it次迭代计算。
经过上面的迭代求解,得到目标函数的最优解,即向量s,利用向量s和预先选取的横向弛豫时间构成的T2向量,即可在对数坐标中得到所求的核磁共振的反演谱。
下面以核磁共振T2谱数据反演为例,进行数值模拟实验,来验证本说明书实施例提供的核磁共振数据反演方法的效果。
如图2所示,是数值模拟实验构造的含有双峰T2谱,两个峰对应的T2值分别为5ms和100ms,然后向正演结果中加入一定的高斯白噪声,模拟得到不同信噪比(SNR)的回波串数据。
图3、图4、图5为正演的回波串数据,其中,图3中回波间隔为0.15ms,回波个数为3000,SNR为150;图4中回波间隔为0.15ms,回波个数为3000,SNR为50;图5中回波间隔为0.15ms,回波个数为3000,SNR为10。将图3、图4、图5中的回波数据个数压缩到7,再分别利用截断奇异值分解方解(TSVD)法、Butler-Reeds-Dawson(BRD)法以及本申请提供的IRLS法处理压缩后的回波数据,得到不同信噪比数据对应的核磁共振T2谱,如图6、图7、图8所示,图6是三种不同反演方法分别处理压缩后的图3中的回波数据得到的核磁共振T2谱对比结果图;图7是三种不同反演方法分别处理压缩后的图4中的回波数据得到的核磁共振T2谱对比结果图;图8是三种不同反演方法分别处理压缩后的图5中的回波数据得到的核磁共振T2谱对比结果图。
由图6、图7、图8可知,随着信噪比的增高,本申请提出的IRLS法反演结果准确性越来越高,且其反演结果优于目前常用的TSVD方法以及BRD法的反演结果,尤其在信噪比很低时,本申请提出的IRLS法反演得到的T2谱幅度更高,与模型更接近。
由本说明书提供的实施例可以看出,本申请可以实现如下技术效果:
在本说明书提供的一个实施例中,建立的目标函数中不含残差L2范数,提高了目标函数的解的稀疏性,从而提高了反演谱的分辨率;
在本说明书提供的一个实施例中,目标函数中不含正则化项,使得反演过程中不需要求取正则化参数,简化了反演步骤。
如图9所示,本申请实施例还提供一种核磁共振数据反演装置,包括:
回波数据获取模块10,用于获取核磁共振回波数据;
目标函数建立模块20,用于基于向量的Lp范数,建立目标函数;其中,所述向量满足预设的约束方程,所述约束方程通过所述核磁共振回波数据确定,所述Lp范数中的参数p的取值范围为大于等于0且小于等于1;
反演谱生成模块30,用于基于所述目标函数的解,得到所述核磁共振回波数据的反演谱。
其中,所述反演谱生成模块30,包括:
参数设置单元301,用于设置参数ε(ε>0)、阈值tol、Lp范数的参数p(0≤p≤1)、最大迭代次数itmax,并将迭代次数it的初始值设置为1;
初始解计算单元302,用于计算向量s(0)=K\y,并确定向量s(0)中的元素个数n;
初始残差计算单元303,用于计算初始残差
Figure BDA0002187902520000071
迭代计算单元304,用于进行第it次迭代计算,所述第it次迭代计算包括:计算第it次迭代后的解s(it),以及第it次迭代后的残差r(it);所述第it次迭代后的解s(it),以及第it次迭代后的残差r(it)通过如下方法得到:
s(it)=WKT/(KWKT)y
Figure BDA0002187902520000072
其中,对角矩阵W=diag(1./w),向量w=[w1,w2,…wi,…wn],
Figure BDA0002187902520000081
在第1次迭代时,
Figure BDA0002187902520000082
为所述向量s(0)中的第i个元素;在迭代次数it的值大于等于2时,
Figure BDA0002187902520000083
为前一次迭代后的解s(it-1)中的第i个元素;
第一确定单元305,用于在|r(it)-r|<tol;或,|r(it)-r|≥tol且it=itmax的情况下,令s=s(it),得到所述向量s。
第二确定单元306,用于在|r(it)-r|≥tol且it<itmax的情况下,令s=s(it)、r=r(it);并将迭代次数it增加1后,进行所述第it次迭代计算。
如图10所示,本申请实施例还提供一种计算机设备,包括处理器以及用于存储处理器可执行指令的存储器,所述处理器执行所述指令时实现上述任意实施例中所述核磁共振数据反演方法的步骤。
本申请实施例还提供一种计算机可读存储介质,其上存储有计算机指令,所述指令被执行时实现上述任意实施例中所述核磁共振数据反演方法的步骤。
在20世纪90年代,对于一个技术的改进可以很明显地区分是硬件上的改进(例如,对二极管、晶体管、开关等电路结构的改进)还是软件上的改进(对于方法流程的改进)。然而,随着技术的发展,当今的很多方法流程的改进已经可以视为硬件电路结构的直接改进。设计人员几乎都通过将改进的方法流程编程到硬件电路中来得到相应的硬件电路结构。因此,不能说一个方法流程的改进就不能用硬件实体模块来实现。例如,可编程逻辑器件(Programmable Logic Device,PLD)(例如现场可编程门阵列(Field Programmable GateArray,FPGA))就是这样一种集成电路,其逻辑功能由用户对器件编程来确定。由设计人员自行编程来把一个数字系统“集成”在一片PLD上,而不需要请芯片制造厂商来设计和制作专用的集成电路芯片。而且,如今,取代手工地制作集成电路芯片,这种编程也多半改用“逻辑编译器(logic compiler)”软件来实现,它与程序开发撰写时所用的软件编译器相类似,而要编译之前的原始代码也得用特定的编程语言来撰写,此称之为硬件描述语言(Hardware Description Language,HDL),而HDL也并非仅有一种,而是有许多种,如ABEL(Advanced Boolean Expression Language)、AHDL(Altera Hardware DescriptionLanguage)、Confluence、CUPL(Cornell University Programming Language)、HDCal、JHDL(Java Hardware Description Language)、Lava、Lola、MyHDL、PALASM、RHDL(RubyHardware Description Language)等,目前最普遍使用的是VHDL(Very-High-SpeedIntegrated Circuit Hardware Description Language)与Verilog。本领域技术人员也应该清楚,只需要将方法流程用上述几种硬件描述语言稍作逻辑编程并编程到集成电路中,就可以很容易得到实现该逻辑方法流程的硬件电路。
本领域技术人员也知道,除了以纯计算机可读程序代码方式实现控制器以外,完全可以通过将方法步骤进行逻辑编程来使得控制器以逻辑门、开关、专用集成电路、可编程逻辑控制器和嵌入微控制器等的形式来实现相同功能。因此这种控制器可以被认为是一种硬件部件,而对其内包括的用于实现各种功能的装置也可以视为硬件部件内的结构。或者甚至,可以将用于实现各种功能的装置视为既可以是实现方法的软件模块又可以是硬件部件内的结构。
上述实施例阐明的装置、模块,具体可以由计算机芯片或实体实现,或者由具有某种功能的产品来实现。
为了描述的方便,描述以上装置时以功能分为各种模块分别描述。当然,在实施本申请时可以把各模块的功能在同一个或多个软件和/或硬件中实现。
通过以上的实施方式的描述可知,本领域的技术人员可以清楚地了解到本申请可借助软件加必需的通用硬件平台的方式来实现。基于这样的理解,本申请的技术方案本质上或者说对现有技术做出贡献的部分可以以软件产品的形式体现出来,在一个典型的配置中,计算设备包括一个或多个处理器(CPU)、输入/输出接口、网络接口和内存。该计算机软件产品可以包括若干指令用以使得一台计算机设备(可以是个人计算机,服务器,或者网络设备等)执行本申请各个实施例或者实施例的某些部分所述的方法。该计算机软件产品可以存储在内存中,内存可能包括计算机可读介质中的非永久性存储器,随机存取存储器(RAM)和/或非易失性内存等形式,如只读存储器(ROM)或闪存(flash RAM)。内存是计算机可读介质的示例。计算机可读介质包括永久性和非永久性、可移动和非可移动媒体可以由任何方法或技术来实现信息存储。信息可以是计算机可读指令、数据结构、程序的模块或其他数据。计算机的存储介质的例子包括,但不限于相变内存(PRAM)、静态随机存取存储器(SRAM)、动态随机存取存储器(DRAM)、其他类型的随机存取存储器(RAM)、只读存储器(ROM)、电可擦除可编程只读存储器(EEPROM)、快闪记忆体或其他内存技术、只读光盘只读存储器(CD-ROM)、数字多功能光盘(DVD)或其他光学存储、磁盒式磁带,磁带磁磁盘存储或其他磁性存储设备或任何其他非传输介质,可用于存储可以被计算设备访问的信息。按照本文中的界定,计算机可读介质不包括短暂电脑可读媒体(transitory media),如调制的数据信号和载波。
本说明书中的各个实施例均采用递进的方式描述,各个实施例之间相同相似的部分互相参见即可,每个实施例重点说明的都是与其他实施例的不同之处。尤其,对于装置实施例而言,由于其基本相似于方法实施例,所以描述的比较简单,相关之处参见方法实施例的部分说明即可。
本申请可用于众多通用或专用的计算机系统环境或配置中。例如:个人计算机、服务器计算机、手持设备或便携式设备、平板型设备、多处理器系统、基于微处理器的系统、置顶盒、可编程的消费电子设备、网络PC、小型计算机、大型计算机、包括以上任何系统或设备的分布式计算环境等等。
本申请可以在由计算机执行的计算机可执行指令的一般上下文中描述,例如程序模块。一般地,程序模块包括执行特定任务或实现特定抽象数据类型的例程、程序、对象、组件、数据结构等等。也可以在分布式计算环境中实践本申请,在这些分布式计算环境中,由通过通信网络而被连接的远程处理设备来执行任务。在分布式计算环境中,程序模块可以位于包括存储设备在内的本地和远程计算机存储介质中。
虽然通过实施例描绘了本申请,本领域普通技术人员知道,本申请有许多变形和变化而不脱离本申请的精神,希望所附的权利要求包括这些变形和变化而不脱离本申请的精神。

Claims (9)

1.一种核磁共振数据反演方法,其特征在于,包括:
获取核磁共振回波数据;
基于向量的Lp范数,建立目标函数;其中,所述向量满足预设的约束方程,所述约束方程通过所述核磁共振回波数据确定,所述Lp范数中的参数p的取值范围为大于等于0且小于等于1;
根据所述目标函数的解,生成所述核磁共振回波数据的反演谱;其中,
所述目标函数为:
Figure FDA0002501870310000011
其中,
Figure FDA0002501870310000012
为向量的Lp范数;向量s为所述目标函数的解,向量s中的所有元素都大于等于0;p为Lp范数的参数;
所述约束方程为:
Ks=y
其中,向量y为不同时刻的回波幅度组成的向量,所述不同时刻的回波幅度通过所述核磁共振回波数据得到;K为核矩阵,所述核矩阵K中的元素为:
Figure FDA0002501870310000013
其中,t为回波幅度衰减时间,T2为横向弛豫时间。
2.根据权利要求1所述的方法,其特征在于,利用以下步骤得到所述目标函数的解:
设置参数ε、阈值tol、Lp范数中的参数p、最大迭代次数itmax,并将迭代次数it的初始值设置为1;所述参数ε和所述阈值tol的取值范围为大于0;
计算向量s(0)=K\y,并确定向量s(0)中的元素个数n;
计算初始残差
Figure FDA0002501870310000014
进行第it次迭代计算,所述第it次迭代计算包括:计算第it次迭代后的解s(it),以及第it次迭代后的残差r(it);所述第it次迭代后的解s(it),以及第it次迭代后的残差r(it)通过如下方法得到:
s(it)=WKT/(KWKT)y
Figure FDA0002501870310000021
其中,对角矩阵W=diag(1./w),向量w=[w1,w2,…wi,…wn],
Figure FDA0002501870310000022
在第1次迭代时,
Figure FDA0002501870310000023
为所述向量s(0)中的第i个元素;在迭代次数it的值大于等于2时,
Figure FDA0002501870310000024
为前一次迭代后的解s(it-1)中的第i个元素;
在|r(it)-r|<tol;或,|r(it)-r|≥tol且it=itmax的情况下,令s=s(it),得到所述向量s。
3.根据权利要求2所述的方法,其特征在于,还包括:
在|r(it)-r|≥tol且it<itmax的情况下,令s=s(it)、r=r(it)
将迭代次数it增加1后,进行第it次迭代计算。
4.根据权利要求2所述的方法,其特征在于,所述令s=s(it),得到向量s,包括:
令s=s(it)
若向量s中含有小于0的元素,将所述小于0的元素设置为0,得到向量s。
5.一种核磁共振数据反演装置,其特征在于,包括:
回波数据获取模块,用于获取核磁共振回波数据;
目标函数建立模块,用于基于向量的Lp范数,建立目标函数;其中,所述向量满足预设的约束方程,所述约束方程通过所述核磁共振回波数据确定,所述Lp范数中的参数p的取值范围为大于等于0且小于等于1;
反演谱生成模块,用于根据所述目标函数的解,生成所述核磁共振回波数据的反演谱;其中,
所述目标函数为:
Figure FDA0002501870310000025
其中,
Figure FDA0002501870310000026
为向量的Lp范数;向量s为所述目标函数的解,向量s中的所有元素都大于等于0;p为Lp范数的参数;
所述约束方程为:
Ks=y
其中,向量y为不同时刻的回波幅度组成的向量,所述不同时刻的回波幅度通过所述核磁共振回波数据得到;K为核矩阵,所述核矩阵K中的元素为:
Figure FDA0002501870310000031
其中,t为回波幅度衰减时间,T2为横向弛豫时间。
6.根据权利要求5所述的装置,其特征在于,所述反演谱生成模块,包括:
参数设置单元,用于设置参数ε、阈值tol、Lp范数中的参数p、最大迭代次数itmax,并将迭代次数it的初始值设置为1;所述参数ε和所述阈值tol的取值范围为大于0;
初始解计算单元,用于计算向量s(0)=K\y,并确定向量s(0)中的元素个数n;
初始残差计算单元,用于计算初始残差
Figure FDA0002501870310000032
迭代计算单元,用于进行第it次迭代计算,所述第it次迭代计算包括:计算第it次迭代后的解s(it),以及第it次迭代后的残差r(it);所述第it次迭代后的解s(it),以及第it次迭代后的残差r(it)通过如下方法得到:
s(it)=WKT/(KWKT)y
Figure FDA0002501870310000033
其中,对角矩阵W=diag(1./w),向量w=[w1,w2,…wi,…wn],
Figure FDA0002501870310000034
在第1次迭代时,
Figure FDA0002501870310000035
为所述向量s(0)中的第i个元素;在迭代次数it的值大于等于2时,
Figure FDA0002501870310000036
为前一次迭代后的解s(it-1)中的第i个元素;
第一确定单元,用于在|r(it)-r|<tol;或,|r(it)-r|≥tol且it=itmax的情况下,令s=s(it),得到所述向量s。
7.根据权利要求6所述的装置,其特征在于,所述反演谱生成模块还包括:
第二确定单元,用于在|r(it)-r|≥tol且it<itmax的情况下,令s=s(it)、r=r(it);并将迭代次数it增加1后,进行所述第it次迭代计算。
8.一种计算机设备,包括处理器以及用于存储处理器可执行指令的存储器,所述处理器执行所述指令时实现权利要求1-4中任意一项所述方法的步骤。
9.一种计算机可读存储介质,其上存储有计算机指令,所述指令被执行时实现权利要求1-4中任意一项所述方法的步骤。
CN201910822315.7A 2019-09-02 2019-09-02 一种核磁共振数据反演方法、装置、存储介质及设备 Active CN110618463B (zh)

Priority Applications (1)

Application Number Priority Date Filing Date Title
CN201910822315.7A CN110618463B (zh) 2019-09-02 2019-09-02 一种核磁共振数据反演方法、装置、存储介质及设备

Applications Claiming Priority (1)

Application Number Priority Date Filing Date Title
CN201910822315.7A CN110618463B (zh) 2019-09-02 2019-09-02 一种核磁共振数据反演方法、装置、存储介质及设备

Publications (2)

Publication Number Publication Date
CN110618463A CN110618463A (zh) 2019-12-27
CN110618463B true CN110618463B (zh) 2020-07-28

Family

ID=68922935

Family Applications (1)

Application Number Title Priority Date Filing Date
CN201910822315.7A Active CN110618463B (zh) 2019-09-02 2019-09-02 一种核磁共振数据反演方法、装置、存储介质及设备

Country Status (1)

Country Link
CN (1) CN110618463B (zh)

Families Citing this family (2)

* Cited by examiner, † Cited by third party
Publication number Priority date Publication date Assignee Title
CN113743682B (zh) * 2021-11-03 2022-02-18 中国科学院精密测量科学与技术创新研究院 一种基于有监督深度神经网络的nmr弛豫时间反演方法
CN116012264B (zh) * 2023-03-27 2023-06-13 山东省工业技术研究院 一种基于稀疏约束的图像恢复方法

Citations (4)

* Cited by examiner, † Cited by third party
Publication number Priority date Publication date Assignee Title
CN105891249A (zh) * 2016-03-29 2016-08-24 上海理工大学 时域核磁共振谱反演的方法
CN106772645A (zh) * 2016-12-15 2017-05-31 中国石油大学(北京) 基于一般先验信息约束的核磁共振数据反演方法和装置
CN108009125A (zh) * 2017-12-15 2018-05-08 中国石油大学(北京) 基于l0正则化的核磁共振回波数据反演方法及装置
CN108038082A (zh) * 2017-11-29 2018-05-15 中国石油大学(北京) 基于TwIST的二维核磁共振谱的反演方法和装置

Family Cites Families (1)

* Cited by examiner, † Cited by third party
Publication number Priority date Publication date Assignee Title
US10359532B2 (en) * 2014-12-10 2019-07-23 Schlumberger Technology Corporation Methods to characterize formation properties

Patent Citations (4)

* Cited by examiner, † Cited by third party
Publication number Priority date Publication date Assignee Title
CN105891249A (zh) * 2016-03-29 2016-08-24 上海理工大学 时域核磁共振谱反演的方法
CN106772645A (zh) * 2016-12-15 2017-05-31 中国石油大学(北京) 基于一般先验信息约束的核磁共振数据反演方法和装置
CN108038082A (zh) * 2017-11-29 2018-05-15 中国石油大学(北京) 基于TwIST的二维核磁共振谱的反演方法和装置
CN108009125A (zh) * 2017-12-15 2018-05-08 中国石油大学(北京) 基于l0正则化的核磁共振回波数据反演方法及装置

Non-Patent Citations (1)

* Cited by examiner, † Cited by third party
Title
基于先验信息约束的核磁共振数据反演新方法;金国文 等;《中国石油大学学报(自然科学版)》;20190430;第43卷(第2期);第53-60页 *

Also Published As

Publication number Publication date
CN110618463A (zh) 2019-12-27

Similar Documents

Publication Publication Date Title
Jia et al. Intelligent interpolation by Monte Carlo machine learning
Leistedt et al. Exact wavelets on the ball
Huder et al. COV-OBS. x2: 180 years of geomagnetic field evolution from ground-based and satellite observations
Venema et al. A stochastic iterative amplitude adjusted Fourier transform algorithm with improved accuracy
CN110618463B (zh) 一种核磁共振数据反演方法、装置、存储介质及设备
Albert et al. Calculation of last closed drift shells for the 2013 GEM radiation belt challenge events
CN109612900B (zh) 一种储层岩石渗透率预测方法、装置及存储介质
CN108038149B (zh) 一种温度场数据重构方法
CN110879412A (zh) 地下横波速度反演方法、装置、计算设备及存储介质
Alamaniotis et al. Kernel-based machine learning for background estimation of NaI low-count gamma-ray spectra
CN110794458B (zh) 基于时频分析的含气性检测方法、装置及存储介质
Qian et al. Unsupervised erratic seismic noise attenuation with robust deep convolutional autoencoders
Morzfeld et al. Variational particle smoothers and their localization
Lan et al. Robust high-dimensional seismic data interpolation based on elastic half norm regularization and tensor dictionary learning
CN110686610B (zh) 基于自适应网格的光学变形测量方法及电子设备
Mahadevan et al. Metrics for Intercomparison of Remapping Algorithms (MIRA) protocol applied to Earth system models
Wang et al. Denoising with weak signal preservation by group-sparsity transform learning
CN114444393A (zh) 基于时间卷积神经网络的测井曲线构建方法及装置
CN104391325A (zh) 不连续非均质地质体检测方法和装置
CN114049530A (zh) 混合精度神经网络量化方法、装置及设备
Nikakhtar et al. Optimal transport reconstruction of biased tracers in redshift space
Arshakian et al. Wavelet-based cross-correlation analysis of structure scaling in turbulent clouds
Shaw et al. Multidimensional method-of-lines transport for atmospheric flows over steep terrain using arbitrary meshes
Revunova Randomization approach to the reconstruction of signals resulted from indirect measurements
Conley et al. Bootstrapping the long run

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