周 冉,宋華珠,劉 翔
(武漢理工大學 計算機科學與技術(shù)學院,武漢 430070)
顆粒粒度檢測方法因能夠測試顆粒粒徑,已成為材料在制作加工或使用過程中控制粒徑大小的重要手段,在材料粒度分析領(lǐng)域有極大的研究與應用價值[1-2]。傳統(tǒng)的粒度測試方法有篩分法和激光法,前者靈活方便,后者精度較高,但同時前者精度難以達到要求,后者依賴于精密儀器導致成本高昂[3-4]。因此,研究人員提出了基于圖像分析顆粒粒度測試方法,該類方法簡單便捷、對儀器要求不高,且能通過算法改良來提高結(jié)果可靠性、精度等。同時,該類方法還擁有非入侵式的特點,能夠在不影響實驗結(jié)果的前提下,反映顆粒的沉降情況,對液體內(nèi)部的顆粒進行識別、跟蹤、粒度分析等。
近年來,許多國內(nèi)外專家致力于研究顆粒圖像的識別與分割。文獻[5]提出一種無參數(shù)回歸方法來為3D圖像信息降噪。文獻[6]改良了傳統(tǒng)非局部均值濾波(Non-local Mean Filtering,NMF)并結(jié)合粒子概率(Particle Probability,PP)還原被泊松-高斯噪音污染過的圖像。上述方法的圖像識別結(jié)果較好,但顆粒與背景區(qū)分度不高,并且灰度模糊。針對該問題,本文通過分析顆粒消失、聚合以及分散運動對成像灰度的影響,對識別圖像進行灰度增強、像素特征轉(zhuǎn)換、像素點連通性判斷等處理,并在此基礎(chǔ)上提出一種基于Haar-PP混合特征的沉降顆粒識別策略。
沉降顆粒的識別、跟蹤是一種顆粒運動的動態(tài)捕捉過程。由于在實際過程中無法有效觀測到沉降顆粒的運動過程,因此本文采用激光散射的方法實現(xiàn)對顆粒重量、尺寸、分布及聚集態(tài)結(jié)構(gòu)等信息的捕捉[7-8]。通過在顆粒運動過程中不斷地捕捉瞬態(tài)圖像,可得到一幀一幀的激光散射像點圖。實現(xiàn)顆粒跟蹤的重點在于完成運動狀態(tài)下顆粒狀態(tài)的有效識別。考慮到顆粒灰度存在的動態(tài)變化且有較多灰度值低(灰度模糊)的沉降顆粒,因此在識別上需要實現(xiàn)動態(tài)判定有效目標。本文制定的沉降顆粒激光散射像點識別策略如圖1所示。

圖1 沉降顆粒激光散射像點識別策略
針對灰度模糊及運動帶來的灰度變化問題,首先需要完成對每一幀中像點的灰度增強,改變圖像前景與背景的灰度對比度,降低噪聲影響;其次為了更好地進行PP特征的判斷,則必須對識別窗口進行規(guī)格大小的確定,并根據(jù)區(qū)域像素連通性計算顆粒面積;最后對特殊處理后所得的Haar特征像素點增加粒子概率屬性,普遍認為粒子概率越高的像素點被識別通過的概率越高。為了進一步提高準確率,本文在粒子概率的基礎(chǔ)上,通過概率加權(quán)的方法得到最終的識別結(jié)果。
對于顆粒特征的提取與識別研究,首先需要完成對顆粒群中不同個體的形態(tài)、運動分析。粒子形態(tài)特征的選取需要根據(jù)粒子本身的狀態(tài)來決定,如圖2(a)所示,在顆粒沉降視頻中,絕大多數(shù)沉降顆粒雖然大小相近,但其周長、面積、極軸、極軸夾角和質(zhì)心位置等各種數(shù)據(jù)特征卻不盡相同[9-10]。除了上述所提及的一般化的形態(tài)特征,也有小部分的沉降顆粒與一般顆粒因本身物理性質(zhì)或者拍攝原因會產(chǎn)生較大差別,主要體現(xiàn)在大顆粒與長尾顆粒上。在顆粒沉降視頻中,一般沉降顆粒的形態(tài)特征與圖2(a)相符,但是也有個別沉降顆粒的大小會大于圖2(a)中顆粒。如圖2(b)所示,沉降顆粒在沉降液體中發(fā)生結(jié)塊現(xiàn)象,因此整個顆粒的大小會大于一般顆粒的大小。由于攝像機清晰度的問題,在拍攝的顆粒沉降視頻中,顆粒可能會帶有較長的尾跡。沉降顆粒中的長尾顆粒如圖2(c)所示。

