基于希爾伯特黃變換的輸氣管道泄漏音波時頻特性分析

劉翠偉，李玉星，孟令雅，孫玉萍

(1.中國石油大學 儲運與建筑工程學院，山東 青島 266555；2.中國石油大學 信息與控制工程學院，山東 青島 266555)

隨天然氣工業的迅速發展，天然氣長輸管道與城市燃氣輸配管道鋪設及使用亦越多。輸氣管道泄漏檢測技術[1-3]研究及應用對輸氣管線系統安全運行起重要作用。音波泄漏檢測技術[4-6]具有靈敏度高、定位精度高、誤報率低、檢測時間短、適應性強等優點；音波法可測量管線流體中微弱動態壓力變化，與管線運行壓力絕對值無關；音波法響應頻率、檢測范圍更寬。雖輸氣管道音波法泄漏檢測技術優勢明顯，但對泄漏音波信號進行處理獲得直觀有效音波特征量成為制約該技術用于實際的瓶頸。

希爾伯特黃變換(HHT)[7]為對時間序列數據先進行經驗模態分解，再對各分量信號作希爾伯特變換。與傳統的信號處理方法相比HHT變換具有諸多優勢，即能分析非線性非平穩信號、具有完全自適應性、不受Heisenberg測不準原理制約、適合突變信號、瞬時頻率用求導獲得等。該方法用于地震波檢測、生物醫學等領域研究[8-10]效果較好。

為提取輸氣管道音波法泄漏檢測技術有效的音波特征量及盡快用于實際工程，搭建實驗室設備進行實驗研究[11]，并對實驗結果進行HHT分析，建立完整的HHT處理音波信號理論體系；選勝利油田孤島集中處理站與集賢集中處理站一段管線進行現場實驗，并對實驗信號進行HHT分析，驗證音波法泄漏檢測技術的現場應用效果。

1 希爾伯特黃變換理論

希爾伯特黃變換為自適應的時頻分析法，用其進行數據分析主要步驟為，①用EMD方法分解信號，提取原信號中固有模態函數分量；對分解所得各固有模態函數分量作Hilbert變換，獲得其Hilbert譜與Hilbert邊際譜，從而進行分析處理。

1.1 EMD基本原理

經驗模態分解(EMD)法能將非平穩、非線性信號分解成一組穩態、線性數據序列集，即本質模態函數(IMF)。以時間序列數據s(t)為例，經驗模態分解算法最后可得：

(1)

式中：ck(t)為IMF分量；rn(t)為平均趨勢分量。

1.2 瞬時頻率及HHT譜

Hilbert-Huang變換方法基于信號局部特征時間尺度，將信號自適應分解為若干IMF分量之和，計算每個IMF分量的瞬時頻率及瞬時幅值。

對式(1)每個固有模態函數ck(t)作Hilbert變換得：

(2)

通過構造解析信號，求出幅值函數、相位函數，進而求出瞬時頻率為

(3)

忽略殘余項rn，其數據可表示為

(4)

式中：頻率ωj(t)、幅值aj(t)為時間的變量，可構成時間、頻率、幅值三維時頻譜圖，此在時頻平面的信號幅度分布稱為Hilbert時頻譜，簡稱Hilbert譜，以H(ω,t)表示。其可精確描述信號幅值在整個頻率段隨時間、頻率的變化規律。

邊際譜h(ω)定義為

(5)

式中：T為序列時間長度。

將三維時頻譜對時間積分，便形成僅有頻率、幅值的二維譜圖。邊際譜可由統計意義表征整組數據每個頻率點的積累幅值分布，反映信號幅值在整個頻率段的變化。

2 實驗室信號HHT分析

