CN114442157A - 地震反射系数的确定方法和装置 - Google Patents

地震反射系数的确定方法和装置 Download PDF

Info

Publication number
CN114442157A
CN114442157A CN202011208807.6A CN202011208807A CN114442157A CN 114442157 A CN114442157 A CN 114442157A CN 202011208807 A CN202011208807 A CN 202011208807A CN 114442157 A CN114442157 A CN 114442157A
Authority
CN
China
Prior art keywords
objective function
seismic
reflection
function
reflection coefficient
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.)
Pending
Application number
CN202011208807.6A
Other languages
English (en)
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 CN202011208807.6A priority Critical patent/CN114442157A/zh
Publication of CN114442157A publication Critical patent/CN114442157A/zh
Pending legal-status Critical Current

Links

Images

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/30—Analysis
    • G01V1/307—Analysis for determining seismic attributes, e.g. amplitude, instantaneous phase or frequency, reflection strength or polarity
    • 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/282—Application of seismic models, synthetic seismograms
    • G—PHYSICS
    • G01—MEASURING; TESTING
    • G01V—GEOPHYSICS; GRAVITATIONAL MEASUREMENTS; DETECTING MASSES OR OBJECTS; TAGS
    • G01V2210/00—Details of seismic processing or analysis
    • G01V2210/40—Transforming data representation
    • G01V2210/45—F-x or F-xy domain
    • G—PHYSICS
    • G01—MEASURING; TESTING
    • G01V—GEOPHYSICS; GRAVITATIONAL MEASUREMENTS; DETECTING MASSES OR OBJECTS; TAGS
    • G01V2210/00—Details of seismic processing or analysis
    • G01V2210/60—Analysis
    • G01V2210/61—Analysis by combining or comparing a seismic data set with other data
    • G01V2210/616—Data from specific type of measurement
    • G01V2210/6161—Seismic or acoustic, e.g. land or sea measurements

Landscapes

  • Engineering & Computer Science (AREA)
  • Remote Sensing (AREA)
  • Physics & Mathematics (AREA)
  • Life Sciences & Earth Sciences (AREA)
  • Acoustics & Sound (AREA)
  • Environmental & Geological Engineering (AREA)
  • Geology (AREA)
  • General Life Sciences & Earth Sciences (AREA)
  • General Physics & Mathematics (AREA)
  • Geophysics (AREA)
  • Geophysics And Detection Of Objects (AREA)

Abstract

本说明书提供了一种地震反射系数的确定方法和装置。基于该方法，在根据叠后地震资料和地震子波褶积模型，构建得到关于目标区域的地震反射系数的第一目标函数之后，可以先构建与上述第一目标函数对应的对偶近似函数，作为第二目标函数，从而可以将直接求解第一目标函数的最优值问题转换为了相应的对偶问题，避免了对第一目标函数的直接求解；再通过采用基于半定规划理论的预设处理方式来处理所述第二目标函数，以确定出反射时间，这样可以不需要对反射系数的稀疏域进行离散化，能够有效地避免出现基失配的问题，从而可以准确地反演得到精度较高的反射系数，解决了现有方法中存在的所确定的反射系数误差大、精度低的技术问题。

Description

