履帶車輛制動器扭振信號瞬時頻率特征提取方法研究

周鋮，羅楊，魏江，曹宏瑞，2，蘭海，張萬昊

(1.西安交通大學 機械工程學院，陜西 西安 710049;2.西安交通大學 機械制造系統工程國家重點實驗室，陜西 西安 710049;3.中國北方車輛研究所，北京 100072)

0 引言

履帶車輛因其良好的機動性能在現代軍事、農業及建筑業都發揮著關鍵性作用［1-2］。隨著我國綜合國力的提升，科學技術以及軍事力量也有極大的發展，為保證履帶車輛在運輸和戰斗過程中具有更高的機動性能，履帶車輛發動機也向高速大功率方向發展［3］，對于工況復雜的履帶車輛而言，其制動器的扭振現象愈加顯著。

發動機是履帶車輛的核心，除提供動力外，也是扭振的激勵源［4-5］。當發動機產生扭振后經過一系列的傳動裝置最終將激勵傳遞到制動器，從而使得制動器發生扭振。制動器發生扭振對其可靠性是一個潛在的極大威脅，在很大程度上扭振會造成制動器制動效果降低甚至失效，從而引發事故。因此對制動器扭振的研究在履帶車輛的研究中尤為重要，對其參數的測試不僅可以為制動器的內部尺寸設計提供參數，降低制動器的扭振幅值，并且能為履帶車輛制動器中制動片的更換提供參考，保證履帶車輛制動器在工作過程中的可靠性。

目前，測量扭振的方法還在不斷發展的過程中，對于扭振信號的獲取精度、有效的數據分析方法和扭振監測的實時性都需進一步得到提高［6-7］，尤其針對處于復雜惡劣、空間狹小結構緊湊的環境中的旋轉部件。本文就某型履帶車輛制動器的測試環境惡劣、工況復雜、干擾強等特點，展開了對履帶車輛制動器扭振瞬時頻率特征提取的研究，采用脈沖時序計數法［8-10］的扭振信號獲取方法，研究基于參數優化變分模態分解(VMD) 的扭振信號瞬時頻率特征提取方法，通過采用優化的VMD 算法實現對扭振信號的分解降噪以及重構，采用零點插值計算扭振信號兩脈沖之間的時間間隔，完成對扭振信號瞬時轉速波動信號的計算，而后采用頻譜分析，提取扭振信號的瞬時頻率特征。

1 扭振信號瞬時轉速獲取方法

1.1 扭振信號的測量原理

脈沖時序計數法扭振測量技術是利用安裝在軸上的等分齒輪結構，通過非接觸式磁電、光等傳感器產生相應的脈沖信號［11］。碼盤或等分齒輪每轉過一定角度就會產生一個脈沖，測得兩個脈沖間隔時間，就可獲得兩脈沖之間的平均角速度，進而獲得扭轉角和角速度差。當轉軸不發生扭振時，轉軸的瞬時速度等于其平均速度，傳感器經過等分齒輪時輸出的脈沖為等間隔的［12］;當轉軸發生扭振時，相當于在轉軸平均速度上疊加了轉速波動，則傳感器輸出的脈沖為不等間隔的。

假設等分齒輪齒數為Z 個，則每個齒的分度Δθ 為

或用弧度制(rad) 表示Δθ 為

根據傳感器輸出波形可采用高頻A/D 采樣法［13］或高頻計數法［14］求出轉過第i 齒時的脈沖信號的周期Ti(s)，則轉過第i 齒(i=1，2，3，…，Z) 時的瞬時轉速ni為

由此可得轉軸轉過一周的平均轉速為

則等分齒輪轉過第i 齒的扭振瞬時轉速波動為

1.2 扭振信號瞬時轉速的提取

由1.1 節計算公式可知，扭振測量的關鍵是獲取扭振信號周期和計算瞬時轉速，瞬時轉速的提取主要采用零點插值法［15-16］。

零點插值法是指通過線性插值的方式找到原始信號波形與零線的交點序列{p1，p2，…，pZ}，如圖1所示。

圖1 零點法求瞬時轉速示意圖Fig.1 Diagram of instantaneous speed by the zero point method

分析圖2 所示的信號片段，可知單個齒經過所產生的信號為pi-1至pi+1段。

