楊丹迪, 路敬祎, 周怡娜, 王 闖, 董宏麗

(東北石油大學(xué) a. 復(fù)雜系統(tǒng)與先進(jìn)控制研究院; b. 黑龍江省網(wǎng)絡(luò)化與智能控制重點(diǎn)實(shí)驗(yàn)室, 黑龍江 大慶 163318)

0 引 言

管道運(yùn)輸因其成本低, 動(dòng)力消耗小等優(yōu)勢(shì), 被廣泛應(yīng)用于石油天然氣的輸送。目前, 由于管道運(yùn)行年限的不斷增加, 老化程度逐漸嚴(yán)重, 而且人為破壞、 盜取行為也嚴(yán)重威脅著油氣管道的安全, 泄漏風(fēng)險(xiǎn)逐漸增大。 而管道泄漏將給社會(huì)造成很大的經(jīng)濟(jì)損失以及危害, 因此人們對(duì)管道安全問(wèn)題非常關(guān)注, 對(duì)管道泄漏檢測(cè)技術(shù)[1-4]的研究已成為國(guó)內(nèi)外學(xué)術(shù)界的研究熱點(diǎn)。目前管道泄漏檢測(cè)方法有很多, 主要有負(fù)壓波檢測(cè)法、 聲波法、 漏磁檢測(cè)法、 光纖檢測(cè)法和人工智能算法等。聲波法[5-6]因其方法簡(jiǎn)單、 維護(hù)費(fèi)用低得到了廣泛應(yīng)用。經(jīng)驗(yàn)?zāi)B(tài)分解(EMD: Empirical Mode Decomposition)[7-8]是一種將時(shí)域信號(hào)按頻率尺度自適應(yīng)分解的數(shù)值算法, 非常適用于非平穩(wěn)信號(hào)的分析, 但在信號(hào)分解時(shí), 容易存在模態(tài)混疊現(xiàn)象, 從而導(dǎo)致分解的IMF分量失真。

Dragomiretskiy等[9]于2014年提出VMD(Variational Mode Decomposition)算法, 可以很好地抑制EMD方法的缺陷[10-11]。文獻(xiàn)[12]將VMD算法與相關(guān)系數(shù)與峭度結(jié)合作為前置濾波處理, 選取合適濾波器長(zhǎng)度及解卷積周期對(duì)重構(gòu)信號(hào)進(jìn)行最大相關(guān)峭度解卷積運(yùn)算以實(shí)現(xiàn)對(duì)軸承故障進(jìn)行準(zhǔn)確分離, 并通過(guò)與EEMD(Ensemble Empirical Mode Decomposition)降噪后解卷積處理的對(duì)比, 驗(yàn)證了VMD前置濾波處理的優(yōu)越性。文獻(xiàn)[13]提出了一種VMD算法與模糊相關(guān)分類(lèi)結(jié)合的自適應(yīng)去噪方法, 能檢測(cè)管道中出現(xiàn)的小泄漏, 與支持向量機(jī)方法相比具有更高的檢測(cè)率, 但VMD算法存在預(yù)設(shè)尺度K難以選取及難以確定分解后有效的IMF分量的問(wèn)題。文獻(xiàn)[14]提出了互信息與VMD結(jié)合的去躁方法在天然氣管道泄漏檢測(cè)中的應(yīng)用, 通過(guò)計(jì)算VMD分解后的IMF分量與原始管道信號(hào)的互信息值, 篩選有效的IMF分量進(jìn)行重構(gòu), 克服了VMD在分解時(shí)選取IMF分量的盲目性。文獻(xiàn)[15]提出VMD與豪斯多夫距離結(jié)合的方法對(duì)信號(hào)進(jìn)行去噪, 先用VMD對(duì)信號(hào)進(jìn)行分解后得到IMF分量, 計(jì)算各分量和信號(hào)的概率密度函數(shù), 最后通過(guò)各分量的概率密度函數(shù)和原始信號(hào)的概率密度函數(shù)之間的豪斯多夫距離判斷出有效IMF分量對(duì)信號(hào)進(jìn)行重構(gòu)。

