周 攀，周東岳，鄧 磊

（國網(wǎng)新源控股有限公司技術(shù)中心，北京市 100073）

基于EMD方法的水力測(cè)量信號(hào)分析

周 攀，周東岳，鄧 磊

（國網(wǎng)新源控股有限公司技術(shù)中心，北京市 100073）

壓力脈動(dòng)具有非平穩(wěn)的特點(diǎn)，尤其在過渡過程中表現(xiàn)出非常強(qiáng)的非時(shí)變特性；應(yīng)用EMD方法進(jìn)行壓力脈動(dòng)信息提取解決了現(xiàn)場(chǎng)實(shí)測(cè)過程中壓力脈動(dòng)特性難以確定的問題，采用基于連續(xù)均方誤差法有效地解決了壓力脈動(dòng)信號(hào)難以重構(gòu)的問題。本文以現(xiàn)場(chǎng)機(jī)組甩負(fù)荷試驗(yàn)數(shù)據(jù)為計(jì)算基礎(chǔ)，利用EMD方法實(shí)現(xiàn)了現(xiàn)場(chǎng)實(shí)測(cè)壓力信號(hào)的分解和脈動(dòng)信息的提取，為分析機(jī)組各部位的壓力和脈動(dòng)信息提供一個(gè)重要手段，證明了該方法對(duì)實(shí)測(cè)信號(hào)的分析具有很好的實(shí)用性。

過渡過程；經(jīng)驗(yàn)?zāi)B(tài)分解；固有模態(tài)函數(shù)；壓力脈動(dòng)；水力量測(cè)

0 引言

壓力脈動(dòng)是作用于被作用對(duì)象上并不均勻的，在某個(gè)部位有較集中的或是較大的壓力，且這種壓力單次持續(xù)的時(shí)間不長，有可能呈現(xiàn)一定的周期性。水力壓力脈動(dòng)是水輪機(jī)運(yùn)行時(shí)不可避免的現(xiàn)象，是導(dǎo)致機(jī)組振動(dòng)的主要原因之一，在一定條件下甚至?xí)鸾Y(jié)構(gòu)的固有振動(dòng)，嚴(yán)重影響到機(jī)組的穩(wěn)定運(yùn)行和過流部件的使用壽命。因此，對(duì)水輪機(jī)運(yùn)行時(shí)開展水力量測(cè)，獲取其壓力脈動(dòng)并進(jìn)行監(jiān)測(cè)分析，研究其變化規(guī)律，對(duì)于機(jī)組的運(yùn)行維護(hù)具有十分重要的意義[1]。

針對(duì)水力測(cè)量信號(hào)非平穩(wěn)性的特點(diǎn)，結(jié)合數(shù)值分析方法，本文提出了基于經(jīng)驗(yàn)?zāi)B(tài)分解（Empirical Mode Decomposition，EMD）的手段，將復(fù)雜水力測(cè)量信號(hào)分解成若干個(gè)本征模式分量函數(shù)（Intrinsic Mode Function，IMF），這些模態(tài)函數(shù)能很好地反映信號(hào)在任何時(shí)間局部的頻率特性[2]?；贓MD方法的水力測(cè)量信號(hào)分析，能夠計(jì)算出水力信號(hào)的壓力均值和脈動(dòng)值，并對(duì)此進(jìn)行Hilbert-Huang變換（簡稱HHT），獲取信號(hào)的內(nèi)在特性，對(duì)機(jī)組調(diào)試和運(yùn)行具有一定的指導(dǎo)意義。

1 經(jīng)驗(yàn)?zāi)B(tài)分析原理介紹

1.1 經(jīng)驗(yàn)?zāi)B(tài)分解條件

EMD方法將復(fù)雜信號(hào)中不同時(shí)間尺度或趨勢(shì)逐級(jí)分解為有限個(gè)IMF，IMF函數(shù)應(yīng)滿足以下兩個(gè)條件：

（1）在整個(gè)數(shù)據(jù)集中極值點(diǎn)的個(gè)數(shù)與穿越零點(diǎn)的個(gè)數(shù)相等或只相差一個(gè)。

（2）在任一點(diǎn)由極大值點(diǎn)形成的包絡(luò)和所有極小值點(diǎn)形成的包絡(luò)的均值等于零。

