基于交替分裂Bregman迭代算法的一次波和多次波分离方法

摘要

本发明属于地震勘探技术中地震信号处理领域，具体公开了一种基于交替分裂Bregman迭代算法的一次波和多次波分离方法。对于基于3D匹配滤波器的多次波自适应相减方法，本发明利用交替分裂Bregman迭代算法求解对一次波施加稀疏约束的优化问题，实现3D匹配滤波器的估计，并利用估计的3D匹配滤波器自适应分离3D数据窗口中的一次波和多次波。相比于传统的迭代重加权最小二乘算法，本发明中的交替分裂Bregman迭代算法在每一个3D数据窗口估计3D匹配滤波器时，只需计算一次矩阵-矩阵相乘和矩阵求逆，能够有效地降低优化问题求解的计算复杂度，提高一次波和多次波自适应分离的计算效率。

说明书

技术领域

本发明属于地震勘探技术中的地震信号处理领域，具体涉及一种基于交替分裂Bregman 迭代算法的一次波和多次波分离方法。

背景技术

SRME(SurfaceRelatedMultipleElimination)是海洋地震勘探中广泛采用的多次波压制 方法。SRME预测的多次波与真实多次波存在时间和空间差异，可以利用匹配滤波器将预测 多次波从原始数据中自适应减去。通常，随着匹配滤波器维数的增加(从1D、2D到3D匹 配滤波器)，多次波能更好地从原始数据中自适应减去，但需要更多的计算时间。基于3D 匹配滤波器的多次波自适应相减方法的数学模型为(Li,Z.,andW.Lu,2013,Adaptivemultiple subtractionbasedon3Dblindseparationofconvolvedmixtures:Geophysics,78,V251-V266)：

v＝d-Hx(1)

其中，v表示估计一次波，d表示原始数据，x表示3D匹配滤波器，H表示预测多次波 的褶积矩阵。

基于3D匹配滤波器的多次波自适应相减方法在相互重叠的3D数据窗口内，通过估计 3D匹配滤波器来自适应分离一次波和多次波。为估计3D匹配滤波器，传统的多次波自适应 相减方法对估计一次波施加能量最小化约束。另外，为避免滤波器估计的不稳定性，滤波器 系数也假设满足能量最小化约束。相应的优化问题为：

$\arg \underset{x}{min}\left[\right||d-Hx|{|}_2^2+\mu \left|\right|x\left|{|}_2^2\right]---\left(2\right)$

其中，μ为正则化参数。方程(2)中的3D匹配滤波器可以采用最小二乘算法进行求解：

x＝(H^TH+μI)^-1H^Td(3)

其中，I为单位矩阵。

最小二乘算法需要一次波和多次波的正交性假设。当一次波和多次波相互重叠或有强一 次波同相轴存在时，最小二乘算法会产生残余多次波或造成一次波的损伤。为克服正交性假 设的缺点，对一次波施加稀疏约束已经引入到多次波自适应相减方法中。另外，假设3D匹 配滤波器系数满足能量最小化约束来确保3D匹配滤波器估计的稳定性，相应的优化问题为：

$\arg \underset{x}{min}\left[\right||d-Hx|{|}_1+\lambda \left|\right|x\left|{|}_2^2\right],---\left(4\right)$

其中，λ为正则化参数。可以采用迭代重加权最小二乘算法(Guitton,A.,andD.J.Verschuur, 2004,AdaptivesubtractionofmultiplesusingtheL₁-norm:GeophysicalProspecting,52,27-38)来估计式(4) 中的3D匹配滤波器。然而，迭代重加权最小二乘算法在每一步迭代均需计算一次矩阵-矩阵 相乘和矩阵求逆，计算复杂度较高。

发明内容

针对现有技术中存在的上述技术问题，本发明提出了一种基于交替分裂Bregman迭代的 一次波和多次波分离方法，该方法能够有效降低对一次波施加稀疏约束的优化问题求解的计 算复杂度，提高一次波和多次波分离的计算效率。

