張海強(qiáng)，王 楠，劉 松，李柄奎

（1.廣州地鐵集團(tuán)有限公司 運(yùn)營(yíng)事業(yè)總部，廣州 510335；2. 南京理工大學(xué) 自動(dòng)化學(xué)院，南京 210094；3.中國(guó)移動(dòng)江蘇公司 揚(yáng)州分公司，揚(yáng)州 225000）

隨著城市軌道交通技術(shù)的快速發(fā)展，地鐵運(yùn)行速度逐步提升，鋼軌表面受到的沖擊愈加嚴(yán)重，導(dǎo)致軌道表面容易產(chǎn)生波磨磨耗。當(dāng)車(chē)輛行駛至波磨路段時(shí)，會(huì)引起強(qiáng)烈的震動(dòng)和沖擊，導(dǎo)致列車(chē)部件受損和乘客舒適度下降，嚴(yán)重時(shí)可能會(huì)發(fā)生列車(chē)脫軌事故。弦測(cè)法[1-2]和慣性基準(zhǔn)法[3-4]是目前軌道波磨檢測(cè)的主要方法，得到了深入研究。孫銳等人[5]以實(shí)測(cè)軌道波磨數(shù)據(jù)為基礎(chǔ)采用ABAQUS 仿真軟件，建立了輪軌關(guān)系有限元模型，研究了不同參數(shù)對(duì)車(chē)輛駛過(guò)波磨時(shí)的動(dòng)態(tài)影響，并論證了不同波磨之間存在共生關(guān)系。張厚貴等人[6]通過(guò)大量的試驗(yàn)和現(xiàn)場(chǎng)調(diào)研，對(duì)波磨形成機(jī)理進(jìn)行了探討，結(jié)果表明，當(dāng)鋼軌模態(tài)頻率因軌下結(jié)構(gòu)與輪對(duì)反共振頻率一致時(shí)，輪軌系統(tǒng)將發(fā)生強(qiáng)烈振動(dòng)，引起軌面磨耗損傷，在反饋- 循環(huán)機(jī)制下，會(huì)生成與軌下模態(tài)頻率相對(duì)應(yīng)的某種波長(zhǎng)的波磨。Oyarzabal 等人[7]應(yīng)用鋼軌波磨萌生和擴(kuò)展（RACING）對(duì)軌道支撐結(jié)構(gòu)中14 個(gè)最為明顯的因素進(jìn)行優(yōu)化，以實(shí)現(xiàn)波磨產(chǎn)生概率的最小化。Valehrach J 等人[8]對(duì)軌道波磨的產(chǎn)生環(huán)境進(jìn)行分析，發(fā)現(xiàn)列車(chē)行駛速度是影響波磨發(fā)展的因素之一。Kurzeck 等人[9]對(duì)小半徑曲線(xiàn)軌道的波磨進(jìn)行分析，驗(yàn)證了輪軌間的摩擦振動(dòng)對(duì)軌道波磨的產(chǎn)生具有很大的相關(guān)性，并且討論了輪軌之間的摩擦自激振動(dòng)需要的條件。Zhai 等人[10]利用SIMPACK 軟件進(jìn)行轉(zhuǎn)向架加速度數(shù)據(jù)仿真，利用連續(xù)小波分析得到的特征頻率和經(jīng)驗(yàn)?zāi)B(tài)分解得到的特征模態(tài)可以對(duì)軌道波磨進(jìn)行診斷識(shí)別。黃文等人[11]研究了小波包能量熵方法，通過(guò)對(duì)波磨仿真數(shù)據(jù)的處理，能夠獲得軌道波磨位置與波長(zhǎng)，并通過(guò)現(xiàn)場(chǎng)數(shù)據(jù)驗(yàn)證了檢測(cè)方法的可行性。王乃珍等人[12]通過(guò)建立灰色區(qū)間預(yù)測(cè)模型對(duì)軌道的不平順狀態(tài)進(jìn)行預(yù)測(cè)，并基于京九線(xiàn)的實(shí)測(cè)數(shù)據(jù)進(jìn)行驗(yàn)證。

