基于spike-angle 圖和Fourier 變換的磨粒圖像形狀特征提取*

魯華杰 郭恒光 劉文彪

（海軍航空大學岸防兵學院 煙臺 264001）

1 引言

磨粒分析技術是一種有效的機械設備狀態監測方法，通過采集機械設備油液中的磨粒，分析磨粒所屬的類型，判斷設備的主要摩擦方式，為機械設備的維護和保養提供決策依據。磨粒的形狀特征是判斷磨粒所屬類型的重要依據，不同類型的磨粒，由于磨損類型和形成機理不同，會表現出不同的形狀特征，目前應用較多的磨粒形狀特征包括磨粒的基本幾何參數［1~2］、傅里葉描繪子［3~6］、分形維數［7~11］和Spike［12~14］參數等。

M.G.Hamblin 和G.W.Stachowiak 在1995 年 提出spike參數，它是在分形維數的基礎上，對磨粒輪廓進行更細致的描述，在每一步的起點與終點之間，引入兩點之間磨粒輪廓上的某一點，這3 個點構成一個“Spike”，然后對Spike 值進行計算，最后將考慮步長大小和數量情況下的Spike平均值作為Spike參數。

Spike 參數反映了該磨粒輪廓的角度變化狀況，該參數越大，磨粒輪廓越尖銳、磨損作用越大。基于Spike 參數的思想，本文提出了spike-angle 參數和spike-angle 圖的概念。首先計算磨粒的輪廓點在不同步長下的spike-angle參數，遍歷磨粒的所有輪廓點，得到磨粒的spike-angle 圖，spike-angle圖的行表示在單一步長下，磨粒所有輪廓點的spike-angle 參數，反映了磨粒輪廓的整體特征，spike-angle 圖的列表示磨粒的單個輪廓點在不同步長時的spike-angle參數，反映了磨粒輪廓的局部特征。然后對spike-angle 圖的行和列分別進行傅里葉變換，并對傅里葉變換系數的幅值進行歸一化。最后對每一階歸一化后的傅里葉變換系數幅值進行直方圖統計，得到反映磨粒輪廓整體特征和局部特征的直方圖，以這兩個直方圖作為磨粒的形狀特征。

2 Spike參數

Spike 參數的計算類似于分形維數的計算，如圖1所示，以輪廓點A為起點，按一定的步長找到輪廓點B 作為終點，其中步長是指起點和終點之間磨粒輪廓像素點的數量。在A 和B 之間的輪廓上尋找一點P，三個點構成一個三角形區域ΔAPB，稱之為“Spike”，用來描述磨粒位于起點和終點之間的輪廓。為了描述ΔAPB 的大小和尖銳程度，定義了Spike值Sv，如下式所示：

圖1 Spike參數計算原理圖

式中：θ為ΔAPB 中PA和PB之間的夾角，h為ΔAPB頂點P到邊AB的垂直距離。

Spike參數SP按下式計算：

式中：Svmax=max(cos(θ2)h)；m為在給定步長時，有效的Sv的數量；n為所使用的不同步長的數量。

3 spike-angle參數和spike-angle圖

由spike 參數的定義可以看出，Spike 參數的實質為通過起點和終點間三角形區域ΔAPB頂點P的角度平均值來描述磨粒輪廓的特征。但是在磨粒類型識別時，若僅使用Spike 參數作為磨粒形狀特征時，識別效果不佳。基于Spike參數的思想，本文提出了基于改進的Spike 參數，稱之為spike-angle參數（spike-angle parameter）。

為了使最后獲得的形狀特征具有尺度不變性，首先對磨粒二值圖像中磨粒大小進行歸一化，磨粒大小歸一化到Q，也就是說磨粒二值圖像強度的和為Q。然后提取磨粒的輪廓，并對磨粒輪廓點進行均勻采樣，得到N個輪廓點，最終的磨粒輪廓為

以第i個輪廓點為起點A，步長為s（最大步長為M），得到終點B，則A和B之間的輪廓點序列為

A點、B點和P點序列構成s個三角形，根據式（1）計算每一個三角形的Spike值。

則第i個輪廓點，在步長為s時的spike-angle 參數的定義為

式中：Svmax=max(cos(θk2)hk)；hmax為Sv最大時，對應的三角形區域的頂點P到邊AB的垂直距離。

依次取磨粒的每一輪廓點為起點，計算其在不同步長下的spike-angle，得到所有輪廓點在不同步長下的spike-angle，所有spike-angle 參數構成一個矩陣，將其定義為spike-angle圖。

