CN106772583A - 一种地震绕射波分离方法和装置 - Google Patents

Abstract

本发明提供了一种地震绕射波分离方法和装置，涉及勘探地震的技术领域，其中，该方法包括获取预设地质区域内携带有地下地质信息的地震炮集数据；将对该地震炮集数据进行预处理得到的单炮预处理数据和预设的偏移速度模型输入至三维单炮角度域成像公式中，对地震炮集数据进行波场反传处理，得到与预设地质区域内的地下成像点一一对应的传播射线的方位角、出射角和振幅信息；根据该传播射线的方位角、出射角和振幅信息，生成三维角度域成像矩阵；通过预设的三维绕射波分离模型，该三维角度域成像矩阵中分离出低秩矩阵分量，将该低秩矩阵分量确定为地震绕射波。本发明提高了分离出的绕射波振幅完整性和波形一致性，进而提高了地质结构成像的分辨率。

Description

一种地震绕射波分离方法和装置

技术领域

本发明涉及勘探地震技术领域，尤其是涉及一种地震绕射波分离方法和装置。

背景技术

碳酸盐岩油气藏已成为油气资源增储上产的主要领域。但是，部分地域的碳酸盐岩层系结构具有很大的特殊性，使得碳酸盐岩储层形成和分布较为复杂，导致岩溶溶洞、裂缝等地质体无法精细成像。

现有技术中，石油工业界地震勘探以反射波为主，但是通过反射波勘探分辨率有限，无法有效识别碳酸盐岩层系结构的地质体。同时，由于该碳酸盐岩层系结构的地质体地震响应为绕射波，因此有效分离绕射波对碳酸盐岩缝洞型油气藏勘探至关重要。常规的绕射波分离方法，大多利用反射波与绕射波运动学特征，通过信号处理方法分离绕射波。然而采集到的三维炮集数据中绕射波与反射波运动学特性相似度较高，仅依靠常规运动学波场分离方法难以有效处理，造成碳酸盐岩层系结构的成像分辨率较低。

针对上述地震绕射波分离方式分离出的绕射波振幅完整性和波形一致性较差的问题，目前尚未提出有效的解决方案。

发明内容

有鉴于此，本发明的目的在于提供一种地震绕射波分离方法和装置，以提高分离出的绕射波振幅完整性和波形一致性。

第一方面，本发明实施例提供了一种地震绕射波分离方法，包括：获取预设地质区域内携带有地下地质信息的地震炮集数据；其中，该地下地质信息包括地质构造信息和地质岩性变化信息；对地震炮集数据进行波场反传处理，得到与预设地质区域内的地下成像点一一对应的传播射线的方位角、出射角和振幅信息；根据该传播射线的方位角、出射角和振幅信息，生成三维角度域成像矩阵；从该三维角度域成像矩阵中分离出低秩矩阵分量，将该低秩矩阵分量确定为地震绕射波。

结合第一方面，本发明实施例提供了第一方面的第一种可能的实施方式，其中，上述对地震炮集数据进行波场反传处理，得到与预设地质区域内的地下成像点一一对应的传播射线的方位角、出射角和振幅信息包括：对地震炮集数据进行预处理，得到单炮预处理数据；其中，该单炮预处理数据为可用于直接成像的地震炮集数据；该预处理包括对地震炮集数据进行噪声去除处理，以及对地震炮集数据与预先保存的历史地震数据进行一一对应；将单炮预处理数据和预设的偏移速度模型输入至三维单炮角度域成像公式中，对地震炮集数据进行波场反传处理，得到与预设地质区域内的地下成像点一一对应的传播射线的方位角、出射角和振幅信息；其中，该三维单炮角度域成像公式中包含三维振幅补偿因子。

结合第一方面的第一种可能的实施方式，本发明实施例提供了第一方面的第二种可能的实施方式，其中，上述三维单炮角度域成像公式包括：

其中，δ为脉冲函数，代表三维角度域成像矩阵；震源s激发的射线经地下空间任意成像点x到达检波点位置r；矢量k_S为震源到成像点的射线参数，矢量k_r为检波点到成像点的射线参数；参数θ为出射角；参数为方位角；矢量k为假定反射界面的法向量；k通过下述方式计算获得：θ_s和分别为k_s的出射角和方位角；θ_r和分别为k_r的出射角和方位角；θ_m,分别为假定反射界面的出射角和方位角；n_x为三维坐标系x方向法向量，且n_x＝(1,0,0)；u(s,r，t)为地震数据，t为地震数据的记录时间；t₀为射线走时；W_3D(s,x,r)为三维振幅补偿因子。

