CN114966861B - 基于Lp伪范数和γ范数稀疏低秩约束的地震去噪方法 - Google Patents

基于Lp伪范数和γ范数稀疏低秩约束的地震去噪方法 Download PDF

Info

Publication number
CN114966861B
CN114966861B CN202210534230.0A CN202210534230A CN114966861B CN 114966861 B CN114966861 B CN 114966861B CN 202210534230 A CN202210534230 A CN 202210534230A CN 114966861 B CN114966861 B CN 114966861B
Authority
CN
China
Prior art keywords
seismic
time
norm
frequency domain
denoising
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
CN202210534230.0A
Other languages
English (en)
Other versions
CN114966861A (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.)
Chengdu Univeristy of Technology
Original Assignee
Chengdu Univeristy of Technology
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 Chengdu Univeristy of Technology filed Critical Chengdu Univeristy of Technology
Priority to CN202210534230.0A priority Critical patent/CN114966861B/zh
Publication of CN114966861A publication Critical patent/CN114966861A/zh
Application granted granted Critical
Publication of CN114966861B publication Critical patent/CN114966861B/zh
Active legal-status Critical Current
Anticipated expiration legal-status Critical

Links

Classifications

    • G—PHYSICS
    • G01—MEASURING; TESTING
    • G01V—GEOPHYSICS; GRAVITATIONAL MEASUREMENTS; DETECTING MASSES OR OBJECTS; TAGS
    • G01V1/00—Seismology; Seismic or acoustic prospecting or detecting
    • G01V1/28—Processing seismic data, e.g. for interpretation or for event detection
    • G01V1/36—Effecting static or dynamic corrections on records, e.g. correcting spread; Correlating seismic signals; Eliminating effects of unwanted energy
    • G—PHYSICS
    • G01—MEASURING; TESTING
    • G01V—GEOPHYSICS; GRAVITATIONAL MEASUREMENTS; DETECTING MASSES OR OBJECTS; TAGS
    • G01V1/00—Seismology; Seismic or acoustic prospecting or detecting
    • G01V1/28—Processing seismic data, e.g. for interpretation or for event detection
    • G01V1/32—Transforming one recording into another or one representation into another
    • 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
    • G—PHYSICS
    • G01—MEASURING; TESTING
    • G01V—GEOPHYSICS; GRAVITATIONAL MEASUREMENTS; DETECTING MASSES OR OBJECTS; TAGS
    • G01V2210/00—Details of seismic processing or analysis
    • G01V2210/30—Noise handling
    • G01V2210/32—Noise reduction

Landscapes

  • Engineering & Computer Science (AREA)
  • Physics & Mathematics (AREA)
  • General Physics & Mathematics (AREA)
  • Mathematical Physics (AREA)
  • Mathematical Analysis (AREA)
  • Remote Sensing (AREA)
  • Life Sciences & Earth Sciences (AREA)
  • Theoretical Computer Science (AREA)
  • Data Mining & Analysis (AREA)
  • Pure & Applied Mathematics (AREA)
  • Computational Mathematics (AREA)
  • Mathematical Optimization (AREA)
  • Environmental & Geological Engineering (AREA)
  • Computing Systems (AREA)
  • Geophysics (AREA)
  • Algebra (AREA)
  • General Life Sciences & Earth Sciences (AREA)
  • Geology (AREA)
  • Databases & Information Systems (AREA)
  • Software Systems (AREA)
  • General Engineering & Computer Science (AREA)
  • Acoustics & Sound (AREA)
  • Geophysics And Detection Of Objects (AREA)

Abstract

本发明提供了一种基于Lp伪范数和γ范数稀疏低秩约束的地震去噪方法，涉及油气勘探及地震资料处理领域，其包括步骤1：对地震资料进行预处理并获取时频域的地震记录步骤2：根据时频域的地震记录构建基于LP伪范数稀疏约束和γ范数低秩约束的目标函数；步骤3：采用交替方向乘子法对目标函数进行迭代计算得到时频域去噪后的地震记录步骤4：将时频域去噪后的地震记录转换到时空域。本发明能够消除地震记录中的随机噪声，并降低地震记录在去噪过程中的损失，提高了地震记录的质量。

Description