地震反射系数的确定方法和装置
技术领域
本说明书属于地球物理勘探技术领域，尤其涉及地震反射系数的确定方法和装置。
背景技术
在地球物理勘探中，在数据较少(例如缺少测井资料等)的情况下，往往会先利用叠后地震资料进行稀疏反演，以得到反射系数的反演结果。
而基于现有方法，在进行具体的稀疏反演时，往往是使用傅里叶基(DFT基)作为反射系数的变换基，采用离散化的方法将连续参数空间简化为有限的网格点集，进行相应的处理。但是，实际上反射系数是在连续参数空间中分布的。因此，通过上述离散化的方法，无论网格点集中的网格点划分得多细密，都无法保证全部的数据信号都处在网格单元的中心，进而会导致出现基失配的问题，影响最终所确定的反射系数的精度。
针对上述技术问题，目前尚未提出有效的解决方案。
发明内容
本说明书提供了一种地震反射系数的确定方法和装置，可以不需要对反射系数的稀疏域进行离散化，能够有效地避免出现基失配的问题，进而能准确地反演得到精度较高的反射系数，解决了现有方法中存在的所确定的反射系数反射时间误差大、稀疏度不足从而导致精度低的技术问题。
本说明书提供的一种地震反射系数的确定方法，包括：
获取目标区域的叠后地震资料；
根据所述叠后地震资料和地震子波褶积模型，构建关于目标区域的地震反射系数的第一目标函数；
根据所述第一目标函数，构建对应的对偶近似函数，作为第二目标函数；
通过基于半定规划理论的预设处理方式处理所述第二目标函数，以确定出反射时间；
根据所述反射时间，确定出目标区域的地震反射系数。
在一个实施例中，通过基于半定规划理论的预设处理方式处理所述第二目标函数，以确定出反射时间，包括：
获取并利用Hermite矩阵，将第二目标函数转换为基于半定规划理论的函数，作为第三目标函数；
根据所述第三目标函数，确定出反射时间。
在一个实施例中，获取并利用Hermite矩阵，将第二目标函数转换为基于半定规划理论的函数，作为第三目标函数，包括：构建以下算式作为第三目标函数：
Figure BDA0002758068070000024
基于
Figure BDA0002758068070000021
其中，Q具体可以表示为Hermite矩阵，c具体可以表示为目标系数，i为矩阵的行编号，j为矩阵的列编号，n为矩阵的大小。
在一个实施例中，在获取并利用Hermite矩阵，将第二目标函数转换为基于半定规划理论的函数，作为第三目标函数之后，所述方法还包括：
根据目标区域的叠后地震资料，确定噪声水平；
根据所述噪声水平，确定出正则化调整参数；
根据所述正则化调整参数和所述第三目标函数，构建第四目标函数；
相应的，根据所述第四目标函数，确定出反射时间。
在一个实施例中，根据所述正则化调整参数和所述第三目标函数，构建第四目标函数，包括：构建以下算式作为第四目标函数：
Figure BDA0002758068070000022
基于
Figure BDA0002758068070000023
其中，λ具体可以表示为正则化调整参数。
在一个实施例中，根据所述第四目标函数，确定出反射时间，包括：
求解所述第四目标函数，确定出目标系数；
根据所述目标系数和约束函数的多项式，计算中间参数；
根据所述中间参数，确定所述反射时间。
在一个实施例中，根据所述反射时间，确定出目标区域的地震反射系数，包括：
将所述反射时间代入基于傅里叶变换的反射系数的表征函数，并使用最小二乘法求解所述表征函数，以得到对应的反射幅值。
在一个实施例中，在根据所述反射时间，确定出目标区域的地震反射系数之后，所述方法还包括：
根据所述目标区域的地震反射系数，指导对目标区域进行地震勘探。
本说明书提供了一种地震反射系数的确定装置，包括：
获取模块，用于获取目标区域的叠后地震资料；
第一构建模块，用于根据所述叠后地震资料和地震子波褶积模型，构建关于目标区域的地震反射系数的第一目标函数；
第二构建模块，用于根据所述第一目标函数，构建对应的对偶近似函数，作为第二目标函数；
处理模块，用于通过基于半定规划理论的预设处理方式处理所述第二目标函数，以确定出反射时间；
确定模块，用于根据所述反射时间，确定出目标区域的地震反射系数。
本说明书提供了一种服务器，包括处理器以及用于存储处理器可执行指令的存储器，所述处理器执行所述指令时实现获取目标区域的叠后地震资料；根据所述叠后地震资料和地震子波褶积模型，构建关于目标区域的地震反射系数的第一目标函数；根据所述第一目标函数，构建对应的对偶近似函数，作为第二目标函数；通过基于半定规划理论的预设处理方式处理所述第二目标函数，以确定出反射时间；根据所述反射时间，确定出目标区域的地震反射系数。
本说明书提供的一种地震反射系数的确定方法和装置，在根据叠后地震资料和地震子波褶积模型，构建得到关于目标区域的地震反射系数的第一目标函数之后，可以先构建出与上述第一目标函数对应的对偶近似函数，作为第二目标函数，从而可以将直接求解第一目标函数的最优值的问题转换为了相应的对偶问题，避免了对第一目标函数的直接求解；再通过采用基于半定规划理论的预设处理方式来处理所述第二目标函数，以确定出反射时间，这样可以不需要对反射系数的稀疏域进行离散化，能够有效地避免出现基失配的问题，从而可以准确地反演得到精度较高的反射系数，解决了现有方法中存在的所确定的反射系数误差大、精度低的技术问题。
附图说明
为了更清楚地说明本说明书实施例，下面将对实施例中所需要使用的附图作简单地介绍，下面描述中的附图仅仅是本说明书中记载的一些实施例，对于本领域普通技术人员来讲，在不付出创造性劳动性的前提下，还可以根据这些附图获得其他的附图。
图1是本说明书的一个实施例提供的地震反射系数的确定方法的流程示意图；
图2是本说明书的一个实施例提供的服务器的结构组成示意图；
图3是本说明书的一个实施例提供的地震反射系数的确定装置的结构组成示意图；
图4是在一个场景示例中，应用本说明书实施例提供的地震反射系数的确定方法的一种实施例的示意图；
图5是在一个场景示例中，应用本说明书实施例提供的地震反射系数的确定方法的一种实施例的示意图；
图6是在一个场景示例中，应用本说明书实施例提供的地震反射系数的确定方法的一种实施例的示意图；
图7是在一个场景示例中，应用本说明书实施例提供的地震反射系数的确定方法的一种实施例的示意图；
图8是在一个场景示例中，应用本说明书实施例提供的地震反射系数的确定方法的一种实施例的示意图；
图9是在一个场景示例中，应用本说明书实施例提供的地震反射系数的确定方法的一种实施例的示意图；
图10是在一个场景示例中，应用本说明书实施例提供的地震反射系数的确定方法的一种实施例的示意图；
图11是在一个场景示例中，应用本说明书实施例提供的地震反射系数的确定方法的一种实施例的示意图。
具体实施方式
为了使本技术领域的人员更好地理解本说明书中的技术方案，下面将结合本说明书实施例中的附图，对本说明书实施例中的技术方案进行清楚、完整地描述，显然，所描述的实施例仅仅是本说明书一部分实施例，而不是全部的实施例。基于本说明书中的实施例，本领域普通技术人员在没有作出创造性劳动前提下所获得的所有其他实施例，都应当属于本说明书保护的范围。
考虑到基于现有的稀疏反演方法，在基于叠后地震资料反演获取反射系数时，往往需要使用例如傅里叶作为变换基，采用离散化的方法将连续参数空间简化为有限的网格点集，在进行具体的数据处理。但是，实际上反射系数是在连续参数空间中分布的。通过上述离散化的方法，无论网格点集中的网格点划分得多细密，都无法保证全部的数据信号都处在网格单元的中心，导致很容易出现基失配的问题。而基失配程度越大，又会使得反演结果的误差越大，进而会影响最终确定反射系数的精度。
可见，基于现有方法，通过使用离散的稀疏基(即傅里叶基)对反射系数的稀疏域进行离散化处理，无法避免出现基失配的问题，导致现有方法具体实施时，往往会存在所确定的反射系数误差大、精度低的技术问题。
针对产生上述问题的根本原因，本说明书考虑在进行稀疏反演的过程中，可以摒弃使用离散的稀疏基对反射系数的稀疏域进行离散化处理的方式，从而可以在源头避免出现基失配的问题。相对的，在根据叠后地震资料和地震子波褶积模型，构建得到关于目标区域的地震反射系数的第一目标函数之后，可以先构建一个与上述第一目标函数对应的对偶近似函数，作为第二目标函数，将直接求解第一目标函数的最优值问题转换为了相应的对偶问题，避免了对第一目标函数的直接求解；进一步，再通过采用基于半定规划理论的预设处理方式来处理所述第二目标函数，以确定出反射时间，这样可以不需要对反射系数的稀疏域进行离散化，从而可以准确地反演得到精度较高的反射系数，解决现有方法中存在的所确定的反射系数误差大、精度低的技术问题。
基于上述思路，参阅图1所示，本说明书实施例提供了一种地震反射系数的确定方法。该方法具体实施时，可以包括以下内容。
S101：获取目标区域的叠后地震资料。
S102：根据所述叠后地震资料和地震子波褶积模型，构建关于目标区域的地震反射系数的第一目标函数。
S103：根据所述第一目标函数，构建对应的对偶近似函数，作为第二目标函数。
S104：通过基于半定规划理论的预设处理方式处理所述第二目标函数，以确定出反射时间。
S105：根据所述反射时间，确定出目标区域的地震反射系数。
在一个实施例中，上述目标区域具体可以理解为待进行地震勘探的区域。上述目标区域的叠后地震资料具体可以是从所述目标区域中采集得到的地震数据。具体的，上述叠后地震资料可以包括地震记录、地震子波等数据。
具体实施时，可以根据叠后地震资料，利用地震记录频谱除以地震子波频谱，得到对应的反射系数的频谱数据(也称反射系数频谱)。例如，可以记为y＝[R(0) R(1) … R(N-1)]T)。其中，上述反射系数的频谱数据(y)中的N可以表示为频域采样点的总数。此外，还可以根据叠后地震资料，通过观察等方式，估计出相应的噪声水平。例如，估计噪声是强还是弱。
在一个实施例中，上述地震褶积模型(convolution model)，具体可以理解为一种制作合成地震记录的模型。基于该地震褶积模型，是假设每道地震记录是由地震子波与地下模型各层的反射函数的褶积所构成的。
在一个实施例中，根据地震褶积模型，地震记录(可以记为s(t))可表示为以下形式：
s(t)＝w(t)*r(t)， (1)
其中，w(t)具体可以表示地震子波，r(t)具体可以表示反射系数序列。
在频率域式(1)可以进一步表示为以下形式：
S(f)＝W(f)×R(f)， (2)
其中，S(f)具体可以表示为地震记录的谱，W(f)具体可以表示为地震子波的谱，R(f)具体可以表示为反射系数序列的谱。
在一个实施例中，如果基于现有方法具体处理时，通常可以先假设地震子波的谱已知，则可以利用地震记录的谱除以地震子波的谱以获得反射系数序列的谱，并结合离散傅里叶变换，得到在频率域进行反射系数反演的以下算式：
Figure BDA0002758068070000061
其中，
Figure BDA0002758068070000062
基于层状介质假设，反射系数序列是稀疏的。进而结合算式(3)，可知r在傅里叶基下也是稀疏的。
又考虑到反射系数序列本质上属于连续时间信号，因而只有在N→∞时，算式(3)才能精确地表示反射系数序列与其谱之间的关系，此时算式(3)可以表示为：
R＝ψ1θ ψ1∈CN×N (4)
其中，ψ1具体可以表示为精确的傅里叶基，θ具体可以表示为真实反射系数序列。然而，实际处理时，N的取值显然不可能无限大，即反演网格不可能无限小，因此势必会出现基失配。于是，算式(3)可以表示为：
R＝ψ0x ψ0∈CN×N (5)
其中，ψ0为常规离散化下的DFT基，x为基于ψ0重构的反射系数序列，此时是不可能使其中各个反射位置精确地落在由ψ0所划分的网格点上的。因此，求得的x往往会存在较多误差。结合算式(4)和(5)，可以得到以下算式：
x＝ψθ， (6)
其中，
Figure BDA0002758068070000063
ψ具体可以表示出ψ0与ψ1之间的失配程度。若ψ0＝ψ1，则ψ＝I，即θ＝x，此时不存在基失配的问题。若ψ0≠ψ1，则令ψ＝I+E，I为单位矩阵，E为扰动项。
基于上述算式可知，当反射系数全部位于预设网格点上时，利用压缩感知原理的常用地震反演方法能够实现精确稀疏反演，但正常情况下，反射系数的真实位置往往偏离预设网格点，存在基失配的问题，导致反演结果会受基失配影响而存在误差。
在一个实施例中，基于本说明书提供的方法，可以先根据上述地震子波褶积模型，结合层状介质假设，将反射系数序列r可表示为以下形式：
Figure BDA0002758068070000071
其中，N具体可以表示为地震信号采样点(时域采样点)的个数，aj具体可以表示为反射系数在编号为j的采样点的振幅(即反射幅值)，tj具体表示为反射系数在编号为j的采样点的时间位置(即反射时间)，δ具体可以表示为一个脉冲函数。根据傅里叶变换，可以将算式(7)转换为以下形式：
Figure BDA0002758068070000072
其中，j具体可以表示为时域采样点的编号，aj具体可以表示为编号为j的时域采样点的反射幅值，tj具体可以表示为反射系数在编号为j的时域采样点处的反射时间。
需要补充的是，基于层状介质假设，可以认为地层的反射系数是稀疏的，褶积模型可以看作地震记录的稀疏表示模型，进而可以利用反射系数的稀疏性约束进行地震反演。
结合上述算式，为了从低频地震记录中恢复出频带较宽的稀疏脉冲反射系数序列，可以利用叠后地震资料中的反射系数的频谱数据，建立反射系数在L1范数极小化意义下的最优化求解问题，从构建得到关于目标区域的地震反射系数的第一目标函数。
在一个实施例中，上述根据所述叠后地震资料和地震子波褶积模型，构建关于目标区域的地震反射系数的第一目标函数，具体实施时，可以按照以下算式构建第一目标函数：
min||x||1基于y＝Fx (9)
其中，x可以表示为反射系数的幅值数据，具体可以表示为x＝[a0 a1 … aN-1]T，y可以表示为反射系数的频谱数据，具体可以表示为y＝[R(0) R(1) … R(N-1)]T，T表示转置，F具体可以表示为傅里叶变换矩阵，||.||1具体可以表示为求解一范数。
在一个实施例中，在具体处理上述第一目标函数时，考虑到求解第一目标函数中的一范数带有绝对值算子，是一个高度非线性的凸优化问题，因此可以先将上述求解第一目标函数的最优值问题先转换成相对应的对偶问题。具体的，可以先根据第一目标函数，构建对应的对偶近似函数(或者称对偶近似方程)作为第二目标函数。再通过处理上述第二目标函数，来处理对偶问题，得到所需要的数据结果，避免对第一目标函数的直接求解，从而可以避免对绝对值算子进行求导，减少求解过程的误差，更加高效地得到较为准确的求解结果。
在一个实施例中，根据所述第一目标函数，构建对应的对偶近似函数，作为第二目标函数，具体实施时，可以按照以下算式构建第二目标函数：
max Re<y，c>，基于||F*c||∞≤1 (10)
其中，c具体可以表示为目标系数，具体可以是一个目标系数向量，Re<y，c>具体可以表示为复向量内积后取实部的计算，F具体可以表示为傅里叶矩阵，*表示共轭，||F*c||∞≤1为约束条件。其中，上述复向量内积可以表示为：
Figure BDA0002758068070000081
上述约束函数可以表示为：
Figure BDA0002758068070000082
其中，k具体可以表示为频域采样点的编号，fc具体可以表示为高截频。
上述处理是考虑到由于算式(9)中的一范数带有绝对值算子，通常无法直接求导。又注意到算式(9)求解一范数极小化的问题可以视为一个高度非线性的优化问题，因此可以利用半定规划理论进行求解。若要利用半定规划理论对其进行求解，又需要先过渡到其对偶函数。通过上述处理，可以将原本求解第一目标函数中的x的一范数最小值的问题，转换为求解对偶近似函数(即算式(10))中的最大值问题。
在一个实施例中，在处理上述第二目标函数的过程中，为了避免使用离散的稀疏基(例如，傅里叶基等)对反射系数的稀疏域进行离散化处理所带来的基失配的问题，可以选择采用基于半定规划理论的预设处理方式来处理所述第二目标函数，以确定出反射时间，避免在处理过程中由于出现基失配的问题所带来的误差。
在一个实施例中，通常对于求一组变量有可能受限于线性不等式，线性对称矩阵不等式，半定约束的线性函数的最小(最大)的一类问题可以称为半定规划问题。
在一个实施例中，上述通过基于半定规划理论的预设处理方式处理所述第二目标函数，以确定出反射时间，具体实施时，可以包括以下内容：
S1：获取并利用Hermite矩阵，将第二目标函数转换为基于半定规划理论的函数，作为第三目标函数；
S2：根据所述第三目标函数，确定出反射时间。
在一个实施例中，基于半定规划理论，在求解第二目标函数(10)时，往往会利用到以下定理(Candès和Fernandez-Granda，2014)：一个三角多项式
Figure BDA0002758068070000083
以1为界，其中
Figure BDA0002758068070000084
有且仅有一个Hermite矩阵
Figure BDA0002758068070000085
满足以下关系：
Figure BDA0002758068070000086
其中，由于半正定约束等价于Q-cc*≥0，因此，存在以下算式：
z*Qz≥|c*z|2。 (13)
其中，Z*具体可以表示为Z的共轭转置。
基于上述定理，设zk＝ei2πkt，其中，t可以表示为反射时间。由算式(12)可以得出以下关系：
z*Qz＝1，|c*z|2＝|(F*c)(t)|2。 (14)
因此，可以得到以下约束条件：
1≥|(F*c)(t)|2。 (15)
通过上述约束条件可以表明F*c是均匀有界的。因此，可以引入并利用Hermite矩阵，将第二目标函数转换为基于半定规划理论的函数，得到第三目标函数。进而可以将关于第二目标函数的对偶问题，进一步转化为关于第三目标函数的半定规划问题。
在一个实施例中，上述获取并利用Hermite矩阵，将第二目标函数转换为基于半定规划理论的函数，作为第三目标函数，具体实施时，可以包括：构建以下算式作为第三目标函数：
Figure BDA0002758068070000091
基于
Figure BDA0002758068070000092
其中，Q为Hermite矩阵，c为目标系数，i为矩阵的行编号，j为矩阵的列编号，n为矩阵的大小。
在一个实施例中，可以先通过求解上述第三目标函数，得到对应的目标系数；再根据上述目标系数，计算出对应的反射时间。
在一个实施例中，为了能够进一步更加精准地确定出反射时间，还可以引入关于稀疏度的正则化调整参数，对所述第三目标函数进行适应性调整，得到更加符合目标区域的反射系数的真实稀疏情况的第四目标函数，进而可以利用上述第四目标函数代替第三目标函数，得到精度相对更高的反射时间，以进一步减少误差。
在一个实施例中，在获取并利用Hermite矩阵，将第二目标函数转换为基于半定规划理论的函数，作为第三目标函数之后，所述方法具体实施时，还可以包括以下内容：根据目标区域的叠后地震资料，确定噪声水平；根据所述噪声水平，确定出正则化调整参数；根据所述正则化调整参数和所述第三目标函数，构建第四目标函数；相应的，根据所述第四目标函数，确定出反射时间。
具体实施时，可以根据目标区域的叠后地震资料，通过观察等方式估计出上述噪声水平。
在一个实施例中，上述根据所述噪声水平，确定出正则化调整参数，具体实施时，可以包括：在确定噪声水平的数值较大的情况下，可以将所述正则化调整参数的数值设置得相对较大，以压制多余的反射，避免反演结果出现噪声导致的虚假反射；在确定噪声水平的数值较小的情况下，可以将所述正则化调整参数的数值设置得相对较小，以更好地保留各层反射系数，从而可以得到相对于第三目标函数更加准确、效果更好的第四目标函数。
在一个实施例中，上述根据所述正则化调整参数和所述第三目标函数，构建第四目标函数，具体实施时，可以包括：构建以下算式作为第四目标函数：
Figure BDA0002758068070000101
基于
Figure BDA0002758068070000102
其中，λ具体可以表示为正则化调整参数。
进一步，可以通过利用上述第四目标函数代替之前所使用第三目标函数，更加准确地确定出对应的反射时间。
在一个实施例中，上述根据所述第四目标函数，确定出反射时间，具体实施时，可以包括以下内容：求解所述第四目标函数，确定出目标系数；根据所述目标系数和约束函数的多项式，计算中间参数(可以记为z)；根据所述中间参数，确定所述反射时间。
在一个实施例中，具体实施时，可以利用Matlab中的cvx程序包具体求解第三目标函数或第四目标函数，以计算出目标系数c。当然，上述所列举的求解第三目标函数或第四目标函数的方式只是一种示意性说明。具体实施时，根据具体情况和精度要求，灵活地采用其他合适的求解方式来求解得到目标系数c。
在一个实施例中，具体计算中间参数时，可以先将上述已经求解得到的目标系数c先代入如下所示的约束函数的多项式中：
Figure BDA0002758068070000103
其中，z具体可以表示为待求的中间参数，上述约束函数的多项式具体可以是根据算式(15)变形得到的，下标中的n等于2fc+1(相当于频率采样点个数，其中fc为高截频)，
Figure BDA0002758068070000104
中的-表示共轭。
进一步，再对上述算式(18)进行求解，得到三角多项式的根值z作为所述中间参数。由于z的模量是1，只需要考虑单位圆上的根。
进而，可以根据中间参数z计算出对应的反射时间。
再将上述反射时间代入算式(8)(可以记为基于傅里叶变换的反射系数的表征函数)，即
Figure BDA0002758068070000105
并对该算式(8)进行求解处理，以确定出对应的反射幅值，从而可以确定出所需要的反射系数。
在一个实施例中，具体求解上述代入了反射时间的算式(8)时，可以采用最小二乘法进行求解处理，从而可以更加高效地求解得到对应的反射幅值。当然，上述所列的最小二乘法只是一种示意性说明。具体实施时，根据具体情况和精度要求，还可以采用其他合适的求解处理方式进行反射幅值的求解。对此，本说明书不作限定。
在一个实施例中，上述根据所述反射时间，确定出目标区域的地震反射系数，具体实施时，可以包括以下内容：将所述反射时间代入基于傅里叶变换的反射系数的表征函数，并使用最小二乘法求解所述表征函数，以得到对应的反射幅值。
在一个实施例中，在根据所述反射时间，确定出目标区域的地震反射系数之后，所述方法具体实施时，还可以包括以下内容：根据所述目标区域的地震反射系数，指导对目标区域进行地震勘探。
在一个实施例中，在对目标区域勘探开发的各个阶段，尤其是地震勘探初期，所拥有的数据资料往往相对有限，例如，可能还没有目标区域的测井资料等。这时，可以先利用本说明书所提供的上述方法，对所获取的叠后地震资料进行处理，得到精度较高、误差较小，且分辨率较好的反射系数。进而可以利用上述反射系数对目标区域进行前期地震勘探，得到前期地震勘探结果。后续可以利用上述前期地震勘探结果作为参考依据，指导进行进一步的更加深入的地震勘探。
由上可见，本说明书实施例提供的地震反射系数的确定方法，由于在根据叠后地震资料和地震子波褶积模型，构建得到关于目标区域的地震反射系数的第一目标函数之后，可以先构建与上述第一目标函数对应的对偶近似函数，作为第二目标函数，从而将直接求解第一目标函数的最优值问题转换为了相应的对偶问题，避免了对第一目标函数的直接求解；再通过采用基于半定规划理论的预设处理方式来处理所述第二目标函数，以确定出反射时间，这样可以不需要对反射系数的稀疏域进行离散化，能够有效地避免出现基失配的问题，从而可以准确地反演得到精度较高的反射系数，解决了现有方法中存在的所确定的反射系数反射时间误差大、稀疏度不足从而导致精度低的技术问题。
本说明书实施例还提供一种服务器，包括处理器以及用于存储处理器可执行指令的存储器，所述处理器具体实施时可以根据指令执行以下步骤：获取目标区域的叠后地震资料；根据所述叠后地震资料和地震子波褶积模型，构建关于目标区域的地震反射系数的第一目标函数；根据所述第一目标函数，构建对应的对偶近似函数，作为第二目标函数；通过基于半定规划理论的预设处理方式处理所述第二目标函数，以确定出反射时间；根据所述反射时间，确定出目标区域的地震反射系数。
为了能够更加准确地完成上述指令，参阅图2所示，本说明书实施例还提供了另一种具体的服务器，其中，所述服务器包括网络通信端口201、处理器202以及存储器203，上述结构通过内部线缆相连，以便各个结构可以进行具体的数据交互。
其中，所述网络通信端口201，具体可以用于获取目标区域的叠后地震资料。
所述处理器202，具体可以用于根据所述叠后地震资料和地震子波褶积模型，构建关于目标区域的地震反射系数的第一目标函数；根据所述第一目标函数，构建对应的对偶近似函数，作为第二目标函数；通过基于半定规划理论的预设处理方式处理所述第二目标函数，以确定出反射时间；根据所述反射时间，确定出目标区域的地震反射系数。
所述存储器203，具体可以用于存储相应的指令程序。
在本实施例中，所述网络通信端口201可以是与不同的通信协议进行绑定，从而可以发送或接收不同数据的虚拟端口。例如，所述网络通信端口可以是负责进行web数据通信的端口，也可以是负责进行FTP数据通信的端口，还可以是负责进行邮件数据通信的端口。此外，所述网络通信端口还可以是实体的通信接口或者通信芯片。例如，其可以为无线移动网络通信芯片，如GSM、CDMA等；其还可以为Wifi芯片；其还可以为蓝牙芯片。
在本实施例中，所述处理器202可以按任何适当的方式实现。例如，处理器可以采取例如微处理器或处理器以及存储可由该(微)处理器执行的计算机可读程序代码(例如软件或固件)的计算机可读介质、逻辑门、开关、专用集成电路(Application SpecificIntegrated Circuit，ASIC)、可编程逻辑控制器和嵌入微控制器的形式等等。本说明书并不作限定。
在本实施例中，所述存储器203可以包括多个层次，在数字系统中，只要能保存二进制数据的都可以是存储器；在集成电路中，一个没有实物形式的具有存储功能的电路也叫存储器，如RAM、FIFO等；在系统中，具有实物形式的存储设备也叫存储器，如内存条、TF卡等。
本说明书实施例还提供了一种基于上述地震反射系数的确定方法的计算机存储介质，所述计算机存储介质存储有计算机程序指令，在所述计算机程序指令被执行时实现：获取目标区域的叠后地震资料；根据所述叠后地震资料和地震子波褶积模型，构建关于目标区域的地震反射系数的第一目标函数；根据所述第一目标函数，构建对应的对偶近似函数，作为第二目标函数；通过基于半定规划理论的预设处理方式处理所述第二目标函数，以确定出反射时间；根据所述反射时间，确定出目标区域的地震反射系数。
在本实施例中，上述存储介质包括但不限于随机存取存储器(Random AccessMemory,RAM)、只读存储器(Read-Only Memory,ROM)、缓存(Cache)、硬盘(Hard DiskDrive,HDD)或者存储卡(Memory Card)。所述存储器可以用于存储计算机程序指令。网络通信单元可以是依照通信协议规定的标准设置的，用于进行网络连接通信的接口。
在本实施例中，该计算机存储介质存储的程序指令具体实现的功能和效果，可以与其它实施方式对照解释，在此不再赘述。
参阅图3所示，在软件层面上，本说明书实施例还提供了一种地震反射系数的确定装置，该装置具体可以包括以下的结构模块。
获取模块301，具体可以用于获取目标区域的叠后地震资料；
第一构建模块302，具体可以用于根据所述叠后地震资料和地震子波褶积模型，构建关于目标区域的地震反射系数的第一目标函数；
第二构建模块303，具体可以用于根据所述第一目标函数，构建对应的对偶近似函数，作为第二目标函数；
处理模块304，具体可以用于通过基于半定规划理论的预设处理方式处理所述第二目标函数，以确定出反射时间；
确定模块305，具体可以用于根据所述反射时间，确定出目标区域的地震反射系数。
在一个实施例中，上述处理模块304具体可以包括以下结构单元：
转换单元，具体可以用于获取并利用Hermite矩阵，将第二目标函数转换为基于半定规划理论的函数，作为第三目标函数；
确定单元，具体可以用于根据所述第三目标函数，确定出反射时间。
在一个实施例中，上述处理模块304具体实施时，可以构建以下算式作为第三目标函数：
Figure BDA0002758068070000134
基于
Figure BDA0002758068070000131
其中，Q为Hermite矩阵，c为目标系数，i为矩阵的行编号，j为矩阵的列编号，n为矩阵的大小。
在一个实施例中，上述处理模块304在获取并利用Hermite矩阵，将第二目标函数转换为基于半定规划理论的函数，作为第三目标函数之后，还用于根据目标区域的叠后地震资料，确定噪声水平；根据所述噪声水平，确定出正则化调整参数；根据所述正则化调整参数和所述第三目标函数，构建第四目标函数。相应的，上述确定模块305可以用于根据所述第四目标函数，确定出反射时间。
在一个实施例中，所述处理模块304还可以构建以下算式作为第四目标函数：
Figure BDA0002758068070000132
基于
Figure BDA0002758068070000133
其中，λ为正则化调整参数。
在一个实施例中，所述确定模块305具体实施时，可以用于求解所述第四目标函数，确定出目标系数；根据所述目标系数和约束函数的多项式，计算中间参数；根据所述中间参数，确定所述反射时间。
在一个实施例中，上述确定模块305在根据所述反射时间，确定出目标区域的地震反射系数时，可以将反射时间代入基于傅里叶变换的反射系数的表征函数，并使用最小二乘法求解所述表征函数，以得到对应的反射幅值，从而确定出目标区域的地震反射系数。
在一个实施例中，上述装置还可以包括勘探模块，具体实施时，可以用于根据所述目标区域的地震反射系数，对目标区域进行前期地震勘探。
需要说明的是，上述实施例阐明的单元、装置或模块等，具体可以由计算机芯片或实体实现，或者由具有某种功能的产品来实现。为了描述的方便，描述以上装置时以功能分为各种模块分别描述。当然，在实施本说明书时可以把各模块的功能在同一个或多个软件和/或硬件中实现，也可以将实现同一功能的模块由多个子模块或子单元的组合实现等。以上所描述的装置实施例仅仅是示意性的，例如，所述单元的划分，仅仅为一种逻辑功能划分，实际实现时可以有另外的划分方式，例如多个单元或组件可以结合或者可以集成到另一个系统，或一些特征可以忽略，或不执行。另一点，所显示或讨论的相互之间的耦合或直接耦合或通信连接可以是通过一些接口，装置或单元的间接耦合或通信连接，可以是电性，机械或其它的形式。
由上可见，本说明书实施例提供的地震反射系数的确定装置，可以不需要对反射系数的稀疏域进行离散化，能够有效地避免出现基失配的问题，从而可以准确地反演得到精度较高的反射系数，解决了现有方法中存在的所确定的反射系数误差大、精度低的技术问题。
在一个具体的场景示例中，可以应用本说明书提供的地震反射系数的确定方法进行反演测试。
通常为了进行上述反演的基本思路是求解反射系数L1范数的极小值，对于这一非线性的优化问题，稀疏反演中常用的基追踪，匹配追踪等算法通常将该问题转化为线性规划问题进行求解。而基于本说明书所提供的地震反射系数的确定方法则可以将半定规划理论应用于地震稀疏反演领域，在求解L1范数最小化的稀疏约束问题时通过求其对应的半定规划问题避免了稀疏域离散化，从而解决稀疏反射系数反演中的基失配问题，获得精度更高的反演结果。
具体的处理流程可以参阅图4所示。具体实施时，可以包括以下流程步骤：
S1：根据原始数据(例如，叠后地震资料)和地震子波褶积矩阵(例如，地震子波褶积模型)建立反射系数反演目标方程(即，第一目标函数)；
S2：根据目标方程列出其对偶方程以及相应的近似方程(即，第二目标函数)；
S3：利用cvx程序求解对偶近似方程的解(即，得到目标系数)；
S4：根据所求方程解计算反射系数的多项式的解，既求得反射系数时间(记为反射时间)；
S5：根据所求反射时间构建准确的傅里叶矩阵，利用最小二乘法求反射系数幅值。
在本场景示例中，具体实施时，可以分别选取二维层状模型(可以参阅图5所示)，以及大小为201*201*251的实际三维数据体(可以参阅图9所示)应用本方法进行反射系数反演测试。其中，图5中的二维层状模型具体可以采用带有楔形结构的二维稀疏反射系数模型，并且使反射系数随机分布。其中，图5、图6、图7、图8中的CDP(Common Depth Point)表示CDP道集，即具有共同深度反射点的相应各记录道组成共深度点(或共反射点)道集。
其中，图6中所展示的二维地震记录可以理解为图5中的二维层状模型与30Hz雷克子波褶积得到的合成地震记录，并加入了10％随机噪声。
为了能够更好地比较出本方法与现有方法的区别，可以先使用现有方法，利用基追踪算法对图5中的二维层状模型进行反演，得到相应的结果，可以参阅图7所示。
同时，在利用本方法，按照上述流程步骤，通过半定规划的方法对图5中的二维层状模型进行反演，得到相应的结果，可以参阅图8所示。
图9展示出了某地区实际叠后地震资料。图10展示出了反演时所用地震子波。图11具体展示出了利用半定规划方法反演图9中的叠后地震资料所得到的结果。
通过对比图7和图8所展示的内容可以发现：虽然用基追踪算法得到的反演结果也呈现出原始模型的地层格架，但是基于半定规划理论得到的结果能够得到更清晰的反射时间和更准确的振幅，并且保持了良好的横向连续性和较强的抗噪性。
根据图11所示的三维数据反演结果可以看出：图9中无法识别的薄层也可以清晰地显示在图11所示的高分辨率剖面中，有效地提高了反演的分辨率。可见基于本说明书提供的方法，经半定规划反演得到的反射系数剖面是连续的，能反映更多的地层信息，保留了原始地震剖面的地层格架。
通过上述场景示例，验证了本说明书提供的地震反射系数的确定方法可以不需要对反射系数的稀疏域进行离散化，能够有效地避免出现基失配的问题，从而可以准确地反演得到精度较高的反射系数，解决了现有方法中存在的所确定的反射系数误差大、精度低的技术问题。
虽然本说明书提供了如实施例或流程图所述的方法操作步骤，但基于常规或者无创造性的手段可以包括更多或者更少的操作步骤。实施例中列举的步骤顺序仅仅为众多步骤执行顺序中的一种方式，不代表唯一的执行顺序。在实际中的装置或客户端产品执行时，可以按照实施例或者附图所示的方法顺序执行或者并行执行(例如并行处理器或者多线程处理的环境，甚至为分布式数据处理环境)。术语“包括”、“包含”或者其任何其他变体意在涵盖非排他性的包含，从而使得包括一系列要素的过程、方法、产品或者设备不仅包括那些要素，而且还包括没有明确列出的其他要素，或者是还包括为这种过程、方法、产品或者设备所固有的要素。在没有更多限制的情况下，并不排除在包括所述要素的过程、方法、产品或者设备中还存在另外的相同或等同要素。第一，第二等词语用来表示名称，而并不表示任何特定的顺序。
本领域技术人员也知道，除了以纯计算机可读程序代码方式实现控制器以外，完全可以通过将方法步骤进行逻辑编程来使得控制器以逻辑门、开关、专用集成电路、可编程逻辑控制器和嵌入微控制器等的形式来实现相同功能。因此这种控制器可以被认为是一种硬件部件，而对其内部包括的用于实现各种功能的装置也可以视为硬件部件内的结构。或者甚至，可以将用于实现各种功能的装置视为既可以是实现方法的软件模块又可以是硬件部件内的结构。
本说明书可以在由计算机执行的计算机可执行指令的一般上下文中描述，例如程序模块。一般地，程序模块包括执行特定任务或实现特定抽象数据类型的例程、程序、对象、组件、数据结构、类等等。也可以在分布式计算环境中实践本说明书，在这些分布式计算环境中，由通过通信网络而被连接的远程处理设备来执行任务。在分布式计算环境中，程序模块可以位于包括存储设备在内的本地和远程计算机存储介质中。
通过以上的实施例的描述可知，本领域的技术人员可以清楚地了解到本说明书可借助软件加必需的通用硬件平台的方式来实现。基于这样的理解，本说明书的技术方案本质上可以以软件产品的形式体现出来，该计算机软件产品可以存储在存储介质中，如ROM/RAM、磁碟、光盘等，包括若干指令用以使得一台计算机设备(可以是个人计算机，移动终端，服务器，或者网络设备等)执行本说明书各个实施例或者实施例的某些部分所述的方法。
本说明书中的各个实施例采用递进的方式描述，各个实施例之间相同或相似的部分互相参见即可，每个实施例重点说明的都是与其他实施例的不同之处。本说明书可用于众多通用或专用的计算机系统环境或配置中。例如：个人计算机、服务器计算机、手持设备或便携式设备、平板型设备、多处理器系统、基于微处理器的系统、置顶盒、可编程的电子设备、网络PC、小型计算机、大型计算机、包括以上任何系统或设备的分布式计算环境等等。
虽然通过实施例描绘了本说明书，本领域普通技术人员知道，本说明书有许多变形和变化而不脱离本说明书的精神，希望所附的权利要求包括这些变形和变化而不脱离本说明书的精神。

