CN110208859A - 基于超声成像的油基泥浆井裂缝定量参数智能计算方法 - Google Patents

基于超声成像的油基泥浆井裂缝定量参数智能计算方法 Download PDF

Info

Publication number
CN110208859A
CN110208859A CN201910375453.5A CN201910375453A CN110208859A CN 110208859 A CN110208859 A CN 110208859A CN 201910375453 A CN201910375453 A CN 201910375453A CN 110208859 A CN110208859 A CN 110208859A
Authority
CN
China
Prior art keywords
crack
fracture
value
point
imaging
Prior art date
Legal status (The legal status is an assumption and is not a legal conclusion. Google has not performed a legal analysis and makes no representation as to the accuracy of the status listed.)
Granted
Application number
CN201910375453.5A
Other languages
English (en)
Other versions
CN110208859B (zh
Inventor
蔡明�
章成广
信毅
韩闯
陈义群
周箩鱼
姜龙
Current Assignee (The listed assignees may be inaccurate. Google has not performed a legal analysis and makes no representation or warranty as to the accuracy of the list.)
Yangtze University
Original Assignee
Yangtze University
Priority date (The priority date is an assumption and is not a legal conclusion. Google has not performed a legal analysis and makes no representation as to the accuracy of the date listed.)
Filing date
Publication date
Application filed by Yangtze University filed Critical Yangtze University
Priority to CN201910375453.5A priority Critical patent/CN110208859B/zh
Publication of CN110208859A publication Critical patent/CN110208859A/zh
Application granted granted Critical
Publication of CN110208859B publication Critical patent/CN110208859B/zh
Active legal-status Critical Current
Anticipated expiration legal-status Critical

Links

Classifications

    • G—PHYSICS
    • G01—MEASURING; TESTING
    • G01V—GEOPHYSICS; GRAVITATIONAL MEASUREMENTS; DETECTING MASSES OR OBJECTS; TAGS
    • G01V1/00—Seismology; Seismic or acoustic prospecting or detecting
    • G01V1/40—Seismology; Seismic or acoustic prospecting or detecting specially adapted for well-logging
    • G01V1/44—Seismology; Seismic or acoustic prospecting or detecting specially adapted for well-logging using generators and receivers in the same well
    • G01V1/48—Processing data

Landscapes

  • Physics & Mathematics (AREA)
  • Life Sciences & Earth Sciences (AREA)
  • Engineering & Computer Science (AREA)
  • Acoustics & Sound (AREA)
  • Environmental & Geological Engineering (AREA)
  • Geology (AREA)
  • Remote Sensing (AREA)
  • General Life Sciences & Earth Sciences (AREA)
  • General Physics & Mathematics (AREA)
  • Geophysics (AREA)
  • Investigating Or Analyzing Materials By The Use Of Ultrasonic Waves (AREA)

Abstract

本发明公开一种基于超声成像的油基泥浆井裂缝定量参数智能计算方法，其在分析超声成像裂缝响应特征的基础上，通过滤波处理突显裂缝响应，提高信噪比，然后通过图像智能识别算法自动拾取裂缝轨迹并计算裂缝长度、密度和倾角，接着沿裂缝轨迹法向方向自动探索并统计计算裂缝宽度，再根据裂缝长度、宽度、密度参数计算裂缝面孔率；同时，利用人工或计算机自动处理岩心滚扫照片获取岩心裂缝参数，并利用岩心裂缝参数对超声成像计算的裂缝参数进行刻度校正，以提高超声成像裂缝参数计算的精度。本发明用于水基和油基泥浆井的裂缝参数定量评价并指导储层评价，提高复杂非常规储层测井评价的准确性和可靠性，进而为制定合理高效的开发方案提供依据。

Description

基于超声成像的油基泥浆井裂缝定量参数智能计算方法
技术领域
本发明涉及应用地球物理学和地球物理测井领域，具体涉及一种基于超声成像的油基泥浆井裂缝定量参数智能计算方法。
背景技术
随着石油天然气工业深入发展，大型整装高孔高渗砂岩油气藏越来越难以找到，许多大公司不得不把注意力转移到裂缝－复杂岩性油气藏的勘探和开发上来。致密裂缝型油气藏是21世纪石油增储上产的重要领域之一，在我国，裂缝型低渗透储层和数量比例更为突出，裂缝型低渗储层油气产量占整个石油天然气产量的一半以上，占未来准备投产的石油天然气储量的三分之二以上。
对于低渗致密裂缝性储层，由于其基质的低孔低渗特性，裂缝则作为其主要渗流途径，在孔道间起到连接作用，提高了储层的渗透率，为储层产能的提高提供了基础。故天然裂缝的识别和特征参数的评价是这类储层评价的一个十分重要的方面。另外，对于页岩气、致密气等非常规储层，往往需要压裂改造，而压裂施工后人造裂缝的发育情况评价同样十分重要。
利用测井资料识别并精细评价裂缝是储层裂缝评价最主要的手段，裂缝测井评价方法主要包括常规测井评价法、成像测井评价法、阵列声波测井评价法和反射声波成像测井评价法。常规测井评价法侧重于裂缝的定性评价。成像测井评价法主要是利用微电阻率成像和超声成像测井资料对穿过井轴的裂缝进行评价，可以获得裂缝密度、裂缝宽度、裂缝倾角、裂缝长度、裂缝面孔率等裂缝定量参数，被认为是目前可信度最高的一种裂缝测井评价技术。在水基泥浆钻井液条件下，微电阻率成像测井已广泛应用于裂缝的识别和裂缝定量参数的计算，并取得了较好的应用效果。近年来，随着超深致密油气藏、页岩油气藏等非常规油气藏的广泛勘探开发，大量的井采用了油基钻井液；而油基泥浆的导电性差，导致基于电学性质的裂缝识别与评价方法应用效果明显变差。超声成像测井不受泥浆电阻率的影响，可以由扫描测量记录到的回波波形处理得到的声阻抗和回波时间参数提供井壁360°方位的成像图，已广泛应用于油基泥浆井裂缝的识别和评价，且一般情况下效果优于微电阻率成像(如图1)。
利用超声成像测井资料识别裂缝并定量计算裂缝参数是基于超声成像的裂缝评价法的关键步骤，直接影响着该方法的应用效果。
目前，国内外已有商业化应用的超声成像测井仪(如CBIL和UXPL等)，但鲜见基于超声成像的裂缝参数定量评价方法相关的文献报道，针对裂缝评价的超声成像资料处理主要是由人工观察拾取超声成像图上可能的裂缝轨迹，并统计裂缝的数目和长度，进而计算裂缝的密度和倾角参数。这种处理方法效率较低，人为因素对处理结果影响较大，且不能实现裂缝宽度、裂缝面孔率等参数的精细计算。利用现代先进信号处理方法对超声成像测井资料进行自动处理，智能识别裂缝，定量计算裂缝参数并利用岩心裂缝参数对超声成像裂缝参数进行刻度以提高裂缝参数计算的精度是有前景的超声成像裂缝评价方法发展方向。
发明内容
本发明的目的在于克服现有技术的不足，提供了一种基于超声成像的油基泥浆井裂缝定量参数智能计算方法，该方法是在分析超声成像裂缝响应特征的基础上，通过滤波处理突显裂缝响应，提高信噪比，然后通过图像智能识别算法自动拾取裂缝轨迹并计算裂缝长度、密度和倾角，接着沿裂缝轨迹法向方向自动探索并统计计算裂缝宽度，再根据裂缝长度、宽度、密度参数计算裂缝面孔率；同时，利用人工或计算机自动处理岩心滚扫照片获取岩心裂缝参数，并利用岩心裂缝参数对超声成像计算的裂缝参数进行刻度校正，以提高超声成像裂缝参数计算的精度。
为实现上述目的，本发明所设计一种基于超声成像的油基泥浆井裂缝定量参数智能计算方法，包括以下步骤：
1)数据解编、加载与图像着色
利用超声波探头对井壁发射超声波脉冲并接收回波信号，得到由不同深度不同方位对应的回波幅度(或回波时间)值组成的二维数据体，即为声成像测井数据(为了便于后续裂缝评价处理，需首先解编加载原始测井数据到处理软件中，并以图像形式显示声成像测井二维数据体)；
2)预处理(由于测量过程中仪器偏心、仪器非匀速运动、仪器旋转等问题会引起成像图上裂缝响应形态和位置发生畸变，在进行图像处理之前，必须要对声成像测井数据进行预处理)
在声成像图图像处理之前，对声成像测井数据进行加速度校正(加速校正是利用测井记录的加速度曲线等信息对测井数据的深度进行修正的处理，是比较成熟的测井资料预处理方法)、方位校正(方位校正是为了校正仪器在测井过程中旋转引起的成像测量点方位偏差)和回波时间偏心校正(超声成像测井通常记录回波幅度(AMP) 和回波时间(TT)数据，回波幅度与井壁的声阻抗和路径有关，而回波时间只与路径有关。当测井过程中仪器不居中时，测量的对称点的回波时间就不相等，利用时间计算的井眼半径就不合理，回波幅度也不能真正反映井壁声阻抗的变化，因此要进行偏心校正)的预处理：
a.方位校正是根据仪器方位曲线数据将各深度点的成像样点方位调整到统一的大地方位坐标系中，具体调整公式为：
其中，N为一周的扫描样点总数，s0[i]为一周扫描的第i个成像样点值，AZ为仪器方位角，％为求余运算，为向下取整运算，s1[i]为方位校正后的第i个成像样点值；
b.偏心校正处理步骤如下：
①对每一深度处回波时间数据进行统计，计算最大值MAX和最小值MIN；
②计算直方图：将统计的最小值到最大值区间段[MIN,MAX]细分为N等分，每一区间范围为：
[MIN+(i-1)*(MAX-MIN)/N,MIN+i*(MAX-MIN)/N](细分后的区间段按从小到大的顺序排列)；
其中，N为时间成像图一周内的样点数，i为区间编号；
然后再对全井段所有数据进行分区统计，得到每一区间的频数；
③修正直方图极值
将直方图的两端各去掉2.5％的点后，可得到直方图新的两端极值，最小值为MIN1，最大值MAX1(MIN1是两端各去掉2.5％样点周的最小值，MAX1是两端各去掉2.5％样点周的最大值)；
④以仪器的1号极板方位作为Y轴(如所示)将位于[MIN1,MAX1]范围内的所有回波时间样点值TT投影到X和Y轴上，则有：
TTX[i]＝TT[i]*sin(α[i])
TTY[i]＝TT[i]*cos(α[i])
α[i]＝i/m*360
其中，i为处理深度点上位于[MIN1,MAX1]范围内的回波时间样点序号，m为处理深度点上位于[MIN1,MAX1]范围内的回波时间样点总数；则X和Y轴上的时间投影之和为：
偏心距：偏心方位：α＝arctan(XX/YY)；
⑤偏心距偏心方位滤波
a)在由各深度偏心距和偏心方位构成的数据序列中选取N个点作为窗长对目标井段所有深度点的偏心距和偏心方位进行滤波处理，其中，N为奇数，选3、5、7、 9……；处理方法为：
计算出窗长内偏心距、偏心方位平均值和方差，则平均值为本点(窗长内的中心点)的偏心距和偏心方位；当偏心距的方差大于2.5时，则本点的偏心距为空值；当偏心方位的方差大于4.5时，则本点的偏心方位为空值；
b)对于偏心距、偏心方位为空值的点，以本点为中心上下各取4个有效点(非空值点)计算其平均值，用计算的平均值取代本空值点；
⑥偏心幅度校正处理(根据回波时间偏心校正结果结合超声波在泥浆中传播的衰减规律即可进行偏心幅度校正；偏心校正是对一周所有成像点进行处理的；这样就消除了仪器偏心对测量结果的影响，得到仪器居中时的正确测量结果)
根据仪器偏心时仪器轴中心点与井筒截面的相对位置几何关系，由仪器偏心时记录的回波时间计算出仪器居中时的理论回波时间，实现仪器偏心引起的测量误差校正；具体校正公式为：
其中，N为测量一周的扫描点数，i为参考方位(即偏心方位)回波时间样值在本深度点所有方位回波时间样值序列中的序号，
j为实际测量方位回波时间样值在本深度点所有方位回波时间样值序列中的序号，
TT[j]为实际测量方位回波时间样值，
α为实际测量方位与参考方位的夹角，
TT1[j]为仪器居中时实际测量方位的理论回波时间样值，
h为TT1[j]在参考方位的投影值，
d为TT[j]在参考方位的投影值与h之差；
3)滤波及裂缝自动拾取
a.提高声成像图的对比度：采用直方图均衡化法提高成像图的对比度。直方图均衡化是通过对原图像进行灰度变换，使变换后的图像直方图能均匀分布，这样就能使原图像中具有相近灰度并占有大量像素点的区域的灰度范围得到展宽，使大区域中的微小灰度变化显现出来，从而使得目标更清晰。相对于一般直方图均衡化算法，限定对比度直方图均衡化法更有利于突出图像边缘，增强图像对比度以及抑制图像噪声)；
b.提取裂缝响应图像；
利用拉普拉斯算子提取裂缝的二值图像(拉普拉斯算子作为一个二阶导数算子，常常用于边缘增强及图像分割，它的优势在于：①它是基于零穿越检测算法，定位比较准确；②对细节有较强的响应，因此，对于对比度不强的裂缝也能提取出来)；
c.提高裂缝响应轨迹的连续性(提取的裂缝响应轨迹二值图像中有大量的孤立噪声点，这是由于拉普拉斯算子裂缝响应轨迹提取方法对孤立噪声点较为敏感引起的。为了进一步降低孤立噪声点的影响，并提高裂缝响应轨迹的连续性，先采用自适应滤波进行去噪处理，再利用形态学处理方法补充裂缝响应轨迹间断。需要说明的是，裂缝响应图像的提取和噪声的消除也可以采用二维小波变换或大数据人工智能等方法实现)
先采用自适应滤波进行去噪处理，再利用形态学处理方法补充裂缝响应轨迹间断；具体步骤如下：
(1)自适应滤波处理(自适应滤波处理的基本依据为小窗口滤波时，图像细节信息保持较好，但当噪声密度较大时，滤波性能较差；相反，较大窗口进行滤波时，去噪性能较好，但细节容易模糊，因此需要结合二者的优点)：
①将图像分成大小相等的子块；
②估计各子块的噪声方差大小；
③依据噪声方差确定滤波窗口，进行滤波
(2)裂缝响应轨迹间断弥补处理
利用下述公式对B作关于原点的映射，再将其连续的在A上移动来进行膨胀运算进而实现裂缝响应轨迹间断弥补处理；
膨胀的运算符为A用B来膨胀写作其定义式为
其中，A为待处理的二值图像，B为卷积模板(有圆形、方形、十字形模板)，x 代表像素点，为B的反射；
d.裂缝自动拾取
采用Hough变换自动拾取裂缝响应轨迹；具体实现步骤为：
①计算正弦曲线的频率：ω＝1/2πr，其中r为井眼半径；
②通过正弦曲线的对称性确定基线位置y0；
③采用二维Hough变换提取该位置所对应的正弦曲线的幅度A和初相位参数
④根据上述参数及正弦曲线函数即在成像图上标示出轨迹；其中，井壁成像的平面展开图上裂缝响应轨迹为一条正弦曲线，其方程表示为：
4)裂缝定量参数计算
拾取裂缝响应轨迹后，进一步定量计算裂缝的产状参数、裂缝密度、裂缝长度、裂缝宽度和裂缝面孔率，其中，
采用沿裂缝轨迹法线方向自动搜索的方式来统计计算裂缝的宽度(裂缝宽度是评价裂缝有效性的关键参数，也是后续计算裂缝面孔率的关键输入参数)，具体实现方式如下：
①在裂缝轨迹上等间距选择N1(N1>10，建议取30，N1越大裂缝宽度计算精度越高)个点；
②以选择的点为起点沿裂缝轨迹法向方向搜索局部像素极小点MIN和极大点MAX1和MAX2；
③计算裂缝边界像素值，计算公式为：
其中c1权重系数，取值范围为(0,1)，一般取0.4～0.8，由实验或经验确定；
④根据裂缝边界像素值EP1和EP2，沿裂缝轨迹法向方向上自动搜索裂缝边界像素点，并计算两个裂缝边界像素点之间的距离，即得裂缝轨迹上所选点处的裂缝宽度；
⑤将计算的裂缝轨迹上全部所选点处的裂缝宽度值从小到大排序，去掉排序后裂缝宽度值序列两端各10％的宽度值，计算剩余80％的裂缝宽度值的平均值，并作为该条裂缝的宽度；
⑥由上述计算的单条裂缝宽度结合裂缝平均宽度定义(单位井段中裂缝轨迹宽度的平均值)即计算出裂缝平均宽度值，计算公式为：
其中，FVAi为第i条裂缝的宽度；
5)岩心裂缝参数提取
采用上述声成像裂缝参数评价方法对岩心滚扫照片进行自动处理提取得到岩心裂缝参数；
或者，利用图像处理软件对岩心滚扫照片进行人工处理；其中，利用Image J图像处理软件提取得到岩心裂缝参数，步骤如下：
①将岩心滚扫彩色照片转换为灰度图像；
②测量全直径岩心柱直径进而可以计算出岩心柱的周长，故利用岩心柱的实际周长和其滚扫照片中的周长之间的关系来刻度滚扫照片的像素尺度；
③测量裂缝长度和裂缝宽度(裂缝宽度测量过程中，裂缝边界的确定是裂缝宽度准确测量的关键。可以根据裂缝处像素灰度值的变化规律由以上成像裂缝宽度计算方法部分所述的裂缝边界判定准则确定裂缝的边界。另外，为了提高裂缝宽度测量的准确性，应选择裂缝比较稳定规则处进行测量，且应多测量几处求取平均值作为此裂缝的宽度)；
④裂缝平均宽度和裂缝面孔率的计算；
统计一定窗长内岩心照片上所有裂缝的宽度值，并计算平均值即得裂缝平均宽度，具体计算公式为：
其中，Wi为第i条裂缝的宽度，Nw为统计窗长内裂缝的条数；
根据上述计算的岩心裂缝的长度和宽度结合裂缝面孔率的定义，计算出统计窗长内的裂缝面孔率，具体计算公式为：
其中，Li为第i条裂缝的长度，Lw为统计窗口长度，Rcore为岩心柱的半径；
6)超声成像裂缝参数刻度校正(岩心观察与描述能提供关于裂缝参数和发育情况、力学性质、充填特征和含油气性等的第一手资料，是最直接、最有效和最可靠的裂缝评价方式，故评价的裂缝参数精度较高。为了提高声成像裂缝参数评价的精度，有必要利用岩心裂缝参数对声成像裂缝参数进行刻度，成像裂缝参数的刻度实质上是通过对比分析成像资料处理得到的裂缝参数和岩心观察与描述得到的裂缝参数，建立二者之间的刻度系数，以用于校正成像资料处理得到的裂缝参数，最终达到精细评价储层裂缝发育情况的目的)
利用岩心裂缝参数对声成像裂缝参数进行刻度，步骤如下：
(1)通过上述岩心观察与描述方法获取岩心裂缝参数，并利用岩心地面自然伽马值把岩心深度归位到常规GR曲线统一的深度尺度上；
(2)利用上述声成像测井资料处理方法处理成像原始数据得到成像图中的裂缝面孔率等参数，之后将处理得到的裂缝面孔率等曲线样值也归位到常规GR曲线统一的深度尺度上；
(3)对照岩心照片裂缝参数和声成像处理得到的成像裂缝参数，分析两者存在的关系，得到两者之间的刻度系数，即得到岩心校正刻度系数(从而达到用岩心照片裂缝面孔率刻度成像面孔率的目的)；
(4)利用岩心校正刻度系数对超声成像资料处理得到的裂缝参数进行刻度校正即可得到更符合实际的裂缝定量参数；具体校正公式为：
其中，COEH为裂缝宽度校正系数，COEA为裂缝面孔率校正系数，FVAH为刻度校正后的裂缝宽度，FVAHUXPL为声成像处理得到的裂缝宽度；FVPA为刻度校正后的裂缝面孔率，FVPAUXPL为声成像处理得到的裂缝面孔率。
进一步地，所述步骤3)第a小步中，提高声成像图的对比度具体处理步骤如下：
①将预处理后的声成像图像分为大小相等的若干个子块；
②统计各个子块的局部直方图；
③计算子块内所有像素数平均分配到的各灰度级值A；
④设定一个剪切系数C，取值范围为[0,1]。则实际剪切值；
V＝A+C(Nx×Ny-A)
其中，Nx，Ny分别为子块水平方向和垂直方向的像素数目；
⑤利用剪切值，对局部直方图各个像素进行剪切，多余的像素重新分配到各直方图的各灰度级中，于是可以得到每个灰度级应该分到的像素数NACP为：NACP＝NCLIP/G
其中G表示为子块中的灰度级，NCLIP表示已剪切的像素总数；
假设原直方图为H(i)，由此得到重新分配后的直方图CH(i)：
⑥对每个子区域的直方图进行均衡化处理，处理后裂缝响应轨迹更清晰。
再进一步地，所述步骤3)第b小步中，
①构造拉普拉斯算子：
②通过拉普拉斯算子提取裂缝响应轨迹的二值图像。
再进一步地，所述步骤4)中，
a.裂缝产状参数计算
根据拾取的裂缝轨迹和方位信息很容易确定裂缝的倾向和走向，裂缝倾角可由如下公式计算得到：
其中，H为裂缝低点到高点的高度，R为井眼半径；
b.裂缝密度和长度计算
由裂缝响应轨迹拾取结果并结合裂缝密度定义(单位长度井壁上所见到的裂缝总条数)容易得到裂缝密度计算公式如下：
其中，Nw为统计窗长内裂缝条数，单位：条；Lw为窗口长度，单位：m；
由裂缝响应轨迹拾取结果利用线积分的方式容易计算出每条裂缝的长度，再根据裂缝长度的定义得到裂缝密度计算公式如下：
其中，Li为线积分方式计算的第i条裂缝的长度，单位：m；
c.裂缝面孔率计算
由上述计算的裂缝长度和宽度参数结合裂缝面孔率的定义，即计算出裂缝面孔率，计算公式为：
本发明的有益效果：
1)本发明可以实现裂缝的自动识别与拾取，处理效率高，人为干扰因素少；
2)本发明可以对现有的超声成像资料进行处理，并定量计算裂缝相关参数；
3)本发明可以利用岩心滚扫照片比较准确地提取岩心裂缝参数，并可进一步用于各种测井计算裂缝参数的刻度；
4)本发明可以实现超声成像裂缝参数的刻度校正，提高成像裂缝参数计算的精度；
5)本发明可以用于水基泥浆井和油基泥浆井的裂缝参数定量评价，并用于指导储层评价，提高复杂非常规储层测井评价的准确性和可靠性，进而为制定合理高效的开发方案提供依据。
附图说明
图1为X油基泥浆井超声成像和电阻率成像对比图；
图2为基于超声成像的油基泥浆井裂缝定量参数智能计算流程图；
图3为探头高速(左)和低速(右)对井壁超声扫描的螺旋型路径；
图4为方位校正前后超声成像对比图；
图5为回波时间样点投影示意图；
图6为仪器偏心示意图；
图7为偏心校正前后的幅度成像对比图；
图8为提高对比度前(左)后(右)的声成像对比图；
图9为提取的裂缝响应轨迹二值图像(右)与原始声成像(左)对比图；
图10为孤立噪声点的消除及裂缝间断弥补处理前(左)后(右)的对比图；
图11为裂缝响应示意图，其中，左图为裂缝穿过井孔示意图，右图为井壁成像平面展开图；
图12为裂缝响应轨迹拾取效果图，其中，第二道为预处理后的成像图，第三道为包含裂缝响应轨迹拾取结果的成像图；
图13为裂缝宽度计算方法示意图；
图14为Image J软件主菜单界面；
图15为像素尺度刻度操作界面；
图16为裂缝长度测量操作界面；
图17为裂缝长度记录结果；
图18为裂缝宽度测量操作界面；
图19为裂缝宽度记录结果；
图20为利用岩心裂缝参数刻度成像裂缝参数流程；
图21为X地区油基泥浆井岩心裂缝宽度与声成像裂缝宽度交会图；
图22为X地区油基泥浆井岩心裂缝面孔率与声成像面孔率交会图。具体实施方式
下面结合具体实施例对本发明作进一步的详细描述，以便本领域技术人员理解。
如图2所示的基于超声成像的油基泥浆井裂缝定量参数智能计算方法，包括以下步骤：
1.数据解编、加载与图像着色
利用超声波探头对井壁发射超声波脉冲并接收回波信号，得到由不同深度不同方位对应的回波幅度(或回波时间)值组成的二维数据体，即为声成像测井数据(为了便于后续裂缝评价处理，需首先解编加载原始测井数据到处理软件中，并以图像形式显示声成像测井二维数据体)；
2.预处理
超声成像测井是仪器在上提的同时超声波探头以一定的转速旋转并向井壁发射超声波脉冲和接收回波信号(超声波探头每旋转一周记录N个方位的井壁回波信号)，探头对井壁超声扫描的路径如图3所示；通过记录在不同深度不同方位处的井壁回声信号幅度和传播时间来反映井壁的声阻抗和井径变化，从而提供高分辨率环井壁 360°全方位的岩石物理二维图像信息，把地层岩性、裂缝、孔洞和层理等地层特征引起的声阻抗差异，转换成图像上不同色标显示，通过对声成像图进行处理即可获得裂缝发育信息。声成像的物理基础是不同介质的声阻抗和表面的粗糙程度不同，对声波的反射能力不同。
由于测量过程中仪器偏心、仪器非匀速运动、仪器旋转等问题会引起成像图上裂缝响应形态和位置发生畸变，在进行图像处理之前，必须要对声成像测井数据进行加速度校正、方位校正和偏心校正的预处理。
(1)加速校正
加速校正是利用测井记录的加速度曲线等信息对测井数据的深度进行修正的处理，是比较成熟的测井资料预处理方法。
(2)方位校正
方位校正是为了校正仪器在测井过程中旋转引起的成像测量点方位偏差，是根据仪器方位曲线数据将各深度点的成像样点方位调整到统一的大地方位坐标系中，具体调整公式为：
其中，N为一周的扫描样点总数，s0[i]为一周扫描的第i个成像样点值，AZ为仪器方位角，％为求余运算，为向下取整运算，s1[i]为方位校正后的第i个成像样点值。
方位校正处理前后超声成像对比图由所示，由图可以看出，方位校正处理后蓝色方框内的暗色裂缝响应轨迹为正弦曲线形态，与理论上的裂缝响应轨迹形态一致，说明本方位校正处理可以修正因仪器旋转引起的裂缝响应轨迹形态畸变。
(3)偏心校正
超声成像测井通常记录回波幅度(AMP)和回波时间(TT)数据，回波幅度与井壁的声阻抗和路径有关，而回波时间只与路径有关。当测井过程中仪器不居中时，测量的对称点的回波时间就不相等，利用时间计算的井眼半径就不合理，回波幅度也不能真正反映井壁声阻抗的变化，因此要进行偏心校正。回波时间偏心校正处理方法和步骤如下：
①对每一深度处回波时间数据进行统计，计算最大值MAX和最小值MIN；
②计算直方图：将统计的最小值到最大值区间段[MIN,MAX]细分为N等分，每一区间范围为
[MIN+(i-1)*(MAX-MIN)/N,MIN+i*(MAX-MIN)/N](细分后的区间段按从小到大的顺序排列，N为时间成像图一周内的样点数，i为区间编号)，然后再对全井段所有数据进行分区统计，就可得到每一区间的频数。
③修正直方图极值
将直方图的两端各去掉2.5％的点后，可得到直方图新的两端极值，最小值为MIN1，最大值MAX1。
④以仪器的1号极板方位作为Y轴(如所示)将位于[MIN1,MAX1]范围内的所有回波时间样点值TT投影到X和Y轴上，则有：
TTX[i]＝TT[i]*sin(α[i])
TTY[i]＝TT[i]*cos(α[i]) (2)
α[i]＝i/m*360
其中，i为处理深度点上位于[MIN1,MAX1]范围内的回波时间样点序号，m为处理深度点上位于[MIN1,MAX1]范围内的回波时间样点总数。则X和Y轴上的时间投影之和为：
则有偏心距：偏心方位：α＝arctan(XX/YY)。
⑤偏心距偏心方位滤波
滤波处理第一步：选9点作为窗长对目标井段所有深度点的偏心距和偏心方位进行滤波处理，处理方法为：计算出窗长内偏心距、偏心方位平均值和方差，则平均值为本点(窗长内的中心点)的偏心距和偏心方位。当偏心距的方差大于2.5时，则本点的偏心距为空值；当偏心方位的方差大于4.5时，则本点的偏心方位为空值。
滤波处理第二步：对于偏心距、偏心方位为空值的点，以本点为中心上下各取4 个有效点(非空值点)计算其平均值，用计算的平均值取代本空值点。
⑥偏心校正处理
图6为仪器偏心示意图，由图可以看出，当仪器不居中时，仪器中心点在B点，这时需要将B点记录的数据校正到井轴中心O点处。假设图中P点(本点在数据上为第i点)为长轴点(偏心方位)，测量点为A点(本点在数据上为第j点)，N为测量一周的扫描点数，测量回波时间为：TT[j](即AB)，偏心校正后的回波时间为： TT1[j](即AO)。则对于测量点A根据几何关系有：
根据式(4)结合原始回波时间TT、计算的偏心方位和偏心距即可求解出TT1[j]。再根据回波时间偏心校正结果结合超声波在泥浆中传播的衰减规律即可进行偏心幅度校正。偏心校正是对一周所有成像点进行处理的。这样就消除了仪器偏心对测量结果的影响，得到仪器居中时的正确测量结果。偏心校正前后的幅度成像对比图如所示，其中第2、3道为原始静态和动态幅度成像图，第5、6道为偏心校正后的静态和动态幅度成像图，由图中蓝色方框部分的对比可以看出，偏心校正处理后裂缝响应轨迹更清晰，连续性更好。
3.滤波及裂缝自动拾取
1)提高声成像图的对比度
采用直方图均衡化法提高成像图的对比度。直方图均衡化是通过对原图像进行灰度变换，使变换后的图像直方图能均匀分布，这样就能使原图像中具有相近灰度并占有大量像素点的区域的灰度范围得到展宽，使大区域中的微小灰度变化显现出来，从而使得目标更清晰。相对于一般直方图均衡化算法，限定对比度直方图均衡化法更有利于突出图像边缘，增强图像对比度以及抑制图像噪声。具体处理步骤如下：
①将预处理后的声成像图像分为大小相等的若干个子块；
②统计各个子块的局部直方图；
③计算子块内所有像素数平均分配到的各灰度级值A；
④设定一个剪切系数C，取值范围为[0,1]。则实际剪切值；
V＝A+C(Nx×Ny-A) (5)
其中，Nx，Ny分别为子块水平方向和垂直方向的像素数目。
⑤利用剪切值，对局部直方图各个像素进行剪切，多余的像素重新分配到各直方图的各灰度级中，于是可以得到每个灰度级应该分到的像素数NACP为：NACP＝NCLIP/G
其中G表示为子块中的灰度级，NCLIP表示已剪切的像素总数。
假设原直方图为H(i)，由此得到重新分配后的直方图CH(i)：
⑥对每个子区域的直方图进行均衡化处理。
提高对比度处理前后的声成像对比图如图8所示，由图可以看出，处理后裂缝响应轨迹更清晰。
2)提取裂缝响应图像
利用拉普拉斯算子提取裂缝的二值图像。拉普拉斯算子作为一个二阶导数算子，常常用于边缘增强及图像分割，它的优势在于：①它是基于零穿越检测算法，定位比较准确；②对细节有较强的响应，因此，对于对比度不强的裂缝也能提取出来。其主要步骤为：
a.构造拉普拉斯算子：
b.通过拉普拉斯算子提取裂缝响应轨迹的二值图像。
提取的裂缝响应轨迹二值图像与原始声成像对比图如图9所示。
3)提高裂缝响应轨迹的连续性
如上图9所示，提取的裂缝响应轨迹二值图像中有大量的孤立噪声点，这是由于拉普拉斯算子裂缝响应轨迹提取方法对孤立噪声点较为敏感引起的。为了进一步降低孤立噪声点的影响，并提高裂缝响应轨迹的连续性，先采用自适应滤波进行去噪处理，再利用形态学处理方法补充裂缝响应轨迹间断。需要说明的是，裂缝响应图像的提取和噪声的消除也可以采用二维小波变换或大数据人工智能等方法实现。
a.自适应滤波处理
自适应滤波处理的基本依据为小窗口滤波时，图像细节信息保持较好，但当噪声密度较大时，滤波性能较差；相反，较大窗口进行滤波时，去噪性能较好，但细节容易模糊，因此需要结合二者的优点。自适应滤波处理的具体流程为：
①将图像分成大小相等的子块；
②估计各子块的噪声方差大小；
③依据噪声方差确定滤波窗口，进行滤波；
b.裂缝响应轨迹间断弥补处理
滤波后去除了大部分的孤立噪声点。然而，由于原始成像中有些部分裂缝响应的像素和背景相差不大，这就使得提取裂缝响应轨迹的时候对应的部分就会存在缺失，形成间断。
形态学中的膨胀具有将裂缝间断桥接起来的功能，因此，利用形态学处理方法可以适当地弥补裂缝间断。其基本原理为：
膨胀的运算符为A用B来膨胀写作其定义式为
其中，A为待处理的二值图像，B为卷积模板(有圆形、方形、十字形模板)，x 代表像素点，为B的反射；
也就是说，如果将B看作一个卷积模板，膨胀运算就是先对B作关于原点的映射，再将其连续的在A上移动来实现的。因此，通过膨胀运算，能补充一些裂缝间断，从而有助于后续的裂缝识别。孤立噪声点的消除及裂缝间断弥补处理前后的对比图如图 10所示，由图可以看出，处理后孤立噪声点明显减少，裂缝响应轨迹连续性更好。
4)裂缝拾取
如图11所示，倾斜裂缝穿过井孔时，井壁成像的平面展开图上裂缝响应轨迹为一条正弦曲线，用方程表示为：
基于上述的裂缝响应轨迹特征，采用Hough变换拾取裂缝响应轨迹。Hough变换的基本原理为：图像空间的直线与描述它的极坐标表达式参数有唯一的对应关系，若以这些参数形成一个空间，则图像空间的一条直线必对应于参数空间阵列中的一个点。建立这种直线和点的关系的运算就称为Hough变换，Hough变换最初用于检测直线。
对于正弦曲线来说，式(8)中所描述的共有4个参数：振幅A、频率ω、初始相角基线位置(偏移量)y0。因此，对于任意可能的正弦曲线的Hough变换对应四维空间中的一个点。而要建立四个参数的Hough变换空间，将消耗十分巨大的资源。为了减少时间和存储空间的需求，提高检测速度，根据测井图像中正弦裂缝的特点，将四维Hough变换转化为二维Hough变换，其过程为：
①计算正弦曲线的频率：ω＝1/2πr；
②通过正弦曲线的对称性确定基线位置y0；
③采用二维Hough变换提取该位置所对应的正弦曲线的幅度和初相位参数。
通过Hough变换拾取的裂缝响应轨迹如图12中黑色虚线所示，由图可以看出，该方法可以较好地拾取裂缝响应轨迹。
4.裂缝定量参数计算
拾取裂缝响应轨迹后就可以进一步定量计算裂缝的产状参数、裂缝密度、裂缝长度、裂缝宽度和裂缝面孔率等。
1)裂缝产状参数计算
根据拾取的裂缝轨迹和方位信息很容易确定裂缝的倾向(如图11所示，裂缝轨迹高点方位到地点方位的方向即为裂缝倾向)和走向(水平面上，与裂缝倾向垂直的方向即为裂缝走向)，裂缝倾角可由如下公式计算得到：
其中H为裂缝低点到高点的高度，R为井眼半径。
2)裂缝密度和长度计算
由裂缝响应轨迹拾取结果并结合裂缝密度定义(单位长度井壁上所见到的裂缝总条数)容易得到裂缝密度计算公式如下：
其中，Nw为统计窗长内裂缝条数，单位：条；Lw为窗口长度，单位：m。
由裂缝响应轨迹拾取结果利用线积分的方式容易计算出每条裂缝的长度，再根据裂缝长度的定义(每平方米井壁上所有裂缝长度之和)容易得到裂缝密度计算公式如下：
其中，Li为线积分方式计算的第i条裂缝的长度，单位：m。
3)裂缝宽度计算
裂缝宽度是评价裂缝有效性的关键参数，也是后续计算裂缝面孔率的关键输入参数。本技术方案采用沿裂缝轨迹法线方向自动搜索的方式来统计计算裂缝的宽度，具体实现方式如下：
①如图13(a)所示，在裂缝轨迹上等间距选择N1(N1>10，建议取30，N1越大裂缝宽度计算精度越高)个点；
②以选择的点为起点沿裂缝轨迹法向方向搜索局部像素极小点MIN和极大点MAX1和MAX2(如图13(b))；
③计算裂缝边界像素值，计算公式为：
其中c1权重系数，取值范围为(0,1)，一般取0.4～0.8，可由实验或经验确定。
④根据裂缝边界像素值EP1和EP2，沿裂缝轨迹法向方向上自动搜索裂缝边界像素点，并计算两个裂缝边界像素点之间的距离，即得裂缝轨迹上所选点处的裂缝宽度；
⑤将计算的裂缝轨迹上全部所选点处的裂缝宽度值从小到大排序，去掉排序后裂缝宽度值序列两端各10％的宽度值，计算剩余80％的裂缝宽度值的平均值，并作为该条裂缝的宽度；
⑥由上述计算的单条裂缝宽度结合裂缝平均宽度定义(单位井段中裂缝轨迹宽度的平均值)即可计算出裂缝平均宽度值，计算公式为：
其中FVAi为第i条裂缝的宽度。
4)裂缝面孔率计算
由上述计算的裂缝长度和宽度参数结合裂缝面孔率的定义(单位井段上裂缝在井壁上所占面积与成像测井覆盖井壁的面积之比)即可计算出裂缝面孔率，计算公式为：
5.岩心裂缝参数提取
岩心裂缝参数提取可以采用上述声成像裂缝参数评价方法对岩心滚扫照片进行自动处理，
或者，利用图像处理软件对岩心滚扫照片进行人工处理。下面以Image J图像处理软件(Image J软件能够显示，编辑，分析，处理，保存，打印8位，16位，32位的图片，还能进行图片的区域和像素统计以及间距和角度计算，距离量取精度为1个像素的长度)为例介绍岩心裂缝参数提取的人工处理步骤。
1)打开Image J软件，主菜单界面如图14所示，在工具栏选择“File”中的“open”选项打开要处理的岩心滚扫彩色照片并将图片转换为灰度图像(便于后续的裂缝宽度量取等处理)
2)像素尺度刻度。全直径岩心柱近似为圆柱体，测量其直径进而可以计算出岩心柱的周长，故可利用岩心柱的实际周长和其滚扫照片中的周长之间的关系来刻度滚扫照片的像素尺度。具体操作步骤为：如图15所示，在工具栏选择“Straight Line”，然后沿着岩心照片横向画一条直线，其长度等于岩心照片的宽度；再选择“Analyze”中的“Set Scale”进行刻度，“Set Scale”界面中“Distance in pixels”为直线的像素长度，“Known distance”为岩心实际周长，“Pixel aspect ratio”为像素比率选择默认，“Unit of length”为“Known distance”的长度单位。
3)裂缝长度的测量。具体操作步骤为：在工具栏选择“Freehand Line”，在此模式下可以画任意曲线；将岩心照片放大到合适的尺寸，根据裂缝形态画出完整的裂缝轨迹(如图16所示)，然后按“M”键记录所画曲线的长度，即得裂缝的长度(如图17所示)。
4)裂缝宽度的测量。具体操作步骤为：在工具栏选择“Straight Line”，然后按“+”键放大照片中裂缝部分达到像素级，沿裂缝轨迹法线方向由裂缝的一个边界到另一个边界画一条直线段代表裂缝的宽度(如图18所示)，然后按“M”键记录所画线段的长度，即得裂缝的宽度(如图19所示)。
裂缝宽度测量过程中，裂缝边界的确定是裂缝宽度准确测量的关键。可以根据裂缝处像素灰度值的变化规律由以上成像裂缝宽度计算方法部分所述的裂缝边界判定准则确定裂缝的边界。另外，为了提高裂缝宽度测量的准确性，应选择裂缝比较稳定规则处进行测量，且应多测量几处求取平均值作为此裂缝的宽度。
5)裂缝平均宽度和裂缝面孔率的计算。
统计一定窗长内岩心照片上所有裂缝的宽度值，并计算平均值即得裂缝平均宽度，具体计算公式为：
其中，Wi为第i条裂缝的宽度，Nw为统计窗长内裂缝的条数。
根据上述计算的岩心裂缝的长度和宽度结合裂缝面孔率的定义可以计算出统计窗长内的裂缝面孔率，具体计算公式为：
其中，Li为第i条裂缝的长度，Lw为统计窗口长度，Rcore为岩心柱的半径。
6.超声成像裂缝参数刻度校正
岩心观察与描述能提供关于裂缝参数和发育情况、力学性质、充填特征和含油气性等的第一手资料，是最直接、最有效和最可靠的裂缝评价方式，故评价的裂缝参数精度较高。为了提高声成像裂缝参数评价的精度，有必要利用岩心裂缝参数对声成像裂缝参数进行刻度。
成像裂缝参数的刻度实质上是通过对比分析成像资料处理得到的裂缝参数和岩心观察与描述得到的裂缝参数，建立二者之间的刻度系数，以用于校正成像资料处理得到的裂缝参数，最终达到精细评价储层裂缝发育情况的目的。成像裂缝参数的刻度实现步骤如下：
(1)通过上述岩心观察与描述方法获取岩心裂缝参数，并利用岩心地面自然伽马值把岩心深度归位到常规GR曲线统一的深度尺度上；
(2)利用上述声成像测井资料处理方法处理成像原始数据得到成像图中的裂缝面孔率等参数，之后将处理得到的裂缝面孔率等曲线样值也归位到常规GR曲线统一的深度尺度上；
(3)对照岩心照片裂缝参数和声成像处理得到的成像裂缝参数，分析两者存在的关系，得到两者之间的刻度系数，从而达到用岩心照片裂缝面孔率刻度成像面孔率的目的；
统计X地区10口油基泥浆井的岩心裂缝参数和声成像资料处理得到的裂缝参数，并分别绘制岩心裂缝参数和成像裂缝参数的交会图如图21和图22所示，图中直线为线性拟合趋势线。由图21可以看出，整体上声成像裂缝宽度约为岩心裂缝宽度的11.345倍；由图22可以看出，整体上声成像裂缝面孔率约为岩心裂缝面孔率的 11.189倍；这两个系数就是根据岩心裂缝参数和声成像裂缝参数之间的关系得到的声成像裂缝参数刻度系数。实际声成像资料处理时利用这个系数对处理得到的裂缝参数进行刻度校正后可以更好地反映裂缝的真实状态和发育情况，具体声成像裂缝参数刻度公式为
其中，FVAH为刻度校正后的裂缝宽度，FVAHUXPL为声成像处理得到的裂缝宽度；FVPA为刻度校正后的裂缝面孔率，FVPAUXPL为声成像处理得到的裂缝面孔率。
其它未详细说明的部分均为现有技术。尽管上述实施例对本发明做出了详尽的描述，但它仅仅是本发明一部分实施例，而不是全部实施例，人们还可以根据本实施例在不经创造性前提下获得其他实施例，这些实施例都属于本发明保护范围。

