一种基于频变属性的储层渗透率评估方法.pdf

本发明提供一种基于频变属性的储层渗透率评估方法，应用数值模拟分析了K1属性对渗透率变化的响应特征,结果表明该属性对于渗透率的变化有着明显的反映。且在不同顶底波阻抗差异下会出现不同的响应特征。较于流度属性只考虑优势低频，K1属性将高频、低频地震资料结合使用。在测井分析和实例应用中，K1属性能够和测井资料吻合并且有效地反映储层流体的流动能力，应用数值模拟分析K1属性对渗透率变化的响应特征，及大地提高了对储层渗透率评价的准确性。

1.一种基于频变属性的储层渗透率评估方法，其特征在于：针对不同地区的地质情况，
选择适当的岩石物理模型进行数值模拟，分析该地区K1属性对渗透率的响应特征；计算叠
后地震资料的时频谱；从时频信号中计算出K1属性；结合数值模拟结果和测井数据评价储
层的渗透率分布。
2.根据权利要求1所述的基于频变属性的储层渗透率评估方法，其特征在于，具体步骤
如下：
(1)K1属性的提取
K1属性是地震频谱中高频信号与低频信号的比值：
$K_1={\mathrm{\Σ}}_{{\mathrm{\ω}}_m}^{{\mathrm{\ω}}_2}\left|S\left(\mathrm{\ω}\right)\right|/{\mathrm{\Σ}}_{{\mathrm{\ω}}_1}^{{\mathrm{\ω}}_m}\left|S\left(\mathrm{\ω}\right)\right|---\left(1\right)$
式中S(ω)为地震频谱，ω1、ωm和ω2是地震信号的低频频率、峰值频率和高频频率；
为了将时域的信号转换成时频信号，利用反褶积短时傅里叶变换提取地震资料的时频
信号,信号x(u)的短时傅立叶变换定义为:
${\mathrm{STFT}}_X(t,f)={\∫}_{-\mathrm{\∞}}^{+\mathrm{\∞}}x\left(u\right)h*(u-t)e^{(-j2\mathrm{\π}fu)}du---\left(8\right)$
其中h(u-t)是窗函数，常用的是高斯窗，其中*是共轭转置；
短时傅立叶变换谱定义为：
${\mathrm{SPEC}}_X(t,f)=|{\∫}_{-\mathrm{\∞}}^{+\mathrm{\∞}}x\left(u\right)h*(u-t)e^{(-j2\mathrm{\π}fu)}du{|}^2---\left(9\right)$
对于一个信号x(t),它的Wigner-Ville分布定义为：
${\mathrm{WVD}}_X(t,f)={\∫}_{-\mathrm{\∞}}^{+\mathrm{\∞}}x(t+\frac{\mathrm{\τ}}{2})x*(t-\frac{\mathrm{\τ}}{2})e^{(-j2\mathrm{\π}f\mathrm{\τ})}d\mathrm{\τ}---\left(10\right)$
Wigner-Ville分布有交叉项，它是由Wigner-Ville分布非线性引起的；
广义S变换的表达式为：
$GST(\mathrm{\τ},f)={\∫}_{-\mathrm{\∞}}^{+\mathrm{\∞}}x\left(t\right)\frac{\mathrm{\λ}\left|f\right|p}{\sqrt{2\mathrm{\π}}}{e}^{\frac{-{\mathrm{\λ}}^2{f}^{2p}{(\mathrm{\τ}-t)}^2}{2}}{e}^{-j2\mathrm{\π}ft}dt,\mathrm{\λ}\>0,p\>0---\left(11\right)$
广义S变换通过引入两个参数λ和p，改造了S变换的高斯
窗函数，灵活地调节高斯窗函数随频率尺度f的变化趋势；
信号x(t)的短时傅立叶变换谱写成以下二维褶积的形式：
$\begin{array}{c}{\mathrm{SPEC}}_X\left(t,f\right)=|{\mathrm{STFT}}_X\left(t,f\right){|}^2\\ ={\∫}_{-\mathrm{\∞}}^{+\mathrm{\∞}}{\∫}_{-\mathrm{\∞}}^{+\mathrm{\∞}}{\mathrm{WVD}}_h\left(u,v\right)\×{\mathrm{WVD}}_x\left(t-u,f-v\right)dudv\end{array}---\left(12\right)$
其中WVDx和WVDh分别为信号x(t)和窗函数h(u)的Wigner-Ville分布，短时傅立叶变换
谱的交叉项在大部分情况为0；
将地震信号带入上式中，得到了时频地震信号，然后将时频信号带入公式(1)中，求得
K1属性。
(2)储层渗透率评价
通过数值模拟的结果了解该储层的K1属性总体趋势，渗透率随K1属性趋势变化而变
化。