基于以上研究, 為更好地選取VMD算法分解后有效IMF分量及解決信號(hào)中的噪聲干擾問(wèn)題, 減少高頻段有用信息的丟失, 筆者提出了一種基于相關(guān)系數(shù)、 峭度和小波變換的VMD改進(jìn)算法, 通過(guò)仿真及試驗(yàn)將筆者算法與其他算法進(jìn)行對(duì)比, 以信噪比和均方誤差作為性能指標(biāo), 驗(yàn)證了筆者提出的算法在信號(hào)去噪上的優(yōu)越性和應(yīng)用在管道泄漏檢測(cè)中的有效性。

1 VMD算法及相關(guān)理論

1.1 VMD算法原理

VMD算法是一種可變尺度的信號(hào)處理方法, 可將一個(gè)實(shí)際信號(hào)f分解成K個(gè)模態(tài)函數(shù), 使每個(gè)模態(tài)的估計(jì)帶寬之和最小, 具體分解步驟如下。

1) 采用Hilbert變換, 計(jì)算每個(gè)uk的解析函數(shù)以獲得其相應(yīng)的單側(cè)頻譜

(1)

2) 對(duì)每個(gè)模態(tài)uk, 通過(guò)與其對(duì)應(yīng)的中心頻率的指數(shù)項(xiàng)混疊, 將每個(gè)模態(tài)uk的頻譜轉(zhuǎn)移到相應(yīng)基帶

(2)

3) 由解調(diào)信號(hào)的高斯平滑度和梯度平方范數(shù)估計(jì)帶寬。

4) 由上述步驟得到的約束變分問(wèn)題為

(3)

其中uk={u1,u2,…,uk}為各模態(tài)函數(shù),ωk={ω1,ω2,…,ωk}為各中心頻率。

通過(guò)引入拉格朗日乘子λ(t)和二次懲罰因子?將約束變分問(wèn)題轉(zhuǎn)換為非約束變分問(wèn)題, 二次懲罰因子是典型的實(shí)現(xiàn)重構(gòu)信號(hào)保真度的方法, 拉格朗日乘子則用于實(shí)現(xiàn)精確重構(gòu)。將兩者結(jié)合得到拓展的拉格朗日表達(dá)式如下

采用乘法算子交替方向法解決上述變分問(wèn)題, 迭代優(yōu)化uk+1、ωk+1、λk+1即可求得擴(kuò)展拉格朗日表達(dá)式的“鞍點(diǎn)”, 變分問(wèn)題的解為

(5)

中心頻率的更新方法為

(6)

VMD算法迭代步驟如下。

1)初始化u1、ω1、λ1、n=0。

2) 令n=n+1, 執(zhí)行整個(gè)循環(huán), 更新uk。

3) 對(duì)所有ω≥0, 使

(7)

4) 重復(fù)步驟2), 步驟3)直到滿足約束條件

(8)

從整體上看, 在迭代求解變分模型的過(guò)程中各個(gè)IMF分量的中心頻率和帶寬不斷更新, 直到滿足迭代終止條件。最終, 根據(jù)信號(hào)的頻域特征得到K個(gè)IMF分量, 完成了信號(hào)的自適應(yīng)分割, 能有效避免模態(tài)混疊現(xiàn)象。

然而, VMD算法在迭代過(guò)程中需要預(yù)先確定IMF分量的個(gè)數(shù), 但K值通常難以確定, 導(dǎo)致不能得到有效的IMF分量。

1.2 相關(guān)系數(shù)

相關(guān)系數(shù)(CC: Correlation Coefficient)是由著名統(tǒng)計(jì)學(xué)家卡爾-皮爾遜提出的, 用以反映變量之間相關(guān)關(guān)系密切程度的統(tǒng)計(jì)指標(biāo)[16], 筆者引用相關(guān)系數(shù)判斷IMF分量和原始信號(hào)之間的相關(guān)程度, 以兩變量與各自平均值的離差為基礎(chǔ), 通過(guò)兩個(gè)離差相乘反映兩變量之間的相關(guān)程度。

相關(guān)系數(shù)的計(jì)算公式為

(9)

其中x,y為兩個(gè)矢量,Ex,Ey為對(duì)應(yīng)矢量的期望。r反應(yīng)了變量之間相關(guān)程度的高低。

1.3 峭 度

峭度(Kurtosis)是描述信號(hào)波峰尖度的無(wú)量綱參數(shù)[18], 定義為信號(hào)的歸一化4階中心矩, 反映隨機(jī)變量分布特性的數(shù)值統(tǒng)計(jì)量, 計(jì)算表達(dá)式為

