基于文氏谱的起伏海面虚反射特征研究❋

高俊杰, 张 进,2,3❋❋, 彭阳阳, 董博艺

(1. 中国海洋大学海洋地球科学学院, 山东 青岛 266100;2. 青岛海洋科学与技术试点国家实验室 海洋矿产资源评价与探测技术功能实验室, 山东 青岛 266237;3.中国海洋大学 海底科学与探测技术教育部重点实验室, 山东 青岛 266100)

在海上地震勘探中,震源和检波器往往沉放在海面以下,这是由海洋地震勘探特殊施工条件所决定的,当震源激发后,作为强反射层的海面或海底,将会把地震波反射回来并形成虚反射。在海上作业时,虚反射的广泛发育严重影响地震资料的质量,需加以消除。常规的资料处理工作中,海面往往被假设为水平且反射系数简单定义为-1,实际上,海面呈现出起伏不平的复杂变化,导致实际地震记录中含有复杂的虚反射,常规的水平海面假设将会在反射系数和延迟时间的计算过程中引入误差,这就导致了虚反射算子估计不准确,影响地下储层的最终成像。

很多学者在模拟起伏海面方面做了大量的工作,以期得较为准确描述起伏海面反射系数。1983年Jovanovich等[1]提出了粗糙海面下的高斯曲面相干系数表达式,并用该相干系数来代替海面反射系数;1988年Macaskill[2]提出了粗糙海面下的平面波二维散射公式;2002年Laws和Kragh[3-4]分析时移地震成像观测结果后认为起伏海面造成的成像误差不可忽略;2005年Saenger[5]提出了随机海面下描述声波散射近似式;Robert等[6-7]随后在2006年又提出了用低频信息重构海面的方法;2015年Egorov等[8]提出了近似计算粗糙海面反射系数的积分表达式;2016年Zhang等[9]将海面反射系数看成是频率和波高的函数;2016年Park等[10]利用矩量法计算随风速变化的海面反射系数,并检验了模型准确性和有效区域;2017年Zhang等[11]通过高斯统计公式计算粗糙海面的反射系数实现粗糙海面虚反射压制;2021年Lv等[12]提出了一种基于环境参数和实时测量数据的浅海海面模型。上述方法重构的海面海浪形态不能完全符合海浪的统计学观测规律,与真实海面形态存在一定误差,从而影响反射系数和延迟时间的计算精度。

为了更真实地模拟起伏海面,本文引入海浪谱概念,基于文圣常等[13]提出并改进的文氏海浪谱(下文简称文氏谱),采用有效波高、海浪成长系数等各项有效海浪参量,建立符合中国海域海浪起伏特征的起伏海面模型;在满足差分条件下,使用有限差分交错网格方法进行正演模拟,通过时频分析方法描述起伏海面虚反射的振幅特征和频率特征。

1 起伏海面条件下地震虚反射

根据传播路径,虚反射可分为直达波虚反射、激发虚反射、接收虚反射、激发-接收虚反射[14]四种,如图1所示。

图1 一次反射与虚反射传播路径示意图

通常情况下,直达波虚反射在数据预处理过程中就被切除,所以实际地震记录关注三种虚反射,分别是激发虚反射、接收虚反射、激发-接收虚反射。

用s(t)表示实际地震记录,则:

s(t)=u(t)+Ru(t-Δts)+Ru(t-Δtr)+R2u(t-Δts-Δtr)。

(1)

式中:u(t)为一次反射波;R为反射系数;Δts、Δtr和(t-Δts-Δtr)分别表示激发、接收和激发接收虚反射延迟时间;u(t-Δts)和u(t-Δtr)为激发虚反射和接收虚反射。通过傅氏变换,式(1)为:

s(t)=U(f)(1+Re-i2πfΔts)(1+Re-i2πfΔtr)=U(f)·G(f)。

(2)

可以发现,含有虚反射的地震记录S在频率域相当于一次反射波U(f)和算子G(f)的乘积,其陷波周期主要与虚反射延迟时间Δt相关。

海面的反射系数是虚反射模拟以及虚反射压制处理中的一个重要参数,Jovanovich等[15]结合经典的散射理论并利用Kirchhoff近似推导出了高斯曲面的相关系数表达式,并提出该相关系数就是起伏海面下的海面反射系数,其表达式如下:

(3)