Claims (10)

1.一种地震反射系数的确定方法，其特征在于，包括：
获取目标区域的叠后地震资料；
根据所述叠后地震资料和地震子波褶积模型，构建关于目标区域的地震反射系数的第一目标函数；
根据所述第一目标函数，构建对应的对偶近似函数，作为第二目标函数；
通过基于半定规划理论的预设处理方式处理所述第二目标函数，以确定出反射时间；
根据所述反射时间，确定出目标区域的地震反射系数。
2.根据权利要求1所述的方法，其特征在于，通过基于半定规划理论的预设处理方式处理所述第二目标函数，以确定出反射时间，包括：
获取并利用Hermite矩阵，将第二目标函数转换为基于半定规划理论的函数，作为第三目标函数；
根据所述第三目标函数，确定出反射时间。
3.根据权利要求2所述的方法，其特征在于，获取并利用Hermite矩阵，将第二目标函数转换为基于半定规划理论的函数，作为第三目标函数，包括：构建以下算式作为第三目标函数：
Figure FDA0002758068060000011
基于
Figure FDA0002758068060000012
其中，Q为Hermite矩阵，c为目标系数，i为矩阵的行编号，j为矩阵的列编号，n为矩阵的大小。
4.根据权利要求3所述的方法，其特征在于，在获取并利用Hermite矩阵，将第二目标函数转换为基于半定规划理论的函数，作为第三目标函数之后，所述方法还包括：
根据目标区域的叠后地震资料，确定噪声水平；
根据所述噪声水平，确定出正则化调整参数；
根据所述正则化调整参数和所述第三目标函数，构建第四目标函数；
相应的，根据所述第四目标函数，确定出反射时间。
5.根据权利要求4所述的方法，其特征在于，根据所述正则化调整参数和所述第三目标函数，构建第四目标函数，包括：构建以下算式作为第四目标函数：
Figure FDA0002758068060000021
基于
Figure FDA0002758068060000022
其中，λ为正则化调整参数。
6.根据权利要求5所述的方法，其特征在于，根据所述第四目标函数，确定出反射时间，包括：
求解所述第四目标函数，确定出目标系数；
根据所述目标系数和约束函数的多项式，计算中间参数；
根据所述中间参数，确定所述反射时间。
7.根据权利要求1所述的方法，其特征在于，根据所述反射时间，确定出目标区域的地震反射系数，包括：
将所述反射时间代入基于傅里叶变换的反射系数的表征函数，并使用最小二乘法求解所述表征函数，以得到对应的反射幅值。
8.根据权利要求1所述的方法，其特征在于，在根据所述反射时间，确定出目标区域的地震反射系数之后，所述方法还包括：
根据所述目标区域的地震反射系数，指导对目标区域进行地震勘探。
9.一种地震反射系数的确定装置，其特征在于，包括：
获取模块，用于获取目标区域的叠后地震资料；
第一构建模块，用于根据所述叠后地震资料和地震子波褶积模型，构建关于目标区域的地震反射系数的第一目标函数；
第二构建模块，用于根据所述第一目标函数，构建对应的对偶近似函数，作为第二目标函数；
处理模块，用于通过基于半定规划理论的预设处理方式处理所述第二目标函数，以确定出反射时间；
确定模块，用于根据所述反射时间，确定出目标区域的地震反射系数。
10.一种服务器，其特征在于，包括处理器以及用于存储处理器可执行指令的存储器，所述处理器执行所述指令时实现权利要求1至8中任一项所述方法的步骤。
CN202011208807.6A 2020-11-03 2020-11-03 地震反射系数的确定方法和装置 Pending CN114442157A (zh)