為研究希爾伯特黃變換在輸氣管道泄漏音波信號處理中的應用及能獲得泄漏檢測有用的音波特征量，搭建實驗管道[11]，采集實驗數據，對實驗數據進行希爾伯特黃變換分析。在距每個泄漏點10 cm處安裝PCB公司M106B音波傳感器，以確保能檢測到泄漏點音波信號，采樣率3 kHz。音波傳感器工作原理及數據采集系統見圖1。音波數據采集系統主要采集沿線布置的音波傳感器數據。

圖1 音波檢測系統

在輸氣管道上安裝音波傳感器時應使受力隔膜置于水平，以便輸氣管線壓力穩定時使隔膜受力平衡，不引起電壓輸出，讀數為0；輸氣管線發生泄漏瞬間，因音波擾動使受力隔膜受聲壓作用，受力隔膜變形使石英晶體內部產生極化并產生符號相反電荷從而產生電壓值，隨聲壓增大受力隔膜變形增大，電壓增大，直到聲壓達到幅值，電荷開始釋放輸出電流信號，隨電流信號輸出電荷電量逐漸減小，從而使輸出電流信號逐漸減小，直至變0；輸氣管線持續泄漏時，聲壓始終作用于受力隔膜，但此時隔膜受力平衡不再發生變形，石英晶體不產生極化及電荷，因此無電流輸出，讀數為0。

音波傳感器采集的信號經自身前置放大后，通過同軸電纜002C30與442B104型信號調理器相連，主要為音波傳感器提供工作電源，去除信號中直流成分；內置有模擬濾波器，可去除高頻無用信號。從信號調理器輸出的信號通過電纜與數據采集卡相連，由于采集的為動態數據，因此實驗室選NI公司PCI-4474動態數據采集卡，該采集卡有4個通道，采樣率高達102.4 kS/s，分辨率24 bits。所采信號最終傳至監控主機，通過數據采集與監控系統處理，實現泄漏檢測及定位。

泄漏檢測系統的數據采集、處理主要通過LabVIEW軟件實現。音波數據采集系統主要采集管線沿線的動態壓力數據，并對數據進行實時存儲。

2.1 穩定信號與泄漏音波信號對比分析

穩定運行時取音波信號及壓力為1 MPa、泄漏孔徑0.6 mm、點1泄漏時音波信號分別進行HHT分析獲得信號幅頻特性、相頻特性、瞬時頻率、EMD分量、HHT譜及邊際譜，通過對比分析，提取用于長距離輸氣管道泄漏音波信號檢測的特征量。

2.1.1 原始信號

由于實驗室條件理想，所含背景噪音較小，見圖2。由圖2看出，穩定時音波信號波動幅值處于0.01 kPa及0.1 kPa數量級；泄漏時音波信號在泄漏時刻存在一突變值，下降值達13.5 kPa。

音波信號屬于突變信號，時域分析只能顯示其波形與幅值隨時間變化關系，用于泄漏檢測的音波特征量為動態壓力幅值，不能全面描述信號特征。為進一步分析音波信號從而提取泄漏檢測可用的特征量對信號進行HHT分析。獲得幅頻特性、相頻特性、瞬時頻率、EMD分量、HHT二維譜能量、HHT三維譜能量及邊際譜特性7個音波特征量，并分析該特征量是否可用于泄漏檢測。

2.1.2 幅頻曲線

由圖3看出，穩定與泄漏時比較知，穩定時幅值在頻域上均勻分布；泄漏時幅值集中在低頻段，尤其在0～100 Hz低頻段區間，泄漏時幅值數量級遠大于穩定時的幅值數量級；穩定時、泄漏時均在頻率1 000～1 500 Hz間存在極值點，但由于此頻段信號不能遠傳，因此，高頻段幅頻值并非有效特征量。低頻段幅頻特性為有效特征量。

2.1.3 相頻曲線

由圖4看出，通過比較穩定信號相頻曲線與泄漏信號相頻曲線，無法獲得能明顯區分泄漏信號結論，因此相頻特性亦非理想中音波特征量。

圖2 實驗室原始信號圖

2.1.4 瞬時頻率

