0 引言

可控源海洋电磁(MCSEM)勘探实现油气储层识别和描述的主要原理在于电磁能量在海底地层中传播时遇到不同电阻率层表现出衰减和相位的差异^{[1, 2, 3, 4]}。油气储层电阻率比围岩高几到几十倍，电磁能量在储层中衰减较慢，因此，储层中的电磁场幅度比围岩中大^{[5, 6]}。然而，在浅水域实施MCSEM勘探时，发射源在空气-海水界面激发的空气波在大偏移距的海底接收器响应信号中占主导，使来自海底地层的有效信号淹没于强空气波信号中，难以区分高阻层地层引起的异常^{[7, 8]}。相应地，在深水域，由于空气波在海水层中的双程衰减，对接收信号的影响可以忽略。因此，MCSEM方法在海洋油气勘探中的成功应用目前主要限于深水域(水深>300 m)^{[9, 10]}。

图 1给出了海洋电磁信号的基本传播路径和模式，其中的“空气波”是受弹性波传播的启发，最初是由Chave和Cox^[11]命名的，主要指由源到海水面垂直扩散的信号产生、在空气层以波的形式传播、通过海水层扩散返回海底的信号，它在空气层传播时几乎没有衰减。与深水域相比，在浅水域时接收器记录的空气波要强得多，这是浅水域海洋电磁勘探的困难所在。

图 1 水平电偶极子源发射低频电磁波在海水表面、水中和海底地层传播示意图 Fig. 1 Schematic diagram of propagation paths of low-frequency electromagnetic wave emitted from a horizontal electric dipole to the sea surface,seawater and undersea stratum

图选项

关于“空气波”性质的研究，早在50年前Wait^[12]就曾论述过，后来Baños^[13]、Bannister^[14]和King等^[15]也已开展过研究，只是他们没有用“空气波”的名称而已，他们将之称为“侧面波”并给出了半空间水层之上空气层中的空气波表达式。近年来，Constable等^[1]、MacGregor等^[4]、Løseth等^[16]基于层状介质的TE(横电)和TM(横磁)模式分解理论，分析了空气波产生的机制和影响大小。Um等^[5]用电流密度剖面分析了不同海水深度上空气波能量占主导的临界偏移距位置，得到了垂直电场分量不受空气波影响的重要认识。Ziolkowski等^[17]提出了时间域的信号分离等多种压制空气波方法。Amundsen等^[8]基于“空气波”的波动性质，将接收的电磁信号分离为上行和下行的电磁波场，通过去除下行波场实现空气波的压制。Lu等^[18]提出了应用数值模拟的方法估计空气波信号，再从接收信号中减去的方法，减弱空气波的影响。国内，殷长春等^[19]、沈金松等^[20]和刘长胜等^[21]基于一维层状介质的电磁场响应，分析了空气波随水深变化对接收器响应的影响和相应的压制方法。

本研究基于电磁场的模式分解理论导出海底接收信号中空气波与空气层、海水层、海底地层的相互作用关系，用空气波产生的物理机制解释空气波与海底地层信号的复杂耦合特征。在Nordskag等^[22]研究的基础上，利用有限水深模型，建立空气波在海底地层和海水界面之间交混回响衰减模式，导出有限水深时空气波分量的近似表达式，构建模拟电磁场相减压制空气波的方法。同时，还依据Amundsen等^[8]的波场分离方法，推导海洋电磁测量数据的上行和下行波场分离关系，讨论对剖面数据实现空气波压制的效果。通过层状背景模型和高阻含油气储层模型的模拟数据比较2种空气波压制方法的效果，分析各自的优缺点。 1 电磁场基本方程

由频率域Maxwell方程，低频水平电偶极子(HED)的电磁场可表示为^[23]

点击浏览公式

