提高分辨率的地震数据处理方法、介质、终端和装置与流程

【技术领域】

本发明涉及地震数据处理领域，尤其涉及一种提高分辨率的地震数据处理方法、介质、终端和装置。

背景技术：

高分辨率是地下浅层空间(0-50米)探测开发追求的永恒目标。由于浅层地裂缝规模和尺度较小，勘探浅表地裂缝需要极宽的频带，而人工激发的频带较窄，低频和高频分量能量较弱，通常采用室内数据处理方法拓宽频带。当前提高浅层地震勘探地裂缝分辨率最常用的方法，如反褶积、反q滤波等拓频技术，可以不同程度的提高分辨率，但由于浅表地层的不均一性，物理分层不明显，反褶积算子和q值不易准确求得，拓宽频带的同时，子波旁瓣增多、增大，降低信噪比和保真度。

技术实现要素：

本发明提供了一种提高分辨率的地震数据处理方法、介质、终端和装置，解决了以上所述的技术问题。

本发明解决上述技术问题的技术方案如下：一种提高分辨率的地震数据处理方法，包括以下步骤：

s1，通过加伯正变换对地震数据的原始时间域信号进行谱分解生成初始时频谱，并从所述初始时频谱中提取至少一个平均q值；

s2，采用内点算法对所述至少一个平均q值进行反演，生成所述地震数据中每层地层的层间q值；

s3，采用所述层间q值对所述初始时频谱进行反q滤波，生成新的时频谱；

s4，通过加伯逆变换将所述新的时频谱转换为反q滤波补偿后的时间域信号。

在一个优选实施方式中，所述s1具体包括以下步骤：

s101，通过加伯正变换对地震数据的原始时间域信号s(t)进行谱分解生成初始二维时频谱s(τ,ω)，其中为τ旅行时间，ω为角频率；

s102，将所述初始二维时频谱s(τ,ω)转化为一维时频谱s(χ)，并对所述一维时频谱s(χ)进行最大幅值归一化后生成振幅补偿算子如下：

其中χ＝ωτ，smax为s(χ)对应的最大幅值；

s103，建立地震波品质因子q和速度v之间经验公式的初始解以及初始解的子空间并采用所述振幅补偿算子建立估计模型如下：

其中，σ²为稳定因子，a0为τ＝0时的振幅，常数q^-1为τ＝0和当前旅行时间τ之间的平均q值；

s104，在所述子空间内利用线性搜索算法生成至少一个不同时刻下所述估计模型对应的平均q值。

在一个优选实施方式中，通过s1得到m个时刻的平均q值，分别记为且所述地震数据包括n个地层，每层地层的层间q值分别记为则所述s2具体包括以下步骤：

s201，建立层间q值的最小化目标函数如下：

s.t.dq^-1/dt＝0(b)

其中，

s202，将所述最小化目标函数转化为以下矩阵形式，并求解生成每个地层的层间q值：

ax＝b

s.t.x＞0

其中，矩阵a和b的元素构成为：

a(m+n-1)×n＝[um×n,λv(n-1)×n]^t

其中，下标表示矩阵的维数，矩阵为m个平均q值，x为待反演的n个层间q值，λ为调谐参数，对所述最小化目标函数和约束条件起加权作用；

所述矩阵a的jacobian矩阵u为由产生的下三角矩阵，约束矩阵v是由元素(1,-1)构成的双对角矩阵。

在一个优选实施方式中，所述s3具体为：在地表波场延拓算法上引入所述层间q值，并对地震信号进行反q滤波生成粘弹介质对应的新时频谱如下：

以及生成弹性介质对应的新时频谱如下：

其中，σ²为稳定因子，ωh为调谐参数，γ＝1/πq。

本发明实施例的第二方面提供了一种计算机可读存储介质，存储有计算机程序，所述计算机程序被处理器执行时，实现以上所述的提高分辨率的地震数据处理方法。

本发明实施例的第三方面提供了一种提高分辨率的地震数据处理终端，包括所述的计算机可读存储介质和处理器，所述处理器执行所述计算机可读存储介质上的计算机程序时实现以上所述提高分辨率的地震数据处理方法的步骤。

本发明实施例的第四方面提供了一种提高分辨率的地震数据处理装置，包括提取模块、反演模块、滤波模块和逆变换模块，

所述提取模块用于通过加伯正变换对地震数据的原始时间域信号进行谱分解生成初始时频谱，并从所述初始时频谱中提取至少一个平均q值；