圖2 線性插值零點法計算瞬時轉速Fig.2 Calculation of instantaneous speed by the zero point linear interpolation method

轉過第i 齒經歷的時間為

式中:T 為采樣卡采樣周期。

第i 齒的瞬時轉速為

2 參數優化VMD

2.1 VMD 算法基本原理

VMD 是一種非遞歸、自適應分解算法，每個本征模式分量(IMF) 分量的帶寬和中心頻率利用迭代確定［17］，由K 個IMF 組成的模態函數uk(t) (k=1，2，…，K) 的具體步驟如下:

1) 通過Hilbert 變換，得到各模態函數uk(t) 的解析信號，獲得uk(t)的單邊頻譜uk(t)，δ(t) 為Dirichlet 函數，* 為卷積符號。

2) 進行頻率混合，對各模態函數uk(t) 的解析信號混合一個預先估計中心頻率{ω1，ω2，…，ωK}表示各IMF 分量的頻率中心。將各模態的頻譜變換到基頻帶

3) 通過高斯平滑解調信號獲得每段帶寬，即L2范數梯度的平方根。構造出約束的變分模型為

式中:{uk}={ u1，u2，…，uK} 表示分解得到的K 個IMF 分量;?t為對函數求時間的偏導數;x(t) 為分解的主信號。

為求解上述約束的變分模型的最優解，引入懲罰因子α 構造增廣Lagrange 函數為

式中:λ 為Lagrange 乘子;α 為懲罰參數;〈·〉表示內積運算;λ(t) 用來保持約束條件的嚴格性。

將上述Lageange 函數從時域轉換到頻域，并進行極值求解，可獲得迭代n +1 次過程中模態分量在頻域中的表達式，即

4) 使用交替方向乘子算法(ADMM)［18］進行約束變分模型的最優解的求取，將原始信號將被分解為K 個窄帶IMF 分量。

具體算法如下:

3) 根據如式(13) 更新λn+1:

式中:τ 為正則系數。

2.2 參數優化VMD

由PSO 算法［19-21］的原理可知，PSO 能夠根據適應度函數的最值作為優化準則來尋找VMD 算法的最優參數組合。利用PSO 搜尋VMD 算法的影響參數時，需確定一個適應度函數，粒子每次更新位置時計算一次適應度值，通過對比新粒子適應度值進行更新。本文基于扭振信號的特點采用能量熵［22］作為PSO 的適應度函數搜尋VMD 最優參數組合［Kopt，αopt］。

能量熵Lp為

式中:Lp(k) 表示第k 個IMF 分量的能量熵;pk表示第k 個單分量IMF 能量的歸一化形式;E(k) 表示第k 個單分量IMF 的能量。

基于PSO 優化VMD 參數組合的主要操作步驟如下:

1) 定義粒子的位置和速度范圍。設粒子群個體的位置對應VMD 算法參數［K，α］的一個候選解;初始化所有粒子位置和速度，并在設定范圍內選取隨機值，通過計算確定目標函數。

2) 確定PSO 算法的適應度函數，在不同位置對扭振信號進行VMD 運算，并采用式(14) 和式(15) 計算各IMF 的能量熵;

3) 根據式(16) 計算最小能量熵為局部極小值min Lp(k)，即VMD 的最優參數組合［Kopt，αopt］的一個候選解;

4) 比較不同位置的適應度函數值局部極小值min Lp(k)，并對個體局部極值和群體全局極值進行更新;

5) 對粒子的速度和位置進行更新。轉至步驟2，直至達到預設迭代次數，輸出最佳適應度函數值min Lp和粒子位置［Kopt，αopt］。

3 扭振信號處理流程

由扭振信號的測量方法可知，扭振的獲取是一種間接測量法，獲取的扭振信號中包含大量的噪聲，故而必須采用一定的手段或方法進行降噪。本文采用VMD 算法對原始扭振信號進行分解降噪，并利用PSO 算法搜尋VMD 算法最優參數組合，由此提出基于PSO-VMD 的扭振信號處理方法，具體實現步驟如下:

1) 初始化PSO 算法的相關參數。PSO 算法設定的各項參數如表1 所示。

表1 PSO 各項參數Table 1 PSO parameters