式中：E为电场强度(V/m)；H为磁场强度(A/m)；J^s_e是水平电偶极子源的电流密度，J^s_e=px，p和x分别是水平电偶极子的偶极矩和x方向的单位矢量；磁阻抗=iωμ；导纳=σ+iωε。假设介质磁导率μ=μ₀=4π×10^-7H/m，介质电导率为σ(S/m)，介质的介电常数ε=ε_rε₀，相对介电常数为ε_r，真空介电常数ε₀=1/(36π)×10^-9F/m。

在层状均匀介质中，方程(1)的电场和磁场可分解为横电和横磁2种模式及相应的矢量位函数表示^{[20, 23]}：

点击浏览公式

点击浏览公式

点击浏览公式

其中：E^TE、E^TM和H^TE、H^TM分别表示横电(TE)和横磁(TM)模式的电场和磁场强度；F和A为2种模式对应的位函数。F和A满足如下非均匀Helmholtz方程^[23]：

点击浏览公式

式中，k²=-。利用式(3)和式(4)的电磁场与电磁位的关系，可得海水层中电磁场表达式。

对于图 2所示的无限水深单界面的情况，若以海水-空气界面为z=0，向下为正，利用电磁场连续性条件，即可求得海水层中HED源激发的空气和海水层中的电磁位^[20]：

点击浏览公式

点击浏览公式

点击浏览公式

点击浏览公式

其中：A₀^TM(ρ,z)、A₁^TM(ρ,z)、F₀^TE(ρ,z)和F₁^TE(ρ,z)分别表示在空气层和海水层中的TM和TE模式的电磁位函数；ρ=(x²+y²)^1/2为径向距离；h为源到空气-海水界面的垂直距离；λ²=k_x²+k_y²，λ为柱坐标下的空间水平方向波数，k_x和k_y分别是直角坐标系下二维傅里叶变换x、y方向的波数；J₀为第一类零阶贝塞尔函数； u_i=(k_x²+k_y²-k_i²)^1/2，i=0，1分别是与空气和海水层参数相关的波数，k_i²=-_ii，i=0，1时k_i分别为空气和海水层的波数；₀=₁=iωμ₀为空气和海水层的磁阻抗；_i=σ_i+iε_iω，i=0，1分别为空气和海水层的导纳；R_0u^TM、R_1d^TM、R_0u^TE和R_1d^TE分别为空气和海水交界面上TM和TE模式电磁场的上行和下行电磁波场的反射和透射系数，且R_0u^TM=，R_1d^TM=，R_0u^TE=，R_1d^TE=，Z_i=，i=0，1分别与空气及海水层的电磁参数有关。 2 海水半空间模型中空气波的渐近表达式

为更好地分析源信号与空气界面、海水层的相互作用，以半空间海水层模型(图 2)的电场径向分量为例对空气波的表达式进行分析。根据式(6)－式(9)，TM和TE两种模式的径向电场分量可表示如下^[20]：

点击浏览公式

点击浏览公式

式中：P^TM(u₁,z)=e^-u₁|z-h|+R_1d^TMe^-u₁(z+h)；P^TE(u₁,z)=e^-u₁|z-h|+R_1d^TEe^-u₁(z+h)；J₁为第一类一阶贝塞尔函数。

图 2 无限水深MCSEM分析介质模型 Fig. 2 Model of infinite water depth for MCSEM analysis

图选项

实际勘探中，源激发的TM和TE两种模式的电磁信号与空气-海水界面的相互作用全部体现在与R_1d^TM和R_1d^TE2个反射系数相关的项中，它与海水深度、空气与水的电导率及激发频率有关。对于TM模式，由式(10)的反射系数R_1d^TM==可知，在MCSEM的激发频率下，_i≈σ_i(i=0,1)，空气电导率σ₀→0，则R_1d^TM=-1。即海水层中的TM模式电磁信号在空气界面只发生相位翻转，返回海水层。因此，TM模式信号不受空气波影响。然而，对于空气和海水层(μ_i=μ₀(i=0,1))中的TE模式，式(11)中的R_1d^TE==≈，当空气电导率σ₀→0时，u₀=(k_x²+k_y²-k₀²)^1/2≈(k_x²+k_y²)^1/2=λ；因此，无论海水电导率如何变化，在空间波数λ不为零时，R_1d^TE都不为1。即通常情况下，海水层中的TE模式电磁信号均在海水-空气界面发生反射，而且，空气-海水界面的反射系数随海水电导率和横向空间波数发生变化，这是大偏移距上空气波影响大的原因之一。由式(11)知道，TE模式中只有反射项与空气波相关。则与空气波相关的TE模式径向电场分量为

