CN106597535A - 一种提高弹性波逆时偏移计算率和空间分辨率的方法 - Google Patents

Abstract

本发明公开了一种提高弹性波逆时偏移计算率和空间分辨率的方法。考虑在震源波场正向传播过程中，计算每个网格点的能量，保存最大能量密度的时刻和相应的波场值；在检波器波场逆时传播过程中，利用成像条件在每个网格点提取最大能量的检波器波场，用保存的最大能量震源波场做归一化，获得角度依赖的反射系数成像剖面。相比于激发时间成像条件，本发明自动校正了水平分量在震源两侧的极性反转，在多炮叠加时避免振幅损失；相比于归一化互相关成像条件，本发明在震源波场正向传播过程中无需存储波场快照，节省大量磁盘空间和I/O吞吐任务，提高了计算效率，产生的低频假象更弱，水平分量的成像能力更优，具有更高的空间分辨率。

Description

一种提高弹性波逆时偏移计算率和空间分辨率的方法

技术领域

本发明涉及地震勘探领域，特别是一种提高弹性波逆时偏移计算率和空间分辨率的方法。

背景技术

在勘探地震学领域，地震偏移成像技术是获取地下构造形态有效的途径，其中，叠前深度域偏移是最为精确的成像手段之一。叠前逆时偏移包括波场外推和成像条件的运用，其中成像条件直接关系到成像的质量和计算效率。目前比较常用的激发时间成像条件，无需存储震源波场，但仅仅利用了初至或者最大振幅波场，且成像剖面在量纲上不是反射系数，图像不够清晰；互相关成像条件需要存储正传的波场快照。尽管互相关成像条件具有实现简单、稳定等特点，但是成像剖面低频假象严重，同时需要硬盘存储震源的波场快照，在利用成像条件时又需要读取震源波场的波场快照，严重降低了逆时偏移的效率。

发明内容

本发明的目的是克服现有技术的不足，提供一种利用修正的成像条件提高弹性波逆时偏移计算率和空间分辨率的方法。

具体步骤如下：

(1)采用一阶速度-应力的二维弹性波动方程模拟地震波场

在波场模拟时，使用交错网格有限差分，把一阶速度-应力的二维弹性波动方程离散成空间为高阶(采用8阶)，时间为二阶精度的交错网格差分,人工边界采用PML吸收边界条件^]来压制来自边界的反射波，在计算上述的能量密度时，需要在时空上进行插值；为了提高计算效率，采用线性插值；为了高效地去除逆时偏移的成像剖面中直达波和震源这部分低频的假象，采用空间高通滤波方法。

(2)成像条件的修正

考虑到归一化互相关成像剖面低频假象严重，同时需要硬盘存储震源的波场快照，在利用成像条件时又需要读取震源波场的波场快照，严重降低了逆时偏移的效率；激发振幅成像条件无需任何硬盘存储，从而避免大量的硬盘读写操作，大大提高了计算效率；相比于激发时间成像条件，由于该成像条件用最大能量震源波场做归一化，自动校正了水平分量在震源两侧的极性反转，在多炮叠加时避免振幅损失，并且得到角度依赖的反射系数剖面；但是，该成像条件目前仅用在声波逆时偏移中，将其直接用在弹性波中还存在一定的不稳定性；为此，提出修正激发振幅成像条件，使逆时偏移更加稳定健壮。

归一化互相关成像条件主要通过计算震源波场和接收波场的零延迟互相关，并进行归一化处理：

其中，S,R分别表示震源波场与检波器波场；分别表示对炮集和时间进行叠加。

激发振幅成像条件通过在震源波场的正传过程中，在每个时刻计算出网格点的能量密度，并保存下相应的最大能量密度对应的走时和波场值；在检波器波场的逆传过程中，利用保存的走时提取出相应的检波器波场的波场值，并用保存的震源波场做归一化，从而得到角度依赖的反射系数剖面，激发振幅成像条件表达形式为：

S,R分别表示最大能量密度对应的震源波场和提取出的检波器波场；表示对炮集进行叠加；不同于声波波场，弹性波波场可能在某些网格点上始终为0，如水平速度分量剖面上通过震源的垂直对称轴；(1)和(2)式则能被一个很小的数除，表现数值不稳定。

为此,我们对其进行修正，并得到如下稳定的激发振幅成像条件：