種群數量m 反映了計算時需要的粒子的數目，m 過小，PSO 計算速度快但易于過早收斂，過大，尋優效率降低，通常設置為［10 -100］，Nmax為算法最大迭代次數。

種群維度D 一般用于反映優化的參數個數，本文主要優化VMD 算法的分解個數K 和懲罰因子α兩個參數，故而設置為2。

學習因子c1和c2反映了粒子間的學習能力，c1表示自我經驗、c2表示種群經驗。一般c1和c2取固定常數，并且取c1和c2相等，早期的研究中通常將取值規定為2［19-21］。

慣性權重可調節PSO 算法局部和全局搜索能力，Shi 發現，當w 在0.9～0.4 線性變化時PSO 算法優化能力更好，即采用線性遞減策略:

2) PSO 搜尋VMD 算法的參數組合。以VMD算法中的分解個數K 和懲罰因子α 目標函數，以能量熵為適應度函數進行優化，得到最優參數組合［Kopt，αopt］。

3) 將最優參數組合［Kopt，αopt］重新代入到VMD 算法中完成扭振信號的分解，得到Kopt個IMF分量;

4) 計算各IMF 與原始扭振信號的相關系數rK，并取相關系數rK最大(或較大) 對應的IMF 進行扭振信號的重構;

相關系數用來描述經VMD 后得到的各IMF 與原始扭振信號的關聯程度，通過計算各IMF 與原始信號自相關系數將其歸一化處理［23］，其表達式為

式中: N0為扭振信號y(t) 對應的原始信號個數;為第k 個IMF 的自相關系數;Ry(k) 為扭振信號y(t) 的自相關系數。rK越大表示第K 個IMF與原始信號的相關性越大，含有的轉頻信息也越豐富。

5) 將重構后的扭振信號進行零點插值，提取扭振瞬時轉速波動，最后對扭振瞬時轉速波動進行頻譜分析，得到扭振瞬時頻譜特征。

以上扭振信號處理流程如圖3 所示。

圖3 算法流程圖Fig.3 Flow chart of the algorithm

4 實驗分析

為了驗證所提方法的有效性，搭建如圖4 所示的扭振頻率可調的扭振實驗樣機進行實驗?，F場測試中采用加工有100 齒的等分齒輪，安裝在實驗樣機的合適位置(見圖4(b))，美國Lion 公司的CPL290 電容式位移傳感器采集扭振數據。通過億恒8 通道ECON MI-7008 數據采集儀記錄，采樣頻率為20 480 Hz，采樣時間為2 s，設定實驗臺轉速為600 r/min、扭振頻率為50 Hz。采集的扭振信號時域波形如圖5 所示。

圖4 扭振實驗樣機Fig.4 Torsional vibration experimental prototype

圖5 扭振信號時域圖Fig.5 Time domain diagram of torsional vibration signals

如圖5 所示為扭振試驗樣機轉速為600 r/min、扭振頻率為50 Hz 的扭振信號時域波形圖以及扭振信號的局部放大圖。從時域圖中可明顯觀察到采集到的等分齒輪上信號具有很強的周期性，經過計算每一個周期內共有100 個脈沖信號(類正弦信號)，對應等分齒輪的齒數，說明采集到的扭振信號是正確的。

圖6 為扭振信號的頻譜分析，從頻譜的低頻部分可得到扭振轉速600 r/min 的轉頻10 Hz，從高頻部分可得到等分齒輪齒的通過頻率1 000 Hz，且存在10 Hz 的邊頻帶，從而說明采集的扭振數據的正確性。但是從頻譜圖中無法得知50 Hz 的來源，說明采用簡單的頻譜分析并不能有效提取出扭振頻率。

圖6 扭振信號頻譜分析圖Fig.6 Spectrum analysis of torsional vibration signals

采用PSO 優化VMD 參數，得到如圖7 所示的適應度函數變化曲線，從圖7 中可得，PSO 迭代1 000 次，粒子群適應度函數值最后達到了收斂，完成收斂時粒子群的迭代次數為331，對應VMD 參數組合為［K，α］=［5，300］，將優化得到的最優參數代入VMD 中對扭振信號進行分解，得到5 個單分量IMF，如圖8 所示。

圖7 適應度函數曲線Fig.7 Curve of fitness function

從圖8 中可以發現經過參數優化后的VMD 能夠將扭振信號完全分離。