Claims (4)

1.一种基于超声成像的油基泥浆井裂缝定量参数智能计算方法，其特征在于：包括以下步骤：
1)数据解编、加载与图像着色
利用超声波探头对井壁发射超声波脉冲并接收回波信号，得到由不同深度不同方位对应的回波幅度值或回波时间值组成的二维数据体，即为声成像测井数据；
2)预处理
在声成像图图像处理之前，对声成像测井数据进行加速度校正、方位校正和回波时间偏心校正的预处理：
a.方位校正是根据仪器方位曲线数据将各深度点的成像样点方位调整到统一的大地方位坐标系中，具体调整公式为：
其中，N为一周的扫描样点总数，s0[i]为一周扫描的第i个成像样点值，AZ为仪器方位角，％为求余运算，为向下取整运算，s1[i]为方位校正后的第i个成像样点值；
b.偏心校正处理步骤如下：
①对每一深度处回波时间数据进行统计，计算最大值MAX和最小值MIN；
②计算直方图：将统计的最小值到最大值区间段[MIN,MAX]细分为N等分，每一区间范围为：
[MIN+(i-1)*(MAX-MIN)/N,MIN+i*(MAX-MIN)/N]；
其中，N为时间成像图一周内的样点数，i为区间编号；
然后再对全井段所有数据进行分区统计，得到每一区间的频数；
③修正直方图极值
将直方图的两端各去掉2.5％的点后，得到直方图新的两端极值，最小值为MIN1，最大值MAX1；
④以仪器的1号极板方位作为Y轴；将位于[MIN1,MAX1]范围内的所有回波时间样点值TT投影到X和Y轴上，则有：
TTX[i]＝TT[i]*sin(α[i])
TTY[i]＝TT[i]*cos(α[i])
α[i]＝i/m*360
其中，i为处理深度点上位于[MIN1,MAX1]范围内的回波时间样点序号，m为处理深度点上位于[MIN1,MAX1]范围内的回波时间样点总数；则X和Y轴上的时间投影之和为：
偏心距：偏心方位：α＝arctan(XX/YY)；
⑤偏心距偏心方位滤波
a.在由各深度偏心距和偏心方位构成的数据序列中选取N个点作为窗长对目标井段所有深度点的偏心距和偏心方位进行滤波处理，其中，N为奇数，选3、5、7、9……；处理方法为：
计算出窗长内偏心距、偏心方位平均值和方差，则平均值为本点的偏心距和偏心方位；当偏心距的方差大于2.5时，则本点的偏心距为空值；当偏心方位的方差大于4.5时，则本点的偏心方位为空值；
b.对于偏心距、偏心方位为空值的点，以本点为中心上下各取4个有效点计算其平均值，用计算的平均值取代本空值点；
⑥偏心校正处理
根据仪器偏心时仪器轴中心点与井筒截面的相对位置几何关系，由仪器偏心时记录的回波时间计算出仪器居中时的理论回波时间，实现仪器偏心引起的测量误差校正；具体校正公式为：
其中，N为测量一周的扫描点数，i为参考方位回波时间样值在本深度点所有方位回波时间样值序列中的序号，
j为实际测量方位回波时间样值在本深度点所有方位回波时间样值序列中的序号，
TT[j]为实际测量方位回波时间样值，
α为实际测量方位与参考方位的夹角，
TT1[j]为仪器居中时实际测量方位的理论回波时间样值，
h为TT1[j]在参考方位的投影值，
d为TT[j]在参考方位的投影值与h之差；
3)滤波及裂缝自动拾取
a.提高声成像图的对比度；
b.提取裂缝响应图像；
利用拉普拉斯算子提取裂缝的二值图像；
c.提高裂缝响应轨迹的连续性
先采用自适应滤波进行去噪处理，再利用形态学处理方法补充裂缝响应轨迹间断；具体步骤如下：
(1)自适应滤波处理：
①将图像分成大小相等的子块；
②估计各子块的噪声方差大小；
③依据噪声方差确定滤波窗口，进行滤波；
(2)裂缝响应轨迹间断弥补处理
利用下述公式对B作关于原点的映射，再将其连续的在A上移动来进行膨胀运算进而实现裂缝响应轨迹间断弥补处理；
膨胀的运算符为A用B来膨胀写作其定义式为
其中，A为待处理的二值图像，B为卷积模板，x代表像素点，为B的反射；
d.裂缝自动拾取
采用Hough变换自动拾取裂缝响应轨迹；具体实现步骤为：
①计算正弦曲线的频率：ω＝1/2πr，其中r为井眼半径；
②通过正弦曲线的对称性确定基线位置y0；
③采用二维Hough变换提取该位置所对应的正弦曲线的幅度A和初相位参
④根据上述参数及正弦曲线函数即在成像图上标示出轨迹；其中，井壁成像的平面展开图上裂缝响应轨迹为一条正弦曲线，其方程表示为：
4)裂缝定量参数计算
拾取裂缝响应轨迹后，进一步定量计算裂缝的产状参数、裂缝密度、裂缝长度、裂缝宽度和裂缝面孔率，其中，
采用沿裂缝轨迹法线方向自动搜索的方式来统计计算裂缝的宽度，具体实现方式如下：
①在裂缝轨迹上等间距选择N1个点；
②以选择的点为起点沿裂缝轨迹法向方向搜索局部像素极小点MIN和极大点MAX1和MAX2；
③计算裂缝边界像素值，计算公式为：
其中，c1权重系数，取值范围为(0,1)，一般取0.4～0.8，由实验或经验确定；
④根据裂缝边界像素值EP1和EP2，沿裂缝轨迹法向方向上自动搜索裂缝边界像素点，并计算两个裂缝边界像素点之间的距离，即得裂缝轨迹上所选点处的裂缝宽度；
⑤将计算的裂缝轨迹上全部所选点处的裂缝宽度值从小到大排序，去掉排序后裂缝宽度值序列两端各10％的宽度值，计算剩余80％的裂缝宽度值的平均值，并作为该条裂缝的宽度；
⑥由上述计算的单条裂缝宽度结合裂缝平均宽度定义即计算出裂缝平均宽度值，计算公式为：
其中，FVAi为第i条裂缝的宽度；
5)岩心裂缝参数提取
采用上述声成像裂缝参数评价方法对岩心滚扫照片进行自动处理提取得到岩心裂缝参数；
或者，利用图像处理软件对岩心滚扫照片进行人工处理；其中，利用Image J图像处理软件提取得到岩心裂缝参数，步骤如下：
①将岩心滚扫彩色照片转换为灰度图像；
②测量全直径岩心柱直径进而计算出岩心柱的周长，故利用岩心柱的实际周长和其滚扫照片中的周长之间的关系来刻度滚扫照片的像素尺度；
③测量裂缝长度和裂缝宽度；
④裂缝平均宽度和裂缝面孔率的计算；
统计窗长内岩心照片上所有裂缝的宽度值，并计算平均值即得裂缝平均宽度，具体计算公式为：
其中，Wi为第i条裂缝的宽度，Nw为统计窗长内裂缝的条数；
根据上述计算的岩心裂缝的长度和宽度结合裂缝面孔率的定义，计算出统计窗长内的裂缝面孔率，具体计算公式为：
其中，Li为第i条裂缝的长度，Lw为统计窗口长度，Rcore为岩心柱的半径；
6)超声成像裂缝参数刻度校正
利用岩心裂缝参数对声成像裂缝参数进行刻度，步骤如下：
(1)通过上述岩心观察与描述方法获取岩心裂缝参数，并利用岩心地面自然伽马值把岩心深度归位到常规GR曲线统一的深度尺度上；
(2)利用上述声成像测井资料处理方法处理成像原始数据得到成像图中的裂缝面孔率等参数，之后将处理得到的裂缝面孔率等曲线样值也归位到常规GR曲线统一的深度尺度上；
(3)对照岩心照片裂缝参数和声成像处理得到的成像裂缝参数，分析两者存在的关系，得到两者之间的刻度系数，即得到岩心校正刻度系数(从而达到用岩心照片裂缝面孔率刻度成像面孔率的目的)；
(4)利用岩心校正刻度系数对超声成像资料处理得到的裂缝参数进行刻度校正即得到更符合实际的裂缝定量参数；具体校正公式为：
其中，COEH为裂缝宽度校正系数，COEA为裂缝面孔率校正系数，FVAH为刻度校正后的裂缝宽度，FVAHUXPL为声成像处理得到的裂缝宽度；FVPA为刻度校正后的裂缝面孔率，FVPAUXPL为声成像处理得到的裂缝面孔率。
2.根据权利要求1所述基于超声成像的油基泥浆井裂缝定量参数智能计算方法，其特征在于：所述步骤3)第a小步中，提高声成像图的对比度具体处理步骤如下：
①将预处理后的声成像图像分为大小相等的若干个子块；
②统计各个子块的局部直方图；
③计算子块内所有像素数平均分配到的各灰度级值A；
④设定一个剪切系数C，取值范围为[0,1]，则实际剪切值；
V＝A+C(Nx×Ny-A)
其中，Nx，Ny分别为子块水平方向和垂直方向的像素数目；
⑤利用剪切值，对局部直方图各个像素进行剪切，多余的像素重新分配到各直方图的各灰度级中，于是得到每个灰度级应该分到的像素数NACP为：NACP＝NCLIP/G
其中，G表示为子块中的灰度级，NCLIP表示已剪切的像素总数；
假设原直方图为H(i)，由此得到重新分配后的直方图CH(i)：
⑥对每个子区域的直方图进行均衡化处理，处理后裂缝响应轨迹更清晰。
3.根据权利要求1所述基于超声成像的油基泥浆井裂缝定量参数智能计算方法，其特征在于：所述步骤3)第b小步中，
①构造拉普拉斯算子：
②通过拉普拉斯算子提取裂缝响应轨迹的二值图像。
4.根据权利要求1所述基于超声成像的油基泥浆井裂缝定量参数智能计算方法，其特征在于：所述步骤4)中，
a.裂缝产状参数计算
根据拾取的裂缝轨迹和方位信息很容易确定裂缝的倾向和走向，裂缝倾角由如下公式计算得到：
其中，H为裂缝低点到高点的高度，R为井眼半径；
b.裂缝密度和长度计算
由裂缝响应轨迹拾取结果并结合裂缝密度定义(单位长度井壁上所见到的裂缝总条数)容易得到裂缝密度计算公式如下：
其中，Nw为统计窗长内裂缝条数，单位：条；Lw为窗口长度，单位：m；
由裂缝响应轨迹拾取结果利用线积分的方式容易计算出每条裂缝的长度，再根据裂缝长度的定义得到裂缝密度计算公式如下：
其中，Li为线积分方式计算的第i条裂缝的长度，单位：m；
c.裂缝面孔率计算
由上述计算的裂缝长度和宽度参数结合裂缝面孔率的定义，即计算出裂缝面孔率，计算公式为：
CN201910375453.5A 2019-05-07 2019-05-07 基于超声成像的油基泥浆井裂缝定量参数智能计算方法 Active CN110208859B (zh)