结合第一方面的第二种可能的实施方式，本发明实施例提供了第一方面的第三种可能的实施方式，其中，上述三维振幅补偿因子W_3D(s,x,r)包括：

其中，v_s为震源位置速度，α_s为震源位置射线入射角，α_r为检波点位置射线出射角，N ₁,N ₂分别为第一支射线和第二支射线的走时对震源位置和检波点位置的混合导数，T表示矩阵转置操作，走时根据三维波前重建法射线追踪计算，该计算考虑多值走时；第一支射线为震源到成像点的射线；第二支射线为检波点到成像点的射线；Σ和Γ是与地震观测方式有关的矩阵，在共炮点观测情况下Σ＝0,Γ＝I；其中I为单位矩阵；i表示复数的虚数单位，κ₁,κ₂分别为第一支射线和第二支射线的焦散点数目，κ₁,κ₂根据三维射线追踪动力学方程计算获得。

结合第一方面，本发明实施例提供了第一方面的第四种可能的实施方式，其中，上述从三维角度域成像矩阵中分离出低秩矩阵分量，将低秩矩阵分量确定为地震绕射波包括：通过预设的三维绕射波分离模型，从三维角度域成像矩阵中分离出低秩矩阵分量，将低秩矩阵分量确定为地震绕射波；其中，该预设的三维绕射波分离模型包括：

其中，为第i个成像点x_i位置处的三维角度域成像矩阵；为三维角度域成像矩阵分解后的低秩矩阵分量；为三维角度域成像矩阵分解后的稀疏矩阵分量；x_i为第i个成像点；参数θ为出射角；参数为方位角。

结合第一方面的第四种可能的实施方式，本发明实施例提供了第一方面的第五种可能的实施方式，其中，上述预设的三维绕射波分离模型还包括：

其中，J(L,S,Y,β)为目标函数，Y为拉格朗日乘子矩阵，T表示矩阵转置操作，λ为正则化参数，β为保真惩罚因子，||·||_*为核范数，即矩阵奇异值之和，||·||₁为l₁范数，即矩阵中各元素绝对值之和，||·||_F为Frobenius范数，Frobenius范数为矩阵全部元素的平方和的平方根；L为三维角度域成像矩阵分解后的低秩矩阵分量；S为三维角度域成像矩阵分解后的稀疏矩阵分量；R为三维角度域成像矩阵。

结合第一方面的第五种可能的实施方式，本发明实施例提供了第一方面的第六种可能的实施方式，其中，上述从三维角度域成像矩阵中分离出低秩矩阵分量，将低秩矩阵分量确定为地震绕射波包括：设置正则化参数λ和预设的最大迭代次数N；其中，λ>0；设置迭代次数初始值k＝1、低秩矩阵分量的初始值L⁰、稀疏矩阵分量的初始值S⁰、拉格朗日乘子初始值Y⁰、以及保真惩罚因子初始值β⁰；将k＝1、L⁰、S⁰、Y⁰和β⁰作为初始值，对三维角度域成像矩阵进行迭代处理，迭代处理包括如下步骤：通过进行奇异值分解计算，得到奇异值对角阵；其中，R为三维角度域成像矩阵；U和V的列为基向量；∑为对角阵；奇异值对角阵的对角上的元素为奇异值；通过对奇异值对角阵中的奇异值a_i进行软门限阈值操作，得到新对角矩阵其中，x为预设的固定值；根据新对角矩阵计算低秩矩阵分量L^k和稀疏矩阵分量S^k；判断L^k和S^k是否满足关系式且k≤N；如果是，更新k＝k+1、拉格朗日乘子Y^k＝Y^k-1+β^k-1(R-L^k-S^k)、以及保真惩罚因子β^k＝ωβ^k-1(ω>0)其中ω为比例因子；并继续进行迭代处理；如果否，确定L^k为分离的地震绕射波。

结合第一方面的第六种可能的实施方式，本发明实施例提供了第一方面的第七种可能的实施方式，其中，上述根据新对角矩阵计算低秩矩阵分量L^k和稀疏矩阵分量S^k包括：根据新对角矩阵计算低秩矩阵分量：计算稀疏矩阵分量：其中，j表示矩阵的第j列，||.||₂表示l₂范数。

