裂缝储层作为一种重要的含油气储层, 其中发育的裂缝既可以作为油气储集空间, 又可以作为渗滤通道。裂缝密度、发育方向、开度及流体饱和等情况, 是指导勘探开发的重要指标, 如何准确预测这些参数一直是难点问题。由于实际地层中裂缝发育情况是复杂的, 一般建立等效模型^[1]来简化裂缝, 便于研究。裂缝的等效化模型大致可以分为3种:线性滑移模型^[2]、Hudson模型^[3-4]以及Thomson模型^[5]; 线性滑移理论是近年最流行的等效理论, 因为它忽略了裂缝的微观结构, 仅仅使用裂缝法向和切向柔度参数就可以描述裂缝, 而且两者的比值可以有助于确定岩石裂缝和主体骨架的岩石物理参数。当地层中发育一组平行高角度的裂缝时, 可以将其等效为横观各向同性(HTI)介质, 通过Ruger近似公式^[6]建立裂缝储层弹性物性参数与纵波方位地震数据之间的定量数学关系, 最终可以利用反演方法^[7-8]得到裂缝信息。叠前方位数据的多参数反演时, 反演方程的系数矩阵条件数大, 反演过程会变得极不稳定, 相比于常规叠后反演, 其固有的病态性更为严重。为了避免求解过程中大量复杂的矩阵运算, Downton等^[9]应用傅里叶级数展开理论, 针对HTI介质, 将地震数据在方位角上展开为多项与入射角度有关的三角函数, 将反演矩阵求解问题转换为方位数据求和问题, 显著提高了计算速度与计算结果的准确性。但实际地层中通常发育着多个方向的裂缝, 而且经常是相互正交的^[10-11], 显然HTI介质无法满足实际要求, 这时需要对复杂介质做更深入研究。笔者在前人研究的基础上, 结合线性滑移理论, 对正交各向异性裂缝介质进行等效模型化, 进一步研究傅里叶级数展开法在正交裂缝参数预测中的应用。

1 旋转对称正交裂缝介质刚度表示

线性滑动理论忽略了裂缝的形状和微结构, 认为裂缝将地层分割成彼此平行的柔性面, 各柔性面在边界处连续, 满足线性滑动边界条件。含有裂缝岩石的柔度S表示为各向同性背景围岩的柔度S_b与裂缝柔度ΔS之和^[12]:

$ \mathit{\boldsymbol{S}} = {\mathit{\boldsymbol{C}}^{ - 1}} = {\mathit{\boldsymbol{S}}_{\rm{b}}} + \sum\limits_{n = 1}^n \Delta {\mathit{\boldsymbol{S}}^n}. $

(1)

式中, n为裂缝组号; C为含有裂缝岩石的刚度矩阵。各向同性背景介质柔度矩阵S_b为

$ {\mathit{\boldsymbol{S}}_{\rm{b}}} = \left[ {\begin{array}{*{20}{c}} {\frac{{{\lambda _{\rm{b}}} + {\mu _{\rm{b}}}}}{{{\mu _{\rm{b}}}\left( {3{\lambda _{\rm{b}}} + 2{\mu _{\rm{b}}}} \right)}}}&{\frac{{ - {\lambda _{\rm{b}}}}}{{2{\mu _{\rm{b}}}\left( {3{\lambda _{\rm{b}}} + 2{\mu _{\rm{b}}}} \right)}}}&{\frac{{ - {\lambda _{\rm{b}}}}}{{2{\mu _{\rm{b}}}\left( {3{\lambda _{\rm{b}}} + 2{\mu _{\rm{b}}}} \right)}}}&0&0&0\\ {\frac{{ - {\lambda _{\rm{b}}}}}{{2{\mu _{\rm{b}}}\left( {3{\lambda _{\rm{b}}} + 2{\mu _{\rm{b}}}} \right)}}}&{\frac{{{\lambda _{\rm{b}}} + {\mu _{\rm{b}}}}}{{{\mu _{\rm{b}}}\left( {3{\lambda _{\rm{b}}} + 2{\mu _{\rm{b}}}} \right)}}}&{\frac{{ - {\lambda _{\rm{b}}}}}{{2{\mu _{\rm{b}}}\left( {3{\lambda _{\rm{b}}} + 2{\mu _{\rm{b}}}} \right)}}}&0&0&0\\ {\frac{{ - {\lambda _{\rm{b}}}}}{{2{\mu _{\rm{b}}}\left( {3{\lambda _{\rm{b}}} + 2{\mu _{\rm{b}}}} \right)}}}&{\frac{{ - {\lambda _{\rm{b}}}}}{{2{\mu _{\rm{b}}}\left( {3{\lambda _{\rm{b}}} + 2{\mu _{\rm{b}}}} \right)}}}&{\frac{{{\lambda _{\rm{b}}} + {\mu _{\rm{b}}}}}{{{\mu _{\rm{b}}}\left( {3{\lambda _{\rm{b}}} + 2{\mu _{\rm{b}}}} \right)}}}&0&0&0\\ 0&0&0&{\frac{1}{{{\mu _{\rm{b}}}}}}&0&0\\ 0&0&0&0&{\frac{1}{{{\mu _{\rm{b}}}}}}&0\\ 0&0&0&0&0&{\frac{1}{{{\mu _{\rm{b}}}}}} \end{array}} \right]. $