式中:f为入射波频率;φ为均方根振幅;θ是入射角;v为海水速度。根据公式(3),绘制出海面反射系数随入射角和频率的关系(见图2),当海面均方根振幅一定时,海面反射系数绝对值随入射角减小、地震波频率增大而减小[15]。在垂直入射条件下,当频率130 Hz时海面反射系数为-0.75,误差可达25%。因此,如果忽略海面起伏不平的复杂变化,普遍假设海水面平坦,将在虚反射反射系数和延迟时间计算中将引入巨大误差,故海面起伏形态与海面反射系数有关,且不可忽略。

图2 起伏海面反射系数随入射角和地震波频率的变化关系

2 基于文氏谱的起伏海面地震波场数值模拟

早在20世纪50年代初,海浪被认为是由许多振幅、频率、方向和相位不同的简单波动的叠加[16],因为这些简单波的振幅或相位被学者们规定为随机量,所以这些简单波叠加所得的结果也是随机量。其中由S(ω)表示角频率为ω的组成波的能量,则描述有限区域内海浪能量的函数S(ω)就被称为海浪谱。

海面起伏形态一般与水域深度、海面风速、地理位置和海浪状态等因素相关,常见的海浪谱主要有Neumann谱(诺依曼谱)、Pierson-Moskowitz谱(P-M谱),ITTC谱和双参数谱[17],而上述海浪谱主要来自北大西洋海域,与中国海域有所差异。为弥补上述海浪谱的缺陷,文圣常等结合中国黄、东、南、渤海资料,对理论风浪频谱进行改进和更广泛的检验,提出更符合中国海域海面起伏特征的海浪谱:

(4)

图3 不同风速下的文氏谱(a)和海面起伏模型(b)

为探究起伏海面条件下的虚反射特征,设置如图4所示的正演模型,其中海底深度为50 m,海水速度为1 500 m/s,地层速度为2 400 m/s,采用左侧单边放炮的激发方式,正演子波为30 Hz雷克子波,时间采样间隔为0.5 ms,空间采样间隔为1 m。根据文氏谱分别设置不同有效波高的海面进行正演模拟。

图4 起伏海面正演模型示意图

本文利用等距差分格式的有限差分方法,建立声波方程的层状速度模型进行数值模拟,并设置PML层,在兼顾精度与计算效率的同时减少边界效应的影响。

根据二维均匀介质中的弹性波方程可推导出声波方程[18]为:

(5)

式中:ρ为地层密度;u为应力分量;v为速度;vx、vz分别表示在x和z方向上的速度分量。

用交错网格有限差分方法(见图5)对式(5)中各项偏导数进行离散[19-21],如式(6):

图5 有限差分方程的交错网格示意图

(6)

式中:Δx、Δz分别表示 X 和 Z 方向上的空间网格步长;i、j分别表示相应的空间网格点;Δt表示时间网格步长;k表示相应的时间网格点;用U、V、W来分别表示应力、速度分量vx、vz的离散形式。

3 起伏海面的虚反射特征分析

本文基于文氏谱起伏海面建模和有限差分正演模拟方法,分别探究起伏海面虚反射对地震记录的改造作用、起伏海面条件下的虚反射时频谱特征。

3.1 起伏海面虚反射对地震记录的改造作用

为探究起伏海面虚反射对地震记录的改造作用,保持震源沉放深度为水下3 m,检波器沉放深度为水下10 m,分别设置海面有效波高H为0 m(水平海面)、2 m(低起伏海面)、3 m(中起伏海面)和4 m(高起伏海面)进行正演模拟。为符合实际海洋地震勘探现状,本次模拟中采用30 Hz雷克子波作为震源子波,时间采样间隔为1 ms。由此观测方式计算的激发、接收虚反射的延迟时间差值小于单位时间间隔,激发虚反射、接收虚反射交混在一起难以分离,故本文分析虚反射的整体响应。

由不同有效波高起伏海面模型的地震单炮正演记录(见图6)发现:

与水平海面模型的地震记录(见图6(a))相比,当海面为起伏形态时,记录中各同相轴粗糙程度明显增加,出现抖动现象。随着海面起伏增大,地震记录同相轴更加粗糙,能量分布不均匀程度增加,同相轴难以分辨(见图6(b)—(d))。