本文提出了一種基于參數(shù)優(yōu)化變分模態(tài)分解（VMD，Variational Mode Decomposition）和平滑偽維格納分布（SPWVD ，Smooth Pseudo Wigner-Ville Distribution）的軌道波磨辨識(shí)方法，對(duì)列車(chē)軸箱垂向振動(dòng)加速度信號(hào)進(jìn)行處理分析，以分解信號(hào)的包絡(luò)熵為特征量，設(shè)定故障閾值，識(shí)別軌道波磨故障，計(jì)算波磨波長(zhǎng)并確定軌道波磨發(fā)生的位置。本方法不僅消除了VMD 分解參數(shù)選取不確定的弊端，而且克服了振動(dòng)信號(hào)時(shí)頻分解過(guò)程中出現(xiàn)的模態(tài)混疊問(wèn)題，具有良好的準(zhǔn)確性，可輔助軌道波磨辨識(shí)，提高軌道維保效率。

1 基本理論

1.1 變步長(zhǎng)最小均方算法

變步長(zhǎng)最小均方（VSSLMS，VariableStep Size Least Mean Square）算法[13-14]具有計(jì)算量小、自適應(yīng)能力高等優(yōu)點(diǎn)，且其收斂速度在算法不同時(shí)期可以采用不同步長(zhǎng)來(lái)控制，因此在信號(hào)處理的各個(gè)領(lǐng)域得到了廣泛的應(yīng)用。

令d(n)表示期望輸出信號(hào)，W(n)表示濾波器的抽頭權(quán)向量，e(n)表示誤差信號(hào)，X(n)表示輸入信號(hào)向量，則它們有如下關(guān)系：

式中，μ(n)為步長(zhǎng)參數(shù)，并且步長(zhǎng)迭代關(guān)系為：

式中，μmax是獲得最大可能收斂速率的步長(zhǎng)；μmin是具有跟蹤能力的最小步長(zhǎng)；α＞0，β＜1，α為步長(zhǎng)遺傳因子，主要作用是確定算法收斂時(shí)的步長(zhǎng)值，β決定步長(zhǎng)受誤差瞬時(shí)能量的影響程度，決定著算法的開(kāi)始和跟蹤速度。

在VSSLMS 的初期階段，e(n)很大，且指數(shù)an衰減很小，因此相對(duì)的步長(zhǎng)μ(n)比較大，算法收斂比較快。在VSSLMS 的收斂末期，e(n)較小，而指數(shù)an衰減較大，所以對(duì)應(yīng)的步長(zhǎng)μ(n)較小。由于VSSLMS 能較好地達(dá)到穩(wěn)態(tài)誤差和收斂速率之間的平衡，且容易實(shí)現(xiàn)，因此本文用該方法對(duì)振動(dòng)信號(hào)進(jìn)行濾波，取得了較好的效果。

1.2 參數(shù)優(yōu)化VMD

1.2.1 VMD

Dragomiretskiy 等人[15]提出了VMD 方法，屬于一類(lèi)非遞歸分解方法，自適應(yīng)同步提取復(fù)雜信號(hào)中的調(diào)幅-調(diào)頻分量[16]。變分模式分解將信號(hào)x(t)分解為一組本征模態(tài)分量（IMF）ck(t)，通過(guò)這些本征模態(tài)分量，可以重構(gòu)信號(hào)，在頻域內(nèi)迭代更新各個(gè)本征模態(tài)分量ck(t)，將本征模態(tài)分量頻譜的重心作為中心頻率ωk。