基于Lp伪范数和γ范数稀疏低秩约束的地震去噪方法
技术领域
本发明涉及油气勘探及地震资料处理领域，尤其是一种基于Lp伪范数和γ范数稀疏低秩约束的地震去噪方法。
背景技术
野外地震记录会受到随机噪声干扰，这对于解释工作带来了不小的麻烦，因此地震记录去噪是油气勘探中非常重要的过程。基于稀疏低秩理论的地震记录去噪技术是地震去噪的重要方法，原理是利用时频域地震记录的稀疏性质和低秩性质进行去噪处理。
基于L1范数稀疏约束和核范数低秩约束的基于Lp伪范数和γ范数稀疏低秩约束的地震去噪方法是基于稀疏低秩理论地震去噪的重要方法。该方法采用L1范数对时频域的地震记录进行稀疏约束，采用核范数对时频域的地震记录进行低秩约束，提高了地震去噪能力。Anvari等人将L1范数和核范数用于时频域的地震去噪中，证明了该方法的可行性。但是L1范数并不能完全地挖掘稀疏性，核范数也不能精确的表示数据矩阵的秩，导致去噪结果噪音仍有残余且有效信息存在损失，因此需要一种去噪方法克服以上的问题。
发明内容
本发明的目的在于：本发明提供了一种基于Lp伪范数和γ范数稀疏低秩约束的地震去噪方法，解决了现有基于L1范数稀疏约束和核范数低秩约束的基于Lp伪范数和γ范数稀疏低秩约束的地震去噪方法中所造成的去噪结果噪音仍有残余且有效信息存在损失的问题，提高去噪效果。
本发明采用的技术方案如下：
根据本发明的第一方面，提供一种基于Lp伪范数和γ范数稀疏低秩约束的地震去噪方法，所述方法包括：根据地震资料得到时频域的地震记录根据时频域的地震记录/>构建基于LP伪范数稀疏约束和γ范数低秩约束的目标函数；采用交替方向乘子法对目标函数进行迭代计算得到时频域去噪后的地震记录/>将时频域去噪后的地震记录/>转换到时空域。
进一步，所述根据地震资料得到时频域的地震记录包括：输入时空域含噪的地震记录S和子波数据w；利用子波数据w确定短时傅里叶变换的窗口长度，滑动距离，并对时空域含噪的地震记录S进行短时傅里叶变换，得到时频域的地震记录/>
进一步，所述根据时频域的地震记录构建基于LP伪范数稀疏约束和γ范数低秩约束的目标函数，包括：
根据时频域的地震记录构建目标函数中的保真项，如公式(1)所示：
其中，表示时频域含噪的地震记录，/>表示时频域去噪的地震记录，/>表示Frobenius范数；
对时频域去噪的地震记录加上LP伪范数进行稀疏约束，如公式(2)所示：
其中，λ1表示LP伪范数稀疏约束的权重系数，表示LP伪范数，其表达式如公式(3)所示：
其中a、b、m和n分别表示时频地震记录的第a行、第b列、总行数和总列数；
对时频域去噪的地震记录加上γ范数进行低秩约束，如公式(4)所示：
其中，λ2表示γ范数低秩约束的权重系数，||||γ表示γ范数，其表达式如公式(5)所示；
其中，表示时频域去噪的地震记录/>经过奇异值分解后的第k个特征值；
引入拉格朗日乘子项R和对偶项C，得到基于LP伪范数和γ范数稀疏低秩约束的地震去噪目标函数，如公式(6)所示：
其中λ3表示对偶项权重系数。
进一步，所述采用交替方向乘子法对目标函数进行迭代计算得到时频域去噪后的地震记录包括：
分别设置稀疏约束项、低秩约束项和对偶项的权重系数λ1、λ2和λ3的值，以及LP伪范数中参数P的值和γ范数中参数γ的值；
设置拉格朗日乘子项初始值R1＝0及其对偶项初始值C1＝0；
根据交替方向乘子算法以及软阈值收缩算法更新时频域待去噪的地震记录计算公式如公式(7)所示：
其中，表示矩阵点乘运算，i+1表示第i+1次更新，sign表示符号函数，表达式如公式(8)所示：
其中，ya,b表示矩阵中的第a行，第b列元素；
根据交替方向乘子算法以及差分凸规划更新拉格朗日乘子项Ri+1，计算公式如公式(9)所示：
Ri+1＝U·diag(σ*)·VT (27)
其中，diag表示对角矩阵，U是m×m阶酉矩阵，σ*是特征值，VT是n×n阶酉矩阵，矩阵U和VT由矩阵(Xi+1+Ci)经过奇异值分解后得到，奇异值分解如公式(10)所示：
SVD(Xi+1+Ci)＝U·diag(σ)·VT (28)
其中，SVD表示奇异值分解，通过如下公式(11)对特征值σ进行更新：
其中，t+1表示第t+1次迭代，经过多次迭代后，σt+1收敛到最优点σ*；
根据交替方向乘子算法以及费马引理更新对偶项Ci+1，如公式(12)所示：
Ci+1＝Ci+Xi+1-Ri+1 (30)
判断若是，则根据公式(7)更新时频域待去噪的地震记录若否，则结束循环。
进一步，所述将时频域去噪后的地震记录转换到时空域，包括：利用子波数据w确定短时傅里叶逆变换的窗口长度，滑动距离；利用短时傅里叶逆变换将时频域去噪的地震记录/>转换到时空域。
根据本发明的第二方面，提供一种基于Lp伪范数和γ范数稀疏低秩约束的地震去噪装置，所述装置包括：预处理模块，被配置为根据地震资料得到时频域的地震记录目标函数构建模块，被配置为根据时频域的地震记录/>构建基于LP伪范数稀疏约束和γ范数低秩约束的目标函数；去噪模块，被配置为采用交替方向乘子法对目标函数进行迭代计算得到时频域去噪后的地震记录/>转换模块，被配置为将时频域去噪后的地震记录/>转换到时空域。
进一步，所述去噪模块被进一步配置为：
分别设置稀疏约束项、低秩约束项和对偶项的权重系数λ1、λ2和λ3的值，以及LP伪范数中参数P的值和γ范数中参数γ的值；
设置拉格朗日乘子项初始值R1＝0及其对偶项初始值C1＝0；
根据交替方向乘子算法以及软阈值收缩算法更新时频域待去噪的地震记录计算公式如公式(7)所示：
其中，表示矩阵点乘运算，i+1表示第i+1次更新，sign表示符号函数，表达式如公式(8)所示：
其中，ya,b表示矩阵中的第a行，第b列元素；
根据交替方向乘子算法以及差分凸规划更新拉格朗日乘子项Ri+1，计算公式如公式(9)所示：
Ri+1＝U·diag(σ*)·VT (33)
其中，diag表示对角矩阵，U是m×m阶酉矩阵，σ*是特征值，VT是n×n阶酉矩阵，矩阵U和VT由矩阵(Xi+1+Ci)经过奇异值分解后得到，奇异值分解如公式(10)所示：
SVD(Xi+1+Ci)＝U·diag(σ)·VT (34)
其中，SVD表示奇异值分解，通过如下公式(11)对特征值σ进行更新：
其中，t+1表示第t+1次迭代，经过多次迭代后，σt+1收敛到最优点σ*；
根据交替方向乘子算法以及费马引理更新对偶项Ci+1，如公式(12)所示：
Ci+1＝Ci+Xi+1-Ri+1 (36)
判断若是，则根据公式(7)更新时频域待去噪的地震记录若否，则结束循环。
根据本发明的第三方面，提供一种电子设备，包括：控制器；存储器，用于存储一个或多个程序，当所述一个或多个程序被所述控制器执行时，使得所述控制器实现本发明各个实施例所述的基于Lp伪范数和γ范数稀疏低秩约束的地震去噪方法。
根据本发明的第四方面，提供一种计算机可读存储介质，其上存储有计算机可读指令，当所述计算机可读指令被计算机的处理器执行时，使计算机执行本发明各个实施例所述的基于Lp伪范数和γ范数稀疏低秩约束的地震去噪方法。
根据本发明的第五方面，提供一种计算机程序产品或计算机程序，该计算机程序产品或计算机程序包括计算机指令，该计算机指令存储在计算机可读存储介质中，其特征在于，计算机设备的处理器从计算机可读存储介质读取该计算机指令，处理器执行该计算机指令，使得该计算机设备执行本发明各个实施例所述的基于Lp伪范数和γ范数稀疏低秩约束的地震去噪方法。
本发明在L1范数稀疏约束和核范数低秩约束的基础上，用LP伪范数和γ范数分别替换L1范数和核范数，并结合交替方向乘子法提出了一种基于LP伪范数稀疏约束和γ范数低秩约束的基于Lp伪范数和γ范数稀疏低秩约束的地震去噪方法，对含噪地震记录采用上述中的去噪方法进行去噪处理，通过反复迭代，输出最优去噪结果。本发明解决了现有L1稀疏约束和核范数低秩约束去噪结果噪音仍有残余且有效信息存在损失的问题，提高去噪效果。
附图说明
此处的附图被并入说明书中并构成本说明书的一部分，示出了符合本申请的实施例，并与说明书一起用于解释本申请的原理。显而易见地，下面描述中的附图仅仅是本申请的一些实施例，对于本领域普通技术者来讲，在不付出创造性劳动的前提下，还可以根据这些附图获得其他的附图。在附图中：
图1为本发明的基于Lp伪范数和γ范数稀疏低秩约束的地震去噪方法流程框图；
图2为本发明的无噪声的地震剖面；
图3为本发明的含噪声的地震剖面；
图4为本发明的去噪结果示意图；
图5为本发明的所去噪声示意图；
图6为本发明的装置结构图；
图7为本发明的电子设备的计算机系统的结构示意图。
具体实施方式
为了使本发明的目的、技术方案及优点更加清楚明白，以下结合附图及实施例，对本发明进行进一步详细说明。应当理解，此处所描述的具体实施例仅用以解释本发明，并不用于限定本发明，即所描述的实施例仅仅是本发明一部分实施例，而不是全部的实施例。通常在此处附图中描述和示出的本发明实施例的组件可以以各种不同的配置来布置和设计。
因此，以下对在附图中提供的本发明的实施例的详细描述并非旨在限制要求保护的本发明的范围，而是仅仅表示本发明的选定实施例。基于本发明的实施例，本领域技术人员在没有做出创造性劳动的前提下所获得的所有其他实施例，都属于本发明保护的范围。
需要说明的是，术语“第一”和“第二”等之类的关系术语仅仅用来将一个实体或者操作与另一个实体或操作区分开来，而不一定要求或者暗示这些实体或操作之间存在任何这种实际的关系或者顺序。而且，术语“包括”、“包含”或者其任何其他变体意在涵盖非排他性的包含，从而使得包括一系列要素的过程、方法、物品或者设备不仅包括那些要素，而且还包括没有明确列出的其他要素，或者是还包括为这种过程、方法、物品或者设备所固有的要素。在没有更多限制的情况下，由语句“包括一个……”限定的要素，并不排除在包括所述要素的过程、方法、物品或者设备中还存在另外的相同要素。
下面结合图1-5对本发明作详细说明。
本发明解决的技术问题：解决了L1范数稀疏约束和核范数低秩约束导致的去噪结果噪音仍有残余且有效信息存在损失的问题，提高去噪效果。
请参阅图1，是本发明实施例的基于Lp伪范数和γ范数稀疏低秩约束的地震去噪方法的整体流程框图。本发明实施例提供一种基于Lp伪范数和γ范数稀疏低秩约束的地震去噪方法，包括如下步骤：
步骤1：对地震资料进行预处理并获取时频域的地震记录
步骤2：根据时频域的地震记录构建基于LP伪范数稀疏约束和γ范数低秩约束的目标函数。
步骤3：采用交替方向乘子法对目标函数进行迭代计算得到时频域去噪后的地震记录
步骤4：将时频域去噪后的地震记录转换到时空域。
在一些实施例中，步骤1包括如下步骤：
步骤1.1：输入时空域含噪的地震记录S和子波数据w；
步骤1.2：利用子波数据w确定短时傅里叶变换的窗口长度，滑动距离，并对时空域含噪的地震记录S进行短时傅里叶变换，得到时频域的地震记录
在一些实施例中，步骤2包括如下步骤：
步骤2.1：根据时频域的地震记录构建目标函数中的保真项，如公式1所示：
其中，表示时频域含噪的地震记录，/>表示时频域去噪的地震记录，/>表示Frobenius范数；
步骤2.2：在公式(1)的基础上，对时频域去噪的地震记录加上LP伪范数进行稀疏约束，如公式2所示：
其中，λ1表示LP伪范数稀疏约束的权重系数，表示LP伪范数，其表达式如公式3所示：
其中a、b、m和n分别表示时频地震记录的第a行、第b列、总行数和总列数；
步骤2.3：在公式(2)的基础上，对时频域去噪的地震记录加上γ范数进行低秩约束，如公式4所示：
其中，λ2表示γ范数低秩约束的权重系数，||||γ表示γ范数，其表达式如公式5所示；
其中，表示时频域去噪的地震记录/>经过奇异值分解后的第k个特征值；
步骤2.4：在公式(4)的基础上，引入拉格朗日乘子项R和对偶项C，得到基于LP伪范数和γ范数稀疏低秩约束的地震去噪目标函数，如公式4所示：
其中λ3表示对偶项权重系数。
在一些实施例中，步骤3包括如下步骤：
步骤3.1：分别设置稀疏约束项、低秩约束项和对偶项的权重系数λ1、λ2和λ3的值，LP伪范数和γ范数中P与γ的值。
步骤3.2：设置拉格朗日乘子项初始值R1＝0及其对偶项初始值C1＝0。
步骤3.3：根据交替方向乘子算法以及软阈值收缩算法更新时频域待去噪的地震记录计算公式如公式5所示：
其中，表示矩阵点乘运算，i+1表示第i+1次更新，sign表示符号函数，表达式如公式6所示：
其中，ya,b表示矩阵中的第a行，第b列元素。
步骤3.4：根据交替方向乘子算法以及差分凸规划更新拉格朗日乘子项Ri+1，计算公式如公式9所示：
Ri+1＝U·diag(σ*)·VT (45)
其中，diag表示对角矩阵，U是m×m阶酉矩阵，σ*是特征值，VT是n×n阶酉矩阵，矩阵U和VT由矩阵(Xi+1+Ci)经过奇异值分解后得到，如公式10所示：
SVD(Xi+1+Ci)＝U·diag(σ)·VT (46)
其中，SVD表示奇异值分解。再对特征值σ进行更新，如公式11所示：
其中，t+1表示第t+1次迭代，经过几次迭代后，σt+1收敛到最优点σ*。
步骤3.5：根据交替方向乘子算法以及费马引理更新对偶项Ci+1，如公式12所示：
Ci+1＝Ci+Xi+1-Ri+1 (48)
步骤3.6：判断如果是，则返回步骤3.3进行循环，如果否，结束循环。
在一些实施例中，步骤4包括如下步骤：
步骤4.1：利用子波数据w确定短时傅里叶逆变换的窗口长度，滑动距离；
步骤4.2：利用短时傅里叶逆变换将时频域去噪的地震记录转换到时空域。
综上所述，本发明在L1范数稀疏约束和核范数低秩约束的基础上，将L1范数和核范数分别用LP伪范数和γ范数替换，解决了L1范数稀疏约束和核范数低秩约束造成的去噪结果噪音仍有残余且有效信息存在损失的问题，提高去噪效果。
下面本发明实施例将结合具体的实验以进一步说明本发明的可行性和进步性。基于如上所述的基于Lp伪范数和γ范数稀疏低秩约束的地震去噪方法，本发明实施例对野外地震记录进行试验。所述野外地震记录是在油气勘探中利用相关地震设备所采集到的地震记录。野外地震记录至少包括时空域含噪的地震记录S和子波数据w。
在利用本发明实施例所述的基于Lp伪范数和γ范数稀疏低秩约束的地震去噪方法对野外地震记录进行去噪处理，结果如图4和图5所示。将图4和图2进行对比，其中图4为本发明的去噪结果示意图，图2为本发明的无噪声的地震剖面，可以发现去噪后的地震记录与真实的无噪声的地震记录非常的接近，并结合图3，图3为本发明的含噪声的地震剖面，可以发现地震记录的去噪结果剖面相较于含噪声的地震记录剖面能更清晰的反映地下的层位和构造，证明了本发明的正确性。图5中可以看到所去除的噪声剖面上没有出现层位结构性的特征，说明去噪后地震记录的损失非常小，进一步证明了该方法的正确性。本发明在L1范数稀疏约束和核范数低秩约束的基础上，用LP伪范数和γ范数分别替换L1范数和核范数，并结合交替方向乘子法提出了一种基于LP伪范数稀疏约束和γ范数低秩约束的基于Lp伪范数和γ范数稀疏低秩约束的地震去噪方法。对含噪地震记录采用上述中的去噪方法进行去噪处理，通过反复迭代，输出最优去噪结果。本发明成功解决了现有L1稀疏约束和核范数低秩约束去噪结果噪音仍有残余且有效信息存在损失的问题，提高去噪效果。
请参阅图6，是本发明的装置结构图。本发明实施例还提供一种基于Lp伪范数和γ范数稀疏低秩约束的地震去噪装置，所述装置600包括：
预处理模块601，被配置为根据地震资料得到时频域的地震记录
目标函数构建模块602，被配置为根据时频域的地震记录构建基于LP伪范数稀疏约束和γ范数低秩约束的目标函数；
去噪模块603，被配置为采用交替方向乘子法对目标函数进行迭代计算得到时频域去噪后的地震记录
转换模块604，被配置为将时频域去噪后的地震记录转换到时空域。
在一些实施例中，所述预处理模块被进一步配置为：输入时空域含噪的地震记录S和子波数据w；利用子波数据w确定短时傅里叶变换的窗口长度，滑动距离，并对时空域含噪的地震记录S进行短时傅里叶变换，得到时频域的地震记录
在一些实施例中，所述目标函数构建模块被进一步配置为：
根据时频域的地震记录构建目标函数中的保真项，如公式(1)所示：
其中，表示时频域含噪的地震记录，/>表示时频域去噪的地震记录，/>表示Frobenius范数；
对时频域去噪的地震记录加上LP伪范数进行稀疏约束，如公式(2)所示：
其中，λ1表示LP伪范数稀疏约束的权重系数，表示LP伪范数，其表达式如公式(3)所示：
其中a、b、m和n分别表示时频地震记录的第a行、第b列、总行数和总列数；
对时频域去噪的地震记录加上γ范数进行低秩约束，如公式(4)所示：
其中，λ2表示γ范数低秩约束的权重系数，||||γ表示γ范数，其表达式如公式(5)所示；
其中，表示时频域去噪的地震记录/>经过奇异值分解后的第k个特征值；
引入拉格朗日乘子项R和对偶项C，得到基于LP伪范数和γ范数稀疏低秩约束的地震去噪目标函数，如公式(6)所示：
其中λ3表示对偶项权重系数。
在一些实施例中，所述去噪模块被进一步配置为：
分别设置稀疏约束项、低秩约束项和对偶项的权重系数λ1、λ2和λ3的值，以及LP伪范数中参数P的值和γ范数中参数γ的值；
设置拉格朗日乘子项初始值R1＝0及其对偶项初始值C1＝0；
根据交替方向乘子算法以及软阈值收缩算法更新时频域待去噪的地震记录计算公式如公式(7)所示：
其中，表示矩阵点乘运算，i+1表示第i+1次更新，sign表示符号函数，表达式如公式(8)所示：
其中，ya,b表示矩阵中的第a行，第b列元素；
根据交替方向乘子算法以及差分凸规划更新拉格朗日乘子项Ri+1，计算公式如公式(9)所示：
Ri+1＝U·diag(σ*)·VT (57)
其中，diag表示对角矩阵，U是m×m阶酉矩阵，σ*是特征值，VT是n×n阶酉矩阵，矩阵U和VT由矩阵(Xi+1+Ci)经过奇异值分解后得到，奇异值分解如公式(10)所示：
SVD(Xi+1+Ci)＝U·diag(σ)·VT (58)
其中，SVD表示奇异值分解，通过如下公式(11)对特征值σ进行更新：
其中，t+1表示第t+1次迭代，经过多次迭代后，σt+1收敛到最优点σ*；
根据交替方向乘子算法以及费马引理更新对偶项Ci+1，如公式(12)所示：
Ci+1＝Ci+Xi+1-Ri+1 (60)
判断若是，则根据公式(7)更新时频域待去噪的地震记录若否，则结束循环。
在一些实施例中，所述转换模块被进一步配置为：利用子波数据w确定短时傅里叶逆变换的窗口长度，滑动距离；利用短时傅里叶逆变换将时频域去噪的地震记录转换到时空域。
需要说明的是，上述实施例所提供的地震去噪装置与前述实施例所提供的基于Lp伪范数和γ范数稀疏低秩约束的地震去噪方法属于同一构思，其中各个模块和单元执行操作的具体方式已经在方法实施例中进行了详细描述，这里不再赘述。
本发明实施例提供了一种电子设备，包括：控制器；存储器，用于存储一个或多个程序，当一个或多个程序被控制器执行时，以执行上述各个实施例中所述的基于Lp伪范数和γ范数稀疏低秩约束的地震去噪方法。
请参阅图7，图7是本发明的一示例性实施例示出的电子设备的计算机系统的结构示意图，其示出了适于用来实现本申请实施例的电子设备的计算机系统的结构示意图。
需要说明的是，图7示出的电子设备的计算机系统700仅是一个示例，不应对本发明实施例的功能和使用范围带来任何限制。
如图7所示，计算机系统700包括中央处理单元(Central Processing Unit，CPU)701，其可以根据存储在只读存储器(Read-Only Memory，ROM)702中的程序或者从存储部分708加载到随机访问存储器(Random Access Memory，RAM)703中的程序而执行各种适当的动作和处理，例如执行上述实施例中的方法。在RAM 703中，还存储有系统操作所需的各种程序和数据。CPU 701、ROM702以及RAM 703通过总线704彼此相连。输入/输出(Input/Output，I/O)接口705也连接至总线704。
以下部件连接至I/O接口705：包括键盘、鼠标等的输入部分706；包括诸如阴极射线管(Cathode Ray Tube，CRT)、液晶显示器(Liquid Crystal Display，LCD)等以及扬声器等的输出部分707；包括硬盘等的存储部分708；以及包括诸如LAN(Local Area Network，局域网)卡、调制解调器等的网络接口卡的通信部分709。通信部分709经由诸如因特网的网络执行通信处理。驱动器710也根据需要连接至I/O接口705。可拆卸介质711，诸如磁盘、光盘、磁光盘、半导体存储器等等，根据需要安装在驱动器710上，以便于从其上读出的计算机程序根据需要被安装入存储部分708。
特别地，根据本发明的实施例，上文所述的基于Lp伪范数和γ范数稀疏低秩约束的地震去噪方法中描述的各个步骤过程可以被实现为计算机软件程序。例如，本发明实施例包括一种计算机程序产品，其包括承载在计算机可读介质上的计算机程序，该计算机程序包含用于执行流程图所示的方法的计算机程序。在这样的实施例中，该计算机程序可以通过通信部分709从网络上被下载和安装，和/或从可拆卸介质711被安装。在该计算机程序被中央处理单元(CPU)701执行时，执行本发明的系统中限定的各种功能。
需要说明的是，本发明实施例所示的计算机可读介质可以是计算机可读信号介质或者计算机可读存储介质或者是上述两者的任意组合。计算机可读存储介质例如可以是电、磁、光、电磁、红外线、或半导体的系统、装置或器件，或者任意以上的组合。计算机可读存储介质的更具体的例子可以包括但不限于：具有一个或多个导线的电连接、便携式计算机磁盘、硬盘、随机访问存储器(RAM)、只读存储器(ROM)、可擦式可编程只读存储器(Erasable Programmable Read Only Memory，EPROM)、闪存、光纤、便携式紧凑磁盘只读存储器(Compact Disc Read-Only Memory，CD-ROM)、光存储器件、磁存储器件、或者上述的任意合适的组合。在本发明中，计算机可读存储介质可以是任何包含或存储程序的有形介质，该程序可以被指令执行系统、装置或者器件使用或者与其结合使用。而在本发明中，计算机可读的信号介质可以包括在基带中或者作为载波一部分传播的数据信号，其中承载了计算机可读的计算机程序。这种传播的数据信号可以采用多种形式，包括但不限于电磁信号、光信号或上述的任意合适的组合。计算机可读的信号介质还可以是计算机可读存储介质以外的任何计算机可读介质，该计算机可读介质可以发送、传播或者传输用于由指令执行系统、装置或者器件使用或者与其结合使用的程序。计算机可读介质上包含的计算机程序可以用任何适当的介质传输，包括但不限于：无线、有线等等，或者上述的任意合适的组合。
描述于本发明实施例中所涉及到的模块可以通过软件的方式实现，也可以通过硬件的方式来实现，所描述的模块也可以设置在处理器中。其中，这些模块的名称在某种情况下并不构成对该模块本身的限定。
本发明实施例还提供了一种计算机可读存储介质，其上存储有计算机程序，该计算机程序被处理器执行时实现如前的识别方法。该计算机可读存储介质可以是上述实施例中描述的电子设备中所包含的，也可以是单独存在，而未装配入该电子设备中。
本发明实施例还提供了一种计算机程序产品或计算机程序，该计算机程序产品或计算机程序包括计算机指令，该计算机指令存储在计算机可读存储介质中。计算机设备的处理器从计算机可读存储介质读取该计算机指令，处理器执行该计算机指令，使得该计算机设备执行上述各个实施例中提供的基于Lp伪范数和γ范数稀疏低秩约束的地震去噪方法。
本发明实施例还提供了一种计算机系统，包括中央处理单元(CentralProcessing Unit，CPU)，其可以根据存储在只读存储器(Read-Only Memory，ROM)中的程序或者从存储部分加载到随机访问存储器(Random Access Memory，RAM)中的程序而执行各种适当的动作和处理，例如执行上述实施例中的基于Lp伪范数和γ范数稀疏低秩约束的地震去噪方法。在RAM中，还存储有系统操作所需的各种程序和数据。CPU、ROM以及RAM通过总线彼此相连。输入/输出(Input/Output，I/O)接口也连接至总线。
以下部件连接至I/O接口：包括键盘、鼠标等的输入部分；包括诸如阴极射线管(Cathode Ray Tube，CRT)、液晶显示器(Liquid Crystal Display，LCD)等以及扬声器等的输出部分；包括硬盘等的存储部分；以及包括诸如LAN(Local Area Network，局域网)卡、调制解调器等的网络接口卡的通信部分。通信部分经由诸如因特网的网络执行通信处理。驱动器也根据需要连接至I/O接口。可拆卸介质，诸如磁盘、光盘、磁光盘、半导体存储器等等，根据需要安装在驱动器上，以便于从其上读出的计算机程序根据需要被安装入存储部分。
上述内容，仅为本发明的较佳示例性实施例，并非用于限制本发明的实施方案，本领域普通技术人员根据本发明的主要构思和精神，可以十分方便地进行相应的变通或修改，故本发明的保护范围应以权利要求书所要求的保护范围为准。

