基于能量衰减系数的弹性波CT反演

(水利部长江勘测技术研究所，武汉 430011)

1 研究背景

弹性波计算机断层(Computed Tomography，CT)成像技术是一项成熟的地球物理探测方法，在地下岩溶、孤石勘察、采空区勘察等多方面得到应用。弹性波CT反演可分为速度反演和衰减系数反演。

弹性波CT速度反演测试值稳定，受测试条件影响小，对小缺陷敏感性较小；弹性波CT衰减反演受测试条件影响较大，对小缺陷敏感[1]。

目前大部分弹性波CT工程采用速度反演，速度反演以直达波初至时刻反演，生成速度剖面，划分异常。此种反演方法在高速体中寻找低速异常时效果较好，如岩溶，破碎带的测试；但在低速体中需找高速体时，则存在异常体分辨率不高，异常界限不清晰，如地下孤石的探测。

弹性波CT衰减系数反演虽然理论成熟，但一直没有得到广泛应用。主要原因有两点：

(1)震源每次激发能量无法精确一致。弹性波CT目前有多种发射震源，如压电陶瓷、超磁材料、电火花、炸药等。虽然压电陶瓷、超磁材料、电火花可通过发射电压来控制发射能量；炸药可以通过质量来控制爆炸能量，但均无法像电磁波CT那样准确控制发射能量。

(2)检波器接到弹性波能量无法精确计算。震源激发的弹性波到达各个检波器，采样得到每道波形数据。如果按全部波形幅值计算弹性波能量，则需要长时间采样。如果以波形最大振幅计算弹性波能量，则存在最大振幅可能为横波或面波，或者为其他转换波的问题。

本文尝试借鉴电磁波CT能量衰减方式进行弹性波CT衰减系数反演。电磁波CT是在固定电磁波发射频率后，固定发射强度，确保电磁波发射能量固定，接收天线收到经过介质吸收后的电磁波强度，进行反演后计算出介质对电磁波的吸收系数。电磁波CT目前采用能量衰减进行反演形成电磁波吸收剖面，划分异常。

2 计算方法

目前弹性波CT处理数据采用读取每炮每道记录的波形初至时刻，代入方程组反演出剖面网格速度。能量衰减系数反演则读取每道记录波形首波振幅值(见图1)。

图1 测试原理及衰减系数反演首波振幅

经过距离校正归一化，衰减系数拟合，计算出各道衰减系数值，最后进行反演计算。

2.1 首波振幅距离校正归一化

弹性波能量与振幅的平方呈正比关系，介质单位面积能量与震源距离平方呈反比。根据以上弹性波能量与振幅、震源距离的关系，将每炮记录首波振幅能量按式(1)进行首波振幅距离校正归一化。

(1)

将记录到的各道首波振幅及炮检距按式(1)计算，得到图2。由图2可知，波形首波振幅由大到小的前5位通道号依次为：2，1，15，3，12。经过距离校正归一化后，能量由大到小的前5位通道号依次为： 15，12，14， 2， 1。

图2 首波振幅及归一化能量

2.2 能量衰减系数拟合

经研究表明，小振幅的弹性波，传播时的衰减呈指数形式[2-3]。

将单炮数据记录首波振幅能量距离校正归一化值En、发射震源与各道检波器距离Ln按弹性波能量衰减公式进行最小二乘法拟合，得到每炮记录衰减系数λ，即

E=E0exp(-λL) 。

(2)

式中：E为震源振动传播到检测点的能量值；E0为震源能量值；λ为能量衰减系数；L为震源到检波器距离。

介质能量衰减系数大，表明介质对弹性波能量吸收强，介质软弱、破碎；介质能量衰减系数小，表明介质对弹性波能量吸收弱，介质坚硬、完整。

将首波振幅校正归一化及拟合值按式(3)求取拟合相对误差(见图3)。

(3)

图3 不同通道的归一化能量、拟合能量、拟合相对误差

2.3 计算每道衰减值

每炮记录衰减系数、每道拟合相对误差代入式(4)计算每炮每道衰减值。

Λn=λ/(1+εn) 。

(4)

式中Λn为第n道记录能量衰减值，单炮记录每道衰减系数见图4。

图4 单炮记录每道衰减系数

图5 速度和衰减系数反演成果对比

通过以上计算，可计算出每炮记录各道的衰减值。

2.4 反演计算

弹性波剖面共M炮

记录，每炮记录为N道波形数据，可以按计算出衰减矩阵Λm×n、炮检距Lm×n(m为炮编号，1≤m≤M；n为检波器道编号，1≤n≤N)。

将剖面分解为J×K正交网格化(竖直方向分解为J行，水平方向分解为K列)，得到剖面衰减网格矩阵Sj×k(1≤j≤J；1≤k≤K)。

CT反演计算为解矩阵方程组，即

LS=Λ。

