2017, Vol. 60
随钻声波测井在近年来发展迅速，已经成为了地球物理勘探和油藏工程应用中的重要技术.四极子波在随钻条件下测量横波速度时具有较好的性质(苏远大, 2006; Tang et al., 2002; Dubinsky et al., 2003)，被广泛用作地层横波速度测量.然而，作为在钻铤存在条件下沿钻井传播的导波，四极子波具有很强的频散效应：四极子声波仅在低频时才以地层横波速度传播，其波速随频率增加大幅下降；在泥浆密度和仪器尺寸较大的情况下，频散效应更为严重(Tang et al., 2003).其后果是，对随钻四极声波测井数据用常规的时间慢度相关法(STC，为一非频散常规处理方法)(Kimball and Marzetta, 1984)处理，得到的慢度大于实际地层慢度(Geerits and Tang, 2003, Kimball, 1998).
为了从随钻四极声波数据中得到准确的地层横波慢度，国内外学者开展了大量的相关研究.其中应用比较普遍的是基于理论模型的频散校正方法(Tang et al., 2006; Zhang et al., 2010; Lee et al., 2016; Su et al., 2016)，这类方法需要假设均匀地层模型计算理论频散曲线，计算时需要输入较多的模型参数(如地层密度、速度及井径等)，当实际地层情况严重偏离理论模型时，难以得到理想结果.另一类方法是数据驱动法.根据导波频散曲线低频截止频率处的慢度为地层慢度这一基本事实(唐晓明和郑传汉，2004)，Huang和Yin(2005)提出了一种数据驱动的频散校正方法.该方法将频率-慢度域中的二维频散数据向慢度轴投影，得到与慢度相关的统计直方图，再通过确定该直方图的边缘的方法来估算实际地层慢度.然而，当该方法应用于随钻数据时，噪声的严重干扰会导致该方法的失效.为弥补噪声干扰的影响，Tang等(2010)提出了一种对声波频散数据进行曲线拟合的方法.通过构造解析函数代替理论频散曲线，降低了个别频率点上频散数据出现奇异值的影响.由于拟合所用的解析函数是从实际频散数据来构造的，该方法也是一种数据驱动法.
本文采用的数据驱动频散处理方法是将解析函数拟合的思想用于慢度直方图的拟合.首先通过慢度-频率相干法计算实际频散曲线，然后计算慢度概率密度函数(SPD, 相当于上述的慢度直方图)，根据概率密度函数形态特征采用与之相适应的拟合函数，再利用统计特征值匹配的方法确定该函数的参数值，最后对该拟合函数寻边来确定地层慢度.对该频散处理方法用有限差分模拟的随钻四极子声波数据进行了验证，并结合现场随钻测井资料的进行了处理.处理结果表明该方法快速有效，能较为准确地确定地层慢度.
2 方法概述本文所述频散波处理方法需要首先从实际测井资料中提取频率-慢度频散数据(曹正良等, 2005).假设有N个接收器的阵列声波测井数据中的导波为Xn(t), (n=1, 2, …, N)，其傅里叶变换后的频谱为Xn(ω), (n=1, 2, …, N).此阵列声波频谱的慢度-频率相关函数为(唐晓明和郑传汉, 2004; Lee et al., 2016)
|
(1) |
其中*为复共轭，接收器间距为d，ω为角频率.对每一频率ω求取相关函数最大值对应的慢度s，便得到实际数据的频散数据Sd(ω).
通过该频散数据可以计算慢度概率密度函数SPD(Slowness Probability Density Function, Tang et al., 2010)，定义如下：
|
(2) |
其中ωH和ωL分别为频散数据Sd(ω)的频率上限与下限，σs为Sd(ω)的均方差开根.与传统的统计直方图(Huang and Yin, 2005)相比较，上述的SPD为一连续函数，可以减少频散数据中异常值的影响.随钻四极声波的频散特征是其慢度在低频时趋于地层横波慢度，即该导波慢度在地层横波慢度附近取值的概率密度增加，使得SPD函数呈现一明显的峰值.因此，确定该峰值的边缘即可提取地层的慢度.然而在随钻环境下，测井数据中的强噪声干扰有时会严重影响SPD函数的形态，使其边缘的确定十分困难.本文采用一种基于函数拟合的计算方法，通过对拟合函数的解析计算来代替对SPD函数直接寻边.
对SPD数据的拟合函数的考虑是基于上述的导波慢度随频率的单调变化特征，这一特征使得SPD函数在峰值边缘处较为陡峭，根据这种形态特征提出一种相适应的拟合函数:
|
(3) |
其中s为慢度变量，该拟合参数共包括四个参数：A、α、β、s0，A确定了拟合函数的幅度值，α、β确定了拟合函数的形态，s0确定了拟合函数的起始值.该函数有一明显的峰值，且其低慢度一侧甚为陡峭.虽然其它函数或函数组合也能满足这种形态，但(3)式的解析性质可以使以下的统计特征值匹配计算变得十分简便.用上述解析函数直接拟合SPD数据，可以确定(3)式中的四个参数.然而这是一个计算量较大的四参数的反演问题，不适合测井数据的快速有效处理.我们采用以下统计特征值匹配的方法来确定这些参数.
根对(2)、(3)式所示的概率密度及其拟合函数，计算二者的期望、均值和方差(蔡德勒等, 2012)等统计特征值I0、sc、σc2，再令计算的二函数的各统计特征值分别相等，即匹配二者相应的特征值，得到
|
(4) |
上式中右端是根据解析函数(3)式求得的各统计特征值，其中Γ为伽马函数；Ω为慢度区间，可以根据实际数据并参考STC处理结果设置合理区间.令这些特征值分别等于SPD数据计算出的相应特征值，便得到确定拟合函数f(s)中四个参数的三个方程，仍需的第四个方程由调节f(s)中慢度起始值s0得到.
为了确定拟合函数的起始慢度s0，以拟合均方差构建如下的目标函数：
|
(5) |
其中f(s; s0)为式(3)中的拟合函数，现以s0为待求参量.拟合函数的参数A、α、β已经通过特征值匹配方法确定，上式就变成了一维函数极小值优化问题，当E(s)最小时，拟合函数与SPD函数拟合最好，此时对应的s0就是拟合函数的起始慢度.
拟合函数确定以后，其寻边问题由函数的解析性质很容易解决.SPD函数在其边缘处增加最快，所以其边缘对应于拟合函数的一阶导函数极值.由(3)式求取该极值处的慢度值，最终得到实际地层的慢度值se:
|
(6) |
上述方法通过求解声波频散数据的慢度概率密度函数的寻边问题，增强了声波数据处理的抗噪性和稳定性，能够快速准确的提取实际地层慢度，适合于现场随钻声波测井资料的快速处理.
3 方法验证为验证本文方法的正确性和有效性，利用有限差分方法模拟在慢地层中传播的随钻四极子波波形(Liu et al., 1996; Wang and Tang, 2003).随钻声波测井仪器模型的剖面见图 1，其中包括声源、接收阵列、以及放置在它们中间位置的刻槽隔声体，模拟参数见表 1.模拟所得阵列波形如图 2a所示，图中突出的振相便是四极声波，其慢度接近于地层横波慢度.利用常规的慢度-时间相干方法(STC)对该四极子波形数据进行处理，处理结果如图 2b的相干图所示，据此获得地层横波慢度为1290 μs·m-1，由于STC方法并没有考虑四极子波的频散效应，获得的慢度值比实际地层横波慢度大不少，所以需要对该波形数据用频散方法处理.
|
图 1 随钻声波测井仪器模型 Fig. 1 The LWD acoustic tool model |
|
表 1 计算随钻多极子波的地层及钻铤参数 Table 1 The parameters of drill collar and formation for calculating the LWD multipole waves |
|
图 2 模拟数据及其处理结果 (a)四极子波形; (b)四极子波的STC处理结果; (c)四极子波的频散曲线; (d) SPD函数(散点)及其拟合函数(红色). Fig. 2 Synthetic LWD data and processing result (a) Synthetic LWD quadrupole waveform data; (b) The STC correlogram of the data; (c) Frequency dispersion of the wave data; (d) SPD (markers) of the dispersion data and its fitting function (red curve). |
下面利用本文所述方法对模拟波形数据进行处理.该方法主要利用四极子波的低频信息，所以在提取频散数据时，设置频率区间为1~3 kHz，慢度区间为1150~1350 μs·m-1，所得频散数据(散点)如图 2c所示.可以看出，实际频散数据显示出明显的频散效应，慢度随频率增加，低频段平缓趋于地层横波慢度.根据式(2)计算SPD函数(散点)，如图 2d所示.由于低频数据较为平坦，在数据在趋于地层慢度时概率密度最大，显示出突出的峰值.利用式(3)所示的拟合函数采用上述的方法对该SPD函数进行拟合，拟合函数(红色)如图 2d所示，可以看出该函数很好地拟合了SPD数据的低慢度部分.通过式(6)计算获得地层横波慢度为1215.5 μs·m-1，与模拟参数符合很好，表明该方法是准确有效的.
4 实际数据处理 4.1 中速地层处理实例利用本文方法对一中速地层的随钻四极声波数据进行处理，井段长度为350 m，处理结果如图 3所示.图中第一道为岩性曲线，给出了自然伽马、井径曲线.第二道为随钻四极阵列数据第一接收器的变密度波形图.处理结果如图中第三道所示，其中给出了提取的横波慢度(红线)以及拟合函数(变密度图)，拟合函数可以作为一种质量监控工具，其边缘对应于实际地层慢度.可以看出，提取的慢度很好地追踪了拟合函数的边缘.图中第四道分别给出了非频散STC方法计算的四极横波慢度(蓝线)和用本文方法提取的横波慢度(红线).由于随钻四极子波的频散，频散与非频散两条慢度曲线的差异就反映了波的频散效应，该差异即为常规STC提取的地层横波慢度所需的频散校正量，如第四道中灰色填充曲线所示.该层段校正量约为7%~8%，说明该频散效应影响甚大不容忽略.作为对比，图中同时给出了单极子横波慢度曲线(黑线)，可以看出频散处理结果与单极横波慢度符合很好.由单极提取的横波慢度一般不受频散的影响，反映地层的横波慢度.这种符合说明本文的频散处理方法能从频散很强的随钻四极波形数据中有效地获取地层横波慢度.
|
图 3 中速地层随钻声波数据处理实例 Fig. 3 Example of LWD data processing for a moderately fast formation |
在慢速地层中，单极声波测井无法得到地层横波慢度，四极声波测井成为地层横波慢度随钻测量的主要技术，随钻四极子波的频散处理方法因此变得十分重要.利用本文方法对一500 m井段的超慢地层的随钻四极声波数据进行处理，处理结果如图 4a所示.从第二道波形图中可以看到四极子波的到时随着深度的增加逐渐提前，四极子波受到的频散效应影响逐渐减弱.第三道为频散处理结果，所提取的地层横波慢度可以很好地追踪拟合函数的边缘.第四道分别给出了常规STC处理结果(蓝线)和频散校正结果(红线)，二者之间的差异为图中的灰色填充曲线所示.在浅部未固结的超慢地层中，四极子波频散效应较强，频散校正量约为10%左右.随着深度的增加，地层横波慢度逐渐减小，四极子波频散效应逐渐减弱，在深部快地层处频散校正量降至2%~3%左右.
|
图 4 超慢地层随钻声波数据实例 Fig. 4 Example of LWD data processing for an unconsolidated slow formation |
为进一步验证频散处理结果的正确性，图 4a第四道还给出根据图 1的理论模型所做频散校正(Su et al., 2016)得到的横波慢度(黑线).在这个实例中，基于理论模型的频散处理与数据驱动的频散处理所得结果十分一致.这两种结果的一致性，除了证明随钻四极声波理论(Tang et al., 2002)与实际数据吻合之处，还说明了本文处理方法的正确性.图 4b还分别给出了深层(X827m)、中层(X559m)、浅层(X401m)的三个有代表性的深度点的实际频散数据.图 4c给出了相应的SPD(散点)及其拟合函数(红线).频散数据的低频部分跳动较大，导致SPD的边缘较难确定.但是本文的拟合函数的方法仍能对SPD进行有效寻边，从而获得较为准确的地层横波慢度(三深度点的地层横波慢度由图 4b中的箭头给出).
5 结论为了解决随钻四极声波测量地层横波慢度中的频散问题，本文采用了并实现了一种新的数据驱动的频散处理方法.该方法将四极子声波低频慢度的提取转化为波的慢度概率密度函数的寻边问题，并给出了一种曲线拟合的快速寻边方法.即首先利用慢度-频率相干法对四极子波波形数据进行处理获得其频散数据，通过该频散数据计算SPD函数，最后对该SPD函数进行寻边获得实际地层慢度.对实际随钻四极声波数据的处理结果表明，该方法可以准确提取地层横波慢度，且不需要复杂的理论模型计算，为随钻四极子横波测井提供了一种快速有效的处理方法.
Cao Z L, Wang K X, Xie R H, et al.
2005. A comparative study on application of three dispersion analysis to array sonic logs. Chinese Journal of Geophysics, 48(6): 1449-1459.
|
|
Dubinsky V, Tang X M, Bolshakov A, et al. 2003. Engineering aspects of LWD quadrupole measurements and field test results.//SPE Annual Technical Conference and Exhibition. Denver, Colorado:Society of Petroleum Engineers, doi:10.2118/84248-MS.
|
|
Geerits T, Tang X M.
2003. Centroid phase slowness as a tool for dispersion correction of dipole acoustic logging data. Geophysics, 68(1): 101-107.
DOI:10.1190/1.1543197 |
|
Huang X J, Yin H Z. 2005. A data-driven approach to extract shear and compressional slowness from dispersive waveform data.//SEG Technical Program Expanded Abstracts 2005. Society of Exploration Geophysicists, 384-387.
|
|
Kimball C V, Marzetta T L.
1984. Semblance processing of borehole acoustic array data. Geophysics, 49(3): 274-281.
DOI:10.1190/1.1441659 |
|
Kimball C V.
1998. Shear slowness measurement by dispersive processing of the borehole flexural mode. Geophysics, 63(2): 337-344.
DOI:10.1190/1.1444333 |
|
Lee S Q, Tang X M, Su Y D, et al.
2016. Model-based dispersive processing of borehole dipole wave data using an equivalent-tool theory. Geophysics, 81(1): D35-D43.
DOI:10.1190/geo2015-0210.1 |
|
Liu Q H, Schoen E, Daube F, et al.
1996. A three-dimensional finite difference simulation of sonic logging. Journal of the Acoustical Society of America, 100(1): 72-79.
DOI:10.1121/1.415869 |
|
Su Y D, Sun J M, Fan Y R, et al.
2006. On Numerical study of the LWD multipole source propagation characteristics in slow formation. Well Logging Technology, 30(3): 205-207.
DOI:10.3969/j.issn.1004-1338.2006.03.005 |
|
Su Y D, Jiang C, Zhuang C X, et al.
2016. Joint inversion of logging-while-drilling multipole acoustic data to determine P-and S-wave velocities in unconsolidated slow formations. Geophysics, 81(5): D553-D560.
DOI:10.1190/geo2016-0133.1 |
|
Tang X M, Dubinsky V, Bolshakov A, et al. 2002. Shear-velocity measurement in the logging-while-drilling environment:Modeling and field evaluations.//Presented at the SPWLA 43rd Annual Logging Symposium. Oiso, Japan:Society of Petrophysicists and Well-Log Analysts.
|
|
Tang X M, Patterson D, Dubinsky V, et al. 2003. Logging-While-Drilling shear and compressional measurements in varying environments.//SPWLA 44th Annual Logging Symposium. Galveston, Texas:Society of Petrophysicists and Well-Log Analysts.
|
|
Tang X M, Dubinsky V, Patterson D J. 2006. Development of a low-frequency LWD quadrupole shear-wave technology to improve quality of formation shear velocity measurement.//Presented at theSPE Annual Technical Conference and Exhibition. San Antonio, Texas, USA:SPE.
|
|
Tang X M, Li C, Patterson D J.
2010. A curve-fitting technique for determining dispersion characteristics of guided elastic waves. Geophysics, 75(3): E153-E160.
DOI:10.1190/1.3420736 |
|
Wang T, Tang X M.
2003. Finite-difference modeling of elastic wave propagation:A nonsplitting perfectly matched layer approach. Geophysics, 68(5): 1749-1755.
DOI:10.1190/1.1620648 |
|
Zhang Z X, Mochida M, Kubota M, et al. 2010. Shear slowness estimation by inversion of LWD borehole quadrupole mode.//80th Annual International Meeting, SEG, Expanded Abstracts. SEG:528-532.
|
|
埃伯哈德·蔡德勒等. 2012. 数学指南: 实用数学手册. 李文林译. 北京. 科学出版社.
|
|
曹正良, 王克协, 谢荣华, 等.
2005. 三种阵列声波测井数据频散分析方法的应用与比较. 地球物理学报, 48(6): 1449–1459.
|
|
唐晓明, 郑传汉. 2004. 定量测井声学. 赵晓敏译. 北京: 石油工业出版社.
|
|
苏远大, 孙建孟, 范宜仁, 等.
2006. 基于软地层随钻多极子声源传播特性的数值研究. 测井技术, 30(3): 205–207.
DOI:10.3969/j.issn.1004-1338.2006.03.005 |
|