第二方面，本发明实施例提供了一种地震绕射波分离装置，包括：数据获取模块，用于获取预设地质区域内携带有地下地质信息的地震炮集数据；其中，该地下地质信息包括地质构造信息和地质岩性变化信息；波场反传处理模块，用于对地震炮集数据进行波场反传处理，得到与预设地质区域内的地下成像点一一对应的传播射线的方位角、出射角和振幅信息；矩阵生成模块，用于根据传播射线的方位角、出射角和振幅信息，生成三维角度域成像矩阵；分离模块，用于从三维角度域成像矩阵中分离出低秩矩阵分量，将低秩矩阵分量确定为地震绕射波。

结合第二方面，本发明实施例提供了第二方面的第一种可能的实施方式，其中，上述波场反传处理模块包括：预处理单元，用于对地震炮集数据进行预处理，得到单炮预处理数据；其中，该单炮预处理数据为可用于直接成像的地震炮集数据；预处理包括对地震炮集数据进行噪声去除处理，以及对地震炮集数据与预先保存的历史地震数据进行一一对应；波场反传处理单元，用于将单炮预处理数据和预设的偏移速度模型输入至三维单炮角度域成像公式中，对地震炮集数据进行波场反传处理，得到与预设地质区域内的地下成像点一一对应的传播射线的方位角、出射角和振幅信息；其中，该三维单炮角度域成像公式中包含三维振幅补偿因子。

本发明实施例带来了以下有益效果：

本发明实施例提供的一种地震绕射波分离方法和装置，通过对预设地质区域内携带有地下地质信息的地震炮集数据进行波场反传处理，可以得到与预设地质区域内的地下成像点一一对应的传播射线的方位角、出射角和振幅信息；根据该传播射线的方位角、出射角和振幅信息，可以生成三维角度域成像矩阵，并从该三维角度域成像矩阵中分离出低秩矩阵分量，进而将低秩矩阵分量确定为地震绕射波；上述方式通过构建三维角度域成像矩阵，并利用分离低秩矩阵分量获得地震绕射波，可以提高分离出的绕射波振幅完整性和波形一致性，进而提高地质结构成像的分辨率。

本发明的其他特征和优点将在随后的说明书中阐述，并且，部分地从说明书中变得显而易见，或者通过实施本发明而了解。本发明的目的和其他优点在说明书、权利要求书以及附图中所特别指出的结构来实现和获得。

为使本发明的上述目的、特征和优点能更明显易懂，下文特举较佳实施例，并配合所附附图，作详细说明如下。

附图说明

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

图1为本发明实施例提供的一种地震绕射波分离方法的流程图；

图2为本发明实施例提供的一种地震绕射波分离方法中，从三维角度域成像矩阵中分离出低秩矩阵分量，将低秩矩阵分量确定为地震绕射波的具体流程图；

图3为本发明实施例提供的一种地震绕射波分离装置的结构示意图；

图4为本发明实施例提供的一种地震绕射波分离装置中，分离模块的具体结构示意图。

具体实施方式

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

考虑到现有的地震绕射波分离方式分离出的绕射波振幅完整性和波形一致性较差的问题，本发明实施例提供了一种地震绕射波分离方法和装置，该技术可以应用于根据绕射波的特性对复杂地质结构与岩性的分析中，还可以应用于根据绕射波进行油气藏勘探中；该技术可以采用相关的软件和硬件实现，下面通过实施例进行描述。

实施例一：

参见图1所示的一种地震绕射波分离方法的流程图，该方法包括如下步骤：

步骤S102，获取预设地质区域内携带有地下地质信息的地震炮集数据；其中，该地下地质信息包括地质构造信息和地质岩性变化信息；具体地，该地下地质信息可以为岩层结构、断层、岩溶洞穴和岩性突变点等信息；

步骤S104，对上述地震炮集数据进行波场反传处理，得到与预设地质区域内的地下成像点一一对应的传播射线的方位角、出射角和振幅信息；

步骤S106，根据传播射线的方位角、出射角和振幅信息，生成三维角度域成像矩阵；具体地，该三维角度域成像矩阵与上述方位角和出射角相关，该三维角度域成像矩阵可以用于绕射波分离；

步骤S108，从三维角度域成像矩阵中分离出低秩矩阵分量，将该低秩矩阵分量确定为地震绕射波。在实际实现时，三维角度域成像矩阵中可以分离出低秩矩阵分量和稀疏矩阵分量；其中，该稀疏矩阵分量可以确定为地震反射波。

本发明实施例提供的一种地震绕射波分离方法，通过对预设地质区域内携带有地下地质信息的地震炮集数据进行波场反传处理，可以得到与预设地质区域内的地下成像点一一对应的传播射线的方位角、出射角和振幅信息；根据该传播射线的方位角、出射角和振幅信息，可以生成三维角度域成像矩阵，并从该三维角度域成像矩阵中分离出低秩矩阵分量，进而将低秩矩阵分量确定为地震绕射波；上述方式通过构建三维角度域成像矩阵，并利用分离低秩矩阵分量获得地震绕射波，可以提高分离出的绕射波振幅完整性和波形一致性，进而提高地质结构成像的分辨率。

