地震数据的处理方法和装置与流程

本发明涉及地震数据处理的技术领域，尤其是涉及一种地震数据的处理方法和装置。

背景技术：

在地震勘探中，往往会由于采集成本、坏道、噪音、地形等诸多因素的影响导致采集的地震数据不满足采样定理，不满足采样定理的地震数据会影响多次波消除、偏移成像等后期处理效果。因此，需要提前对地震数据进行规则化处理。

现有技术中最常用的数据规则化方法是基于各种变换的方法。具体地，假设地震数据在变化域是稀疏的，此时可以利用压缩感知理论，实现地震数据的稀疏分解及重构，例如，傅里叶变换、radon变化、curvelet变化等。

但是，地震数据的稀疏表达问题是一个l0范数最优化问题，直接求解十分困难，这是一个np完全问题。因此，可以通过匹配追踪法进行求解，但这是一个贪婪的算法，每次迭代时选取一个局部最优解逐步逼近原始信号。之后发展了正交匹配追踪法、正则化匹配追踪法等，虽然这类方法执行简单，但是计算速度较慢，且不是全局最优。更加常用的方法是把l0范数最优化问题松弛为l1范数的凸优化问题，这样就可以用常见的连续凸优化算法求解，例如，基追踪算法、内点法、投影梯度算法等，但这类算法的解不是最稀疏的。

技术实现要素：

本发明的目的在于提供一种地震数据的处理方法和装置，以缓解通过传统的正则化处理方法对地震数据进行处理之后处理效果较差，以及处理效率较慢的技术问题。

本发明实施例的一个方面，提供了一种地震数据的处理方法，包括：将目标地震数据表示为高斯束叠加形式的地震数据，其中，所述高斯束叠加形式的地震数据中包括所述目标地震数据的波形信息；将所述高斯束叠加形式的地震数据进行稀疏处理，处理得到波形函数矩阵；基于处理之后得到的所述波形函数矩阵进行地震数据的重构，得到重构之后的地震数据。

进一步地，将所述高斯束叠加形式的地震数据进行稀疏处理，处理得到波形函数矩阵包括：基于波形信息矩阵建立稀疏矩阵优化模型，其中，所述波形信息矩阵中包括所述波形信息；对所述稀疏矩阵优化模型中的波形信息矩阵进行稀疏处理，得到处理之后的所述稀疏矩阵优化模型，以使稀疏处理之后的所述波形信息矩阵中的每个数值均不等于预设数值；对处理之后的所述稀疏矩阵优化模型进行稀疏求解，求解得到所述波形函数矩阵。

进一步地，对所述稀疏矩阵优化模型中的波形信息矩阵进行稀疏处理包括：对所述目标地震数据中的高斯束束中心相邻道上进行相关性度量分析，得到地震子波，其中，所述地震子波中包含所述地震子波所处的时间和倾角方向；根据所述地震子波所处的时间和所述倾角方向确定所述波形信息矩阵中数值等于所述预设数值的向量；在所述波形信息矩阵中删除数值等于所述预设数值的向量，得到处理之后的所述稀疏矩阵优化模型。

进一步地，对处理之后的所述稀疏矩阵优化模型进行稀疏求解，求解得到所述波形函数矩阵包括：反复通过第一公式和第二公式对稀疏矩阵优化模型进行迭代计算，求解得到所述波形函数矩阵，直至满足第三公式时，确定所述第二公式中的输出结果为所述波形函数矩阵，其中，α1，α2和αk为对角线矩阵λ中的数值，gk为迭代变量，p为稀疏度参数，k为迭代指数，为矩阵算子，d为数据矩阵，为稀疏处理之后的所述波形信息矩阵，im表示m维恒等算子，为矩阵的第k行向量，

进一步地，基于处理之后得到的所述波形函数矩阵进行地震数据的重构，得到重构之后的地震数据包括：通过公式d(x,t)＝∑γcγψγ(x,t)进行地震数据的重构，得到重构之后的地震数据，其中，为所述波形函数矩阵，d(x,t)为所述重构之后的地震数据。