(2)

式中, λ_b、μ_b为背景围岩的拉梅参数。对于正交裂缝, 具有两组不同的裂缝, 则式(1)可以展开为

$ \mathit{\boldsymbol{S}} = {\mathit{\boldsymbol{S}}_{\rm{b}}} + \Delta {\mathit{\boldsymbol{S}}_{{\rm{fl}}}} + \Delta {\mathit{\boldsymbol{S}}_{{\rm{f}}2}}. $

(3)

若以x₁轴平行于第一组裂缝的对称轴方向, x₂轴平行于第二组裂缝对称轴方向建立坐标系, 并假设裂缝具有旋转对称轴, 且其在不同方向上没有相互作用^[12-14], 则两组裂缝的附加柔度矩阵ΔS_f1、ΔS_f2可以表示为

$ \Delta {\mathit{\boldsymbol{S}}_{{\rm{fl}}}} = {\mathop{\rm diag}\nolimits} \left( {{S_{{\rm{NI}}}}\quad 0\quad 0\quad 0\quad {S_{{\rm{T1}}}}\quad {S_{{\rm{T1}}}}} \right), $

(4)

$ \Delta {\mathit{\boldsymbol{S}}_{12}} = {\mathop{\rm diag}\nolimits} \left( {0\quad {S_{{\rm{N}}2}}\quad 0\quad {S_{{\rm{T}}2}}\quad 0\quad {S_{{\rm{T}}2}}} \right). $

(5)

式中, S_N、S_T分别为裂缝附加柔度矩阵ΔS中的法向柔度和切向柔度, 下标1、2分别对应第一组裂缝和第二组裂缝。

引入Hsu等^[15]提出的无量纲参数弱度项:

$ {\delta _{{\rm{N}}i}} = \frac{{\left( {{\lambda _{\rm{b}}} + 2{\mu _{\rm{b}}}} \right){S_{{\rm{N}}i}}}}{{1 + \left( {{\lambda _{\rm{b}}} + 2{\mu _{\rm{b}}}} \right){S_{{\rm{N}}i}}}},{\delta _{{\rm{T}}i}} = \frac{{{\mu _{\rm{b}}}{S_{{\rm{T}}i}}}}{{1 + {\mu _{\rm{b}}}{S_{{\rm{T}}i}}}},i = 1,2. $

(6)

式中, δ_N1和δ_T1分别为第一组裂缝的法向弱度和切向弱度; δ_N2和δ_T2分别为第二组裂缝的法向弱度和切向弱度。

进一步对正交裂缝介质的柔度矩阵S(式(3))求逆, 并将裂缝弱度(式(6))带入化简, 得到正交裂缝介质的刚度矩阵C为

$ \mathit{\boldsymbol{C}} = \left[ {\begin{array}{*{20}{c}} {{C_{11}}}&{{C_{12}}}&{{C_{13}}}&0&0&0\\ {{C_{12}}}&{{C_{22}}}&{{C_{23}}}&0&0&0\\ {{C_{13}}}&{{C_{23}}}&{{C_{33}}}&0&0&0\\ 0&0&0&{{C_{44}}}&0&0\\ 0&0&0&0&{{C_{55}}}&0\\ 0&0&0&0&0&{{C_{66}}} \end{array}} \right]. $

(7)

其中刚度矩阵C的各个系数为