一种基于频变属性的储层渗透率评估方法

技术领域

本发明涉及地质技术领域，尤其涉及一种基于频变属性的储层渗透率评估方法。

背景技术

对于石油工程师来说，渗透率无疑是一项必须加以重点关注的地层参数。它是确
定一口井是否应当完井和投产的依据。同时，它也是钻进油层保护、完井射孔方案选择、最
佳排液位置和生产速率以及三次采油措施的制定的基础。

受现有地震勘探技术的影响，难以从地震资料中直接反演渗透率。目前流度属性
已经被证实能够间接反映储层渗透率，但是受储层厚度和优势频率的不确定性影响，用其
评价储层渗透率存在不确定性和多解性。

发明内容

本发明的目的在于解决上述现有技术存在的缺陷，提供一种基于频变属性的储层
渗透率评估方法，应用数值模拟分析K1属性对渗透率变化的响应特征，及大地提高了对储
层渗透率评价的准确性。

为实现上述目的本发明采用如下方案：

一种基于频变属性的储层渗透率评估方法，针对不同地区的地质情况，选择 适当
的岩石物理模型进行数值模拟，分析该地区K1属性对渗透率的响应特征；计算叠后地震资
料的时频谱；从时频信号中计算出K1属性；结合数值模拟结果和测井数据评价储层的渗透
率分布。

进一步，具体步骤如下：

(1)K1属性的提取

K1属性是地震频谱中高频信号与低频信号的比值：

$K_1={\mathrm{\Σ}}_{{\mathrm{\ω}}_m}^{{\mathrm{\ω}}_2}\left|S\left(\mathrm{\ω}\right)\right|/{\mathrm{\Σ}}_{{\mathrm{\ω}}_1}^{{\mathrm{\ω}}_m}\left|S\left(\mathrm{\ω}\right)\right|---\left(1\right)$

式中S(ω)为地震频谱，ω1、ωm和ω2是地震信号的低频频率、峰值频率和高频频
率；

为了将时域的信号转换成时频信号，利用反褶积短时傅里叶变换提取地震资料的
时频信号,信号x(u)的短时傅立叶变换定义为:

${\mathrm{STFT}}_X(t,f)={\∫}_{-\mathrm{\∞}}^{+\mathrm{\∞}}x\left(u\right)h*(u-t)e^{(-j2\mathrm{\π}fu)}du---\left(8\right)$

其中h(u-t)是窗函数，常用的是高斯窗，其中*是共轭转置；

短时傅立叶变换谱定义为：

${\mathrm{SPEC}}_X(t,f)={|{\∫}_{-\mathrm{\∞}}^{+\mathrm{\∞}}x\left(u\right)h*(u-t)e^{(-j2\mathrm{\π}fu)}du|}^2---\left(9\right)$

对于一个信号x(t),它的Wigner-Ville分布定义为：

${\mathrm{WVD}}_X(t,f)={\∫}_{-\mathrm{\∞}}^{+\mathrm{\∞}}x(t+\frac{\mathrm{\τ}}{2})x*(t-\frac{\mathrm{\τ}}{2})e^{(-j2\mathrm{\π}f\mathrm{\τ})}d\mathrm{\τ}---\left(10\right)$

Wigner-Ville分布有交叉项，它是由Wigner-Ville分布非线性引起的；

广义S变换的表达式为：