假設(shè)本征模態(tài)分量ck(t)在頻域內(nèi)緊致聚集在中心頻率ωk周?chē)?，即具有稀疏性質(zhì)。每個(gè)本征模態(tài)分量的稀疏性質(zhì)通過(guò)它的帶寬來(lái)描述：（1）應(yīng)用Hibert變換構(gòu)造其解析信號(hào)，以便獲得非負(fù)頻率的單邊頻譜；（2）通過(guò)和頻率調(diào)整至各自中心頻率處的指數(shù)諧波相乘，將頻譜平移至基帶；（3）通過(guò)梯度的平方l2范數(shù)估計(jì)帶寬。

為了將約束變分優(yōu)化問(wèn)題轉(zhuǎn)化為無(wú)約束形式，增加懲罰項(xiàng)和Lagrance 乘子，以便快速收斂并增強(qiáng)約束，所以，將最小化的目標(biāo)函數(shù)轉(zhuǎn)換為增廣Lagrance 函數(shù)：

式中，λ(t)為L(zhǎng)agrance 乘子，α 為數(shù)據(jù)真實(shí)性約束的平衡參數(shù)表示內(nèi)積。

求解與式（3）等價(jià)的最小化問(wèn)題，應(yīng)用最優(yōu)解更新每個(gè)本征模態(tài)分量ck(t)：

將式（4）頻域內(nèi)的解取Fourier 逆變換并取實(shí)部，獲得時(shí)域里的本征模態(tài)分量。

求解和式（3）等價(jià)的最小化問(wèn)題，應(yīng)用最優(yōu)解更新每個(gè)本征模態(tài)分量ck(t)的中心頻率ωk：

變分模式分解的完整算法如下。

（2）對(duì)于k=1，2，…，K 及所有ω(ω ≥0)，更新每個(gè)本征模態(tài)分量ck(t)及其中心頻率ωk；

（3）對(duì)于所有ω ≥0，更新Lagrance 乘子；

（4）檢查收斂條件：

1.2.2 VMD參數(shù)優(yōu)化

雖然變分模態(tài)分解具有非遞歸同步分解本質(zhì)，從而有效避免了遞歸分解算法的缺點(diǎn)，如對(duì)噪聲的采樣敏感性，以及在極值插值過(guò)程中產(chǎn)生的上下包絡(luò)過(guò)擬合和欠擬合等問(wèn)題，但是對(duì)振動(dòng)信號(hào)分解結(jié)果影響比較大的參數(shù)分解模態(tài)數(shù)K 與懲罰參數(shù)α 卻不能自適應(yīng)選取，這很容易對(duì)分解結(jié)果造成較大影響。因此，提出了一種基于變分模態(tài)分解信號(hào)中IMF 分量最大峭度值參數(shù)優(yōu)化選擇方法。參數(shù)優(yōu)化VMD 算法的流程如圖1 所示，詳細(xì)步驟如下：

（1）設(shè)定參數(shù)α 的搜索區(qū)間為[200 3 000]，搜索步長(zhǎng)為10，參數(shù)K 的搜索區(qū)間設(shè)為[2 10]，搜索步長(zhǎng)為1，在（K，α）的坐標(biāo)平面上構(gòu)建了一個(gè)二維網(wǎng)格；

（2）選?。↘，α）平面網(wǎng)格點(diǎn)上的參數(shù)值，對(duì)原始信號(hào)實(shí)施變分模態(tài)分解；

（3）求出分解后的原始信號(hào)各IMF 分量中最大峭度值并保存；

（4）遍歷網(wǎng)格中所有的點(diǎn)，并完成步驟（2）、步驟（3）；

（5）從保存的峭度值中篩選出最大的，這個(gè)最大峭度值所采用的參數(shù)K 值與α 值即為最優(yōu)。

圖1 參數(shù)優(yōu)化VMD算法流程

1.3 SPWVD

維格納分布（WVD，Wigner-Ville-Distribution）是非平穩(wěn)信號(hào)分析中常見(jiàn)的二次型時(shí)頻分布，廣泛應(yīng)用在各種領(lǐng)域[17]，但WVD 可能會(huì)產(chǎn)生交叉項(xiàng)，從而輸出虛假分量。因此，為避免信號(hào)分析時(shí)的交叉項(xiàng)干擾，獲得良好的信號(hào)時(shí)頻分析效果，本文采用時(shí)域加窗的WVD 對(duì)信號(hào)的WVD 在時(shí)域和頻域做平滑處理的方法進(jìn)行時(shí)頻分析[18]。

