劉名楊，楊 進(jìn)，鄭 偉，樊恩東

(重慶大學(xué) 光電工程學(xué)院 光電技術(shù)及系統(tǒng)教育部重點(diǎn)實(shí)驗(yàn)室，重慶 400044)

供水管道泄漏不僅帶來(lái)經(jīng)濟(jì)損失、環(huán)境污染的問(wèn)題，甚至對(duì)城市居民的健康安全構(gòu)成威脅，更重要的是造成了水資源浪費(fèi)問(wèn)題。造成管道泄漏的原因主要包括管道老化、化學(xué)腐蝕、外力損傷以及疏于管理維護(hù)。據(jù)2020年《2020—2026年中國(guó)城市供水行業(yè)深度調(diào)研及投資前景分析報(bào)告》指出[1]，中國(guó)大陸地區(qū)平均每年城市管網(wǎng)泄漏率達(dá)15.3%，部分城市甚至超過(guò)25%，這高于其他許多發(fā)達(dá)國(guó)家標(biāo)準(zhǔn)。

在管道泄漏檢測(cè)領(lǐng)域中，聲發(fā)射檢測(cè)技術(shù)以其成本低、非開挖、高效率的特點(diǎn)被得到廣泛應(yīng)用。經(jīng)典的管道聲發(fā)射檢測(cè)如圖1所示，聲學(xué)或振動(dòng)傳感器緊貼于泄漏管道兩端獲取泄漏信號(hào)[2]。但在實(shí)際的泄漏檢測(cè)中，泄漏聲發(fā)射信號(hào)可能會(huì)受到各種各樣的管道外部和管道內(nèi)部噪聲干擾，外部噪聲有汽車噪聲、機(jī)器噪聲、人為噪聲等等，內(nèi)部噪聲以管道自身噪聲為主，如支管流致噪聲等其他非直管道接頭引起的噪聲。

圖1 經(jīng)典管道聲發(fā)射檢測(cè)示意圖

為了有效去除噪聲、提取源泄漏信號(hào)，目前的研究大多采用小波分析[3-4]、EMD(empirical mode decomposition)類[5-6]、VMD(variational mode decomposition)類[7-8]等算法來(lái)去除或抑制噪聲。Ting等[9]提出一種改進(jìn)的雙樹復(fù)小波變換技術(shù)(dual-tree complex wavelet transform, DTCWT)將含噪泄漏信號(hào)分解成若干個(gè)子信號(hào)，利用軟閾值去除子信號(hào)中的噪聲，再重建泄漏信號(hào)；郭晨城等[10]提出改進(jìn)的 EMD 信號(hào)增強(qiáng)方法對(duì)含噪泄漏信號(hào)進(jìn)行信號(hào)增強(qiáng)處理，該方法雖然不需要掌握泄漏信號(hào)和噪聲信號(hào)的先驗(yàn)信息，但要對(duì)其特征進(jìn)行數(shù)學(xué)假設(shè)解析；Diao等[11]將粒子群算法和最大熵法結(jié)合起來(lái)，優(yōu)化了VMD的懲罰因子和模態(tài)數(shù)選擇，隨后對(duì)含噪泄漏信號(hào)進(jìn)行分解重構(gòu)，去除背景噪聲以提取源泄漏信號(hào)。