考虑到获取到的地震炮集数据可处理性较差，上述对地震炮集数据进行波场反传处理，得到与预设地质区域内的地下成像点一一对应的传播射线的方位角、出射角和振幅信息，包括如下步骤：

(1)对地震炮集数据进行预处理，得到单炮预处理数据；其中，该单炮预处理数据为可用于直接成像的地震炮集数据；上述预处理包括对地震炮集数据进行噪声去除处理，以及对地震炮集数据与预先保存的历史地震数据进行一一对应；进一步地，上述预处理还可以包括观测系统加载；

(2)将上述单炮预处理数据和预设的偏移速度模型输入至三维单炮角度域成像公式中，对地震炮集数据进行波场反传处理，得到与预设地质区域内的地下成像点一一对应的传播射线的方位角、出射角和振幅信息；其中，三维单炮角度域成像公式中包含三维振幅补偿因子。

具体地，上述步骤(2)也可以通过下述方式完成：根据上述单炮预处理数据和输入的偏移速度模型，由三维单炮角度域成像公式完成该单炮预处理数据波场反传，得到对应于地下任意成像点位置的传播射线的方位角、出射角和振幅信息。

上述方法通过将获取到的地震炮集数据进行预处理，提高了数据后续的可处理性。

进一步地，上述三维单炮角度域成像公式可以表示为：

其中，δ为脉冲函数，代表三维角度域成像矩阵；震源s激发的射线经地下空间任意成像点x到达检波点位置r；矢量k_S为震源到成像点的射线参数，矢量k_r为检波点到成像点的射线参数；参数θ为出射角；参数为方位角；矢量k为假定反射界面的法向量；k通过下述方式计算获得：θ_s和分别为k_s的出射角和方位角；θ_r和分别为k_r的出射角和方位角；θ_m,分别为假定反射界面的出射角和方位角；n_x为三维坐标系x方向法向量，且n_x＝(1,0,0)；u(s,r,t)为地震数据，t为地震数据的记录时间；t₀为射线走时；W_3D(s,x,r)为三维振幅补偿因子。

上述三维振幅补偿因子W_3D(s,x,r)可以具体表示为：

其中，v_s为震源位置速度，α_s为震源位置射线入射角，α_r为检波点位置射线出射角，N ₁,N ₂分别为第一支射线和第二支射线的走时对震源位置和检波点位置的混合导数，T表示矩阵转置操作，走时根据三维波前重建法射线追踪计算，该计算考虑多值走时；第一支射线为震源到成像点的射线；第二支射线为检波点到成像点的射线；Σ和Γ是与地震观测方式有关的矩阵，在共炮点观测情况下Σ＝0,Γ＝I；其中I为单位矩阵；i表示复数的虚数单位，κ₁,κ₂分别为第一支射线和第二支射线的焦散点数目，κ₁,κ₂根据三维射线追踪动力学方程计算获得。

为了准确高效地分离地震绕射波，上述从三维角度域成像矩阵中分离出低秩矩阵分量，将低秩矩阵分量确定为地震绕射波，可以通过下述方式实现：通过预设的三维绕射波分离模型，从三维角度域成像矩阵中分离出低秩矩阵分量，将该低秩矩阵分量确定为地震绕射波；其中，该预设的三维绕射波分离模型包括：

其中，为第i个成像点x_i位置处的三维角度域成像矩阵；为三维角度域成像矩阵分解后的低秩矩阵分量；为三维角度域成像矩阵分解后的稀疏矩阵分量；x_i为第i个成像点；参数θ为出射角；参数为方位角。

上述三维角度域成像矩阵中每个元素的数值为振幅值，包含了反射波与绕射波信息；根据斯奈尔定理，反射波的入射线和出射线位于同一平面内，且出射角等于入射角；反射波分布在三维角度域成像矩阵特定的方位角和出射角位置，具有稀疏性；根据惠更斯原理，绕射波以球面波形式传播，因此分布在角度域成像矩阵的各个方位角和出射角位置，具有低秩性；其中，上述稀疏性是指矩阵中的大部分元素为零，可用于表征三维角度域成像矩阵中的反射波特征；上述低秩性是指矩阵中的元素呈重复或近似均匀分布，可用于表征三维角度域成像矩阵中的绕射波特征。