所述反演模块用于采用内点算法对所述至少一个平均q值进行反演，生成所述地震数据中每层地层的层间q值；

所述滤波模块用于采用所述层间q值对所述初始时频谱进行反q滤波，生成新的时频谱；

所述逆变换模块用于通过加伯逆变换将所述新的时频谱转换为反q滤波补偿后的时间域信号。

在一个优选实施方式中，所述提取模块具体包括：

第一变换单元，用于通过加伯正变换对地震数据的原始时间域信号s(t)进行谱分解生成初始二维时频谱s(τ,ω)，其中为τ旅行时间，ω为角频率；

第二变换单元，用于将所述初始二维时频谱s(τ,ω)转化为一维时频谱s(χ)，并对所述一维时频谱s(χ)进行最大幅值归一化后生成振幅补偿算子如下：其中χ＝ωτ，smax为s(χ)对应的最大幅值；

第一模型建立单元，用于建立地震波品质因子q和速度v之间经验公式的初始解以及初始解的子空间并采用所述振幅补偿算子建立估计模型如下：

其中，σ²为稳定因子，a0为τ＝0时的振幅，常数q^-1为τ＝0和当前旅行时间τ之间的平均q值；

第一解算单元，用于在所述子空间内利用线性搜索算法生成至少一个不同时刻下所述估计模型对应的平均q值。

在一个优选实施方式中，所述第一解算单元生成m个时刻的平均q值，分别记为且所述地震数据包括n个地层，每层地层的层间q值分别记为则所述反演模块具体包括：

第二模型建立单元，用于建立层间q值的最小化目标函数如下：

s.t.dq^-1/dt＝0(b)

其中，

第二解算单元，用于将所述最小化目标函数转化为以下矩阵形式，并求解生成每个地层的层间q值：

其中，矩阵a和b的元素构成为：

其中，下标表示矩阵的维数，矩阵为m个平均q值，x为待反演的n个层间q值，x＝[q1^-1,q2^-1,…,qn^-1]^t，λ为调谐参数，对所述最小化目标函数和约束条件起加权作用；

所述矩阵a的jacobian矩阵u为由产生的下三角矩阵，约束矩阵v是由元素(1,-1)构成的双对角矩阵。

在一个优选实施方式中，所述滤波模块具体用于在地表波场延拓算法上引入所述层间q值，并对地震信号进行反q滤波生成粘弹介质对应的新时频谱如下：

以及用于生成弹性介质对应的新时频谱如下：

其中，σ²为稳定因子，ωh为调谐参数，γ＝1/πq。

本发明提出了一种提高分辨率的地震数据处理方法、介质、终端和装置，具有有益效果如下：

(1)叠前和叠后数据均可适用；

(2)信噪比保持性好，可拓宽低频和高频的弱反射信号；

(3)相对振幅关系及振幅能量级别保持不变。

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

【附图说明】

为了更清楚地说明本发明实施例的技术方案，下面将对实施例中所需要使用的附图作简单地介绍，应当理解，以下附图仅示出了本发明的某些实施例，因此不应被看作是对范围的限定，对于本领域普通技术人员来讲，在不付出创造性劳动的前提下，还可以根据这些附图获得其他相关的附图。

图1是实施例1提供的提高分辨率的地震数据处理方法的流程示意图；

图2是采用实施例1提供的地震数据处理方法处理前后的地震数据图；

图3是实施例2提供的提高分辨率的地震数据处理装置的结构示意图；

图4是实施例3提供的提高分辨率的地震数据处理终端的结构示意图。

【具体实施方式】

为了使本发明的目的、技术方案和有益技术效果更加清晰明白，以下结合附图和具体实施方式，对本发明进行进一步详细说明。应当理解的是，本说明书中描述的具体实施方式仅仅是为了解释本发明，并不是为了限定本发明。

图1是实施例1提供的提高分辨率的地震数据处理方法的流程示意图，如图1所示，包括以下步骤：

s1，通过加伯正变换对地震数据的原始时间域信号进行谱分解生成初始时频谱，并从初始时频谱中提取至少一个平均q值。作为地震波的动力学参数之一，品质因子q值和速度v都表示为岩石的一种物理性质，尤其品质因子q值可以表征地层的衰减程度。当前已有许多针对两者统计关系的研究，并建立了二者之间的经验公式。如waters(1978)和李庆忠(1993)分别针对不同地区、不同岩性进行统计建立了q-v之间的经验公式，形式如下：

