摘 要

　　　　　　　　　——试论从地震波形剖面到波阻抗剖面的转换摘　要：地震记录终极地质解释的目标包含两重含义：①能够展示地震记录中客观存在的全部地质信息；②所提

　　　　　　　　　——试论从地震波形剖面到波阻抗剖面的转换

摘　要：地震记录终极地质解释的目标包含两重含义：①能够展示地震记录中客观存在的全部地质信息；②所提供的各种信息的精度应该是客观上能够达到的最高水平。为了实现上述目标，探讨了从地震波形剖面到波阻抗剖面转换的思路和方法：从地震记录道中提取一个时变的反射子波，与地震记录道按褶积的逆过程作反褶积运算，可以求得反射波的振幅—时间序列；在把反射波振幅系列标定为反射系数系列的同时，推算出各反射层的波阻抗，从而把地震波形剖面转换为波阻抗剖面，以达到充分分离地震记录道中客观上所蕴含的各种地层地质信息的目的，同时获得更为精确的各类地震信息数据，用以拓展地震资料的解释应用领域。该文也可以被看做是编制相关转换软件的一个简要设计，对转换中涉及的重要参数进行了详细分析。经转换所获得的地震波阻抗剖面具有很高的实用价值，并最终服务于提升布井的针对性和钻井的成功率。

关键词：地震记录道  地震子波  反射系数  子波提取  地震波形剖面  波阻抗剖面  转换  勘探成功率

The ultimate geological interpretation from the seismic data：On the transformation from seismic wave section to acoustic impedance section

Abstract：The goal of ultimate geological interpretation from seismic data actually contains two meanillgs：First，complete and objective geological data would be obtained front the seismic record；second，the accuracy of all the provided data should be objectively up to the highest level．To this end，the idea and methods were discussed of the transformation from seismic wave section to acoustic impedance section.First，a time-varying reflected wavelet was extracted from the seismic recording channel．Then，the deconvoiution of seismic data was done according to the inverse process of convolution，thus an amplitude time series was obtained for the reflected wavelet．Meanwhile the amplitude series of a reflected wavelet was set to the reflection coefficient series，the acoustic impedance of each reflected layer was thus calculated．Therefore，the seismic wave section was transformed into acoustic impedance section，which helps completely strip all kinds of geological data of each layer objectively reflected from the seismic recording channel，meanwhile all types of seismic data were most exactly achieved，as a result，the domain of seismic data interpretation and application was thereby extended to a new area．On the other hand，this study was also regarded as a brief design of programming such a data-transformation software,the involved important parameters in which was also discussed in detail．The transformed acoustic impedance section is of great value in serving the targeted well layout and improving the drilling success ratio．

Key words：seismic recording channel，seismic wavelet，reflection coefficient，wavelet extraction，seismic wave section，acoustic impedance section，transformation，exploration Success ratio

从地震记录道中提取一个时变的反射子波，与地震记录道按褶积的逆过程作反褶积运算，可以求得反射波的振幅一时间序列。在把反射波振幅系列标定为反射系数系列的同时，推算出各反射层的波阻抗，从而把地震波形剖面转换为波阻抗剖面，以达到充分分离地震记录道中客观上所蕴含的各种地层地质信息的目的，同时获得更为精确的各类地震信息数据，用以拓展地震资料的解释应用领域。本文可以被看做是编制这一转换软件的一个简要设计。

1　地震记录道的基本组成元素

无干扰的地震记录道是地震子波与反射系数系列褶积的结果，即：

w(t)*R(t)＝s(t)　　     　　　　    (1)

式中w(t)表示地震子波；R(t)表示反射系数系列；s(t)表示地震记录道。

因此，地震记录中客观上所可能包含的地质信息，只会存在于地震子波或者反射系数系列这两个仅有的构造元素或载体中。

观察地震测井道(或VSP)初至直达波的变化，地震子波(地震单元波)在地下非均匀层状介质的扩散传播中，存在明显的变化。总的趋势是子波振幅包络变长，包络峰值滞后，包络内的振荡次数增加，振荡周期加长，即频率变低。有一个从所谓最小相位到混合相位再到最大相位的变化，本文参考文献[1]对形成这一变化的机制进行过初步的讨论。