进一步地，将目标地震数据表示为高斯束叠加形式的地震数据包括：获取高斯束基函数的半宽w和所述高斯束中任意相邻的两个束中心之间的束中心间距a；根据所述半宽和所述束中心间距计算所述高斯束的束中心位置xm，其中，xm＝am，m＝1,2,…,m，m为束中心的个数；通通过公式确定所述目标地震数据中每道地震数据与所述束中心相对应的基函数，其中，tk为所述高斯束基函数的束中心所处的时间位置，px为水平慢度，v为地震速度常数，ωr为相关低频；通过所述目标地震数据中每道地震数据与所述束中心相对应的基函数，将所述目标地震数据表示为所述高斯束叠加形式的地震数据。

本发明实施例的另一个方面，还提供了一种地震数据的处理装置，包括：第一处理单元，用于将目标地震数据表示为高斯束叠加形式的地震数据，其中，所述高斯束叠加形式的地震数据中包括所述目标地震数据的波形信息；第二处理单元，用于将所述高斯束叠加形式的地震数据进行稀疏处理，处理得到波形函数矩阵；重构单元，用于基于处理之后得到的所述波形函数矩阵进行地震数据的重构，得到重构之后的地震数据。

进一步地，所述第二处理单元包括：建立模块，用于基于波形信息矩阵建立稀疏矩阵优化模型，其中，所述波形信息矩阵中包括所述波形信息；稀疏处理模块，用于对所述稀疏矩阵优化模型中的波形信息矩阵进行稀疏处理，得到处理之后的所述稀疏矩阵优化模型，以使稀疏处理之后的所述波形信息矩阵中的每个数值均不等于预设数值；求解模块，用于对处理之后的所述稀疏矩阵优化模型进行稀疏求解，求解得到所述波形函数矩阵。

进一步地，所述稀疏处理模块包括：分析子模块，用于对所述目标地震数据中的高斯束束中心相邻道上进行相关性度量分析，得到地震子波，其中，所述地震子波中包含所述地震子波所处的时间和倾角方向；确定子模块，用于根据所述地震子波所处的时间和所述倾角方向确定所述波形信息矩阵中数值等于所述预设数值的向量；删除子模块，用于在所述波形信息矩阵中删除数值等于所述预设数值的向量，得到处理之后的所述稀疏矩阵优化模型。

进一步地，所述求解模块包括：计算子模块，用于反复通过第一公式和第二公式对稀疏矩阵优化模型进行迭代计算，求解得到所述波形函数矩阵，直至满足第三公式时，确定所述第二公式中的输出结果为所述波形函数矩阵，其中，α1，α2和αk为对角线矩阵λ中的数值，gk为迭代变量，p为稀疏度参数，k为迭代指数，为矩阵算子，d为数据矩阵，为稀疏处理之后的所述波形信息矩阵，im表示m维恒等算子，为矩阵的第k行向量，

在本发明实施例中，首先将目标地震数据表示为高斯束叠加形式的地震数据，其中，高斯束叠加形式的地震数据包括目标地震数据包括波形信息；接下来，将高斯束叠加形式的地震数据进行稀疏处理，处理得到波形函数矩阵；最后，基于波形函数矩阵进行地震数据的重构，以完成地震数据的正则化处理。通过波形函数矩阵对地震数据进行重构，能够将空间随机采样的地震数据差值到规则网格上，以消除随机噪音和特定倾角的噪音，进而缓解了通过传统的正则化处理方法对地震数据进行处理之后处理效果较差，以及处理效率较慢的技术问题，从而实现了提高了地震数据的正则化处理效果的技术效果。

附图说明

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

图1是根据本发明实施例的一种地震数据的处理方法的流程图；

图2是根据本发明实施例的一种地震数据的处理装置的示意图；

图3是根据本发明实施例的另一种可选地地震数据的处理装置的示意图；

图4是根据本发明实施例的另一种可选地地震数据的处理装置的示意图；

图5是根据本发明实施例的另一种可选地地震数据的处理装置的示意图。