q＝c·vⁿ(1)

式中c为比例常数，n为指数因子。由于在地震数据处理中要做大量的速度分析，因此很容易得到速度v，进而得到用于q值迭代的初始解，并可建立初始解的子空间然后在时频域进行子空间γ内的q值分析。具体包括以下步骤：

s101，对于一个给定的地震数据的原始时间域信号s(t)，利用加伯正变换，即gabor正变换产生初始二维时频谱s(τ,ω)，式中τ为旅行时间，ω为角频率。具体实施时，可采用短时窗gabor变换。

s102，然后令χ＝ωτ，将初始二维时频谱s(τ,ω)转化为一维时频谱s(χ)，其物理意义即为原始时间域信号s(t)随χ的一维振幅衰减度量，按最大幅值归一化得其中smax为s(χ)对应的最大幅值。因此得到振幅补偿算子为：

式中，λd为振幅补偿算子，σ²为稳定因子，χa是衰减系数a在χ处的取值。

s103，考虑衰减时，理论振幅衰减表达式为：

式中，a0为τ＝0时的振幅，常数q^-1为τ＝0和当前旅行时间τ之间的平均q值，因此振幅衰减补偿算子为：

最终子空间γ内q值估计转化为最小化问题，即建立估计模型如下：

s104，在子空间内利用线性搜索算法生成至少一个不同时刻下估计模型对应的平均q值，记为

类似于速度扫描，通过上述步骤可以得到一系列不同时刻tm下的平均q值记为然后可以执行步骤s2，采用内点算法对该m个平均q值进行层间q值反演,如图1所示。内点法(interiorpointmethod)是一种用于求解带约束的优化命题的方法，其特点是将所构造的无约束目标函数——惩罚函数定义在可行域内，并在可行域内求惩罚函数的极值点，从而在进行层间q值反演过程中不依赖初始解同时保证反演得到的层间q值准确。步骤s2具体包括以下步骤：

s201，假设地层按照层间q值分为n层，记则由层间q值计算平均q值的计算公式为：

式中，表示在时刻tm计算的平均q值，假定层厚为δt，则有：

tm＝(n-1)δt+cδt(7)

其中cδt(c∈(0,1])为假定层厚所对应的时间，则式(6)变为：

因此，层间q值的计算转化为最小化问题，即建立最小化目标函数：

式中方程(9(b))为约束条件，为通过前面步骤得到的平均q值。

s202，对于该最小化问题的求解，首先将方程(9)转化为矩阵相乘形式为：

矩阵a和b的元素构成为：

其中，下标表示矩阵的维数，矩阵为m个平均q值，x为待反演的n个层间q值，x＝[q1^-1,q2^-1,…,qn^-1]^t，λ为调谐参数，对最小化目标函数和约束条件起加权作用。构成矩阵a的jacobian矩阵u为由方程(8)产生的下三角矩阵，约束矩阵v是由元素(1,-1)构成的双对角矩阵。

然后执行s3，采用层间q值对初始时频谱进行反q滤波，生成新的时频谱。具体实现过程中，平面波在粘弹介质中传播其单程波解析解为：

r(z+δzω)＝p(z，ω)exp[-ik(ω)δz](12)

式中，i为虚数单位，ω为角频率，k(ω)为波数。地层的q效应通过引入复波数表示为：

其中，qr，vr为参考频率的品质因子和相速度，ωh为调谐参数，γ＝1/πq。

因此，反q滤波表达为：

其中两个指数项分别补偿速度频散引起的相位影响和能量衰减引起的振幅影响。

对于q值随旅行时间连续变化时，将地表波场延拓到时间深度τ有：

其中，γ(τ)＝1/πq(τ)，进一步变形方程(15)得：

式中：

λ(τ，ω)为振幅补偿算子，β(τ，ω)为振幅衰减算子，σ²为稳定因子，通过经验公式σ²＝exp[-(0.23glim+1.63)]得到。参数glim为规定的增益限制，一般为0～70db。

由方程(15)对于弹性介质有：

则方程(15)变为：

方程(19)为基于gabor变换的提高分辨率方程，对于不同时间深度的波场，对地表地震记录进行gabor变换，再结合振幅补偿算子及相位补偿算子在时频域进行因子调节，最后对所得的新的时频谱p(τ,ω)利用逆gabor变换可得到经由提高分辨率后处理的地震记录p(t)。图2是采用以上地震数据处理方法处理前后的地震数据图，如图2所示，处理前数据的分辨率较低，处理后有效明显提高了各频段弱信号，综合提高了数据的分辨率，同时振幅保持较好，明显地增强了地震数据的可解释性。