(5)

式(5)为大型稀疏矩阵，求解方法较多，目前应用较多的迭代方法有代数重建法(Algebraic Reconstruction Technique，ART)、联合迭代重建法(Simultaneous Iterative Reconstruction Technique，SIRT)、最小二乘正交分解法(Least Squares QR decomposition，LSQR)。本文采取基于最小二乘法准则的联合迭代重建法[4]。

3 工程实例

在武汉地铁机场线孤石弹性波CT勘察中，采用直达弹性波初至时刻反演，生成弹性波速度剖面。根据孤石波速高于周围粘土波速，进行孤石的判别。对发现的孤石进行地下破碎处理后，盾构机顺利通过，目前武汉地铁机场线已建成运行。

在后期进行技术总结发现，弹性波CT速度反演探测孤石总体效果很好，但还是存在部分问题：

(1)实际工程中，孤石与周围粘土波速差异大，为波速突变。而在速度反演剖面中，孤石与周围黏土的波速为渐变，不利于孤石边界的划分。

(2)由于尺度效应，射线穿过黏土、孤石，反演出现黏土波速提高，孤石波速降低的情况。

速度和衰减系数反演成果对比如图5所示。 图5中GP15和GP17为2个钻孔，图中反映了2个钻孔之间区域的反演成果。

图5(a)中的速度反演，即存在上述2个问题。波速反演最小波速1 227 m/s，最大波速2 240 m/s，平均波速1 654 m/s，波速相对极差0.61；黏土区波速1 200～2 000 m/s，孤石区波速2 000～2 400 m/s。

极差是样本数据的最大值与最小值的差，可以评价样本数据离散程度；相对极差是极差除以样本数据的平均值。弹性波CT成果一般以等值线图形式体现，等值线背景值即样本平均值，等值线有效范围与精度与相对极差成正相关关系。相对极差大则等值线精度高，有利于异常识别及异常范围界定。

鉴于波速反演中问题，对剖面数据进行弹性波CT衰减系数反演，如图5(b)。在衰减系数反演中，最小衰减系数0.12，最大衰减系数0.74，平均衰减系数0.27，衰减系数相对极差2.30，参数的相对差异增大。在波速反演中2处波速高值孤石异常区在衰减反演中同样部位出现衰减低值孤石异常区，异常位置、形态对应较好。在衰减系数反演图中，在剖面中右部，孔深8～11 m，存在倾斜椭圆状衰减低值异常区，在波速反演中，该部位无异常，经钻孔验证为破碎孤石区，说明衰减系数反演对于小缺陷的灵敏性高于速度反演。

4 衰减系数反演与波速反演差异

通过对比和验证，弹性波CT能量衰减系数反演与弹性波CT速度反演均能较好适用于黏土中孤石探测。弹性波CT能量衰减系数反演与波速反演有以下3个不同之处。

(1)对检波器一致性要求提高。在弹性波CT速度反演中需读取首波初至时刻，对检波器时间一致性有严格要求。在弹性波CT衰减系数反演中需读取直达波第一个波峰幅值，对检波器的波形一致性也提出严格要求。需要在工作前对检波器进行幅值一致性试验，对于幅值偏差超过允许值的检波器需要进行幅值校正。

(2)弹性波能量的确定。检波器接收到从震源激发出弹性波的能量应该为全部波形的能量或最大振幅的能量。检波器接收到的波形包括直达波及转换波等，波形成分复杂，各种波难以分离，全波波形能量或最大振幅的能量来衡量波形能量均有不恰当之处，故在衰减系数反演中选取直达波第一个波峰幅值来计算弹性波能量。波速反演只与直达波初至时间有关，与波形能量无关，故无此要求。

(3)弹性波路径。弹性波在非均匀介质传播时，直达波在震源与检波器间路径并不是直线传播，而是遵循费马原理，路径为弯曲射线[5]，故目前弹性波CT速度反演多采用弯曲射线追踪[6-8]。本文中弹性波CT 衰减系数反演，能量计算采取直达波第1个波峰幅值，故反演路径也应该为弯曲射线。直达波能量无法计算路径，本文衰减系数反演中按直线传播反演。

能量衰减系数反演如能与弹性波CT速度反演的弯曲路径相结合，应该能够进一步提高反演结果精度。

5 结 论

(1)本文通过将单炮记录各道的直达波第一个波峰幅值经过距离校正后，进行归一化处理，拟合震源到检波器间扇形区域的整体衰减系数，根据各道拟合的相对误差，计算出震源到各道检波器间直线衰减系数，进行弹性波CT衰减系数反演，得到衰减系数剖面，作为划分异常依据，取得较好效果。

(2)衰减系数反演对于小缺陷的灵敏性高于速度反演。

(3)衰减系数反演得到的衰减系数相对极差大于速度反演得到的速度相对极差，更利于异常划分。