Claims (7)

1.一种基于Lp伪范数和γ范数稀疏低秩约束的地震去噪方法，其特征在于：所述方法包括：
根据地震资料得到时频域的地震记录
根据时频域的地震记录构建基于LP伪范数稀疏约束和γ范数低秩约束的目标函数；
采用交替方向乘子法对目标函数进行迭代计算得到时频域去噪后的地震记录
将时频域去噪后的地震记录转换到时空域；
所述根据地震资料得到时频域的地震记录包括：
输入时空域含噪的地震记录S和子波数据w；
利用子波数据w确定短时傅里叶变换的窗口长度、滑动距离，并对时空域含噪的地震记录S进行短时傅里叶变换，得到时频域的地震记录
所述根据时频域的地震记录构建基于LP伪范数稀疏约束和γ范数低秩约束的目标函数，包括：
根据时频域的地震记录构建目标函数中的保真项，如公式(1)所示：
其中，表示时频域含噪的地震记录，/>表示时频域去噪的地震记录，/>表示Frobenius范数；
对时频域去噪的地震记录加上LP伪范数进行稀疏约束，如公式(2)所示：
其中，λ1表示LP伪范数稀疏约束的权重系数，表示LP伪范数，其表达式如公式(3)所示：
其中a、b、m和n分别表示时频地震记录的第a行、第b列、总行数和总列数；
对时频域去噪的地震记录加上γ范数进行低秩约束，如公式(4)所示：
其中，λ2表示γ范数低秩约束的权重系数，|| ||γ表示γ范数，其表达式如公式(5)所示；
其中，表示时频域去噪的地震记录/>经过奇异值分解后的第k个特征值；
引入拉格朗日乘子项R和对偶项C，得到基于LP伪范数和γ范数稀疏低秩约束的地震去噪目标函数，如公式(6)所示：
其中λ3表示对偶项权重系数。
2.根据权利要求1所述的基于Lp伪范数和γ范数稀疏低秩约束的地震去噪方法，其特征在于：所述采用交替方向乘子法对目标函数进行迭代计算得到时频域去噪后的地震记录包括：
分别设置稀疏约束项、低秩约束项和对偶项的权重系数λ1、λ2和λ3的值，以及LP伪范数中参数P的值和γ范数中参数γ的值；
设置拉格朗日乘子项初始值R1＝0及其对偶项初始值C1＝0；
根据交替方向乘子算法以及软阈值收缩算法更新时频域待去噪的地震记录计算公式如公式(7)所示：
其中，表示矩阵点乘运算，i+1表示第i+1次更新，sign表示符号函数，表达式如公式(8)所示：
其中，ya,b表示矩阵中的第a行，第b列元素；
根据交替方向乘子算法以及差分凸规划更新拉格朗日乘子项Ri+1，计算公式如公式(9)所示：
Ri+1＝U·diag(σ*)·VT (9)
其中，diag表示对角矩阵，U是m×m阶酉矩阵，σ*是特征值，VT是n×n阶酉矩阵，矩阵U和VT由矩阵(Xi+1+Ci)经过奇异值分解后得到，奇异值分解如公式(10)所示：
SVD(Xi+1+Ci)＝U·diag(σ)·VT (10)
其中，SVD表示奇异值分解，通过如下公式(11)对特征值σ进行更新：
其中，t+1表示第t+1次迭代，经过多次迭代后，σt+1收敛到最优点σ*；
根据交替方向乘子算法以及费马引理更新对偶项Ci+1，如公式(12)所示：
Ci+1＝Ci+Xi+1-Ri+1 (12)
判断若是，则根据公式(7)更新时频域待去噪的地震记录/>若否，则结束循环，其中，tol为设定的误差值。
3.根据权利要求2所述的基于Lp伪范数和γ范数稀疏低秩约束的地震去噪方法，其特征在于：所述将时频域去噪后的地震记录转换到时空域，包括：
利用子波数据w确定短时傅里叶逆变换的窗口长度、滑动距离；
利用短时傅里叶逆变换将时频域去噪的地震记录转换到时空域。
4.一种基于Lp伪范数和γ范数稀疏低秩约束的地震去噪装置，其特征在于：所述装置包括：
预处理模块，被配置为根据地震资料得到时频域的地震记录
目标函数构建模块，被配置为根据时频域的地震记录构建基于LP伪范数稀疏约束和γ范数低秩约束的目标函数；
去噪模块，被配置为采用交替方向乘子法对目标函数进行迭代计算得到时频域去噪后的地震记录
转换模块，被配置为将时频域去噪后的地震记录转换到时空域；
所述根据地震资料得到时频域的地震记录包括：
输入时空域含噪的地震记录S和子波数据w；
利用子波数据w确定短时傅里叶变换的窗口长度、滑动距离，并对时空域含噪的地震记录S进行短时傅里叶变换，得到时频域的地震记录
所述根据时频域的地震记录构建基于LP伪范数稀疏约束和γ范数低秩约束的目标函数，包括：
根据时频域的地震记录构建目标函数中的保真项，如公式(1)所示：
其中，表示时频域含噪的地震记录，/>表示时频域去噪的地震记录，/>表示Frobenius范数；
对时频域去噪的地震记录加上LP伪范数进行稀疏约束，如公式(2)所示：
其中，λ1表示LP伪范数稀疏约束的权重系数，表示LP伪范数，其表达式如公式(3)所示：
其中a、b、m和n分别表示时频地震记录的第a行、第b列、总行数和总列数；
对时频域去噪的地震记录加上γ范数进行低秩约束，如公式(4)所示：
其中，λ2表示γ范数低秩约束的权重系数，|| ||γ表示γ范数，其表达式如公式(5)所示；
其中，表示时频域去噪的地震记录/>经过奇异值分解后的第k个特征值；
引入拉格朗日乘子项R和对偶项C，得到基于LP伪范数和γ范数稀疏低秩约束的地震去噪目标函数，如公式(6)所示：
其中λ3表示对偶项权重系数。
5.根据权利要求4所述的地震去噪装置，其特征在于：所述去噪模块被进一步配置为：
分别设置稀疏约束项、低秩约束项和对偶项的权重系数λ1、λ2和λ3的值，以及LP伪范数中参数P的值和γ范数中参数γ的值；
设置拉格朗日乘子项初始值R1＝0及其对偶项初始值C1＝0；
根据交替方向乘子算法以及软阈值收缩算法更新时频域待去噪的地震记录计算公式如公式(7)所示：
其中，表示矩阵点乘运算，i+1表示第i+1次更新，sign表示符号函数，表达式如公式(8)所示：
其中，ya,b表示矩阵中的第a行，第b列元素；
根据交替方向乘子算法以及差分凸规划更新拉格朗日乘子项Ri+1，计算公式如公式(9)所示：
Ri+1＝U·diag(σ*)·VT (9)
其中，diag表示对角矩阵，U是m×m阶酉矩阵，σ*是特征值，VT是n×n阶酉矩阵，矩阵U和VT由矩阵(Xi+1+Ci)经过奇异值分解后得到，奇异值分解如公式(10)所示：
SVD(Xi+1+Ci)＝U·diag(σ)·VT (10)
其中，SVD表示奇异值分解，通过如下公式(11)对特征值σ进行更新：
其中，t+1表示第t+1次迭代，经过多次迭代后，σt+1收敛到最优点σ*；
根据交替方向乘子算法以及费马引理更新对偶项Ci+1，如公式(12)所示：
Ci+1＝Ci+Xi+1-Ri+1 (12)
判断若是，则根据公式(7)更新时频域待去噪的地震记录/>若否，则结束循环，其中，tol为设定的误差值。
6.一种电子设备，包括：控制器；存储器，用于存储一个或多个程序，其特征在于：当所述一个或多个程序被所述控制器执行时，使得所述控制器实现权利要求1至3所述的基于Lp伪范数和γ范数稀疏低秩约束的地震去噪方法。
7.一种计算机可读存储介质，其特征在于：其上存储有计算机可读指令，当所述计算机可读指令被计算机的处理器执行时，使计算机执行权利要求1至3中任一项所述的基于Lp伪范数和γ范数稀疏低秩约束的地震去噪方法。
CN202210534230.0A 2022-05-17 2022-05-17 基于Lp伪范数和γ范数稀疏低秩约束的地震去噪方法 Active CN114966861B (zh)

