陳道云，孫守光，李 強(qiáng)

（北京交通大學(xué)　機(jī)械與電子控制工程學(xué)院，北京100044）

基于經(jīng)驗(yàn)?zāi)B(tài)分解及時(shí)域積分的軌道高低不平順獲得方法

陳道云，孫守光，李 強(qiáng)

（北京交通大學(xué)機(jī)械與電子控制工程學(xué)院，北京100044）

分析了信號(hào)中的趨勢(shì)項(xiàng)產(chǎn)生原因，并推導(dǎo)了時(shí)域積分過(guò)程中趨勢(shì)項(xiàng)的產(chǎn)生過(guò)程；利用經(jīng)驗(yàn)?zāi)B(tài)分解的方法提取并去除了實(shí)測(cè)列車軸箱垂向加速度信號(hào)的趨勢(shì)項(xiàng)，然后對(duì)去趨勢(shì)項(xiàng)的加速度信號(hào)在時(shí)域內(nèi)連續(xù)兩次積分得到軌道高低不平順數(shù)據(jù)；分別用經(jīng)驗(yàn)?zāi)B(tài)分解法和最小二乘法去除了積分過(guò)程中產(chǎn)生的趨勢(shì)項(xiàng)，并將兩種方法的處理結(jié)果進(jìn)行了對(duì)比分析。分析結(jié)果表明，利用最小二乘法去除信號(hào)及積分過(guò)程中產(chǎn)生的趨勢(shì)項(xiàng)存在較大的計(jì)算誤差，而利用經(jīng)驗(yàn)?zāi)B(tài)分解的方法則可以有效去除信號(hào)及積分過(guò)程產(chǎn)生的趨勢(shì)項(xiàng)，從而得到軌道高低不平順數(shù)據(jù)。

加速度；經(jīng)驗(yàn)?zāi)B(tài)分解；最小二乘法原理；趨勢(shì)項(xiàng)；時(shí)域積分；軌道高低不平順

軌道不平順對(duì)機(jī)車車輛系統(tǒng)是一種外部激擾，是產(chǎn)生機(jī)車車輛系統(tǒng)振動(dòng)的主要根源。按其對(duì)車輛激擾作用的方向可分為垂向軌道不平順和橫向軌道不平順。左右軌垂向、橫向不平順又可組合成常說(shuō)的高低、水平、軌向及軌距4種不平順。

目前軌道不平順的獲取主要依靠軌檢車，系統(tǒng)很復(fù)雜。為了簡(jiǎn)單地獲得軌道不平順數(shù)據(jù)，本文采用了一種基于經(jīng)驗(yàn)?zāi)B(tài)分解及時(shí)域內(nèi)軸箱加速度積分的方法，該方法易于實(shí)現(xiàn)且精度較高，比較適合對(duì)軌道進(jìn)行簡(jiǎn)易評(píng)估。

對(duì)于加速度信號(hào)的積分，之前有大量的科研工作者對(duì)其進(jìn)行了不同形式的研究，而對(duì)于趨勢(shì)項(xiàng)的消除方法，普遍采用的是最小二乘法原理［1-2］，這種方法對(duì)處理平穩(wěn)且數(shù)據(jù)量小的信號(hào)具有很大的優(yōu)勢(shì)。對(duì)于本文中研究的非平穩(wěn)隨機(jī)的加速度信號(hào)，采用這種方法得到的結(jié)果往往是不準(zhǔn)確的。為此，本文利用經(jīng)驗(yàn)?zāi)B(tài)分解技術(shù)對(duì)加速度信號(hào)及時(shí)域積分產(chǎn)生的趨勢(shì)項(xiàng)進(jìn)行了去除，得到了軌道高低不平順數(shù)據(jù)。

1　趨勢(shì)項(xiàng)

1.1信號(hào)本身的趨勢(shì)項(xiàng)