由圖5看出，通過比較穩定信號瞬時頻率與泄漏信號瞬時頻率，仍無法獲得能明顯區分泄漏信號結論，因此瞬時頻率亦非理想音波特征量。

2.1.5 EMD分量

由圖6看出，穩定時EMD分量即各IMF無明顯特征的平穩波動，且波動范圍較小；泄漏時EMD分量即各IMF均在泄漏時刻有明顯振蕩特性，其它時刻為平穩波動，因此EMD分量中IMF可作為音波特征量，具體哪個IMF可作為音波特征量需進一步研究。

通過以上對穩定工況信號及泄漏工況音波信號的時域特征、頻域特征分析，獲得頻域內泄漏信號特征。為獲得信號的頻域特征，須用一段時間內信號時域內所有信息，且當信號在某時刻小鄰域內發生變化時，信號整個頻譜會受影響，而對頻譜變化無法標定發生變化時間及劇烈程度，即無法獲得泄漏信號的產生時刻及該時刻對應幅度。實驗室因外界干擾噪聲較少，因此可通過觀測時域波形內幅值對應時刻判斷泄漏發生時刻，但現場存在較多噪聲影響時，時域內奇異點判斷不再準確，且無法獲得泄漏時刻頻域特征。因此可采用時頻聯合分析方法將信號轉換到時間-頻率二維平面或時間-頻率-能量三維平面，使處理后信號不僅含頻率信息及對應時間信息、能量或幅值信息。由此獲得到HHT二維譜能量、HHT三維譜能量等特征量。

2.1.6 HHT譜二維圖

由圖7可知，或穩定或泄漏時，音波信號能量在頻域上集中在低頻段，穩定時能量在時間域分布均勻，泄漏時能量在時間域分布集中在泄漏時刻；而穩定時信號在10～20 Hz之間也存在散亂分布，而泄漏信號集中在0～3 Hz頻段，且其它頻段分布非常均勻；而HHT譜二維圖僅能看出信號時頻分布，信號能量分布無法表示，因此二維HHT譜非理想音波特征量。

圖5 實驗室信號瞬時頻率圖

2.1.7 HHT譜三維圖

由圖8的HHT三維譜看出，穩定時音波信號能量值在頻域上集中在低頻段，而在時域上均勻分布，能量范圍處于0.01 kPa數量級；泄漏時音波信號能量值在頻域上仍集中在低頻段，在時域上存在突變點且位于泄漏時刻，泄漏發生前音波信號能量在時域上分布均勻，泄漏發生后音波信號能量在時域上同樣分布均勻，突變處能量范圍處于10 kPa數量級。因此，HHT三維譜將音波信號表現在時間-頻率-能量三維平面上，使處理后的信號不僅含頻率信息，也對應時間信息及能量或幅值信息，可直觀檢測泄漏信號的音波特征量。

2.1.8 邊際譜

由圖9看出，低頻段區間，泄漏時幅值數量級遠大于穩定時幅值數量級，因此，低頻段邊際譜為有效特征量。

由于實驗采集的動態壓力信號在穩定工況或泄漏工況，均非由單一頻率構成，而由不同頻率諧波信號疊加而成，該諧波信號幅值變化與相位變化均會最終作用于動態壓力信號引起信號特征變化，因此需研究動態壓力信號頻譜特征，明確泄漏引起的動態壓力信號所處能量占優頻段及幅度能量。綜上，音波特征量經分析提取有用音波特征量為低頻段幅頻特性、EMD分量、HHT三維譜及低頻段邊際譜，其中低頻段幅頻特性、低頻段邊際譜雖能明確區分是否發生泄漏，但無法確定時域及泄漏發生時刻，EMD分量具體選取哪個分量需進一步研究分析且可用信息復雜多變，而HHT三維譜在時域上直觀明確顯示泄漏發生時刻，在頻域上表現可提取用于遠傳的低頻段，在能量域也能使泄漏、穩定工況表現出明顯差異，因此HHT三維譜為最直觀有效的音波特征量，故選HHT三維譜進行現場實驗信號HHT分析。