根據(jù)WVD 的特性，兩個(gè)時(shí)域信號(hào)乘積的WVD與各自WVD 在頻域的卷積相等，因此可以表示成：

根據(jù)上式能夠得知，在時(shí)域加窗之后的結(jié)果與在頻域進(jìn)行低通平滑濾波效果等同。如果再在頻域加一個(gè)窗函數(shù)，等同于在時(shí)域做卷積。兩個(gè)時(shí)域信號(hào)卷積后的WVD 等于各自WVD 在時(shí)間軸上的卷積，相當(dāng)于對(duì)WVD 變換在時(shí)域做平滑。經(jīng)過(guò)時(shí)域與頻域加窗處理后的WVD 稱(chēng)為平滑偽維格納分布（SPWVD）。

式中，h(t)和g(t)是奇數(shù)長(zhǎng)度的窗函數(shù)，滿(mǎn)足h(0)=1 和G(0)=1，G(f)為g(t)的傅里葉變換。

2 算法流程

本文提出的基于參數(shù)優(yōu)化VMD 和SPWVD 的波磨辨識(shí)算法如圖2 所示，步驟描述如下：

圖2 軌道波磨辨識(shí)算法流程

（1）采集軸箱垂向振動(dòng)信號(hào)，進(jìn)行去除干擾及噪聲等濾波處理；

（2）將濾波后的振動(dòng)信號(hào)進(jìn)行參數(shù)優(yōu)化的變分模態(tài)分解；

（3）選取信號(hào)的本征模態(tài)分量，計(jì)算波磨的特征值包絡(luò)譜熵；

（4）根據(jù)閾值判斷是否存在波磨；

（5）若采集的垂向振動(dòng)信號(hào)，存在波磨，則對(duì)分解后的振動(dòng)信號(hào)進(jìn)行SPWVD 時(shí)頻分析，確定波磨中心頻率及波磨位置，計(jì)算波磨波長(zhǎng)。

3 仿真分析

3.1 波磨仿真

軌道波磨能夠用一個(gè)近似余弦的信號(hào)來(lái)描述，為了更好地對(duì)軌道波磨進(jìn)行仿真，本文進(jìn)行4 種不同的波磨仿真實(shí)驗(yàn)。鋼軌長(zhǎng)度為40 m，采樣頻率為1 kHz，車(chē)速為10 m/s，一側(cè)軌道添加美國(guó)五級(jí)軌道垂向不平順（AAR5），另一側(cè)軌道在特定位置添加不同波磨，詳細(xì)的參數(shù)設(shè)置如表1 所示。

表1 軌道波磨參數(shù)設(shè)置

按照表1 所述參數(shù)，將軌道不平順數(shù)據(jù)制作成tre 格式的文件，導(dǎo)入SIMPACK 多體動(dòng)力學(xué)軟件，并建立相應(yīng)的軌道激勵(lì)函數(shù)，添加至列車(chē)仿真模型，然后分別完成速度、采集頻率、時(shí)間等一系列參數(shù)設(shè)置后，進(jìn)行離線(xiàn)分析，保存軸箱振動(dòng)的仿真結(jié)果數(shù)據(jù)，得到振動(dòng)信號(hào)時(shí)域圖如圖3、圖4 所示，由圖可見(jiàn)，當(dāng)添加的軌道波磨波深相同時(shí)，隨著波長(zhǎng)的增加，軸箱垂向振動(dòng)信號(hào)幅值呈現(xiàn)下降的趨勢(shì)。

圖3 僅添加AAR5的軸箱垂向振動(dòng)時(shí)域

3.2 波磨故障辨識(shí)