實(shí)際采集的加速度振動(dòng)信號(hào)，由于放大器隨著溫度變化產(chǎn)生的零點(diǎn)漂移、傳感器頻率范圍外低頻性能的不穩(wěn)定以及傳感器周圍的環(huán)境干擾，往往會(huì)偏離基線，甚至偏離基線的大小還會(huì)隨著時(shí)間變化。偏離基線隨時(shí)間變化的過(guò)程稱為信號(hào)的趨勢(shì)項(xiàng)，它的存在會(huì)直接影響到信號(hào)的正確性，因此必須將其去除。

1.2時(shí)域積分及其趨勢(shì)項(xiàng)

在時(shí)域內(nèi)對(duì)軸箱加速度信號(hào)兩次積分是獲得軌道高低不平順的最直接方法。數(shù)值積分的方法有很多，如梯形公式、Sim pson公式、Newton-Cotes公式、復(fù)化梯形公式、復(fù)化Sim pson公式等［3］。為了提高積分的精度且方便運(yùn)算，本文采用復(fù)化Sim pson公式對(duì)加速度信號(hào)進(jìn)行兩次積分，分別得到速度信號(hào)和位移信號(hào)：

式中

然而，直接通過(guò)時(shí)域積分的方法獲得的軌道高低不平順數(shù)據(jù)是不準(zhǔn)確的，這種偏差的原因正是由于趨勢(shì)項(xiàng)的存在。

假設(shè)加速度信號(hào)中含有微小直流分量，即

一次積分得：

二次積分得：

式中δ，η分別為ε一次、二次積分后產(chǎn)生的積分常量。

由式（3）～式（5）可以看出，在一次、二次積分后分別產(chǎn)生了一次、二次趨勢(shì)項(xiàng)，因此應(yīng)予以去除。

1.3趨勢(shì)項(xiàng)的去除

（1）基于經(jīng)驗(yàn)?zāi)B(tài)分解的趨勢(shì)項(xiàng)去除

1998年美國(guó)宇航局的N orden E H uang等［4］提出了一種新的信號(hào)處理方法，該方法主要由兩部分組成：首先采用經(jīng)驗(yàn)?zāi)B(tài)分解（E m pirical M ode Deco m position，以下簡(jiǎn)稱E M D）將原始時(shí)間信號(hào)分解成一組本征模函數(shù)（Intrinsic M ode Function，以下簡(jiǎn)稱I M F），然后再對(duì)每一個(gè)I M F進(jìn)行Hilbert變換。

本文利用其中的E M D方法提取并去除信號(hào)中的趨勢(shì)項(xiàng)。算法如下：

（1）找出信號(hào)x（t）的所有極大點(diǎn)和極小點(diǎn)，用3次樣條曲線分別擬合為原數(shù)據(jù)序列的上包絡(luò)線U（t）和下包絡(luò)線L（t），上、下包絡(luò)線的均值

將原數(shù)據(jù)序列x（t）減去m1（t）可得到一個(gè)去掉低頻的新數(shù)據(jù)序列h1（t），其式為

判斷h1（t）是否為I M F的條件有兩個(gè)：

①對(duì)于一列數(shù)據(jù)，極值點(diǎn)和過(guò)零點(diǎn)數(shù)目必須相等或至多相差一點(diǎn)；

②在任意點(diǎn)，由局部極大點(diǎn)和極小點(diǎn)構(gòu)成的兩條包絡(luò)線的平均值為零。

如果h1（t）不滿足上述條件，那么將h1（t）看成，則

重復(fù)以上過(guò)程m次，直到所得的h1m（t）滿足I M F所必需的條件，此時(shí)h1m（t）就是第一個(gè)I M F，記為I M F1（t），它表示信號(hào)中的最高頻成分。

（2）用x（t）減去I M F1（t）得到一個(gè)去掉高頻成分的新數(shù)據(jù)序列r1（t），其式為

將r1（t）看作是原數(shù)據(jù)序列重復(fù)步驟（1），可得到一系列I M Fi（t）和最后一個(gè)不可分解的序列rn（t）。由此，原數(shù)據(jù)序列x（t）可表示為一組I M F分量和一個(gè)殘余分量的和：

