CN104597489A - 一种震源子波优化设置方法和装置 - Google Patents

Abstract

本发明提供一种震源子波优化设置方法和装置，其中，该方法包括以下步骤：从地震资料中提取地震子波作为目标子波和初始震源子波；利用目标子波和初始震源子波，采用正则化反演算法更新初始震源子波；当满足预设收敛条件时输出优化的震源子波。本发明采用了正则化算法反演得到优化的震源子波，可以将从地震资料中提取的初始地震子波直接作为震源子波输入到二维偏移和反演程序中，从而可以为声波和弹性波数值模拟、偏移成像、全波形反演提供可靠的震源子波。

Description

一种震源子波优化设置方法和装置

技术领域

本发明涉及地震勘探技术领域，特别涉及一种震源子波优化设置方法和装置。

背景技术

随着勘探对象的日趋复杂，地震成像技术迫切地需要一种适用于任意复杂构造介质成像的精确偏移方法，而逆时偏移(Reverse Time Migration，RTM)方法的出现给地球物理学家带来了希望。基于双程波动方程的逆时偏移方法具备理论简单、成像精确且无倾角限制、适应任意复杂速度模型等诸多优点。目前，逆时偏移方法已被广泛应用于复杂介质的成像处理中。

另外，地球物理学的基本问题就是用地面或井中的各种观测资料定性或定量地对地下的地质结构和矿产资源做出判断，即解地球物理反演问题。地震反演的主要任务就是利用地震资料，反推地下的波阻抗或速度信息，进行储层参数估算、储层预测和油藏描述，为油气勘探提供可靠的基础资料。地震反演按不同的标准可以分为不同的类型。按照所用的信息不同分为基于射线理论的走时类反演、基于AVO(AmplitudeVersus Offset)理论的振幅类反演和基于波动方程理论的全波形反演。全波形反演充分利用地震波运动学和动力学信息，能得到精度更高的岩性参数，具有揭示复杂地质构造和岩性细节信息的潜力。

逆时偏移和全波形反演的关键步骤都是波动方程数值模拟。常用的数值模拟方法包括有限元素法，有限差分法和伪谱法等。其中，有限差分法计算量小，易于实现，能够较精确地模拟任意非均匀介质中的地震波场。由于超大的计算量和存储需求，全三维逆时偏移和全波形反演算法很难在短时间内被大规模应用于实际生产中。对于常见的水平层状地层，可以采用一维偏移和反演算法；而对于有一定延伸的地质体，在垂直于其走向的剖面，可以采用二维偏移和算法。然而，当采用二维偏移和反演方法处理实际资料时，存在数据不匹配的问题。野外实际观测到的资料是由三维点源激发产生的，而一维和二维算法分别对应面源和线源。由于不同维度间存在着波形和相位的畸变，从观测数据中提取的地震子波不能直接作为震源输入到二维偏移和反演程序中。另外，在一维和二维数值模拟中，地震记录的波形和加入的震源也是有差异的。这些都会大大降低成像和反演精度。近年来，多分量资料在地震勘探中被广泛应用。因此，发展适用于声波和弹性波正演模拟的震源优化算法具有较好应用前景。然而，现有技术中还没有提出一种方法可以为声波和弹性波数值模拟、偏移成像、全波形反演提供可靠的震源子波。

发明内容

本发明实施例提供了一种震源子波优化方法，用以为声波和弹性波数值模拟、偏移成像、全波形反演提供可靠的震源子波，包括：

从地震资料中提取地震子波作为目标子波和初始震源子波；

利用目标子波和初始震源子波，采用正则化反演算法更新初始震源子波；

当满足预设收敛条件时输出优化的震源子波。

在一个实施例中，所述利用目标子波和初始震源子波，采用正则化反演算法更新初始震源子波，包括：利用目标子波和初始震源子波，采用正则化反演算法更新用于二维声波有限差分数值模拟的初始震源子波；利用目标子波和初始震源子波，采用正则化反演算法更新用于弹性波有限差分数值模拟的初始震源子波。