具体实施方式

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

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

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

实施例一

根据本发明实施例，提供了一种地震数据的处理方法的实施例，需要说明的是，在附图的流程图示出的步骤可以在诸如一组计算机可执行指令的计算机系统中执行，并且，虽然在流程图中示出了逻辑顺序，但是在某些情况下，可以以不同于此处的顺序执行所示出或描述的步骤。

图1是根据本发明实施例的一种地震数据的处理方法的流程图，如图1所示，该方法包括如下步骤：

步骤s102，将目标地震数据表示为高斯束叠加形式的地震数据，其中，高斯束叠加形式的地震数据中包括目标地震数据的波形信息。

在远场假设条件下，地震场波可以认为是由局部平面波组成，通过高斯束的叠加可以得到平面波，因此，地震数据可以表示为高斯束叠加的形式。

步骤s104，将高斯束叠加形式的地震数据进行稀疏处理，处理得到波形函数矩阵。

步骤s106，基于处理之后得到的波形函数矩阵进行地震数据的重构，得到重构之后的地震数据。

在本发明实施例中，首先将目标地震数据表示为高斯束叠加形式的地震数据，其中，高斯束叠加形式的地震数据包括目标地震数据包括波形信息；接下来，将高斯束叠加形式的地震数据进行稀疏处理，处理得到波形函数矩阵；最后，基于波形函数矩阵进行地震数据的重构，以完成地震数据的正则化处理。通过波形函数矩阵对地震数据进行重构，能够将空间随机采样的地震数据差值到规则网格上，以消除随机噪音和特定倾角的噪音，进而缓解了通过传统的正则化处理方法对地震数据进行处理之后处理效果较差，以及处理效率较慢的技术问题，从而实现了提高了地震数据的正则化处理效果的技术效果。

在本发明实施例的一个可选实施方式中，将目标地震叠加数据表示为高斯束叠加形式的地震数据包括如下步骤：

步骤s1021，获取高斯束基函数的半宽w和高斯束中任意相邻的两个束中心之间的束中心间距a。

首先，确定高斯束基函数的半宽w和高斯束中任意相邻的两个束中心之间的束中心间距a。

步骤s1022，根据半宽和束中心间距计算高斯束的束中心位置xm，其中，xm＝am，m＝1,2,…,m，m为束中心的个数。

接下来，就可以根据公式xm＝am确定高斯束的束中心位置。需要说明的是，为了加快计算效率，可以对地震数据进行分块处理，处理得到多个子区块。其中，为了减少误差，相邻子区块有一定的重合。

步骤s1023，通过下述公式确定目标地震数据中每道地震数据与束中心相对应的基函数，其中，tk为高斯束基函数的束中心所处的时间位置，px为水平慢度，v为地震速度常数，ωr为相关低频，该公式表示为：

在对地震数据进行分块处理之后，就可以在多个子区块中任意选取一个子区块进行处理。为了提高扫描的稳定性，可以首先通过傅里叶变换对子区块中的每道地震数据进行插值加密；然后，计算每道地震数据到高斯束束中心的距离；接下来，使用下述计算公式计算与高斯束束中心向对应的基函数，即，通过逆傅里叶变换的形式得到时间-空间域的高斯束基函数。

需要说明的是，px满足下述公式px＝sinθ/v，v为一个地震速度常数；且ω满足下述条件：0<ωl≤|ω|≤ωh<2πfnyquist，ωl和ωh分别为地震数据的最低和最高波动频率，fnyquist为尼奎斯特频率。基函数的半宽w和相邻束中心的间距a通过如下公式计算，va为地下介质速度的平均值。

步骤s1024，通过目标地震数据中每道地震数据与束中心相对应的基函数，将目标地震数据表示为所述高斯束形式的地震数据。

综上，地震数据d(x,t)可以表示为：d(x,t)＝∑γcγψγ(x,t)，其中，ψγ(x,t)为基函数ψ(x,t；xm,tk,px)的简写，cγ是基函数的系数，γ表示参数(xm,tk,px)的集合。