优选地，根据增广拉格朗日函数定义，上述预设的三维绕射波分离模型还可以表示为：

其中，J(L,S,Y,β)为目标函数，Y为拉格朗日乘子矩阵，T表示矩阵转置操作，λ为正则化参数，β为保真惩罚因子，||·||_*为核范数，即矩阵奇异值之和，||·||₁为l₁范数，即矩阵中各元素绝对值之和，||·||_F为Frobenius范数，Frobenius范数为矩阵全部元素的平方和的平方根；L为三维角度域成像矩阵分解后的低秩矩阵分量；S为三维角度域成像矩阵分解后的稀疏矩阵分量；R为三维角度域成像矩阵。

参见图2所示的一种地震绕射波分离方法中，从三维角度域成像矩阵中分离出低秩矩阵分量，将低秩矩阵分量确定为地震绕射波的具体流程图；该方法可以通过上述三维绕射波分离模型执行；该方法包括如下步骤：

步骤S202，设置正则化参数λ和预设的最大迭代次数N；其中，λ>0；在实际实现时，λ的数值可以根据经验取值，通过该λ的取值范围为0<λ<1；

步骤S204，设置迭代次数初始值k＝1、低秩矩阵分量的初始值L⁰、稀疏矩阵分量的初始值S⁰、拉格朗日乘子初始值Y⁰、以及保真惩罚因子初始值β⁰；

步骤S206，通过进行奇异值分解计算，得到奇异值对角阵；其中，R为三维角度域成像矩阵；U和V的列为基向量；∑为对角阵；奇异值对角阵的对角上的元素为奇异值；

步骤S208，通过对奇异值对角阵中的奇异值a_i进行软门限阈值操作，得到新对角矩阵其中，x为预设的固定值；

步骤S210，根据新对角矩阵计算低秩矩阵分量L^k和稀疏矩阵分量S^k；

步骤S212，判断L^k和S^k是否满足关系式且k≤N；如果是，执行步骤S214；如果否，执行步骤S216；

步骤S214，更新k＝k+1、拉格朗日乘子Y^k＝Y^k-1+β^k-1(R-L^k-S^k)、以及保真惩罚因子β^k＝ωβ^k-1(ω>0)，其中ω为比例因子；执行步骤S208；

步骤S216，确定L^k为分离的地震绕射波。

上述方式通过迭代的方法，可以高效地获得振幅完整性和波形一致性均很好的地震绕射波。

进一步地，上述根据新对角矩阵计算低秩矩阵分量L^k和稀疏矩阵分量S^k包括：根据新对角矩阵计算低秩矩阵分量：计算稀疏矩阵分量：其中，j表示矩阵的第j列，||.||₂表示l₂范数。

实施例二：

对应于上述方法实施例，参见图3所示的一种地震绕射波分离装置的结构示意图，该装置包括如下部分：

数据获取模块302，用于获取预设地质区域内携带有地下地质信息的地震炮集数据；其中，该地下地质信息包括地质构造信息和地质岩性变化信息；

波场反传处理模块304，用于对地震炮集数据进行波场反传处理，得到与预设地质区域内的地下成像点一一对应的传播射线的方位角、出射角和振幅信息；

矩阵生成模块306，用于根据传播射线的方位角、出射角和振幅信息，生成三维角度域成像矩阵；

分离模块308，用于从三维角度域成像矩阵中分离出低秩矩阵分量，将低秩矩阵分量确定为地震绕射波。

本发明实施例提供的一种地震绕射波分离装置，通过对预设地质区域内携带有地下地质信息的地震炮集数据进行波场反传处理，可以得到与预设地质区域内的地下成像点一一对应的传播射线的方位角、出射角和振幅信息；根据该传播射线的方位角、出射角和振幅信息，可以生成三维角度域成像矩阵，并从该三维角度域成像矩阵中分离出低秩矩阵分量，进而将低秩矩阵分量确定为地震绕射波；上述方式通过构建三维角度域成像矩阵，并利用分离低秩矩阵分量获得地震绕射波，可以提高分离出的绕射波振幅完整性和波形一致性，进而提高地质结构成像的分辨率。

考虑到获取到的地震炮集数据可处理性较差，上述波场反传处理模块包括：(1)预处理单元，用于对地震炮集数据进行预处理，得到单炮预处理数据；其中，该单炮预处理数据为可用于直接成像的地震炮集数据；上述预处理包括对地震炮集数据进行噪声去除处理，以及对地震炮集数据与预先保存的历史地震数据进行一一对应；(2)波场反传处理单元，用于将单炮预处理数据和预设的偏移速度模型输入至三维单炮角度域成像公式中，对地震炮集数据进行波场反传处理，得到与预设地质区域内的地下成像点一一对应的传播射线的方位角、出射角和振幅信息；其中，上述三维单炮角度域成像公式中包含三维振幅补偿因子。上述方法通过将获取到的地震炮集数据进行预处理，提高了数据后续的可处理性。