由于仿真信號(hào)的振動(dòng)數(shù)據(jù)不含有噪聲信號(hào)，故不再應(yīng)用VSSLMS 方法進(jìn)行濾波，直接采用參數(shù)優(yōu)化VMD 方法進(jìn)行模態(tài)分解，根據(jù)軸箱垂向振動(dòng)信號(hào)的模態(tài)分解結(jié)果，計(jì)算分解后信號(hào)的包絡(luò)熵值，K、α參數(shù)值及包絡(luò)熵值，如表2 所示。

圖4 不同波長(zhǎng)0.09 mm波深軸箱垂向振動(dòng)時(shí)域

表2 軸箱垂向振動(dòng)參數(shù)優(yōu)化VMD的分解參數(shù)及信號(hào)包絡(luò)熵

由表2 中結(jié)果可知，僅含美國(guó)五級(jí)軌道譜的第一組包絡(luò)熵值為0.920 7，而其他添加軌道不平順的包絡(luò)熵?cái)?shù)據(jù)遠(yuǎn)大于第一組的值，因此，包絡(luò)熵可以作為軌道是否發(fā)生波磨的特征值，參考僅含美國(guó)五級(jí)軌道譜的包絡(luò)熵值，為留有一定的工程余量設(shè)置閾值為0.95，大于閾值判定為存在軌道波磨故障，反之則無(wú)軌道波磨故障。接下來(lái)為了進(jìn)一步分析軌道波磨發(fā)生的位置及波長(zhǎng)范圍，本文選取第4 組數(shù)據(jù)，用SPWVD 方法對(duì)分解后的垂向振動(dòng)信號(hào)進(jìn)行時(shí)頻分析，分析結(jié)果如圖5 所示。

圖5 200 mm波長(zhǎng)0.09 mm波深時(shí)頻

由圖5 可知，在A 處1.41 s ～1.69 s 時(shí)發(fā)生頻率集中，頻率范圍為[44.91 Hz 56.63 Hz]。時(shí)頻圖中零點(diǎn)起始位置也有頻率集中現(xiàn)象，但經(jīng)過(guò)分析可知，此位置是由美國(guó)五級(jí)軌道譜所引起，并非添加軌道波磨所致。由于采樣頻率設(shè)置為1 kHz，根據(jù)波磨波長(zhǎng)λ與振動(dòng)頻率f的關(guān)系式可知：

式中，v為列車(chē)速度。

為了準(zhǔn)確計(jì)算波磨波長(zhǎng)，f取頻率集中區(qū)的中心頻率，計(jì)算可知，圖5 中添加軌道波磨的位置在14.1 m ～16.9 m，這與在15 m ～16 m 處的添加波磨位置基本穩(wěn)合，波長(zhǎng)為199.63 mm，與仿真200 mm波磨波長(zhǎng)非常接近。

4 實(shí)例分析

為了進(jìn)一步驗(yàn)證波磨故障辨識(shí)算法的有效性，本節(jié)利用軌道波磨在線(xiàn)檢測(cè)系統(tǒng)設(shè)備，對(duì)廣州地鐵8號(hào)線(xiàn)的軌道線(xiàn)路進(jìn)行檢測(cè)試驗(yàn)，檢測(cè)列車(chē)運(yùn)行速度為72 km/h，采樣頻率是10 kHz。本文截取一段1 km長(zhǎng)度包含異常的軸箱垂向振動(dòng)信號(hào)數(shù)據(jù)進(jìn)行分析，其振動(dòng)信號(hào)波形如圖6 所示。由圖6 可知，實(shí)測(cè)軸箱垂向振動(dòng)信號(hào)中含有大量噪聲干擾，因此需要對(duì)采集的振動(dòng)信號(hào)進(jìn)行濾波處理，本文采用VSSLMS方法濾波，濾波后的波形如圖7 所示。

圖6 軸箱垂向振動(dòng)信號(hào)波形（實(shí)測(cè)值）

圖7 軸箱垂向振動(dòng)信號(hào)波形（濾波后）

