許永強(qiáng)，姜志宏，蔡改貧，郭進(jìn)山，熊 洋

（江西理工大學(xué) 機(jī)電工程學(xué)院，江西 贛州 341000）

基于AR模型的磨機(jī)振動(dòng)信號(hào)特征提取方法

許永強(qiáng)，姜志宏，蔡改貧，郭進(jìn)山，熊 洋

（江西理工大學(xué) 機(jī)電工程學(xué)院，江西 贛州 341000）

磨礦過(guò)程中的球磨機(jī)筒體內(nèi)部負(fù)荷（填充率、料球比）研究是選礦設(shè)備節(jié)能降耗的重要內(nèi)容。以試驗(yàn)球磨機(jī)為對(duì)象，通過(guò)采集軸承座振動(dòng)信號(hào)，采用AR模型對(duì)振動(dòng)信號(hào)進(jìn)行特征提取和進(jìn)行功率譜估計(jì)，研究了5種充填率條件下的磨機(jī)負(fù)荷參數(shù)與信號(hào)時(shí)域特征的相關(guān)性，得出隨著磨礦過(guò)程中筒體中鋼球、物料的變化，低頻段、高頻段的頻譜能量值曲線的相應(yīng)變化規(guī)律。研究表明，有量綱時(shí)域特征參數(shù)可以很好地表征特征信號(hào)與負(fù)荷狀態(tài)參數(shù)的相關(guān)性；無(wú)量綱時(shí)域特征參數(shù)可以很好地解釋振動(dòng)時(shí)域信號(hào)的波形分布以及振動(dòng)沖擊特征，根據(jù)頻譜能量分布變化規(guī)律，可提取能夠表征其對(duì)應(yīng)磨機(jī)負(fù)荷狀態(tài)的振動(dòng)特征，為磨機(jī)負(fù)荷預(yù)測(cè)提供依據(jù)。

振動(dòng)信號(hào)；磨機(jī)負(fù)荷；AR模型；頻域特征

磨礦是選礦生產(chǎn)流程中一個(gè)很重要工序，也是選礦的耗能大戶[1]?？茖W(xué)、準(zhǔn)確地預(yù)測(cè)球磨機(jī)內(nèi)部負(fù)荷并開(kāi)發(fā)磨礦優(yōu)化控制技術(shù)是實(shí)現(xiàn)選礦廠節(jié)能降耗、提質(zhì)提產(chǎn)的根本任務(wù)之一[2]。目前，常采用的間接檢測(cè)方式是振動(dòng)法，振動(dòng)法是通過(guò)檢測(cè)磨機(jī)筒體、軸承等部位的振動(dòng)信號(hào)，采用合適的振動(dòng)信號(hào)處理方法，建立磨機(jī)內(nèi)部負(fù)荷與振動(dòng)信號(hào)的相關(guān)性，來(lái)預(yù)測(cè)磨機(jī)負(fù)荷，此方法具有靈敏度高、干擾信號(hào)小等優(yōu)點(diǎn)[3]。在信號(hào)處理方面，可分為經(jīng)典功率譜估計(jì)和現(xiàn)代功率譜估計(jì)，而現(xiàn)代功率譜估計(jì)方法中以參數(shù)模型功率譜估計(jì)為代表。其中的自回歸（Auto Regressive，AR）模型，是利用前期若干時(shí)刻的隨機(jī)變量來(lái)描述以后某時(shí)刻隨機(jī)變量的線性回歸模型，對(duì)非平穩(wěn)、非線性信號(hào)處理有較精確的譜估計(jì)及較優(yōu)良的譜分辨率，廣泛應(yīng)用于航行軌跡、特征提取、故障識(shí)別、負(fù)荷預(yù)測(cè)等工程實(shí)際問(wèn)題[4]。因此本文選擇小波及AR模型對(duì)軸承振動(dòng)信號(hào)進(jìn)行分析處理[5]。

1 試驗(yàn)設(shè)計(jì)

1.1 試驗(yàn)設(shè)備

