CN107861918A - 基于m‑稀疏算法的核磁共振回波数据反演方法与装置 - Google Patents

基于m‑稀疏算法的核磁共振回波数据反演方法与装置 Download PDF

Info

Publication number
CN107861918A
CN107861918A CN201711260011.3A CN201711260011A CN107861918A CN 107861918 A CN107861918 A CN 107861918A CN 201711260011 A CN201711260011 A CN 201711260011A CN 107861918 A CN107861918 A CN 107861918A
Authority
CN
China
Prior art keywords
object function
magnetic resonance
nuclear magnetic
echo data
formula
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
Application number
CN201711260011.3A
Other languages
English (en)
Other versions
CN107861918B (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 CN201711260011.3A priority Critical patent/CN107861918B/zh
Publication of CN107861918A publication Critical patent/CN107861918A/zh
Application granted granted Critical
Publication of CN107861918B publication Critical patent/CN107861918B/zh
Active legal-status Critical Current
Anticipated expiration legal-status Critical

Links

Classifications

    • 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/15—Correlation function computation including computation of convolution operations
    • 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/16—Matrix or vector computation, e.g. matrix-matrix or matrix-vector multiplication, matrix factorization

Landscapes

  • Engineering & Computer Science (AREA)
  • Physics & Mathematics (AREA)
  • General Physics & Mathematics (AREA)
  • Mathematical Physics (AREA)
  • Pure & Applied Mathematics (AREA)
  • Mathematical Analysis (AREA)
  • Mathematical Optimization (AREA)
  • Computational Mathematics (AREA)
  • Data Mining & Analysis (AREA)
  • Theoretical Computer Science (AREA)
  • Computing Systems (AREA)
  • Algebra (AREA)
  • Databases & Information Systems (AREA)
  • Software Systems (AREA)
  • General Engineering & Computer Science (AREA)
  • Magnetic Resonance Imaging Apparatus (AREA)

Abstract

本发明实施例提供一种基于M‑稀疏算法的核磁共振回波数据反演方法与装置。该方法包括：基于反演核矩阵和采集的核磁共振回波数据构建目标函数，其中，所述目标函数的解中非零元素的个数小于预设阈值，所述目标函数中不包括正则化项；使用M‑稀疏算法求解所述目标函数，获取所述核磁共振回波数据的反演谱。本实施例构建的目标函数中不包括正则化项，这样在反演过程中不需要求解正则化参数，简化了反演流程，提高了反演的速度。

Description

基于M-稀疏算法的核磁共振回波数据反演方法与装置
技术领域
本发明实施例涉及核磁共振技术，尤其涉及一种基于M-稀疏算法的核磁共振回波数据反演方法与装置。
背景技术
核磁共振测井是核磁共振技术在石油工业油气勘探与开发中的拓展与应用，该技术的物理基础是利用氢原子核自身的磁性及其与外加磁场的相互作用。具体是通过测量地层孔隙流体中氢核的核磁共振弛豫性质来探测地层孔隙特性和识别流体性质，已经成为世界各大油气田勘探开发中的一项重要测量技术。
核磁共振测井采集的原始数据是由成百上千个自旋回波组成的弛豫衰减曲线，需要通过反演才能得到反映地层信息的核磁共振谱。CPMG(Carr-Purcell-Meiboom-Gill)脉冲序列是核磁共振测井中最常用的测量一维T2谱的序列之一。研究人员通常利用CPMG脉冲序列测量的原始回波数据和横向弛豫时间的分布特点，进行核磁共振数据反演的相关研究。
目前，国内外研究人员使用的核磁共振回波数据反演方法，以增加罚函数的正则化方法为主。正则化方法主要对解的某些性质进行限制，以得到合理的反演结果。但是，正则化方法的关键在于如何选择准确的正则化参数，进而使得反演过程复杂，且不准确的正则化参数会导致反演结果差。
发明内容
本发明实施例提供一种基于M-稀疏算法的核磁共振回波数据反演方法与装置，以解决现有的方法反演过程复杂，且反演结果差的问题。
第一方面，本发明实施例提供一种基于M-稀疏算法的核磁共振回波数据反演方法，包括：
基于M-稀疏算法的核磁共振回波数据反演方法，其特征在于，包括：
基于反演核矩阵和采集的核磁共振回波数据构建目标函数，其中，所述目标函数的解中非零元素的个数小于预设阈值，所述目标函数中不包括正则化项；
使用M-稀疏算法求解所述目标函数，获取所述核磁共振回波数据的反演谱。
第一方面的一种可能的实现方式，所述基于所述反演核矩阵构建目标函数，具体包括：
构建如下式所述的目标函数：
其中，所述f≥0表示向量f中的所有元素均大于或等于0，||f||0表示向量中非零元素的个数，所述N为预设的常数。
第一方面的一种可能的实现方式，所述使用M-稀疏算法求解所述目标函数，获取所述核磁共振回波数据的反演谱，具体包括：
A、根据公式和公式确定所述Γi+1；
B、判断所述Γi+1与Γi是否相等；
C、若是，则令执行D；
D、计算并判断ξ＜tol是否成立；若是则执行F，若否则执行E；
E、判断当前迭代次数是否大于预设的最大迭代次数；若是则执行F，若否令i＝i+1，返回执行A；
F、若所述fi+1中有负值的项，将所述fi+1中负值的项置0，生成新的fi+1，并将所述新的fi+1确定为所述目标函数的目标解；
其中，所述i为迭代次数，所述gi＝AT(b-Afi)，所述所述为向量g只保留对应Γi中元素位置的子向量，所述是矩阵A只保留对应所述Γi中元素位置的列的子矩阵，所述所述是向量|k|中每个元素按从大到小排列后的第M个元素，所述tol为预设阈值，f1＝0，Γ1＝supp(HM(ATb))，所述supp(p)代表向量p的支撑集。
第一方面的一种可能的实现方式，所述判断所述Γi+1与Γi是否相等之后，所述方法还包括：
G、当所述Γi+1与所述Γi不相等时，判断是否成立；若是，则令执行D；若否，则执行H；
H、根据公式μi＝μi/(k(1-c))更新所述μi，根据公式更新所述直到所述成立为止；
I、令和返回执行所述D；
其中，所述k为常数，且0＜k＜1，所述c为常数。
可选的，所述c为0.01。
第二方面，本发明实施例还提供一种基于M-稀疏算法的核磁共振回波数据反演装置，包括：
构建模块，用于基于反演核矩阵和采集的核磁共振回波数据构建目标函数，其中，所述目标函数的解中非零元素的个数小于预设阈值，所述目标函数中不包括正则化项；
获取模块，用于使用M-稀疏算法求解所述目标函数，获取所述核磁共振回波数据的反演谱。
第二方面的一种可能的实现方式，所述构建模块，具体用于构建如下式所述的目标函数：
其中，所述f≥0表示向量f中的所有元素均大于或等于0，||f||0表示向量中非零元素的个数，所述N为预设的常数。
第二方面的一种可能的实现方式，所述获取模块包括确定单元、判断单元和执行单元，其中，
所述确定单元，用于根据公式和公式确定所述Γi+1；
所述判断单元，用于判断所述Γi+1与Γi是否相等；
所述执行单元，用于当所述判断单元判断所述Γi+1与Γi相等时，令且计算
所述判断单元，还用于判断ξ＜tol是否成立；
所述执行单元，还用于当所述判断单元判断ξ＜tol时，将所述fi+1中负值的项置0，生成新的fi+1，并将所述新的fi+1确定为所述目标函数的目标解；当所述判断单元判断ξ＜tol不成立时，使得所述判断单元判断当前迭代次数是否大于预设的最大迭代次数，若是则将所述fi+1中负值的项置0，生成新的fi+1，并将所述新的fi+1确定为所述目标函数的目标解，若否则令i＝i+1，且使得所述确定模块根据公式和公式确定所述Γi+1。
其中，所述i为迭代次数，所述gi＝AT(b-Afi)，所述所述为向量g只保留对应Γi中元素位置的子向量，所述是矩阵A只保留对应所述Γi中元素位置的列的子矩阵，所述所述是向量|k|中每个元素按从大到小排列后的第M个元素，所述tol为预设阈值，f1＝0，Γ1＝supp(HM(ATb))，所述supp(p)代表向量p的支撑集。
第二方面的一种可能的实现方式，所述判断单元，还用于在判断Γi+1与所述Γi不相等时，判断是否成立；
所述执行单元，还用于在所述判断模块判断成立时，计算
所述执行单元，还用于在所述判断模块判断不成立时，根据公式μi＝μi/(k(1-c))更新所述μi，根据公式更新所述直到所述成立为止，此时令和计算
本发明实施例提供的基于M-稀疏算法的核磁共振回波数据反演方法与装置，基于反演核矩阵和采集的核磁共振回波数据构建目标函数，其中，目标函数的解中非零元素的个数小于预设阈值，目标函数中不包括正则化项；接着使用M-稀疏算法求解所述目标函数，获取所述核磁共振回波数据的反演谱。即本实施例构建的目标函数中不包括正则化项，这样在反演过程中不需要求解正则化参数，简化了反演流程，提高了反演的速度。
附图说明
为了更清楚地说明本发明实施例或现有技术中的技术方案，下面将对实施例或现有技术描述中所需要使用的附图作一简单地介绍，显而易见地，下面描述中的附图是本发明的一些实施例，对于本领域普通技术人员来讲，在不付出创造性劳动性的前提下，还可以根据这些附图获得其他的附图。
图1为本发明提供的基于M-稀疏算法的核磁共振回波数据反演方法实施例一的流程图；
图2为本发明提供的基于M-稀疏算法的核磁共振回波数据反演方法实施例二的流程图；
图3是数值模拟实验构造的含有双峰的T2谱示意图；
图4(a)中回波间隔为0.2ms、信噪比SNR为10的回波串数据示意图；
图4(b)中回波间隔为0.2ms、信噪比SNR为20的回波串数据示意图；
图4(c)中回波间隔为0.2ms、信噪比SNR为40的回波串数据示意图；
图5(a)是三种不同反演方法分别处理图4(a)中的数据得到的核磁共振T2谱与模型的对比结果图；
图5(b)是三种不同反演方法分别处理图4(b)中的数据得到的核磁共振T2谱与模型的对比结果图；
图5(c)是三种不同反演方法分别处理图4(c)中的数据得到的核磁共振T2谱与模型的对比结果图；
图6为本发明提供的基于M-稀疏算法的核磁共振回波数据反演装置实施例一的示意图；
图7为本发明提供的基于M-稀疏算法的核磁共振回波数据反演装置实施例二的示意图。
具体实施方式
为使本发明实施例的目的、技术方案和优点更加清楚，下面将结合本发明实施例中的附图，对本发明实施例中的技术方案进行清楚、完整地描述，显然，所描述的实施例是本发明一部分实施例，而不是全部的实施例。基于本发明中的实施例，本领域普通技术人员在没有作出创造性劳动前提下所获得的所有其他实施例，都属于本发明保护的范围。
本发明实施例提供的基于M-稀疏算法的核磁共振回波数据反演方法，可以应用于石油、天然气等的核磁共振测井数据处理过程。具体通过核磁共振采样参数确定反演核矩阵，并根据反演核矩阵构建目标函数，由于该目标函数不包括正则化项。这样在反演过程中不需要求解正则化参数，简化了反演流程，计算精度高，且在不同信噪比数据情况下，能得到稳定的反演结果。
下面以具体地实施例对本发明的技术方案进行详细说明。下面这几个具体的实施例可以相互结合，对于相同或相似的概念或过程可能在某些实施例不再赘述。
图1为本发明提供的基于M-稀疏算法的核磁共振回波数据反演方法实施例一的流程图，如图1所示，本实施例的方法可以包括：
S101、基于反演核矩阵和采集的核磁共振回波数据构建目标函数，其中，所述目标函数的解中非零元素的个数小于预设阈值，所述目标函数中不包括正则化项。
本实施例的执行主体是具有数据处理能力的终端设备，例如计算机、智能手机、上位机或测试机等。
首先，获取核磁共振回波数据，核磁共振回波数据为施加了噪声的回波数据，通常可以表示为一个向量，该向量的维度即核磁共振回波数据的个数。
由于一维核磁共振CPMG脉冲序列测量数据的反演问题，实质为求解如公式(1)所示的积分方程：
其中，t为所述核磁共振回波数据的采样时刻，b(t)表示t采样时刻的回波幅度，f(T2)表示核磁共振T2谱。
上述公式(1)可以用矩阵形式来表示，如公式(2)所示：
b＝Af，(2)；
其中，b为b(t)的向量，f为f(T2)的向量，A是反演核矩阵。
接着，基于上述反演核矩阵和采集的核磁共振回波数据构建目标函数。
其中，本实施例在构建目标函数时，使用非零元素个数来约束目标函数的解，此时构建的目标函数不含正则化项。这样，本实施例的方法与现有的正则化方法相比，在求解目标函数的解时，不需要求解正则化参数，进而简化了反演流程，提高了反演速度。
基于上述步骤可知，当使用矩阵b＝Af表示一维核磁共振CPMG脉冲序列测量数据的反演问题时，本实施例的目标函数可以如(3a)所示：
可选的，本实施例的目标函数还可以如(3b)所示：
可选的，本实施例的目标函数还可以如(3c)所示：
其中，所述k为常数，所述f≥0表示向量f中的所有元素均大于或等于0，||f||0表示向量f中非零元素的个数，所述N为预设的常数。
S102、使用M-稀疏算法求解所述目标函数，获取所述核磁共振回波数据的反演谱。
储层岩石通常包含多种不同大小的孔隙，并且常常含有多种流体成分，此时孔隙中存在多种弛豫组分，即横向弛豫时间(T2)不是单值，而是一个分布的T2谱。
基于此，本实施例根据上述步骤构建目标函数，对该目标函数进行求解，并将目标函数的解作为核磁共振回波数据的反演谱。本实施例构建的目标函数中不包括正则化项，这样在反演过程中不需要求解正则化参数，简化了反演流程，提高了反演的速度。
本实施例可以使用迭代法对目标函数进行求解，获得目标函数的解。具体是，使用M-稀疏算法求解目标函数。
获得目标函数的解之后，结合布点好的弛豫时间，进而可以输出核磁共振回波数据的反演谱。
本发明实施例提供的基于M-稀疏算法的核磁共振回波数据反演方法，基于反演核矩阵和采集的核磁共振回波数据构建目标函数，其中，目标函数的解中非零元素的个数小于预设阈值，目标函数中不包括正则化项；接着使用M-稀疏算法求解所述目标函数，获取所述核磁共振回波数据的反演谱。即本实施例构建的目标函数中不包括正则化项，这样在反演过程中不需要求解正则化参数，简化了反演流程，提高了反演的速度。
图2为本发明提供的基于M-稀疏算法的核磁共振回波数据反演方法实施例二的流程图，在上述实施例的基础上，本实施例涉及的是使用M-稀疏算法求解所述目标函数，获取所述核磁共振回波数据的反演谱的具体过程。上述S102具体可以包括：
S201、根据公式和公式确定所述Γi+1。
本实施例使用M-稀疏算法构建目标函数解f的迭代函数的具体过程为：
设定gi＝AT(b-Afi)，(4)
接着，设定
其中，i为迭代次数，为向量g只保留对应Γi中元素位置的子向量，是矩阵A只保留对应所述Γi中元素位置的列的子矩阵， 是向量|k|中每个元素按从大到小排列后的第M个元素。
然后，根据公式(4)和公式(5)，确定如公式(6)和公式(7)所示的迭代函数：
其中，supp(p)代表向量p的支撑集。
在进行迭代之前，首先需要设定迭代函数的初始值，例如，设定f1＝0，Γ1＝supp(HM(ATb))。
S202、判断所述Γi+1与Γi是否相等。
若是，则令执行S203。
根据上述设定的初始条件和迭代函数，进行迭代，获得Γi+1后，将Γi+1与之前的Γi进行比较，判断Γi+1与Γi是否相等。例如，第一次迭代时，i＝1，根据上述迭代公式获得的Γ2，接着，判断Γ2与Γ1是否相同，若相同，则执行下述S203。
若Γi+1与Γi不相等时，可以继续进行迭代，也可以对迭代函数的参数gi和参数μi进行更新，进而更新目标函数，并使用更新后的目标函数进行再次迭代。
S203、计算并判断ξ＜tol是否成立。
若是则执行S205，若否则执行S204。
其中，tol为预设阈值。
S204、判断当前迭代次数是否大于预设的最大迭代次数。
若是则执行S205，若否则令i＝i+1，并返回执行S201。
其中预设的最大迭代次数根据实际需要进行设定，实施例对此不做限制。
S205、若所述fi+1中有负值的项，将所述fi+1中负值的项置0，生成新的fi+1，并将所述新的fi+1确定为所述目标函数的目标解。
即本实施例，经过上述迭代，得到最优的f，最后，利用T2和f两组数据输出核磁共振T2谱。
继续参照图2所示，在本实施例的一种可能的实现方式中，为了进一步提高目标函数解的精确性，当S202在判断Γi+1与Γi不相等时，本实施例的方法还可以包括：
S206、判断是否成立。
若是，则令执行S203；若否，则执行S207。
S207、根据公式μi＝μi/(k(1-c))更新所述μi，根据公式更新所述直到所述成立为止。
S208、令和
并返回执行所述S203。
其中，所述k为常数，且0＜k＜1，所述c为常数。
可选的，所述c为0.01。
本发明实施例提供的基于M-稀疏算法的核磁共振回波数据反演方法，通过M-稀疏算法对目标函数进行求解，可以获得最优解，进而基于该最优解获得的核磁共振T2谱更加的精确，进而提高了核磁共振回波数据反演的准确性。
下面以构造出的核磁共振回波数据的目标函数为例，对本发明提供的核磁共振回波数据反演方法进行详细的说明。
图3是数值模拟实验构造的含有双峰的T2谱示意图，两个峰对应的T2值分别为5ms和150ms。
接着，向正演结果中加入一定的高斯白噪声，模拟得到不同信噪比(SNR)的回波串数据，如图4(a)、图4(b)和图4(c)所示。其中，图4(a)中回波间隔为0.2ms、SNR为10的回波串数据示意图。图4(b)中回波间隔为0.2ms、SNR为20的回波串数据示意图。图4(c)中回波间隔为0.2ms、SNR为40的回波串数据示意图。
然后，使用三种不同反演方法，即截断奇异值分解方法TSVD、Tikhonov正则化方法，以及本申请的M-稀疏算法，分别处理图4(a)、图4(b)和图4(c)所示的不同信噪比数据，得到的核磁共振T2谱与模型的对比结果如图5(a)、图5(b)和图5(c)所示。
图5(a)是三种不同反演方法分别处理图4(a)中的数据得到的核磁共振T2谱与模型的对比结果图。图5(b)是三种不同反演方法分别处理图4(b)中的数据得到的核磁共振T2谱与模型的对比结果图。图5(c)是三种不同反演方法分别处理图4(c)中的数据得到的核磁共振T2谱与模型的对比结果图。
由图5(a)、图5(b)和图5(c)可知，本申请的反演结果优于目前常用的TSVD方法以及Tikhonov正则化方法的反演结果，尤其在信噪比很低时，本申请提出的方法反演得到的T2谱分辨率更高。且随着信噪比的增高，本申请的反演结果准确性越来越高。
为了进一步验证本发明方法的优越性，表1列举了信噪比分别为10、20和40时的三种方法反演结果与模型的相对误差，相对误差的计算公式如下：
其中f是反演结果，fmodel是构造的模型。
由表1可知，相比于TSVD方法以及Tikhonov正则化方法，本发明的方法反演的结果与模型的相对误差更小，反演结果更贴近真实的模拟数据。
表1
本发明实施例提供的基于M-稀疏算法的核磁共振回波数据反演方法，通过M-稀疏算法对目标函数进行求解，可以获得目标函数的精确解，进而基于该精确解获得的核磁共振T2谱更加的精确，进而提高了核磁共振回波数据反演的准确性。并且本实施例的方法，在信噪比很低时，可以得到分辨率高的T2谱。
图6为本发明提供的基于M-稀疏算法的核磁共振回波数据反演装置实施例一的示意图。如图6所示，本实施例的核磁共振回波数据反演装置100可以包括：构建模块110和获取模块120，其中，
构建模块110，用于基于反演核矩阵和采集的核磁共振回波数据构建目标函数，其中，所述目标函数的解中非零元素的个数小于预设阈值，所述目标函数中不包括正则化项；
获取模块120，用于使用M-稀疏算法求解所述目标函数，获取所述核磁共振回波数据的反演谱。
本实施例的装置，可以用于执行上述方法实施例的技术方案，其实现原理和技术效果类似，此处不再赘述。
在本实施例的一种可能的实现方式中，上述构建模块110，具体用于构建目标函数：
其中，所述f≥0表示向量f中的所有元素均大于或等于0，||f||0表示向量中非零元素的个数，所述N为预设的常数。
图7为本发明提供的基于M-稀疏算法的核磁共振回波数据反演装置实施例二的示意图。在上述实施例的基础上，如图7所示，本实施例的获取模块120包括确定单元121、判断单元122和执行单元123，其中，
确定单元121，用于根据公式和公式确定所述Γi+1；
判断单元122，用于判断所述Γi+1与Γi是否相等；
执行单元123，用于当判断单元122判断所述Γi+1与Γi相等时，令且计算
判断单元122，还用于判断ξ＜tol是否成立；
执行单元123，还用于当判断单元122判断ξ＜tol时，将所述fi+1中的负值的项置0，生成新的fi+1，并将所述新的fi+1确定为所述目标函数的目标解；当判断单元122判断ξ＜tol不成立时，使得判断单元122判断当前迭代次数是否大于预设的最大迭代次数，若是则将所述fi+1中负值的项置0，生成新的fi+1，并将所述新的fi+1确定为所述目标函数的目标解，若否则令i＝i+1，且使得所述确定模块根据公式和公式确定所述Γi+1。
其中，所述i为迭代次数，所述gi＝AT(b-Afi)，所述所述为向量g只保留对应Γi中元素位置的子向量，所述是矩阵A只保留对应所述Γi中元素位置的列的子矩阵，所述所述是向量|k|中每个元素按从大到小排列后的第M个元素，所述tol为预设阈值，f1＝0，Γ1＝supp(HM(ATb))，所述supp(p)代表向量p的支撑集。
在本实施例的一种可能的实现方式中，
判断单元122，还用于在判断Γi+1与所述Γi不相等时，判断是否成立；
执行单元123，还用于在所述判断模块122判断成立时，计算
执行单元123，还用于在所述判断模块122判断不成立时，根据公式μi＝μi/(k(1-c))更新所述μi，根据公式更新所述直到所述成立为止，此时令和计算
其中，所述k为常数，且0＜k＜1，所述c为常数。
本实施例的装置，可以用于执行上述所示方法实施例的技术方案，其实现原理和技术效果类似，此处不再赘述。
本领域普通技术人员可以理解：实现上述各方法实施例的全部或部分步骤可以通过程序指令相关的硬件来完成。前述的程序可以存储于一计算机可读取存储介质中。该程序在执行时，执行包括上述各方法实施例的步骤；而前述的存储介质包括：ROM、RAM、磁碟或者光盘等各种可以存储程序代码的介质。
最后应说明的是：以上各实施例仅用以说明本发明的技术方案，而非对其限制；尽管参照前述各实施例对本发明进行了详细的说明，本领域的普通技术人员应当理解：其依然可以对前述各实施例所记载的技术方案进行修改，或者对其中部分或者全部技术特征进行等同替换；而这些修改或者替换，并不使相应技术方案的本质脱离本发明各实施例技术方案的范围。

Claims (10)

1.一种基于M-稀疏算法的核磁共振回波数据反演方法，其特征在于，包括：
基于反演核矩阵和采集的核磁共振回波数据构建目标函数，其中，所述目标函数的解中非零元素的个数小于预设阈值，所述目标函数中不包括正则化项；
使用M-稀疏算法求解所述目标函数，获取所述核磁共振回波数据的反演谱。
2.根据权利要求1所述的方法，其特征在于，所述基于反演核矩阵和采集的核磁共振回波数据构建目标函数，具体包括：
构建如下式所述的目标函数：
<mrow> <munder> <munder> <mi>argmin</mi> <mrow> <mi>f</mi> <mo>&amp;GreaterEqual;</mo> <mn>0</mn> </mrow> </munder> <mrow> <mo>|</mo> <mo>|</mo> <mi>f</mi> <mo>|</mo> <msub> <mo>|</mo> <mn>0</mn> </msub> <mo>&amp;le;</mo> <mi>N</mi> </mrow> </munder> <mo>|</mo> <mo>|</mo> <mi>A</mi> <mi>f</mi> <mo>-</mo> <mi>b</mi> <mo>|</mo> <msubsup> <mo>|</mo> <mn>2</mn> <mn>2</mn> </msubsup> <mo>;</mo> </mrow>
其中，所述f≥0表示向量f中的所有元素均大于或等于0，||f||0表示向量中非零元素的个数，所述N为预设的常数。
3.根据权利要求2所述的方法，其特征在于，所述使用M-稀疏算法求解所述目标函数，获取所述核磁共振回波数据的反演谱，具体包括：
A、根据公式和公式确定所述Γi+1；
B、判断所述Γi+1与Γi是否相等；
C、若是，则令执行D；
D、计算并判断ξ＜tol是否成立；若是则执行F，若否则执行E；
E、判断当前迭代次数是否大于预设的最大迭代次数；若是则执行F，若否令i＝i+1，返回执行A；
F、若所述fi+1中有负值的项，将所述fi+1中负值的项置0，生成新的fi+1，并将所述新的fi +1确定为所述目标函数的目标解；
其中，所述i为迭代次数，所述gi＝AT(b-Afi)，所述所述为向量g只保留对应Γi中元素位置的子向量，所述是矩阵A只保留对应所述Γi中元素位置的列的子矩阵，所述所述是向量|k|中每个元素按从大到小排列后的第M个元素，所述tol为预设阈值，f1＝0，Γ1＝supp(HM(ATb))，所述supp(p)代表向量p的支撑集。
4.根据权利要求3所述的方法，其特征在于，所述判断所述Γi+1与Γi是否相等之后，所述方法还包括：
G、当所述Γi+1与所述Γi不相等时，判断是否成立；若是，则令执行D；若否，则执行H；
H、根据公式μi＝μi/(k(1-c))更新所述μi，根据公式更新所述直到所述成立为止；
I、令和返回执行所述D；
其中，所述k为常数，且0＜k＜1，所述c为常数。
5.根据权利要求4所述的方法，其特征在于，所述c为0.01。
6.一种基于M-稀疏算法的核磁共振回波数据反演装置，其特征在于，包括：
构建模块，用于基于反演核矩阵和采集的核磁共振回波数据构建目标函数，其中，所述目标函数的解中非零元素的个数小于预设阈值，所述目标函数中不包括正则化项；
获取模块，用于使用M-稀疏算法求解所述目标函数，获取所述核磁共振回波数据的反演谱。
7.根据权利要求6所述的装置，其特征在于，所述构建模块，具体用于构建如下式所述的目标函数：
<mrow> <munder> <munder> <mi>argmin</mi> <mrow> <mi>f</mi> <mo>&amp;GreaterEqual;</mo> <mn>0</mn> </mrow> </munder> <mrow> <mo>|</mo> <mo>|</mo> <mi>f</mi> <mo>|</mo> <msub> <mo>|</mo> <mn>0</mn> </msub> <mo>&amp;le;</mo> <mi>N</mi> </mrow> </munder> <mo>|</mo> <mo>|</mo> <mi>A</mi> <mi>f</mi> <mo>-</mo> <mi>b</mi> <mo>|</mo> <msubsup> <mo>|</mo> <mn>2</mn> <mn>2</mn> </msubsup> <mo>;</mo> </mrow>
其中，所述f≥0表示向量f中的所有元素均大于或等于0，||f||0表示向量中非零元素的个数，所述N为预设的常数。
8.根据权利要求7所述的装置，其特征在于，所述获取模块包括确定单元、判断单元和执行单元；
所述确定单元，用于根据公式和公式确定所述Γi+1；
所述判断单元，用于判断所述Γi+1与Γi是否相等；
所述执行单元，用于当判断单元判断所述Γi+1与Γi相等时，令且计算
所述判断单元，还用于判断ξ＜tol是否成立；
所述执行单元，还用于当所述判断单元判断ξ＜tol时，将所述fi+1中负值的项置0，生成新的fi+1，并将所述新的fi+1确定为所述目标函数的目标解；当所述判断单元判断ξ＜tol不成立时，使得所述判断单元判断当前迭代次数是否大于预设的最大迭代次数，若是则将所述fi+1中负值的项置0，生成新的fi+1，并将所述新的fi+1确定为所述目标函数的目标解，若否则令i＝i+1，且使得所述确定模块根据公式和公式确定所述Γi+1。
其中，所述i为迭代次数，所述gi＝AT(b-Afi)，所述所述为向量g只保留对应Γi中元素位置的子向量，所述是矩阵A只保留对应所述Γi中元素位置的列的子矩阵，所述所述是向量|k|中每个元素按从大到小排列后的第M个元素，所述tol为预设阈值，所述supp(p)代表向量p的支撑集。
9.根据权利要求8所述的装置，其特征在于，
所述判断单元，还用于在判断Γi+1与所述Γi不相等时，判断是否成立；
所述执行单元，还用于在所述判断模块判断成立时，计算
所述执行单元，还用于在所述判断模块判断不成立时，根据公式μi＝μi/(k(1-c))更新所述μi，根据公式更新所述直到所述成立为止，并令和计算
其中，所述k为常数，且0＜k＜1，所述c为常数。
10.根据权利要求9所述的装置，其特征在于，所述c为0.01。
CN201711260011.3A 2017-12-04 2017-12-04 基于m-稀疏算法的核磁共振回波数据反演方法与装置 Active CN107861918B (zh)

Priority Applications (1)

Application Number Priority Date Filing Date Title
CN201711260011.3A CN107861918B (zh) 2017-12-04 2017-12-04 基于m-稀疏算法的核磁共振回波数据反演方法与装置

Applications Claiming Priority (1)

Application Number Priority Date Filing Date Title
CN201711260011.3A CN107861918B (zh) 2017-12-04 2017-12-04 基于m-稀疏算法的核磁共振回波数据反演方法与装置

Publications (2)

Publication Number Publication Date
CN107861918A true CN107861918A (zh) 2018-03-30
CN107861918B CN107861918B (zh) 2020-01-03

Family

ID=61704985

Family Applications (1)

Application Number Title Priority Date Filing Date
CN201711260011.3A Active CN107861918B (zh) 2017-12-04 2017-12-04 基于m-稀疏算法的核磁共振回波数据反演方法与装置

Country Status (1)

Country Link
CN (1) CN107861918B (zh)

Cited By (1)

* Cited by examiner, † Cited by third party
Publication number Priority date Publication date Assignee Title
CN109282841A (zh) * 2018-09-19 2019-01-29 湘潭大学 一种无线无源声表面波传感器的超分辨率测量方法

Citations (4)

* Cited by examiner, † Cited by third party
Publication number Priority date Publication date Assignee Title
CN102565865A (zh) * 2011-12-07 2012-07-11 中国石油大学(北京) 降噪核磁共振测井回波信号的获得方法及装置
US20160061986A1 (en) * 2014-08-27 2016-03-03 Schlumberger Technology Corporation Formation Property Characteristic Determination Methods
CN106201991A (zh) * 2016-07-07 2016-12-07 中国石油大学(北京) 核磁共振回波数据反演方法及装置
CN106707353A (zh) * 2016-12-01 2017-05-24 中国科学技术大学 一种低信噪比下考虑基线偏移的核磁共振t2谱反演方法

Patent Citations (4)

* Cited by examiner, † Cited by third party
Publication number Priority date Publication date Assignee Title
CN102565865A (zh) * 2011-12-07 2012-07-11 中国石油大学(北京) 降噪核磁共振测井回波信号的获得方法及装置
US20160061986A1 (en) * 2014-08-27 2016-03-03 Schlumberger Technology Corporation Formation Property Characteristic Determination Methods
CN106201991A (zh) * 2016-07-07 2016-12-07 中国石油大学(北京) 核磁共振回波数据反演方法及装置
CN106707353A (zh) * 2016-12-01 2017-05-24 中国科学技术大学 一种低信噪比下考虑基线偏移的核磁共振t2谱反演方法

Non-Patent Citations (2)

* Cited by examiner, † Cited by third party
Title
GUANQUN SU ET AL.: "An inversion method of 2D NMR relaxation spectra in low fields based on LSQR and L-curve", 《JOURNAL OF MAGNETIC RESONANCE》 *
蒋川东 等: "基于L1范数的低场核磁共振T2谱稀疏反演方法", 《物理学报》 *

Cited By (2)

* Cited by examiner, † Cited by third party
Publication number Priority date Publication date Assignee Title
CN109282841A (zh) * 2018-09-19 2019-01-29 湘潭大学 一种无线无源声表面波传感器的超分辨率测量方法
CN109282841B (zh) * 2018-09-19 2021-03-02 湘潭大学 一种无线无源声表面波传感器的超分辨率测量方法

Also Published As

Publication number Publication date
CN107861918B (zh) 2020-01-03

Similar Documents

Publication Publication Date Title
Mohamed et al. History matching and uncertainty quantification: multiobjective particle swarm optimisation approach
JP5797262B2 (ja) 成熟炭化水素産地をシミュレートするための生産シミュレータ
CN1693896B (zh) 确定地层流体性质的方法
CN103336305B (zh) 一种基于灰色理论划分致密砂岩储层岩石物理相的方法
CN111027882A (zh) 一种基于高阶神经网络利用常规测井资料评价脆性指数的方法
CN114662414B (zh) 一种基于图小波神经网络模型的油藏生产预测方法
Parker et al. Assigning uncertainties in the inversion of NMR relaxation data
CN106201991A (zh) 核磁共振回波数据反演方法及装置
CN112098449B (zh) 一种储层毛管压力曲线的获取方法和装置
CN114154427A (zh) 基于深度学习的体积压裂裂缝扩展预测方法和系统
CN106093083A (zh) 一种采用数字岩心模拟建立岩电关系的方法
CN107808032A (zh) 页岩气储层含水饱和度计算方法及系统
Zhang et al. Reservoir geomechanical parameters identification based on ground surface movements
CN113050189B (zh) 一种测井曲线的重构方法、装置、设备及存储介质
Song et al. Potential for vertical heterogeneity prediction in reservoir basing on machine learning methods
WO2004086280A2 (fr) Methode pour former rapidement un modele stochastique representatif de la distribution d’une grandeur physique dans un milieu heterogene par une selection appropriee de realisations geostatistiques
CN110029990B (zh) 一种核磁共振测井方法和装置
CN108009125A (zh) 基于l0正则化的核磁共振回波数据反演方法及装置
CN110320227A (zh) 一种二维核磁共振d-t2谱反演方法及装置
Zhou et al. Determination of pore size distribution in tight gas sandstones based on Bayesian regularization neural network with MICP, NMR and petrophysical logs
Zhang et al. Efficient surrogate modeling based on improved vision transformer neural network for history matching
CN107861918A (zh) 基于m‑稀疏算法的核磁共振回波数据反演方法与装置
CN115130529A (zh) 一种叠前地震反演方法和装置
Shaochuan Bayesian multiple changepoint detection for stochastic models in continuous time
CN109655394B (zh) 一种孔喉结构参数约束下的核磁共振t2谱渗透率计算方法

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