为了实现上述目的，本发明采用如下技术方案：

基于交替分裂Bregman迭代算法的一次波和多次波分离方法，包括如下步骤：

a设置变量初始值，需要设置初始值的变量包括3D数据窗口的时间长度T₀，3D数据窗 口的空间长度X₀，3D数据窗口的道集个数Y₀，3D匹配滤波器的时间长度K，3D匹配滤波 器的空间长度R，3D匹配滤波器的道集个数G，一次波阈值s_β，阻尼因子β和迭代次数

b输入原始数据道集中一个待处理的3D数据窗口中的数据d，然后利用预测多次波数据 以及3D数据窗口的参数T₀、X₀和Y₀、3D匹配滤波器的参数K、R和G构造褶积矩阵H，并 采用克莱斯基分解计算逆矩阵

c基于交替分裂Bregman迭代算法，利用步骤b得到的逆矩阵对地震道集的所有3D 数据窗口的数据d逐一进行处理；

d判断原始数据道集中所有3D数据窗口内的数据d是否全部处理完毕；如果否，返回步 骤c；如果全部处理完毕，则首先采用3D汉宁窗将每一个3D数据窗口中估计的一次波进行 加权，并融合为一个数据体然后采用同样的方式将3D汉宁窗融合为一个数据体最 终的一次波估计结果表示为：其中，/表示逐个元素的相除操作。

优选地，所述步骤c具体包括：

c1设置迭代数m＝1，利用步骤b得到的逆矩阵求取一次波的初始估计值v⁽⁰⁾： $v^{\left(0\right)}=d-H\left(\bar{H}\right(H^{T}d\left)\right),$并令b₁⁽⁰⁾＝0，u₁⁽⁰⁾＝0；

c2计算向量y^(m)＝v^(m-1)+b₁^(m-1)；

c3对y^(m)利用如下的距离算子prox计算向量u₁^(m)：

${u_1}^{\left(m\right)}={\mathrm{prox}}_s^1\left(y^{\left(m\right)}\right)=\left\{{\left(\right|y_{i,j,n}^{\left(m\right)}|-s_{\beta }{C}^{\left(m\right)})}_{+}\mathrm{sgn}\left({y_{i,j,n}}^{\left(m\right)}\right)\right\},$

其中，i＝1,2,…,T₀,j＝1,2,…,X₀,n＝1,2,…,Y₀，s＝s_βC^(m)，0＜s_β＜1；y^(m)＝{y_i,j,n^(m)}， y_i,j,n^(m)表示向量y^(m)中下标为(i,j,n)的元素，C^(m)＝max(|y_i,j,n^(m)|)，${\left(b\right)}_{+}=\left(\begin{array}{cc}b,& b\ge 0\\ 0,& b<0\end{array}\right),$$\mathrm{sgn}\left(b\right)=\left(\begin{array}{cc}1,& b>0\\ 0,& b=0\\ -1,& b<0\end{array}\right);$

c4计算向量b₁^(m)＝b₁^(m-1)-[u₁^(m)-v^(m-1)]；

c5计算3D匹配滤波器$X^{\left(m\right)}=\bar{H}\left(H^T\right(d+{b_1}^{\left(m\right)}-{u_1}^{\left(m\right)}\left)\right);$

c6估计一次波v^(m)＝d-Hx^(m)；

c7令m＝m+1，如果返回到步骤c2；如果输出当前3D数据窗口估计 的一次波结果。

本发明具有如下优点：

对于基于3D匹配滤波器的多次波自适应相减方法，本发明利用交替分裂Bregman迭代 算法求解对一次波施加稀疏约束的优化问题，实现3D匹配滤波器的估计，并利用估计的3D 匹配滤波器自适应分离3D数据窗口中的一次波和多次波。相比于传统的迭代重加权最小二 乘算法，本发明中的交替分裂Bregman迭代算法在每一个3D数据窗口估计3D匹配滤波器时， 只需计算一次矩阵-矩阵相乘和矩阵求逆，能有效地降低优化问题求解的计算复杂度，提高一 次波和多次波自适应分离的计算效率。