应理解，上述实施例中各步骤的序号的大小并不意味着执行顺序的先后，各过程的执行顺序应以其功能和内在逻辑确定，而不应对本发明实施例的实施过程构成任何限定。

本发明实施例还提供了一种计算机可读存储介质，存储有计算机程序，计算机程序被处理器执行时，实现以上的提高分辨率的地震数据处理方法。图3是实施例2提供的提高分辨率的地震数据处理装置的结构示意图，如图3所示，包括提取模块100、反演模块200、滤波模块300和逆变换模块400，

提取模块100用于通过加伯正变换对地震数据的原始时间域信号进行谱分解生成初始时频谱，并从初始时频谱中提取至少一个平均q值；

反演模块200用于采用内点算法对至少一个平均q值进行反演，生成地震数据中每层地层的层间q值；

滤波模块300用于采用层间q值对初始时频谱进行反q滤波，生成新的时频谱；

逆变换模块400用于通过加伯逆变换将新的时频谱转换为反q滤波补偿后的时间域信号。

优选实施例中，提取模块100具体包括：

第一变换单元101，用于通过加伯正变换对地震数据的原始时间域信号s(t)进行谱分解生成初始二维时频谱s(τ,ω)，其中为τ旅行时间，ω为角频率；

第二变换单元102，用于将初始二维时频谱s(τ,ω)转化为一维时频谱s(χ)，并对一维时频谱s(χ)进行最大幅值归一化后生成振幅补偿算子如下：其中χ＝ωτ，smax为s(χ)对应的最大幅值；

第一模型建立单元103，用于建立地震波品质因子q和速度v之间经验公式的初始解以及初始解的子空间并采用振幅补偿算子建立估计模型如下：

其中，σ²为稳定因子，a0为τ＝0时的振幅，常数q^-1为τ＝0和当前旅行时间τ之间的平均q值；

第一解算单元104，用于在子空间内利用线性搜索算法生成至少一个不同时刻下估计模型对应的平均q值。

另一优选实施例中，第一解算单元104生成m个时刻的平均q值，分别记为且地震数据包括n个地层，每层地层的层间q值分别记为{qi^-1,i＝1,2,…n}，则反演模块200具体包括：

第二模型建立单元201，用于建立层间q值的最小化目标函数如下：

其中，

第二解算单元202，用于将最小化目标函数转化为以下矩阵形式，并求解生成每个地层的层间q值：

其中，矩阵a和b的元素构成为：

其中，下标表示矩阵的维数，矩阵为m个平均q值，x为待反演的n个层间q值，x＝[q1^-1,q2^-1,…,qn^-1]^t，λ为调谐参数，对最小化目标函数和约束条件起加权作用；

矩阵a的jacobian矩阵u为由产生的下三角矩阵，约束矩阵v是由元素(1,-1)构成的双对角矩阵。

优选实施例中，滤波模块300具体用于在地表波场延拓算法上引入层间q值，并对地震信号进行反q滤波生成粘弹介质对应的新时频谱如下：

以及用于生成弹性介质对应的新时频谱如下：

其中，

σ²为稳定因子，ωh为调谐参数，γ＝1/πq。

本发明实施例还提供了一种提高分辨率的地震数据处理终端，包括的计算机可读存储介质和处理器，处理器执行计算机可读存储介质上的计算机程序时实现以上提高分辨率的地震数据处理方法的步骤。图4是本发明实施例3提供的提高分辨率的地震数据处理终端的结构示意图，如图4所示，该实施例的提高分辨率的地震数据处理终端8包括：处理器80、可读存储介质81以及存储在可读存储介质81中并可在处理器80上运行的计算机程序82。处理器80执行计算机程序82时实现上述各个方法实施例中的步骤，例如图1所示的s1至s4。或者，处理器80执行计算机程序82时实现上述各装置实施例中各模块的功能，例如图3所示模块100至400的功能。

示例性的，计算机程序82可以被分割成一个或多个模块，一个或者多个模块被存储在可读存储介质81中，并由处理器80执行，以完成本发明。一个或多个模块可以是能够完成特定功能的一系列计算机程序指令段，该指令段用于描述计算机程序82在提高分辨率的地震数据处理终端8中的执行过程。