圖8 實驗室信號HHT譜三維圖

圖9 實驗室信號邊際譜

2.2 不同工況泄漏音波信號分析

通過不同泄漏工況實驗，獲得同一泄漏點發生泄漏時泄漏孔徑0.45 mm，0.6 mm，0.7 mm，0.8 mm，運行壓力1 MPa，2 MPa，3 MPa，4 MPa，5 MPa，距泄漏點0.1 m，48.09 m，108.78 m，159.51 m處HHT三維譜特征值為HHT三維譜在某低頻段范圍內的泄漏幅值。

由圖10看出，①同一泄漏孔徑下，隨傳播距離增大HHT三維譜特征量值逐漸減小，減小趨勢不隨距離呈單一規律性，且不同壓力下變化趨勢基本一致。②同一壓力下，隨傳播距離的增大HHT三維譜特征量值逐漸減小；隨泄漏孔增大HHT三維譜特征量值逐漸增大。③同一傳播距離下，隨泄漏孔徑增大HHT三維譜特征量值逐漸增大；隨壓力增大HHT三維譜特征量值逐漸增大。

3 現場實驗信號HHT分析

據實驗管道所得實驗結果確定HHT三維譜為音波特征量，選勝利油田孤島集中處理站與集賢集中處理站一段管線進行現場實驗，探究HHT三維譜作為音波特征量的現場應用效果。先選集賢集中處理站內一段管線進行實驗，獲得同站實驗結果后選孤島站至集賢站間一段長約28～30 km管線進行實驗，檢驗HHT三維譜作為音波特征量效果。

3.1 同站實驗

3.1.1 現場實驗工況

勝利油田孤島站來氣經分離器分離，穩定后進入計量站計量，后經管道輸送東營，泄漏通過泄漏點處球閥打開實現且泄漏處閥門內徑與管線內徑均10 mm，泄漏點位置見圖11，在距泄漏點下游約65 cm處接收音波信號。管線運行工況為壓力1 MPa，瞬時流量15 720 m3/h，管線規格Φ426×6 mm，溫度27.5°。

圖11 集賢站所選管線流程圖

同站實驗為：① 記錄管線不泄漏時音波信號，采樣頻率1 000 Hz；② 打開球閥，記錄泄漏時音波信號，采樣頻率1 000 Hz，將泄漏信號及穩定信號進行HHT求解、對比分析；③ 記錄先打開閥門，后關閉閥門過程的音波信號并處理分析，閥門開啟與關閉時間差20 s，采樣頻率3 000 Hz。具體處理過程如下。

3.1.2 泄漏音波信號的HHT分析

3.1.2.1 泄漏信號與穩定信號對比分析

(1) 原始信號。管線穩定時音波信號在0值附近波動且波動較小，但由于現場工況背景噪聲強烈，使壓力波動值較實驗室工況壓力穩定時壓力波動值大；管線發生泄漏瞬間，管內壓力平衡被破壞，音波信號有急劇下降沿，下降值達0.23 kPa，突變到一定幅值后迅速恢復至0附近，產生較小波動后，動態壓力值基本趨于0恢復初始狀態。

圖12 現場實驗原始信號

(2) HHT譜三維圖。將所得音波信號進行HHT求解，獲得穩定及泄漏時音波信號的HHT譜三維圖，見圖13。由圖13可知，穩定時音波信號能量值集中于頻域的低頻段，在時域上均勻分布，能量范圍屬0.01數量級；泄漏時，音波信號能量值仍集中于頻域低頻段，在時域上存在突變點且位于泄漏時刻，泄漏發生前音波信號能量在時域上分布均勻，泄漏發生后，音波信號能量在時域上亦分布均勻，能量范圍屬0.1數量級，突變處能量值達0.16 kPa；因此，HHT三維譜不僅可在時域上直觀檢測泄漏時刻，且可提取能量值占優的遠傳頻段。

