绕射波场获取方法和装置与流程

本发明涉及地震勘探的建模技术领域，尤其是涉及一种绕射波场获取方法和装置。

背景技术：

地层中存在断层、小断块、裂缝和陷落柱等非均匀不连续地质体，这些地质体附近往往可能存储有煤炭、石油等矿产资源，因此精确定位这些地质体是矿产资源高效开采的有力保障。然而，由于断层的断楞和裂缝等这些突变点在空间展布上通常小于地震子波波长，因此，地震波在传播过程中如遇到上述突变点则会引起绕射波。此时，绕射波将作为由小尺度地质异常体引起的地震响应而成为非均匀不连续地质体的可靠信息载体，其中，绕射波具有着高分辨探测能力，对定位非均匀不连续地质体具有重要意义。

基于绕射波的优点，相关技术人员开启了对绕射波的研究，其中，在对绕射波进行研究时，主要是在野外采集地震信号，从而对地震信号进行分析，得到绕射波。通常情况下，野外作业过程中检波器接收到的地震信号不仅包含绕射波，还包含反射波。反射波成像条件是假定无限大光滑界面，故一般只能反映大尺度地质背景。此外，绕射波属于异常波中的一种，能量衰减较快，地质体产生的绕射能量约为反射能量0.1-0.01倍，为弱信号。所以绕射波搭载的信号常淹没在反射波成像结果中，因此需要去除反射波。

现有技术中，平面波分解是最常用的去除反射波的方法之一。平面波分解是利用平面波解构滤波器对反射波场进行预测滤波，从而达到压制反射波场、提取绕射波场的目的。其中，平面波分解方法仅仅依赖于反射波场局部倾角值，局部倾角的预测直接决定了所提取的绕射波场的精度。由于平面波分解方法没有考虑到反射波倾角场稀疏的特性，导致所预测的反射波倾角场精度降低，使反射波场不能得到较好的压制，绕射波场精度也有所降低。

技术实现要素：

本发明实施例提供了一种绕射波场获取方法和装置，缓解了现有技术中在确定绕射波场时，精确度较差的技术问题。

根据本发明实施例的一个方面，提供了一种绕射波场获取方法，包括：获取目标区域内的共偏移距地震道集数据，其中，所述共偏移距地震道集数据为携带有地层界面信息的数据；根据所述共偏移距地震道集数据，在平面波分解技术的基础上，通过L1范数正则化算法构建目标函数；采用预设迭代算法对所述目标函数达到最小值的可行解进行迭代求解，并将所述可行解作为目标反射波倾角场；根据所述目标反射波倾角场和所述共偏移距地震道集数据确定绕射波场。

进一步地，所述根据所述共偏移距地震道集数据，在平面波分解技术的基础上，通过L1范数正则化算法构建目标函数包括：利用平面波分解技术，根据所述共偏移距地震道集数据建立数据拟合项；对反射波倾角场通过稀疏变换方法变换之后，将稀疏变换系数中大于所述预设系数分界值的系数作为高频系数予以保留，再通过L1范数对所述高频系数进行正则化约束而建立正则化项；将所述数据拟合项和所述正则化项相加构建所述目标函数。

进一步地，所述稀疏变换方法为小波变换方法，所述目标函数构建如下：

其中，E(σ)表示所述目标函数，E₁为所述数据拟合项的滤波结果，E₂为所述正则化项的数值，d为所述共偏移距地震道集数据，σ为所述反射波倾角场，α为尺度化调节因子，||||₂表示L2范数，||||₁表示L1范数，Φ(σ)表示所述高频系数，C为利用平面波分解方法构建的平面波解构滤波器。

进一步地，所述采用预设迭代算法对所述目标函数达到最小值的可行解进行迭代求解，并将所述可行解作为目标反射波倾角场包括：

变换步骤，通过Split Bregman算法对所述目标函数达到最小值的方程进行变换，得到变换方程min表示最小化；

分解步骤，通过Split Bregman算法对所述变换方程进行分解，分解得到第一变量的迭代公式第二变量的迭代公式和第三变量的迭代公式b^k+1＝b^k+Φ(σ^k+1)-p^k+1，λ为所述Split Bregman算法中的参数；

第一计算步骤，通过共轭梯度算法对所述第一变量的迭代公式进行求解，得到所述反射波倾角场的迭代值；