点击浏览公式

式中：E_1ρ^air(ρ,z)为海水层中空气波水平分量。若发射源位置h=z_s，接收器位置z=z_r，基于Baños^[13]的远场近似结果，式(12)可近似为与发射器水深、空气界面中的传播距离及接收器水深相关的3个乘积项：

点击浏览公式

式中：k₀≈0 m^-1；k₁=(iωμ₀σ₁)^1/2；ω是圆频率；θ是源的方位角；z_s和z_r分别为源和接收器的深度。 3 有限水深模型中空气波的近似表达与压制 3.1 有限水深模型中空气波的多次鸣震近似

基于远场近似关系式(13)知道，在半空间海水层的情况下，空气波主要与4个因素相关，即海水层的电导率σ₁、源与接收器的深度z_s和z_r，以及源与接收器之间的偏移距ρ。对于有限水深，Bannister^[14]对侧面电磁波的研究表明，水平电偶极子源的TE模式分量基本由源垂直向上传播，在海水与空气界面感应产生空气波，在海底地层和海水界面之间形成一系列快速衰减的交混回响信号。同样，由源激发的垂直向下传播的TE模式信号到达海底也形成向上反射的交混回响衰减信号。类似地，在接收器一侧，最初垂直向下的TE模式信号在海底地层和海水面之间也发生交混回响。然而，由于海水层电导率较大，其趋肤深度较小，在海水层中对空气波分量有较大贡献的交混回响次数不会太大。因此，对于有限水深的空气波传播，笔者从发射器一侧、海水空气界面和接收器一侧3个路径上电磁信号的传播和交混回响来分析。

图 3给出了源激发的空气波在海水和海底界面之间反射的示意图。海水的电导率为σ₁，深度为z_b处海底地层的电导率为σ₂。图 3a给出由z_s源处激发向上扩散的信号到达z=0处的海水面。由式(13)第一个指数项，发射器一侧的影响可以近似地表示为^[24]

点击浏览公式

式中，复波数k₁决定了信号的衰减。

a.来自源的初始上行波场；b,c.由初始上行波场感应的多次交混回响；d.初始下行波场和海底反射在海水面感应的空气波；e,f.由下行波场诱发的多次交混回响。 图 3 源侧反射和交混回响信号 Fig. 3 Signals of the reflection and reverberation near the source side

图选项

图 3d为源产生的信号向下扩散再在海底反射向上传播，并在海水面诱发空气波。参考式(14)在z=0处的空气波可用下式近似：

点击浏览公式

其中：R=(-)/(+)，是TE模式垂直传播的平面波反射系数；z_b为有限水深时的海底深度。

图 3b和图 3c给出了图 3a中描述的向上传播的源信号在空气波产生前所发生的一次或二次鸣震，对应的近似关系分别为^[24]：

点击浏览公式

点击浏览公式

同样，图 3e和图 3f给出了图 3d中最初向下传播的源信号在产生空气波前发生了一次或二次鸣震，对应的关系分别为^[24]：

点击浏览公式

点击浏览公式

空气波交混回响过程的无限持续所形成的总信号为^{[7, 24]}

点击浏览公式

其中，

S是一个滤波器，它表示源一侧产生的所有空气波分量。两个花括号中的项依次代表来自源的初始信号和海底反射信号。

类似地，图 4中给出了空气波分量由海底反射和鸣震在接收器一侧的改造作用。在图 4a中，空气波分量到达z_r处的接收器；然后，它在海底反射又回到接收器。这一过程可以近似地描述为^[24]

点击浏览公式