造成地震记录同相轴“抖动”的原因为:对于直达波,当炮检距较大时,直达波的入射和反射平行于海面,故与虚反射发生相对规则的干涉,海面起伏的影响较小;对于反射波,当海面为起伏状态时,相对直达波其入射方向于海面的夹角角度更大,故一次反射波和虚反射发生干涉的不规则程度相对增加,受海面起伏的影响更大。

为探究不同海面起伏程度对地下反射层成像效果的影响,将多个单炮记录进行动校正、叠加等处理,获得成像剖面(见图7)如下:

((a)水平海面;(b)低起伏海面(H=2 m);(c)中起伏海面(H=3 m);(d)高起伏海面(H=4 m)。 (a) Horizontal sea surface; (b) Slightly undulating sea surface (H=2 m); (c) Moderately undulating sea surface (H=3 m); (d) Heavily undulating sea surface (H=4 m).)

对比水平、低、中、高起伏海面条件地下反射层成像结果可见:相对于水平海面,起伏海面记录有着较为粗糙的同相轴,能量分布不均匀,且有可能在实际记录中出现虚假同相轴。此例中,当起伏海面H=2 m时,地层成像结果(见图7(b))与水平海面条件下的成像结果(见图7(a))基本一致;当起伏海面H=3 m时,0.55 s、1.17 s处出现了较为明显的虚反射同相轴(见图7(c)),在图7(d)中H=4 m时,起伏海面导致的虚反射更加严重,不仅强烈扭曲了海底一次波反射同相轴,且与一次波同相轴发生了强烈的干涉作用,一次波位置和能量发生变化(0.95 s),成像剖面出现假同相轴(0.8 s、1.25 s处)。

3.2 起伏海面条件下的虚反射时频谱特征

为进一步探究起伏海面条件下的虚反射时频谱特征,分别抽取上述四种不同波高起伏海面的单道虚反射记录(见图8(a),(c),(e),(g)),通过短时傅里叶变换方法获得时频谱(见图8(b),(d),(f),(h))。

((a)水平海面(H=0 m)单道记录;(b)水平海面时频谱;(c)低起伏海面(H=2 m)单道记录;(d)低起伏海面时频谱;(e)中起伏海面(H=3 m)单道记录;(f)中起伏海面时频谱;(g)高起伏海面(H=4 m)单道记录;(h)高起伏海面时频谱。(a) Single-trace record of horizontal sea level (H=0 m); (b) The spectrum of (a); (c) Single-trace record of slightly undulating sea surface (H=2 m); (d) The spectrum of (c); (e) Single-trace record of moderately undulating sea surface (H=3 m); (f) The spectrum of (e); (g) Single-trace record of heavily undulating sea surface (H=4 m); (h) The spectrum of (g).)

通过对比发现:与水平海面单道记录的时频分析结果相比,直达波虚反射能量强度随着海面有效波高的增大而相对减小。随着起伏海面有效波高的增加,有效波虚反射能量逐渐增强,与直达波虚反射响应时间差变化剧烈;多次波能量强度增大且分布紊乱程度加剧。之所以造成上述现象,除了受海绵起伏影响较大的一次反射波和虚反射的不规则干涉,还因为海面起伏造成的不均匀照明,加剧了反射能量的空间不均匀分布。

通过在模型四周设置PML层,可以得到不含虚反射的有效波记录。在同一个地层模型、相同的震检组合下,得到含有虚反射的水平海面正演结果(见图9(a)红)与不含虚反射的水平海面有效波正演结果(见图9(a)黑),并分别选取偏移距为30、35和40的单道虚反射记录进行频谱分析。同理,得到含有虚反射的有效波高2m的起伏海面正演结果(见图9(b)红)与不含虚反射的水平海面有效波正演结果(见图9(b)黑),并对不同偏移距的单道记录进行频谱分析,探究不同波高起伏海面条件下的虚反射信号频谱特征。结果如图9所示。

通过对比发现:不同海面的正演结果频谱中都存在陷波点,水平海面含虚反射记录频谱中陷波点分别在17、55 Hz的位置;起伏海面含虚反射记录频谱中陷波点分别在19、52 Hz的位置。相较于水平海面,起伏海面频谱陷波效应明显增强,低频和高频端陷波点处能量明显增强,主频幅值相对降低。随着偏移距的增加,水平、起伏海面条件下陷波点呈现相似规律,即主频幅值随着偏移距增加而逐渐减小,主频位置向高频方向移动;但起伏海面变化更为剧烈。

