张其虎，张 希，孙进辉，范天林

(1.昆明理工大学 国土资源工程学院，昆明 650093；2.玉溪矿业有限公司，云南 玉溪 653405；3.云南铜业股份有限公司矿山研究院，云南 玉溪 653405)

凿岩爆破是地下矿山掘进过程中的重要施工手段，但其产生的爆轰波不可避免地对地表临近建筑物的承载能力及稳定性带来不利影响，矿山爆破的声音、震动、爆破飞屑物等会给周围环境带来较大的影响，甚至会对周边居民的日常生活造成严重的困扰[1-3]；而且爆破作业受到技术、资金等建设条件的限制以及环境、地质方面的更高要求，所以人们对爆破技术的要求越加苛刻[4-6]。又因为爆轰波是在水平面位置处以波的形式做往复运动，本身持续时间短以及波形传播介质和赋存条件复杂多变，从而导致爆破地震的复杂性[7]。众多研究者分别从理论研究、试验验证、破坏准则等方面研究爆破振动规律。雷明等[8]基于对爆轰波峰值速度的理论研究，分析了井下爆炸产生的应力波对围岩区域的影响因素，从而对岩巷稳定性进行了评估。张家斌[9]利用统筹理论建立了东川岭矿的爆轰波运动轨迹模型，以此来研究振动波的传播规律。李宇星等[10]从爆破振动三要素的角度出发，基于爆破振动及波形传播规律对周围建筑物进行稳定性评估分析。邱贤阳等[11]将实测数据和数值模拟结合，得出具有现实意义的理论数值依据，用作周围建筑物承载能力和稳定性的评判标准。梁升等[12]运用ANSYS/LS-DYNA数值模拟手段，分析爆破振动规律，为建构筑物提供安全稳定性的理论依据。结合国内外学者的研究成果[13-17]，笔者以金沙矿业官房矿段地下爆破作业为依托，对现场爆破作业过程中的实测数据信号进行采集，运用EMD-HHT分析法从爆破能量分布规律、爆破振动波形的频谱特征、地震波衰减分析理论和爆破振动速度等方面进行研究分析。在爆破振动波形图中可以得到质点的峰值振速及主振频率等参量的数值大小，两者都是与爆破振动效应相关的重要因素，峰值振速反映爆破荷载的大小，主频率反映振动物的动态响应情况特征。基于对上述物理量对爆破振动波形的研究分析得到的相关结论，可为今后的地下浅孔爆破提供更为精确的研究思路和方法，指导地下矿山浅孔爆破作业。

1 EMD-HHT分析法

EMD-HHT分析法是根据实测数据本身的时辰维度进行频谱分解，不再需要固定的先遣波形基底，利用Matlab对爆破振动波形数据信号进行EMD分解，通过编辑相应的程序语言在Matlab软件上对波形数据运行处理，运行后会得到有限个不同的IMF分量，对去除噪音分量的IMF进行重构，将主要IMF分量进行Hilbert变换，得到能量分布图谱即Hilbert谱，再进行Hilbert变换，最后提取出不同IMF分量的包络幅值，并将频率幅值分解为时频谱、边际谱、瞬时能量谱等，用以识别实际微差延期时间[18-20]。Hilbert时频图的表达形式为：

(1)

式中：Re—残余分量r的实部；i—IMF的个数，i=1，2，3，…，n；ai、wi—常数；ai(t)—幅值函数；ωi(t)—频率函数。

将H(ω,t)对时间积分，就得到Hilbert边际谱：

(2)

再对Hilbert谱h(ω)的幅值利用积分计算方式对时间这个参量进行积分，最后通过求和可得Hilbert瞬时能量谱Es，其表达形式为：

(3)

对Es在频域内积分可得到信号的Hilbert能量，其表达式为：

(4)

采用EMD方法对给定信号S(t)进行分解，最终原始信号S(t)可分解为若干IMF分量，每个IMF分量均为时间的函数ci(t)和一个余项rn(t)的和，即：

(5)

从式(5)中可以看出，将时间的函数ci(t)进行i=1到n的求和，即可得IMF分量c(t)，并将分量c(t)再一次进行Hilbert变换，即：

(6)

式中：PV—柯西主值，因此可以得到构造解析信号z(t)：

z(t)=c(t)+jH[c(t)]=a(t)ejΦ(T)

(7)

式中：a(t)—幅值函数；j—相关系数。a(t)为z(t)的幅值，也称为信号包络，其幅值函数为：

(8)

Φ(t)为相位函数：

(9)

2 爆破振动监测与信号分析

2.1 工程地质及爆破参数