进一步地，公式d(x,t)＝∑γcγψγ(x,t)可以表示成矩阵-向量方程组的形式dj＝acj，其中，为带有波形信息系数向量(也即，上述波形信息矩阵中的波形信息)，j＝1,2,…,l。矩阵a是由高斯束基函数形成的有限秩算子，c＝[c1,c2,…,cl]是系数向量cj所组成的系数矩阵(也即，上述波形信息矩阵)，d＝[d1,d2,…,dl]为矩阵形式所表示的地震数据。

在将目标地震数据表示为高斯束叠加形式的地震数据之后，就可以将高斯束叠加形式的地震数据进行稀疏处理，处理得到波形函数矩阵，具体地，包括如下步骤：

步骤s1041，基于波形信息矩阵建立稀疏矩阵优化模型，其中，波形信息矩阵中包括目标地震数据中的波形信息；

通过上述描述可知，已将目标地震数据转化为如公式dj＝acj的描述形式。此时，就需要计算向量cj。为了求解向量cj，在本发明实施例中建立了lq-lp(q>0,p>0)矩阵化模型：

其中，是稳定因子，αj是正则化参数。当0<p≤1时，模型(1)将获得系数向量cj中非零个数最少的解。但是，发明人发现，地震数据具有列相关性，即地震信号在矩阵表达上横向之间具有联系，基于该认识，新建立约束极小化矩阵优化模型，即l2,q-l2,p模型q>0,p>0，其中，l2,p范数定义做p∈(0,1]，c^k为矩阵c的第k行，‖c^k‖2表示欧几里得范数，表示一个正则化矩阵，αk>0代表矩阵c的第k行的正则化参数。可以把矩阵优化问题转化为一个迹优化问题，mincj(c)＝tr((ac-d)^t(ac-d))+tr(c^thm)，其中，tr(·)表示矩阵的迹运算,h是一个对角矩阵，定义作：c^k(k＝1,2,…,k)为矩阵c的第k行。下面的核心问题是快速求解迹优化问题。

但是，直接求解向量cj的非零值并且要求非零值的个数最少是不可行的。为了克服这种np完全问题，需要充分结合地震数据的特点。因此，在本发明实施例中，需要结合地震数据的特点对公式dj＝acj中的波形信息矩阵c进行稀疏处理，以降低稀疏分解的计算量。

步骤s1042，对稀疏矩阵优化模型中的波形信息矩阵进行稀疏处理，得到处理之后的稀疏矩阵优化模型，其中，稀疏处理之后的波形信息矩阵中的每个数值均不等于预设数值。

进一步地，对稀疏矩阵优化模型中的波形信息矩阵进行稀疏处理包括如下步骤：

步骤s11，对目标地震数据中的高斯束束中心相邻道上进行相关性度量分析，得到地震子波，其中，地震子波中包含地震子波所处的时间和倾角方向；

步骤s12，根据地震子波所处的时间和倾角方向确定波形信息矩阵中数值等于预设数值的向量；

步骤s13，在波形信息矩阵中删除数值等于预设数值的向量，得到处理之后的稀疏矩阵优化模型。

具体地，从地下反射(或绕射)而来的地震子波，具有时间和方向特性，在局部上表现为同相轴的形式。因此可以通过倾角扫描的方式，提取地震子波所处的时间和倾角方向，得到向量cj中非零值的个数和所处的位置，从而降低了稀疏分解的计算量。

其中，由于地震子波在局部上表现为同相轴的形式，在时间上连续，倾角方向相同，因此可以通过下述公式对地震数据进行相关性度量的分析，以提取地震子波：

其中，n为时间窗的半宽，xm为束中心的位置。

通过倾角扫描的方式，以得到波形向量cj中非零值的个数和位置，以实现进而达到稀疏的目的。其中，非零值即为不等于预设数值的向量，也就是说，在本发明实施例中，零即为上述预设数值。