在一个实施例中，所述利用目标子波和初始震源子波，采用正则化反演算法更新用于二维声波有限差分数值模拟的初始震源子波，按照如下公式计算：

s₁＝||J₁Δm-Δp||²+r||LΔm||²；

其中，

s₁为迭代过程中每一步的极小化目标函数值；

Δp为二维声波有限差分数值模拟值与实际观测值的残量，l维；

p^obs为实际观测到的波场值，l维；

p^cal为迭代过程中每次的正演模拟波场值，l维；

J₁为雅克比矩阵，l行q列；

Δm为每次的震源子波修正量，q维；

r为正则化系数。

在一个实施例中，所述收敛条件为Δp的模小于预设的极小值。

在一个实施例中，所述利用目标子波和初始震源子波，采用正则化反演算法更新用于弹性波有限差分数值模拟的初始震源子波，按照如下公式计算：

s₂＝||J₂Δm-Δv||²+r||LΔm||²；

其中，

s₂为迭代过程中每一步的极小化目标函数值；

Δv为弹性波有限差分数值模拟值与实际观测值的残量，包含水平分量和垂直分量，2l维；

v^obs为实际观测到的波场值，包含水平分量和垂直分量，2l维；

v^cal为迭代过程中每次的正演模拟波场值，包含水平分量和垂直分量，2l维；

J₂为雅克比矩阵，2l行q列；

Δm为每次的震源子波修正量，q维；

r为正则化系数。

在一个实施例中，所述收敛条件为Δv的模小于预设的极小值。

在一个实施例中，根据得到的优化的震源子波，通过声波和弹性波的有限差分方程数值模拟得到波形；

对比模拟得到的波形和目标子波的波形，对优化的震源子波进行验证。

本发明实施例还提供了一种震源子波优化设置装置，用以为声波和弹性波数值模拟、偏移成像、全波形反演提供可靠的震源子波，包括：

子波选取模块，用于从地震资料中提取地震子波作为目标子波和初始震源子波；

初始震源子波更新模块，用于利用目标子波和初始震源子波，采用正则化反演算法更新初始震源子波；

优化震源子波输出模块，用于当满足预设收敛条件时输出优化的震源子波。

在一个实施例中，所述初始震源子波更新模块具体用于，利用目标子波和初始震源子波，采用正则化反演算法更新用于二维声波有限差分数值模拟的初始震源子波；利用目标子波和初始震源子波，采用正则化反演算法更新用于弹性波有限差分数值模拟的初始震源子波。

在一个实施例中，该装置还包括：模拟波形获得模块，用于根据得到的优化的震源子波，通过声波和弹性波的有限差分方程数值模拟得到波形；

验证模块，用于对比模拟得到的波形和目标子波的波形，对优化的震源子波进行验证。

本发明实施例中，采用了正则化算法反演得到优化的震源子波，可以将从地震资料中提取的初始地震子波直接作为震源子波输入到二维偏移和反演程序中，从而可以为声波和弹性波数值模拟、偏移成像、全波形反演提供可靠的震源子波。

附图说明

此处所说明的附图用来提供对本发明的进一步理解，构成本申请的一部分，并不构成对本发明的限定。在附图中：

图1是本发明实施例提供的一种震源子波优化设置方法流程图；

图2是本发明实施例提供的一种震源子波优化设置装置示意图；

图3是本发明实施例提供的一种30HZ的雷克子波波形图；

图4是本发明实施例提供的一种采用30HZ的雷克子波作为初始震源子波情况下，一维声波方程有限差分数值模拟获得的波形图；

图5是本发明实施例提供的一种采用30HZ的雷克子波作为初始震源子波情况下，二维声波方程有限差分数值模拟获得的波形图；

图6是本发明实施例提供的一种采用30HZ的雷克子波作为初始震源子波情况下，三维声波方程有限差分数值模拟获得的波形图；

图7是本发明实施例提供的一种一维声波情况下以30HZ的雷克子波为目标子波反演得到的震源子波波形图；