Priority Applications (1)

Application Number Priority Date Filing Date Title
CN201910375453.5A CN110208859B (zh) 2019-05-07 2019-05-07 基于超声成像的油基泥浆井裂缝定量参数智能计算方法

Applications Claiming Priority (1)

Application Number Priority Date Filing Date Title
CN201910375453.5A CN110208859B (zh) 2019-05-07 2019-05-07 基于超声成像的油基泥浆井裂缝定量参数智能计算方法

Publications (2)

Publication Number Publication Date
CN110208859A true CN110208859A (zh) 2019-09-06
CN110208859B CN110208859B (zh) 2020-10-02

Family

ID=67785579

Family Applications (1)

Application Number Title Priority Date Filing Date
CN201910375453.5A Active CN110208859B (zh) 2019-05-07 2019-05-07 基于超声成像的油基泥浆井裂缝定量参数智能计算方法

Country Status (1)

Country Link
CN (1) CN110208859B (zh)

Cited By (9)

* Cited by examiner, † Cited by third party
Publication number Priority date Publication date Assignee Title
CN111353396A (zh) * 2020-02-20 2020-06-30 江苏东印智慧工程技术研究院有限公司 一种基于SCSEOCUnet的混凝土裂缝分割方法
CN112835124A (zh) * 2021-03-10 2021-05-25 长江大学 基于成像测井和阵列声波测井资料的裂缝有效性评价方法
CN113050168A (zh) * 2021-03-10 2021-06-29 长江大学 基于阵列声波测井和声波远探测测井资料的裂缝有效性评价方法
CN113075748A (zh) * 2021-03-10 2021-07-06 中石化石油工程技术服务有限公司 基于成像测井和声波远探测测井资料的裂缝有效性评价方法
CN113392924A (zh) * 2021-06-29 2021-09-14 中海油田服务股份有限公司 声电成像测井图的识别方法及相关设备
CN113565495A (zh) * 2020-04-29 2021-10-29 中国石油天然气股份有限公司 储层有效性等级评价方法
CN113689453A (zh) * 2021-08-24 2021-11-23 中石化石油工程技术服务有限公司 一种测井图像裂缝自动识别方法、装置、设备及存储介质
CN113808190A (zh) * 2021-09-23 2021-12-17 西南石油大学 一种基于电成像测井图像的页岩纹层信息定量提取方法
CN115829990A (zh) * 2022-12-15 2023-03-21 吉林大学 一种基于成像测井图像处理的天然裂缝识别方法