IMF表征了數(shù)據(jù)內(nèi)在振動(dòng)模式。由該定義可知由過零點(diǎn)所定義的固有模態(tài)函數(shù)的每一個(gè)振動(dòng)周期，只有一個(gè)振動(dòng)模式，沒有其他復(fù)雜的奇波。IMF沒有約束為是一個(gè)窄帶信號(hào)，并且可以是頻率和振幅的調(diào)制，還可以是非穩(wěn)態(tài)的；單由頻率或者單由振幅調(diào)試的信號(hào)也可成為固有模態(tài)函數(shù)。

EMD分解方法基于以下假設(shè)條件：

（1）數(shù)據(jù)至少有兩個(gè)極值點(diǎn)，一個(gè)最大值和一個(gè)最小值。

（2）數(shù)據(jù)的局部時(shí)域特性是由極值點(diǎn)間的時(shí)間尺度唯一確定的。

（3）如果數(shù)據(jù)沒有極值點(diǎn)而單有拐點(diǎn)，則可以通過數(shù)據(jù)微分一次或多次求得極值，然后再通過積分來獲得分解結(jié)果。

1.2 經(jīng)驗(yàn)?zāi)B(tài)分解步驟

這種方法的本質(zhì)是通過數(shù)據(jù)的特征時(shí)間尺度來獲得固有波動(dòng)模式，然后分解數(shù)據(jù)。這種分解過程也被稱為“篩過程”，具體步驟如下：

設(shè)壓力波動(dòng)信號(hào)x（t）為待分解的原始信號(hào)，首先找出x（t）中所有極大值點(diǎn)，通過3次樣條函數(shù)擬合出極大值包絡(luò)線e+（t）；同理找出x（t）的極小值點(diǎn)，以同樣的方法擬合出極小值包絡(luò)線e-（t）。將上下包絡(luò)線間的均值定義為m1（t），結(jié)果見圖1。

圖1 信號(hào)上、下包絡(luò)線及均值

則原始信號(hào)x（t）去除均值后的第一分量為：

通常經(jīng)過第一次篩分所得分量h1（t）并不是IMF函數(shù)，需經(jīng)過多次篩分。將h1（t）作為新的待篩分樣板，經(jīng)過k次篩分直到h1k（t）滿足IMF的定義，則原始信號(hào)x（t）的一階段IMF分量為：

c1（t）中包含了壓力信號(hào)最短周期的分量，將c1（t）從x（t）中分離出，則得到一個(gè)去掉高頻分量的新信號(hào)r1（t）：

因?yàn)閞1（t）中仍然包含較長周期的分量，因此將r1（t）作為新的帶篩分樣本，經(jīng)相同的篩分過程獲得第二個(gè)IMF分量c2（t），繼續(xù)EMD篩分過程直到余量變得很小或者為一個(gè)單調(diào)函數(shù)rn（t），產(chǎn)生的余量為：

式中：rn（t）為趨勢(shì)項(xiàng)，代表信號(hào)的平均趨勢(shì)或均值。

1.3 經(jīng)驗(yàn)?zāi)B(tài)分解篩選判據(jù)

在實(shí)際情況中，上下包絡(luò)線的均值無法為零，法國學(xué)者Rilling[3]提出IMF的判據(jù)如下：

式中：e+（t）和e-（t）分別是極大值和極小值包絡(luò)線。設(shè)定3個(gè)閾值θ1、θ2、和α，相應(yīng)地滿足以下兩個(gè)條件，即可認(rèn)為分解得到的IMF滿足要求：

（1）σ（t）中小于θ1的比例不低于1-α，即：

式中：T是信號(hào)的總時(shí)間長度，#A表示在集合A中元素的個(gè)數(shù)。

（2）對(duì)于整個(gè)時(shí)間范圍內(nèi)的時(shí)刻t，要滿足：

以上兩個(gè)條件θ1、θ2和α的值默認(rèn)為θ1=0.05，θ2=0.5，α=0.95。

2 基于連續(xù)均方誤差法的信號(hào)重構(gòu)