图8是本发明实施例提供的一种二维声波情况下以30HZ的雷克子波为目标子波反演得到的震源子波波形图；

图9是本发明实施例提供的一种三维声波情况下以30HZ的雷克子波为目标子波反演得到的震源子波波形图；

图10是本发明实施例提供的一种20HZ的雷克子波波形图；

图11是本发明实施例提供的一种采用20HZ的雷克子波作为初始震源子波情况下，二维弹性波方程有限差分数值模拟获得的波形图(包括水平分量和垂直分量)；

图12是本发明实施例提供的一种二维弹性波情况下以20HZ的雷克子波为目标子波反演得到的震源子波波形图。

具体实施方式

为使本发明的目的、技术方案和优点更加清楚明白，下面结合实施方式和附图，对本发明做进一步详细说明。在此，本发明的示意性实施方式及其说明用于解释本发明，但并不作为对本发明的限定。

图1是本发明实施例提供的一种震源子波优化设置方法流程图，如图1所示，该方法包括：

步骤101：从地震资料中提取地震子波作为目标子波和初始震源子波；

步骤102：利用目标子波和初始震源子波，采用正则化反演算法更新初始震源子波；

步骤103：当满足预设收敛条件时输出优化的震源子波。

实际中，获得地震子波的方法有很多(如：直接观测法、自相关法、多项式求根法、测井资料求取子波法和对数分解法等)，本发明方法中不考虑地震子波的获取方法。

对于如何选择初始震源子波，本发明方法中将从地震资料中提取地震子波即作为目标子波，又作为初始震源子波。还可以选择其他的初始震源子波。好的反演算法应该对初始模型依赖性较小。本发明方法中我们采用的初始震源子波初值与真实震源子波相差较大。

具体实施时，我们利用目标子波和初始震源子波，采用正则化反演算法更新用于二维声波和弹性波有限差分数值模拟的初始震源子波。

二维声波方程如下：

$\frac{1}{v^2}\frac{{\∂}^{2}p}{\∂t^2}=\frac{{\∂}^{2}p}{\∂x^2}+\frac{{\∂}^{2}p}{\∂z^2}+s(x,z)---\left(1\right)$

其中：p为波场值，v为速度，s(x,z)为初始震源子波。

正演模拟中采用有限差分法，差分格式如下：

$\frac{{\∂}^{2}p}{\∂t^2}\≈\frac{1}{\mathrm{\Δ}{t}^2}[p_{i,j}^{n+1}-2p_{i,j}^n+p_{i,j}^{n-1}]---\left(2\right)$

$\frac{{\∂}^{2}p}{\∂x^2}\≈\frac{1}{h^2}[a_0{p}_{i,j}^n+\underset{m=1}{\overset{M}{\mathrm{\Σ}}}{a}_m(p_{i-m,j}^n+p_{i+m,j}^n)]---\left(3\right)$

$\frac{{\∂}^{2}p}{\∂z^2}\≈\frac{1}{h^2}[a_0{p}_{i,j}^n+\underset{m=1}{\overset{M}{\mathrm{\Σ}}}{a}_m(p_{i,j-m}^n+p_{i,j+m}^n)]---\left(4\right)$

其中，Δt是时间采样间隔，h是空间采样间隔，a_m为差分系数，M为算子长度。

已知模型、模拟参数及初始震源子波s(x,z)，通过方程(2)至(4)可以获得模型中任意一点的记录。相反，我们通过模型中某一点的波形(目标子波)也可以反推初始震源子波s(x,z)，反演算法如下。

利用Taylor展开把非线性问题线性化，忽略高阶项得：

式(5)中，Δp为二维声波有限差分数值模拟值与实际观测值的残量，l维；p^obs为实际观测到的波场值(目标子波的波场值)，l维；p^cal为迭代过程中每次的正演模拟波场值，l维；J₁为雅克比矩阵，l行q列，其元素为偏导数算子，可以通过数值离散求得；Δm为每次的震源子波修正量，q维。