$GST\left(\mathrm{\τ},f\right)={\∫}_{-\mathrm{\∞}}^{+\mathrm{\∞}}x\left(t\right)\frac{\mathrm{\λ}\left|f\right|p}{\sqrt{2\mathrm{\π}}}{e}^{-\frac{-{\mathrm{\λ}}^2{f}^{2\mathrm{\ρ}}{\left(\mathrm{\τ}-t\right)}^2}{2}}{e}^{-j2\mathrm{\π}ft}dt,\mathrm{\λ}>0,p>0---\left(11\right)$

广义S变换通过引入两个参数λ和p，改造了S变换的高斯窗函数，灵活地调节高斯
窗函数随频率尺度f的变化趋势；

信号x(t)的短时傅立叶变换谱写成以下二维褶积的形式：

$\begin{array}{c}{\mathrm{SPEC}}_X(t,f)={\left|{\mathrm{STFT}}_X(t,f)\right|}^2\\ ={\∫}_{-\mathrm{\∞}}^{+\mathrm{\∞}}{\∫}_{-\mathrm{\∞}}^{+\mathrm{\∞}}{\mathrm{WVD}}_h(u,v)\×{\mathrm{WVD}}_x(t-u,f-v)dudv\end{array}---\left(12\right)$

其中WVDx和WVDh分别为信号x(t)和窗函数h(u)的Wigner-Ville分布，短时傅立叶
变换谱的交叉项在大部分情况为0；

将地震信号带入上式中，得到了时频地震信号，然后将时频信号带入公式(1)中，
求得K1属性。

(2)储层渗透率评价

通过数值模拟的结果了解该储层的K1属性总体趋势，渗透率随K1属性趋势变化而
变化。

本发明基于双重孔隙介质模型，应用数值模拟分析了K1属性对渗透率变化的响应
特征,结果表明该属性对于渗透率的变化有着明显的反映。且在不同顶底波阻抗差异下会
出现不同的响应特征。较于流度属性只考虑优势低频，K1属性将高频、低频地震资料结合使
用。在测井分析和实例应用中，K1属性能够和测井资料吻合并且有效地反映储层流体的流
动能力，应用数值模拟分 析K1属性对渗透率变化的响应特征，及大地提高了对储层渗透率
评价的准确性。

附图说明

图1 K1属性随渗透率的变化曲线；

图2海上某区过井A地震剖面；

图3过井A的地震合成记录与层位标定；

图4井A岩心实验孔隙度与渗透率关系；

图5过井A的K1属性剖面；

图6顶层下延15ms的K1属性切片；

图7顶层下延15ms的流度属性切片。

具体实施方式

为使本发明的目的、技术方案和优点更加清楚，下面本发明中的技术方案进行清
楚、完整地描述，显然，所描述的实施例是本发明一部分实施例，而不是全部的实施例。基于
本发明中的实施例，本领域普通技术人员在没有作出创造性劳动前提下所获得的所有其他
实施例，都属于本发明保护的范围。

本发明提出了一套利用K1属性的储层渗透率评价方法，该属性从高频地震信号和
低频地震信号中提取出的。该方法的步骤为：首先针对不同地区的地质情况，选择适当的岩
石物理模型进行数值模拟，分析该地区K1属性对渗透率的响 应特征；计算叠后地震资料的
时频谱；从时频信号中计算出K1属性；结合数值模拟结果和测井数据评价储层的渗透率分
布，以下是实现过程：

步骤一：首先针对不同地区的地质情况，选择适当的岩石物理模型进行数值模拟，
分析该地区K1属性对渗透率的响应特征；

步骤二：计算叠后地震资料的时频谱；

步骤三：从时频信号中计算出K1属性；

步骤四：结合数值模拟结果和测井数据评价储层的渗透率分布。

以下是实现过程：

1 数值模拟

K1属性是地震频谱中高频信号与低频信号的比值：

$K_1={\mathrm{\Σ}}_{{\mathrm{\ω}}_m}^{{\mathrm{\ω}}_2}\left|S\left(\mathrm{\ω}\right)\right|/{\mathrm{\Σ}}_{{\mathrm{\ω}}_1}^{{\mathrm{\ω}}_m}\left|S\left(\mathrm{\ω}\right)\right|---\left(1\right)$