3.1.2.2 泄漏發生與閥門關閉音波信號分析

(1)原始信號。由圖14看出，穩定時信號在0值附近均勻波動，泄漏發生時信號產生一下降沿，閥門關閉時信號產生一上升沿，且閥門關閉點與泄漏點間時間差為19.2 s，與設定的20 s接近，誤差源于手動操作及時域計時的不精確。

圖13 現場實驗HHT譜三維圖

圖14 閥門開啟與關閉時音波信號原始信號圖

圖15 閥門開啟與關閉時音波信號HHT譜三維圖

(2)HHT譜三維圖。由圖15看出，該三維圖存在兩突變點(30 590，1，0.160 6)，(88 670，1，0.137 9)，即音波信號能量存在兩突變值且較閥門關閉時刻相比泄漏時刻突變值更大，兩突變值間時間差為19.36 s，與所設20 s更接近，較時域的突變值時間差更準確，提取可用于遠傳頻段為0～1 Hz。

3.2 長輸實驗

3.2.1 現場實驗工況

選孤島站至集賢站間28～30 km管線，在孤島站出口管線放氣模擬泄漏，打開泄漏點處球閥，泄漏處閥門內徑與管線內徑均為10 mm。在集賢站進口采集信號，采樣頻率1 000 Hz，布置見圖16。管線運行工況：孤島站管線運行壓力1.075 MPa，集賢站運行壓力1 MPa，瞬時流量15 720 m3/h，管線規格Φ426×6，運行溫度為27.5°。記錄閥門先開后關過程的音波信號并處理分析，閥門開、關時間差30s。具體處理過程見圖16。

圖16 孤島站至集賢站長輸實驗

3.2.2 泄漏音波信號HHT分析

(1) 原始信號。由于音波信號傳播28～30 km，使泄漏時刻產生的音波信號隱藏于背景噪聲中，僅從時域觀察是否發生泄漏變困難，見圖17。

圖17 長輸管道實驗原始信號

圖18 長輸管道實驗音波信號HHT譜三維圖

(2) HHT譜三維圖。將所得音波信號進行HHT求解，獲得音波信號的HHT譜三維圖，見圖18。由圖18看出，該圖存在一突變點(56 400，1，0.174 4)，即音波信號能量存在一突變值，考慮同站實驗時音波信號HHT譜三維圖的泄漏時刻音波信號能量突變值強于閥門關閉時刻音波信號能量突變值，長輸管道對音波信號能量衰減作用相同，音波信號能量突變值較大時刻即泄漏時刻。此時閥門關閉時信號難以檢測，此因閥門關閉前，由于氣體可壓縮性及管道彈性等因素，使泄漏后管道內動態壓力未及時歸零而在某負值附近波動，即閥門關閉前音波信號已占有能量，因此在閥門關閉后壓力上升引起的能量變化并不明顯，無法觀測閥門關閉時刻。

4 結 論

通過適當運用希爾伯特黃變換信號處理方法，對實驗室及現場實驗所得音波信號進行HHT分析，結論如下：

(1) 綜合分析HHT變換獲得音波特征量，提取可用音波特征量為EMD分量、低頻段幅頻特性、HHT三維譜及低頻段邊際譜，其中HHT三維譜為最直觀有效的音波特征量，選HHT三維譜進行現場實驗信號HHT分析。將泄漏音波信號通過HHT分析轉換到時間-頻率二維平面或時間-頻率-能量三維平面，使處理后信號不僅含頻率信息，亦對應時間信息及能量或幅值信息。

(2) 通過同站實驗中泄漏信號與穩定信號對比分析，HHT三維譜可直觀檢測泄漏信號；通過同站實驗泄漏發生與閥門關閉音波信號分析，HHT三維譜能量突變點明顯，能提高低頻段時頻分辨率。