圖2 沉降顆粒的不同形態(tài)
顆粒運動的軌跡是無規(guī)則的,但在平面運動上可分解為水平與垂直方向。在大量的激光散射采樣中,通過對個體運動狀態(tài)的分析得出一般情況下的顆粒運動過程,如圖3(a)所示。圖3(a)為沉降顆粒在150幀內(nèi)發(fā)生的位移圖,實線末端表示顆粒的始末位置,虛線表示垂直移動距離。2條實線的差距代表沉降顆粒在150幀之后,水平方向移動距離。可以看出一般情況下,沉降顆粒一直做沉降運動并伴隨著水平方向上的微小位移。但如果受到外界擾動,則顆粒群會表現(xiàn)出明顯的類圓運動。譬如沉降環(huán)境中加入了新的沉降物質(zhì),顆粒的沉降運動就會明顯波動,如圖3(b)所示。圖3(b)左圖的沉降顆粒的總體沉降軌跡類似于圓弧。這2種現(xiàn)象可能同時發(fā)生,都屬于類圓運動的范疇。

圖3 顆粒的不同運動狀態(tài)
沉降顆粒的實際運動過程為空間三維運動,而表現(xiàn)在圖像識別中卻只能為二維運動,但這并不影響圖像處理對三維運動的判斷。在成像過程中,顆粒的遠近程度實質(zhì)上表現(xiàn)在像點的灰度值上。因此,沉降顆粒的散射像點圖中的灰度值是動態(tài)的,具體可以表現(xiàn)在2個方面:同一幀圖像中不同顆粒呈現(xiàn)不同灰度值,不同幀圖像中同一顆粒呈現(xiàn)不同灰度值。
對于單幀圖像中的顆粒群而言,由于遠近程度、大小、角度等因素,其灰度值是不同的。圖4(a)為圖3(b)中隨機2個沉降顆粒及其像素分布圖,可以看出左邊的顆粒灰度值較小,右邊較大,具體表現(xiàn)在對應的灰度直方圖中的曲線也明顯不同。本文同時進行了對同一顆粒不同幀的灰度分析,如圖4(b)所示。圖4(b)中為2種顆粒在不同4幀圖像中的灰度表現(xiàn),結(jié)果表明顆粒會隨著運動而表現(xiàn)出不同的灰度值。

圖4 2種形式的灰度特征
圖像預處理步驟如下:
1)灰度增強階段
顆粒識別前的圖像預處理策略,是為了使得識別環(huán)境更加友好,降低圖像噪音,增強顆粒信號,以達到改變前景與背景之間的對比度、凸顯沉降顆粒等目的[11]。混合重力沉降系統(tǒng)環(huán)境簡單,噪音干擾較少,因此只在識別前進行灰度增強處理,使得:
(1)

2)像素特征轉(zhuǎn)換
在像素轉(zhuǎn)換階段,核心內(nèi)容是確定識別窗口的大小。為了防止單一識別窗口效果不佳,本文以7×7為基礎(chǔ),每次增強1.5倍并設(shè)立多個識別窗口。共選擇7×7,11×11,17×17,25×25和37×37,5個識別窗口來綜合計算Haar特征信息,圖5為窗口確定示意圖。考慮到沉降顆粒的背景因素,在綜合運用多規(guī)格識別窗口的Harr特征信息時,將對5個窗口的Haar特征值進行加權(quán)處理。

圖5 以像素點d為中心的3個識別窗口
3)PP特征轉(zhuǎn)換
在PP特征轉(zhuǎn)換階段,根據(jù)像素概率理論,以像素點d為中心的局部區(qū)域內(nèi),將連通區(qū)域與局部區(qū)域面積之比作為概率值,因此需要對連通區(qū)域面積進行合理的算法設(shè)計。該算法從像素點d開始,向上下左右4個方向延伸,判斷鄰近像素點e是否符合連通條件,符合則以e為基礎(chǔ)繼續(xù)朝4個方向搜索,直至沒有其他像素點符合條件,或達到圖像邊界。該算法時間復雜度為O(n),n為連通區(qū)域內(nèi)像素個數(shù)。圖6為連通性示意圖。

圖6 像素點連通性判斷
算法1連通區(qū)域搜索算法
輸入輸入沉降顆粒灰度圖像I,輸入像素點d坐標(x,y)
輸出連通區(qū)域Ω
ConnectivitySearch(I,d)
Step 1:x0←x+1得到新像素點d0
if fh(d0)≥λh
then
d0計入連通區(qū)域Ω
Ω←Ω+ConnectivitySearch(I,d0)
else return φ
Step 2:x1←x-1得到新像素點d1
同Step 1
Step 3:y0←y-1得到新像素點d3
同Step 1
Step 4:y2←y+1得到新像素點d4
同Step 1
Haar特征是一種基于矩形區(qū)域灰度強度差來進行識別的特征描述算子[12-15]。本文建立3個Haar特征作為識別顆粒的條件之一,如圖7所示。以d為中心的識別窗口被劃分成黑色和白色2個子區(qū)域,白色區(qū)域被標記為A區(qū)域,黑色區(qū)域被標記為B區(qū)域。