spike-angle圖的行表示在單一步長下，磨粒所有輪廓點的spike-angle參數，反映了磨粒輪廓的整體特征，spike-angle圖的列表示磨粒的單個輪廓點在不同步長時的spike-angle參數，反映了磨粒輪廓的局部特征，同時由于計算spike-angle時采用了不同的步長，spike-angle圖又反映了磨粒輪廓的多尺度信息。因此，spike-angle圖包含了磨粒的整體形狀特征、局部形狀特征以及多尺度形狀特征。Spike-angle參數不受磨粒圖像旋轉和平移的影響，而且由于對磨粒圖像強度進行了歸一化處理，所以spike-angle圖具有平移、尺度和旋轉不變性。

圖2（a）、圖3（a）和圖4（a）分別為三種典型的磨粒：切削磨粒、嚴重滑動磨損磨粒和疲勞磨損磨粒，圖2（b）、圖3（b）和圖4（b）所示分別為三種典型的磨粒對應的spike-angle 圖。由圖2（b）、圖3（b）和圖4（b）可以看出，不同類型的磨粒，對應的spike-angle圖區別很顯著，而且spike-angle圖包含了磨粒的整體形狀特征、局部形狀特征、以及多尺度形狀特征，因此根據磨粒的spike-angle圖對磨粒形狀進行識別是可行的。spike-angle 圖的數據量大，直接用于磨粒形狀識別，效率低，需要耗費大量時間，這就需要研究從spike-angle 圖中提取信息，作為磨粒的形狀特征。

圖2 切削磨粒及其spike-angle圖

圖3 嚴重滑動磨損磨粒及其spike-angle圖

圖4 疲勞磨損磨粒及其spike-angle圖

4 磨粒圖像形狀特征提取

在形狀的輪廓描述和形狀識別中，廣泛使用傅里葉描述子，實現使用少數低價系數得到高質量的形狀描述的效果。本文對spike-angle 圖的行和列分別進行傅里葉變換，提取spike-angle圖中提取信息，進而得到磨粒的形狀特征。

對spike-angle 圖的第s行做傅里葉變換，相當于對磨粒的整體輪廓在步長為s時的spike-angle值序列做傅里葉變換，得到的傅里葉變換系數Global_as，n為

對spike-angle 圖的第i列做傅里葉變換，相當于對磨粒的第i個輪廓點在所有步長下的spike-angle值序列做傅里葉變換，得到的傅里葉變換系數Local_ai，n為

為了使得最后的形狀特征具有平移、尺度和旋轉不變性，取傅里葉變換系數的幅值并除以第一個系數的幅值，得到歸一化的系數幅值。因為少數低階系數可對形狀進行高質量的表示，因此可采用前m個歸一化后的系數幅值，作為磨粒形狀的特征，即：

由于傅里葉變換的第一個系數的幅值遠大于其他系數的幅值，因此歸一化后的傅里葉系數Global_FSAPs，j和Local_FSAPi，j的 取 值 區 間 為(0，1)。

圖5 為歸一化的傅里葉變換系數幅值曲線，其中實線、點線和虛線對應的分別為切削磨粒、嚴重滑動磨損磨粒和疲勞磨損磨粒。圖5（a）~（d）分別為對spike-angle圖的行做傅里葉變換之后，得到的前4 階歸一化系數幅值，圖5（e）~（h）分別為對spike-angle 圖的列做傅里葉變換之后，得到的前4階歸一化系數幅值。由圖5 中的每個圖中的曲線對比都可以看出，不同磨粒的歸一化系數幅值差別是比較明顯的，這也從一個角度證明了本文提出的spike-angle參數和spike-angle圖是可行的，并且包含了磨粒形狀的信息，可以從中提取磨粒的形狀特征，用于磨粒類型識別。

圖5 歸一化的傅里葉變換系數幅值

雖然只采用前m個歸一化后的傅里葉系數幅值，但是得到的系數矩陣仍然較大，也不能直接用作磨粒的形狀特征，下面就對每一階歸一化后的傅里葉變換系數幅度進行直方圖統計，以得到的直方圖作為磨粒的形狀特征。

由于每一階歸一化后的傅里葉系數幅值取值區間差別比較大，因此無法采用統一的區間劃分對所有的系數幅值進行直方圖統計。本文采用了以下的方法：采集大量的不同類型的磨粒，并用本文的方法得到所有磨粒的歸一化的系數矩陣Global_FSAP和Local_FSAP，確定Global_FSAP矩陣和Local_FSAP矩陣的每一個列的最大值，即每一階歸一化系數幅值的最大值。然后將每一階系數幅值的取值區間劃分為8 個小區間，即將0 到每一階系數幅值的最大值之間的區間劃分為8 個小區間，根據相應的劃分區間，分別對每一階歸一化后的傅里葉系數幅值進行直方圖統計，并對直方圖作歸一化處理，則得到的直方圖為