PD207矿段为构造角砾岩型铅锌矿，矿体厚度约为1.0 m，其岩石物理力学性质如表1所示。选择连续两次的爆破作业振动监测结果进行分析。因矿体厚度较薄，采用浅孔爆破，采场回采作业孔深为2.5 m，炮孔直径40 mm，孔距为1.0 m，排距0.8 m，炮孔不做填塞，炮孔采用梅花形布孔。第一次爆破总药量29.2 kg，最大单响药量4.8 kg；第二次爆破总药量24 kg，最大单响药量4.8 kg。采用磁电雷管加毫秒延期导爆管雷管的起爆方式，共布置了5个爆破振动监测点，分别放置4套TC-4850测振仪和1套BMView-8016测振仪。

表1 岩石物理力学性质

2.2 基于EMD-HHT法的信号处理及分析

选取BMView-8016所测得的数据来进行EMD-HHT分析，在分析过程中通常截取爆破振动波的垂直振动速度波形整体中的最大峰值振动段波形作为研究的对象，因为最大峰值处对应的最大瞬时能量即为最大能量荷载区域[21-23]。分别截取两次典型的垂直方向通道的峰值振动段波形的两组信号作为研究对象，记做A和B。

1)EMD-HHT分解

将已选择好的两次爆破作业振动监测结果的数据运用配套软件将波形数据导出，在导出的波形图中选择波形峰值最大处，并将波形最大处所对应的数据转为Excel文件格式导出，除去波形前后的数据和负延时所对应的信号数据。进而运用Matlab2018b数据处理软件对导出的数据进行处理，然后运行出数据波形，再对得出的波形进行研究分析。从波形图中可以得到12个分量以及残余分量r，不同的IMF分量包含了不同的信号特征，从图中可以看出，信号频率是由高频到低频顺序排列，在这些排列规律中可以发现A和B的各个IMF分量，如图1和图2所示。

图1 A的IMF分量Fig.1 IMF component of A

图2 B的IMF分量Fig.2 IMF component of B

从图1和图2中分析得出IMF1是速度和时间的关系，随时间的变化速度变化频率非常大，可以看出速度随时间的变化较为明显，从振动频率幅值可以看出频率变化最大，但占有的能量最小，从中可以分析出IMF1是进行振动监测时引入的高频振动幅值。IMF2是频率和振幅的关系，随频率的增加，振幅在逐渐增大，但将振幅值扩大时可以看出随频率增加振幅从变大到变小，中间会出现峰值振幅。IMF3又是速度随时间的变化，将速度值变小时可以看到0.6 s左右时会有速度变化，将速度变大会看到速度为正弦波形变化。IMF4可以看出随频率的增加振幅逐渐趋于零。IMF5可以看出随时间的增加速度在不断变化。在上面的IMF信号分量中，IMF2～IMF5是进行研究的主要对象，在这其中含有信号的能量占有总能量中的绝大部分，属于爆破振动信号中的优势振动频率。而最后一个分量IMF6则可能是仪器的飘零或信号本身最终的趋势结果，没有实际研究的价值和意义，可以忽略。因重构后的信号波形更符合爆破振动波的实际情况，将A和B的信号数据进行EMD重构，A和B的EMD重构信号如图3和图4所示。

图3 A的EMD重构信号Fig.3 EMD reconstruction signal of A

图4 B的EMD重构信号Fig.4 EMD reconstruction signal of B

从EMD重构信号中可以得出速度随时间的变化规律，A的速度约为B的速度的两倍，A在0.33 s左右时候速度达到最大0.035 cm/s，B在0.6 s时候速度达到最大0.070 cm/s。

2)边际谱分析

对经过EMD重构后的A和B数据信号，再进行边际谱分析。从数据信号的边际谱中可以看出能量和频率的相互关系，能量在频率带上的集中程度和爆破振动信号在不同频率上的变化幅度[24]。A和B信号的边际谱如图5和图6所示。

图5 A的边际谱Fig.5 Marginal spectrum of A

图6 B的边际谱Fig.6 Marginal spectrum of B

从图5和图6可以看出，信号A的主振频率主要分布在0～100 Hz频带范围内，最大能量约为3.3×10-6，80 Hz以下的低频范围内分布集中，包含多个频带，主振带分别为30和40 Hz。截取的信号B的频率主要分布在0～80 Hz频带范围内，最大能量约为1.0×10-5，60 Hz以下的低频范围内分布集中，包含多个频带，主振带分别为30和45 Hz。信号A和B的频率主要分布区域相差不大，但是在峰值处的能量B约为A的三倍。

3)瞬时能量谱分析

边际谱反映的是能量随频率的变化规律，但瞬时能量谱能够直观地看出瞬时能量随时间的分布规律[25]。截取的A和B信号的瞬时能量谱如图7和图8所示。

图7 A的瞬时能量谱Fig.7 Instantaneous energy spectrum of A

图8 B的瞬时能量谱Fig.8 Instantaneous energy spectrum of B