圖8 VMD 各IMFFig.8 Decomposition of each IMF by VMD

計算VMD 的5 個IMF 與扭振信號的相關系數rK，結果如表2 所示，選擇相關系數最大的IMF2重構扭振信號，IMF2時域圖和頻譜分析如圖9 所示。從圖9 中可以得到等分齒輪的工頻以及以轉軸轉頻10 Hz 的邊頻帶。對重構的位移信號采用式(6) 的零點插值法求取單齒經過的時間間隔，并用式(7)求解瞬時轉速，進而利用式(5) 求解得到如圖10 所示的瞬時轉速波動數據。對瞬時轉速波動信號做頻譜分析得到如圖11 所示的轉速波動頻率即為扭振頻率。

圖9 IMF2時域波形及頻譜Fig.9 IMF2 Time domain waveform and spectrum

圖10 瞬時轉速波動Fig.10 Instantaneous speed fluctuation

表2 各IMF 分量相關系數Table 2 Correlation coefficients of IMF components

從圖11 頻譜分析局部放大圖中可以清楚地看到50.02 Hz 的扭振頻率及其倍頻，同時在扭振倍頻兩側存在10 Hz 的邊頻帶，從而證明所提算法能夠實現對扭振信號頻率的提取，且在考慮頻譜分辨率的情況下可以認為扭振提取誤差應小于1%。

圖11 瞬時轉速波動的頻譜分析Fig.11 Spectrum analysis of instantaneous speed fluctuation

為了進一步驗證本算法的優越性，現對未經優化的扭振信號和經集合經驗模態分解(EEMD) 分解［24］的扭振信號進行分析，分析結果分別如圖12和圖13 所示。

從圖12(a) 可看出扭振瞬時轉速波動信號毫無規律，從圖12(b) 扭振瞬時轉速波動信號的頻譜局部放大可看出，在頻譜圖中并未找到轉軸轉頻10 Hz分量，雖然找到了50.02 Hz 成分，但是在其二倍頻處并未找到以10 Hz 為間隔的邊頻帶，并不能說明此50.05 Hz 是否就是所要尋找的扭振瞬時頻率特征。

圖12 未經優化的扭振信號分析結果Fig.12 Results of unoptimized torsional vibration signals

采用EEMD 對扭振信號進行分解，經過分解發現，前5 個IMF 基本能反映扭振信號的模態，而后求解前5 個IMF 與原始信號的相關系數，選擇相關系數最大的IMF2重構扭振信號。

從圖13(d) 扭振瞬時轉速波動信號的頻譜局部放大可看出，在頻譜圖中雖然找到了50.024 4 Hz 的扭振頻率，但是并未找到轉軸轉頻10 Hz 分量，同時也無法找到扭振頻率的倍頻及以10 Hz 為邊頻帶的成分，從而證明了所提的算法的分解效果更佳。同時又分別對經VMD 分解得到的其他IMF 進行了同樣的分析，經過分析發現均無法得到扭振頻率特征，充分證明了選擇IMF2進行分析的合理性，同時進一步驗證了所提算法的有效性。而對于各瞬時轉速波動圖中所產生的突變成分可能是由于被測轉軸質量的不平衡或傳感器不嚴格對中所導致，將在后續的工作中作進一步分析。

圖13 經EEMD 的扭振信號分析結果Fig.13 Analysis results of torsional vibration signal decomposed by EEMD

5 結論

本文提出了一種使用粒子群參數優化的VMD和零點線性插值法的扭振特征頻率提取方法。得出主要結論如下:

1) 所提的方法關鍵在于設置信號能量熵為粒子群參數優化的適應度函數，為VMD 算法提供分解個數 和懲罰因子的最優參數組合。在最優參數組合下，VMD 算法能夠準確提取原信號中轉頻信息最豐富的分量，因而使所提方法具有獨特優越性。

2) 對最優分量信號使用零點線性插值法能夠更為準確地還原軸系瞬時轉速波動信號，即扭振信號，其頻譜能夠準確體現扭振特征頻率。

3) 慣性載荷扭振實驗結果表明: 在轉速為600 r/min、扭振激勵頻率為50 Hz 的條件下，對采樣頻率為20 480 Hz 的位移信號使用所提方法能精確提取扭振瞬時頻率特征，誤差小于1%。