为了准确高效地分离地震绕射波，上述分离模块还用于：通过预设的三维绕射波分离模型，从三维角度域成像矩阵中分离出低秩矩阵分量，将该低秩矩阵分量确定为地震绕射波；其中，该预设的三维绕射波分离模型包括：

其中，为第i个成像点x_i位置处的三维角度域成像矩阵；为三维角度域成像矩阵分解后的低秩矩阵分量；为三维角度域成像矩阵分解后的稀疏矩阵分量；x_i为第i个成像点；参数θ为出射角；参数为方位角。

参见图4所示的一种地震绕射波分离装置中，分离模块的具体结构示意图，该装置包括如下部分：

第一设置模块402，用于设置正则化参数λ和预设的最大迭代次数N；其中，λ>0；

第二设置模块404，用于设置迭代次数初始值k＝1、低秩矩阵分量的初始值L⁰、稀疏矩阵分量的初始值S⁰、拉格朗日乘子初始值Y⁰、以及保真惩罚因子初始值β⁰；

分解计算模块406，用于通过进行奇异值分解计算，得到奇异值对角阵；其中，R为三维角度域成像矩阵；U和V的列为基向量；∑为对角阵；上述奇异值对角阵的对角上的元素为奇异值；

阈值操作模块408，用于通过对奇异值对角阵中的奇异值a_i进行软门限阈值操作，得到新对角矩阵其中，x为预设的固定值；

计算模块410，用于根据新对角矩阵计算低秩矩阵分量L^k和稀疏矩阵分量S^k；

判断模块412，用于判断L^k和S^k是否满足关系式且k≤N；

更新模块414，用于如果L^k和S^k满足关系式且k≤N，更新k＝k+1、拉格朗日乘子Y^k＝Y^k-1+β^k-1(R-L^k-S^k)、以及保真惩罚因子β^k＝ωβ^k-1(ω>0)，其中ω为比例因子；并继续进行迭代处理；

确定模块416，用于如果L^k和S^k不满足关系式且k≤N，确定L^k为分离的地震绕射波。

上述方式通过迭代的方法，可以高效地获得振幅完整性和波形一致性均很好的地震绕射波。

进一步地，上述计算模块410包括：

第一计算单元，用于根据新对角矩阵计算低秩矩阵分量：

第二计算单元，用于计算稀疏矩阵分量：

其中，j表示矩阵的第j列，||.||₂表示l₂范数。

对比而言，现有技术中关于绕射波分离的研究中，Harlan等(1984)利用Radon变换去除反射波，并借鉴统计学原理分离绕射波。Bansal和Imhof(2005)通过信号处理方法，研究地震处理流程中的标准模块，分析了不同反射波去除方法。Taner等(2006)利用平面波分解方法压制反射波，分离绕射波。通过研究倾角域绕射波与反射波几何形态差异，Landa和Fomel(2008)提出一种基于平面波滤波的倾角域绕射波分离方法。Khaidukov等(2004)提出聚焦-切除-反聚焦方法，实现叠前时间域绕射波成像，该方法对速度模型依懒性强，反射波去除难度大，在实际应用中具有局限性。Figueiredo等(2013)利用模式识别技术研究绕射波自动成像方法。

上述常规绕射波分离方法，大多利用反射波与绕射波运动学特征，通过信号处理方法分离绕射波，且没有针对三维炮集数据开展研究。三维炮集数据中绕射波与反射波运动学特性相似，仅依靠常规运动学波场分离方法难以有效处理，但炮集数据中绕射波由同一震源激发波形一致性强，有利于高分辨率绕射波成像。因此本发明通过研究斯奈尔定理和惠更斯原理，构建了三维角度域成像矩阵，用于三维炮集绕射波分离，该技术利用低秩与稀疏优化分解方法，通过捕捉地震数据动力学特征，分离绕射波，能够保证绕射波分离完整性和波形特征的一致性，有利于高分辨率成像，在碳酸盐岩缝洞型油气藏勘探开发中具有重要应用价值。

本发明实施例所提供的一种地震绕射波分离方法和装置的计算机程序产品，包括存储了程序代码的计算机可读存储介质，所述程序代码包括的指令可用于执行前面方法实施例中所述的方法，具体实现可参见方法实施例，在此不再赘述。