由于我们不可能从地震记录中分离出每一个纯净的瞬时反射子波，以便从一个单层(例如研究目的层)顶、底反射子波间的变化来得到该层的某种地质信息，加之由单层引起的变化量十分微弱且影响因素十分复杂。因此，通过地震子波获得地层地质信息的想法是不现实的。

反射系数系列取决于各个地层界面上、下介质的波阻抗差，即由地层的波阻抗(波速与密度的乘积)所决定。单个界面的反射系数(R)为：

 

式中r表示密度；口表示波的传播速度；下标1表示界面上层，下标2表示界面下层。

地层的波阻抗由地层的岩性(包括岩石成分和结构等)和物性(包括岩石弹性、孔隙性和所含流体性质等)所决定。因此，可以说地震子波是地震记录中的物理元素，反射系数系列则是地震记录中的地质元素。通过从地震记录中提取反射系数系列是获取地层地质信息的必然和可能的途径，而利用得到的反射系数系列把地震波形剖面转换为波阻抗剖面将全面精确地展示出地震记录巾所蕴含的绝大部分地质信息。

2　关于前人“伪速度测井"方法的讨论

从式(1)可知，在地震子波、反射系数系列及地震记录道这3个元素中，原则上只要已知其中的两个，就可以求出第三者。

西方地球物理学家首先想到的是利用钻井得到的声速和密度测井曲线求取波阻抗曲线(常常假定密度为常数，而用声速测井代替波阻抗)，从而推算出反射系数系列，通过后者和井旁地震记录道进行运算，求得井旁地震子波的反褶积算子，也称为反子波。运算的原则是要求所求解的未知量(反褶积算子)与地震记录道褶积的期望输出为已知的反射系数系列。最后利用这一算子与剖面上各地震道记录褶积求得各接收道上的反射系数系列并最终推出波阻抗(或声速)曲线。经过如此运算得到的井旁道声速曲线，可以与井下实测曲线十分逼近。

该方法被命名为“伪速度测井”，或者“合成速度测井”“地震约束测井”等名称。

但上述方法的思路显然存在以下问题：

1)首先假定在简单的地震波动曲线中，蕴含着如声速测井这般细密丰富的地层信息完全足一个主观的认定。

2)井下声速测井与地震记录道的精度(或分辨率)是大相径庭的。前者是在深度域上记录，后者则是在时间域上进行的。如果假定的前提存在，那么要求相匹配的时深转换精度是目前无法想象的。

3)当人们从地震记录道上选取一个时窗去对应声速测井的某一段曲线时，时窗内的地震记录中存在截断效应，即时窗的前部混有时窗开始以前到达的某些波动的后续振荡，而时窗尾部的某些波动的后续能量却又被截除，这必然给反褶积算子的计算带来畸变。同时，时窗起点、尾点与其所对应的声速测井曲线相关联的首、尾是否称得上对应，也值得怀疑。

4)地震道上记录的振幅变化是一个非线性变化的量(所谓的真振幅恢复不可能改变这一特征)，而测井曲线上的读数则是一个线性变化的量，将两个不同性质变化的量通过数学运算进行约束或者逼近，仍然是一个未经论证的问题。

5)如上所述，地震子波是时变的，甚至可能是随时空变化的。伪速度测井将一口井上得到的反子波按时空不变地推广到一个区域，同样值得商榷。

某些研究者的实践证明了这些疑点，他们得到的结论是：“方法推不出去”，指的是从一口井上得到的反褶积算子运用到邻近的另一口井上，伪速度测井曲线与声速测井曲线就不再逼近了。而另一些研究者则错误地把无关的某个地震道记录与测井的声速曲线进行“约束”或“合成”，却仍然可以得到不错的逼近效果。

3　关于地震子波的讨论

综上所述，合理而且客观的方法应该是从地震记录中提取子波，再利用提取到的子波与地震记录道作运算，进而求出反射系数系列。较之反射系数系列，地震子波毕竟是一个相对简单、有章可循、容易认识的未知量。

为了从地震记录道中提取子波，笔者在本文参考文献[1]中推出了一个子波的数学表达式，即：

 

式中A为子波的瞬时(t)振幅，K为单元子波的振幅倍数，对于反射子波，如果入射子波为一个单元子波，那么K就是反射系数，这时K是一个分数，T是振荡的主周期，(e^-at2-e^-bt2)为子波包络，其中t是从包络起点t0开始计算的时间，t0可以是从激发开始计算的时间，也就是记录上的波至时间，子波包络还可以写成：