這些方法一般都是從直的泄漏管采集泄漏信號(hào)，且認(rèn)為混合泄漏信號(hào)是等效于噪聲同源泄漏信號(hào)的疊加，它們對(duì)噪聲的處理過(guò)程是一個(gè)去除或抑制加性噪聲的過(guò)程。但是，管內(nèi)噪聲和泄漏聲信號(hào)耦合疊加以后，以管道作為傳播信道進(jìn)行傳播，其傳播過(guò)程是一個(gè)卷積過(guò)程，因此，對(duì)含管內(nèi)噪聲的泄漏信號(hào)進(jìn)行分解處理是一個(gè)解卷積的過(guò)程。此外，包含支管流致噪聲在內(nèi)，所有的管道內(nèi)部噪聲都屬于一種后驗(yàn)噪聲，其時(shí)域、頻域和概率密度特征都是未知的。但前述幾種方法都是基于外部噪聲的功率譜特征、時(shí)頻特征等先驗(yàn)信息可知或假設(shè)可知，以實(shí)現(xiàn)對(duì)聲發(fā)射信號(hào)進(jìn)行處理。因此，上述方法技術(shù)不能直接被運(yùn)用于對(duì)管內(nèi)噪聲的抑制或去除。在支管泄漏檢測(cè)中，泄漏信號(hào)和支管流致噪聲具有少量的先驗(yàn)知識(shí)，且支管流致噪聲和泄漏信號(hào)耦合以后與傳播信道(管道)之間是卷積關(guān)系，因此，從含支管流致噪聲的泄漏信號(hào)中分離出源泄漏信號(hào)來(lái)定位管道泄漏是一個(gè)典型的盲解卷積分離(blind deconvolution separation, BDS)問(wèn)題[12]。

在盲分離中，快速獨(dú)立主成分分析(fast independent component analysis，F(xiàn)astICA)因收斂速度快，且分離性能好，被得到廣泛應(yīng)用。與此同時(shí)，一些改進(jìn)的FastICA也被得到廣泛應(yīng)用。Jiménez-González等[13]將時(shí)空ICA(space and time independent component analysis, ST-ICA)與FastICA相結(jié)合提出ST-FastICA，實(shí)現(xiàn)對(duì)胎兒心電圖信號(hào)的高質(zhì)量分離；Dong等[14]將四階矩累積量的聯(lián)合對(duì)角化運(yùn)用于FastICA算法中，大大降低FastICA的運(yùn)算量；Xie等[15]提出一種改進(jìn)的復(fù)數(shù)域FastICA算法，該方法在分離矩陣的更新迭代過(guò)程中引入噪聲因子，實(shí)現(xiàn)對(duì)低信噪比混合信號(hào)的分離。

但上述主成分分析(independent component analysis，ICA)類盲分離算法大多未考慮待分解信號(hào)的概率分布特征，固定的選擇非線性函數(shù)g來(lái)構(gòu)建算法的代價(jià)函數(shù)和迭代學(xué)習(xí)規(guī)則，這樣做有兩大缺點(diǎn)：其一是算法根據(jù)觀測(cè)樣本的概率密度特征去估計(jì)得分函數(shù)(同樣由g構(gòu)建)有一定的難度；其二是固定選擇的非線性函數(shù)g在觀測(cè)樣本空間內(nèi)不一定是連續(xù)可導(dǎo)的。這使得算法對(duì)混合信號(hào)的分離效率和質(zhì)量受到一定影響，進(jìn)一步使得管道泄漏定位精確度受限。

為了從含有支管流致噪聲的泄漏信號(hào)中高效地分離提取源泄漏信號(hào)，文中提出一種新穎的復(fù)數(shù)域盲分離處理方法，即復(fù)值域高效快速獨(dú)立主成分分析(complex efficient fast independent component analysis，C-EFastICA)方法。該方法基于現(xiàn)有的高效快速獨(dú)立主成分分析(efficient fast independent component analysis，EFastICA)算法[16-17]，EFastICA作為FastICA的一種改進(jìn)算法，它與FastICA相比：① 其計(jì)算復(fù)雜度與FastICA保持相當(dāng)，但收斂速度略快于FastICA；② EFastICA分解精度更高，EFastICA根據(jù)信號(hào)的廣義高斯特征自適應(yīng)地選擇非線性函數(shù)g來(lái)建立代價(jià)函數(shù)和迭代學(xué)習(xí)規(guī)則，對(duì)信號(hào)的分解程度能到達(dá)克拉美羅下界。EFastICA經(jīng)復(fù)值化擴(kuò)展后得到C-FastICA(complex FastICA)，C-FastICA同樣具有分解程度高、分解速度快的優(yōu)點(diǎn)。文中通過(guò)對(duì)泄漏信號(hào)和噪聲信號(hào)進(jìn)行廣義高斯特征計(jì)算，發(fā)現(xiàn)泄漏信號(hào)和噪聲信號(hào)均是超高斯信號(hào)，C-EFastICA根據(jù)其超高斯特性選擇非線性函數(shù)構(gòu)建代價(jià)函數(shù)和迭代學(xué)習(xí)規(guī)則?；旌闲孤┞曅盘?hào)經(jīng)C-EFastICA分解處理以后，再使用互相關(guān)技術(shù)得到時(shí)延估計(jì)值，結(jié)合已知的泄漏聲傳播速度完成泄漏定位工作。