第二计算步骤，通过shrink方法对所述第二变量的迭代公式进行求解，计算出所述高频系数的迭代值；

第三计算步骤，将所述反射波倾角场的迭代值和所述高频系数的迭代值代入至所述第三变量的迭代公式中，计算得到所述第三变量b的迭代值；

判断步骤，判断当前迭代终止条件是否满足预设迭代终止条件，所述预设迭代终止条件包括：迭代次数等于预设次数，或者，所述共轭梯度算法中的梯度值小于预设梯度值，如果判断出不满足所述预设迭代终止条件，返回执行所述第一计算步骤，如果判断出满足所述预设迭代终止条件，则终止迭代；

平滑处理步骤，采用三角平滑算法对所述反射波倾角场的迭代值进行平滑处理，并将平滑处理之后的所述反射波倾角场的迭代值作为所述目标反射波倾角场。

进一步地，根据所述目标反射波倾角场和所述共偏移距地震道集数据确定绕射波场包括：将所述目标反射波倾角场代入至公式C(σ^*)d中进行滤波运算，其中，d为所述共偏移距地震道集数据，σ^*为所述目标反射波倾角场，C为利用平面波分解方法构建的平面波解构滤波器。

根据本发明实施例的另一方面，还提供了一种绕射波场获取装置，包括：获取单元，用于获取目标区域内的共偏移距地震道集数据；构建单元，用于根据所述共偏移距地震道集数据，在平面波分解技术的基础上，通过L1范数正则化算法构建目标函数；求解单元，用于采用预设迭代算法对所述目标函数达到最小值的可行解进行迭代求解，并将所述可行解作为目标反射波倾角场；确定单元，用于根据所述目标反射波倾角场和所述共偏移距地震道集数据确定绕射波场。

进一步地，所述根据所述共偏移距地震道集数据，在平面波分解技术的基础上，通过L1范数正则化算法构建目标函数包括：利用平面波分解技术，根据所述共偏移距地震道集数据建立数据拟合项；对反射波倾角场通过稀疏变换方法变换之后，将稀疏变换系数中大于所述预设系数分界值的系数作为高频系数予以保留，再通过L1范数对所述高频系数进行正则化约束而建立正则化项；将所述数据拟合项和所述正则化项相加，并将相加的结果作为所述目标函数。

进一步地，所述稀疏变换方法为小波变换方法，所述构建单元构建的目标函数如下：

其中，E(σ)表示所述目标函数，E₁为所述数据拟合项的滤波结果，E₂为所述正则化项的数值，d为所述共偏移距地震道集数据，σ为所述反射波倾角场，α为尺度化调节因子，||||₂表示L2范数，||||₁表示L1范数，Ф(σ)表示所述高频系数，C为利用平面波分解方法构建的平面波解构滤波器。

进一步地，所述求解单元包括：

变换模块，通过Split Bregman算法对所述目标函数达到最小值的方程进行变换，得到变换方程min表示最小化；

分解模块，通过Split Bregman算法对所述变换方程进行分解，分解得到第一变量的迭代公式第二变量的迭代公式和第三变量的迭代公式b^k+1＝b^k+Φ(σ^k+1)-p^k+1，λ为所述Split Bregman算法中的参数；

第一计算模块，通过共轭梯度算法对所述第一变量的迭代公式进行求解，得到所述反射波倾角场的迭代值；

第二计算模块，通过shrink方法对所述第二变量的迭代公式进行求解，计算出所述高频系数的迭代值；

第三计算模块，将所述反射波倾角场的迭代值和所述高频系数的迭代值代入至所述第三变量的迭代公式中，计算得到所述第三变量b的迭代值；

判断模块，判断当前迭代终止条件是否满足预设迭代终止条件，所述预设迭代终止条件包括：迭代次数等于预设次数，或者，所述共轭梯度算法中的梯度值小于预设梯度值，如果判断出不满足所述预设迭代终止条件，返回执行所述第一计算步骤，如果判断出满足所述预设迭代终止条件，则终止迭代；

平滑处理模块，采用三角平滑算法对所述反射波倾角场的迭代值进行平滑处理，并将平滑处理之后的所述反射波倾角场的迭代值作为所述目标反射波倾角场。