E M D分解的收斂特性表明，分解得到的殘余信號(hào)分量rn（t）為單調(diào)函數(shù)，其中包含了測(cè)試信號(hào)中頻率最低的成分，其周期大于采樣信號(hào)的長(zhǎng)度，所以rn（t）即為測(cè)試信號(hào)中包含的趨勢(shì)項(xiàng)［5］。

（2）基于最小二乘法原理的趨勢(shì)項(xiàng)去除

信號(hào)測(cè)試的采樣頻率記為fs，總采樣時(shí)間記為T，則采樣時(shí)間間隔Δt=1/fs，實(shí)測(cè)的加速度信號(hào)記為a（ti），其中：

構(gòu)造一個(gè)多項(xiàng)式函數(shù)

其中Pk為多項(xiàng)式系數(shù)；φ為所有次數(shù)不超過(guò)n（N≤T/ Δt-1）的多項(xiàng)式構(gòu)成的函數(shù)類，使得函數(shù)λn（ti）與離散數(shù)據(jù)a（ti）的誤差平方和最小，即：

若想確定擬合趨勢(shì)項(xiàng)，只需確定某一組Pk的值，使得E最小。由多元函數(shù)求極值的必要條件：

其中，j=0，1，2，…，n因而可以得到：

用矩陣表示為

可以證明，方程組的系數(shù)矩陣是一個(gè)對(duì)稱正定矩陣，故存在惟一解。

2　加速度信號(hào)測(cè)試及去除趨勢(shì)項(xiàng)處理

2.1測(cè)試情況

2014年11月12日對(duì)尚未開(kāi)通的滬昆高鐵線路上運(yùn)行的CRH3 A型動(dòng)車組的軸箱垂向加速度等參數(shù)進(jìn)行了線路實(shí)測(cè)，測(cè)試區(qū)間為南昌西至玉山南。

圖1　CRH3 A型動(dòng)車組

測(cè)試選用IC壓電式單向加速度傳感器，其靈敏度為25.1 m V/g，頻率范圍為0.7～11 000 Hz，量程為100g，布置在1-2L號(hào)軸箱與一系減振器之間。采用德國(guó)IM C數(shù)據(jù)采集系統(tǒng)對(duì)加速度信號(hào)進(jìn)行采集，如圖2所示。

2.2加速度信號(hào)去除趨勢(shì)項(xiàng)

在車輛運(yùn)行平穩(wěn)的區(qū)段，截取1 s時(shí)長(zhǎng)的加速度測(cè)試數(shù)據(jù)，觀察加速度值隨時(shí)間的變化情況，如圖3所示。

圖2　加速度傳感器安裝位置及數(shù)據(jù)采集系統(tǒng)

圖3　含有趨勢(shì)項(xiàng)的加速度信號(hào)

現(xiàn)有的消除趨勢(shì)項(xiàng)的方法大多需要預(yù)先假定信號(hào)中趨勢(shì)項(xiàng)的類型，這就需要對(duì)測(cè)試信號(hào)中包含趨勢(shì)項(xiàng)的特征具有一定的預(yù)先驗(yàn)證知識(shí)，因而不適用于處理隨機(jī)信號(hào)。

本文利用E M D算法，通過(guò)M atlab編程對(duì)上述加速度信號(hào)進(jìn)行分解，所得的I M F及其對(duì)應(yīng)的功率譜（PSD）分析如圖4所示。

對(duì)圖4進(jìn)行分析可以發(fā)現(xiàn)：

（1）軸箱加速度信號(hào)被成功分解為I M F1～I(xiàn) M F10及Trend共11個(gè)分量。I M F1頻率最高，波長(zhǎng)最短，在隨后的分量中，頻率逐漸變低，波長(zhǎng)變長(zhǎng)，其中Trend的頻率最低。各個(gè)I M F分別以不同的分辨率體現(xiàn)了信號(hào)中該頻段范圍內(nèi)的振動(dòng)模態(tài)。

