巖礦薄片偏光序列圖像的礦物顆粒分割方法①

吳耀坤,滕奇志,何海波,彭志偉

1(四川大學(xué) 電子信息學(xué)院 圖像信息研究所,成都 610065)

2(成都西圖科技有限公司,成都 610024)

圖像處理技術(shù)在巖礦薄片鑒定分析中的應(yīng)用已經(jīng)成為地質(zhì)行業(yè)的常用手段.礦物顆粒分割是薄片鑒定和礦物成分識(shí)別的基礎(chǔ),其效果直接影響后續(xù)鑒定分析的準(zhǔn)確性.利用巖礦薄片在固定視域下的正交偏光序列圖像進(jìn)行礦物顆粒的提取是一種新的研究方法,已取得一定成果.在實(shí)際生產(chǎn)應(yīng)用中主要采用基于灰度閾值的分割算法[1]、基于統(tǒng)計(jì)區(qū)域合并(SRM)的算法[2]以及基于邊緣流[3]的算法.

超像素分割算法[4]的最大優(yōu)勢(shì)在于邊緣精準(zhǔn)定位以及算法復(fù)雜度低.本文使用基于熵率的超像素算法[5,6]對(duì)偏光序列融合圖像進(jìn)行初步分割,通過(guò)快速區(qū)域合并算法[7,8]合并具有相似性的過(guò)分割區(qū)域,同時(shí),基于礦物顆粒在正交偏光序列圖像中的變化規(guī)律,對(duì)結(jié)果進(jìn)行二次合并,改善超像素算法過(guò)分割嚴(yán)重的問(wèn)題,實(shí)現(xiàn)了對(duì)巖礦薄片圖像中顆粒的分割提取.本文分割算法的流程圖如圖1所示.

圖1 本文算法流程圖

1 基于偏光序列圖像的超像素分割

超像素[9]是指具備相似顏色、灰度、紋理等特征的相鄰像素組成的圖像塊,超像素分割實(shí)質(zhì)上就是依據(jù)各種特征將相似或相同的像素聚集到一起的聚類問(wèn)題.超像素分割算法在國(guó)內(nèi)外已廣泛應(yīng)用,主要分為基于圖論和梯度下降的兩類方法,本文采用基于圖論的熵率超像素分割算法對(duì)巖礦薄片偏光序列圖進(jìn)行礦物顆粒分割.

1.1 熵率超像素分割

熵率法超像素[5]提出了將圖上隨機(jī)游走熵率和平衡項(xiàng)相結(jié)合的目標(biāo)函數(shù),并在聚類過(guò)程中不斷優(yōu)化目標(biāo)函數(shù)實(shí)現(xiàn)對(duì)圖像的分割.該方法分割出的超像素塊相對(duì)規(guī)則、均勻,其主要思想分為以下四個(gè)方面.

(1)圖的構(gòu)造:采取映射的方式,將圖片映射至無(wú)向圖G=(V,E),其中V表示像素點(diǎn)集,E為邊集,并將像素間的相似性定義為邊的權(quán)值 ω .目標(biāo)是尋找到集合A?E,使生成的無(wú)向圖G=(V,A)中包括K個(gè)連通的子圖.

(2)隨機(jī)游走熵率:熵率用來(lái)度量隨機(jī)過(guò)程的不確定性,可以很好的表征信息的冗余程度、內(nèi)部結(jié)構(gòu)的緊湊性,隨機(jī)游走模型的加入能夠獲取緊湊、均勻的超像素區(qū)域塊.無(wú)向圖G=(V,A)上的隨機(jī)游走熵率計(jì)算公式如下:

(3)平衡條件:在目標(biāo)函數(shù)中添加平衡條件,可以約束限制熵率的隨機(jī)性.令A(yù)為 已選擇邊集,劃分出A中超像素塊集為為G=(V,A)中連通子集的數(shù)目,ZA為子圖聚類的概率分布:

由此,定義平衡函數(shù)為:

其中,ZA的熵H(ZA)可以使聚類結(jié)果有相似的尺寸,同時(shí)NA約束集群的數(shù)目.圖2中用不同粗細(xì)的線條連接代表不同的集群,可見圖2(a)中的平衡項(xiàng)大,集群分布更均勻.

圖2 平衡條件對(duì)集群的作用

(4)目標(biāo)函數(shù):組合隨機(jī)游走熵率和平衡條件,得到聚類最優(yōu)化的目標(biāo)函數(shù)如式(5)所示,式中集合滿足A?E,且 λ為正實(shí)數(shù).聚類過(guò)程中通過(guò)最大化目標(biāo)函數(shù)得到圖像的最佳分割結(jié)果.

1.2 基于偏光序列圖像的熵率超像素分割

巖礦薄片在單偏光下獲取的圖像中,礦物顆粒之間粘連比較嚴(yán)重、缺少清晰的邊界,當(dāng)顆粒顏色相差不大時(shí),很難辨認(rèn)出具體是幾個(gè)顆粒,因此,給顆粒分割帶來(lái)極大的困難.在正交偏光鏡下,礦物顆粒在視域下會(huì)表現(xiàn)出一系列特有的光學(xué)特性[10],如消光特性、干涉色等.正交偏光下的消光特性[10]在非均質(zhì)礦物中普遍存在,當(dāng)正交偏光鏡旋轉(zhuǎn)一周,礦物共會(huì)出現(xiàn)四次消光現(xiàn)象,并且,不同礦物顆粒在不同偏光角度下消光狀態(tài)也可能不一樣,加上礦物顆粒在正交偏光下的干涉色,由此顯現(xiàn)出明顯的邊緣.本文所需的序列圖像由自主研制的偏光圖像采集與拼接系統(tǒng)獲取,先在單偏光下采集視域的一張圖像,然后保持顯微鏡載物臺(tái)固定,每旋轉(zhuǎn) 10°正交偏光鏡采集一張圖像,共旋轉(zhuǎn) 120°,共采集13張固定視域下的正交偏光圖像,圖3為部分序列圖像.

圖3 偏光序列圖像示意圖

考慮到礦物顆粒本身的復(fù)雜性,以及單個(gè)正交偏光角度下礦物顆粒的不完整性給分割帶來(lái)的影響,首先選取一個(gè)消光周期(90°)內(nèi)的正交偏光序列圖像,即拍攝的前十張正交偏光序列圖像,依據(jù)式(6)進(jìn)行圖像差分[11]融合,形成一張擁有視域中全部顆粒的正交偏光圖像,如圖4(a)所示.

其中,g(x,y)為融合后圖像,fi(x,y)是序列圖中第i張正交偏光圖像,且g(x,y)的R、G、B分量gR(x,y)、gG(x,y)、gB(x,y)均已映射至 0 ～255.

為了減弱R、G、B分量間的相互影響,提高礦物顆粒分割的效果,本文將融合后的圖像g(x,y)轉(zhuǎn)化至Lab空間.對(duì)Lab空間下的融合圖像采用1.1節(jié)中的算法進(jìn)行圖像分割,其分割結(jié)果如圖4(b)所示.

圖4 熵率超像素分割

2 礦物顆粒目標(biāo)的區(qū)域合并

2.1 基于融合圖像與單偏光圖像的快速區(qū)域合并

由圖4(b)的結(jié)果可以看出,利用超像素分割算法可以完整的保留礦物顆粒的邊緣信息,但是使礦物顆粒出現(xiàn)過(guò)分割現(xiàn)象.一般而言,屬于同一個(gè)礦物顆粒的區(qū)域具有相似性,利用其相似性將同質(zhì)區(qū)域合并.本文首先采用文獻(xiàn)[7]中提出的基于RAG(區(qū)域鄰接圖)和NNG(最近鄰接圖)數(shù)據(jù)結(jié)構(gòu)的快速區(qū)域合并算法對(duì)同質(zhì)區(qū)域進(jìn)行合并,其中合并判定準(zhǔn)則結(jié)合了目標(biāo)區(qū)域融合圖像的差分直方圖和單偏光圖像的灰度直方圖信息.