采用阻尼最小二乘法，求s₁＝||J₁Δm-Δp||²的极小值。迭代过程中，解m_k+1＝m_k+Δm。ε²为阻尼因子，m_k为第k次迭代的震源子波。

为了提高反演算法的稳定性，我们采用正则化算法，即求s₁＝||J₁Δm-Δp||²+r||LΔm||²极小。其中，s₁为迭代过程中每一步的极小化目标函数值；r为正则化系数。迭代过程中，解m_k+1＝m_k+Δm。

正则化中，在残量Δp中加入了表示光滑度的量，这里使用参变量的三阶导数使原来的不适定问题正则化：

通过上述步骤进行正则化反演优化之后，当满足预设收敛条件时输出优化的震源子波，即当残量Δp的模小于某个极小值时，跳出循环，输出反演得到的优化的震源子波。

当获得优化的震源子波之后，我们还需要对其进行验证，说明本发明方法可以为声波和弹性波数值模拟、偏移成像、全波形反演提供可靠的震源子波。具体的，用优化的震源子波和方程(2)至(4)数值模拟得到某点的波形，对比模拟得到的波形和目标波形的差异。如果二者比较一致，则得到的震源子波是正确的、可靠的。

二维弹性介质速度-应力方程如下：

$\mathrm{\ρ}\frac{\∂v_x}{\∂t}=\frac{\∂{\mathrm{\τ}}_{\mathrm{xx}}}{\∂x}+\frac{\∂{\mathrm{\τ}}_{\mathrm{xz}}}{\∂z}---\left(7\right)$

$\mathrm{\ρ}\frac{\∂v_z}{\∂t}=\frac{\∂{\mathrm{\τ}}_{\mathrm{xz}}}{\∂x}+\frac{\∂{\mathrm{\τ}}_{\mathrm{zz}}}{\∂z}---\left(8\right)$

$\frac{\∂{\mathrm{\τ}}_{\mathrm{xx}}}{\∂t}=(\mathrm{\λ}+2\mathrm{\μ})\frac{\∂v_x}{\∂x}+\mathrm{\λ}\frac{\∂v_z}{\∂z}+s(x,z)---\left(9\right)$

$\frac{\∂{\mathrm{\τ}}_{\mathrm{zz}}}{\∂t}=\mathrm{\λ}\frac{\∂v_x}{\∂x}+(\mathrm{\λ}+2\mathrm{\μ})\mathrm{\λ}\frac{\∂v_z}{\∂z}+s(x,z)---\left(10\right)$

$\frac{\∂{\mathrm{\τ}}_{\mathrm{xz}}}{\∂t}=\mathrm{\μ}(\frac{\∂v_z}{\∂x}+\frac{\∂v_x}{\∂z})---\left(11\right)$

式中，v_x是水平速度矢量，v_z是垂直速度矢量，τ_xx是正应力在水平方向的分量，τ_zz是正应力在垂直方向的分量，τ_xz是剪切应力，λ(x,z)和μ(x,z)是拉梅系数。

正演模拟中采用有限差分法，差分格式如下：

$\mathrm{\ρ}\frac{v_{\mathrm{xi},j}^{n+1}-v_{\mathrm{xi},j}^n}{\mathrm{\Δt}}=\frac{1}{h}\underset{m=1}{\overset{M}{\mathrm{\Σ}}}{a}_m\left[\begin{array}{c}({\mathrm{\τ}}_{\mathrm{xxi}+m-1/2,j}^{n+1/2}-{\mathrm{\τ}}_{\mathrm{xxi}-m+1/2,j}^{n+1/2})\\ +({\mathrm{\τ}}_{\mathrm{xzi},j+m-1/2}^{n+1/2}-{\mathrm{\τ}}_{\mathrm{xzi},j-m+1/2}^{n+1/2})\end{array}\right]---\left(12\right)$