（2）軸箱加速度信號(hào)的頻譜比較豐富，大部分信號(hào)分量集中在1 000 H z以下。I M F1～I(xiàn) M F6屬于測(cè)試信號(hào)的優(yōu)勢(shì)頻段，體現(xiàn)了軸箱垂向振動(dòng)的主要時(shí)頻特征。I M F7～I(xiàn) M F10屬于測(cè)試信號(hào)中包含的低頻成分。余量Trend頻率很低且所占能量很大，其為單調(diào)函數(shù)，周期大于采樣信號(hào)的長(zhǎng)度，因此屬于低頻趨勢(shì)項(xiàng)，須將其去除。

3　軌道高低不平順數(shù)據(jù)的獲取

在時(shí)域內(nèi)，對(duì)去除趨勢(shì)項(xiàng)的軸箱加速度信號(hào)連續(xù)兩次積分便得到軌道的高低不平順，如圖5所示。

可見(jiàn)時(shí)域內(nèi)積分獲得的位移激擾含有明顯的趨勢(shì)項(xiàng)，這是因?yàn)樵谌コ思铀俣刃盘?hào)的趨勢(shì)項(xiàng)之后，新的加速度信號(hào)中不可避免地殘留有微小的直流分量，在積分過(guò)程中，這種微小的直流分量被逐漸放大，以至于偏離正確值，因而需對(duì)其進(jìn)行去除趨勢(shì)項(xiàng)的處理?，F(xiàn)分別利用E M D方法和最小二乘法，編寫M atlab程序?qū)Ψe分得到的位移信號(hào)去除趨勢(shì)項(xiàng)。

3.1E M D方法去除積分趨勢(shì)項(xiàng)

利用E M D原理編寫M atlab程序?qū)ι鲜鑫灰菩盘?hào)進(jìn)行分解，所得的I M F 分量及其對(duì)應(yīng)的功率譜（PSD）分析如圖6所示。

對(duì)圖6進(jìn)行分析，發(fā)現(xiàn)經(jīng)E M D分解得到的4個(gè)分量中，Trend的頻率最低且所占能量較大，其為單調(diào)函數(shù)，周期大于采樣信號(hào)的長(zhǎng)度，因此為趨勢(shì)項(xiàng)，須將其去除。去除趨勢(shì)項(xiàng)后的位移信號(hào)如圖7所示。

3.2最小二乘法去除積分趨勢(shì)項(xiàng)

利用前面提到的最小二乘法也可以去除由時(shí)域積分產(chǎn)生的趨勢(shì)項(xiàng)，去除趨勢(shì)項(xiàng)后的位移數(shù)據(jù)如圖8所示。

3.3兩種去除趨勢(shì)項(xiàng)方法的對(duì)比

首先將E M D得到的趨勢(shì)項(xiàng)和最小二乘法擬合的二次趨勢(shì)項(xiàng)進(jìn)行對(duì)比，如圖9所示。

可以看出，兩條擬合趨勢(shì)曲線的大致走向是一致的，相比最小二乘法擬合趨勢(shì)項(xiàng)，E M D分解得到的趨勢(shì)項(xiàng)波動(dòng)更明顯一些，這與E M D分解的本質(zhì)特征有關(guān)，因?yàn)檫@種分解完全依據(jù)數(shù)據(jù)自身的時(shí)間尺度特征來(lái)進(jìn)行信號(hào)分解，具有自適應(yīng)特征，因而提取趨勢(shì)項(xiàng)的精度非常高，而利用最小二乘法擬合的趨勢(shì)項(xiàng)則是一條標(biāo)準(zhǔn)的二次函數(shù)曲線，與E M D相比，精度偏低。

然后將兩種方法去除趨勢(shì)項(xiàng)后得到的位移數(shù)據(jù)進(jìn)行對(duì)比，如圖10所示。

觀察圖10可以發(fā)現(xiàn)，通過(guò)最小二乘法去除趨勢(shì)項(xiàng)得到的位移信號(hào)波動(dòng)幅度要明顯大于通過(guò)E M D方法去除趨勢(shì)項(xiàng)得到的位移信號(hào)，但是二者的波形走勢(shì)卻基本相同。導(dǎo)致前者波動(dòng)幅度大的原因是利用最小二乘法原理去除趨勢(shì)項(xiàng)本身存在一定缺陷。