图 4b和图 4c分别考虑海水层中的一次和二次鸣震，其数学表达为^[24]

点击浏览公式

点击浏览公式

a.接收器记录的来自空气波的下行波场和海底反射；b,c.接收器一侧经过多次鸣震改造的空气波下行波场和海底反射。 图 4 接收器一侧反射和交混回响信号 Fig. 4 Signals of the reflection and reverberation near the receiver side

图选项

最后，考虑无限多次鸣震的衰减，接收器一侧的反射和鸣震级数可表示为^[24]

点击浏览公式

式中，，代表海底反射及鸣震在接收器一侧对空气波的滤波作用。

对于有限水深模型的空气波，除了式(13)中的项外，还包括了源一侧和接收器一侧的反射和鸣震：

点击浏览公式

方程(25)假设了HED源激发的TE模式场最初是在海水层中垂直向上和向下传播的，且也在海水层中激发一系列垂直传播并以指数形式衰减的反射和鸣震级数，每个TE模式分量在海水表面诱发幅度逐渐衰减的空气波分量。在接收器一侧，也发生与源一侧类似的反射和鸣震。 3.2 崎岖海底和海水电导率横向变化模型中空气波的近似

在导出一般空气波表达式(25)时，假设海水电导率为常数。然而，实践中，由于海水中的含盐度和温度梯度随深度变化，海水电导率随深度发生微小变化。若电导率的垂向梯度已知，那么在空气波的关系中可以将k₁z项用∫₀^zk₁(z′)dz′来代替。方程(25)空气波近似关系对MCSEM的长偏移距响应分析仍然适用。对于源和接收器深度存在较大变化的崎岖海底(图 5)，前文海水深度不变的假设可以放宽，因为前文鸣震近似的分析中，源一侧和接收器一侧反射和鸣震的影响是解耦的。

设源一侧的水深为z_bs，海底的反射系数为 R_s：

式中：σ_1s为源位置海水的电导率；σ_2s为z_bs处海底地层的电导率。将R_s代入式(20)中的S^{[7, 22, 24]}：

点击浏览公式

同样，在接收器一侧，若海水深为z_br及海底反射系数为R_r：

式中：σ_1r为接收器位置海水的电导率；σ_2r为z_br处海底地层的电导率。用z_br和R_r代替方程(24)中的z_b和R，可以得到新的^[24]：

点击浏览公式

利用式(26)和式(27)的修正关系，可在源和接收器附近为局部平坦崎岖海底(图 5)实现空气波压制。它在除海底反射外局部崎岖区域不发生较大散射时，可以较好地模拟崎岖海底的空气波。需要注意的是，前文分析的鸣震近似方法主要考虑了发射和接收机附近电性分布变化的影响，同时在式(26)和式(27)中也考虑了源和接收器处的海水深度变化，但没有考虑发射接收间距较大情况下电性变化和海底起伏较大时对空气波的影响，这也是本方法有待改进的致命缺陷之一。

图 5 源和接收器在深度存在较大差异的崎岖海底 Fig. 5 Rugged seabed model with different depth levels of source and receiver positions

图选项

为分析前述空气波近似关系(式(27))的空气波压制效果，利用有限水深模型(图 6a)计算的近似空气波的可探测性与无限水深模型(图 6b)进行对比，分析交混回响近似。

a.100 m水深五层和三层模型；b.15 000 m水深五层和三层模型。 图 6 空气波压制效果模型示意图 Fig. 6 Models for indication of the air wave mitigation effect

图选项

图 6a中左右两模型的差异为储层的有无：左图是五层模型，存在一层电导率为0.02 S/m的储层；右图是三层模型，无储层。图 6b中两模型与图 6a的差异为海水层的厚度，图 6b中两模型海水深度都为15 000 m，用于近似模拟无限水深模型，计算无空气波时的电场幅度。

