基于CEEMDAN+二次时延估计的管道泄漏定位技术

时间：2024-07-28

杨岳鹏徐晨姚世聪岳景寅李志强王海林

渤海钻探第四钻井工程分公司

油气管道泄漏往往造成巨大的经济损失和环境污染，严重影响企业的社会形象，因此对油气管道进行泄漏定位检测显得尤为重要[1]。目前，管道泄漏检测的方法主要有流量平衡法、瞬态模型法、神经网络法、负压波法和音波法[2]，其中负压波法和音波法的应用最为广泛。

1 油气管道泄漏检测方法及筛选

负压波法主要依靠压力传感器对泄漏信号进行捕捉，频率范围在60～100 Hz，泄漏中出现的大幅压降在液体管道上具有较高的灵敏性，但对于可压缩的气体管道适应性较差；音波法按照音波频率不同分为次声波（频率范围0.000 1～20 Hz）和超声波（频率范围16 000～1012Hz）法等，BRODETSKY 等[3]通过研究音波的衰减情况，证明了频率越大，音波信号的衰减越大，故次声波较其余音波在介质中传播距离远且能量损失低。次声波法检测的关键在于信号奇异点的获取和次声波到达上、下游传感器的时延估计，但管道在泄漏的过程中，系统本体噪声、管道运行噪声和外部随机噪声等均会对次声波信号造成影响，使信号湮灭在噪声中，尤其对缓慢泄漏和微小泄漏的检测成功率不高。近年来，国内外学者针对信号降噪和时延估计进行了大量研究，主要采用小波分析[4]、经验模态分解（EMD）[5-7]和局部均值分解（LMD）[8]等方法。但上述方法在信号处理中仍存在一定缺陷，如小波分析会造成重构信号振动失真，丢失有效信息；EMD 存在模态混叠和残余噪声问题；LMD 存在端点效应、平滑次数较多的现象，且平滑时步长无法最优化处理，计算速度较慢。

针对以上不足，以次声波泄漏检测原理为基础，采用改进的自适应白噪声完全集合经验模态分解（CEEMDAN）对原始信号进行降噪分解，利用相关系数选取与原始信号相关度较大的固有模态分量（IMF），去除噪声和虚假信息，运用二次时延估计实现时差计算[9]，对泄漏点进行定位。通过现场试验对算法进行验证，以期为管道泄漏定位提供理论依据和实际参考。

2 泄漏定位原理

2.1 次声波泄漏定位原理

当输气管道发生泄漏时，气体从泄漏处喷出，管内气体的流动参数发生改变，泄漏处气体密度减小，压力降低，泄漏点相邻两侧的气体在压差的作用下向泄漏点处补充，进而使泄漏点相邻两侧的气体密度进一步减小，压力进一步降低，依次类推向更远处传播。输气管道泄漏产生的次声波本质上属于湍流脉动引起的四极子声源和偶极子声源，声波衰减遵循指数规律，安装在管道上、下游两端的声波传感器监听并捕捉传来的声波信号，根据式（1）对泄漏点进行定位。

式中：L1为泄漏点距离上游传感器的距离，m；a为次声波声速，m/s；V为流体速度，m/s；L2为上、下游传感器的距离；Δt为上、下游接收到次声波的时间差。

一般a远大于V，故式（1）可简化为公式（2）：

可见，泄漏定位的关键在于Δt的计算精度，而信号中含有的大量噪声会对时延估计造成较大影响，故需对信号进行降噪处理[10-11]。

2.2 CEEMDAN信号降噪原理

CEEMDAN 是在EMD、EEMD（集合平均经验模态分解）算法基础上的改进，其重构误差几乎为0，解决了不同信号添加噪声后出现不同模态的问题。算法步骤如下：

（1）设采用EMD算法分解后得到的第i个固有模态分量为imfi()，在原始信号x(t)中添加了标准差为kε、均值为0 的高斯白噪声n(t)，得到混合信号x(t)+kεn(t)，其中ε为原始信号标准差，k为白噪声幅度。对混合信号进行分解，按时间尺度由大到小得到M次分解时的不同模态分量，将M个模态分量求总体平均后，得到CEEMDAN算法的第一个固有模态分量IMF1(t)，如公式（3）所示：

其中，残差r1(t)=x(t)-IMF1(t)。