Priority Applications (1)

Application Number Priority Date Filing Date Title
CN202210534230.0A CN114966861B (zh) 2022-05-17 2022-05-17 基于Lp伪范数和γ范数稀疏低秩约束的地震去噪方法

Applications Claiming Priority (1)

Application Number Priority Date Filing Date Title
CN202210534230.0A CN114966861B (zh) 2022-05-17 2022-05-17 基于Lp伪范数和γ范数稀疏低秩约束的地震去噪方法

Publications (2)

Publication Number Publication Date
CN114966861A CN114966861A (zh) 2022-08-30
CN114966861B true CN114966861B (zh) 2024-03-26

Family

ID=82983526

Family Applications (1)

Application Number Title Priority Date Filing Date
CN202210534230.0A Active CN114966861B (zh) 2022-05-17 2022-05-17 基于Lp伪范数和γ范数稀疏低秩约束的地震去噪方法

Country Status (1)

Country Link
CN (1) CN114966861B (zh)

Families Citing this family (2)

* Cited by examiner, † Cited by third party
Publication number Priority date Publication date Assignee Title
CN115469364B (zh) * 2022-10-21 2023-04-28 成都理工大学 一种基于平稳Framelet变换的地震去噪方法及系统
CN116736378B (zh) * 2023-05-05 2024-02-02 广东海洋大学 一种基于Lq范数约束的用于地震勘探的速度分析方法

