基于Hilbert变换的复合波场弹性波逆时偏移成像方法

基于hilbert变换的复合波场弹性波逆时偏移成像方法
技术领域
1.本发明涉及基于hilbert变换的复合波场弹性波逆时偏移成像方法，属于石油地球物理勘探技术领域。


背景技术：

2.弹性波逆时偏移(ertm)算法基于弹性波动方程，能有效解决地下弹性介质的成像问题，利用转换横波波长更短的优势，提供比声波逆时偏移分辨率更高的地下介质成像。传统弹性波逆时偏移成像在成像结果中存在强烈的干扰噪音和假象。从行波分离的角度出发，可以剔除部分噪音及假象信息并保留有效成像信息，以达到高质量高分辨率的弹性波逆时偏移成像。基于波场分离的弹性波逆时偏移在充分利用弹性波在成像精度上优势的同时保证了对于复杂结构的成像质量。
3.弹性波逆时偏移技术作为对地下复杂构造成像的一大利器，自20世纪80年代被提出以来迅速发展。在是针对弹性波逆时偏移算法计算效率和成像精度方面不断优化。随着当代超算以及并行技术的突破，弹性波逆时偏移算法的计算效率得到了极大的提升，使其满足了对庞大的地震数据的储存以及计算的可行性，计算成本不再是限制弹性波逆时偏移算法的主要方面。这也就进一步提高了对弹性波逆时偏移算法可靠性和精确性的要求。传统弹性波逆时偏移成像在成像质量存在的问题，其一在偏移速度模型梯度较大的区域，出现回转波互相关，致使其成像点在真实界面的上方，出现假象；其二：当参与成像的两个波矢量夹角为超大角(接近180度)时会出现噪音。针对这两个问题，早期一些学者提出了采取行波分离的方式，将上下行波分离，只进行正传下行波与反传上行波的互相关成像，这样可以同时去除回转波假象和大的波矢量夹角产生的低频噪声。这样的处理在传统弹性波逆时偏移算法的基础上的确有着极大的优化，但随之也带来了两个急需解决的问题。
4.问题一，如果采用传统的频率波速域分离上下行波，则会产生较大的计算量和储存量，相比于此方法，利用fourier变换和hilbert变换的关系，对原始炮检波场依次进行时间方向和z方向的hilbert变换，再与原始炮检波场相加减即可得到分离后的上下行波。这种方法可以避免存储整个波场，提高了计算效率和资源利用率。问题二，对原始炮检波场分离后的行波，进行目的性的选取和剔除，只选择不含低频噪声和假象的正传下行波与反传上行波互相关成像，这样不仅去除了回转波假象和偏移低频噪声，也连带着丢失了弹性波逆时偏移算法对于地下高陡复杂结构的有效成像，这本该是弹性波逆时偏移成像的一大优势。
5.针对以上两个问题，如何去建立一种基于行波分离且保留有效回转波成像的弹性波逆时偏移算法，成为了之后我们主要讨论的问题。大致思路可以确定为：利用hilbert变换在行波分离上高计算效率和资源利用率这一优势，去建立一个复合的正传波场，其包括了完整的下行波和部分小角度的上行波，与完整的反传上行波进行互相关成像。这样既可以有效去除回转波假象又可以保留弹性波逆时偏移在高陡结构的棱柱波成像优势。同时相对于频率波数域的行波分离也大大提高了计算效率和资源利用率。
6.现有的基于hilbert变换的弹性波逆时偏移成像方法，均是对地震波场的一些特殊方向(如全波场分离)进行行波分离，以得到完整波场按照多个均匀角度的分割的分离波场。再一一对这些波场进行互相关成像，最后达到去噪去假象的效果。本发明在此基础之上，对重建后的波场数据网格进行一个方向的hilbert变换，再对波场数据网格进行还原以得到对特定方向进行hilbert变换后的波场数据，为组建复合正传p波，反传p波和反传s波场以及得到高精度的pp，ps弹性波逆时偏移成像提供了先决条件。


技术实现要素：