圖7 3個不同的用于沉降顆粒識別的Haar特征
約定Haar特征函數(shù)為fh(·),fh(d,l)為像素點d在特征l下的Haar特征值,有
fh(d,l)=Bl-Al
(2)
其中,d為像素中心,l=1,2,3,代表圖7所示的7×7規(guī)格的沉降顆粒識別窗口(Haar特征識別窗口),Al代表白色區(qū)域的灰度值之和,Bl代表黑色區(qū)域灰度值之和。
最終的Haar特征值為3個不同大小Haar特征值中最大的一個,即:
(3)
當Haar特征值fh(d)高于所定閾值,則代表該像素點位于沉降顆粒中。
當算法遍歷完圖像中所有像素點,原圖像就轉(zhuǎn)化為Haar特征圖像,則有:
(4)

(5)
其中,thh為用于區(qū)分沉降顆粒的閾值。像素點d粒子概率的大小取決于識別窗口本身以及周圍像素點Haar特征值的大小。如果周圍超過Haar特征值閾值λh的像素越多,且與像素點d以區(qū)域Ω鄰近,那么像素點d越可能是沉降顆粒的一部分。如圖8所示,經(jīng)過計算可得fpp(d)>fpp(e),從圖8可以看出,像素點d的識別窗口內(nèi)基本是沉降顆粒,而像素點e的識別窗口內(nèi)只有很少量的顆粒是沉降顆粒,因此判斷像素點d為沉降顆粒可能性更大。

圖8 像素點d與e之間的粒子概率對比

(6)
(7)
沉降顆粒識別算法如算法2所示,該算法的時間復雜度為O(Height×Width)。其中,Height為沉降顆粒圖像高度,Width為寬度。
算法2沉降顆粒識別
輸入沉降顆粒灰度增強圖像I,二值化閾值thpw
輸出概率加權(quán)圖像Ipw
Do while l←1:3
Do while h←1:Height
Do while w←1:Width
Step1:得到對應像素點d,
Step2:計算
end
end
end
Do while h←1:Height
Do while w←1:Width
Step3:得到對應像素點d
Step4:計算像素點d的連通區(qū)域
Step5:計算fpp(d)

end
end
實驗數(shù)據(jù)來源于“異相顆粒共沉降裝置”中的顆粒沉降實驗。顆粒沉降實驗在混合重力沉降系統(tǒng)上完成,系統(tǒng)內(nèi)的激光與沉降顆粒碰撞形成反射,通過系統(tǒng)外部的攝像機在側(cè)面拍攝得到沉降顆粒激光散射像點視頻。實驗所得視頻共2個,視頻1長27 min 37 s,視頻2長20 min 7 s,包含25 frame/s,分辨率為576像素×960像素。對該視頻的解析以及對顆粒沉降圖像的分析借助于OpenCV3的視頻解析模塊與圖像分析模塊。在OpenCV3模塊的基礎(chǔ)上將本文算法形成代碼,識別沉降顆粒,跟蹤顆粒沉降過程,具體實驗環(huán)境如表1所示。

表1 實驗環(huán)境
為了能識別到盡可能多的沉降顆粒,將Haar特征閾值λh確定為典型粒子區(qū)域與典型背景區(qū)域灰度平均值的20%,經(jīng)過計算,thh=21。
通過Haar-pp混合特征策略得到的最終識別結(jié)果如圖9所示。圖9第2行和第3行為原圖像經(jīng)過灰度增強、Haar特征轉(zhuǎn)換、PP特征轉(zhuǎn)換后對比度增強的結(jié)果。圖9第1行和第4行分別為經(jīng)過Haar特征圖像和PP特征轉(zhuǎn)換得到概率加權(quán)圖像。從圖9第2行第3個圖像可以發(fā)現(xiàn),圖像從灰度增強圖像中識別出更多灰度值低的沉降顆粒,但是同時帶來的問題是噪音的增多。通過PP特征轉(zhuǎn)換來降低噪音帶來的干擾,處理結(jié)果如圖10所示。

圖9 沉降顆粒識別結(jié)果

圖10 轉(zhuǎn)換過程中的降噪結(jié)果
圖11描述各方法真正類率(True Positive Rate,TPR)和假正類率(False Positive Rute,FPR)的接受者操作特性曲線(Receiver Operating Characteristic Curve,ROC)。從圖11可以看出,較高的TPR代表本文算法從圖像中選出真的沉降顆粒的能力較高,且本文算法和文獻[6]算法優(yōu)于高斯混合模型與Ostu閾值分割算法。

圖11 不同方法的ROC曲線
本文算法在計算Haar特征方面做了簡化,需要的計算量為Harr特征轉(zhuǎn)換后識別的粒子的0.2倍,較文獻[6]算法提升22.3 s,且特征提取效果略優(yōu)于文獻[6]算法。表2描述了各種方法的平均計算時間。

表2 各算法平均計算時間對比 s
本文通過研究灰度變化與顆粒追蹤的關(guān)系,分析沉降顆粒形態(tài)、運動對其灰度的影響,并運用識別圖像進行灰度增強、像素特征轉(zhuǎn)換、像素點連通性判斷等消除灰度模糊的方法,提出一種基于Haar-PP混合特征的沉降顆粒識別算法。實驗結(jié)果表明,該算法能有效降低灰度模糊帶來的誤差,為粒子追蹤方面的研究提供理論依據(jù)。