附图说明

图1为本发明中基于交替分裂Bregman迭代算法的一次波和多次波分离方法的流程图；

图2a为原始数据的共偏移距道集图；

图2b为预测多次波的共偏移距道集图；

图3a为基于最小二乘算法的一次波估计结果的共偏移距道集图；

图3b为基于迭代重加权最小二乘算法的一次波估计结果的共偏移距道集图；

图3c为基于交替分裂Bregman迭代算法的一次波估计结果的共偏移距道集图；

图4a为基于最小二乘算法所去除多次波的共偏移距道集图；

图4b为基于迭代重加权最小二乘算法所去除多次波的共偏移距道集图；

图4c为基于交替分裂Bregman迭代算法所去除多次波的共偏移距道集图；

图5a为原始数据的放大显示结果图(对应于图2a中黑色方框)；

图5b为预测多次波的放大显示结果图；

图5c为基于最小二乘算法所估计一次波的放大显示结果图；

图5d为基于迭代重加权最小二乘算法所估计一次波的放大显示结果图；

图5e为基于交替分裂Bregman迭代算法所估计一次波的放大显示结果图。

具体实施方式

本发明的基本思想是：逐个3D数据窗口地分离一次波和多次波，构建对一次波施加稀 疏约束的优化问题：其中，d为原始数据，x表示3D匹配滤波 器，H表示预测多次波的褶积矩阵，λ为正则化因子。求解上式中的优化问题来估计3D匹配 滤波器，并采用估计的3D匹配滤波器来对3D数据窗口中一次波和多次波进行分离，最后将 所有3D数据窗口内的一次波估计结果进行合并，得到最终的一次波估计结果。

下面结合附图以及具体实施方式对本发明作进一步详细说明：

如图1所示，基于交替分裂Bregman迭代算法的一次波和多次波分离方法，包括步骤：

a设置变量初始值，需要设置初始值的变量包括3D数据窗口的时间长度T₀，3D数据窗 口的空间长度X₀，3D数据窗口的道集个数Y₀，3D匹配滤波器的时间长度K，3D匹配滤波 器的空间长度R，3D匹配滤波器的道集个数G，一次波阈值s_β，阻尼因子β和迭代次数

b输入原始数据道集中一个待处理的3D数据窗口中的数据d，然后利用预测多次波数据 以及3D数据窗口的参数T₀、X₀和Y₀、3D匹配滤波器的参数K、R和G构造褶积矩阵H，矩 阵H的行数为T₀X₀Y₀，列数为KRG，并采用克莱斯基分解计算逆矩阵

c利用步骤b得到的逆矩阵H对地震道集的3D数据窗口的数据d进行处理；

本发明采用交替分裂Bregman迭代算法求解3D匹配滤波器，实现一次波和多次波的自 适应分离，即在每一步迭代，交替分裂Bregman迭代算法交替求解对一次波施加稀疏约束的 优化问题和对3D匹配滤波器施加能量最小化约束的优化问题；具体地讲，在每一步迭代， 交替分裂Bregman迭代算法采用距离算子求解对一次波施加稀疏约束的优化问题，采用最小 二乘算法求解对3D匹配滤波器施加能量最小化约束的优化问题；同时假设正则化参数间的 倍数关系对迭代步骤进行简化；最后利用估计的3D匹配滤波器自适应分离一次波和多次波。

其具体处理过程为：

c1设置迭代数m＝1，利用步骤b得到的逆矩阵求取一次波的初始估计值v⁽⁰⁾： $v^{\left(0\right)}=d-H\left(\bar{H}\right(H^{T}d\left)\right),$并令b₁⁽⁰⁾＝0，u₁⁽⁰⁾＝0；

c2计算向量y^(m)＝v^(m-1)+b₁^(m-1)；

c3对y^(m)利用如下的距离算子prox计算向量u₁^(m)：