式中S(ω)为地震频谱，ω1、ωm和ω2是地震信号的低频频率、峰值频率和高频频
率。地震波在含流体的双孔介质中传播时，其主要的衰减是因为地震波驱动流体的移动。不
同频率的地震波在渗透率的影响下会出现不同的衰减，高频时地震波穿过双孔介质时，流
体来不及反映，所以高频时地震波在双孔介质中传播的衰减会比低频时地震波的衰减小许
多。

针对实际工区含油气储层为碳酸盐岩，为了分析该储层渗透率对K1属性的影响，
首先利用Kozlov的粗糙裂缝面的双重孔隙介质模型。根据储层与盖层的阻抗差异，将储层
分为两类，模型1是随渗透率值变大，盖层与储层的阻抗值差异变小；模型2是随渗透率的变
大，盖层与储层的阻抗值差异也变 大。根据计划，设计了两类模型的岩石骨架参数以及裂
缝参数和盖层与底层的速度与密度。假设岩石中的裂缝是定向排列的，利用Kozlov运用
Hertz的理论表示岩石中的裂缝面通过裂缝面上的凹凸点接触，并结合Thomsen流体影响因
子推导了由裂缝存在造成的岩石额外柔度表达式：

$Z_N=\frac{2}{3}Bce{\left(\frac{c}{K^{2}p}\right)}^{\frac{1}{3}}D---\left(2a\right)$

$Z_T=(1-0.5\mathrm{\ν})ce{\left(\frac{2c}{3(1-\mathrm{\ν}){\mathrm{\μ}}^{2}p}\right)}^{\frac{1}{3}}---\left(2b\right)$

$D={\[1+\frac{2}{3}\left(\frac{K_f}{\mathrm{\φ}K}\right)\left(\frac{K_S}{K_S-K_f}\right)BFce\sqrt[3]{\frac{cK}{3p}}\]}^{-1}---\left(2c\right)$

式中ZN和ZT分别是由于裂缝的存在造成的法向额外柔量和切向额外柔量。D是流体
因子，B＝[(1+ν)/(1-ν)]2/3，ν为泊松比，e为裂缝的密度，c是裂缝粗糙度，K是干燥岩石的体
积模量，P为有效压力，Kf为流体的体积模量，KS为岩石基质的体积模量，F是Hudson推导的流
体连接孔隙和裂隙的流体影响因子的表达式：

$F={\{1+\[\frac{3}{2b}(1-i)\]\[\frac{{\mathrm{\φK}}_f\mathrm{\κ}}{2{\mathrm{\η}}_f\mathrm{\ω}}\]\}}^{-1}---3)$

式中b为裂缝的平均开放性，φ为开孔的孔隙度，κ为渗透率，ω为角频率，ηf为流
体粘滞系数。

对于地震波的复速度计算式为式(4)，式中Z为干燥岩石的体变柔量，式(5)为逆品
质因子的表达式。

$\tilde{V}={\{\[(Z+\widetilde{Z_N})+\frac{4}{3}\mathrm{\μ}\]/\mathrm{\ρ}\}}^2---4)$

$Q^{-1}=\frac{\mathrm{Im}\left(\widetilde{V^2}\right)}{\mathrm{Re}\left(\widetilde{V^2}\right)}---5)$

两个模型的物性参数代入以上公式中，将有效压力设置为5Mpa得到了纵波相速度
带入Ursin所推导的在粘弹性层状各向同性介质传播的地震波反射系数计算方法，其方程
见式(6)，式中τ为储层时间厚度。V2、ρ2分别是储层的速度和密度，ρ1、ρ3和V3分别是盖层和底
层的密度和速度。

$R(\mathrm{\ω},\mathrm{\τ})=\frac{\widetilde{r_1}+\widetilde{r_2}\exp (-i\mathrm{\ω}\mathrm{\τ}+Q_2^{-1}/2)}{1+\widetilde{r_1}\widetilde{r_2}\exp (-i\mathrm{\ω}\mathrm{\τ}+Q_2^{-1}/2)}---6a)$