采用由北京礦冶研究總院機(jī)械研究所生產(chǎn)的型號(hào)為 330mm×330mm Bond功指數(shù)球磨機(jī)進(jìn)行磨礦試驗(yàn)[5]，如圖1（a）所示。為獲取球磨機(jī)軸承座振動(dòng)信號(hào)，在球磨機(jī)軸承座上安裝由東華公司生產(chǎn)的測(cè)量振動(dòng)信號(hào)的DH131加速度傳感器，如圖1（b）所示。

圖1 球磨機(jī)及傳感器Fig.1 Ballm illand sensor

1.2 試驗(yàn)方案

變頻器頻率選擇為50 Hz（轉(zhuǎn)速率為0.8），入料粒度1～9mm，鋼球直徑分別為20mm、30mm、40mm。按料球比均在0.6條件下，設(shè)計(jì)了10%，20%，30%，40%，50%五種不同填充率的試驗(yàn)方案，根據(jù)文獻(xiàn)[6]公式，可以計(jì)算出所填充的入料粒級(jí)分布的各級(jí)鋼球的重量，具體如表1所示。

表1 試驗(yàn)方案Tab.1 Experimentalscheme

1.3 試驗(yàn)結(jié)果

為了評(píng)價(jià)球磨機(jī)在不同負(fù)荷狀態(tài)時(shí)的磨礦效果，根據(jù)表1中的方案進(jìn)行磨礦試驗(yàn)，以能耗和-0.074μm粒度產(chǎn)品的產(chǎn)率為指標(biāo)，得到如表2試驗(yàn)結(jié)果、負(fù)荷狀態(tài)形式及預(yù)測(cè)指標(biāo)。

表2 球磨機(jī)負(fù)荷狀態(tài)形式及預(yù)測(cè)指標(biāo)Tab.2 State form and p rediction index of ballm ill load

2 振動(dòng)信號(hào)時(shí)域統(tǒng)計(jì)特征的提取與分析

選取振動(dòng)信號(hào)有量綱量（峰峰值、均值、標(biāo)準(zhǔn)差）和無(wú)量綱量（偏度、峭度）提取振動(dòng)信號(hào)的時(shí)域特征。根據(jù)文獻(xiàn)[7-11]，選擇10 s的振動(dòng)數(shù)據(jù)用于特征提取。時(shí)域特征有量綱量曲線如圖2所示，時(shí)域特征有量綱量樣本平均值如表3所示。

從圖2和表3可知，峰峰值的平均值隨著填充率的增大而增大。由圖2（b）可知，隨著填充率從10%增大至50%，均值都大于0，填充率為40%時(shí)，均值最大，其余負(fù)荷參數(shù)間有部分交集，區(qū)分度相對(duì)較小。由圖2（c）可知，標(biāo)準(zhǔn)差變化的趨勢(shì)整體與峰峰值大致相同，但填充率為10%與填充率為20%時(shí)，標(biāo)準(zhǔn)差較接近。時(shí)域特征參數(shù)峰峰值、均值、標(biāo)準(zhǔn)差均能很好地表征不同負(fù)荷狀態(tài)下振動(dòng)信號(hào)的變化和總體趨勢(shì)。

圖2 有量綱時(shí)域特征曲線Fig.2 Dimensional tim edomain featureof sam plemean

表3 時(shí)域特征有量綱量樣本平均值Tab.3 Dimensionless tim edom ain characteristicscurve

為了分析不同時(shí)域統(tǒng)計(jì)特征與磨機(jī)負(fù)荷的相關(guān)性，振動(dòng)信號(hào)的無(wú)量綱時(shí)域特征曲線，如圖3所示。

由圖3（a）知，當(dāng)填充率為20%、30%時(shí)，偏度值在0附近變化，其時(shí)域波形大致呈正態(tài)分布形式；填充率10%時(shí)，偏度值小于0，時(shí)域波形呈左偏；填充率40%時(shí)，對(duì)應(yīng)的偏度值大于0，時(shí)域波形呈右偏；填充率為50%時(shí)，對(duì)應(yīng)的偏度增加幅度較大，偏度值在0.5附近，振動(dòng)信號(hào)時(shí)域波形呈右偏。由圖2（b）知，不同負(fù)荷狀態(tài)下振動(dòng)信號(hào)的峭度曲線呈現(xiàn)“鋸齒”狀，周期約為2 s。當(dāng)填充率為20%時(shí)，峭度值最大值約為20。在球磨機(jī)的磨礦過(guò)程中，內(nèi)部的負(fù)荷介質(zhì)隨球磨機(jī)筒體的轉(zhuǎn)動(dòng)而進(jìn)行復(fù)雜的運(yùn)動(dòng)，拋落區(qū)的鋼球會(huì)以一定的周期對(duì)底部區(qū)域的鋼球、物料或筒體產(chǎn)生沖擊。峭度曲線可描述球磨機(jī)拋落區(qū)介質(zhì)的運(yùn)動(dòng)狀態(tài)，表征磨機(jī)振動(dòng)信號(hào)的沖擊特性。