Citations (6)

* Cited by examiner, † Cited by third party
Publication number Priority date Publication date Assignee Title
CN101430270A (zh) * 2007-11-08 2009-05-13 中国石油天然气股份有限公司 一种高蒙皂石砂岩孔隙度和渗透率的分析方法
CN101776769A (zh) * 2010-01-28 2010-07-14 中国科学院地质与地球物理研究所 应用石油包裹体进行油气运移追踪和古今油层判识的数据校正方法
US20120163123A1 (en) * 2010-12-27 2012-06-28 Baker Hughes Incorporated Stress in formations from azimuthal variation in acoustic and other properties
CN104076405A (zh) * 2014-07-18 2014-10-01 北京中石曙光科技有限公司 一种地层电阻率成像方法和微电阻率成像扫描仪
CN106014397A (zh) * 2016-05-13 2016-10-12 中国石油天然气集团公司 一种成像测井缝洞储层参数提取方法
CN106443802A (zh) * 2016-09-20 2017-02-22 中石化石油工程技术服务有限公司 一种基于形态学滤波的电成像储层缝洞体的定量表征方法及系统

Patent Citations (6)

* Cited by examiner, † Cited by third party
Publication number Priority date Publication date Assignee Title
CN101430270A (zh) * 2007-11-08 2009-05-13 中国石油天然气股份有限公司 一种高蒙皂石砂岩孔隙度和渗透率的分析方法
CN101776769A (zh) * 2010-01-28 2010-07-14 中国科学院地质与地球物理研究所 应用石油包裹体进行油气运移追踪和古今油层判识的数据校正方法
US20120163123A1 (en) * 2010-12-27 2012-06-28 Baker Hughes Incorporated Stress in formations from azimuthal variation in acoustic and other properties
CN104076405A (zh) * 2014-07-18 2014-10-01 北京中石曙光科技有限公司 一种地层电阻率成像方法和微电阻率成像扫描仪
CN106014397A (zh) * 2016-05-13 2016-10-12 中国石油天然气集团公司 一种成像测井缝洞储层参数提取方法
CN106443802A (zh) * 2016-09-20 2017-02-22 中石化石油工程技术服务有限公司 一种基于形态学滤波的电成像储层缝洞体的定量表征方法及系统