1 支管流致噪聲與含噪泄漏信號(hào)頻譜特征研究

本章采用管徑100 mm的供水管道采集支管流致噪聲和含支管流致噪聲的泄漏聲信號(hào)。100 mm管徑的管道系統(tǒng)實(shí)物圖如圖2所示，包含支管接頭和模擬泄漏點(diǎn)，其中模擬泄漏點(diǎn)為DN15泄漏閥(泄漏孔直徑為15 mm)。由圖2可知，100 mm管徑管道的模擬泄漏點(diǎn)距離支管接頭為60 m。我們將管內(nèi)壓強(qiáng)控制在0.6 MPa左右，利用集成有加速度計(jì)的數(shù)據(jù)采樣設(shè)備以采樣率10 kHz對(duì)信號(hào)進(jìn)行采樣。

圖2 測(cè)試管道系統(tǒng)示意圖

本試驗(yàn)的主要過(guò)程如下：

步驟1打開模擬泄漏點(diǎn)，放置傳感器在泄漏點(diǎn)附近，獲取含有支管流致噪聲的泄漏信號(hào)；

步驟2關(guān)閉模擬泄漏點(diǎn)，放置傳感器在支管接頭上，僅獲取支管流致噪聲；

步驟3分別計(jì)算含支管流致噪聲的泄漏信號(hào)和支管流致噪聲的頻譜。

100 mm管徑的含噪泄漏信號(hào)和支管流致噪聲頻譜，如圖3所示。從圖3中可知，支管流致噪聲的主譜峰位于100～2 000 Hz頻帶內(nèi)；而含噪泄漏信號(hào)的主頻帶位于120～3 300 Hz，其中500～1 500 Hz頻帶內(nèi)(小矩形框內(nèi))具有一個(gè)顯著的譜峰，這個(gè)顯著的譜峰則是由支管流致噪聲引起的。上述結(jié)果表明：① 泄漏信號(hào)的主頻帶要比支管流致噪聲的主頻帶更寬，進(jìn)一步說(shuō)明泄漏信號(hào)的隨機(jī)性要比支管流致噪聲的隨機(jī)性更強(qiáng)；② 當(dāng)管道存在支管接頭時(shí)，支管流致噪聲對(duì)于泄漏檢測(cè)而言是一種不可忽略的噪聲干擾。支管流致噪聲和泄漏信號(hào)的概率密度分布圖，如圖4所示。圖4中支管流致噪聲和泄漏信號(hào)的概率密度均呈現(xiàn)中心部分又窄又高，而尾部拖得很長(zhǎng)的特征，這是典型的超高斯信號(hào)概率密度特征。為此我們引入一個(gè)四階矩來(lái)度量信號(hào)的高斯性

(a)支管流致噪聲頻譜

圖4 支管流致噪聲信號(hào)、含噪泄漏聲信號(hào)概率分布

(1)

(2)

圖5 泄漏聲信號(hào)和支管流致噪聲的值計(jì)算結(jié)果

2 算法理論背景

根據(jù)泄漏特征，我們建立泄漏時(shí)域卷積模型如下

(3)

式中：s(n)=sleak(n)+sbranch(n)，n1(n)和n2(n)為背景噪聲；a(n)和b(n)為系統(tǒng)混合矩陣，為了使計(jì)算更方便快捷，我們將時(shí)域卷積的混合模型變換到頻域中的乘積模型，再進(jìn)行頻域內(nèi)的盲分離。頻域混合模型如下

(4)

通常將z省略，則有

(5)

2.1 C-EFastICA算法