$\mathrm{\ρ}\frac{v_{\mathrm{zi}+1/2,j+1/2}^{n+1}-v_{\mathrm{zi}+1/2,j+1/2}^n}{\mathrm{\Δt}}=\frac{1}{h}\underset{m=1}{\overset{M}{\mathrm{\Σ}}}{a}_m\left[\begin{array}{c}({\mathrm{\τ}}_{\mathrm{xzi}+m,j+1/2}^{n+1/2}-{\mathrm{\τ}}_{\mathrm{xzi}-m+1,j+1/2}^{n+1/2})\\ +({\mathrm{\τ}}_{\mathrm{zzi}+1/2,j+m}^{n+1/2}-{\mathrm{\τ}}_{\mathrm{zzi}+1/2,j-m+1}^{n+1/2})\end{array}\right]---\left(13\right)$

$\frac{{\mathrm{\τ}}_{\mathrm{xxi}+1/2,j}^{n+3/2}-{\mathrm{\τ}}_{\mathrm{xxi}+1/2,j}^{n+1/2}}{\mathrm{\Δt}}=\frac{1}{h}\underset{m=1}{\overset{M}{\mathrm{\Σ}}}{a}_m\left[\begin{array}{c}(\mathrm{\λ}+2\mathrm{\μ})(v_{\mathrm{xi}+m,j}^{n+1}-v_{\mathrm{xi}-m+1,j}^{n+1})\\ +\mathrm{\λ}(v_{\mathrm{zi}+1/2,j+m-1/2}^{n+1}-v_{\mathrm{zi}+1/2,j-m+1/2}^{n+1})\end{array}\right]---\left(14\right)$

$\frac{{\mathrm{\τ}}_{\mathrm{zzi}+1/2,j}^{n+3/2}-{\mathrm{\τ}}_{\mathrm{zzi}+1/2,j}^{n+1/2}}{\mathrm{\Δt}}=\frac{1}{h}\underset{m=1}{\overset{M}{\mathrm{\Σ}}}{a}_m\left[\begin{array}{c}\mathrm{\λ}(v_{\mathrm{xi}+m,j}^{n+1}-v_{\mathrm{xi}-m+1,j}^{n+1})\\ +(\mathrm{\λ}+2\mathrm{\μ})(v_{\mathrm{zi}+1/2,j+m-1/2}^{n+1}-v_{\mathrm{zi}+1/2,j-m+1/2}^{n+1})\end{array}\right]---\left(15\right)$

$\frac{{\mathrm{\τ}}_{\mathrm{xzi},j+1/2}^{n+3/2}-{\mathrm{\τ}}_{\mathrm{xzi},j+1/2}^{n+1/2}}{\mathrm{\Δt}}=\frac{1}{h}\underset{m=1}{\overset{M}{\mathrm{\Σ}}}{a}_m\left[\begin{array}{c}\mathrm{\μ}(v_{\mathrm{zi}+m-1/2,j+1/2}^{n+1}-v_{\mathrm{zi}-m+1/2,j+1/2}^{n+1})\\ +\mathrm{\μ}(v_{\mathrm{xi},j+m}^{n+1}-v_{\mathrm{xi},j-m+1}^{n+1})\end{array}\right]---\left(16\right)$

通过方程(7)至(11)可以实现弹性波波动方程正演。弹性波数值模拟震源子波优化方法上述声波情况基本一致。

弹性波模拟得到的地震记录是两分量的(v_x,v_z)，同时采用水平分量和垂直分量反演震源。有如下形式：

式(17)中，Δv为弹性波有限差分数值模拟值与实际观测值的残量，包含水平分量和垂直分量，2l维；v^obs为实际观测到的波场值，包含水平分量和垂直分量，2l维；v^cal为迭代过程中每次的正演模拟波场值，包含水平分量和垂直分量，2l维；J₂为雅克比矩阵，2l行q列；Δm为每次的震源子波修正量，q维。

我们采用正则化算法，即求s₂＝||J₂Δm-Δv||²+r||LΔm||²极小。其中，s₂为迭代过程中每一步的极小化目标函数值；r为正则化系数。迭代过程中，解 $(J_1^T{J}_1+{\mathrm{\ϵ}}^2{L}^{T}L)\mathrm{\Δm}=J_1^T\mathrm{\Δv},$ m_k+1＝m_k+Δm。