$\widetilde{r_1}=\frac{{\mathrm{\ρ}}_2\mathrm{Re}\left(\widetilde{V_2}\right)-{\mathrm{\ρ}}_1{V}_1(1+{\mathrm{iQ}}_2^{-1}/2)}{{\mathrm{\ρ}}_2\mathrm{Re}\left(\widetilde{V_2}\right)+{\mathrm{\ρ}}_1{V}_1(1-{\mathrm{iQ}}_2^{-1}/2)}---6b)$

$\widetilde{r_2}=\frac{{\mathrm{\ρ}}_3{V}_3(1-{\mathrm{iQ}}_2^{-1}/2)-{\mathrm{\ρ}}_2\mathrm{Re}\left(\widetilde{V_2}\right)}{{\mathrm{\ρ}}_2\mathrm{Re}\left(\widetilde{V_2}\right)+{\mathrm{\ρ}}_3{V}_3(1-{\mathrm{iQ}}_2^{-1}/2)}---6c)$

得到反射系数后，根据式(7)计算地震频谱。式中W(ω)为子波频谱，其主频为
30Hz。

S(ω,τ)＝R(ω,τ)W(ω) (7)

然后将S(ω,τ)带入式(1)中，计算出模型1和模型2的K1属性。

图1中显示的是K1随渗透率的变化曲线。当模型1层厚度为25ms时，K1随渗透率增大
而减小；当厚度为55ms和85ms时，K1随渗透率增大而增大。模型2层厚度为25ms时，K1随渗透
率增大而增大；厚度为55ms时，在渗透率为0.001至0.02mD段K1随渗透率增大而增大，在
0.02mD至1000mD时K1随渗透率增大而减小；层厚度为85ms，K1变化趋势与厚度为55ms时大致
一致。从以上分析可以看出K1属性在渗透率大于0.01mD且层厚度时间较高时，总体趋势在
模型1中随渗透率增加而增加，在模型2中随渗透率增加而减小。

2 K1属性的提取

为了将时域的信号转换成时频信号，利用反褶积短时傅里叶变换提

取地震资料的时频信号,信号x(u)的短时傅立叶变换定义为:

${\mathrm{STFT}}_X(t,f)={\∫}_{-\mathrm{\∞}}^{+\mathrm{\∞}}x\left(u\right)h*(u-t)e^{(-j2\mathrm{\π}fu)}du---\left(8\right)$

其中h(u-t)是窗函数，常用的是高斯窗，其中*是共轭转置。

短时傅立叶变换谱定义为：

${\mathrm{SPEC}}_X(t,f)={|{\∫}_{-\mathrm{\∞}}^{+\mathrm{\∞}}x\left(u\right)h*(u-t)e^{(-j2\mathrm{\π}fu)}du|}^2---\left(9\right)$

对于一个信号x(t),它的Wigner-Ville分布定义为：

${\mathrm{WVD}}_X(t,f)={\∫}_{-\mathrm{\∞}}^{+\mathrm{\∞}}x(t+\frac{\mathrm{\τ}}{2})x*(t-\frac{\mathrm{\τ}}{2})e^{(-j2\mathrm{\π}f\mathrm{\τ})}d\mathrm{\τ}---\left(10\right)$

Wigner-Ville分布有交叉项，它是由Wigner-Ville分布非线性引起的。广义S变换
的表达式为：

$GST(\mathrm{\τ},f)={\∫}_{-\mathrm{\∞}}^{+\mathrm{\∞}}x\left(t\right)\frac{\mathrm{\λ}\left|f\right|p}{\sqrt{2\mathrm{\π}}}{e}^{\frac{-{\mathrm{\λ}}^2{f}^{2p}{(\mathrm{\τ}-t)}^2}{2}}{e}^{-j2\mathrm{\π}ft}dt,\mathrm{\λ}>0,p>0---\left(11\right)$

广义S变换通过引入两个参数λ和p，改造了S变换的高斯窗函数，灵活地调节高斯
窗函数随频率尺度f的变化趋势。

信号x(t)的短时傅立叶变换谱可以写成以下二维褶积的形式：