进一步地，所述确定模块根据所述目标反射波倾角场和所述共偏移距地震道集数据确定绕射波场包括：将所述目标反射波倾角场代入至公式C(σ^*)d中进行滤波运算，其中，d为所述共偏移距地震道集数据，σ^*为所述目标反射波倾角场，C为利用平面波分解方法构建的平面波解构滤波器。

在本发明实施例中，在获取目标区域内的共偏移距地震道集数据后，根据所述共偏移距地震道集数据，在平面波分解技术的基础上，通过L1范数正则化算法构建目标函数，继而采用预设迭代算法对所述目标函数达到最小值的可行解进行迭代求解，并将所述可行解作为目标反射波倾角场，最后根据所述目标反射波倾角场和所述共偏移距地震道集数据确定绕射波场。本发明提供的绕射波场获取方法中，在平面波分解技术的基础上，通过L1范数正则化算法构建目标函数，L1范数正则化具有特征选择能力，可以产生稀疏解，从而使得反射波倾角场中一些低频分量的输出权值趋近于零，即限定了反射波倾角场的稀疏性，使得求解出的目标反射波倾角场更为符合实际地质情况，进而使得根据所述目标反射波倾角场和所述共偏移距地震道集数据确定的绕射波场更加精确。

附图说明

为了更清楚地说明本发明具体实施方式或现有技术中的技术方案，下面将对具体实施方式或现有技术描述中所需要使用的附图作简单地介绍，显而易见地，下面描述中的附图是本发明的一些实施方式，对于本领域普通技术人员来讲，在不付出创造性劳动的前提下，还可以根据这些附图获得其他的附图。

图1为本发明一种绕射波场获取方法的流程图；

图2为本发明采用预设迭代算法对所述目标函数达到最小值的可行解进行迭代求解的流程图；

图3为本发明一种绕射波场获取装置的结构示意图。

图标：302-获取单元；304-构建单元；306-求解单元；308-确定单元。

具体实施方式

为了使本领域的人员更好地理解本发明方案，下面将结合本发明实施例中的附图对本发明的技术方案进行清楚、完整地描述，显然，所描述的实施例是本发明一部分实施例，而不是全部的实施例。基于本发明中的实施例，本领域普通技术人员在没有做出创造性劳动前提下所获得的所有其他实施例，都属于本发明保护的范围。

如图1所示，一种绕射波场获取方法，包括以下步骤：

步骤S102，获取目标区域内的共偏移距地震道集数据，其中，共偏移距地震道集数据为携带有地层界面信息的数据。

在本发明实施例中，共偏移距地震道集数据为通过获取单元302检测到的地震波，具体过程为：采用人工方法在炮点激发了地震波，地震波向四面八方传播，并且在遇到地下不同岩性的分解面时产生了反射波和绕射波，这些反射波和绕射波又返回地面而引起地面的振动，然后沿与炮点等距的地方设置获取单元302来检测引起地面振动的地震波。

需要说明的是，由于检测到的地震波是受到了地下地层介质的改造，因此地层界面信息至少包括地质构造和地层岩性，携带有地层界面信息的数据包括时间、速度、频率。

步骤S104，根据共偏移距地震道集数据，在平面波分解技术的基础上，通过L1范数正则化算法构建目标函数。

在本发明实施例中，目标函数是基于平面波分解方法，将反射波倾角场作为因变量，通过L1范数正则化算法而构建的，L1范数正则化算法对反射波倾角场进行了约束，由于平面波分解方法中所依赖的反射波倾角场具有稀疏的特性，L1范数正则化算法限定了反射波倾角场的稀疏性，从而使得反射波倾角场中一些低频分量的输出权值趋近于零，即使得反射波倾角场得到了更加精确的表达。

步骤S106，采用预设迭代算法对目标函数达到最小值的可行解进行迭代求解，并将可行解作为目标反射波倾角场。

在本发明实施例中，由于反射波倾角场是目标函数的因变量，因此目标函数达到最小值的可行解即为反射波倾角场的一个值，将反射波倾角场的这个值作为目标反射波倾角场。

步骤S108，根据目标反射波倾角场和共偏移距地震道集数据确定绕射波场。

在本发明实施例中，将目标反射波倾角场和共偏移距地震道集数据代入平面波分解方法的计算公式中即可确定绕射波场。