圖4　I M F及相應(yīng)功率譜

圖5　含趨勢(shì)項(xiàng)的位移信號(hào)

3.4最小二乘法去除積分趨勢(shì)項(xiàng)的缺陷

在時(shí)域內(nèi)對(duì)積分后的數(shù)據(jù)去除趨勢(shì)項(xiàng)時(shí)，當(dāng)截取起點(diǎn)相同而長(zhǎng)度不同的幾組數(shù)據(jù)時(shí)，會(huì)導(dǎo)致在數(shù)據(jù)相同部分得到的位移修正曲線隨著截取數(shù)據(jù)長(zhǎng)度的增加而逐漸發(fā)生偏移。這是因?yàn)楫?dāng)截取的數(shù)據(jù)長(zhǎng)度不同時(shí)，利用最小二乘法得到的擬合多項(xiàng)式系數(shù)便會(huì)不同，從而導(dǎo)致擬合多項(xiàng)式的值也發(fā)生變化。

圖6　I M F分量及相應(yīng)功率譜

圖7　EMD去除趨勢(shì)項(xiàng)得到的位移

圖8　最小二乘法原理去除趨勢(shì)項(xiàng)得到的位移

圖9　兩種方法得到的趨勢(shì)項(xiàng)對(duì)比

圖10　兩種方法去除趨勢(shì)項(xiàng)后得到的位移對(duì)比

例如若要去除一階趨勢(shì)項(xiàng)，則根據(jù)式（15）得擬合多項(xiàng)式系數(shù)

可見(jiàn)，當(dāng)截取的數(shù)據(jù)點(diǎn)數(shù)不同時(shí)，T/Δt的值便會(huì)發(fā)生變化，因而擬合的一次式系數(shù)P0和P1也不同，在起點(diǎn)相同終點(diǎn)不同的各組積分區(qū)間內(nèi)，數(shù)據(jù)相同部分得到的位移修正曲線便會(huì)不同。這也是利用最小二乘法去除趨勢(shì)項(xiàng)的一個(gè)很大的缺陷，這種缺陷正是算法本身的缺陷造成的。

為了用試驗(yàn)數(shù)據(jù)驗(yàn)證這種積分誤差，分別取兩個(gè)時(shí)間區(qū)間［0，0.2］，［0，1］對(duì)加速度數(shù)據(jù)兩次積分并去除二階趨勢(shì)項(xiàng)，得到兩條位移曲線，為了便于觀察，取相同時(shí)間區(qū)間［0，0.2］作圖如圖11。

圖11　不同長(zhǎng)度積分區(qū)間對(duì)相同時(shí)間段積分結(jié)果影響

從圖中可以看到，隨著積分區(qū)間長(zhǎng)度的增大，在相同區(qū)間段［0，0.2］的位移值發(fā)生偏移，這正是因?yàn)樵黾恿藬?shù)據(jù)點(diǎn)數(shù)而導(dǎo)致對(duì)前面的擬合帶來(lái)了影響。為了消除這種影響，可以采用動(dòng)態(tài)遞推法［6］，步驟如下：

采用某段長(zhǎng)度的加速度信號(hào)進(jìn)行積分和擬合修正，將修正位移結(jié)果第一個(gè)點(diǎn)的值作為最終位移的第一個(gè)值，將該段長(zhǎng)度的數(shù)據(jù)點(diǎn)向后遞推一個(gè)采樣點(diǎn)進(jìn)入下一次循環(huán)處理，各次循環(huán)處理得到的修正位移的第一個(gè)值組成最終位移值。

雖然動(dòng)態(tài)遞推法可以有效消除后面數(shù)據(jù)對(duì)前面擬合影響帶來(lái)的誤差，但是該方法的關(guān)鍵在于找到一段積分結(jié)果準(zhǔn)確的數(shù)據(jù)區(qū)間，這樣才能利用動(dòng)態(tài)遞推的原理消除因不同數(shù)據(jù)區(qū)間長(zhǎng)度而引起的積分結(jié)果偏移。對(duì)于軸箱加速度信號(hào)，如何找到這一個(gè)積分結(jié)果準(zhǔn)確的基準(zhǔn)積分區(qū)間還有待進(jìn)一步研究。