根據(jù)式(5)，暫不考慮背景噪聲N1和N2，即認(rèn)為背景噪聲在進(jìn)行盲分離以前已經(jīng)通過(guò)其他方法技術(shù)降低到可以忽略的程度[18-19]，此時(shí)式(5)簡(jiǎn)化為

(6)

略去下標(biāo)，將混合矩陣A與B用混合矩陣C替換整理得到

x=CS

(7)

信號(hào)源盲分離估計(jì)就是指在混合矩陣C和源信號(hào)S未知的條件下，求一個(gè)分離矩陣W，使得W對(duì)觀測(cè)信號(hào)x有如下線性變換

u=Wx

(8)

設(shè)觀測(cè)信號(hào)和分離矩陣、輸出信號(hào)均為復(fù)數(shù)形式

x=xr+jxi

(9)

W=Wr+jWi

(10)

u=ur+jui

(11)

令xH為x的復(fù)共軛轉(zhuǎn)置，wH為w的復(fù)共軛轉(zhuǎn)置，且它們的實(shí)部和虛部不相關(guān)，則有E(wwH)=I, 設(shè)此時(shí)的輸出信號(hào)模型為u=wHx，N個(gè)向量wi按列組成了W。

本算法第一步是對(duì)信號(hào)去均值處理和預(yù)白化處理

(12)

第二步，根據(jù)信號(hào)的廣義高斯特征選擇非線性g函數(shù)。非線性函數(shù)有以下3種可供選擇

(13)

由于泄漏信號(hào)和支管流致噪聲均是超高斯信號(hào)，則此處選擇非線性g函數(shù)為g(x)=xexp(-η3|x|)，此函數(shù)在超高斯樣本空間內(nèi)連續(xù)可導(dǎo)，令參數(shù)η3為3.348。

(14)

(15)

最終整理得到

(16)

(17)

下一步是計(jì)算

(18)

(19)

(20)

(21)

(22)

式中，k=1,…,d。計(jì)算

(23)

式中，k，l=1,…,d。接下來(lái)計(jì)算

(24)

(25)

(26)

(27)

算法1Complex-Efficient FastICA(C-EFastICA)

初始化w=randn+j·randn;

根據(jù)g(x)=x/(1+x2);

Foriter≥0

計(jì)算

收斂、停機(jī)！

3 泄漏定位檢測(cè)試驗(yàn)

3.1 C-EFastICA算法運(yùn)用于管道泄漏定位

本節(jié)將利用加速度傳感器采集兩路含支管流致噪聲的泄漏信號(hào)，再利用C-EFastICA算法對(duì)含支管流致噪聲的泄漏信號(hào)進(jìn)行分離，主要的泄漏檢測(cè)步驟如圖6所示。將根據(jù)該流程圖的流程，利用實(shí)際的100 mm管徑供水管道進(jìn)行試驗(yàn)研究。

圖6 泄漏定位檢測(cè)流程圖

3.2 試驗(yàn)條件

試驗(yàn)選擇直徑為100 mm、長(zhǎng)度70 m的鑄鐵供水管道采集泄漏信號(hào)，管道內(nèi)部壓力范圍為0.55～0.60 MPa，入口流速約為2 m/s，模擬泄漏點(diǎn)以DN15小口閥門(泄漏孔徑15 mm)代替，且模擬泄漏點(diǎn)距離傳感器1有42 m，信號(hào)采集設(shè)備采用自主研制的泄漏信號(hào)數(shù)據(jù)采集設(shè)備，數(shù)據(jù)采樣率設(shè)置為10 kHz，詳細(xì)參數(shù)設(shè)置如表1所示。泄漏檢測(cè)系統(tǒng)示意圖，如圖7所示。兩路傳感器被放置在管道兩端盡頭，采集到的泄漏信號(hào)經(jīng)無(wú)線傳輸?shù)街鳈C(jī)。兩路傳感器和支管接頭實(shí)物圖，如圖8所示。

圖7 泄漏檢測(cè)系統(tǒng)示意圖(m)

(a)

表1 試驗(yàn)中的詳細(xì)參數(shù)條件

3.3 泄漏定位過(guò)程