圖3 無(wú)量綱時(shí)域特征曲線圖Fig.3 Dimensionless timedomain characteristicscurve

由圖2、圖3的振動(dòng)信號(hào)的時(shí)域特征可知，有量綱時(shí)域特征參數(shù)可以很好地表征特征信號(hào)與負(fù)荷狀態(tài)參數(shù)的相關(guān)性；無(wú)量綱時(shí)域特征參數(shù)可以很好地解釋振動(dòng)時(shí)域信號(hào)的波形分布以及振動(dòng)沖擊特征。通過(guò)提取不同負(fù)荷狀態(tài)下的振動(dòng)信號(hào)的時(shí)域特征參數(shù)，根據(jù)時(shí)域特征參數(shù)與負(fù)荷狀態(tài)之間的關(guān)系，可為描述不同振動(dòng)信號(hào)的物理特征提供依據(jù)。

3 基于AR模型的信號(hào)頻域特征提取

為了更加全面地提取振動(dòng)的特征信號(hào)，采用AR模型對(duì)不同負(fù)荷狀態(tài)下的振動(dòng)信號(hào)進(jìn)行功率譜估計(jì)，從統(tǒng)計(jì)的角度描述隨機(jī)信號(hào)的頻域特性[12]，尋找能夠反映磨機(jī)負(fù)荷的頻域特征。

3.1 AR模型的建立

觀測(cè)數(shù)據(jù)的離散序列x（n）可用線性差分方程表示，如式（1）所示。

式中：w（n-i）表示白噪聲序列。

對(duì)式（1）進(jìn)行z變換可得式（2）：

假設(shè)輸入白噪聲的譜密度Pw（z）=σw2，變換可得式（3），即：

式（4）即為AR模型，p為AR模型的階數(shù)，aj為AR模型參數(shù)，j=1，2，...，p。

3.2 AR模型參數(shù)求解

AR模型參數(shù)可通過(guò)從觀測(cè)數(shù)據(jù)序列的自相關(guān)序列中計(jì)算出固定階次AR模型參數(shù)，從而得到觀測(cè)數(shù)據(jù)序列的功率譜估計(jì)。

根據(jù)自相關(guān)序列定義，n時(shí)刻的自相關(guān)函數(shù)，如式（5）所示。

令bi=0，b0=1，則系統(tǒng)差分方程可表示為：

式中：m=0，1，…，p。

將式（5）按照矩陣形式展開(kāi)，并根據(jù)自相關(guān)函數(shù)的偶對(duì)稱性質(zhì)，整理可得式（6）：

3.3 AR模型最佳階數(shù)的確定

AR模型的階次對(duì)功率譜估計(jì)有重要影響，AR模型的最佳階次范圍p在N/2～N/3之間[13]。根據(jù)不同信號(hào)的特征，計(jì)算得到最佳的AR模型階次。采用信息論準(zhǔn)則（AIC）法，選擇AR模型的階次[11]，即

式中：p為模型的階次；N為數(shù)據(jù)個(gè)數(shù)；σp2為不同階次模型的預(yù)測(cè)誤差。

根據(jù)表1磨機(jī)負(fù)荷狀態(tài)的試驗(yàn)參數(shù)，其振動(dòng)信號(hào)的AIC值與AR模型階次的變化情況，如圖4所示。

圖4 AIC值隨模型階數(shù)p變化曲線圖Fig.4 Curveof AIC valuewith thevariation ofm odelorder p

