如何控制爆破振动，减少爆破施工对周围环境和人员的危害，一直是爆破工程安全技术人员研究的重要课题^[1-2]。目前，爆破领域衡量爆破振动强度的物理量主要是速度和加速度^[3]。虽然理论上已经有许多测试爆破振动速度和加速度的仪器与设备，但从实施难易程度看，爆破振动速度传感器的测量容易受到环境或传感器安装条件的限制，有时无法直接测量。和爆破振动速度传感器相比，爆破振动加速度传感器的安装要容易，这是因为爆破振动加速度传感器具有体积小、重量轻、安装方便、灵敏度高等优点。因此，在实际工程中爆破振动加速度测试得到了广泛的应用^[4-5]。

理论上，对监测到的实测爆破振动加速度信号进行一次数值积分即可得到爆破振动速度信号。对实测爆破振动加速度信号进行数值积分的算法，按数据处理过程可以分为频域积分法和时域积分法两种。频域积分法主要是利用傅里叶正、反变换，积分在频域内以傅里叶系数来表示，这就不可避免地存在由于频率截断误差而引起的能量泄漏^[6]。时域积分法通常采用梯形公式，受零点漂移或灵敏度漂移影响较大，需要对爆破振动信号进行剔除趋势项或滤波处理，具有代表性的算法主要有差分法、低通滤波法、最小二乘法、小波变换法或固有模态分解法等^[7-9]。然而，采用上述方法仍不能将爆破振动信号中的时域积分趋势项予以去除，残余的微小误差在时域积分过程中也会积累而逐渐放大，甚至会使爆破振动信号发生非线性失真，而时域积分算法本身也会产生趋势项，从而严重降低爆破振动信号频谱分析的精度。

为了更精确地提取爆破振动信号中振动峰值和各频率区间内能量等重要特征，必须对爆破振动加速度信号时域积分中的趋势项予以去除。本文中以此为切入点，在模式自适应小波原理的基础上^[10]，通过对实测爆破振动加速度信号进行梯形数值积分，构建将时域积分后的爆破振动速度信号作为小波基函数的方法，用构造模式自适应爆破振动小波基函数的方法来剔除时域积分后爆破振动速度信号中的趋势项，从而提高爆破振动信号时域和频域分辨率。

1.   模式自适应小波理论

1.1   连续小波变换

设函数ψ(t)∈L²(R)，其连续傅里叶变换为$\hat{\psi }(\omega )$^[11]。当$\hat{\psi }(\omega )$满足下列容许性条件：

时，则称函数ψ(t)为一个基本小波。若将该基本小波ψ(t)通过伸缩与平移等运算处理后，就可得到一组连续小波序列，称之为连续小波基函数，即：

式中：a为尺度因子(a≠0)，b为平移因子。

对于任意平方可积的函数或信号f(t)，即f(t)∈L²(R)，其关于小波基函数ψ_{a, b}(t)的连续小波变换为：

式中：＜f(t), ψ_{a, b}(t)＞为函数或信号f(t)与小波基函数ψ_{a, b}(t)的内积，$\overline{\psi \left( \frac{t-b}{a} \right)}$为$\psi \left( \frac{t-b}{a} \right)$的共轭函数。

由式(2)可以看出，连续小波变换的实质是将函数或信号f(t)以小波基函数序列ψ_{a, b}(t)为基底进行展开，由此可以求得函数或信号f(t)在小波基函数序列ψ_{a, b}(t)上的投影，即采用不同时间、不同尺度的小波基函数序列对被分析函数或信号f(t)进行相似性分析。若要使被分析函数或信号f(t)的时频特征在小波变换系数上更好地呈现，则需要将被分析函数或信号f(t)在小波基函数序列ψ_{a, b}(t)上的投影系数尽可能大。换一句话说，小波基函数序列ψ_{a, b}(t)的时域波形应与被分析函数或信号f(t)的特征波形具有良好的局部形状相似性。

1.2   模式自适应连续小波

模式自适应连续小波是利用最佳平方逼近方法来设计一个与给定函数(或信号波形)的局部特征相似度高的小波^[12]。具体的实现方法有多项式逼近法和常数正交函数空间投影法2种。若给定的函数(或信号波形)比较简单，则可以选择多项式逼近法；若给定的函数(或信号波形)比较复杂，则可以选择常数正交函数空间投影法。参数规则度定义了模式自适应连续小波在紧支撑区间上的边界约束，可以是“无”或“连续”或“可微”。若对给定函数(或信号波形)进行多项式逼近，参数规则度选择连续约束条件时，多项式拟合阶次应大于等于3；参数规则度选择可微约束条件时，多项式拟合阶次应大于等于5。

基于上述构造模式自适应连续小波的思想^[13-14]，其过程大致可以分为3个步骤：