其中，直方圖的前8 個分量對應的是系數矩陣的第1 列系數的歸一化直方圖，直方圖的第9 個分量到第16 分量對應的是系數矩陣的第2 列系數的歸一化直方圖，直方圖的第8m-1 個分量到第m分量對應的是系數矩陣的第m列系數的歸一化直方圖。最后以所得到的兩個直方圖Global_H和Local_H作為磨粒的形狀特征，分別稱之為整體直方圖特征和局部直方圖特征。

圖6（a）、7（a）和8（a）分別為三種典型磨粒的整體直方圖特征，圖6（b）、7（b）和8（b）分別為三種典型磨粒的局部直方圖特征，提取直方圖特征時，對歸一化傅里葉系數幅值矩陣的前4 階系數進行直方圖統計。由各種典型磨粒的整體直方圖特征和局部直方圖特征的對比可以看出，不同類型磨粒的直方圖特征區別是很明顯的，各階歸一化傅里葉系數幅值出現的區間以及頻率相差也比較大，因此采用整體直方圖特征和局部直方圖特征作為磨粒的形狀特征是可行的。

圖6 切削磨損磨粒直方圖特征

圖7 嚴重滑動磨粒直方圖特征

圖8 疲勞磨損磨粒直方圖特征

5 磨粒圖像類型識別

以某型導彈發射車發動機為監測對象，定期采集其發動機潤滑油，并分離出其中的磨粒，然后使用光學顯微鏡拍照獲取磨粒圖像。根據研究需要，本文對其中的正常磨粒、切削磨粒、嚴重滑動磨損磨粒和疲勞磨損磨粒進行類型識別研究。

每種類型的磨粒分別獲得30個樣本，將120個樣本隨機分為兩組，其中1組70個樣本用于訓練算法，50 個樣本用于算法測試。第一步根據本文提出spike-angle參數計算方法，得到每個磨粒圖像的spike-angle圖；然后對spike-angle圖的行和列分別進行傅里葉變換，并對傅里葉變換系數幅值進行歸一化；最后對得到歸一化傅里葉系數幅值矩陣按列進行直方圖統計，獲得磨粒圖像的形狀特征。根據磨粒圖像的形狀特征，采用支持向量機的方法［15］，對磨粒類型進行識別，并與其它形狀特征提取方法做比較。

由于SVM 是一種二元分類器，只能對兩種類型進行識別，當對多種磨粒類型進行識別時，必須構造一種多元分類器才能進行多種磨粒類型的識別。本文采用1 對1 的方法構造多元分類器，對任意兩個類型磨粒，訓練一個SVM 將二者分開。因此，對于N元分類問題，需要構造N(N-1) 2 個SVM。本文對4 種磨粒類型進行識別時，需要建立和訓練的SVM為6個。

使用磨粒樣本對多元分類器訓練后，在不改變模型參數的情況下，采用另一組的50 個樣本對分類器進行測試。測試結果如表1 和表2 所示，其中表1 為訓練實驗結果，表2 為測試結果。由表1 和表2 可以看出，本文提出的基于spike-angle 圖的磨粒形狀特征提取方法，優于其他幾種常用的形狀特征。其中，Spike 參數和簡單的幾何參數對應的識別率較低，這兩個形狀特征不具有平移、尺度和旋轉不變性，而且只利用簡單的數值來表征磨粒形狀，不能充分反映磨粒的形狀特征。Zernike 矩和傅里葉描繪子對應的識別率要高于Spike參數和簡單的幾何參數，這兩個形狀特征都具有平移、尺度和旋轉不變性，反映了磨粒形狀的整體特征，也是應用較多的形狀特征。本文提出的磨粒形狀特征，包含了磨粒的整體形狀特征、局部形狀特征和多尺度形狀特征，能夠充分利用磨粒的形狀信息，因此在磨粒類型識別中能夠取得好的效果。

表1 磨粒類型識別訓練實驗結果

表2 磨粒類型識別測試實驗結果

6 結語

本文根據磨粒spike 參數的計算思想，提出了spike-angle值和spike-angle圖的概念，并且研究了基于spike-angle 圖的磨粒形狀特征提取方法。該方法首先計算磨粒的所有輪廓點在不同步長時的spike-angle 值，得到磨粒的spike-angle 圖；然后對spike-angle圖的行和列分別進行傅里葉變換，并對得到的傅里葉變換系數幅值進行歸一化處理；最后對每一階歸一化的傅里葉變換系數幅值進行直方圖統計，以獲得的整體直方圖特征和局部直方圖特征作為磨粒的形狀特征。本文所提的磨粒形狀特征提取方法，包含了磨粒形狀的整體信息、局部信息和多尺度信息，磨粒類型識別訓練實驗結果和磨粒類型識別測試實驗結果也表明了本文所提的方法的有效性。