7.为了克服现有技术中的问题，本发明提供基于hilbert变换的复合波场弹性波逆时偏移成像方法，本发明利用hilbert变换实现非特殊角度的行波分离，并组建部分上行波+下行波的正传复合波场和部分下行波+上行波的反传复合波场，让其分别与反传上行波和正传下行波互相关，得到保留高陡结构的弹性波逆时偏移成像结果。
8.本发明解决上述技术问题所提供的技术方案是：基于hilbert变换的复合波场弹性波逆时偏移成像方法，包括以下步骤：
9.步骤s1、将原始炮检波场数据进行时间方向上的hilbert变换得到h
t
；
10.步骤s2、将h
t
沿着β方向组织数据并插值；
11.步骤s3、对插值后的波场进行β方向的hilbert变换后归位，即得到h
βht
项；
12.步骤s4、将h
t
数据水平方向180度旋转后，沿着β方向组织数据并插值，再进行β方向的hilbert变换后归位，即得到对称方向的h-βht
项；
13.步骤s5、将原始炮检波场数据直接沿着β方向组织数据并插值，再进行β方向的hilbert变换后归位，即得到h
β
；
14.步骤s6、将h
β
沿着β方向组织数据并插值，再进行β方向的hilbert变换后归位，即得到h
β
h-β
项；
15.步骤s7、利用任意β方向的全波场分离方法将上述各项带入斜方向的hilbert全波场分离公式即可得到原始炮检波场数据的四个波场分量；
16.步骤s8、通过正传或反传波场减去上波场分量或下波场分量得到正传复合波场和反传复合波场；
17.步骤s9、利用正传p波复合波场分别与反传p波，反传s波的上行波波场互相关，正传p波下行波波场与反传复合p波，反传s波复合波场分别进行互相关，得到四个通过复合波场来保留高陡构造的弹性波逆时偏移结果。
18.进一步的技术方案是，所述步骤s7中斜方向的hilbert全波场分离公式为：
19.sa(x,t)＝s-h
βht
(s)-h-βht
(s)-h
β
h-β
(s)
20.sb(x,t)＝s+h
βht
(s)-h-βht
(s)+h
β
h-β
(s)
21.sc(x,t)＝s-h
βht
(s)+h-βht
(s)+h
β
h-β
(s)
22.sd(x,t)＝s+h
βht
(s)+h-βht
(s)-h
β
h-β
(s),
23.ra(x,t)＝r-h
βht
(r)-h-βht
(r)-h
β
h-β
(r)
24.rb(x,t)＝r+h
βht
(r)-h-βht
(r)+h
β
h-β
(r)
25.rc(x,t)＝r-h
βht
(r)+h-βht
(r)+h
β
h-β
(r)
26.rd(x,t)＝r+h
βht
(r)+h-βht
(r)-h
β
h-β
(r).
27.式中：x为地下成像点位置；t为炮检点传播到成像点的时间；sa(x,t)、sb(x,t)、sc(x,t)、sd(x,t)分别为时空域震源波场的四个分量；ra(x,t)、rb(x,t)、rc(x,t)、rd(x,t)分别为时空域检波点波场的四个分量。
28.进一步的技术方案是，所述步骤s8中正传复合波场的计算公式为：
29.s
fd
(x,t)＝[3*s+h
βht
(s)-h-βht
(s)-h
β
h-β
(s)]/4
[0030]
式中：s
fd
(x,t)为正传复合波场。
[0031]
进一步的技术方案是，所述步骤s8中反传复合波场的计算公式为：
[0032]rfu
(x,t)＝[3*r-h
βht
(r)+h-βht
(r)-h
β
h-β
(r)]/4
[0033]
式中：r
fu
(x,t)为反传复合波场。
[0034]
进一步的技术方案是，所述步骤s9中pp成像结果的计算公式为：
[0035][0036]
式中：i
pp
为pp成像结果；x为地下成像点位置；t为炮检点传播到成像点的时间；t
max
为最大记录时间；src和rcv分别为震源波场和检波点波场；u和fd分别为反传上行波和正传复合下行波。
[0037]
进一步的技术方案是，所述步骤s9中ps成像结果的计算公式为：
[0038][0039]
式中：i
ps
为ps成像结果；x为地下成像点位置；t为炮检点传播到成像点的时间；t
max
为最大记录时间；src和rcv分别为震源波场和检波点波场；d和fu分别为正传下行波和反传复合上行波。
[0040]
进一步的技术方案是，所述步骤s3、s4、s5、s6中任意可变的β。
[0041]
本发明具有以下有益效果：本发明立足于弹性波逆时偏移成像，结合传统弹性波逆时偏移成像存在的一些缺点，从优化其计算效率以及成像结果精度的角度出发，在现有的弹性波逆时偏移成像技术的基础上，提出两个创新点：
[0042]
一、建立了基于hilbert变换的非特殊角度行波分离算法，解决了hilbert在行波分离中存在分离角度特殊单一的局限性，为后续组建成像信息更为丰富的复合波场提供了支撑，保留了其相对于频率-波数域分离方法更为高效的优势。使计算效率更高效的hilbert变换法更加灵活多用。
[0043]
二、本发明利用组建复合波场，在保证消除回转波假象以及削弱偏移低频噪声的同时，保留了弹性波逆时偏移对高陡结构的成像优势，并且对比与传统弹性波逆时偏移算法，其计算量仅增加了一倍左右；为弹性波逆时偏移算法在对地下储层复杂结构成像方面，提供了新的方法思路。
附图说明
[0044]
图1为复合波场快照图；
[0045]
图2为垂直方向行波分离数据输入示意图；
[0046]
图3为插值示意图；
[0047]
图4为特定斜方向行波分离数据输入示意图；
[0048]
图5为阶梯速度模型图；
[0049]
图6为阶梯模型t＝4s波场快照图；(a)传统正传下行波场；(b)复合正传下行波场；
[0050]
图7为传统逆时偏移结果，(a)pp；(b)ps；
[0051]
图8为震源下行波与检波点上行波互相关成像结果，(a)pp；(b)ps；
[0052]
图9为震源下行复合波场与检波点上行波互相关成像结果，(a)pp；(b)ps；
[0053]
图10为震源下行波与检波点复合上行波互相关成像结果，(a)pp；(b)ps。
具体实施方式
[0054]
下面将结合附图对本发明的技术方案进行清楚、完整地描述，显然，所描述的实施例是本发明一部分实施例，而不是全部的实施例。基于本发明中的实施例，本领域普通技术人员在没有做出创造性劳动前提下所获得的所有其他实施例，都属于本发明保护的范围。
[0055]
本发明的基于hilbert变换的复合波场弹性波逆时偏移成像方法，包括以下步骤：
[0056]
步骤s1、将原始炮检波场数据进行时间方向上的hilbert变换得到h
t
；
[0057]
步骤s2、将h
t
沿着β方向组织数据(如图2所示)并插值(如图3所示)；
[0058]
步骤s3、对插值后的波场进行β方向的hilbert变换后(如图4所示)归位，即得到h
βht
项；
[0059]
步骤s4、将h
t
数据水平方向180度旋转后，沿着β方向组织数据并插值，再进行β方向的hilbert变换后归位，即得到对称方向的h-βht
项；
[0060]
步骤s5、将原始炮检波场数据直接沿着β方向组织数据并插值，再进行β方向的hilbert变换后归位，即得到h
β
；
[0061]
步骤s6、将h
β
水平方向180度旋转后，继续沿着β方向组织数据并插值，再进行β方向的hilbert变换后归位，即得到h
β
h-β
项；
[0062]
步骤s7、利用任意β方向的全波场分离方法将上述各项带入斜方向的hilbert全波场分离公式即可得到原始炮检波场数据的四个波场分量；
[0063]
sa(x,t)＝s-h
βht
(s)-h-βht
(s)-h
β
h-β
(s)
[0064]
sb(x,t)＝s+h
βht
(s)-h-βht
(s)+h
β
h-β
(s)
[0065]
sc(x,t)＝s-h
βht
(s)+h-βht
(s)+h
β
h-β
(s)
[0066]
sd(x,t)＝s+h
βht
(s)+h-βht
(s)-h
β
h-β
(s)，
[0067]
ra(x,t)＝r-h
βht
(r)-h-βht
(r)-h
β
h-β
(r)
[0068]
rb(x,t)＝r+h
βht
(r)-h-βht
(r)+h
β
h-β
(r)
[0069]
rc(x,t)＝r-h
βht
(r)+h-βht
(r)+h
β
h-β
(r)
[0070]
rd(x,t)＝r+h
βht
(r)+h-βht
(r)-h
β
h-β
(r).
[0071]
步骤s8、通过正传或反传波场减去上波场分量或下波场分量得到正传复合波场和反传复合波场；
[0072]sfd
(x,t)＝[3*s+h
βht
(s)-h-βht
(s)-h
β
h-β
(s)]/4
[0073]rfu
(x,t)＝[3*r-h
βht
(r)+h-βht
(r)-h
β
h-β
(r)]/4
[0074]
步骤s9、利用正传p波复合波场分别与反传p波，反传s波的上行波波场互相关，正传p波下行波波场与反传复合p波，反传s波复合波场分别进行互相关，得到四个通过复合波场来保留高陡构造的弹性波逆时偏移结果。
[0075][0076][0077]
式中：i
pp
和i
ps
分别表示pp和ps成像结果；x为地下成像点位置；t为炮检点传播到成像点的时间；t
max
为最大记录时间；src和rcv分别表示震源波场和检波点波场；u、d和fu、fd分别表示上行波、下行波和复合上行波、复合下行波。
[0078]
上述β是可变的，在选取更大β的情况下可以引入更丰富的成像信息，但同时也会导致更多的低频噪声。通过实验，β＝arctan1/3可以达到保留高陡结构和去噪的最佳折衷效果。
[0079]
实施例
[0080]
为验证该方法的有效性，选择如图5所示的平滑后模型进行偏移，采用主频为25hz的雷克子波作为震源子波，大小为10000m
×
3000m的阶梯速度模型，垂直构造的位置在水平7km。背景p波，s波与密度数值如下：
[0081]
图6所展示的为t＝4s时刻的正传复合下行波场和传统下行波波场的波场快照；
[0082]
如图6中，在黑色和白色箭头所示位置，复合波场保留着部分小角度的上行波场，丰富了波场信息的同时滤除了大角度的上行波场。通过分析反射路径可以知道，黑色箭头所指位置部分上行波是由于在高陡构造发生二次反射所产生，故推测复合波场互相关对于高陡构造成像是具有实质性效果的。
[0083]
对图5所示的模型进行平滑后偏移成像；
[0084]
如图7所示，传统弹性波逆时偏移成像对于高陡倾斜构造成像效果明显，但在速度剧烈变化区域会出现回转波假象。众所周知，该假象是由传统弹性波逆时偏移中震源上行波和检波点下行波产生；图8所示为dpertm方法的成像结果，该方法即通过hilbert变化实现上下行波分离，从而只保留震源下行波与检波点上行波的互相关成像，从成像结果可以看出回转波假象和对高陡构造的有效成像同时被去除。
[0085]
图9、10分别为震源下行复合波场与检波点上行波互相关和震源下行波与检波点复合上行波互相关成像结果，由于该成像条件中去除震源上行和检波点下行互相关项的同时保留了部分小角度的同方向波场互相关信息，可以在成像结果中看出在去除假象的同时，对于高陡构造的成像信息也被保留了下来。
[0086]
从波场快照以及成像结果的对比中，本发明方法的有效性得到了验证，及通过基于hilbert的全波场分离思想来组建复合波场，以保证在去除假象的同时保留弹性波逆时偏移算法对于高陡构造成像的有效性。
[0087]
本发明的优点，第一：对于dpertm方法，本发明在同样去除低频噪声假象的同时，也保留了部分正传上行波与反传上行波和正传下行波和反传下行波的成像信息，增强了对高陡结构的成像效果；第二：对于fpdprtm方法，本发明使用任意β方向的hilbert变换以及全波场分离思想结合，使得波场分离和复合波场的组建同时实现，避免了后续更繁琐的空间域滤波计算，极大的提高了计算效率和资源利用率。并且将构建复合波场的思想引入弹性波逆时偏移成像当中，能够提取出更多的有效成像信息，进一步提升了弹性波逆时偏移对于复杂构造的成像优势。
[0088]
以上所述，并非对本发明作任何形式上的限制，虽然本发明已通过上述实施例揭
示，然而并非用以限定本发明，任何熟悉本专业的技术人员，在不脱离本发明技术方案范围内，可利用上述揭示的技术内容作出些变动或修饰为等同变化的等效实施例，但凡是未脱离本发明技术方案的内容，依据本发明的技术实质对以上实施例所作的任何简单修改、等同变化与修饰，均仍属于本发明技术方案的范围内。