在本发明实施例中，在获取目标区域内的共偏移距地震道集数据后，根据共偏移距地震道集数据，在平面波分解技术的基础上，通过L1范数正则化算法构建目标函数，继而采用预设迭代算法对目标函数达到最小值的可行解进行迭代求解，并将可行解作为目标反射波倾角场，最后根据目标反射波倾角场和共偏移距地震道集数据确定绕射波场。本发明提供的绕射波场获取方法中，在平面波分解技术的基础上，通过L1范数正则化算法构建目标函数，L1范数正则化具有特征选择能力，可以产生稀疏解，从而使得反射波倾角场中一些低频分量的输出权值趋近于零，即限定了反射波倾角场的稀疏性，使得求解出的目标反射波倾角场更为符合实际地质情况，进而达到了根据目标反射波倾角场和共偏移距地震道集数据精确确定的绕射波场的目的，并缓解了现有技术中在确定绕射波场时，精确度较差的技术问题。

在本发明的一个可选实施方式中，根据共偏移距地震道集数据，在平面波分解技术的基础上，通过L1范数正则化算法构建目标函数包括：利用平面波分解技术，根据共偏移距地震道集数据建立数据拟合项；对反射波倾角场通过稀疏变换方法变换之后，将稀疏变换系数中大于所述预设系数分界值的系数作为高频系数予以保留，再通过L1范数对所述高频系数进行正则化约束而建立正则化项，这里通过预设系数分界值是先对稀疏系数进行了初步的高低频分界然后提取了高频系数；最后将数据拟合项和正则化项相加，并将相加的结果作为目标函数。

在本发明的另一个可选实施方式中，具体地，稀疏变换方法为小波变换方法，构建目标函数如下：

其中，E(σ)表示目标函数，E₁为数据拟合项的滤波结果，E₂为正则化项的数值，d为共偏移距地震道集数据，σ为反射波倾角场，α为尺度化调节因子，||||₂表示L2范数，||||₁表示L1范数，Φ(σ)表示所述高频系数，C为利用平面波分解方法构建的平面波解构滤波器。在小波变换方法中，小波变换系数的表达式为：

其中，a为伸缩因子，b为平移因子，f(t)为要进行小波变换的原始信号，为小波基函数。此外，C的表达式为：

其中，I表示单位算子，P_i,j为从i道至j道的预测误差算子。

如图2所示，在本发明的另一个可选实施方式中，采用预设迭代算法对目标函数达到最小值的可行解进行迭代求解的一种执行方案包括：

S201，变换步骤，通过Split Bregman算法对目标函数达到最小值的方程进行变换，得到变换方程min表示最小化；

S203，分解步骤，为提高迭代效率，引入变量b，S201中变换方程的求解等价于：

通过Split Bregman算法，上式分解为：

第一变量的迭代公式：

第二变量的迭代公式：

第三变量的迭代公式：b^k+1＝b^k+Φ(σ^k+1)-p^k+1；

其中，λ为Split Bregman算法中涉及的参数，具体地，所引入的变量b还获取初始值b⁰；

S205，第一计算步骤，通过共轭梯度算法对第一变量的迭代公式进行求解，得到反射波倾角场的迭代值，具体地，计算第一变量的迭代公式的等式右侧对反射波倾角场的梯度，共轭梯度法中反射波倾角场迭代更新方向为：

其中，g表示第一迭代公式对反射波倾角场的梯度，q为更新方向，β为共轭调节系数，σ⁰为满足共轭梯度算法的需要而获取的反射波倾角场的初始值，令η为步长更新系数，反射波倾角场更新公式为：

σ^k+1＝σ^k+η^kq^k；

S207，第二计算步骤，通过shrink方法对第二变量的迭代公式进行求解，计算出高频系数的迭代值，若定义shrinkage函数如下：

则利用shrink方法求解公式如下：

第一次执行第二计算步骤时，通过b⁰和σ¹求得p¹。

S209，第三计算步骤，将反射波倾角场的迭代值和高频系数的迭代值代入至第三变量的迭代公式中，计算得到第三变量b的迭代值。