Priority Applications (1)

Application Number Priority Date Filing Date Title
CN202011208807.6A CN114442157A (zh) 2020-11-03 2020-11-03 地震反射系数的确定方法和装置

Applications Claiming Priority (1)

Application Number Priority Date Filing Date Title
CN202011208807.6A CN114442157A (zh) 2020-11-03 2020-11-03 地震反射系数的确定方法和装置

Publications (1)

Publication Number Publication Date
CN114442157A true CN114442157A (zh) 2022-05-06

Family

ID=81361801

Family Applications (1)

Application Number Title Priority Date Filing Date
CN202011208807.6A Pending CN114442157A (zh) 2020-11-03 2020-11-03 地震反射系数的确定方法和装置

Country Status (1)

Country Link
CN (1) CN114442157A (zh)

Similar Documents

Publication Publication Date Title
Ng et al. Chaos based analytical techniques for daily extreme hydrological observations
US10436924B2 (en) Denoising seismic data
Ghaderpour et al. LSWAVE: A MATLAB software for the least-squares wavelet and cross-wavelet analyses
CN110954950B (zh) 地下横波速度反演方法、装置、计算设备及存储介质
Rosenberger Real-time ground-motion analysis: Distinguishing P and S arrivals in a noisy environment
CN112882100B (zh) 一种储层参数确定方法、装置、电子设备和存储介质
CN112213782B (zh) 分相位地震数据的处理方法、装置和服务器
Meng et al. Estimation of chirp signals with time-varying amplitudes
CN110879412A (zh) 地下横波速度反演方法、装置、计算设备及存储介质
CN115795282A (zh) 激波管动态压力重构方法、装置、电子设备及存储介质
CN116204768A (zh) 故障诊断方法、装置、计算机设备和存储介质
Yablokov et al. Uncertainty quantification of multimodal surface wave inversion using artificial neural networks
CN114442157A (zh) 地震反射系数的确定方法和装置
CN113970787B (zh) 物性参数反演方法、装置、计算机设备和存储介质
CN115980847A (zh) 一种含气性预测方法、装置、计算机设备及存储介质
CN111337972A (zh) 用于地震数据重建的规则化域选取方法、系统
CN111077573A (zh) 一种确定地层弹性参数的方法、装置及系统
CN113805225A (zh) 一种相位约束的高分辨率地震反演方法及计算机设备
CN108599773A (zh) 一种基于确定性测量矩阵的振动信号数据压缩采集方法
CN116299701B (zh) 断溶体储层识别方法、装置及电子设备
CN117849875B (zh) 一种地震信号分析方法、系统、装置及存储介质
Lv Sparse dictionary learning for noise attenuation in the exactly flattened dimension
CN115542394A (zh) 一种用谱反演提高地震分辨率方法
CN111948706B (zh) 正交各向异性介质地震成像方法及装置
Su et al. Adaptive approach based on curve fitting and interpolation for boundary effects reduction

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