由圖7 可知，經(jīng)過(guò)濾波后的軸箱振動(dòng)信號(hào)在200 m ～300 m 及390 m ～530 m 處具有比較大的幅值波動(dòng)，特別是在里程390 m ～530 m 處，不僅波幅大而且里程較長(zhǎng)，并且從圖7 還可以看到較多的振動(dòng)信號(hào)尖峰，分析可知，這是由于城軌車(chē)輛運(yùn)行到鋼軌接頭處產(chǎn)生的沖擊力變化所致。

為了計(jì)算方便，本文截取350 m ～550 m 里程處的濾波振動(dòng)信號(hào)進(jìn)行數(shù)據(jù)分析，并且對(duì)此段數(shù)據(jù)5倍降頻處理，采樣頻率降為2 kHz，如圖8 所示。

圖8 350 m～550 m處濾波后的實(shí)測(cè)軸箱垂向振動(dòng)信號(hào)

在350 m ～550 m 里程處濾波信號(hào)經(jīng)過(guò)參數(shù)優(yōu)化的變分模態(tài)分解后，計(jì)算各IMF 分量的包絡(luò)熵為13.936，而這段軌道正常軸向振動(dòng)的信號(hào)的包絡(luò)熵值為1.917，故可以判斷出該段軌道存在波磨。

為了進(jìn)一步驗(yàn)證檢測(cè)方法的有效性，在廣州地鐵工程人員的配合下，對(duì)440 m ～480 m 處軌道進(jìn)行檢查，發(fā)現(xiàn)該段線(xiàn)路發(fā)生了如圖9 所示的波磨故障，從圖中可以明顯看到鋼軌表面發(fā)生了波浪形的磨耗，測(cè)量后得到該處軌道波磨的波長(zhǎng)大致在76 mm 左右，與本方法檢測(cè)結(jié)果誤差最大僅為3.92%，從而驗(yàn)證了本文所述方法對(duì)波磨故障檢測(cè)的有效性。

圖9 現(xiàn)場(chǎng)存在軌道波磨路段

5 結(jié)束語(yǔ)

本文通過(guò)采集環(huán)境激勵(lì)下列車(chē)軸箱上振動(dòng)加速度傳感器的響應(yīng)信號(hào)，利用變步長(zhǎng)的LMS 自適應(yīng)濾波方法對(duì)采集的加速度信號(hào)過(guò)濾作降噪處理，采用參數(shù)優(yōu)化VMD 方法對(duì)信號(hào)進(jìn)行分解，以分解信號(hào)的包絡(luò)熵為故障辨識(shí)指標(biāo)，根據(jù)正常軌道的包絡(luò)熵設(shè)定閾值，識(shí)別軌道波磨，最后利用SPWVD 方法對(duì)含有波磨的振動(dòng)信號(hào)進(jìn)行時(shí)頻分析。具體結(jié)論如下：

（1）參數(shù)優(yōu)化VMD 算法不僅克服了VMD 分解參數(shù)不能自適應(yīng)選取的弊端，并且通過(guò)對(duì)比振動(dòng)信號(hào)VMD 分解圖及各IMF 分量的時(shí)頻圖，充分說(shuō)明了VMD 分解在關(guān)于抑制模態(tài)混疊方面的優(yōu)越性；

（2）根據(jù)不同故障振動(dòng)信號(hào)分解時(shí)包絡(luò)譜的差異，提出了以分解信號(hào)的包絡(luò)熵為軌道波磨辨識(shí)指標(biāo)，通過(guò)設(shè)定故障閾值，識(shí)別軌道波磨；

（3）基于SPWVD 方法含軌道波磨振動(dòng)信號(hào)進(jìn)行時(shí)頻分析，根據(jù)分析結(jié)果，可以判斷出軌道波磨位置，并且利用波磨波長(zhǎng)λ與振動(dòng)頻率f的關(guān)系，能夠計(jì)算出波磨波長(zhǎng)。