Citations (5)

* Cited by examiner, † Cited by third party
Publication number Priority date Publication date Assignee Title
CN106646612A (zh) * 2016-12-20 2017-05-10 中国地质大学（北京） 基于矩阵降秩的地震数据重建方法
CN108594304A (zh) * 2018-07-25 2018-09-28 中国石油化工股份有限公司胜利油田分公司勘探开发研究院 基于多极子库利用线性规划求解l1范数的阻抗反演方法
CN109143341A (zh) * 2018-08-20 2019-01-04 成都理工大学 基于Hampel范数的降秩滤波方法
CN110568486A (zh) * 2019-09-17 2019-12-13 电子科技大学 基于同步稀疏低秩张量补全模型的地震信号补全方法
CN111158051A (zh) * 2020-01-07 2020-05-15 自然资源部第一海洋研究所 一种基于稀疏正则化的联合约束随机噪声压制方法

Family Cites Families (3)

* Cited by examiner, † Cited by third party
Publication number Priority date Publication date Assignee Title
US20140078860A1 (en) * 2012-09-19 2014-03-20 Cgg Services Sa Interference noise attenuation method and apparatus
GB2558630A (en) * 2017-01-12 2018-07-18 Seismic Apparition Gmbh Method for dealiasing data
AU2018231078A1 (en) * 2017-03-08 2019-10-10 Saudi Arabian Oil Company Automated system and methods for adaptive robust denoising of large-scale seismic data sets