其中，表示对最大能量密度对应的震源波场值求取绝对值后再取的平均值，N为二维网格的节点总数。sign[S(x,z)]表示对S(x,z)取符号。相比于激发时间成像条件，由于该成像条件用最大能量震源波场做归一化，自动校正了水平分量在震源两侧的极性反转，且在多炮叠加时避免由于极性反转引起的振幅损失。

(3)实例验证

基于典型的倾斜模型和断层模型利用交错网格有限差分模拟地震波场，继而通过归一化互相关成像条件和激发振幅成像条件分别进行逆时偏移，分析两种成像条件下的模型重构图形，结果充分验证了修正的激发振幅成像条件极大地提高了逆时偏移的计算率和空间分辨率。

(4)计算效率对比

为了对比上述两种成像条件的计算效率，分别统计上述两个模型的计算时间和对硬盘的需求量。

本发明充分利用了空间为高阶、时间为二阶精度的交错网格差分，人工边界采用PML吸收边界条件来压制来自边界的反射波以及采用空间高通滤波进行去除逆时偏移的成像剖面中直达波和震源低频的假象，并利用修正的激发振幅成像条件提高了弹性波逆时偏移计算率和空间分辨率，为解决弹性波逆时偏移计算量和空间分辨率低的问题提供了一个有效的实施方案。

附图说明

图1为本发明实施例倾斜模型示意图。

图2为本发明实施例第24炮的P波速度模型(a)，及最大能量走时等值线(b)，x分量的震源波场值(c)，z分量的震源波场值(d)；虚线表示界面的位置。

图3为本发明实施例两种不同成像条件下的第24炮成像剖面，虚线表示界面位置；归一化互相关成像剖面：(a)水平分量；(b)垂直分量；稳定激发振幅成像剖面：(c)水平分量；(d)垂直分量。

图4为本发明实施例偏移剖面的空间高通滤波结果。

图5为本发明实施例两种不同成像条件下的46炮叠加成像剖面；归一化互相关成像剖面：(a)水平分量；(b)垂直分量；稳定激发振幅成像剖面：(c)水平分量；(d)垂直分量。

图6为本发明实施例断层模型。

图7为本发明实施例两种不同成像条件下的81炮叠加成像剖面(未滤波)；归一化互相关成像剖面：(a)水平分量；(b)垂直分量；稳定激发振幅成像剖面：(c)水平分量；(d)垂直分量。

具体实施方式

实施例：

利用修正的成像条件提高弹性波逆时偏移计算率和空间分辨率的方法步骤如下：

(1)采用一阶速度-应力的二维弹性波动方程模拟地震波场

在波场模拟时，使用交错网格有限差分，把一阶速度-应力的二维弹性波动方程离散成空间为高阶(采用8阶)，时间为二阶精度的交错网格差分,人工边界采用PML吸收边界条件^]来压制来自边界的反射波，在计算上述的能量密度时，需要在时空上进行插值。为了提高计算效率，采用线性插值。为了高效地去除逆时偏移的成像剖面中直达波和震源这部分低频的假象，采用空间高通滤波方法。

(2)成像条件的修正

考虑到归一化互相关成像剖面低频假象严重，同时需要硬盘存储震源的波场快照，在利用成像条件时又需要读取震源波场的波场快照，严重降低了逆时偏移的效率。激发振幅成像条件无需任何硬盘存储，从而避免大量的硬盘读写操作，大大提高了计算效率；相比于激发时间成像条件，由于该成像条件用最大能量震源波场做归一化，自动校正了水平分量在震源两侧的极性反转，在多炮叠加时避免振幅损失，并且得到角度依赖的反射系数剖面。但是，该成像条件目前仅用在声波逆时偏移中，将其直接用在弹性波中还存在一定的不稳定性。为此，提出修正激发振幅成像条件，使逆时偏移更加稳定健壮。

归一化互相关成像条件主要通过计算震源波场和接收波场的零延迟互相关，并进行归一化处理：

其中，S,R分别表示震源波场与检波器波场；分别表示对炮集和时间进行叠加。

激发振幅成像条件通过在震源波场的正传过程中，在每个时刻计算出网格点的能量密度，并保存下相应的最大能量密度对应的走时和波场值；在检波器波场的逆传过程中，利用保存的走时提取出相应的检波器波场的波场值，并用保存的震源波场做归一化，从而得到角度依赖的反射系数剖面，激发振幅成像条件表达形式为：