由圖4可知，不同磨機(jī)負(fù)荷對(duì)應(yīng)的振動(dòng)信號(hào)的信息論函數(shù)AIC值，隨著模型階次p先減小后增大，并相應(yīng)取得極小值，且隨著磨機(jī)負(fù)荷增大，對(duì)應(yīng)振動(dòng)信號(hào)的信息論函數(shù)AIC值呈遞增趨勢(shì)。

AIC值取得極小值時(shí)，對(duì)應(yīng)AR模型取得最佳階次，當(dāng)填充率為10%時(shí)，其振動(dòng)信號(hào)模型階次p1= 594，對(duì)應(yīng)的AIC值最大為7.03，即AR模型最佳階次是594。負(fù)荷參數(shù)料球比均在0.6，填充率分別為10%、20%、30%、40%、50%對(duì)應(yīng)的AR模型功率譜估計(jì)的最佳階次和對(duì)應(yīng)的AIC值，如表4所示。

由表4可知，不同負(fù)荷狀態(tài)下的振動(dòng)信號(hào)AR模型的最佳階次和信息論函數(shù)AIC值各不相同。針對(duì)不同振動(dòng)信號(hào)的特征，利用信息論函數(shù)AIC得到其最佳階次，避免了不同振動(dòng)信號(hào)使用相同的AR模型階次進(jìn)行功率譜估計(jì)帶來(lái)的誤差。

表4 不同負(fù)荷AR模型的最佳階次和AIC值Tab.4 Optimalorder and AIC values for different load ARmodels

圖5 正常負(fù)荷狀態(tài)振動(dòng)信號(hào)的AR模型功率譜估計(jì)能量分布Fig.5 ARmodelpower spectrum estimation energy distribution of the vibration signalof thenorm al load state

4 AR模型功率譜估計(jì)及分析

根據(jù)表4中得到AR模型的最佳階次，分別對(duì)不同負(fù)荷參數(shù)所對(duì)應(yīng)的振動(dòng)信號(hào)進(jìn)行AR模型功率譜估計(jì)，得到其在頻域的能量分布。以表2磨機(jī)負(fù)荷狀態(tài)為試驗(yàn)參數(shù)，采用AR模型對(duì)其振動(dòng)信號(hào)分別進(jìn)行功率譜估計(jì)。其中，正常負(fù)荷狀態(tài)的振動(dòng)信號(hào)的功率譜估計(jì)能量分布，如圖5所示。

從圖5（a）、圖5（b）、圖5（c）可以看出，將正常負(fù)荷狀態(tài)下對(duì)應(yīng)的振動(dòng)信號(hào)做AR功率譜估計(jì)，頻譜曲線主要分布在頻率段2 000～3 635 Hz和3 800～8 500Hz之間。通過(guò)對(duì)比圖5（a）、圖5（b）、圖5（c）可知，低頻段2 000～3 635 Hz的頻譜曲線圖，頻譜大致以3 065 Hz為分界，并隨著負(fù)荷參數(shù)填充率的增大，低頻段的左邊頻段2 000～3 065 Hz在2 000～3 635Hz內(nèi)，能量占比逐漸較小。由圖5（c）可知，當(dāng)負(fù)荷參數(shù)填充率為40%時(shí)，頻譜大致以3065Hz呈對(duì)稱狀態(tài)，而高頻率段3800～8500Hz頻譜能量整體增大。表1中，非正常負(fù)荷參數(shù)的振動(dòng)信號(hào)的功率譜能量分布，如圖6所示。