S211，判断步骤，判断当前迭代终止条件是否满足预设迭代终止条件，所述预设迭代终止条件包括：迭代次数等于预设次数，或者，所述共轭梯度算法中的梯度值小于预设梯度值，预设次数和预设剃度值都是预先获取的已知值，如果判断出不满足所述预设迭代终止条件，返回执行所述第一计算步骤，此次迭代过程中求得的反射波倾角场的迭代值、高频系数的迭代值和第三变量b的迭代值代入第一计算步骤计算下次迭代过程中的反射波倾角场的迭代值，如果判断出满足所述预设迭代终止条件，则终止迭代；

S213，平滑处理步骤，采用三角平滑算法对最后一次迭代过程中第一计算步骤求得的所述反射波倾角场的迭代值进行平滑处理，并将平滑处理之后的所述反射波倾角场的迭代值作为所述目标反射波倾角场。

具体地，平滑处理步骤的实施过程为：首先，对所述反射波倾角场的迭代值实施矩形窗平滑B_k(Z)；接着，实施矩形窗平滑B_k(1/Z)，从而达到三角平滑的目的。其中，B_k(Z)和B_k(1/Z)为矩形窗平滑在Z变换中的表达式，矩形窗平滑在Z变换中的具体定义为而三角平滑定义为T_k(Z)＝B_k(1/Z)B_k(Z)，X(Z)表示需要平滑处理的信号，k表示滤波器长度。

在本发明的另一个可选实施方式中，根据目标反射波倾角场和共偏移距地震道集数据确定绕射波场包括：将目标反射波倾角场代入至公式C(σ^*)d中进行滤波运算，其中，d为共偏移距地震道集数据，σ^*为目标反射波倾角场，C为利用平面波分解方法构建的平面波解构滤波器。

本发明实施例还提供了一种绕射波场获取装置，该绕射波场获取装置主要用于执行本发明实施例上述内容所提供的绕射波场获取方法，以下对本发明实施例提供的绕射波场获取装置做具体介绍。

图3是本发明一种可选实施例提供的绕射波场获取装置的示意图，如图3所示，该绕射波场获取装置主要包括获取单元302、构建单元304、求解单元306和确定单元308，其中：

获取单元302，用于获取目标区域内的共偏移距地震道集数据；

构建单元304，用于根据共偏移距地震道集数据，在平面波分解技术的基础上，通过L1范数正则化算法构建目标函数；

求解单元306，用于采用预设迭代算法对目标函数达到最小值的可行解进行迭代求解，并将可行解作为目标反射波倾角场；

确定单元308，用于根据目标反射波倾角场和共偏移距地震道集数据确定绕射波场。

在本发明实施例中，在获取单元302获取目标区域内的共偏移距地震道集数据后，构建单元304根据共偏移距地震道集数据，在平面波分解技术的基础上，通过L1范数正则化算法构建目标函数，求解单元306继而采用预设迭代算法对目标函数达到最小值的可行解进行迭代求解，并将可行解作为目标反射波倾角场，最后确定单元308根据目标反射波倾角场和共偏移距地震道集数据确定绕射波场。本发明提供的绕射波场获取装置中，构建单元304在平面波分解技术的基础上，通过L1范数正则化算法构建目标函数，L1范数正则化具有特征选择能力，可以产生稀疏解，从而使得反射波倾角场中一些低频分量的输出权值趋近于零，即限定了反射波倾角场的稀疏性，使得求解出的目标反射波倾角场更为符合实际地质情况，进而达到了根据目标反射波倾角场和共偏移距地震道集数据精确确定的绕射波场的目的，并缓解了现有技术中在确定绕射波场时，精确度较差的技术问题。

可选地，所述构建单元304根据所述共偏移距地震道集数据，在平面波分解技术的基础上，通过L1范数正则化算法构建目标函数包括：利用平面波分解技术，根据所述共偏移距地震道集数据建立数据拟合项；对反射波倾角场通过稀疏变换方法变换之后，将稀疏变换系数中大于所述预设系数分界值的系数作为高频系数予以保留，再通过L1范数对所述高频系数进行正则化约束而建立正则化项；将所述数据拟合项和所述正则化项相加，并将相加的结果作为所述目标函数。

可选地，稀疏变换方法为小波变换方法，构建单元304构建的目标函数如下：

其中，E(σ)表示目标函数，E₁为数据拟合项的滤波结果，E₂为正则化项的数值，d为共偏移距地震道集数据，σ为反射波倾角场，α为尺度化调节因子，||||₂表示L2范数，||||₁表示L1范数，Ф(σ)表示所述高频系数，C为利用平面波分解方法构建的平面波解构滤波器。