兩路傳感器1、傳感器2采集得到的含支管流致噪聲的泄漏信號(hào)頻域形式如圖9所示。

(a)

為了證明C-EFastICA算法的分離性能，我們利用已知的C-FastICA算法[20]對(duì)含噪泄漏信號(hào)進(jìn)行分解處理。在對(duì)含噪泄漏信號(hào)進(jìn)行分解前，必須對(duì)含噪泄漏信號(hào)進(jìn)行預(yù)處理，抑制加性噪聲n1和n2對(duì)信號(hào)分解的影響，更多關(guān)于抑制加性噪聲的細(xì)節(jié)詳見(jiàn)朱建斌等和張艷艷等的研究。C-EFastICA算法、C-FastICA算法對(duì)含噪泄漏信號(hào)(以傳感器1采集的源數(shù)據(jù)為例進(jìn)行分解處理)的分解結(jié)果分別如圖10和圖11所示。

(a)

(a)

為了評(píng)估算法的收斂速度，碼間串?dāng)_(inter-symbol interference, ISI)的概念被引入，其定義如下

(28)

式中：cij為全局系統(tǒng)矩陣C中的第(i,j)個(gè)元素，當(dāng)每個(gè)通道的ISI輸出曲線越接近0的時(shí)候，算法越收斂。C-EFastICA和C-FastICA的兩路通道的ISI輸出曲線，如圖12所示。從圖12可知，C-EFastICA在300步左右就收斂了，而C-FastICA則要在400步左右時(shí)才收斂，由此可見(jiàn)，C-EFastICA的收斂速度要快于C-FastICA。

(a)

此處引入串音誤差PI來(lái)比較C-EFastICA算法和C-FastICA算法的分離性能，PI值越接近0說(shuō)明算法對(duì)混合信號(hào)的分離程度越高，PI定義如下

(29)

式中：M為源信號(hào)的數(shù)量；cik為C=W·A的元素，A為混合矩陣，W為算法分離得到的分離矩陣。將C-EFastICA算法和C-FastICA算法分別對(duì)含支管流致噪聲的泄漏信號(hào)再進(jìn)行50次分離試驗(yàn)，得到試驗(yàn)結(jié)果如圖13所示。

圖13 C-EFastICA與C-FastICA分別進(jìn)行50次盲分離試驗(yàn)的PI值

從圖13可知，C-EFastICA的PI系數(shù)均小于C-FastICA的PI系數(shù)，C-EFastICA的PI系數(shù)中位數(shù)在0.125左右，而C-FastICA的PI系數(shù)中位數(shù)在0.225左右，即可表明C-EFastICA對(duì)含噪信號(hào)的分離程度越徹底，分離出的源泄漏信號(hào)所含泄漏特征更多。

u1的譜成分主要分布于1 600～4 500 Hz的有效頻帶內(nèi)，而u2的譜成分主要位于1 000～2 300 Hz區(qū)間，見(jiàn)圖10。圖11中，分離信號(hào)u1和u2相比，u1的頻譜特征位于1 000～4 500 Hz，u2的譜特征位于1 000～2 200 Hz。為了從分離結(jié)果中區(qū)分泄漏源信號(hào)和噪聲源信號(hào)，我們對(duì)分離結(jié)果計(jì)算了樣本熵[21]，以熵值的大小差別區(qū)分泄漏源信號(hào)和噪聲源信號(hào)。這是由于源泄漏聲信號(hào)較支管流致噪聲往往具有更大的隨機(jī)性，即源泄漏信號(hào)的有效主頻帶比支管流致噪聲的主頻帶更寬。從圖3也可看出泄漏聲信號(hào)的頻率成分比支管流致噪聲也更豐富，所以泄漏信號(hào)的熵值應(yīng)大于支管流致噪聲。計(jì)算C-EFastICA算法、C-FastICA算法對(duì)含噪泄漏信號(hào)分解得到的u1和u2的樣本熵(SampEn)值如表2所示。在表2中，u1的樣本熵均大于1.2，而u2均不超過(guò)0.9，因此我們認(rèn)為u1是泄漏源信號(hào)，u2是支管流致噪聲。