(1) 根据给定函数(或信号波形)的检测模式x(t)，设置逼近方法和边界条件，构建模式自适应连续小波ψ_x(t)，不失一般性，设模式自适应连续小波的紧支撑区间为[0, 1]，则有$\int_{0}^{1}{{{\psi }_{x}}(t)\text{d}t=0}$以及${{\left\| {{\psi }_{x}}(t) \right\|}_{2}}=1$；

(2) 检测给定函数(或信号波形)上所有的虚警，对任意b和a>0的情况，搜索函数(或信号波形)的局部极大值点；

(3) 检测和丢弃给定函数(或信号波形)上所有错误的虚警，这一规则必须应用于确定每个虚警是否错误。

2.   爆破振动信号趋势项的去除

2.1   爆破振动加速度信号的测试

爆破振动加速度信号测试时，除了考虑选择性能优良的仪器、设备与测试系统外，还需要根据爆破振动加速度信号测试的目的来标定测振系统与设置测振仪量程、采样频率、触发与延迟时间等参数^[15]。

通过对某个矿床进行大量的爆破振动加速度测试，从中抽取1个测点所监测到的爆破振动加速度信号进行分析，实测爆破振动加速度信号a时程曲线及其相应的傅里叶功率谱密度P见图 1。

图  1  爆破振动加速度信号及功率谱密度

Figure  1.  Time history and power spectrum density of blast vibration acceleration signal

下载: 全尺寸图片 幻灯片

2.2   爆破振动加速度信号的积分变换

将爆破振动加速度信号转化为速度信号时，常采用梯形求积的数值积分法。设实测爆破振动加速度信号为{x(k)} (k=0, 1, 2, …N; N为爆破振动速度信号序列中包含离散采样点的个数)，则可以根据梯形求积公式，求得爆破振动速度信号为

根据式(3)，采用MATLAB编程语言对实测爆破振动加速度信号进行梯形数值积分，结果如图 2(a)所示。图 2(b)为时域积分后爆破振动速度信号的傅里叶功率谱密度。从图 2可以看出，爆破振动速度信号时程曲线存在明显的零点漂移现象；爆破振动速度信号功率谱密度中的优势频率主要集中在0~100 Hz频率范围内，再对功率谱密度低频部分进行局部放大，可以清楚地观察到0~5 Hz频率范围内具有较大幅值的低频直流分量。因此，爆破振动速度信号中趋势项的存在严重影响了爆破振动速度信号时程曲线波动特征的准确性和功率谱分析的分辨率，必须予以去除。

图  2  时域积分后的爆破振动速度信号及功率谱密度

Figure  2.  Time history and power spectrum density of blast vibration velocity signal after time integration

下载: 全尺寸图片 幻灯片

2.3   去均值化

爆破振动速度信号的均值是爆破振动速度信号离散值{y(k)} (k=0, 1, 2, …N)在整个时间坐标轴上的积分平均，其物理含义为该爆破振动速度信号变化的中心趋势，或称为零点漂移。由测振仪或传感器的低频直流分量引起，因此，可以通过去均值公式来实现：

式中：y′(k)为去均值后的爆破振动速度信号；y(k)为爆破振动速度信号；μ_y为爆破振动速度信号的算术平均值。

2.4   模式自适应连续小波法去除爆破振动信号趋势项

利用MATLAB软件平台上的小波工具箱(Wavelet Toolbox)，对图 2(a)所示的爆破振动速度信号进行去均值化处理，然后对其进行模式自适应波形匹配，构造出满足小波容许性条件的模式自适应小波基，用构造模式自适应小波基的方法剔除去均值化后爆破振动速度信号中的趋势项。图 3(a)和3(b)分别为去除趋势项后爆破振动速度信号时程曲线及其相应的傅里叶功率谱密度。从图 3(a)可以看出，去除趋势项后爆破振动速度信号第一个记录点和最后一个记录点的速度值均为零，与图 2(a)中的时程曲线进行比较，模式自适应连续小波法已成功去除了去均值化后爆破振动速度信号中趋势项的影响。对比分析图 3(b)与图 2(b)可以看出，模式自适应连续小波法已成功去除了低频直流分量的干扰，同时可知去除趋势项后的爆破振动速度信号功率谱(见图 3(b))与实测爆破振动加速度信号功率谱(见图 1(b))的优势频率均为92 Hz，且两者的相关系数高达0.953 8，由此可以表明用模式自适应连续小波法去除爆破振动信号中的趋势项是完全可行的。

图  3  去除趋势项后的爆破振动速度信号及功率谱密度

Figure  3.  Time history and power spectrum density of blast vibration velocity signal after removal trend

下载: 全尺寸图片 幻灯片

3.   爆破振动信号各频率区间的能量分布

3.1   爆破振动信号各频率区间的划分