${u_1}^{\left(m\right)}={\mathrm{prox}}_s^1\left(y^{\left(m\right)}\right)=\left\{{\left(\right|{y_{i,j,n}}^{\left(m\right)}|-s_{\beta }{C}^{\left(m\right)})}_{+}\mathrm{sgn}{{(}_{i,j,n}}^{\left(m\right)}\right)\},$

其中，i＝1,2,…,T₀,j＝1,2,…,X₀,n＝1,2,…,Y₀，s＝s_βC^(m)，0＜s_β＜1；y^(m)＝{y_i,j,n^(m)}， y_i,j,n^(m)表示向量y^(m)中下标为(i,j,n)的元素，C^(m)＝max(|y_i,j,n^(m)|)，${\left(b\right)}_{+}=\left(\begin{array}{cc}b,& b\ge 0\\ 0,& b<0\end{array}\right),$$\mathrm{sgn}\left(b\right)=\left(\begin{array}{cc}1,& b>0\\ 0,& b=0\\ -1,& b<0\end{array}\right);$

c4计算向量b₁^(m)＝b₁^(m-1)-[u₁^(m)-v^(m-1)]；

c5计算3D匹配滤波器$X^{\left(m\right)}=\bar{H}\left(H^T\right(d+{b_1}^{\left(m\right)}-{u_1}^{\left(m\right)}\left)\right);$

c6估计一次波v^(m)＝d-Hx^(m)；

c7令m＝m+1，如果返回到步骤c2；如果输出当前数据窗口估计的 一次波结果。

d判断原始数据道集中所有数据窗口内的数据d是否全部处理完毕；如果否，返回步骤c； 如果全部处理完毕，则首先采用3D汉宁窗将每一个3D数据窗口中估计的一次波进行加权， 并融合为一个数据体然后采用同样的方式将3D汉宁窗融合为一个数据体最终的一 次波估计结果表示为：其中，/表示逐个元素的相除操作。

此处，同样的方式是指把相同时间、空间和道集位置处的数据进行相加。

下面利用2D实际数据验证本发明方法的有效性：

图2a为原始数据的共偏移距道集图。图2b为预测多次波的共偏移距道集图。图3a为基 于最小二乘算法的一次波估计结果的共偏移距道集图。图3b为基于迭代重加权最小二乘算法 的一次波估计结果的共偏移距道集图。图3c为基于交替分裂Bregman迭代算法的一次波估计 结果的共偏移距道集图。图4a为基于最小二乘算法所去除多次波的共偏移距道集图。图4b 为基于迭代重加权最小二乘算法所去除多次波的共偏移距道集图。图4c为基于交替分裂 Bregman迭代算法所去除多次波的共偏移距道集图。图3a、图3b和图3c中得到相近的一次 波估计结果，计算时间分别为8.88秒、37.10秒和9.63秒。相比于迭代重加权最小二乘算法， 本发明中的交替分裂Bregman迭代算法在保持一次波估计精度的同时，能有效地提高计算效 率。另外，图5a为原始数据的放大显示结果图，对应于图2a中的黑色方框。图5b为预测多 次波的放大显示结果图。图5c为基于最小二乘算法所估计一次波的放大显示结果图。图5d 为基于迭代重加权最小二乘算法所估计一次波的放大显示结果图。图5e为基于交替分裂 Bregman迭代算法所估计一次波的放大显示结果图。其中，图5a至图5e中的白色箭头表明， 相对于传统的最小二乘算法，本发明中的交替分裂Bregman迭代算法在一次波估计结果中能 更好地去除残余多次波。其中，图2a至2b、图3a至3c、图4a至4c及图5a至5e中横坐标 TraceNumber表示道号，纵坐标Time表示时间，单位为毫秒(ms)。

当然，以上说明仅仅为本发明的较佳实施例，本发明并不限于列举上述实施例，应当说 明的是，任何熟悉本领域的技术人员在本说明书的教导下，所做出的所有等同替代、明显变 形形式，均落在本说明书的实质范围之内，理应受到本发明的保护。