4　結(jié) 論

（1）趨勢(shì)項(xiàng)普遍存在于采集得到的信號(hào)以及時(shí)域積分過(guò)程中，運(yùn)用E M D法可有效地消除實(shí)測(cè)信號(hào)以及時(shí)域積分產(chǎn)生的趨勢(shì)項(xiàng)，使通過(guò)軸箱加速度信號(hào)積分得到的軌道高低不平順數(shù)據(jù)更加準(zhǔn)確。

（2）利用最小二乘法在時(shí)域內(nèi)對(duì)積分后的數(shù)據(jù)去除趨勢(shì)項(xiàng)時(shí)，當(dāng)截取起點(diǎn)相同而長(zhǎng)度不同的幾組數(shù)據(jù)時(shí)，會(huì)導(dǎo)致在相同時(shí)間區(qū)間部分得到的位移修正曲線隨著截取數(shù)據(jù)長(zhǎng)度的增加而逐漸發(fā)生偏移。雖然采用動(dòng)態(tài)遞推法可以消除這種偏移，但是如何選擇合理的基準(zhǔn)積分區(qū)間有待研究。

［1］ 陳為真，汪秉文，胡曉婭.基于時(shí)域積分的加速度信號(hào)處理［J］.華中科技大學(xué)學(xué)報(bào)，2010，38（1）：2-4.

［2］ 李東文，熊曉燕，李 博.振動(dòng)加速度信號(hào)處理探討［J］.機(jī)電工程技術(shù)，2008，37（9）：50-52.

［3］ 王兵團(tuán)，張作泉，趙平福.數(shù)值分析簡(jiǎn)明教程［M］.北京：清華大學(xué)出版社，北京交通大學(xué)出版社，2012.

［4］ H uang N E，Shen Z，Long S R，et al.The e-mpirical mode decomposition and the Hi lbe-rt spectrum for nonl inear and non-stationa-ry time series analysis［J］.Proceedings of the Royal Society of London，Series A，1998，454：903-995.

［5］ Rilling G，F(xiàn)landrin P.O n em pirical m ode deco m position and its algorith ms［C］∥Grado.IE E E-E U R A SIPW orkshop on Nonlinear Signal and Image Processing.N SIP-03，June 2003.

［6］ 馬建倉(cāng)，劉 琦，程存虎，等.依據(jù)振動(dòng)加速度信號(hào)繪制航空發(fā)動(dòng)機(jī)轉(zhuǎn)子軸心軌跡的方法［J］.航空計(jì)算技術(shù)，2009，39（4）：10-14.

An Acquired M ethod of Track Vertical Profile Irregularity Based on E M D and Time Domain Integral

C H E N D aoyun，S U N Shouguang，LI Qiang
（School of M echanical，Electronic and Control Engineering，Beijing Jiaotong U niversity，Beijing 100044，China）

This paper analyzed the cause of trend item in signals and derived the generation process of trend item in time do main integral. Then it extracted the trend item of actual measured vertical acceleration signals of the axle box and rem oved it.After that，the new acceleration signals were treated by twice-integration in time do main，w hich finally turned into track vertical profile irregularity data.This paper rem oved the trend item produced in the process of time do main integral by using E m pirical M ode Deco m position（E M D）and least square principle and contrasted the treatment results of two methods.The results showed that big calculation errors will happen by using least square principle to rem ove trend item in signals and process of time do main integral.H owever，E M D can rem ove trend item effectively and guarantee the accuracy of track vertical profile irregularity data.

acceleration；E M D；least square principle；trend term；time do main integral；track vertical profile irregularity

U26011+1

A

10.3969/j.issn.1008-7842.2015.05.03

1008-7842（2015）05-0009-06

陳道云（1988—）男，博士研究生（2015-04-27）