提高分辨率的地震数据处理终端8可包括，但不仅限于，处理器80、可读存储介质81。本领域技术人员可以理解，图3仅仅是提高分辨率的地震数据处理终端8的示例，并不构成对提高分辨率的地震数据处理终端8的限定，可以包括比图示更多或更少的部件，或者组合某些部件，或者不同的部件，例如提高分辨率的地震数据处理终端还可以包括电源管理模块、运算处理模块、输入输出设备、网络接入设备、总线等。

所称处理器80可以是中央处理单元(centralprocessingunit，cpu)，还可以是其他通用处理器、数字信号处理器(digitalsignalprocessor，dsp)、专用集成电路(applicationspecificintegratedcircuit，asic)、现成可编程门阵列(field-programmablegatearray，fpga)或者其他可编程逻辑器件、分立门或者晶体管逻辑器件、分立硬件组件等。通用处理器可以是微处理器或者该处理器也可以是任何常规的处理器等。

可读存储介质81可以是提高分辨率的地震数据处理终端8的内部存储单元，例如提高分辨率的地震数据处理终端8的硬盘或内存。可读存储介质81也可以是提高分辨率的地震数据处理终端8的外部存储设备，例如提高分辨率的地震数据处理终端8上配备的插接式硬盘，智能存储卡(smartmediacard，smc)，安全数字(securedigital，sd)卡，闪存卡(flashcard)等。进一步地，可读存储介质81还可以既包括提高分辨率的地震数据处理终端8的内部存储单元也包括外部存储设备。可读存储介质81用于存储计算机程序以及提高分辨率的地震数据处理终端所需的其他程序和数据。可读存储介质81还可以用于暂时地存储已经输出或者将要输出的数据。

所属领域的技术人员可以清楚地了解到，为了描述的方便和简洁，仅以上述各功能单元、模块的划分进行举例说明，实际应用中，可以根据需要而将上述功能分配由不同的功能单元、模块完成，即将装置的内部结构划分成不同的功能单元或模块，以完成以上描述的全部或者部分功能。实施例中的各功能单元、模块可以集成在一个处理单元中，也可以是各个单元单独物理存在，也可以两个或两个以上单元集成在一个单元中，上述集成的单元既可以采用硬件的形式实现，也可以采用软件功能单元的形式实现。另外，各功能单元、模块的具体名称也只是为了便于相互区分，并不用于限制本申请的保护范围。上述系统中单元、模块的具体工作过程，可以参考前述方法实施例中的对应过程，在此不再赘述。

在上述实施例中，对各个实施例的描述都各有侧重，某个实施例中没有详述或记载的部分，可以参见其它实施例的相关描述。

本领域普通技术人员可以意识到，结合本文中所公开的实施例描述的各示例的单元及方法步骤，能够以电子硬件、或者计算机软件和电子硬件的结合来实现。这些功能究竟以硬件还是软件方式来执行，取决于技术方案的特定应用和设计约束条件。专业技术人员可以对每个特定的应用来使用不同方法来实现所描述的功能，但是这种实现不应认为超出本发明的范围。

在本发明所提供的实施例中，应该理解到，所揭露的装置/终端设备和方法，可以通过其它的方式实现。例如，以上所描述的装置/终端设备实施例仅仅是示意性的，例如，模块或单元的划分，仅仅为一种逻辑功能划分，实际实现时可以有另外的划分方式，例如多个单元或组件可以结合或者可以集成到另一个系统，或一些特征可以忽略，或不执行。另一点，所显示或讨论的相互之间的耦合或直接耦合或通讯连接可以是通过一些接口，装置或单元的间接耦合或通讯连接，可以是电性，机械或其它的形式。

作为分离部件说明的单元可以是或者也可以不是物理上分开的，作为单元显示的部件可以是或者也可以不是物理单元，即可以位于一个地方，或者也可以分布到多个网络单元上。可以根据实际的需要选择其中的部分或者全部单元来实现本实施例方案的目的。

另外，在本发明各个实施例中的各功能单元可以集成在一个处理单元中，也可以是各个单元单独物理存在，也可以两个或两个以上单元集成在一个单元中。上述集成的单元既可以采用硬件的形式实现，也可以采用软件功能单元的形式实现。

本发明并不仅仅限于说明书和实施方式中所描述，因此对于熟悉领域的人员而言可容易地实现另外的优点和修改，故在不背离权利要求及等同范围所限定的一般概念的精神和范围的情况下，本发明并不限于特定的细节、代表性的设备和这里示出与描述的图示示例。