所属领域的技术人员可以清楚地了解到，为描述的方便和简洁，上述描述的装置的具体工作过程，可以参考前述方法实施例中的对应过程，在此不再赘述。

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

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

Claims (10)

1.一种地震绕射波分离方法，其特征在于，包括：

获取预设地质区域内携带有地下地质信息的地震炮集数据；其中，所述地下地质信息包括地质构造信息和地质岩性变化信息；

对所述地震炮集数据进行波场反传处理，得到与预设地质区域内的地下成像点一一对应的传播射线的方位角、出射角和振幅信息；

根据所述传播射线的方位角、出射角和振幅信息，生成三维角度域成像矩阵；

从所述三维角度域成像矩阵中分离出低秩矩阵分量，将所述低秩矩阵分量确定为地震绕射波。

2.根据权利要求1所述的方法，其特征在于，所述对所述地震炮集数据进行波场反传处理，得到与预设地质区域内的地下成像点一一对应的传播射线的方位角、出射角和振幅信息包括：

对所述地震炮集数据进行预处理，得到单炮预处理数据；其中，所述单炮预处理数据为可用于直接成像的地震炮集数据；所述预处理包括对所述地震炮集数据进行噪声去除处理，以及对所述地震炮集数据与预先保存的历史地震数据进行一一对应；

将所述单炮预处理数据和预设的偏移速度模型输入至三维单炮角度域成像公式中，对所述地震炮集数据进行波场反传处理，得到与预设地质区域内的地下成像点一一对应的传播射线的方位角、出射角和振幅信息；其中，所述三维单炮角度域成像公式中包含三维振幅补偿因子。

3.根据权利要求2所述的方法，其特征在于，所述三维单炮角度域成像公式包括：

其中，δ为脉冲函数，代表三维角度域成像矩阵；震源s激发的射线经地下空间任意成像点x到达检波点位置r；矢量k_S为震源到成像点的射线参数，矢量k_r为检波点到成像点的射线参数；参数θ为出射角；参数为方位角；矢量k为假定反射界面的法向量；所述k通过下述方式计算获得：θ_s和分别为所述k_s的出射角和方位角；θ_r和分别为所述k_r的出射角和方位角；分别为假定反射界面的出射角和方位角；n_x为三维坐标系x方向法向量，且n_x＝(1,0,0)；u(s,r,t)为地震数据，t为地震数据的记录时间；t₀为射线走时；W_3D(s,x,r)为三维振幅补偿因子。

4.根据权利要求3所述的方法，其特征在于，所述三维振幅补偿因子W_3D(s,x,r)包括：

$W_{3D}(s,x,r)=\frac{1}{v_s}\sqrt{{\mathrm{cos\&alpha;}}_s{\mathrm{cos\&alpha;}}_r}\frac{\left|\det ({\underset{\&OverBar;}{N}}_1^T\underset{\&OverBar;}{\&Sigma;}+{\underset{\&OverBar;}{N}}_2^T\mathrm{\&Gamma;})\right|}{\sqrt{\left|\det {\underset{\&OverBar;}{N}}_1\right|\left|\det {\underset{\&OverBar;}{N}}_2\right|}}{e}^{-i\frac{\mathrm{\&pi;}}{2}({\mathrm{\&kappa;}}_1+{\mathrm{\&kappa;}}_2)}$

其中，v_s为震源位置速度，α_s为震源位置射线入射角，α_r为检波点位置射线出射角，N ₁,N ₂分别为第一支射线和第二支射线的走时对震源位置和检波点位置的混合导数，T表示矩阵转置操作，所述走时根据三维波前重建法射线追踪计算，该计算考虑多值走时；所述第一支射线为震源到成像点的射线；所述第二支射线为检波点到成像点的射线；Σ和Γ是与地震观测方式有关的矩阵，在共炮点观测情况下Σ＝0,Γ＝I；其中I为单位矩阵；i表示复数的虚数单位，κ₁,κ₂分别为第一支射线和第二支射线的焦散点数目，κ₁,κ₂根据三维射线追踪动力学方程计算获得。

5.根据权利要求1所述的方法，其特征在于，所述从所述三维角度域成像矩阵中分离出低秩矩阵分量，将所述低秩矩阵分量确定为地震绕射波包括：

通过预设的三维绕射波分离模型，从所述三维角度域成像矩阵中分离出低秩矩阵分量，将所述低秩矩阵分量确定为地震绕射波；其中，所述预设的三维绕射波分离模型包括：