從圖6（a）、圖6（b）可知，將非正常負(fù)荷狀態(tài)下對(duì)應(yīng)的振動(dòng)信號(hào)進(jìn)行AR功率譜估計(jì)，頻譜曲線主要分布在頻率段2 000～3 635 Hz和3 800～8 500Hz之間。由圖6（a）可知，當(dāng)負(fù)荷參數(shù)填充率為10%時(shí)，其功率譜曲線的能量分布在頻率段2 000～3 635Hz，以3 065 Hz為分界，呈現(xiàn)“左大右小”分布，在高頻段3 800～8 500Hz以5 830Hz呈對(duì)稱狀態(tài)。由圖6（b）可知，當(dāng)負(fù)荷參數(shù)填充率為50%時(shí)，其功率譜曲線的能量分布在頻率段2 000～3 635Hz，以3 065Hz為分界，呈現(xiàn)“左小右大”分布，在高頻段3 800～8 500 Hz內(nèi)3 800～5 830 Hz的能量值增大，而頻率段5 830～8 500Hz的能量值減小到0。從圖6不同負(fù)荷狀態(tài)振動(dòng)信號(hào)的AR模型功率譜估計(jì)能量分布可知，隨著充填率從10%增大到50%，即填充率的不斷增大，振動(dòng)信號(hào)在頻域的特征頻段主要在2 000～8 500Hz之間變化。從不同負(fù)荷狀態(tài)的振動(dòng)信號(hào)頻譜分布可知，隨著磨礦過(guò)程中筒體中鋼球、物料的變化，特征頻段主要分成低頻段2000～3635Hz，高頻段3 635～8 500 Hz的頻譜能量值曲線會(huì)相應(yīng)的變化。根據(jù)頻譜能量分布變化規(guī)律，就可提取能夠表征其對(duì)應(yīng)磨機(jī)負(fù)荷狀態(tài)的振動(dòng)特征。

圖6 非正常負(fù)荷狀態(tài)振動(dòng)信號(hào)的AR模型功率譜估計(jì)能量分布Fig.6 ARm odelpower spectrum estimation energy distribution of thevibration signalof abnormal load state

5 結(jié)論

（1）采用AR模型對(duì)振動(dòng)信號(hào)進(jìn)行特征提取和進(jìn)行功率譜估計(jì)，能較好地預(yù)測(cè)球磨機(jī)內(nèi)部負(fù)荷特征。

（2）在通過(guò)有量綱量（峰峰值、均值、標(biāo)準(zhǔn)差）和無(wú)量綱量（偏度、峭度）提取振動(dòng)信號(hào)的時(shí)域特征時(shí)，有量綱時(shí)域特征參數(shù)可以很好地表征特征信號(hào)與負(fù)荷狀態(tài)參數(shù)的相關(guān)性；無(wú)量綱時(shí)域特征參數(shù)可以很好地解釋振動(dòng)時(shí)域信號(hào)的波形分布以及振動(dòng)沖擊特征。

（3）隨著磨礦過(guò)程中筒體中鋼球、物料的變化，特征頻段主要分成低頻段2 000～3 635 Hz，高頻段3 635～8 500Hz的頻譜能量值曲線會(huì)相應(yīng)的變化。根據(jù)不同負(fù)荷狀態(tài)振動(dòng)信號(hào)的AR模型功率譜估計(jì)，可提取能夠表征其對(duì)應(yīng)磨機(jī)負(fù)荷狀態(tài)的振動(dòng)特征。參考文獻(xiàn)：

[1] 蔡改貧，熊 洋，姜志宏，等.基于預(yù)磨機(jī)入料粒級(jí)分布建模和參數(shù)辨識(shí)研究[J].礦山機(jī)械，2014，42（9）：71-74. CAIGaipin，XIONG Yang，JIANG Zhihong，et al.Research of parameter identification and into the material graded distribution modelbasedonpre-mill[J].MiningMachinery，2014，42（9）：71-74.

[2] 秦 虎，劉志紅，黃宋魏.碎礦磨礦及浮選自動(dòng)化發(fā)展趨勢(shì)[J].云南冶金，2013，39（3）：13-16. QIN Hu，LIU Zhihong，HUANG Songwei.Tendency of automation developments and broken ore grinding[J].Yunnan Metallurgy，2013，39（3）：13-16.

[3] BEHERA B，MISHRA BK，MURTY CV R.Experimental analysis of charge dynamicsin tumbling millsbyvibration signature technique[J].Minerals Engineering，2007，20（1）：84-91.

[4] 蔡改貧，許 琴，熊 洋，等.球磨機(jī)磨礦系統(tǒng)多變量辨識(shí)建模的研究[J]，中國(guó)鎢業(yè)，2015，30（5）：68-72. CAIGaipin，XU Qin，XIONGYang，et al.Research of identification ofmultivariatemodeling of ballmill[J].China Tungsten Industry，2015，30（5）：68-72.

[5] 羅小燕，盧小江，熊 洋，等.小波分析球磨機(jī)軸承振動(dòng)信號(hào)特征提取方法[J].噪聲與振動(dòng)控制，2016（1）：148-152. LUO Xiaoyan，LU Xiaojiang，XIONG Yang，et al.Extraction of features of Bearing vibration signal of ballmill based on wavelet analysis[J].Noiseand Vibration Control，2016（1）：148-152.