$ \left\{ \begin{array}{l} {C_{11}} = \frac{{\left( {{\lambda _{\rm{b}}} + 2{\mu _{\rm{b}}}} \right)\left( {1 - {\delta _{{\rm{NI}}}}} \right)\left( {1 - {r^2}{\delta _{{\rm{NI}}}}} \right)}}{{1 - {r^2}{\delta _{{\rm{N}}1}}{\delta _{{\rm{N}}2}}}},\\ {C_{12}} = \frac{{{\lambda _{\rm{b}}}\left( {1 - {\delta _{{\rm{N}}1}}} \right)\left( {1 - {\delta _{{\rm{N}}2}}} \right)}}{{1 - {r^2}{\delta _{{\rm{N}}1}}{\delta _{{\rm{N}}2}}}},\\ {C_{13}} = \frac{{{\lambda _{\rm{b}}}\left( {1 - {\delta _{{\rm{N1}}}}} \right)\left( {1 - r{\delta _{{\rm{N}}2}}} \right)}}{{1 - {r^2}{\delta _{{\rm{N}}1}}{\delta _{{\rm{N}}2}}}},\\ {C_{22}} = \frac{{\left( {{\lambda _{\rm{b}}} + 2{\mu _{\rm{b}}}} \right)\left( {1 - {r^2}{\delta _{{\rm{NI}}}}} \right)\left( {1 - {\delta _{{\rm{N}}2}}} \right)}}{{1 - {r^2}{\delta _{{\rm{N}}1}}{\delta _{{\rm{N}}2}}}},\\ {C_{23}} = \frac{{{\lambda _{\rm{b}}}\left( {1 - r{\delta _{{\rm{NI}}}}} \right)\left( {1 - {\delta _{{\rm{N}}2}}} \right)}}{{1 - {r^2}{\delta _{{\rm{NI}}}}{\delta _{{\rm{N}}2}}}},\\ {C_{33}} = \frac{{\left( {{\lambda _{\rm{b}}} + 2{\mu _{\rm{b}}}} \right)\left[ {\left( {1 - {r^2}{\delta _{{\rm{N1}}}}} \right)\left( {1 - {r^2}{\delta _{{\rm{N}}2}}} \right) - 4{r^2}{g^2}{\delta _{{\rm{N}}1}}{\delta _{{\rm{N}}2}}} \right]}}{{1 - {r^2}{\delta _{{\rm{N}}1}}{\delta _{{\rm{N}}2}}}},\\ {C_{44}} = \mu \left( {1 - {\delta _{{\rm{T}}2}}} \right),\\ {C_{55}} = \mu \left( {1 - {\delta _{{\rm{T}}1}}} \right),\\ {C_{66}} = \frac{{\mu \left( {1 - {\delta _{{\rm{T1}}}}} \right)\left( {1 - {\delta _{{\rm{T2}}}}} \right)}}{{1 - {\delta _{{\rm{N}}1}}{\delta _{{\rm{N}}2}}}}. \end{array} \right. $

(8)

其中

$ g=\mu_{\mathrm{b}} /\left(\lambda_{\mathrm{b}}+2 \mu_{\mathrm{b}}\right), $

$ r = 1 - 2g. $

因为裂缝法向弱度δ_N1、δ_N2数量级较小, 所以δ_N1δ_N2乘积近似为零, 式(8)化简为

$ \left\{ \begin{array}{l} {C_{11}} = \left( {{\lambda _{\rm{b}}} + 2{\mu _{\rm{b}}}} \right)\left( {1 - {\delta _{{\rm{N}}1}} - {r^2}{\delta _{{\rm{N}}2}}} \right),\\ {C_{12}} = {\lambda _{\rm{b}}}\left( {1 - {\delta _{{\rm{N}}1}} - {\delta _{{\rm{N}}2}}} \right),\\ {C_{13}} = {\lambda _{\rm{b}}}\left( {1 - {\delta _{{\rm{N}}1}} - r{\delta _{{\rm{N}}2}}} \right),\\ {C_{22}} = \left( {{\lambda _{\rm{b}}} + 2{\mu _{\rm{b}}}} \right)\left( {1 - {\delta _{{\rm{N}}2}} - {r^2}{\delta _{{\rm{N}}1}}} \right),\\ {C_{23}} = {\lambda _{\rm{b}}}\left( {1 - {\delta _{{\rm{N}}2}} - r{\delta _{{\rm{N}}1}}} \right),\\ {C_{33}} = \left( {{\lambda _{\rm{b}}} + 2{\mu _{\rm{b}}}} \right)\left( {1 - {r^2}{\delta _{{\rm{N}}1}} - {r^2}{\delta _{{\rm{N}}2}}} \right),\\ {C_{44}} = {\mu _{\rm{b}}}\left( {1 - {\delta _{{\rm{T}}2}}} \right),\\ {C_{55}} = {\mu _{\rm{b}}}\left( {1 - {\delta _{{\rm{T}}1}}} \right),\\ {C_{66}} = {\mu _{\rm{b}}}\left( {1 - {\delta _{{\rm{T}}1}} - {\delta _{{\rm{T}}2}}} \right). \end{array} \right. $

(9)

2 正交各向异性介质方位反射系数推导

对于弱各向异性介质, 任意对称轴纵波反射系数方程^[16]表示为

$ \begin{array}{l} {R_{{\rm{PP}}}}(\theta ,\varphi ) = R_{{\rm{PP}}}^{{\rm{iso}}}(\theta ) + \frac{1}{2}\Delta {\varepsilon _z} + \frac{1}{2}\left[ {\left( {\Delta {\delta _x} - 8\frac{{{{\bar \beta }^2}}}{{{{\bar \alpha }^2}}}\Delta {\gamma _x}} \right){{\cos }^2}\varphi + } \right.\\ \left( {\Delta {\delta _y} - 8\frac{{{{\bar \beta }^2}}}{{{{\bar \alpha }^2}}}\Delta {\gamma _y}} \right){\sin ^2}\varphi + 2\left( {\Delta {\chi _z} - 4\frac{{{{\bar \beta }^2}}}{{{{\bar \alpha }^2}}}\Delta {\varepsilon _{45}}} \right)\cos \varphi \sin \varphi - \\ \Delta {\varepsilon _z}]{\sin ^2}\theta + \frac{1}{2}\left[ {\Delta {\varepsilon _x}{{\cos }^4}\varphi + \Delta {\varepsilon _y}{{\sin }^4}\varphi + \Delta {\delta _z}{{\cos }^2}\varphi {{\sin }^2}\varphi + } \right.\\ 2\left( {\Delta {\varepsilon _{16}}{{\cos }^2}\varphi + \Delta {\varepsilon _{26}}{{\sin }^2}\varphi } \right)\cos \varphi \sin \varphi ]{\sin ^2}\theta {\tan ^2}\theta . \end{array} $

(10)

式中, α和β分别为纵波、横波速度, 上标横线代表上下两层介质参数的平均值; θ为入射角; φ为方位角; Δ为界面两侧参数的相对变化量; ε_x、ε_y、ε_z、ε₁₆、ε₂₆、ε₄₅、δ_x、δ_y、δ_z、χ_z、γ_x和γ_y为扰动参数。扰动参数ε_x、ε_y、ε_z、ε₁₆、ε₂₆、ε₄₅、δ_x、δ_y、δ_z、χ_z、γ_x和γ_y定义如下:

$ \left\{ {\begin{array}{*{20}{l}} {{\varepsilon _x} = \frac{{{A_{11}} - {\alpha ^2}}}{{2{\alpha ^2}}},{\varepsilon _y} = \frac{{{A_{22}} - {\alpha ^2}}}{{2{\alpha ^2}}},{\varepsilon _z} = \frac{{{A_{33}} - {\alpha ^2}}}{{2{\alpha ^2}}},}\\ {{\varepsilon _{16}} = \frac{{{A_{16}}}}{{{\alpha ^2}}},{\varepsilon _{26}} = \frac{{{A_{26}}}}{{{\alpha ^2}}},}\\ {{\delta _x} = \frac{{{A_{13}} + 2{A_{55}} - {\alpha ^2}}}{{{\alpha ^2}}},{\delta _y} = \frac{{{A_{23}} + 2{A_{44}} - {\alpha ^2}}}{{{\alpha ^2}}},{\delta _z} = \frac{{{A_{12}} + 2{A_{66}} - {\alpha ^2}}}{{{\alpha ^2}}},}\\ {{\chi _z} = \frac{{{A_{36}} + 2{A_{45}}}}{{{\alpha ^2}}},{\gamma _x} = \frac{{{A_{55}} - {\beta ^2}}}{{2{\beta ^2}}},{\gamma _y} = \frac{{{A_{44}} - {\beta ^2}}}{{2{\beta ^2}}},{\varepsilon _{45}} = \frac{{{A_{45}}}}{{{\beta ^2}}}.} \end{array}} \right. $

(11)

其中

$ A_{i j}=C_{i j} / \rho. $

式中, A_ij为密度ρ标准化的刚度参数。

各向同性介质的纵波反射系数R_PP ^iso (θ)可以表示为

$ R_{\mathrm{PP}}^{\mathrm{iso}}(\theta)=A_{\mathrm{iso}}+B_{\mathrm{iso}} \sin ^{2} \theta+C_{\mathrm{iso}} \sin ^{2} \theta \tan ^{2} \theta, $

(12)

$ \left\{ {\begin{array}{*{20}{l}} {{A_{{\rm{iso}}}} = \frac{1}{2}\frac{{\Delta Z}}{{\bar Z}},}\\ {{B_{{\rm{iso}}}} = \frac{1}{2}\frac{{\Delta \alpha }}{{2\bar \alpha }} - 4\frac{{{{\bar \beta }^2}}}{{{{\bar \alpha }^2}}}\left[ {\frac{{\Delta \beta }}{\beta } + \frac{{\Delta \rho }}{{2\bar \rho }}} \right],}\\ {{C_{{\rm{iso}}}} = \frac{{\Delta \alpha }}{{2\bar \alpha }},Z = \rho \alpha .} \end{array}} \right. $

(13)

进一步将正交裂缝的刚度矩阵C(式(1))带入式(11), 得到正交裂缝反射系数的扰动参数为

$ \left\{ \begin{array}{l} {\chi _z} = {\varepsilon _{16}} = {\varepsilon _{26}} = {\varepsilon _{45}} = 0,\\ {\delta _x} = - r{\delta _{{\rm{N}}1}} - {r^2}{\delta _{{\rm{N}}2}} - 2g{\delta _{{\rm{T}}1}},\\ {\delta _y} = - r{\delta _{{\rm{N}}2}} - {r^2}{\delta _{{\rm{N}}1}} - 2g{\delta _{{\rm{T}}2}},\\ {\delta _z} = - r\left( {{\delta _{{\rm{N}}1}} + {\delta _{{\rm{N}}2}}} \right) - 2g\left( {{\delta _{{\rm{T}}1}} + {\delta _{{\rm{T}}2}}} \right),\\ {\varepsilon _x} = - \frac{1}{2}\left( {{r^2}{\delta _{{\rm{N}}2}} + {\delta _{{\rm{N}}1}}} \right),\\ {\gamma _x} = - \frac{1}{2}g{\delta _{{\rm{T}}1}},\\ {\varepsilon _y} = - \frac{1}{2}\left( {{r^2}{\delta _{{\rm{N}}1}} + {\delta _{{\rm{N}}2}}} \right),\\ {\varepsilon _z} = - \frac{1}{2}{r^2}\left( {{\delta _{{\rm{N}}1}} + {\delta _{{\rm{N}}2}}} \right),\\ {\gamma _y} = - \frac{1}{2}g{\delta _{{\rm{T}}2}}. \end{array} \right. $

(14)

将式(14)带入式(10), 最终得到正交介质的反射系数表达式为

$ \begin{array}{l} {R_{{\rm{PP}}}}(\theta ) = \frac{1}{2}\frac{{\Delta Z}}{{\bar Z}} + \frac{1}{2}\left( {\frac{{\Delta \alpha }}{{\bar \alpha }} - 4\frac{{{{\bar \beta }^2}}}{{{{\bar \alpha }^2}}}\frac{{\Delta \mu }}{\mu }} \right){\sin ^2}\theta + \\ \frac{1}{2}\frac{{\Delta \alpha }}{{\bar \alpha }}{\sin ^2}\theta {\tan ^2}\theta - \frac{1}{4}{r^2}\left( {\Delta {\delta _{{\rm{N}}1}} + \Delta {\delta _{{\rm{N}}2}}} \right) - \\ \left[ {\frac{1}{4}r(1 + 2g)\left( {\Delta {\delta _{{\rm{N}}1}} + \Delta {\delta _{{\rm{N}}2}}} \right) - } \right.\\ gr\left( {{{\sin }^2}\varphi \Delta {\delta _{{\rm{NI}}}} + {{\cos }^2}\varphi \Delta {\delta _{{\rm{N}}2}}} \right) - \\ g\left( {\Delta {\delta _{{\rm{T}}2}}{{\sin }^2}\varphi + \Delta {\delta _{{\rm{T}}1}}{{\cos }^2}\varphi } \right)]{\sin ^2}\theta - \\ \frac{1}{4}\left[ {{{\left( {1 - 2g{{\sin }^2}\varphi } \right)}^2}\Delta {\delta _{{\rm{NI}}}} + {{\left( {1 - 2g{{\cos }^2}\varphi } \right)}^2}\Delta {\delta _{{\rm{N}}2}} + } \right.\\ 4g\left( {\Delta {\delta _{{\rm{TI}}}} + \Delta {\delta _{{\rm{T}}2}}} \right){\sin ^2}\varphi {\cos ^2}\varphi ]{\sin ^2}\theta {\tan ^2}\theta . \end{array} $

(15)

3 正交裂缝方位反射系数傅里叶级数表示推导

研究表明任意弱各向异性介质的反射系数可以视为以方位角为周期的函数^[9], 表示成三角函数之和的傅里叶级数形式为

$ R(\varphi ,\theta ) = \sum\limits_{n = 0}^\infty {\left( {{u_n}(\theta )\cos (n\varphi ) + {v_n}(\theta )\sin (n\varphi )} \right)} . $

(16)

对于具有N个采样点的方位数据, 余弦和正弦函数前的系数u_n(θ)和v_n(θ)定义为

$ u_{n}(\theta)=\frac{1}{\pi} \sum\limits_{k=1}^{N} R_{k}(\varphi, \theta) \cos (n \varphi) \mathrm{d} \varphi, $

(17)

$ v_{n}(\theta)=\frac{1}{\pi} \sum\limits_{k=1}^{N} R_{k}(\varphi, \theta) \sin (n \varphi) \mathrm{d} \varphi. $

(18)

对于弱各向异性介质反射系数的近似, 使用四阶傅里叶级数展开精度就可以, 因为更高阶项的值很小, 且无法分清是否是噪音, 可以忽略, 因此反射系数的四阶傅里叶级数展开表示如下:

$ \begin{array}{*{20}{l}} {R(\varphi ,\theta ) = {u_0}(\theta ) + {v_2}(\theta )\sin (2\varphi ) + {u_2}(\theta )\cos (2\varphi ) + }\\ {{v_4}(\theta )\sin (4\varphi ) + {u_4}(\theta )\cos (4\varphi ).} \end{array} $

(19)

其中各阶三角函数前的系数代表的是任意对称轴介质反射系数公式(10)的线性组合, 具体表示为

$ \left\{ {\begin{array}{*{20}{l}} {{u_0}(\theta ) = {\omega _{00}} + {\omega _{10}}{{\sin }^2}\theta + {\omega _{20}}{{\sin }^2}\theta {{\tan }^2}\theta ,}\\ {{v_2}(\theta ) = {\omega _{11}}{{\sin }^2}\theta + {\omega _{21}}{{\sin }^2}\theta {{\tan }^2}\theta ,}\\ {{u_2}(\theta ) = {\omega _{12}}{{\sin }^2}\theta + {\omega _{22}}{{\sin }^2}\theta {{\tan }^2}\theta ,}\\ {{v_4}(\theta ) = {\omega _{23}}{{\sin }^2}\theta {{\tan }^2}\theta ,}\\ {{u_4}(\theta ) = {\omega _{24}}{{\sin }^2}\theta {{\tan }^2}\theta .} \end{array}} \right. $

(20)

参数ω_ij与扰动参数(式(11))之间的关系为

$ \left\{ \begin{array}{l} {\omega _{00}} = {A_{{\rm{iso}}}} + 1/\left( {2{\alpha ^2}} \right)\Delta \left( {{\alpha ^2}{\varepsilon _z}} \right),\\ {\omega _{10}} = {B_{{\rm{iso}}}} + 1/\left( {4{\alpha ^2}} \right)\left[ {\Delta \left( {{\alpha ^2}{\delta _x}} \right) - 8\Delta \left( {{\beta ^2}{\gamma _x}} \right)} \right] + \\ 1/\left( {4{\alpha ^2}} \right)\left[ {\Delta \left( {{\alpha ^2}{\delta _y}} \right) - 8\Delta \left( {{\beta ^2}{\gamma _y}} \right)} \right] - 1/\left( {2{\alpha ^2}} \right)\Delta \left( {{\alpha ^2}{\varepsilon _z}} \right),\\ {\omega _{20}} = {C_{{\rm{iso}}}} + 1/\left( {16{\alpha ^2}} \right)\left[ {3\Delta \left( {{\alpha ^2}{\varepsilon _x}} \right) + 3\Delta \left( {{\alpha ^2}{\varepsilon _y}} \right) + } \right.\\ \Delta \left( {{\alpha ^2}{\delta _z}} \right)],\\ {\omega _{11}} = 1/\left( {2{\alpha ^2}} \right)\left[ {\Delta \left( {{\alpha ^2}{\chi _z}} \right) - 4\Delta \left( {{\beta ^2}{\varepsilon _{45}}} \right)} \right],\\ {\omega _{12}} = 1/\left( {4{\alpha ^2}} \right)\left[ {\Delta \left( {{\alpha ^2}{\delta _x}} \right) - 8\Delta \left( {{\beta ^2}{\gamma _x}} \right) - \left( {\Delta \left( {{\alpha ^2}{\delta _y}} \right) - } \right.} \right.\\ 8\Delta \left( {{\beta ^2}{\gamma _y}} \right))],\\ {\omega _{21}} = 1/\left( {4{\alpha ^2}} \right)\left[ {\Delta \left( {{\alpha ^2}{\varepsilon _{16}}} \right) - \Delta \left( {{\alpha ^2}{\varepsilon _{26}}} \right)} \right],\\ {\omega _{22}} = 1/\left( {4{\alpha ^2}} \right)\left[ {\Delta \left( {{\alpha ^2}{{\cal E}_x}} \right) - \Delta \left( {{\alpha ^2}{\varepsilon _y}} \right)} \right],\\ {\omega _{23}} = 1/\left( {8{\alpha ^2}} \right)\left[ {\Delta \left( {{\alpha ^2}{\varepsilon _{16}}} \right) - \Delta \left( {{\alpha ^2}{\varepsilon _{26}}} \right)} \right],\\ {\omega _{24}} = 1/\left( {16{\alpha ^2}} \right)\left[ {\Delta \left( {{\alpha ^2}{\varepsilon _x}} \right) - \Delta \left( {{\alpha ^2}{\varepsilon _y}} \right) - \Delta \left( {{\alpha ^2}{\delta _z}} \right)} \right]. \end{array} \right. $

(21)

将得到的正交裂缝扰动参数式(14)带入式(21), 化简后可以得到任意观测坐标系下正交介质傅里叶反射系数的具体表达式:

$ \begin{array}{l}{R(\varphi, \theta)=R_{0}(\theta)+R_{2}(\theta) \cos \left(2\left(\varphi-\varphi_{0}\right)\right)+} \\ {R_{4}(\theta) \cos \left(4\left(\varphi-\varphi_{0}\right)\right)}\end{array}. $

(22)

其中φ₀为观测系统坐标系与推导方程所用坐标系的夹角。将余弦函数展开后可以得到式(19)的形式, 而各级数振幅项R_n(θ)和偏移相位φ₀与合并前的傅里叶级数式(16)的关系为

$ \left\{ {\begin{array}{*{20}{l}} {{R_n}(\theta ) = \sqrt {u_n^2(\theta ) + v_n^2(\theta )} ,}\\ {{\varphi _0} = \frac{1}{n}\arctan \left( {\frac{{{v_n}(\theta )}}{{{u_n}(\theta )}}} \right).} \end{array}} \right. $

(23)

最后推导得到的正交裂缝方位反射系数傅里叶各阶级数项的具体表达式为

$ \left\{ {\begin{array}{*{20}{l}} {{R_0}(\theta ) = {A_0} + {B_0}{{\sin }^2}\theta + {C_0}{{\sin }^2}\theta {{\tan }^2}\theta ,}\\ {{R_2}(\theta ) = \frac{1}{2}g\left[ {\Delta \left( {{\delta _{{\rm{T}}1}} - {\delta _{{\rm{T}}2}}} \right) - r\Delta \left( {{\delta _{{\rm{N}}1}} - {\delta _{{\rm{N}}2}}} \right)} \right]{{\sin }^2}\theta + }\\ {\frac{1}{2}g(g - 1)\Delta \left( {{\delta _{{\rm{NI}}}} - {\delta _{{\rm{N}}2}}} \right){{\sin }^2}\theta {{\tan }^2}\theta ,}\\ {{R_4}(\theta ) = \frac{1}{8}g\left[ {\Delta \left( {{\delta _{{\rm{T}}1}} + {\delta _{{\rm{T}}2}}} \right) - g\Delta \left( {{\delta _{{\rm{NI}}}} + {\delta _{{\rm{N}}2}}} \right)} \right]{{\sin }^2}\theta {{\tan }^2}\theta .} \end{array}} \right. $

(24)

其中参数A₀、B₀、C₀与各向同性参数和各向异性参数之间的关系如下:

$ \left\{ {\begin{array}{*{20}{l}} {{A_0} = {A_{{\rm{iso}}}} - \frac{1}{4}{r^2}\Delta \left( {{\delta _{{\rm{NI}}}} + {\delta _{{\rm{N}}2}}} \right),}\\ {{B_0} = {B_{{\rm{iso}}}} + \frac{1}{2}g\Delta \left( {{\delta _{{\rm{T}}1}} + {\delta _{{\rm{T}}2}}} \right) - \frac{1}{4}r\Delta \left( {{\delta _{{\rm{N}}1}} + {\delta _{{\rm{N}}2}}} \right),}\\ {{C_0} = {C_{{\rm{iso}}}} - \frac{1}{8}g\Delta \left( {{\delta _{{\rm{TI}}}} + {\delta _{{\rm{T}}2}}} \right) - \frac{1}{8}\left( {3{g^2} - 4g + 2} \right)\Delta \left( {{\delta _{{\rm{N}}1}} + {\delta _{{\rm{N}}2}}} \right).} \end{array}} \right. $

(25)

当各向同性介质中只发育一组垂直裂缝时, 即裂缝弱度Δδ_N1、Δδ_T1等于零, 或者Δδ_N2、Δδ_T2等于零, 各阶傅里叶反射系数的表达式变为

$ \left\{ \begin{array}{l} {R_0}(\theta ) = {A_{{\rm{iso}}}} - \frac{1}{4}{r^2}\Delta {\delta _{{\rm{N}}1}} + \left( {{B_{{\rm{iso}}}} + \frac{1}{2}g\Delta {\delta _{{\rm{T}}1}} - \frac{1}{4}r\Delta {\delta _{{\rm{N}}1}}} \right){\sin ^2}\theta + \\ \left( {{C_{{\rm{iso}}}} - \frac{1}{8}g\Delta {\delta _{{\rm{T}}1}} - \frac{1}{8}\left( {3{g^2} - 4g + 2} \right)\Delta {\delta _{{\rm{N}}1}}} \right){\sin ^2}\theta {\tan ^2}\theta ,\\ {R_2}(\theta ) = \frac{1}{2}g\left( {\Delta {\delta _{{\rm{TI}}}} - r\Delta {\delta _{{\rm{NI}}}}} \right){\sin ^2}\theta + \frac{1}{2}g(g - \\ 1)\Delta {\delta _{{\rm{NI}}}}{\sin ^2}\theta {\tan ^2}\theta ,\\ {R_4}(\theta ) = \frac{1}{8}g\left( {\Delta {\delta _{{\rm{T}}1}} - g\Delta {\delta _{{\rm{N}}1}}} \right){\sin ^2}\theta {\tan ^2}\theta . \end{array} \right. $

(26)

式(26)即为HTI介质傅里叶反射级数表达式, 和Downtown^[9]给出的形式相同。可见HTI介质傅里叶反射系数公式是正交介质傅里叶反射系数公式的一种特殊情况。

4 正交裂缝密度与傅里叶反射级数项之间的关系

Hudson^[3-4]给出了含有薄硬币形状裂缝的等效模量刚度表达式, 对比线性滑移模型, 可以得到裂缝法向弱度δ_N和切向弱度δ_T与Hudson参数之间的关系:

$ \delta_{\mathrm{N}}=\frac{\lambda_{\mathrm{b}}+2 \mu_{\mathrm{b}}}{\mu_{\mathrm{b}}} U_{11} e, \delta_{\mathrm{T}}=U_{33} e. $

(27)

其中e为裂缝密度, U₁₁、U₃₃为Hudson裂缝参数, 将Hudson裂缝参数带入式(27)化简可以得到^[12]:

$ {\delta _{\rm{N}}} = \frac{{4e}}{{3g(1 - g)\left[ {1 + \frac{1}{{\pi g(1 - g)}}\frac{{{K_{\rm{f}}} + (4/3){\mu _{\rm{f}}}}}{{{\mu _{\rm{b}}}}}\frac{a}{c}} \right]}}, $

(28)

$ {\delta _{\rm{T}}} = \frac{{16e}}{{3(3 - 2g)\left[ {1 + \frac{4}{{{\rm{ \mathsf{ π} }}(3 - 2g)}}\frac{{{\mu _{\rm{f}}}}}{{{\mu _{\rm{b}}}}}\frac{a}{c}} \right]}}. $

(29)

式中, K_f和μ_f分别为裂缝中流体的体积模量和剪切模量; a/c为裂缝高宽比。如果裂缝不含流体或饱含气体, K_f、μ_f近似等于零, 此时, 裂缝的法向和切向弱度为

$ \delta_{\mathrm{Nd}}=\frac{4 e}{3 g(1-g)}, \delta_{\mathrm{Td}}=\frac{16 e}{3(3-2 g)}. $

(30)

对应的反射系数傅里叶级数R₂(θ)和R₄(θ)近似等于:

$ R_{2}(\theta) \approx \frac{2}{3}\left(\frac{4 g}{3-2 g}-\frac{1-2 g}{1-g}\right)\left(\Delta e_{1}-\Delta e_{2}\right) \sin ^{2} \theta, $

(31)

$ R_{4}(\theta) \approx \frac{1}{6} g\left(\frac{4}{3-2 g}-\frac{1}{1-g}\right)\left(\Delta e_{1}+\Delta e_{2}\right) \sin ^{2} \theta \tan ^{2} \theta. $

(32)

可以看出, 对于给定入射角度, 通过联立(31)和(32)两式, 就可以分别得到两个正交方向裂缝的发育密度。

对于饱含流体的薄裂缝, 其宽高比c/a很小, μ_f等于零, 裂缝法向弱度趋于零:

$ \delta_{\mathrm{Ns}}=0, \delta_{\mathrm{Ts}}=\frac{16 e}{3(3-2 g)}. $

(33)

对应的反射系数傅里叶级数R₂(θ)和R₄(θ)近似等于:

$ R_{2}(\theta) \approx \frac{8 g}{3(3-2 g)}\left(\Delta e_{1}-\Delta e_{2}\right) \sin ^{2} \theta, $

(34)

$ R_{4}(\theta) \approx \frac{2 g}{3(3-2 g)}\left(\Delta e_{1}+\Delta e_{2}\right) \sin ^{2} \theta \tan ^{2} \theta. $

(35)

同样, 对于给定入射角度, 通过联立式(34)和(35), 就可以分别得到两个正交方向裂缝的发育密度。

5 模型测试 5.1 单点模型测试

建立表 1所示的模型, 上层为各向同性介质, 纵、横波速度、密度分别为3 800 m/s、1 900 m/s、2 450 kg/m³; 下层为裂缝介质, 只发育一组平行裂缝, 裂缝饱含气体, 裂缝的法向弱度δ_N1为0.09, 切向弱度δ_T1为0.17, 裂缝密度e₁为0.08, 纵、横波速度、密度分别为4 200 m/s、2 100 m/s、2 550 kg/m³。其他参数为零。

按表 1给出的HTI模型参数, 利用正交介质傅里叶反射系数展开式(22), 得到方位各向异性反射系数随入射角度和方位变化响应如图 1所示。可以看出随入射角度增大, 反射系数振幅变化越来越明显, 即方位各向异性特征越来越明显; 随方位角变化, 反射系数在裂缝法向0°取得极大值, 而在裂缝走向±90°取得极小值。符合HTI介质的方位各向异性特征。

进一步研究正交裂缝, 建立如表 2所示的模型, 上层为各向同性介质, 纵横波速度、密度分别为3 800 m/s、1 900 m/s、2 450 kg/m³; 下层为正交裂缝介质, 纵横波速度、密度分别为4 200 m/s、2 100 m/s、2 550 kg/m³。其中第一组裂缝, 裂缝法向平行于30°方向, 裂缝法向弱度δ_N1为0.09, 切向弱度δ_T1为0.17, 裂缝密度e₁为0.08;而第二组裂缝垂直于上一组的裂缝, 即裂缝法向平行于120°方向, 裂缝法向弱度δ_N2为0.05, 切向弱度δ_T2为0.1, 裂缝密度e₂为0.05。由此可以得出正交裂缝偏移相位30°, 正交裂缝方位反射系数随方位角和入射角变化如图 2所示。可以看出, 随入射角度增大, 反射系数振幅变化越来越明显, 方位各向异性明显; 随方位角变化, 同一入射角度情况下, 30°方向为反射系数取得最大值方向, 对应第一组裂缝密度较大的裂缝组法向; -60°方向取得局部极大值, 对应第二组密度较小的裂缝组法向。图 3是正交裂缝方位反射系数傅里叶二阶级数R₂项, 可以看出极大值对应的方位角正是正交裂缝偏移相位角。

5.2 单道模型裂缝密度预测

选取Marmous模型中的一道作为模型数据, 并给出正交裂缝密度, 假设裂缝为含气或干裂缝, 裂缝位置为裂缝密度不为零的地方。模型参数如图 4所示。然后根据褶积模型, 用雷克子波合成不同方位角和入射角的具有不同信噪比的地震记录, 进一步利用傅里叶级数展开原理, 分别得到不同入射角度的傅里叶2、4阶级数, 最后通过联立式(31)和(32)解得正交方向的裂缝密度。

图 5为应用傅里叶级数方法对不同信噪比地震记录进行裂缝密度e₁预测的结果。可以看出当信噪比大于5时, 用该方法预测得到的结果准确, 信噪比过小时, 预测结果与实际裂缝密度差别较大, 因此实际应用中应采用信噪比大于5的地震数据。

6 结束语

通过推导得出弱各向异性条件下正交裂缝方位反射系数的傅里叶级数展开方程, 并根据Hudson裂缝模型, 分别得到了干裂缝和流体饱含裂缝的傅里叶2、4阶级数与裂缝密度的关系。模型测试结果表明, 该公式能够正确描述正交裂缝介质方位各向异性与反射系数之间的关系, 而且应用该方法进行裂缝密度预测时, 在较高信噪比条件下能够准确预测不同方向的裂缝密度。需要说明的是, 该方程应用条件为弱各向异性介质, 且方位采样数据应足够多, 保证傅里叶级数展开方法要求的采样间隔满足Nyquist采样定理; 另外本文中只是进行了模型测试, 实际数据应用中需要对纵横波速度比的影响做进一步研究。