a=e^-at2-e^-bt2=(1-e^-gt2)e^-at2=a1·a2            (4)

显然a1=1-e^-gt2

a2=e^-at2

如图l所示，子波包络描述了振荡质点总的能量变化过程。包络前段上升的半枝反映了质点在开始受压力后能量不断增长积蓄的过程，主要受口值控制；后段下降的半枝，显示着质点能量消耗衰减的过程，主要由a值所决定；包络内的质点振荡曲线是大地介质在激发能量启动后固有的弹性谐振反应，描述着质点能量形式的转换。当质点位于平衡点(即振幅零线)上时，动能达到最大，而势能为零。当质点达到包络线位置时，势能达到最大而动能为零，而处在曲线的任意点位置时，该点到振幅零线的距离表示势能部分，到包络线的距离则表示动能部分，所以质点的波动记录描述的是势能随时间的变化情况。

 

于是在地震记录道上确定(或者提取)一个子波，就是确定(或者提取)a、b、T、K和t0这5个子波参数，其中前3个参数决定子波的波形，K决定子波的强度，t0为波至时间。研究子波随时间的变化，就成为研究子波的前4个参数随第5个参数t0的变化。

子波的a和b值愈大，受其控制的子波包络下降或上升半枝的时间愈短，当a和b值的匹配使得上升和下降的两个半枝相对称，且工值恰好相当子波包络内容纳一个半振荡时，式(3)就表示一个雷克子波。因此，式(3)是包含了雷克子波的更为广泛的子波表达式。雷克子波是在理想的无限均匀介质中，一个球形瞬时扩张和收缩激发下形成的理想子波。在地震测井的初至波记录中，我们只能在最浅层记录道的初至波(即距离激发点不远处)中观察到接近雷克子波的波形，随着波动在非均质的层状介质中传播，子波波形很快就偏离雷克子波。在我们感兴趣的时间段内，雷克子波是不存在的。

4　地震波形剖面到波阻抗剖面的转换

把地震波形剖面转换为波阻抗剖面的运算可以分为3个部分，也就是需要编制3个计算程序：

1)从地震记录道中提取时变的子波，即确定a、b、T这3个参数随时间的变化关系。

2)运用提取的时变子波从地震记录道中求取反射波的振幅系列，即确定K值及其相应的t0(时间)系列。

3)在把反射波振幅系列标定为反射系数系列的同时，推算出各反射层的波阻抗，得到波阻抗的相对变化剖面。

5　子波的提取

如前所述，a、b、T这3个子波参数随时间(或深度)是有变化的，而且其变化量不容忽视。其中了、随时间近于直线型地增大，a和b随时间则是近似双曲线型(以相互垂直的两根坐标线为渐近线)地减小，从一定的深度以下(如l500m左右)，也可以看作直线型的下降关系(本文参考文献[1]中的图8～10)。因此从地震记录中提取子波就被归结为确定子波参数值随时间的变化曲线。这一计算过程可以运用各种方式的优选法，以下提出一个可能的方案。

5．1　T的计算

确定地震曲线在振幅零线(时间轴)上的交点，可以得到一系列的视半周期数据，用最小二乘法按误差最小的原则可以求得半周期随时间变化的趋势线，然后除去与趋势线上、下偏离最大的两个数据再求趋势线，如此反复进行直到前后得到的两条线的差别小于某个认可的误差值时为止。从结果趋势线上可以求得任意时间上的T值。

5．2　a的计算

从图l上可以看出，子波包络在时间tm上有极大值`am(顶杆表示绝对值)，将包络曲线分成前部上升和后部下降的两个半枝。对后半枝而言，当t稍大于tm时，包含b值的e<sup>-gt2</sup>一项很快接近于零[式(4)]，这时`A≈Ke^-at2。

在地震记录道上取下所有峰值(包括波峰及波谷)的绝对值及其对应的时间，组成`Am-t数组(下标m表示峰值)。只要存在3个连续减小的峰值(表示振幅包络衰减的一段)，就可以求得一对a-t0的值。如果有n个连续减小的峰值，就可以求得n-2对。

如果连续减小的3个峰值为`Am1、`Am2、`Am3，对应的峰至时间为t1、t2、t3，那么就有t1=t1-t0、t2=t2-t0和t3=t3-tm。用Dtl2表示t2-tl和Dtl3表示t3-tl，于是：