通过上述步骤进行正则化反演优化之后，当满足预设收敛条件时输出优化的震源子波，即当残量Δv的模小于某个极小值时，跳出循环，输出反演得到的优化的震源子波。

具体实施时，本发明还提供了一种震源子波优化设置装置，如图2所示，该装置包括：

子波选取模块201，用于从地震资料中提取地震子波作为目标子波和初始震源子波；

初始震源子波更新模块202，用于利用目标子波和初始震源子波，采用正则化反演算法更新初始震源子波；

初始震源子波更新模块202具体用于，利用目标子波和初始震源子波，采用正则化反演算法更新用于二维声波有限差分数值模拟的初始震源子波；利用目标子波和初始震源子波，采用正则化反演算法更新用于弹性波有限差分数值模拟的初始震源子波。

优化震源子波输出模块203，用于当满足预设收敛条件时输出优化的震源子波。

该装置还包括：模拟波形获得模块204，用于根据得到的优化的震源子波，通过声波和弹性波的有限差分方程数值模拟得到波形；

验证模块205，用于对比模拟得到的波形和目标子波的波形，对优化的震源子波进行验证。

实施例1

以一个均匀声学介质为例来说明本发明的优势。速度为3000m/s，时间采样间隔为1ms，空间采用间隔为h＝10m，算子长度为M＝5。图3为30HZ的雷克子波波形图；将图3中的雷克子波作为初始震源子波，用于声波方程有限差分数值模拟。当用于一维声波方程有限差分数值模拟，获得的波形图如图4所示；当用于二维声波方程有限差分数值模拟时，获得的波形图如图5所示；当用于三维声波方程有限差分数值模拟时，获得的波形图如图6所示。由图3和图4、5、6对比可知，只有在三维情况下，加入的初始震源子波和模拟的波形的形态是相同的，而一维和二维情况下，加入的初始震源子波和模拟的波形的形态是相差较大。因此，实际中如果直接将野外提取的地震子波作为初始震源子波进行二维数值模拟、逆时偏移、全波形反演必然会降低精度。

采用本文提出的正则化反演算法来优化初始震源子波，将如图3所示的30HZ的雷克子波作为目标子波，通过正则化反演获得一维、二维和三维声波方程有限差分数值模拟中所需要加入的初始震源子波，如图7、8和9所示。想要得到与雷克子波形态相似的波形，在一维、二维和三维声波方程有限差分数值模拟中需要分别加入图7、8和9所示的初始震源子波。由此可知，在低维情况下，初始震源子波的优化设置是很有必要的。反演中，我们采用的初始震源子波为雷克子波(与真实模型差异较大)，仍然能够得到较好的结果，表明本发明中的震源子波优化算法是适用的。

实施例2

以一个均匀弹性介质为例来说明本发明的优势。纵波速度为3200m/s，横波速度为1600m/s，时间采样间隔为1ms，空间采用间隔为h＝10m，算子长度为M＝5。图10为20HZ的雷克子波波形图；将图10中的雷克子波作为初始震源子波，用于二维弹性波方程有限差分数值模拟，获得的波形图(包括水平分量和垂直分量)，如图11所示。由图10和11可知，加入的初始震源子波和模拟得到的波形形态差异较大。因此，弹性波情况下进行震源子波优化也是很有必要的。

我们仍然采用本文提出的正则化反演算法来优化震源子波，将如图10所示的20HZ的雷克子波作为目标子波，通过正则化反演获得弹性波方程有限差分数值模拟中所需要加入的初始震源子波，如图12所示。由图10和图12可知，反演得到的初始震源子波与雷克子波的形态相差甚远。用反演得到的初始震源子波再进行弹性波数值模拟得到的波形与目标子波的形态完全相同，进一步证明了本发明中针对弹性波方程震源优化算法的有效性。