采用小波包分析技术对爆破振动信号进行各频率区间分解，分解的有效深度实际上应根据爆破振动信号分析的目的以及所采用爆破振动记录仪的最小工作频率而定。本次现场实验所选用的采样频率为1 024 Hz，由Shannon采样定理可知，则其Nyquist频率为512 Hz。由于爆破振动记录仪的最小工作频率为5 Hz，根据小波包分解原理以及试验分析的需要，可以将各测点爆破振动信号分解到第7层，则爆破地震波按照从低到高划分为如下128个频率区间：0~4，4~8，8~12，…，508~512 Hz。

3.2   各频率区间信号的能量表征

选取爆破振动信号的全频带[0 Hz, 512 Hz]作为频率区间，对爆破振动信号进行小波包分解，然后由小波包重构公式求取爆破振动信号各频率区间内的重构子信号，则原始爆破振动信号完全可以由这些重构子信号来表示^[16]，即：

式中：S_{7, j}(t) (j=0, 1, 2, …, 2⁷-1)为原始爆破振动信号f(t)第j个节点上的重构子信号。

爆破振动信号各频率区间内的重构子信号能量可以表示为：

式中：x_{j, k}(t) (j=0, 1, 2, …, 2⁷-1;k=1, 2, …, m)为重构子信号S_{7, j}(t)离散采样点的幅值；m为爆破振动信号离散采样点的个数。

由式(6)可知，设原始爆破振动信号的总能量为E₀，则有：

爆破振动信号各频率区间内的重构子信号能量占被分析信号总能量的比例为

3.3   各频率区间信号的能量分布规律

根据式(6)~(8)，用MATLAB语言编写程序，然后分别对图 1(a)所示的爆破振动加速度信号和图 3(a)所示的去除趋势项后的爆破振动速度信号进行小波包分解，得到各频率区间内信号能量分布的计算结果，图 4(a)和4(b)分别为爆破振动加速度信号和去除趋势项后爆破振动速度信号在不同频率范围内的能量分布。通过分别对照图 4(a)与图 1(b)以及图 4(b)与图 3(b)可以看出，与建立在传统Fourier变换基础上的频谱分析相比，基于小波变换的能量分析具有更精细的频率分辨率。从图 4可以看出，爆破振动信号绝大部分能量都集中在40~140 Hz频率范围内，分别占爆破振动信号总能量的89.08%和92.89%，但两者的相关系数仅为0.546 6；实测爆破振动加速度信号与爆破振动速度信号在不同频率范围内能量分布的方差分别为7.214 7和5.955 2，由此表明爆破振动速度信号比实测爆破振动加速度信号在不同频率范围内能量分布更均匀。若将爆破振动信号各频率区间范围划分得稍宽一点，即划分成如下8个频率区间：0~20、20~40、40~60、60~80、80~100、100~120、120~140和140~512 Hz，则爆破振动信号各频率区间内能量占被分析信号总能量的百分比见表 1。由表 1可以计算出，实测爆破振动加速度信号与爆破振动速度信号各频率区间内能量分布的相关系数高达0.968 1。由此表明，各频率区间范围划分越宽，实测爆破振动加速度信号与爆破振动速度信号各频率区间内能量分布的相关程度越高。

图  4  爆破振动信号的能量分布

Figure  4.  Energy distribution of blast vibration signal in different frequency ranges

下载: 全尺寸图片 幻灯片

表  1  爆破振动信号各频率区间内能量分数

Table  1.  Percentage of energy of blast vibration signals in different frequency ranges

加速度信号			速度信号
f/Hz	能量分数/%		f/Hz	能量分数/%
0~20	0.32		0~20	0.65
20~40	1.28		20~40	1.41
40~60	8.95		40~60	7.18
60~80	4.76		60~80	8.66
80~100	25.49		80~100	25.06
100~120	34.92		100~120	41.38
120~140	14.96		120~140	10.61
140~512	9.32		140~512	5.05

下载: 导出CSV 

| 显示表格

3.   结论

(1) 在模式自适应小波原理的基础上，提出了将实测爆破振动加速度信号梯形数值积分后的爆破振动速度信号作为一个模式自适应小波基的方法，采用该方法可成功去除时域积分后爆破振动速度信号中趋势项的干扰，并通过实测爆破振动加速度信号功率谱与爆破振动速度信号功率谱的比较，验证了该方法是切实可行的；

(2) 通过对去除趋势项后的爆破振动速度信号进行能量特征分析，并将此法与爆破振动速度信号傅里叶谱分析作了比较：基于小波变换的能量分析比傅里叶谱分析具有更精细的频率分辨率，更适合于对频率分辨率要求更高的爆破振动信号进行分析；

(3) 将实测爆破振动加速度信号与去除趋势项后爆破振动速度信号的各频率区间内能量特征作了比较，得到频率区间范围划分越宽，实测爆破振动加速度信号与爆破振动速度信号在不同频率范围内能量分布的相关程度越高；反之，相关程度越低。