$\begin{array}{c}{\mathrm{SPEC}}_X(t,f)={\left|{\mathrm{STFT}}_X(t,f)\right|}^2\\ ={\∫}_{-\mathrm{\∞}}^{+\mathrm{\∞}}{\∫}_{-\mathrm{\∞}}^{+\mathrm{\∞}}{\mathrm{WVD}}_h(u,v)\×{\mathrm{WVD}}_x(t-u,f-v)dudv\end{array}---\left(12\right)$

其中WVDx和WVDh分别为信号x(t)和窗函数h(u)的Wigner-Ville分布。短时傅立叶
变换谱的交叉项在大部分情况为0。

(5-15)式也可写成：

SX＝WX**Wh (13)

其中SX是信号x(u)的短时傅立叶变换谱，Wh是信号x(u)的Wigner-Ville分布，Wh是
窗函数h(u)的Wigner-Ville分布，**代表二维褶积。

我们希望反褶积结果Wx-有与Wigner-Ville分布Wx相近的时频分辨率，但是由于短
时傅立叶变换谱而减少了交叉项。

如果知道窗函数h(u)，就可以知道Wigner-Ville分布Wh。从短时傅立叶变换谱获
取Wx-是一个反褶积问题。将地震信号带入上式中，得到了时频地震信号，然后将时频信号带
入公式(1)中，求得K1属性。

3 储层渗透率评价

以实际的海上三维地震资料为例，图2为海上某区过井A地震剖面，绿色实线是目
的层，井A在碳酸盐岩储层中钻遇工业流油，图2中绿色箭头所指示段含85米的油层，亮蓝色
的测井曲线为深侧向电阻率，可以看出高渗高孔的含油气段的电阻率会明显升高，渗透率
可达到700mD，孔隙度达到31％。图3为过井A的地震合成记录与层位标定，第一条测井曲线
是声波阻抗，第二条测井曲线为深侧向电阻率。从图4(孔隙度与渗透率的交汇图)中可以看
出，该储层的孔渗关系复杂，虽然总体呈现孔隙度升高渗透率也升高的趋势，但是同一孔隙
度下的渗透率值也是多变的。

因为储层是碳酸盐岩，盖层是碎屑岩，所以该储层符合文中模型1的情况，即随渗
透率的变大，盖层与储层的波阻抗差异减小，通过数值模拟的结果可以了解该储层的K1属
性总体趋势应该随渗透率的增大而增大。

利用反褶积傅里叶变换计算出了该地震资料的时频谱，拾取了优势频段后带入公
式1中计算出K1属性。图5为过井A的K1属性剖面图，图中的暖色区域和测井资料中高渗高孔
的含油层相吻合。图6是K1属性在顶层下延15ms的切片(灰色标注是无地震数据区域)。切片
中的两个虚线区域的钻井在碳酸盐岩目的层都钻遇工业流油。黑色虚线是井A所在区域，该
井的日产油量可达783.1方；深蓝色虚线是井B所在区域，该井油层厚度达62米厚，渗透 率
达到700mD，孔隙度达到28％，日产油量为94.9方。通过测井资料可以了解两个井区的碳酸
盐岩储层流体活动能力强，渗透率大。可以看出在黑色虚线和蓝色虚线内是K1属性都和K1属
性与测井数据吻合，展示出高渗透率优质碳酸盐岩储层的分布。而在图7的流度属性切片中
蓝色虚线区域内没有展示出井B区内储层较强的流动能力

根据上述技术流程可以看出数值模拟的结果中显示出K1属性能反映渗透率的变
化，K1属性在模型1中总体随渗透率的变大而变大，在模型2中总体随渗透率的变大而变小。
因为K1属性同时利用了高、低频的地震信号，较流度属性更有优势，更能反映储层渗透率的
分布。

最后应说明的是：以上实施例仅用以说明本发明的技术方案，而非对其限制；尽管
参照前述实施例对本发明进行了详细的说明，本领域的普通技术人员应当理解：其依然可
以对前述各实施例所记载的技术方案进行修改，或者对其中部分技术特征进行等同替换；
而这些修改或者替换，并不使相应技术方案的本质脱离本发明各实施例技术方案的精神和
范围。