S,R分别表示最大能量密度对应的震源波场和提取出的检波器波场；表示对炮集进行叠加。不同于声波波场，弹性波波场可能在某些网格点上始终为0，如水平速度分量剖面上通过震源的垂直对称轴。(1)和(2)式则可能被一个很小的数除，表现数值不稳定。

为此,我们对其进行修正，并得到如下稳定的激发振幅成像条件：

其中，表示对最大能量密度对应的震源波场值求取绝对值后再取的平均值，N为二维网格的节点总数。sign[S(x,z)]表示对S(x,z)取符号。相比于激发时间成像条件，由于该成像条件用最大能量震源波场做归一化，自动校正了水平分量在震源两侧的极性反转，且在多炮叠加时避免由于极性反转引起的振幅损失。

(3)实例验证

为了验证修正的激发振幅成像条件在弹性波逆时偏移中的有效性，分别在倾斜模型和断层模型上进行了数值试验。其中倾斜模型(图1)的网格数为1201×601，网格间距10m，上层的倾角是9°；各层的密度及其纵横波速度详见图1。检波器间距为10m；采用双边接收观测系统，两侧各布置150道检波器；震源布置在地表，炮间距为200m，共46炮，布置范围为1500m～10500m；震源为10Hz的雷克子波的爆炸震源；波场模拟时间间隔为1ms，传播至6.0s。

图2为第24炮，采用(3)式计算得到的走时和相应的震源波场值。其中，图2(a)为第24炮的偏移孔径内的速度模型，图2(b)为计算的走时；图中虚线表示两个界面的位置。

图3为对第24炮进行空间高通滤波之前的水平(3a,3c)和垂直(3b,3d)分量的成像剖面。其中，图3(a)和3(b)为归一化互相关成像条件的成像剖面；图3(c)和3(d)为稳定激发振幅成像条件的成像剖面。由于长时程计算，来自边界的虚假反射、以及层间反射波等导致部分虚假反射同相轴(图3中箭头标注的位置)。从图中可以看到，归一化互相关成像条件的成像剖面存在严重的低频假象，如图中椭圆标注范围，垂直分量尤为明显；而激发振幅成像剖面的低频假象比前者弱。从分辨率上看，激发振幅成像条件的分辨率高于归一化互相关成像条件的，尤其表现在垂直分量上；从成像能力方面来看，激发振幅成像条件的水平分量优于互相关成像条件的；但是两种成像条件产生的垂直分量剖面都存在严重的低频假象，这与直达波的切除和边界条件有着密切的关系。

图4是对图3进行空间高通滤波后的结果。从图中可见，空间高通滤波能高效的去除成像剖面上的低频假象。稳定激发振幅成像条件的水平分量的假象比归一化互相关成像条件的更为严重(如箭头所示)，但是前者的分辨率高于后者；并且两者成像条件的垂直分量更获得更加清晰地图像。

图5为倾斜模型的46炮的叠加结果。从图中可见,两种成像条件的垂直分量剖面的界面非常清晰，由于两端叠加次数仅为1，两端的图像不如中间部分(多次叠加)清楚。图5(c)中的第二个界面比图5(a)更加清晰(上述4幅图使用了相同的增益)。

为了验证该成像条件在复杂模型中的有效性，我们在断层模型(图6)上进行了数值试验。该断层模型的网格为1201×601，网格间距为10m；检波器间距为10m；采用双边接收观测系统，两侧各布置200道检波器；震源布置在地表，炮间距为100m，共81炮，布置范围为2000m～10000m；震源为10Hz的雷克子波，为爆炸震源；波场模拟时间间隔为1ms；波场传播计算到6.0s。纵、横波速度和密度见图6所示。

图7为断层模型的81炮叠加结果(未滤波)。从图中可见，稳定激发振幅成像条件也能够获得与归一化互相关成像条件相当的图像，尤其是断层模型的底部的两个界面。总体上来讲，归一化成像条件得到的偏移剖面低频假象更为严重，垂直分量能够获得比水平分量更加清晰的图像。

(4)计算效率对比

为了对比上述两种成像条件的计算效率，分别统计了上述两个模型的计算时间和对硬盘的需求量(表1、2)。倾斜模型和断层模型分别在Lenovo ThinkCentre M8300t和DellPrecision T5600 WorkStation(linux系统)上运行。程序均为有限差分串行程序。