1n(`Am1／`Am2)=a(t2²-t1²)用LAl2表示

ln(`Am1／`Am3)=a(t3²-t1²)用LA13表示

取LAl2与LA13的比值，经过整理化简，可得

 

于是

 

沿地震记录道求得一系列的a-t值(即a-t0)就可以用上述类似的方法优选出a与t的关系曲线。

5．3　b的计算

对地震记录道求取振幅包络，可以得到一系列包络的极大值`Am及其相应的时间tm。在每个峰值的后续包络线(下降的半枝)上适当间隔地取3个包络值A1、A2、A3，按式(5)可求出t0，于是tm-t0=tm。根据包络公式(4)对t取导数，有

 

对于值`Am来说，有¶`A／¶t0＝0。于是可得

tm＝[(1na-lnb)／(a-b)]^1/2

或者1nb-bt²m=1na-at²m          (7)

式(7)右端为一个已知确定的数，左端只有b为未知数，从而可求得b值。

同样通过求得一系列b-t(即b-t0)值，就不难得到b与t的关系曲线。

至此，我们就可以从一道地震记录上，建立起一个时变的子波模型，包括a(t)、b(t)和t(f)3条关系曲线，把这一模型子波视为实际子波的拟合子波，在工区而积内若干控制点卜求得子波模型，可以研究子波参数在面积内的变化，以确定子波是否应被视为时空变化的。该方法具有建立时空变化子波模型的可能性。

6　振幅系列K(t)的建立

在求得拟合子波的基础上，按褶积运算的逆过程去完成反褶积运算，称这种运算为消去反褶积^[2]。

可以采取的运算步骤很多，举例说明，将单元拟合子波振幅包络上升的半枝沿地震记录道从前往后逐(时间)点与记录道的振幅包络作比较，或者说用拟合子波的前端部分沿地震记录道进行搜索。寻找地震记录道振幅包络上与拟合子波的振幅包络形态最为近似的t0时刻。从这一时刻开始，在拟合子波包络的上升半枝区间，通过拟合子波与记录道各对应时间点(t对于拟合子波及t0+t对于记录道)振幅的比值(或反比值)，取算术平均值作为记录道上子波的K值(K可以是负数)。所谓形态最为近似也就是各点的比值最为接近。然后从地震道的“时刻开始减去拟合子波Ka，同时在开始是零的反射系数系列中建立一个反射系数K(t0)。如此逐个消减，或者说是剥离下去，各波至时间上的反射系数逐个形成。当被连续剥离后的地震记录道(可称为道剩余)上剩余的波动能量成为原始记录道总能量的一个很小的比例为止，得到最终的K(t)系列。

7　振幅系列向波阻抗曲线的转换

对一个足够浅的层位和一个足够深的层位给出各自合理的波阻抗值(人们已经知道各种岩类的速度与密度值的可能区间，也知道它们随深度变化的大致关系)。如果在某个位置上有相关的资料，如地震、声速、密度测井或岩样的实验室测定数据，参考这些资料可以使这种设定更有依据，但这并不是绝时必需的。

先暂定一个振幅K与反射系数的转换系数，以便从初始的浅层按式(2)沿振幅系列逐点(层)递推出波阻抗值，用迭代逼近的方法调整振幅与反射系数的变换系数，直到递推占的深层波阻抗值达到前面所给定的数据为止。

反射系数中的正数总和应大于负数总和(指绝对值)，因为深层的波阻抗总是大于浅层的。波阻抗曲线总的趋势是增大的。如果情况相反，则可以认定地震记录道的极性标定是被反转的。

这种固定地震道深浅两头波阻抗的做法，无形中消除了各地震道总能量因地表激发、接收条件不同，包括不同批次资料因接收系统灵敏度的不同而带来的差异，实际上也就是固定了剖面上各道所接收的总反射能量，但这并不妨碍剖面上实际存在的单层或局部波阻抗变化得到正确的保留和显示。

8　波阻抗成果剖面的相对性

这样，我们可以把一个地震波形剖面转换成一个显示波阻抗沿纵向和横向相对变化的剖面。在地震波的采集系统中，存在两个振幅增益控制：①是接近线型的半自动增益控制，用来补偿地震波能量因球面扩张而形成的振幅衰减，以便在地震道上使深浅层反射都得以均衡地表现；②完全非直线型的自动增益控制，它具有压制强反射和放大弱反射的功能，使得强弱反射能同时记录和显示在地震剖面每一道不宽的纵向条带内。这是两个有益和合理的对真振幅的改变。

不同厂家，甚至相同厂家不同批次生产的检波器的机电(机械振动到电压变化)转换系数可以是不同的，而且常常为非线型的，最后，由于我们并不知道到达每一反射界面的入射波振幅，因而反射波振幅在转换成反射系数时也只能是相对的。

所谓的真振幅恢复是不可能的，也是不足取的。因此我们强调所得到的只是一个相对变化的波阻抗剖面。

关于波阻抗剖面的显示方式，既要展示各波阻抗层的构造格局，又要展示研究目的层段的波阻抗变化，这就需要妥善综合地应用曲线和色标两种显示方式，以得到有关信息的最佳表现。

9　地震记录终极地质解释的目标和含义

地震记录终极地质解释的目标应该包含两重含义：①这种解释能够展示地震记录中客观存在的全部地质信息，反过来说，这种解释未能揭示的，就是记录中客观上根本不包含的；②这种解释所提供的表达各种信息在数量上的精度，应该是客观上能够达到的最高水平。

把地震剖面转换为相对波阻抗剖面达到这两点在理论上具有必然性和可能性，但不等于用本文所述的转换方案，就可以做到上述两点。这并非一蹴而就的事情。在整个转换运算中，存在着许多可以优化的环节。

从已经作过的初步尝试中^[1-4]，已经证实或者可以预见到地震波阻抗剖面较之地震波形剖面，在解释和应用上能够现实地达到以下实质性的进展。

1)从地震勘测所提供的基本信息元素而论，在上述的转换运算中，反射时间是通过分离的单个子波的半枝包络从记录上求到的波至初始时间，而不再是波至以后的某个强相位的波峰时间。求到的振幅值是单个子波振幅包络的极大振幅值，而不再是某个强相位的峰值。这两个信息元素都清除了前后波干涉所带来的误差，使得信息数值的定义更为明确，精度更为提高。振幅信息在转换为波阻抗的运算中，自然地减除了因地面激发接收条件不均一和不同批次资料因接收系统灵敏度的差异所带来的影响。这是放心地使用振幅信息所必需的。同时在消去反褶积的过程中所剥离的每个子波，或者说得到的每一个振幅值K，有明确的正负极性(确定反射的硬、软性质)。这是一个十分重要、而在波形剖面中所没能提供的新的信息元素。这些基本信息元素更确切的定义、更准确的量值和更多的内容是转换后的波阻抗剖面提升和扩展其实用价值的根本基础。

2)地震波阻抗曲线较之相对单调的地震波形曲线，形态上更为丰富，更具特征，特别是直接显示了地层的硬(高波阻抗)软(低波阻抗)变化。在岩性剖面已知地点(例如钻井井位旁)，可以不用速度曲线进行正确的地质层位标定。在剖面上对反射层作追踪对比时，特别是在地质构造复杂地区可以大大增加对比的确定性，减少误判，更准确地认识复杂地层构造结构上的性质。

3)在把地震波形记录分解为振幅系列的过程中，实际上极大程度地提高了反射波形的分辨率。分离后的前后两个反射振幅(代表两个反射子波)的间隔可以小于lms。这种分辨率已经不是从波形剖面上读数或识别的分辨率所可比拟的。因此波阻抗剖面只具有展示总体形态和地层变化趋势的功能，而在数据应用上，应该以电子存储器中的数据文件为主，通过各种软件进行管控、显示，在人机联作的环境中提取和应用。本文参考文献[4]中，展示了曾经在四川盆地川东忠县向斜中作过的一次尝试。该向斜中位于上、下页岩中的石炭系白云岩产层白西向东由60m左右的厚度被剥蚀减薄直至尖灭，在东北方向存在一个侵蚀窗。在侵蚀窗内，先后钻的2口井都未钻遇该产层而失败。从常规地震波形剖面上，看不到仟何反射层减薄尖灭的形迹。通过消去反褶积分解出对应石炭系顶底的前后一正一负2个反射系数清楚地由西向东靠拢，直至合并消失。最后编制了该地区的石炭系残厚图(等残厚线距10m)，侵蚀窗的边界被清晰确定。

4)精确的波至时间结合钻井剖面的分层数据，无需地震测井或VSP，可以计算出比较准确的平均速度和层速度，得到更丰富的速度信息，提高对工区速度变化规律的认知程度，给精确研究低幅度平缓构造带来了可能性，有利于波阻抗曲线自身准确的时深转换，更多地掌握速度与其他地质凶素的关联性，为波阻抗曲线的解释和应用建立更明确的依据。

5)岩层的波阻抗由多种因素所决定，包括地层的岩性(岩石成分和结构)和物性(岩石弹性、孔隙性和所含流体性质)等。在具体地区，研究目的层的多个因素往往是稳定和微变的，而明显变化的常常只是其中的某一个因素，于是波阻抗的变化就可以反映发生变化的那一个因素。在常规的波形剖面上，我们可以识别和描绘出一个特殊地下地质体的外形，而难以了解其内部结构。而运用波阻抗剖面，我们则有可能突破这一局限。例如对一个三角洲冲积扇体，从扇体内众多蚯蚓状的反射同相轴图像上，我们可以定性地推测砂泥岩交互沉积的存在，但仅此而已；而利用波阻抗剖面，我们就有可能解释出具体的砂泥岩结构。

6)巧妙和综合地应用波阻抗剖面，可以使许多地震资料的多解性变得单一·和明确。例如地震剖面上相邻两枝同相轴的分离和靠拢，可以是层厚也可以是速度的变化，参考波阻抗的相应变化，就可能作出单一的判断。在早期利用地震剖面振幅信息的亮点技术中，一段局部的强相位同相轴常被认定为下伏产气层的变化，从而因事实上的上覆泥岩盖层的变化而失误，这在波阻抗剖面的解释卜就不难避免。解释人员曾经把泥岩中的煤层反射亮点误认为气藏。如果波阻抗的平面异常与构造形态(常常是低幅度构造)相匹配，那么是气藏的可能性大，而煤层的分布则与构造形态多半是不配套的。

7)我们可以对一个注水采气的气田作出一种大胆的设想，如果通过理论模型计算(已知气层的岩性、厚度和孔隙度)，气层水侵部分与含气部分在波阻抗值和穿过产层的地震波旅行时间两项指标上的差异可以被我们所得到波阻抗剖面所察觉(即这种差异明显大于可能的误差)，那么就可以通过地面地震勘测来确定水线的推进位置，而更加科学地调整注水采气方案。对于被注水所切割的某些尚存的含气区块的测定，有利于更加合理的部署井位而提高采收率。这种设想可能性的客观依据，可以从亮点技术的成功例子中得到推断。诸如此类的针对具体地区的具体设想如果获得成功，将大大地拓展地震勘测的实际应用领域。

8)最后，整个转换处理是地震资料的自我改变，不受其他资料条件的限制，可进入常规处理范畴。

10　展望

可以预见，一旦一个考虑周密而又运用方便的波形到波阻抗的处理软件问世，在各种地质条件地区的应用实践中，地震波阻抗剖面的解释和应用，应将被演绎成为一门内容丰富而实用价值很高的学问，并最终服务于提升布井的针对性和钻井成功率。

 

参考文献

[1]徐中英．地震子波剖析[J]．石油地球物理勘探，l982，17(3)：1-15．

XU Zhongying．The analysis of seismic wavelet[J]．Oil Geophysical Prospecting，1982，17(3)：1-15．

[2]徐中英．地震波列分解[J]．石油地球物理勘探，l983，l8(2)：99-114．

XU Zhongying．The decomposition of seismic wave train[J]．Oil Geophysical Prospecting，1983，18(2)：99-114．

[3]徐中英，韦进平．换子波剖面[J]．石油地球物理勘探，l993，28(4)：439-446．

XU Zhongying，WEI Jinping．Converted wavelet section[J]．Oil Geophysical Prospecting，l993，28(4)：439-446．

[4]唐泽尧，徐中英．天然气开采工程丛书(第一册：气田开发地质)[M]．北京：石油工业出版社，l997：194-l97．

TANG Zeyao，XU Zhongying．Series of natural gas production engineering：Part I Development geology of gas fields[M]．Beijing：Petroleum Industry Press，l997：l94-197．

 

本文作者：徐中英

作者单位：原四川石油管理局地质勘探开发研究院