综上所述，本发明具有以下优点：1.正则化算法可以大大提高反演过程的稳定性，得到的震源子波精度较高。2.计算量小，且容易实现。3.已知空间内某一点的波形即可，对波的相位没有要求(不需要零相位或最小相位的假设)。

以上所述仅为本发明的优选实施例而已，并不用于限制本发明，对于本领域的技术人员来说，本发明实施例可以有各种更改和变化。凡在本发明的精神和原则之内，所作的任何修改、等同替换、改进等，均应包含在本发明的保护范围之内。

Claims (10)

1.一种震源子波优化设置方法，其特征在于，包括：

从地震资料中提取地震子波作为目标子波和初始震源子波；

利用目标子波和初始震源子波，采用正则化反演算法更新初始震源子波；

当满足预设收敛条件时输出优化的震源子波。

2.如权利要求1所述的方法，其特征在于，所述利用目标子波和初始震源子波，采用正则化反演算法更新初始震源子波，包括：利用目标子波和初始震源子波，采用正则化反演算法更新用于二维声波有限差分数值模拟的初始震源子波；利用目标子波和初始震源子波，采用正则化反演算法更新用于弹性波有限差分数值模拟的初始震源子波。

3.如权利要求2所述的方法，其特征在于，所述利用目标子波和初始震源子波，采用正则化反演算法更新用于二维声波有限差分数值模拟的初始震源子波，按照如下公式计算：

s₁＝||J₁Δm-Δp||²+r||LΔm||²；

其中，

s₁为迭代过程中每一步的极小化目标函数值；

Δp为二维声波有限差分数值模拟值与实际观测值的残量，l维；

p^obs为实际观测到的波场值，l维；

p^cal为迭代过程中每次的正演模拟波场值，l维；

J₁为雅克比矩阵，l行q列；

Δm为每次的震源子波修正量，q维；

r为正则化系数。

4.如权利要求3所述的方法，其特征在于，所述收敛条件为Δp的模小于预设的极小值。

5.如权利要求2所述的方法，其特征在于，所述利用目标子波和初始震源子波，采用正则化反演算法更新用于弹性波有限差分数值模拟的初始震源子波，按照如下公式计算：

s₂＝||J₂Δm-Δv||²+r||LΔm||²；

其中，

s₂为迭代过程中每一步的极小化目标函数值；

Δv为弹性波有限差分数值模拟值与实际观测值的残量，包含水平分量和垂直分量，2l维；

v^obs为实际观测到的波场值，包含水平分量和垂直分量，2l维；

v^cal为迭代过程中每次的正演模拟波场值，包含水平分量和垂直分量，2l维；

J₂为雅克比矩阵，2l行q列；

Δm为每次的震源子波修正量，q维；

r为正则化系数。

6.如权利要求5所述的方法，其特征在于，所述收敛条件为Δv的模小于预设的极小值。

7.如权利要求1所述的方法，其特征在于，该方法还包括：

根据得到的优化的震源子波，通过声波和弹性波的有限差分方程数值模拟得到波形；

对比模拟得到的波形和目标子波的波形，对优化的震源子波进行验证。

8.一种震源子波优化设置装置，其特征在于，包括：

子波选取模块，用于从地震资料中提取地震子波作为目标子波和初始震源子波；

初始震源子波更新模块，用于利用目标子波和初始震源子波，采用正则化反演算法更新初始震源子波；

优化震源子波输出模块，用于当满足预设收敛条件时输出优化的震源子波。

9.如权利要求8所述的装置，其特征在于，所述初始震源子波更新模块具体用于，利用目标子波和初始震源子波，采用正则化反演算法更新用于二维声波有限差分数值模拟的初始震源子波；利用目标子波和初始震源子波，采用正则化反演算法更新用于弹性波有限差分数值模拟的初始震源子波。

10.如权利要求8所述的装置，其特征在于，该装置还包括：

模拟波形获得模块，用于根据得到的优化的震源子波，通过声波和弹性波的有限差分方程数值模拟得到波形；

验证模块，用于对比模拟得到的波形和目标子波的波形，对优化的震源子波进行验证。