其中，为第i个成像点x_i位置处的三维角度域成像矩阵；为所述三维角度域成像矩阵分解后的低秩矩阵分量；为所述三维角度域成像矩阵分解后的稀疏矩阵分量；x_i为第i个成像点；参数θ为出射角；参数为方位角。

6.根据权利要求5所述的方法，其特征在于，所述预设的三维绕射波分离模型还包括：

$J(L,S,Y,\mathrm{\&beta;})=\left|\right|L|{|}_{*}+\mathrm{\&lambda;}|\left|S\right|{|}_1+Y^T(R-L-S)+\frac{\mathrm{\&beta;}}{2}\left|\right|R-L-S|{|}_F$

其中，J(L,S,Y,β)为目标函数，Y为拉格朗日乘子矩阵，T表示矩阵转置操作，λ为正则化参数，β为保真惩罚因子，||·||_*为核范数，即矩阵奇异值之和，||·||₁为l₁范数，即矩阵中各元素绝对值之和，||·||_F为Frobenius范数，所述Frobenius范数为矩阵全部元素的平方和的平方根；L为所述三维角度域成像矩阵分解后的低秩矩阵分量；S为所述三维角度域成像矩阵分解后的稀疏矩阵分量；R为三维角度域成像矩阵。

7.根据权利要求6所述的方法，其特征在于，所述从所述三维角度域成像矩阵中分离出低秩矩阵分量，将所述低秩矩阵分量确定为地震绕射波包括：

设置正则化参数λ和预设的最大迭代次数N；其中，所述λ>0；

设置迭代次数初始值k＝1、所述低秩矩阵分量的初始值L⁰、所述稀疏矩阵分量的初始值S⁰、拉格朗日乘子初始值Y⁰、以及保真惩罚因子初始值β⁰；

将所述k＝1、所述L⁰、所述S⁰、所述Y⁰和所述β⁰作为初始值，对所述三维角度域成像矩阵进行迭代处理，所述迭代处理包括如下步骤：

通过进行奇异值分解计算，得到奇异值对角阵；其中，R为三维角度域成像矩阵；U和V的列为基向量；Σ为对角阵；所述奇异值对角阵的对角上的元素为奇异值；

通过对所述奇异值对角阵中的奇异值a_i进行软门限阈值操作，得到新对角矩阵其中，x为预设的固定值；

根据所述新对角矩阵计算低秩矩阵分量L^k和稀疏矩阵分量S^k；

判断所述L^k和S^k是否满足关系式且k≤N；

如果是，更新k＝k+1、拉格朗日乘子Y^k＝Y^k-1+β^k-1(R-L^k-S^k)、以及保真惩罚因子β^k＝ωβ^k-1(ω>0)，其中ω为比例因子；并继续进行所述迭代处理；

如果否，确定所述L^k为分离的地震绕射波。

8.根据权利要求7所述的方法，其特征在于，所述根据所述新对角矩阵计算低秩矩阵分量L^k和稀疏矩阵分量S^k包括：

根据所述新对角矩阵计算低秩矩阵分量：

计算稀疏矩阵分量：其中，j表示矩阵的第j列，||.||₂表示l₂范数。

9.一种地震绕射波分离装置，其特征在于，包括：

数据获取模块，用于获取预设地质区域内携带有地下地质信息的地震炮集数据；其中，所述地下地质信息包括地质构造信息和地质岩性变化信息；

波场反传处理模块，用于对所述地震炮集数据进行波场反传处理，得到与预设地质区域内的地下成像点一一对应的传播射线的方位角、出射角和振幅信息；

矩阵生成模块，用于根据所述传播射线的方位角、出射角和振幅信息，生成三维角度域成像矩阵；

分离模块，用于从所述三维角度域成像矩阵中分离出低秩矩阵分量，将所述低秩矩阵分量确定为地震绕射波。

10.根据权利要求9所述的装置，其特征在于，所述波场反传处理模块包括：

预处理单元，用于对所述地震炮集数据进行预处理，得到单炮预处理数据；其中，所述单炮预处理数据为可用于直接成像的地震炮集数据；所述预处理包括对所述地震炮集数据进行噪声去除处理，以及对所述地震炮集数据与预先保存的历史地震数据进行一一对应；

波场反传处理单元，用于将所述单炮预处理数据和预设的偏移速度模型输入至三维单炮角度域成像公式中，对所述地震炮集数据进行波场反传处理，得到与预设地质区域内的地下成像点一一对应的传播射线的方位角、出射角和振幅信息；其中，所述三维单炮角度域成像公式中包含三维振幅补偿因子。