(1)相似性合并準(zhǔn)則的確定

① 融合圖像的差分直方圖:差分直方圖反映了正交偏光下礦物顆粒圖像在一定方向、一定距離上相鄰點(diǎn)灰度之間差的概率分布[12,13].在正交偏光的融合圖像中,礦物顆粒本身的形態(tài)、紋理比較復(fù)雜,較顏色直方圖而言,差分直方圖具有顆粒目標(biāo)區(qū)域圖像的綜合紋理信息,因此,提取偏光序列融合圖像的差分直方圖作為相似度主要判定準(zhǔn)則.

設(shè)一幅具有M(M=1,2,···,L)灰度級(jí)的圖像F,區(qū)域D內(nèi)的兩個(gè)像素f1=f(x,y)、f2=f(x+dx,y+dy),(dx,dy)∈D,對(duì)于偏移(dx,dy)的像素差分按式(7)計(jì)算如下:

由此可得歸一化的差分直方圖如下:

式中,hd(i)為區(qū)域D中d(x,y)的分布,N為區(qū)域D的像素總數(shù).在計(jì)算每一個(gè)區(qū)域的差分直方圖矩陣時(shí),本文在每個(gè)像素的8鄰域中選取四個(gè)方向(dx,dy)=(d,0),(d,d),(0,d),(-d,d)進(jìn)行合并計(jì)算.

統(tǒng)計(jì)學(xué)中,可以用Bhattacharyya系數(shù)[14]計(jì)算兩個(gè)統(tǒng)計(jì)樣本的重疊量,從而進(jìn)行兩組樣本的相關(guān)性測(cè)量,其計(jì)算公式如式(9)所示.在進(jìn)行直方圖的相似性計(jì)算時(shí),該計(jì)算方法可以得到最好效果.

式中,p(x)、q(x)為兩組不同的概率分布,且0＜BC＜1 .

因此,根據(jù)式(9),定義融合圖像差分直方圖的相似度合并準(zhǔn)則hC為:

其中,Pdp、Pdq分別為兩不同區(qū)域p、q的差分直方圖.

② 單偏光灰度直方圖:單偏光圖像中,同一礦物顆粒所在區(qū)域的灰度值差異變化不大,可將區(qū)域的灰度值作為次要判定準(zhǔn)則.直方圖統(tǒng)計(jì)是提取區(qū)域灰度特征的常用方法,根據(jù)式(9)可得單偏光下灰度相似性合并準(zhǔn)則hA為:

其中,Hp、Hq分別為兩不同區(qū)域p、q的灰度直方圖.

綜合①、②定義本文方法中使用的相似性合并準(zhǔn)則如式(12)所示:

(2)快速區(qū)域合并算法

區(qū)域鄰接圖(RAG)的數(shù)據(jù)結(jié)構(gòu)可以表示圖像中全部區(qū)域之間的鄰近關(guān)系,利用該數(shù)據(jù)結(jié)構(gòu)并依據(jù)(1)中討論的相似性準(zhǔn)則實(shí)現(xiàn)區(qū)域合并,同時(shí),使用最近鄰接圖(NNG)實(shí)現(xiàn)基于RAG區(qū)域合并算法的加速優(yōu)化.具體操作過(guò)程如下:

① 初始化數(shù)據(jù)結(jié)構(gòu):計(jì)算所有超像素區(qū)域下對(duì)應(yīng)融合圖像的差分直方圖以及單偏光圖像的灰度直方圖.將全部的超像素區(qū)域作為單獨(dú)的節(jié)點(diǎn),按照區(qū)域鄰接圖(RAG)的規(guī)則生成一個(gè)無(wú)向鄰接圖,并據(jù)此生成對(duì)應(yīng)的最近鄰接圖(NNG);

② 搜索鄰接情況:采用基于NNG的搜索方法,對(duì)于每一個(gè)區(qū)域確定其所有的鄰接區(qū)域;