K=E(x-μ)4/σ4

(10)

其中σ和μ分別為信號(hào)x的標(biāo)準(zhǔn)差和均值；E為信號(hào)x的期望值。

峭度指標(biāo)對(duì)沖擊信號(hào)較為敏感, 能較好地反映管道泄漏信號(hào)中的沖擊成分, 信號(hào)中的沖擊成分所占比重越大則峭度指標(biāo)的絕對(duì)值也越大。而VMD算法分解出的高頻模態(tài)分量包含噪聲較多, 其沖擊成分占比大、 峭度值高, 從而通過(guò)對(duì)信號(hào)的峭度分析可以實(shí)現(xiàn)VMD算法分解后有效模態(tài)分量的選取。

通過(guò)觀察多組管道工況信號(hào)可知, 當(dāng)VMD算法分解后的模態(tài)分量中無(wú)噪聲或噪聲成分比較少時(shí), IMF分量的峭度值一般小于4, 且IMF分量從低頻到高頻分布, 其峭度值的整體趨勢(shì)明顯增大, 說(shuō)明分量中的噪聲含量所占比重越大峭度值越大。

1.4 小波變換

小波變換(Wavelet Transform)[19]是在短時(shí)傅里葉變換的基礎(chǔ)上發(fā)展起來(lái)的一種信號(hào)時(shí)頻分析方法, 適合分析非平穩(wěn)信號(hào)。小波去噪的本質(zhì)是根據(jù)不同頻帶上的信號(hào)和噪聲的小波分解系數(shù)強(qiáng)度不同的特點(diǎn), 保留原信號(hào)的小波分解系數(shù), 去除噪聲在各頻帶上的小波系數(shù), 最后再將處理后的系數(shù)進(jìn)行小波重構(gòu), 得到去噪后的重構(gòu)信號(hào)。

目前, 小波變換主要分為以下幾種方法：

1) 模極大值小波去噪方法；

2) 小波閾值去噪法[19-20]；

3) 平移不變量法；

4) 基于各尺度下小波系數(shù)相關(guān)性去噪(屏蔽去噪法)。

筆者選用小波自適應(yīng)閾值去噪, 選擇合適的小波基函數(shù)和分解層數(shù)是小波閾值去噪的關(guān)鍵, 選取Daubechies系列中的db5小波基進(jìn)行5層分解, 然后對(duì)分解后的系數(shù)進(jìn)行閾值處理, 最后將閾值處理后的小波系數(shù)進(jìn)行重構(gòu)得到去噪信號(hào)。小波閾值去噪的優(yōu)點(diǎn)是計(jì)算簡(jiǎn)單且快速, 算法易實(shí)現(xiàn), 去噪效果較好, 在實(shí)際中有較為廣泛的應(yīng)用。

在實(shí)際信號(hào)處理中, 有用信息一般包含在低頻信號(hào)中, 而噪聲信號(hào)多為高頻信號(hào), 利用相關(guān)系數(shù)和峭度值篩選出VMD算法分解出的高頻模態(tài), 將其進(jìn)行小波變換就可起到信號(hào)去噪的作用。

1.5 改進(jìn)的VMD算法

通過(guò)對(duì)以上幾種方法的研究, 筆者將VMD算法與相關(guān)系數(shù)和峭度結(jié)合選取有效分量, 經(jīng)過(guò)比較其比結(jié)合單一方法更加準(zhǔn)確。雖然大多數(shù)聲信號(hào)的主要特征集中在低頻段, 但將相關(guān)系數(shù)較小或峭度較高的高頻分量直接濾除會(huì)損失掉高頻部分的有效信息。故筆者提出基于相關(guān)系數(shù)、 峭度和小波變換的VMD改進(jìn)算法,算法流程如圖1所示。

圖1 VMD改進(jìn)算法流程圖Fig.1 Flow chart of VMD improved algorithm

2 仿真試驗(yàn)與結(jié)果分析

為驗(yàn)證筆者提出算法的有效性, 選取各分量為2 Hz、576 Hz、6 912 Hz 3個(gè)頻段的余弦信號(hào)作為原始信號(hào)并加入噪聲強(qiáng)度為0.1的高斯白噪聲進(jìn)行仿真分析, 即

fn(t)=cos(4πt)+0.5cos(1 152πt)+0.05cos(13 824πt)+μ