Patent Citations (5)

* Cited by examiner, † Cited by third party
Publication number Priority date Publication date Assignee Title
CN106646612A (zh) * 2016-12-20 2017-05-10 中国地质大学（北京） 基于矩阵降秩的地震数据重建方法
CN108594304A (zh) * 2018-07-25 2018-09-28 中国石油化工股份有限公司胜利油田分公司勘探开发研究院 基于多极子库利用线性规划求解l1范数的阻抗反演方法
CN109143341A (zh) * 2018-08-20 2019-01-04 成都理工大学 基于Hampel范数的降秩滤波方法
CN110568486A (zh) * 2019-09-17 2019-12-13 电子科技大学 基于同步稀疏低秩张量补全模型的地震信号补全方法
CN111158051A (zh) * 2020-01-07 2020-05-15 自然资源部第一海洋研究所 一种基于稀疏正则化的联合约束随机噪声压制方法

Non-Patent Citations (3)

* Cited by examiner, † Cited by third party
Title
An iterative curvelet thresholding algorithm for seismic random noise attenuation;Wang De-Li et al.;《APPLIED GEOPHYSICS》;20101231;第7卷(第4期);第315-324页 *
低秩矩阵恢复模型改进及其在石油测井中的应用;王艳伟等;《科学技术与工程》;20170131;第17卷(第02期);第74-80页 *
随钻地震参考信号稀疏分解去噪;秦显科等;《中国科技论文》;20161130;第11卷(第21期);第2473-2478页 *