(3) 通過長輸管線實驗驗證表明HHT三維譜為有效音波特征量的可行性及較高自適應性。

[1] 張新年．天然氣管道泄漏檢測技術綜述[J]．內蒙古石油化工,2009(22):106.

ZHANG Xin-nian. Summarization oninspection technique for leakage in natural gas pipelines[J]. Inner Mongulia Petrochemical Industry,2009(22):106.

[2] 蘇欣,袁宗明,范小霞,等．油氣長輸管道檢漏技術綜述[J]．石油化工安全技術,2005,21(4):15-16.

SU Xin，YUAN Zong-ming，FAN Xiao-xia, et al. Review of leak hunting technology for long distance oil-gas pipeline[J]. Petrochemical Safety Technology,2005,21(4):15-16.

[3] 楊向蓮．天然氣管道泄漏檢測技術評價及預防措施[J]．能源技術,2005,26(6):248-251.

YANG Xiang-lian. The leakage detection technique of natural gas pipeline and countermeasures[J]. Energy Technology,2005，26(6):248-251.

[4] 李玉星,彭紅偉,唐建峰,等．天然氣長輸管道泄漏檢測方案對比[J]. 天然氣工業,2008,28(9):101-104.

LI Yu-xing，PENG Hong-wei，TANG Jian-feng,et al. Comparative study on leakage detecting techniques for long-distance natural gas pipelines[J]. Natural Gas Industry,2008,28(9):101-104.

[5] 孟令雅,李玉星,宋立群,等. 輸氣管道泄漏音波傳播特性及監測定位[J]. 天然氣工業,2010,30(11):74-79.

MENG Ling-ya,LI Yu-xing,SONG Li-qun, et al. Acoustic propagation characteristics, position monitoring and locating of gas transmission pipeline leakage[J]. Natural Gas Industry,2010,30(11):74-79.

[6] 劉翠偉,李玉星,李雪潔.基于CFD模擬的輸氣管道閥門流噪聲仿真[J].油氣儲運,2012, 31(9):657-662.

LIU Cui-wei，LI Yu-xing, LI Xue-jie, et al. Simulation of flow noise of valves in gas pipelines on CFD software[J]. Oil & Gas Storage and Transportation,2012,31(9):657-662.

[7] Huang N E，Shen Z，Long S R，et al．The empirical mode decomposition and the Hilbert spectrum for nonlinear and nonstationary time series analysis[J]. Proc. R. Soc. Lond,1998,A: 903-995.

[8] 楊世錫,胡勁松,吳昭同,等. 旋轉機械振動信號基于EMD的希爾伯特變換和小波變換時頻分析比較[J]. 中國電機工程學報,2003,23(6):102-107.

YANG Shi-xi,HU Jin-song,WU Zhao-tong, et al. The comparison of vibration signals’ time-frequency analysis between EMD-based ht and wt method in rotating machinery[J]. Proceedings of the Csee,2003,23(6):102-107.

[9] 侯斌,桂志先,胡敏,等. 基于希爾伯特-黃變換的地震信號時頻分析[J]. 勘探地球物理進展,2009,32(4):248-251.

HOU Bin，GUI Zhi-xian，HU Min, et al. Time-frequency spectral analysis of seismic data based on Hilbert-Huang transform[J]. Progress in Exploration Geophysics, 2009,32(4):248-251.

[10] 王黎黎. 基于希爾伯特-黃變換的時頻分析算法研究[D]. 西安：西安電子科技大學,2009.

[11] 李玉星,唐建峰,王武昌,等.輸氣管道聲波泄漏檢測試驗裝置的設計與構建[J].實驗技術與管理,2010,27(8):63-67.

LI Yu-xing,TANG Jian-feng,WANG Wu-chang, et al. Design and construction of an experimental device for acoustic leak detection in gas pipelines[J]. Experimental Technology and Management,2010,27(8):63-67.