也就是说，通过上述处理方式，就能够将上述迹优化问题mincj(c)＝tr((ac-d)^t(ac-d))+tr(c^thm)转化为求解一个降维的极小化问题，即其中，为去掉其向量c中全零值行后的系数矩阵，为去掉矩阵a中对应列后的矩阵，矩阵的规模比矩阵a小很多，大大降低了反演的计算量。

在上述步骤s1042中，对稀疏矩阵优化模型中的波形信息矩阵进行稀疏处理，并得到处理之后的稀疏矩阵优化模型之后，就可以通过步骤s1043对处理之后的稀疏矩阵优化模型进行稀疏求解，求解得到波形函数矩阵。

具体地，对处理之后的稀疏矩阵优化模型进行稀疏求解，求解得到波形函数矩阵包括如下步骤：

步骤s21，反复通过第一公式和第二公式对稀疏矩阵优化模型进行迭代计算，求解得到波形函数矩阵，直至满足第三公式时，确定第二公式中的输出结果为波形函数矩阵，其中，α1，α2和αk为对角线矩阵λ中的数值，gk为迭代变量，p为稀疏度参数，k为迭代指数，为矩阵算子，d为数据矩阵，为稀疏处理之后的波形信息矩阵，im表示m维恒等算子，为矩阵的第k行向量，

具体地，在对稀疏处理之后的稀疏矩阵优化模型进行稀疏求解时，首先计算函数j(c)关于稀疏处理之前波形信息矩阵c的导数，得到kkt方程，该方程表示如下：

然后，定义令计算可得其中，为矩阵的第k行。

接下来，利用sherman-morrison-woodbury公式可以得到问题的极小解，该极小解满足下述公式：其中，im表示m维恒等算子。

最后，构造的迭代公式，具体构造过程如下：

假设在第k次迭代过程中，生成迭代变量gk，则第k+1次迭代可以利用下述公式得到：

需要说明的是，k为迭代指数，其中，迭代指数k的初值等于0，其中，k＝0时，对应着

也就是说，在本发明实施例中，k＝1,2,…，执行for循环：

其中，约束条件为当ρk≤∈时，终止迭代，输出波形矩阵∈为大于零的常数，即为上述波形函数矩阵。

在将目标地震数据进行稀疏处理，处理得到波形函数矩阵之后，就可以基于分解之后得到的波形函数矩阵进行地震数据的重构，得到重构之后的地震数据。

作为一个可选实施方式中，可以通过公式d(x,t)＝∑γcγψγ(x,t)进行地震数据的重构，得到重构之后的地震数据，其中，为波形函数矩阵，d(x,t)为重构之后的地震数据。需要说明的是，重构之后的地震数据位于规则网格上，消除了随机噪音，并且可以通过控制倾角，消除特定倾角的噪音。

综上，在本发明实施例中，使用高斯束基函数稀疏分解及重构地震数据，并建立l2,q-l2,p(q>0,0<p≤1)范数极小化模型。其中，在求解范数极小化的过程中，为了克服np完全问题，通过倾角扫描的方式得到解的先验信息，并使用sherman-morrison-woodbury公式构造迭代算法计算波形函数矩阵，加快计算速度，实现稀疏分解。本方法可以很自然的应用于随机采样的地震数据。

实施例二

本发明实施例还提供了一种地震数据的处理装置，该地震数据的处理装置主要用于执行本发明实施例上述内容所提供的地震数据的处理方法，以下对本发明实施例提供的地震数据的处理装置做具体介绍。

图2是根据本发明实施例的一种地震数据的处理装置的示意图，如图2所示，该地震数据的处理装置主要包括：第一处理单元21，第二处理单元22和重构单元23，其中：

第一处理单元21，用于将目标地震数据表示为高斯束叠加形式的地震数据，其中，所述高斯束叠加形式的地震数据中包括所述目标地震数据的波形信息；

第二处理单元22，用于将所述高斯束叠加形式的地震数据进行稀疏处理，处理得到波形函数矩阵；

重构单元23，用于基于处理之后得到的波形函数矩阵进行地震数据的重构，得到重构之后的地震数据。