图 7a给出了100 m水深时五层模型和三层模型的原始电场径向分量曲线与鸣震近似关系压制空气波后的电场径向分量曲线。由图 7a看到：2种模型测量的原始数据(未压制)只在4 000~6 000 m的偏移距上出现微小的幅度差，显示了油气储层的存在；用鸣震近似关系压制后，电场幅度在3 000~5 000 m的偏移距上得到微小提升。

a.水深100 m；b.五层模型。 图 7 鸣震近似压制空气波结果对比 Fig. 7 Comparison of the air wave mitigation results by using reverberation approximate expression

图选项

图 7b给出了五层模型在15 000、100 m水深时的原始电场响应与鸣震近似关系压制空气波后100 m水深的标准化电场幅度曲线(分别用各自的三层模型标准化)。从图 7b中看出：由于空气波的影响，水深100 m时，原始电场曲线(未压制)难以显示油气储层，而经过鸣震近似关系压制空气波后，在3 000~5 400 m的偏移距上显示了小幅度提升，但在5 400~10 000 m又呈现小幅度下降，该空气波压制方法效果仍有待改善；将水深100 m鸣震近似后的标准化幅度曲线与水深15 000 m的原始电场响应(未压制)的标准化幅度曲线对比知道，前文实现的鸣震近似空气波压制方法仍然难以较好地消除空气波的影响，这主要是由于前文的鸣震近似方法未考虑占比较大的非垂直传播散射场的影响。 4 上行波与下行波分离及空气波压制方法 4.1 散射电磁波的上行与下行波场分离

前文的鸣震近似压制空气波的方法主要将空气波的传播近似看成海水层中垂直传播和空气层中的测量传播的波场，这是将散射场近似成波场的考虑方法。类似地，将散射场近似看成向上和向下衰减的波场，Løseth等^[25]和Andréis等^[26]提出了将接收器记录的散射场分解成上行和下行衰减场，进而压制空气波的方法。下面以他们^{[25, 26]}的波场分解理论为基础，分析上行和下行散射场分解压制空气波的方法。

假设接收器记录的电磁响应是在某个电导率、磁导率和介电常数为已知常数的平面上获得的，设x=(x₁,x₂,x₃)=(x,y,z)，记深度轴向下为正，与式(1)相似，各向同性无源介质中电磁场方程(1)退化为^{[27, 28]}

点击浏览公式

其中：电场E=(E₁,E₂,E₃)^T；磁场H=(H₁,H₂,H₃)^T；=ε-iσ/ω。对方程(28)作空间水平坐标的傅里叶变换，且用代数运算消去E₃和H₃，得到常微分方程^[25]:

点击浏览公式

其中：

式中：s为水平慢度，s_i=k_i/ω，i=1,2,s²=s₁²+s₂²；q为垂直慢度，q=，为介质的复速度，且垂直慢度的2个特例为q_i=，i=1,2。

在MCSEM测量的低频情况下，复介电常数和复速度为=-iσ/ω，。为了书写方便,下面将各量对频率、波数、深度等显性的依赖关系中相应的参数省略掉。例如：E=E(k₁,k₂,x₃,ω)。对于3D均匀介质中的电偶极子源，电磁场方程存在4个独立的解。其中：2个向上传播的电磁场分量，在源上方；2个向下传播，在源下方^[8]。源的性质决定了这4个分量的幅度。一般地，电场和磁场都是由上行u和下行d传播的分量构成，因此，

点击浏览公式

为了导出电场和磁场之间的关系以及它们与上行和下行分量之间的关系，需要找出系统矩阵A的特征函数，即如下方程的解：

点击浏览公式

其中：λ_N是特征值；而l_N是对应λ_N的特征向量。4个特征值之间有如下关系：

点击浏览公式

它们可以组合成4×4的对角阵：

点击浏览公式

类似地，将特征向量作为列向量也可以形成4×4的矩阵：

点击浏览公式

那么方程(31)的特征解可表示为

点击浏览公式

在介质空间的均匀区域，系统矩阵A是常数，从而特征矢量矩阵L与x₃无关。对于电磁场矢量b的一阶微分方程(29)可以转化成下面的一阶微分方程：

点击浏览公式