其中μ為加性高斯白噪聲。

圖2為不含白噪聲的原始輸入信號(hào), 圖3為加入高斯白噪聲后的復(fù)合信號(hào)。

圖2 不含白噪聲的原始輸入信號(hào) 圖3 加入噪聲的復(fù)合信號(hào) Fig.2 Original input signal without white noise Fig.3 Composite signal with noise

根據(jù)VMD分解經(jīng)驗(yàn), 選取K=4, VMD分解后的IMF分量如圖4a～圖4d所示。

圖4 VMD分解得到的IMF分量Fig.4 IMF components derived from VMD decomposition

由式(9)計(jì)算VMD分解后得到的各個(gè)IMF分量的相關(guān)系數(shù)值, 如表1所示。

表1 各IMF分量的相關(guān)系數(shù)

表2 剩余IMF分量峭度值

由表2判斷出IMF2峭度值小于設(shè)定的峭度閾值4, 故保留IMF2分量作為峭度值篩選出的有效分量； 接著將峭度值大于4的IMF3分量進(jìn)行小波閾值去噪處理, 選取db5小波基進(jìn)行5層分解, 對(duì)分解后的系數(shù)進(jìn)行閾值處理, 將處理后的小波系數(shù)進(jìn)行重構(gòu)濾除噪聲, 最后將去噪后的分量與根據(jù)相關(guān)系數(shù)和峭度選出的有效分量進(jìn)行重構(gòu)。計(jì)算筆者改進(jìn)算法和VMD與相關(guān)系數(shù)、 VMD與小波、 峭度、 豪斯多夫距離、 互信息、 小波去噪的重構(gòu)信號(hào)的信噪比和均方誤差進(jìn)行對(duì)比, 結(jié)果如圖5和圖6所示。

由圖5可以看出, VMD與豪斯多夫距離、 互信息、 小波去噪算法對(duì)信號(hào)去噪后的信噪比分別為6.870 2 dB,6.811 5 dB,6.700 1 dB, 而筆者提出的改進(jìn)VMD算法可以比較明顯提高信號(hào)的信噪比, 達(dá)到了20.084 9 dB； VMD與相關(guān)系數(shù)、 峭度和小波變換單獨(dú)結(jié)合的算法對(duì)信號(hào)進(jìn)行去噪, 信號(hào)信噪比分別為18.452 0 dB,18.953 4 dB和19.783 9 dB, 均低于筆者改進(jìn)算法。

由圖6可以看出, VMD與豪斯多夫距離、 互信息、 小波去噪算法對(duì)信號(hào)去噪后計(jì)算的均方誤差結(jié)果分別為0.359 6,0.362 1,0.366 8, 筆者改進(jìn)算法的均方誤差最小, 為0.078 5； 將VMD與相關(guān)系數(shù)、 峭度和小波變換單獨(dú)結(jié)合的算法對(duì)信號(hào)進(jìn)行去噪計(jì)算出的均方誤差分別為0.094 8,0.089 5,0.080 4, 均大于筆者算法計(jì)算的均方誤差。根據(jù)圖5,圖6的對(duì)比結(jié)果, 可以看出筆者改進(jìn)的VMD算法結(jié)合了3者的優(yōu)點(diǎn), 得到了較為理想的去噪效果。

圖5 SSNR對(duì)比圖 圖6 MMSE對(duì)比圖 Fig.5 SSNR contrast diagram Fig.6 MMSE contrast diagram

3 管道泄漏信號(hào)檢測(cè)與分析

筆者研究所采用的管道數(shù)據(jù)是從天然氣管道泄漏檢測(cè)實(shí)驗(yàn)室采集的, 其管道總長(zhǎng)160 m, 管徑為DN50, 管壁厚為20 mm, 管道內(nèi)可以實(shí)現(xiàn)氣體和液體的運(yùn)輸。管道上設(shè)有多個(gè)泄漏點(diǎn)用于仿真現(xiàn)場(chǎng)管道的泄漏, 并可通過(guò)監(jiān)控臺(tái)對(duì)管道的相關(guān)參數(shù)進(jìn)行監(jiān)控, 實(shí)驗(yàn)室管道泄漏檢測(cè)系統(tǒng)如圖7所示。試驗(yàn)利用Labview編程環(huán)境, 借助NI公司的數(shù)據(jù)采集板卡, 分別利用采樣頻率1 kHz,5 kHz,10 kHz進(jìn)行數(shù)據(jù)采集, 將數(shù)據(jù)采集卡分別與兩個(gè)壓電式聲波傳感器相連, 并利用計(jì)算機(jī)完成數(shù)據(jù)采集。圖8為實(shí)驗(yàn)室采集到的泄漏數(shù)據(jù), 本次實(shí)驗(yàn)截取的長(zhǎng)度為2 048個(gè)采樣點(diǎn)。