可选地，求解单元306包括：

变换模块，通过Split Bregman算法对所述目标函数达到最小值的方程进行变换，得到变换方程min表示最小化；

分解模块，通过Split Bregman算法对所述变换方程进行分解，分解得到第一变量的迭代公式第二变量的迭代公式和第三变量的迭代公式b^k+1＝b^k+Φ(σ^k+1)-p^k+1，λ为所述Split Bregman算法中的参数；

第一计算模块，通过共轭梯度算法对所述第一变量的迭代公式进行求解，得到所述反射波倾角场的迭代值；

第二计算模块，通过shrink方法对所述第二变量的迭代公式进行求解，计算出所述高频系数的迭代值；

第三计算模块，将所述反射波倾角场的迭代值和所述高频系数的迭代值代入至所述第三变量的迭代公式中，计算得到所述第三变量b的迭代值；

判断模块，判断当前迭代终止条件是否满足预设迭代终止条件，所述预设迭代终止条件包括：迭代次数等于预设次数，或者，所述共轭梯度算法中的梯度值小于预设梯度值，如果判断出不满足所述预设迭代终止条件，返回执行所述第一计算步骤，如果判断出满足所述预设迭代终止条件，则终止迭代；

平滑处理模块，采用三角平滑算法对所述反射波倾角场的迭代值进行平滑处理，并将平滑处理之后的所述反射波倾角场的迭代值作为所述目标反射波倾角场。

可选地，所述确定模块308根据所述目标反射波倾角场和所述共偏移距地震道集数据确定绕射波场包括：将所述目标反射波倾角场代入至公式C(σ^*)d中进行滤波运算，其中，d为所述共偏移距地震道集数据，σ^*为所述目标反射波倾角场，C为利用平面波分解方法构建的平面波解构滤波器。

在本发明的描述中，需要说明的是，除非另有明确的规定和限定，术语“安装”、“相连”、“连接”应做广义理解，例如，可以是固定连接，也可以是可拆卸连接，或一体地连接；可以是机械连接，也可以是电连接；可以是直接相连，也可以通过中间媒介间接相连，可以是两个元件内部的连通。对于本领域的普通技术人员而言，可以具体情况理解上述术语在本发明中的具体含义。

所述功能如果以软件功能单元的形式实现并作为独立的产品销售或使用时，可以存储在一个计算机可读取存储介质中。基于这样的理解，本发明的技术方案本质上或者说对现有技术做出贡献的部分或者该技术方案的部分可以以软件产品的形式体现出来，该计算机软件产品存储在一个存储介质中，包括若干指令用以使得一台计算机设备(可以是个人计算机，服务器，或者网络设备等)执行本发明各个实施例所述方法的全部或部分步骤。而前述的存储介质包括：U盘、移动硬盘、只读存储器(ROM，Read-Only Memory)、随机存取存储器(RAM，Random Access Memory)、磁碟或者光盘等各种可以存储程序代码的介质。

在发明的描述中，需要说明的是，术语“中心”、“上”、“下”、“左”、“右”、“竖直”、“水平”、“内”、“外”等指示的方位或位置关系为基于附图所示的方位或位置关系，仅是为了便于描述本发明和简化描述，而不是指示或暗示所指的装置或元件必须具有特定的方位、以特定的方位构造和操作，因此不能理解为对本发明的限制。此外，术语“第一”、“第二”、“第三”仅用于描述目的，而不能理解为指示或暗示相对重要性。

最后应说明的是：以上所述实施例，仅为本发明的具体实施方式，用以说明本发明的技术方案，而非对其限制，本发明的保护范围并不局限于此，尽管参照前述实施例对本发明进行了详细的说明，本领域的普通技术人员应当理解：任何熟悉本技术领域的技术人员在本发明揭露的技术范围内，其依然可以对前述实施例所记载的技术方案进行修改或可轻易想到变化，或者对其中部分技术特征进行等同替换；而这些修改、变化或者替换，并不使相应技术方案的本质脱离本发明实施例技术方案的精神和范围，都应涵盖在本发明的保护范围之内。因此，本发明的保护范围应所述以权利要求的保护范围为准。