Non-Patent Citations (2)

* Cited by examiner, † Cited by third party
Title
姜龙，等: "油基泥浆井致密砂岩超声成像裂缝面孔率的校正刻度方法研究", 《当代化工》 *
袁龙，等: "库车前陆盆地强挤压应力条件下的测井电阻率校正方法", 《石油地球物理勘探》 *

Cited By (16)

* Cited by examiner, † Cited by third party
Publication number Priority date Publication date Assignee Title
CN111353396A (zh) * 2020-02-20 2020-06-30 江苏东印智慧工程技术研究院有限公司 一种基于SCSEOCUnet的混凝土裂缝分割方法
CN113565495A (zh) * 2020-04-29 2021-10-29 中国石油天然气股份有限公司 储层有效性等级评价方法
CN113565495B (zh) * 2020-04-29 2024-01-30 中国石油天然气股份有限公司 储层有效性等级评价方法
CN113050168B (zh) * 2021-03-10 2024-01-26 长江大学 基于阵列声波测井和声波远探测测井资料的裂缝有效性评价方法
CN112835124A (zh) * 2021-03-10 2021-05-25 长江大学 基于成像测井和阵列声波测井资料的裂缝有效性评价方法
CN113050168A (zh) * 2021-03-10 2021-06-29 长江大学 基于阵列声波测井和声波远探测测井资料的裂缝有效性评价方法
CN113075748A (zh) * 2021-03-10 2021-07-06 中石化石油工程技术服务有限公司 基于成像测井和声波远探测测井资料的裂缝有效性评价方法
CN113075748B (zh) * 2021-03-10 2024-04-16 中石化石油工程技术服务有限公司 基于成像测井和声波远探测测井资料的裂缝有效性评价方法
CN112835124B (zh) * 2021-03-10 2024-01-26 长江大学 基于成像测井和阵列声波测井资料的裂缝有效性评价方法
CN113392924A (zh) * 2021-06-29 2021-09-14 中海油田服务股份有限公司 声电成像测井图的识别方法及相关设备
CN113392924B (zh) * 2021-06-29 2023-05-02 中海油田服务股份有限公司 声电成像测井图的识别方法及相关设备
CN113689453A (zh) * 2021-08-24 2021-11-23 中石化石油工程技术服务有限公司 一种测井图像裂缝自动识别方法、装置、设备及存储介质
CN113808190B (zh) * 2021-09-23 2023-07-28 西南石油大学 一种基于电成像测井图像的页岩纹层信息定量提取方法
CN113808190A (zh) * 2021-09-23 2021-12-17 西南石油大学 一种基于电成像测井图像的页岩纹层信息定量提取方法
CN115829990B (zh) * 2022-12-15 2023-10-10 吉林大学 一种基于成像测井图像处理的天然裂缝识别方法
CN115829990A (zh) * 2022-12-15 2023-03-21 吉林大学 一种基于成像测井图像处理的天然裂缝识别方法