从图7和图8中分析得出，信号A和B都是由导爆管起爆引起的结果，在图中波形的峰值处表示瞬时能量的极值，A的波形范围在0.1～0.5 s区域内，B的波形范围在0.2～0.8 s区域内，即表示A和B的瞬时能量分别集中在0.1～0.5 s范围内和0.2～0.8 s范围内，A信号达到峰值瞬时能量在0.32 s左右，峰值瞬时能量约为1.38×10-3；B信号达到峰值瞬时能量在0.58 s左右，峰值瞬时能量约为5.00×10-3。从数据中可以看出，A信号的峰值瞬时能量约为B信号的峰值瞬时能量的3.6倍，从装药量入手，知A与B的单孔最大装药量一致，但测得的振动数据中，瞬时能量差异较大，从正交分析角度入手，说明爆破总用药量也是影响最大瞬时能量的因素之一，振动信号的瞬时能量会随着总装药量的增加而增大。瞬时能量代表着爆破的最大荷载，当瞬时能量过大，则爆破产生的爆破荷载也很大，一定强度范围内的爆破振动，会让爆破区周围建(构)筑物的安全性及稳定性得不到保障，很容易超越建筑本身承受能力，会出现一定的损伤甚至倒塌、爆区地表边坡的滑移及地下巷道冒落等各种损伤破坏问题。为了确保周围建(构)筑物安全，控制爆破地震效应和改善岩石破碎效果，采取孔内间隔装药、控制炸药单耗、降低总装药量、逐孔起爆、预裂爆破等技术以满足降震要求和保证建构筑物的安全。使用高精度高强度数码电子雷管起爆，精确控制延期时间，采用精确延时逐孔起爆的爆破技术手段，再结合上述的降振措施，在工人爆破技术都较高的操作水平下，爆破地震效应将会得到控制，岩石破碎效果将会得到改善，控制爆破技术又上一个新的水平。

4)三维谱分析

三维谱是基于三轴分析理论，用以研究爆破振动信号的一种波形分析理论，基于短时傅里叶变换形成的时间、频率和瞬时能量三维谱图，采用数字图像处理的方法提取信号的特征，可以从空间直接反映能量随时间和频率变化的分布情况。以波形图的颜色深浅作为判别的标准，颜色越鲜艳明亮，则代表能量越高，其爆破荷载也越大。在三维谱中，时间、频率分别代表着X和Y轴，瞬时能量代表Z轴，三者相互作用，以空间立体图反映了A的能量主要集中在0～0.8 s、0～50 Hz范围内，B的能量主要集中在0～1.0 s、0～50 Hz范围内，且可以通过颜色的深浅很好地观察出来，以判别能量大小范围。三维谱得到的数据结果与边际谱和瞬时能量谱分析的结果相一致，说明此三者可以相互检验，互相验证结果，以使得检测分析更加精确。A和B信号的三维谱如图9和图10所示。

图9 A的三维谱Fig.9 Three-dimensional spectrum of A

图10 B的三维谱Fig.10 Three-dimensional spectrum of B

5)Hilbert能量谱

对IMF分量最后进行等高线能量图分析，得到希尔伯特能量谱，该图谱也基于颜色深度和集中度情况来区分，可以进行更加精确的时频分析，同时以时间-频率的方式展现振动信号原始的能量分布，避免了窗函数导致的能量损失以及小波基选取不合适而造成的干扰。信号A和B的等高线能量图如图11和图12所示。

图11 A的等高线能量图Fig.11 Contour energy map of A

图12 B的等高线能量图Fig.12 Contour energy map of B

图表明了不同IMF分量以频率-时间-振幅的分布方式，可以看出，A振动能量主要分布在0～100 Hz以内，并且集中在50 Hz以下的低频区域，其余少部分能量分布在100 Hz以上的频段，B振动能量主要分布在0～60 Hz以内，并且集中在50 Hz以下的低频区域，其余少部分能量分布在60 Hz以上的频段，并且随着质点振速的增加，爆破振动信号频带会朝低频化发展。

3 结论

1)采用EMD-HHT分析法从爆破能量分布规律、爆破振动波形的频谱特征、地震波衰减分析理论和爆破振动速度等方面进行研究分析。装药量的变化会影响振动频率发生变化，药量的增加会导致频率从高频往低频发展，与此同时，最大爆破荷载也会随着药量的变化而变化。

2)在爆破振动波形图中可以得到质点的峰值振速及主频率等参量大小，两者都是与爆破振动效应相关重要因素，峰值振速反映爆破荷载的大小，主频率反映振动物的动态响应情况特征。基于上述物理量对爆破振动波形的研究分析得到的相关结论，可为今后的地下浅孔爆破提供更为精确的研究思路和方法。

3)基于边际谱、瞬时能量谱和三维谱的研究分析，在考虑浅孔爆破振动对建构筑物安全性及稳定性影响的同时，应该重点考虑10～50 Hz低频带，浅孔爆破其主振频率一般为40～100 Hz。