表1 CPU耗时统计

表2硬盘需求情况统计

表1和表2显示：由于归一化成像条件需要存储波场快照，需要巨大的硬盘空间和I/O吞吐操作；而稳定激发振幅成像条件仅需在每个时间步计算能量密度，所需要的额外计算量很小。

综上来看，与归一化互相关成像条件相比，修正的激发振幅成像条件具有计算量小的特点，逆时偏移更稳定健壮，同时水平分量和垂直分量偏移剖面的分辨率均较高，低频假象弱，硬盘需求量为零，避免了I/O吞吐开销，提高了计算效率和空间分辨率。

Claims (1)

1.一种提高弹性波逆时偏移计算率和空间分辨率的方法，其特征在于具体步骤为：

(1)采用一阶速度-应力的二维弹性波动方程模拟地震波场

在波场模拟时，使用交错网格有限差分，把一阶速度-应力的二维弹性波动方程离散成空间为高阶,采用8阶，时间为二阶精度的交错网格差分,人工边界采用PML吸收边界条件^]来压制来自边界的反射波，在计算上述的能量密度时，需要在时空上进行插值；为了提高计算效率，采用线性插值；为了高效地去除逆时偏移的成像剖面中直达波和震源这部分低频的假象，采用空间高通滤波方法；

(2)成像条件的修正

考虑到归一化互相关成像剖面低频假象严重，同时需要硬盘存储震源的波场快照，在利用成像条件时又需要读取震源波场的波场快照，严重降低了逆时偏移的效率；激发振幅成像条件无需任何硬盘存储，从而避免大量的硬盘读写操作，大大提高了计算效率；相比于激发时间成像条件，由于该成像条件用最大能量震源波场做归一化，自动校正了水平分量在震源两侧的极性反转，在多炮叠加时避免振幅损失，并且得到角度依赖的反射系数剖面；但是，该成像条件目前仅用在声波逆时偏移中，将其直接用在弹性波中还存在一定的不稳定性；为此，提出修正激发振幅成像条件，使逆时偏移更加稳定健壮；

归一化互相关成像条件主要通过计算震源波场和接收波场的零延迟互相关，并进行归一化处理：

$I(x,z)=\underset{shot}{\Σ}\frac{\underset{t}{\Σ}S(x,z,t)R(x,z,t)}{\underset{t}{\Σ}S(x,z,t)S(x,z,t)}---\left(1\right)$

其中，S,R分别表示震源波场与检波器波场；分别表示对炮集和时间进行叠加；

激发振幅成像条件通过在震源波场的正传过程中，在每个时刻计算出网格点的能量密度，并保存下相应的最大能量密度对应的走时和波场值；在检波器波场的逆传过程中，利用保存的走时提取出相应的检波器波场的波场值，并用保存的震源波场做归一化，从而得到角度依赖的反射系数剖面，激发振幅成像条件表达形式为：

$I(x,z)=\underset{shot}{\Σ}\frac{R(x,z)}{S(x,z)}---\left(2\right)$

S,R分别表示最大能量密度对应的震源波场和提取出的检波器波场；表示对炮集进行叠加；不同于声波波场，弹性波波场可能在某些网格点上始终为0，如水平速度分量剖面上通过震源的垂直对称轴；(1)和(2)式则能被一个很小的数除，表现数值不稳定；

为此,我们对其进行修正，并得到如下稳定的激发振幅成像条件：

其中，表示对最大能量密度对应的震源波场值求取绝对值后再取的平均值，N为二维网格的节点总数；sign[S(x,z)]表示对S(x,z)取符号；相比于激发时间成像条件，由于该成像条件用最大能量震源波场做归一化，自动校正了水平分量在震源两侧的极性反转，且在多炮叠加时避免由于极性反转引起的振幅损失；

(3)实例验证

基于典型的倾斜模型和断层模型利用交错网格有限差分模拟地震波场，继而通过归一化互相关成像条件和激发振幅成像条件分别进行逆时偏移，分析两种成像条件下的模型重构图形，结果充分验证了修正的激发振幅成像条件极大地提高了逆时偏移的计算率和空间分辨率；

(4)计算效率对比

为了对比上述两种成像条件的计算效率，分别统计上述两个模型的计算时间和对硬盘的需求量。