引入4×1的列向量w=L^-1b，方程(36)对每个w元素变为对上行波场和下行波场解耦的微分方程，分别代表 2个上行波场和2个下行波场的微分方程。因此,w中的元素代表了上行波和下行波场。对特征向量作如下标准化可以得到

点击浏览公式

其中，L^-1的作用是将电磁场矢量分解成上行波和下行波分量，因此，L^-1称为分解矩阵。

事实上，分解只是物理变量的简单线性变换，它是一个稳定的计算过程。由b=Lw知道，L的作用是将上行和下行衰减场分量组合成电磁场矢量，L称为组合矩阵。在波动理论中，λ_N代表相位慢度，l_N表示极化矢量^[29]，将电磁场矢量分解成上行和下行衰减场的具体实现过程见文献^[8]。 4.2 电磁波场分离压制空气波的效果

为了考察电磁波场分离压制空气波的效果并与前文的鸣震近似压制空气波效果对比，这里仍然沿用图 6a所示的层状介质模型。偶极子源发射0.25 Hz的正弦波。

2种模型下，利用一维频域模拟算法计算得到的E₁幅度曲线如图 8a所示。从图 8a看到:2种模型测量的原始数据(未压制)只在3 000~6 500 m的偏移距上出现微小的幅度差，显示了油气储层存在；上、下行波场分离与空气波压制后，电场幅度在空气波较强的大偏移距上得到一定提升。

a.水深100 m；b.五层模型。 图 8 上行波和下行波分离压制空气波结果对比 Fig. 8 Comparison of the air wave mitigation results by using decomposition of up going and down going wave fields

图选项

图 8b给出了五层模型在15 000、100 m水深的原始电场响应(未压制)与水深100 m时上下行波分离压制后的标准化电场幅度曲线。从图 8b中看出：水深100 m时，由于空气波的影响，原始电场曲线(未压制)难以显示油气储层，而经过下行波分离和压制后，上行波场幅度有所提升;在3 000~8 500 m的偏移距上显示了小幅度异常，显示了油气储层的存在。然而，将水深100 m时上下行波分离压制后的标准化幅度曲线与水深15 000 m时原始电场(未压制)的标准化幅度曲线对比可知，前文描述的上行、下行波场分离的空气波压制方法仍然难以完全压制空气波的影响。这主要是由于前文的分离方法难以将与地下地层信号相互耦合的空气波影响消除，而且这些影响在空气波响应中占较大比重。 5 结论和认识

基于Maxwell方程和模式分解理论导出了无限水深时空气波的完整表达式，据此分析了空气波的作用机理和近似表达式。由无限水深空气波的完整表达式知道，在MCSEM测量中空气波分量主要来自于HED源场的TE模式，它几乎从源垂直向上及向下传播；在海水面，上行波信号产生空气波，在海水面和海底之间激发一系列快速衰减的鸣震信号，其中每个信号在海水面又激发额外的空气波信号。从激发源一侧和接收器一侧空气波反射和鸣震的近似数学表达出发，建立了层状介质和崎岖海底的空气波近似耦合关系，给出了可用于崎岖海底空气波近似分析的方法。数值模拟结果显示，鸣震近似空气波模拟方法，只是模拟了垂直入射和反射的空气波部分，它们只占空气波中的极小比例，因此，空气波压制效果不甚理想，若期望进一步改善效果，只有考虑所有非垂直入射和与地下地层信号耦合的多次空气波才可能实现。另外，基于电磁场上行波和下行波分解的方法，虽然从理论上可以较好地压制空气波，但数值模拟结果显示，该方法对空气波压制及提升有效信号的比例仍然有限。这主要是由于上行波场中仍然包含很大部分空气波与地层信号相互干涉的多次波场。对比文中2种压制空气波方法的应用效果发现，对于实际测量的情况，研究如何分离及消除空气波与地层信号相互干涉的多次波场，而尽可能不损伤地层有效信息的压制空气波的方法，才是解决浅水域MCSEM勘探的根本出路。