信號(hào)經(jīng)EMD分解后得到了有限個(gè)頻率從高到低的IMF分量，其中階數(shù)小的IMF分量對(duì)應(yīng)于信號(hào)的高頻成分，一般包含的是信號(hào)尖銳部分或噪聲；階數(shù)大的IMF對(duì)應(yīng)于信號(hào)的低頻成分，一般認(rèn)為其中的噪聲影響小，能夠表征信號(hào)的本質(zhì)特性?；诖藢毫?shí)測(cè)信號(hào)進(jìn)行EMD分解為若干個(gè)IMF分量后，一定存在某個(gè)IMFk分量，使得對(duì)應(yīng)該分量之后的IMFk+1等信號(hào)為主導(dǎo)模態(tài)，而其前述k個(gè)IMF信號(hào)為脈動(dòng)模態(tài)，因此基于EDM的水力測(cè)量信號(hào)分解后的信號(hào)重構(gòu)的關(guān)鍵點(diǎn)就在于尋找到這個(gè)IMFk。

Does medical device product become the new concept or upgrade? 11 8

基于Boudraa[4]提出的基于連續(xù)均方誤差法則可以有效地實(shí)現(xiàn)信號(hào)重構(gòu)，即：

該方法首先找到上式對(duì)應(yīng)的脈動(dòng)能量突變的第k個(gè)IMF分量，然后用從第k+1個(gè)開始的所有IMF對(duì)信號(hào)進(jìn)行重建得到測(cè)量壓力的均值：

其中，分界點(diǎn)k定義為：

因此，經(jīng)過EMD分解后可以得到壓力測(cè)量的均值，然后通過實(shí)測(cè)信號(hào)與均值信號(hào)的差值得到當(dāng)前測(cè)量的壓力脈動(dòng)信號(hào)。

3 現(xiàn)場(chǎng)實(shí)測(cè)信號(hào)分析

在現(xiàn)場(chǎng)試驗(yàn)時(shí)一般測(cè)量的水力信號(hào)包括鋼管壓力、蝸殼壓力、轉(zhuǎn)輪與導(dǎo)葉間壓力、轉(zhuǎn)輪與泄流環(huán)間壓力、尾水管進(jìn)口壓力、尾水管出口壓力。圖2為現(xiàn)場(chǎng)機(jī)組甩100%額定負(fù)荷實(shí)測(cè)曲線。

對(duì)實(shí)測(cè)壓力信號(hào)進(jìn)行EMD分解，再實(shí)施信號(hào)重構(gòu)得到壓力脈動(dòng)信號(hào)，在此基礎(chǔ)利用HHT分解，獲取壓力脈動(dòng)信號(hào)的內(nèi)在特性?，F(xiàn)以蝸殼壓力為例進(jìn)行信號(hào)分析。

3.1 蝸殼壓力信號(hào)分解

第一步：對(duì)機(jī)組甩100%負(fù)荷時(shí)蝸殼壓力進(jìn)行EMD分解，由于數(shù)據(jù)采集精度教高為200Hz，因此進(jìn)行EMD分解時(shí)得到的IMF分量較多有14個(gè)，顯示各IMF分量如圖3所示。

第二步：利用連續(xù)均方誤差法對(duì)信號(hào)進(jìn)行重構(gòu)。即找到使得分量能量突變的IMF階數(shù)。

根據(jù)實(shí)測(cè)結(jié)果k=6；即imf1～imf6為壓力脈動(dòng)的主要成分，而imf7～imf14為均值壓力的主要成分。

第三步：根據(jù)上一步的結(jié)論對(duì)imf7～imf14進(jìn)行信號(hào)重構(gòu)，得到蝸殼均值壓力，如圖4所示。

第四步：通過實(shí)測(cè)信號(hào)與均值信號(hào)的差值得到當(dāng)前測(cè)量的壓力脈動(dòng)信號(hào)，如圖5所示。

3.2 基于FFT的壓力信號(hào)分析

對(duì)機(jī)組甩100%負(fù)荷時(shí)蝸殼壓力進(jìn)行EMD分解得到各IMF分量，對(duì)各IMF分量進(jìn)行FFT分析得到其頻譜如圖6所示。

統(tǒng)計(jì)各階IMF的頻率范圍和最大幅值如表1所示。