为探究不同观测系统条件下的虚反射时频谱特征,设置有效波高为3 m的起伏海面,保持震源沉放深度水下3 m不变,分别设置检波器沉放深度为水下4、6、8和10 m进行正演模拟。抽取单道虚反射记录(见图10(a),(c),(e),(g)),通过短时傅里叶变换方法获得时频谱(见图10(b),(d),(f),(h))。

((a)沉放4 m单道记录;(b)沉放4 m时频谱;(c)沉放6 m单道记录;(d)沉放6 m时频谱;(e)沉放8 m单道记录;(f)沉放8 m时频谱;(g)沉放10 m单道记录;(h)沉放10 m时频谱。(a) Single-trace record of geophone at depth of 4 m; (b) The spectrum of (a); (c) Single-trace record of geophone at depth of 6 m; (d) The spectrum of (c); (e) Single-trace record of geophone at depth of 8 m; (f) The spectrum of (e);(g) Single-trace record of geophone at depth of 10 m; (h) The spectrum of (g).)

由不同检波器深度单道记录的时频分析结果发现:与检波器沉放深度为4 m的虚反射单道记录的时频分析结果(见图10(b))相比,随着检波器沉放深度逐渐增加,深检波器沉放条件下直达波虚反射记录的能量强度逐渐减小,而有效波虚反射的能量强度逐渐增加。

在同一个地层模型、相同的偏移距下,正演模拟得到含有虚反射的水平海面正演结果(见图11(a)红)与不含虚反射的有效波正演结果(见图11(a)黑),分别对检波器沉放深度为4、6和8 m的单道虚反射记录进行频谱分析。同理,正演模拟得到含有虚反射的有效波高3 m的起伏海面正演结果(见图11(b)红)和不含虚反射的水平海面有效波正演结果(见图11(b)黑),并对不同沉放深度的单道记录进行频谱分析,探究不同观测系统条件下的虚反射频谱特征,结果如图11所示。

((a)水平海面含虚反射记录(红)、水平海面无虚反射记录(黑);(b)起伏海面含虚反射记录(红)、水平海面无虚反射记录(黑)。(a) under condition of horizontal sea surface with ghost wave (red), and without ghost wave (black); (b) under condition of undulating sea surface with ghost wave (red), and without ghost wave (black).)

通过对比发现:不同海面的正演结果频谱中都存在陷波点,水平海面含虚反射记录频谱中陷波点分别在13、64 Hz的位置;起伏海面含虚反射记录频谱中陷波点分别在17、64 Hz的位置。相较水平海面,起伏海面频谱的陷波效应明显增强,且频宽变窄。随着检波器沉放深度的增加,起伏海面条件下频谱受改造程度加剧,幅值波动明显增加。根据震源、检波器空间互易原理,激发虚反射与接收虚反射规律一致。

本文设置海浪状态为充分成长型海浪,该状态下由于其内部各组成波之间能量的非线性输送平衡,即使随着时间的变换,海浪能量的频率分布仍具有很好的相似性[22]。经过实验,不同时间、相同状态的海面正演记录所获结论一致。

4 结论

本文利用文氏谱,建立了更符合真实海面形态的起伏海面,并基于交错网格有限差分方法,正演模拟得到了不同有效波高起伏海面地震记录,通过短时傅里叶变换方法进一步研究了起伏海面条件下的虚反射特征,在研究过程中,得到了以下的结论和认识:

(1)当海面起伏时,虚反射常常与地下一次波发生干涉、叠加。当海面起伏程度较大时,海面虚反射对地震记录的改造作用强烈,不仅强烈扭曲了海底一次波反射同相轴,还导致地层成像结果中出现了很多假同相轴,严重影响了地震资料的信噪比和分辨率。

(2)当海面起伏时,虚反射和多次波能量强度增大、分布紊乱程度增加且无法忽略。地震记录受陷波效应逐渐增强,主频幅值逐渐减小,虚反射能量逐渐增加。

(3)随着偏移距的增加,地震记录陷波点幅值变化丰富,主频幅值逐渐减小,主频位置向高频方向移动。

(4)随着震源或检波器沉放深度逐渐增加,海面起伏对信号陷波点位置和能量的改造明显。起伏海面有效波虚反射的能量强度逐渐增加,陷波作用加剧,主频幅值逐渐减小。