在本发明实施例中，首先将目标地震数据表示为高斯束叠加形式的地震数据，其中，高斯束叠加形式的地震数据包括目标地震数据包括波形信息；接下来，将高斯束叠加形式的地震数据进行稀疏处理，处理得到波形函数矩阵；最后，基于波形函数矩阵进行地震数据的重构，以完成地震数据的正则化处理。通过波形函数矩阵对地震数据进行重构，能够将空间随机采样的地震数据差值到规则网格上，以消除随机噪音和特定倾角的噪音，进而缓解了通过传统的正则化处理方法对地震数据进行处理之后处理效果较差，以及处理效率较慢的技术问题，从而实现了提高了地震数据的正则化处理效果的技术效果。

图3是根据本发明实施例的另一种可选地地震数据的处理装置的示意图，如图3所示，第二处理单元22包括：建立模块31，稀疏处理模块32和求解模块33，其中，

建立模块31用于基于波形信息矩阵建立稀疏矩阵优化模型，其中，波形信息矩阵中包括波形信息；稀疏处理模块32用于对稀疏矩阵优化模型中的波形信息矩阵进行稀疏处理，得到处理之后的稀疏矩阵优化模型，以使稀疏处理之后的波形信息矩阵中的每个数值均不等于预设数值；求解模块33用于对处理之后的稀疏矩阵优化模型进行稀疏求解，求解得到波形函数矩阵。

可选地，稀疏处理模块32包括：分析子模块，用于对所述目标地震数据中的高斯束束中心相邻道上进行相关性度量分析，得到地震子波，其中，所述地震子波中包含所述地震子波所处的时间和倾角方向；确定子模块，用于根据地震子波所处的时间和倾角方向确定波形信息矩阵中数值等于预设数值的向量；以及删除子模块，用于在波形信息矩阵中删除数值等于预设数值的向量，得到处理之后的稀疏矩阵优化模型。

可选地，求解模块33包括：计算子模块，用于反复通过第一公式和第二公式对稀疏矩阵优化模型进行迭代计算，求解得到波形函数矩阵，直至满足第三公式时，确定第二公式中的输出结果为波形函数矩阵，其中，α1，α2和αk为对角线矩阵λ中的数值，gk为迭代变量，p为稀疏度参数，k为迭代指数，为矩阵算子，d为数据矩阵，为稀疏处理之后的波形信息矩阵，im表示m维恒等算子，为矩阵的第k行向量，

图4是根据本发明实施例的另一种可选地地震数据的处理装置的示意图，如图4所示，重构单元23包括：重构模块41，重构模块用于通过公式d(x,t)＝∑γcγψγ(x,t)进行地震数据的重构，得到重构之后的地震数据，其中，为波形函数矩阵，d(x,t)为重构之后的地震数据。

图5是根据本发明实施例的另一种可选地地震数据的处理装置的示意图，如图5所示，第一处理单元21包括：获取模块51，计算模块52，第一确定模块53和第二确定模块54，其中，

获取模块51用于获取高斯束基函数的半宽w和高斯束中任意相邻的两个束中心之间的束中心间距a；计算模块52用于根据半宽和束中心间距计算高斯束的束中心位置xm，其中，xm＝am，m＝1,2,…,m，m为束中心的个数；第一确定模块53用于通过公式确定目标地震数据中每道地震数据与束中心相对应的基函数，其中，tk为高斯束基函数的束中心所处的时间位置，px为水平慢度，v为地震速度常数，ωr为相关低频；第二确定模块54用于通过所述目标地震数据中每道地震数据与所述束中心相对应的基函数，将目标地震数据表示为高斯束叠加形式的地震数据。

最后应说明的是：以上各实施例仅用以说明本发明的技术方案，而非对其限制；尽管参照前述各实施例对本发明进行了详细的说明，本领域的普通技术人员应当理解：其依然可以对前述各实施例所记载的技术方案进行修改，或者对其中部分或者全部技术特征进行等同替换；而这些修改或者替换，并不使相应技术方案的本质脱离本发明各实施例技术方案的范围。