（2）将残差与白噪声经EMD 算法分解后的模态分量继续组成混合信号r1(t)+imf1[kεnj(t)]，按照公式（3）分解后，得到CEEMDAN 算法的第一个固有模态分量IMF2(t)，如公式（4）所示：

其中，残差r2(t)=r1(t)-IMF2(t)。

（3）重复以上步骤，当残差信号为单调函数时，停止计算，得到n+1 个固有模态分量和残差，如公式（5）、（6）所示：

CEEMDAN算法加入了由EMD算法分解后的噪声信号，而不是将噪声直接添加到原始信号中，解决了白噪声从高频到低频转移传递的问题。但仍存在白噪声无法完全中和的现象，故采用以正负成对的形式将白噪声添加到信号中，则每次分解得到的模态分量中，一部分残差分量含有正噪声，另一部分残差分量含有负噪声，两者平均后正负消除，使算法更准确。

原始信号分解后得到噪声分量、冗余分量和有效分量，通过相关函数分析频域内两个信号的相关程度，得到原始信号中的有效分量，相关系数公式如公式（7）所示：

式中：R(i)为第i个固有模态分量与原始信号x(t)的相关系数；IMFi(t)为固有模态分量的平均值；为原始信号的平均值。

目前，关于相关系数的阈值尚没有统一认定，但通常认为R(i)的绝对值小于0.2 为微相关，0.2～0.5 为中相关，0.5～0.8 为高相关，0.8～1.0 为显著相关。

2.3 二次时延估计

时延估计算法常采用广义互相关算法，将经相关分析后的IMF分量叠加后进行重构，得到降噪后的信号，将上、下游信号进行互相关函数计算，通过相关峰值的位置得到从泄漏点传播到两端的时间差。设上、下游降噪后的信号分别为x1(t)、x2(t)，计算其互相关函数，如公式（8）所示：

式中：τ即为两个信号的时间差。

求τ即为求R12(τ)的最大值，如公式（9）所示。定义信号x1(t)、x2(t)的互功率谱为G12，则R12(τ)可用G12表示，如公式（10）所示：

式中：φ12(w)为广义互加权函数，当取1时即为基本窗函数，此时R12(τ)表示基本互相关傅里叶变换；X1(w)为信号x1(t)的傅里叶变换；X2*(w)为信号x2(t)的傅里叶变换共轭。

广义互相关分析中假设声源信号与噪声信号互不相关，但在实际中，噪声信号并不完全独立，干扰项依然存在，故在此采用二次时延估计相关算法。将信号x1(t)进行自相关后得到R11(τ)，将R11(τ)与R12(τ)互相关后得到R(τ)，信号经过一次自相关后可提高信噪比，信号噪声的多次相关也抑制了噪声对时差结果的影响。

3 实例分析

以两个集气站间管道为例进行实验，验证算法准确性。管输介质为净化天然气，管径Φ273 mm×6 mm，全长10.154 km，运行压力5.5～6.5 MPa，在首末站设置声波传感器，利用中间5.231 km处的阀室放空阀作为泄漏点，采用实际放空的方式模拟泄漏进行测试，阀门开度10%。一个完整的泄漏信号波形图中共采集104个数据点，信号时长为20 s，次声波信号波形图见图1，纵坐标是归一化后的信号幅值。当管道处于正常工况时，声波传感器捕捉到的信号为背景噪声，当管道一旦发生泄漏，产生的泄漏信号和正常工况下的背景噪声会一同被传感器捕捉。图1中未泄漏时信号相对平稳，泄漏后信号幅值有较大波动，先下降后上升，符合放空阀门启闭的特征，同时下游信号的幅值明显高于上游信号，这是由于泄漏点距离下游传感器更近，信号衰减小，与实际情况相符。但信号波形中均含有不同程度的噪声，需对信号降噪处理。

图1 次声波信号波形图Fig.1 Waveform diagram of infrasonic signal

3.1 泄漏信号分解