圖2 機(jī)組甩100%額定負(fù)荷現(xiàn)場(chǎng)錄波曲線

圖3 蝸殼壓力EMD分解各IMF分量

圖4 蝸殼均值壓力重構(gòu)

圖5 蝸殼壓力脈動(dòng)

圖6 各IMF分量FFT變換圖

表1 蝸殼壓力各IMF分量頻率范圍及幅值表

由表1可知，蝸殼測(cè)點(diǎn)壓力的最大能量出現(xiàn)在45～55Hz，與其他頻率段相比，該頻率段的幅值陡然增大，并向低頻逐步減低。

對(duì)機(jī)組甩100%負(fù)荷時(shí)蝸殼壓力脈動(dòng)信號(hào)進(jìn)行FFT分析得到其頻譜如圖7所示。

經(jīng)FFT變換得到機(jī)組甩100%負(fù)荷時(shí)蝸殼壓力脈動(dòng)信號(hào)的主頻為90Hz和4.5Hz，因此在現(xiàn)場(chǎng)試驗(yàn)時(shí)要特別注意該頻段信號(hào)的脈動(dòng)信息。

圖7 甩100%負(fù)荷蝸殼壓力脈動(dòng)FFT變換圖

4 結(jié)束語

水力測(cè)量信號(hào)經(jīng)EMD分解后得到各IMF分量是從實(shí)測(cè)信號(hào)數(shù)據(jù)中分離出來的，相較小波分析更簡單，無需事先確認(rèn)分解次和分解基[5]；通過客觀、自適應(yīng)的分解過程，使得各IMF分量能夠反映原始信號(hào)的頻帶信息。

現(xiàn)場(chǎng)實(shí)測(cè)壓力信號(hào)在EMD分解的基礎(chǔ)上采用基于連續(xù)均方誤差方法進(jìn)行信號(hào)重構(gòu)，創(chuàng)造性地得到了均值壓力和脈動(dòng)壓力，為解析現(xiàn)場(chǎng)實(shí)測(cè)信號(hào)提供了技術(shù)手段，為分析機(jī)組各部位的壓力和脈動(dòng)信息提供了數(shù)據(jù)基礎(chǔ)，為控制現(xiàn)場(chǎng)試驗(yàn)安全提供了有力保障。

[1]張偉，吳玉林，陳乃祥，等.蓄能機(jī)組壓力脈動(dòng)研究[J].水力發(fā)電學(xué)報(bào)，2001，（72）：72-78.

[2]Patrick Flandrin，Gabriel Rilling，Paulo Gon?alvés.Empirical Mode Decomposition as a Filter Bank [J].IEEE Signal Process Letters，2004，11（2）：112-114.

[3]G.Rilling，P.Flandrin，P.Gon?alvés.On empirical mode decomposition and its algorithms. IEEE-IURASIP workshop on nonlinear signal and image processing，NSIP-03，2003.

[4]Abdel Ouahab Boudraa，Jean-christophe Cexus. EMDBased Signal Filtering.IEEE Transaction on Instrumentation and Measurement. 2007，56（6）：2196-2202.

[5]王婷. EMD算法研究及其在信號(hào)去噪中的應(yīng)用[D].哈爾濱：哈爾濱工業(yè)大學(xué)，2010.

周 攀（1983—），男，工程師，主要從事水電站水機(jī)方向研究和工程調(diào)試工作。E-mail：zp1983101@163.com.cn

Analyze on Hydraulic Metering Signals Based on EMD Method

ZHOU Pan，ZHOU Dongyue，DENG Lei
（State Grid Xinyuan Company Technology Center，Beijing 100073，China）

The pressure fluctuation has the character of nonstationary and time invariant especially in transient process，while using the EMD method can perfectly solve the difficulty of determining the pressure fluctuation in the field testing.It effectively realizes rebuilding of the pressure fluctuation signal basing on root mean squared error. Basing on field load rejection test data，this paper achieves the extracting of the pressure fluctuation by EMD method，provides a magnificent means to analyze the pressure and fluctuation of the unit in each parts，and proves the validity of this method in measuring and analyzing the actual signal.

transient process；empirical mode decomposition；intrinsic mode function；pressure fluctuation；hydraulic measurement signals