a 天然氣管道 b 控制臺(tái)圖7 實(shí)驗(yàn)室管道泄漏檢測(cè)系統(tǒng)Fig.7 Pipeline leakage detection system in Laboratory

圖8 管道原始信號(hào)Fig.8 Pipeline original signal

表3 各個(gè)IMF分量的相關(guān)系數(shù)

表4 剩余IMF分量的峭度值

由表4可以看出, IMF2、 IMF3和IMF5分量的峭度小于設(shè)定的峭度閾值4, 將它們保留作為峭度值篩選出的有效分量； 將IMF4、 IMF6和IMF7分量進(jìn)行小波閾值去噪處理, 將去噪后的分量與根據(jù)相關(guān)系數(shù)和峭度選出的有效分量進(jìn)行重構(gòu), 并計(jì)算重構(gòu)信號(hào)的信噪比和均方誤差作為性能評(píng)價(jià)指標(biāo)。

計(jì)算筆者改進(jìn)算法和VMD與相關(guān)系數(shù)、 小波、 峭度、 豪斯多夫距離、 互信息以及小波去噪的管道重構(gòu)信號(hào)的信噪比和均方誤差進(jìn)行對(duì)比, 結(jié)果如圖9和圖10所示。

圖9 SSNR對(duì)比圖 圖10 MMSE對(duì)比圖 Fig.9 SSNR contrast diagram Fig.10 MMSE contrast diagram

由圖9可以看出, VMD與小波結(jié)合算法、 峭度結(jié)合算法、 小波去噪的重構(gòu)信號(hào)信噪比分別為17.999 6 dB,17.154 0 dB,11.024 6 dB,筆者改進(jìn)的VMD算法重構(gòu)信號(hào)信噪比較高, 為18.290 2 dB, 而VMD與相關(guān)系數(shù)結(jié)合算法、 豪斯多夫距離結(jié)合算法、 互信息結(jié)合算法選取的有效模態(tài)均為IMF1, 故重構(gòu)信號(hào)的信噪比均為13.511 0 dB。由圖10可以看出, VMD與小波結(jié)合算法、 峭度結(jié)合算法、 小波去噪的重構(gòu)信號(hào)的均方誤差分別為0.036 2,0.039 9,0.080 7, 均高于筆者算法的0.035 0, VMD與相關(guān)系數(shù)結(jié)合算法、 豪斯多夫距離結(jié)合算法和互信息結(jié)合算法選取模態(tài)相同, 均方誤差均為0.060 6, 均高于筆者算法。

根據(jù)圖9、 圖10的SNR和MSE的對(duì)比結(jié)果實(shí)驗(yàn), 驗(yàn)證了筆者提出的VMD改進(jìn)算法能提高信號(hào)的信噪比, 降低均方誤差, 能較好地處理復(fù)雜的管道信號(hào), 驗(yàn)證了在管道泄漏檢測(cè)中的有效性。

4 結(jié) 語(yǔ)

鑒于VMD算法分解后的有效模態(tài)分量選擇比較困難的問(wèn)題, 筆者提出了一種改進(jìn)的VMD算法, 將VMD算法與相關(guān)系數(shù)和峭度結(jié)合選出有效模態(tài)分量, 通過(guò)小波變換去除剩余分量的噪聲, 將有效模態(tài)分量和去除噪聲后的剩余模態(tài)分量進(jìn)行重構(gòu), 可以減少高頻段有用信息的丟失, 并能較大程度地去除信號(hào)在管道傳播中混入的噪聲。通過(guò)仿真信號(hào)和實(shí)驗(yàn)室管道采集信號(hào)進(jìn)行試驗(yàn), 與其他算法相比, 筆者提出的算法在信號(hào)去噪方面有著較好的優(yōu)越性, 并且在管道泄漏檢測(cè)方面有著較好的前景。