Also Published As

Publication number Publication date
CN110208859B (zh) 2020-10-02

Similar Documents

Publication Publication Date Title
CN110208859A (zh) 基于超声成像的油基泥浆井裂缝定量参数智能计算方法
Liu et al. Multiple-point simulation integrating wells, three-dimensional seismic data, and geology
CN109085663A (zh) 一种致密砂岩储层层理缝识别方法
CN109782360A (zh) 一种低阻覆盖区深部砂体探测方法
CN105093306A (zh) 一种地球物理勘探中储层自动解释与厚度求取方法
Wang et al. The automatic interpretation of structural plane parameters in borehole camera images from drilling engineering
CN108693562B (zh) 一种用电成像测井测定古流向的方法
CN114609675A (zh) 基于高频旋回对碳酸盐岩地层沉积微地貌的定量恢复方法
Koyan et al. 3D ground-penetrating radar attributes to generate classified facies models: A case study from a dune island
CN109031423B (zh) 基于梯度共生矩阵的叠前地震纹理分析方法
CN108375794A (zh) 基于对称观测的vsp缝洞绕射成像技术方法
You et al. Automatic facies classification from acoustic image logs using deep neural networks
CN111025397B (zh) 地震数据反射波与散射波联合求取深度域速度模型的方法
CN110471105B (zh) 一种层位约束下基于地震解析道的相干增强断层识别方法
CN109283577A (zh) 一种地震层位标定方法
CN113514884A (zh) 一种致密砂岩储层预测方法
CN110244359A (zh) 一种基于改进地震切片技术的异常体边缘检测计算方法
CN111045083A (zh) 一种储层含气性检测方法
CN112906242B (zh) 一种基于朴素贝叶斯法与邻近分类法相结合的地球物理建模方法
CN114428370B (zh) 一种分析古河道形态和内部结构的方法
CN113325474B (zh) 生物礁判别方法
CN113805246B (zh) 一种碳酸盐岩储层连通性评价图版及其生成方法和应用
Liu et al. Automated Detection Of Beddings From Borehole Images: A New Improvement And First Comparative Study
Li et al. Automatic identification of soil layer from borehole digital optical image and GPR based on color features
Nie et al. Research on Localization and three-dimensional Reconstruction of Borehole Wall Fractures Based on Ultrasonic Imaging Data

Legal Events

Date Code Title Description
PB01 Publication
PB01 Publication
SE01 Entry into force of request for substantive examination
SE01 Entry into force of request for substantive examination
GR01 Patent grant
GR01 Patent grant