按照2.2 段的流程对泄漏信号进行自适应，以下游次声波信号为例，得到10 组IMF 分量和残差分量（图2）。其中，前3阶IMF分量呈现宽、高频特性，是噪声分量的主要集中区域；4～7 阶IMF 分量的幅值与原始信号保持一致，属于有效泄漏信号分量；8～10阶IMF分量的信号特征逐渐消失，呈单调函数，属于冗余分量。为验证分解结果的有效性，计算各阶IMF分量与原始信号的相关系数，用于提取有效信号（图3）。前3阶的相关系数均在0.2以下，属于微相关；4-7阶的相关系数较大，均在0.5 以上，属于高相关，其中IMF6 的相关系数最大，为0.712；第8 阶的相关系数突然降低至0.2 左右，9～10阶的IMF分量与原始信号的相关系数几乎为0，说明9～10阶的IMF分量与原始信号没有相关性。综上所述，选取4～7 阶的IMF 分量和残差分量，将其叠加进行重构形成降噪后的信号（图4）。

图2 泄漏信号分解的各阶IMF分量和残差分量（下游信号）Fig.2 IMF components and residual components of leakage signal decomposition(downstream signals)

图3 各阶IMF分量与原始信号的相关系数（下游信号）Fig.3 Correlation coefficients between IMF components of various orders and original signals(downstream signals)

图4 原始信号与降噪信号波形图Fig.4 Waveform diagram of original signal and denoising signal

信号降噪后，波形曲线更加光滑，保留了原始信号的波形趋势，泄漏时刻的波形更加明显，且其余未泄漏时刻的幅值波动明显减小。为评价CEEMDAN 算法的降噪效果，对比了EMD、EEMD算法分解重构后的信号，采用信噪比(SNR)和均方根误差(RMSE)衡量去噪效果，SNR反映有效分量与噪声分量之间的比值，RMSE反映有效分量与原始信号之间的差异程度。因此在滤波后，应预期提高信噪比，减小均方根误差，即SNR越大、RMSE越小，降噪效果越好（表1）。其中，CEEMDAN算法的SNR最大且RMSE最小，说明该算法在改善模态混叠上效果显著，自适应地区分了噪声和有效信号，保证了分解的完整性，信号重构后的误差更小。

表1 原始信号降噪效果对比Tab.1 Comparison of noise reduction effects of original signals

3.2 泄漏定位

分别采用广义互相关算法和二次时延互相关算法进行定位计算，互相关时域分析见图5。广义互相关算法得到的延时采样点数为281，则泄漏时差为0.562 s；二次时延互相关算法得到的延时采样点数为326，则泄漏时差为0.652 s。两种算法的泄漏时差均为正值，这是由于泄漏点距离上游更远，因此传播至上游传感器处花费的时间更长，时差量为正。次声波波速a与管输介质的密度、温度、压力和质点速度等参数相关，采用公式（11）计算如下：

图5 互相关时域分析Fig.5 Cross-correlation time domain analysis

式中：k为介质压缩系数，无量纲，根据管输压力和天然气压缩因子计算；ρ为气体密度，kg/m3；D为管道直径，m；E为管材弹性模量，Pa；e为管道壁厚，m；C为与管道约束条件相关的修正系数，按照第二种约束条件计算，C=1-μ2，μ为泊松系数取0.3。

根据不同工况多次计算取平均值，得到a=405.2 m/s。最终，按照公式（4）计算广义互相关算法的定位结果为5.190 km，定位误差41 m，相对误差0.78%；二次时延互相关的定位结果为5.209 km，定位误差22 m，相对误差0.42%。

3.3 对比验证

为验证CEEMDAN+二次时延估计算法的可靠性，将其与文献[8]中的LMD+小波分析、文献[4]中的小波变换+领域差值两种算法进行对比，并选取阀门开度10%、5%和1%进行测试（表2）。随着泄漏点孔径的减小，互相关分析的峰值点与真实峰值点之间的差距逐渐增大，定位误差也逐渐增大，其中当阀门开度为1%时，只有CEEMDAN+二次时延估计算法对泄漏进行了定位，而其余两种算法由于信号下降的拐点不易识别，故未识别出报警事件，说明CEEMDAN+二次时延估计算法未出现误报或漏报，对微小泄漏的检测具有良好的适应性。当阀门开度为10%时，CEEMDAN+二次时延估计算法的定位相对误差较其余两种算法分别降低了71.03%、60%；当阀门开度为5%时，定位相对误差较其余两种算法分别降低了84.21%、73.38%，对10 km 左右的管线总体定位误差可控制在50 m以内。

表2 不同定位算法结果对比Tab.2 Comparison of results of different positioning methods