[6] 曹雪麗.介質(zhì)充填率和料球比對(duì)球磨機(jī)磨礦效果的影響[D].昆明：昆明理工大學(xué)，2011. CAO Xueli.Effectofmedium change ratio and the ratio ofmaterial tomedia on ballmill[D].Kunming：Kunming University of Science and Technology，2011.

[7] 譚仁鵬.球磨機(jī)振動(dòng)信號(hào)采集與無(wú)線傳輸裝置的設(shè)計(jì)與開(kāi)發(fā)[D].沈陽(yáng)：東北大學(xué)，2009. TANRenpeng.Designand developmentofwirelessdata transmission device and vibration signals collecting of ballmill[D].Shenyang：Northeastern University，2009.

[8] 梁朝霞.球磨機(jī)負(fù)荷檢測(cè)軟測(cè)量建模方法的研究[D].沈陽(yáng)：東北大學(xué)，2008. LIANGChaoxue.Research of softmeasurementmodelingmethod of loadonballmill[D].Shenyang：NortheasternUniversity，2008.

[9] 馮 雪，趙立杰，鄭瀚洋.基于KPCA頻譜特征提取的球磨機(jī)負(fù)荷檢測(cè)方法[J].沈陽(yáng)化工大學(xué)學(xué)報(bào)2014，28（2）：175-180. FENG Xue，ZHAO Lijie，ZHENG Hanyang.Detection method of ballmill load based on KPCA[J].Journal of Shenyang Chemical University，2014，28（2）：175-180.

[10]TANG J，ZHAO L J，YUW，etal.Soft sensormodeling of ball bill load viaprincipal component analysis and support vectormachines [M].Advances in Neural Network Research and Applications，Belin：Springer Belin Heidelberg，2010，67（67）：803-810.

[11]關(guān) 山.在線金屬切削刀具磨損狀態(tài)監(jiān)測(cè)研究的回顧與展望Ⅱ：信號(hào)特征的提取[J].機(jī)床與液壓，2010，38（17）：121-125. GUANShan.Retrospectand prospectof cutting toolwear condition monitoringonmetal removal：Theextraction ofsignalcharacteristics [J].Machine Tool&Hydraulics，2010，38（17）：121-125.

[12]HERNANDEZFE，MEDINACO.Theapplication ofbispectrum on diagnosis of rolling elementbearin-gs：a theoretical approach[J]. Mechanical Systemsand Signal Processing，2008，22（3）：588-596.

ARM odel-based Extraction for Ball-m ill Vibration Signals

XUYongqiang,JIANGZhihong,CAIGaipin,GUO Jinshan,XIONGYang

（SchoolofMechanicaland ElectricalEngineering,JiangXiUniversityofScienceand Technology,Ganzhou 34100,Jiangxi,China)

To study the fill-up rate andmaterial-to-ball ratio of the ball-mill in grinding process,ARmodel for feature extraction of the bearing's vibration signals of themillwas proposed.With a testingmill as the object,the ball-millbearing vibration signalswere collected.By analyzing the correlation of time domain signal characteristics and parameters of themill load identify the type ofmill load (underload,normal,overload);to get the frequency domain character reflect themill load,the vibration signal of differentmill load state were estimated into power spectrum based on ARmodel,which provide the evidence for themill load forecasting.

vibration signal;mill load;ARmodel;frequency domain character

TD4；TH133．3

A

10.3969/j.issn.1009－0622.2017.02.013

2016-11-11

國(guó)家自然基金項(xiàng)目（51464017）；江西省高等學(xué)校科技落地計(jì)劃項(xiàng)目立項(xiàng)（KJLD13045）；江西省研究生創(chuàng)新專項(xiàng)基金項(xiàng)目（YC2015-S284）

許永強(qiáng)（1992-），男，江蘇鹽城人，碩士研究生，研究方向：礦山機(jī)械和物料高效破碎技術(shù)。

姜志宏（1977-），男，新疆奎屯人，副教授，主要從事固體力學(xué)、物料高效破碎技術(shù)研究。