表2 不同算法分離結(jié)果的樣本熵值

以同樣方式對(duì)傳感器2采集得到的含噪泄漏信號(hào)進(jìn)行處理，得到另一路源泄漏信號(hào)，對(duì)于得到的兩路源泄漏信號(hào)，利用互相關(guān)技術(shù)對(duì)它們進(jìn)行處理。如圖7所示，傳感器1距離泄漏點(diǎn)的距離L1=42 m，已知泄漏聲速在鑄鐵管道中的傳播速度為c=1 495 m/s，時(shí)延估計(jì)τ由互相關(guān)估計(jì)得到，因此泄漏點(diǎn)距離傳感器1的距離由下式計(jì)算得到

L1=(L-cτ)/2

(30)

式中，L=L1+L2+L3=70 m。C-EFastICA、C-FastICA算法分解得到的泄漏源信號(hào)的互相關(guān)估計(jì)定位結(jié)果圖，如圖14(a)、圖14(b)所示。兩路混合信號(hào)經(jīng)傳統(tǒng)的直接互相關(guān)算法(direct cross-correlation, DCC)的定位結(jié)果圖，如圖14(c)所示。

(a)C-EFastICA算法的互相關(guān)估計(jì)定位結(jié)果

由圖14(a)、圖14(b)可知，只有在圖14(a)中的峰值最接近于真實(shí)泄漏點(diǎn)距離。對(duì)比圖14(a)和圖14(b)，在圖14(c)中的定位結(jié)果誤差最大，這是由于未對(duì)泄漏信號(hào)中的支管流致噪聲進(jìn)行處理，支管流致噪聲對(duì)互相關(guān)估計(jì)產(chǎn)生了較大的影響。為了進(jìn)一步驗(yàn)證C-EFastICA對(duì)混合泄漏信號(hào)的分離性能優(yōu)于C-FastICA，本試驗(yàn)被重復(fù)數(shù)次，定位結(jié)果如表3所示。

表3中，δ為相對(duì)誤差，C-EFastICA的定位相對(duì)誤差均不超過(guò)12%，而C-FastICA算法的定位誤差達(dá)到15%，所以C-EFastICA算法的定位精度高于C-FastICA算法。而直接互相關(guān)定位法由于受噪聲影響嚴(yán)重，定位相對(duì)誤差達(dá)到60%以上，甚至無(wú)法完成定位。

表3 不同算法對(duì)分支管道泄漏點(diǎn)的定位相對(duì)誤差

4 結(jié) 論

針對(duì)分支泄漏管道中存在支管流致噪聲對(duì)漏點(diǎn)定位干擾的問(wèn)題，文中提出使用基于EFastICA算法的C-EFastICA對(duì)含支管流致噪聲的泄漏信號(hào)進(jìn)行分離計(jì)算，得到源泄漏信號(hào)，結(jié)合先驗(yàn)聲速和管道長(zhǎng)度完成漏點(diǎn)定位，通過(guò)試驗(yàn)證明該方法的可行性，得到以下結(jié)論：

(1)當(dāng)管道中存在分支接頭時(shí)，即泄漏信號(hào)中混有支管流致噪聲時(shí)，采用傳統(tǒng)的直接互相關(guān)方法對(duì)管道進(jìn)行漏點(diǎn)定位誤差較大，甚至無(wú)法完成定位。說(shuō)明了分支管道泄漏定位中，支管流致噪聲對(duì)定位的影響不可忽視。

(2)試驗(yàn)證明，泄漏信號(hào)和支管流致噪聲均是超高斯信號(hào)。與其他ICA類算法固定選擇非線性函數(shù)不同，C-EFastICA能根據(jù)信號(hào)的廣義高斯特性自適應(yīng)地選擇非線性函數(shù)構(gòu)建迭代學(xué)習(xí)規(guī)則，即C-EFastICA算法是基于混合泄漏信號(hào)的物理特征進(jìn)行分解運(yùn)算的。

(3)所提出的C-EFastICA算法相較于同類型的C-FastICA具有更快的收斂速度，在對(duì)混合信號(hào)的分解程度上也具有更高的水平。