③ 合并:依據(jù)(1)中的相似性合并準(zhǔn)則計(jì)算兩鄰接區(qū)域的相似性,若達(dá)到閾值,將兩鄰接區(qū)域合并成一個(gè)區(qū)域;

④ 更新數(shù)據(jù)結(jié)構(gòu):合并之后計(jì)算該區(qū)域相應(yīng)的差分直方圖以及灰度直方圖數(shù)據(jù),并更新RAG和NNG;

⑤ 重復(fù)②～④,直至沒有區(qū)域需要合并為止.

圖5為合并前后對(duì)比圖.可以看出,使用該方法基本可以實(shí)現(xiàn)大部分礦物顆粒的過(guò)分割合并,但由于礦物顆粒本身的復(fù)雜性,仍有部分礦物顆粒的過(guò)分割區(qū)域未實(shí)現(xiàn)合并.本文在此結(jié)果基礎(chǔ)上,利用礦物顆粒在正交偏光序列圖像中的變化特性對(duì)未實(shí)現(xiàn)合并的礦物顆粒再次合并,最終達(dá)到較為理想的結(jié)果.

圖5 初步合并效果

2.2 基于正交偏光序列圖像的二次區(qū)域合并

2.2.1 顆粒目標(biāo)圖層的初步提取

采用熵率超像素算法可以得到良好的顆粒邊界,其中碎屑巖顆粒的邊界也較為理想.由于最終需要得到顆粒目標(biāo)圖形層進(jìn)行粒度分析、成分識(shí)別,因此,利用礦物顆粒在序列圖中的亮度特性采用如下步驟初步獲取顆粒目標(biāo)圖層:

① 對(duì)2.1節(jié)的邊緣結(jié)果進(jìn)行目標(biāo)反轉(zhuǎn)操作得到全部目標(biāo)區(qū)域;

② 依據(jù)式(13)對(duì)①中的所有區(qū)域進(jìn)行消光性篩選,將亮度差低于閾值T的區(qū)域置為背景色,獲得礦物顆粒區(qū)域;

其中,Imax、Imin分別為目標(biāo)區(qū)域在正交偏光的前后序列圖中的最大亮度均值與最小亮度均值,T為亮度篩選閾值條件.經(jīng)大量實(shí)驗(yàn),T取25篩選效果最佳;

③ 對(duì)②中的結(jié)果進(jìn)行去噪、填孔、數(shù)字形態(tài)學(xué)、邊界平滑等操作.

根據(jù)上述步驟初步得到礦物顆粒目標(biāo)的圖形層圖像如圖6所示.

圖6 亮度篩選結(jié)果

2.2.2 顆粒目標(biāo)圖層的消光性合并

經(jīng)過(guò)區(qū)域亮度篩選可基本獲取顆粒目標(biāo)區(qū)域,但仍存在部分顆粒被分為若干個(gè),因此,基于圖6的結(jié)果,利用礦物顆粒目標(biāo)在序列圖中的消光變化特性繼續(xù)對(duì)仍存在的過(guò)分割顆粒目標(biāo)區(qū)域進(jìn)行相似性合并.

對(duì)于圖7(a)中的完整顆粒區(qū)域W以及其在圖7(b)中的過(guò)分割區(qū)域A、B,每個(gè)目標(biāo)區(qū)域在正交偏光序列圖像下的亮度均值變化趨勢(shì)如圖8所示.

圖7 過(guò)分割顆粒示意圖

圖8 同一礦物顆粒的亮度變化趨勢(shì)

雖然礦物顆粒出現(xiàn)了過(guò)分割現(xiàn)象,由于其存在共同的消光特性,亮度均值在序列圖像中具有相同的變化趨勢(shì).在統(tǒng)計(jì)學(xué)中,相關(guān)系數(shù)用來(lái)描述兩個(gè)序列之間的相關(guān)程度,設(shè)兩個(gè)亮度均值序列為a=(a1,a2,···,aN),b=(b1,b2,···,bN),則以式(14)計(jì)算序列的a、b的相關(guān)系數(shù) βab:

其中,cov(a,b)為兩個(gè)序列的協(xié)方差,σa、σb分別為序列a、b的標(biāo)準(zhǔn)差.

因此,本文以兩區(qū)域在正交偏光序列圖下亮度均值序列的相關(guān)系數(shù)作為區(qū)域消光特性的合并準(zhǔn)則,當(dāng)兩區(qū)域的相關(guān)系數(shù) β小于閾值 βT時(shí),判定兩區(qū)域?yàn)橥活w粒目標(biāo),將兩目標(biāo)區(qū)域進(jìn)行合并.對(duì)初步提取的顆粒圖層中所有顆粒目標(biāo)進(jìn)行合并預(yù)測(cè),直至無(wú)顆粒目標(biāo)需要合并為止,得到最終的分割結(jié)果.

3 實(shí)驗(yàn)結(jié)果分析

本文基于熵率超像素算法,針對(duì)礦物顆粒在正交偏光前后序列圖像中的變化特性完成礦物顆粒圖層的提取.本文實(shí)驗(yàn)程序使用C++語(yǔ)言,在VS2015平臺(tái)下基于MFC框架編寫實(shí)現(xiàn).

經(jīng)過(guò)大量薄片圖像的測(cè)試,本文的分割算法能夠取得較好的分割效果.以圖3中礦物薄片的原始偏光序列圖像為例,將本文算法與文獻(xiàn)[2]、文獻(xiàn)[3]中的礦物顆粒分割算法的礦物顆粒提取結(jié)果進(jìn)行對(duì)比分析,如圖9所示.其中,圖9(b)為本文算法得到的顆粒圖層,圖9(c)、圖9(d)分別為SRM和邊緣流算法得到的顆粒目標(biāo)圖層.

圖9 本文算法與其他算法分割結(jié)果圖對(duì)比

結(jié)合該薄片的原始偏光序列圖像(圖3)以及原始融合圖像(圖9(a))可以看出,文獻(xiàn)[2]中利用SRM算法得到的分割結(jié)果容易出現(xiàn)大量的欠分割現(xiàn)象,使很多礦物顆粒目標(biāo)粘連在一起,并且一些礦物顆粒的邊緣定位不夠準(zhǔn)確;文獻(xiàn)[3]中利用邊緣流算法得到的分割結(jié)果大部分比較理想,由于邊緣流算法對(duì)邊界定位的敏感性,對(duì)于表面紋理特征較為復(fù)雜的堿性長(zhǎng)石易出現(xiàn)過(guò)分割情況且不易融合,使提取區(qū)域不夠完整,如圖9(d)中標(biāo)注出的顆粒目標(biāo).本文的分割算法可以較為完整地提取視域中的礦物顆粒,改善了堿性長(zhǎng)石的過(guò)分割現(xiàn)象,同時(shí),對(duì)于粒徑較大的碎屑巖顆粒也有了較為完整的提取結(jié)果,符合實(shí)際生產(chǎn)需求.

4 結(jié)論

本文針對(duì)固定視域下巖礦薄片的偏光序列圖像,提出一種礦物顆粒目標(biāo)提取分割的方法:首先將巖礦薄片在正交偏光下的序列圖像進(jìn)行差分融合后再基于熵率超像素算法分割礦物顆粒,然后提取分割區(qū)域在融合圖像下的差分直方圖以及單偏光圖像下的灰度直方圖作為相似性特征,對(duì)過(guò)分割區(qū)域進(jìn)行快速區(qū)域合并,最后通過(guò)礦物顆粒在正交偏光下的消光特性對(duì)分割區(qū)域篩選后再次進(jìn)行顆粒融合,從而得到最終的顆粒分割結(jié)果.

實(shí)驗(yàn)結(jié)果表明,該算法綜合了視域下偏光序列圖的特征信息,能夠較為準(zhǔn)確的提取礦物顆粒目標(biāo),并且,對(duì)碎屑巖顆粒的提取效果有一定的提升,為礦物顆粒的分割提取提供了一個(gè)新的方向.