Also Published As

Publication number Publication date
CN114966861A (zh) 2022-08-30

Similar Documents

Publication Publication Date Title
CN114966861B (zh) 基于Lp伪范数和γ范数稀疏低秩约束的地震去噪方法
US11386327B2 (en) Block-diagonal hessian-free optimization for recurrent and convolutional neural networks
Ray On Haar wavelet operational matrix of general order and its application for the numerical solution of fractional Bagley Torvik equation
Yu et al. Monte Carlo data-driven tight frame for seismic data recovery
Ratto et al. State dependent parameter metamodelling and sensitivity analysis
Yang et al. Seismic wave propagation and inversion with neural operators
DelSole et al. Confidence intervals in optimal fingerprinting
Messner et al. An accelerated symmetric time-domain boundary element formulation for elasticity
CN107798156B (zh) 一种频率域2.5维粘弹性波数值模拟方法及装置
Huang et al. Erratic noise suppression using iterative structure‐oriented space‐varying median filtering with sparsity constraint
Sharbati et al. Detection and extraction of velocity pulses of near-fault ground motions using asymmetric Gaussian chirplet model
CN112767959A (zh) 语音增强方法、装置、设备及介质
Bindi et al. Seismic input motion determined from a surface–downhole pair of sensors: A constrained deconvolution approach
Verma et al. Solution for generalized fuzzy fractional Kortewege-de Varies equation using a robust fuzzy double parametric approach
CN110646841A (zh) 时变稀疏反褶积方法及系统
CN113640871B (zh) 一种基于重加权l1范数稀疏约束的地震波阻抗反演方法
JP7168639B2 (ja) 計算機システム及びデータ処理方法
US20220027722A1 (en) Deep Relational Factorization Machine Techniques for Content Usage Prediction via Multiple Interaction Types
Liu et al. A sparse surrogate model for structural reliability analysis based on the generalized polynomial chaos expansion
CN112485826B (zh) 绝对波阻抗反演成像方法、装置、设备及存储介质
CN110850473B (zh) 一种基于稀疏变换学习的地震波阻抗反演方法
Geetha et al. Microseismic signal denoising based on variational mode decomposition with adaptive non-local means filtering
Wang et al. Total variation regularized self-supervised Bayesian deep learning for seismic random noise attenuation
Khomrutai Nonlocal equations with regular varying decay solutions
Wang et al. A frequency-dependent ground-motion spatial correlation model of within-event residuals for Fourier amplitude spectra

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