本發(fā)明涉及地質(zhì)勘探技術(shù)領(lǐng)域，尤其涉及一種基于掃描圖像判斷碳酸鹽巖孔隙類型的方法。

背景技術(shù)：

碳酸鹽巖受巖相古地理、沉積構(gòu)造演化或者后期成巖作用等多種因素的共同影響，具有更為復(fù)雜的儲(chǔ)集空間和油氣滲流運(yùn)移規(guī)律，這導(dǎo)致碳酸鹽巖油藏開發(fā)的風(fēng)險(xiǎn)和難度較大。由此，對(duì)于正確選擇碳酸鹽巖的開發(fā)模式和開發(fā)方案，判定碳酸鹽巖儲(chǔ)層的類型顯得尤為重要。

判定碳酸鹽巖儲(chǔ)層類型的傳統(tǒng)模式著重突出沉積成因以及水動(dòng)力學(xué)能量，通常以巖石礦物顆粒骨架為主要分類依據(jù)。這使得油氣運(yùn)移規(guī)律的研究在這種分類方法中難以應(yīng)用。在一種常見的碳酸鹽巖孔隙類型分類方法中，“組構(gòu)選擇性”的關(guān)鍵概念將碳酸鹽巖孔隙分為組構(gòu)選擇性孔隙(包括粒間孔、粒內(nèi)孔、晶間孔和鑄?？椎?、非組構(gòu)選擇性孔隙(包括裂隙、溶溝、溶孔和溶孔等)和中間性組構(gòu)選擇性孔隙(角礫孔、生物潛穴和收縮孔)。這種方法在油氣勘探與開發(fā)過程中得到廣泛地應(yīng)用(choquettepw，praylc.geologicalnomenclatureandclassificationofporosityinsedimentarycarbonates.aapgbulletin，1970，54卷2期，207-250)。而考慮巖石組構(gòu)與孔隙空間之間的關(guān)系的分類方法主要將碳酸鹽巖孔隙空間劃分為粒間孔隙(包括顆粒間孔隙和晶間孔隙)和溶孔，其中粒間孔隙根據(jù)顆粒的大小分別劃分為粗顆粒間孔、中—粉晶顆粒間孔和粉—泥微晶孔，溶孔根據(jù)連通性劃分為孤立溶孔(包括鑄?？住⒘?nèi)孔、晶內(nèi)孔和化石內(nèi)孔)和接觸溶孔(裂縫、溶蝕擴(kuò)大的裂縫、溶孔、角礫孔)(luciafj.petrophysicalparameterestimatedfromvisualdescriptionofcarbonaterocks:afieldclassificationofporespace,journalofpetroleumtechnology，1983，35卷第3期，626-637)。

恒速壓汞實(shí)驗(yàn)通過極低的恒定速度向巖樣喉道及孔隙內(nèi)進(jìn)汞，保證了進(jìn)汞過程在準(zhǔn)靜態(tài)下進(jìn)行，該方法能將喉道和孔隙分開，實(shí)現(xiàn)喉道和孔隙數(shù)量和大小的準(zhǔn)確測(cè)量(李珊、孫衛(wèi)、王力等，恒速壓汞技術(shù)在儲(chǔ)層孔隙結(jié)構(gòu)研究中的應(yīng)用，斷塊油氣田，2013，20卷第4期，485-487)。近年來，隨著核磁共振在巖心分析中的廣泛應(yīng)用，通過巖心核磁共振t2譜分布描述其孔隙分布及孔隙結(jié)構(gòu)特征參數(shù)成為巖心物性分析的新方法。通過壓汞曲線和t2譜分布相態(tài)對(duì)比，可以得到二者之間的轉(zhuǎn)換系數(shù)(王勝，用核磁共振分析巖石孔隙結(jié)構(gòu)特征，新疆石油地質(zhì)，2009，30(6):768-770)。另外，專利《一種獲取碳酸鹽巖巖心孔洞特征參數(shù)的方法及裝置》(cn103325118a)提供了一種利用鉆取巖心圖像獲取全直徑巖心尺度上孔洞地質(zhì)參數(shù)的方法，該方法對(duì)碳酸鹽巖巖心表面孔洞進(jìn)行特征提取和宏觀微觀的分析，能夠?yàn)閮?chǔ)層分布情況預(yù)測(cè)提供良好技術(shù)支持。

以上介紹的孔隙類型分類方法中，考慮組構(gòu)的分類方法主要還是依靠肉眼識(shí)別和經(jīng)驗(yàn)判斷，未能建立對(duì)碳酸鹽巖孔隙類型進(jìn)行定量評(píng)價(jià)的標(biāo)準(zhǔn)體系。壓汞法主要依靠進(jìn)汞飽和度判斷孔徑分布，不能直觀地觀察各種孔隙類型分布特征。核磁方法成像速度慢，價(jià)格昂貴，難以用于批量的巖心掃描。而在碳酸鹽巖研究過程中涉及到圖像處理的相關(guān)技術(shù)，主要是用來定量計(jì)算孔隙和孔洞的參數(shù)，并沒有涉及裂縫的識(shí)別問題。

技術(shù)實(shí)現(xiàn)要素：

針對(duì)上述問題，本發(fā)明提出一種基于ct掃描圖像對(duì)碳酸鹽巖巖心孔隙進(jìn)行定量表征并判斷碳酸鹽巖孔隙類型的方法。該方法利用孔隙參數(shù)定量地表征全巖心孔隙類型(孔、洞、裂縫)的發(fā)育情況，擺脫了現(xiàn)有技術(shù)中依靠研究人員經(jīng)驗(yàn)判斷碳酸鹽巖孔隙類型的傳統(tǒng)模式，為批量判定儲(chǔ)層巖心孔隙類型提供了技術(shù)支持。

一種基于掃描圖像判斷碳酸鹽巖孔隙類型的方法，包括以下步驟：

s100，對(duì)碳酸鹽巖的巖心掃描圖像進(jìn)行包括彩色轉(zhuǎn)灰度和提高信噪比的預(yù)處理；

s200，對(duì)步驟s100處理后的圖像進(jìn)行分割，達(dá)到區(qū)分圖像中的孔隙區(qū)域與基質(zhì)區(qū)域的效果；

s300，對(duì)步驟s200分離出的孔隙區(qū)域提取孔隙參數(shù)，并識(shí)別孔隙區(qū)域中的裂縫和溶孔；

s400，根據(jù)步驟s300獲得的孔隙參數(shù)計(jì)算孔隙分類系數(shù)，根據(jù)孔隙分類系數(shù)的大小判斷巖心的孔隙類型，所述孔隙類型包括裂縫型、溶孔型和基質(zhì)型；

其中，所述孔隙分類系數(shù)是用于評(píng)價(jià)孔隙區(qū)域中裂縫型孔隙與溶孔型孔隙之間相對(duì)大小的參數(shù)，其值tc1為

式中，elc為橫切面裂縫連通系數(shù)；ell為縱切面裂縫發(fā)育系數(shù)；pcr為裂縫和孔洞相對(duì)發(fā)育程度；為巖樣孔隙度；且

式中，fl為橫截面中最大裂縫長度；cl為巖心長度；efl為縱截面中有效長度大于指定長度的裂縫長度總和；cr為巖心直徑；fp為裂縫面積占孔隙區(qū)域百分?jǐn)?shù)；hp為孔洞面積占孔隙區(qū)域百分?jǐn)?shù)。

根據(jù)本發(fā)明的實(shí)施例，上述步驟s100中，可以通過中值濾波、自適應(yīng)中值濾波、高斯濾波、低通濾波、高通濾波或維納濾波來提高圖像的信噪比。

根據(jù)本發(fā)明的實(shí)施例，上述步驟s200中，可以采用基于正態(tài)分布的孔隙度約束算法對(duì)圖像進(jìn)行分割，該算法包括以下步驟：

1)統(tǒng)計(jì)同一巖心不同ct截面處的灰度概率密度曲線，利用正太分布曲線對(duì)灰度概率密度曲線進(jìn)行擬合，求出正態(tài)分布曲線中的均值μ、方差σ，

式中，a為峰值伸展系數(shù)，b為基線偏移量，μ為正態(tài)分布均值，σ為正態(tài)分布方差；

2)在給定的方差距離范圍內(nèi)[λmin，λmax]令λ＝λ0，λ0為初始方差距離；

3)分別計(jì)算出該方差距離下的不同ct截面處所對(duì)應(yīng)的圖像分割閥值μi-λσi，i＝1,2,……n，n為掃描圖像個(gè)數(shù)，然后分別計(jì)算獲得二值圖像中黑色區(qū)域面積，統(tǒng)計(jì)出不同截面下的面孔率spi；

4)計(jì)算不同截面下的面孔率的平均值spave；

5)判斷面孔率的平均值spave與孔隙度φ的差值的絕對(duì)值是否小于預(yù)先給定的條件參數(shù)ε：

若|spave-φ|＜ε，條件中止，且

如果spave＜φ，則λ＝λmax，結(jié)束；

如果spave＞φ，則λ＝λmin，結(jié)束；

若|spave-φ|≥ε，則)，重復(fù)執(zhí)行步驟3)～5)，直到滿足條件中止；

6)基于確定的恒定系數(shù)λ，獲取不同圖像的分割閥值，閥值大小為μi-λσi，i＝1,2,……n，n為掃描圖像個(gè)數(shù)，最終根據(jù)同一測(cè)試巖心不同截面位置處的ct圖像區(qū)分出孔隙區(qū)域與基質(zhì)區(qū)域。

根據(jù)本發(fā)明的實(shí)施例，上述步驟s200中，可以采用以基質(zhì)區(qū)域?yàn)樯L點(diǎn)的種子生長算法對(duì)圖像進(jìn)行分割。

根據(jù)本發(fā)明的實(shí)施例，上述步驟s200中，可以采用孔隙區(qū)域邊緣提取算法對(duì)圖像進(jìn)行分割。

根據(jù)本發(fā)明的實(shí)施例，上述步驟s300中，提取孔隙參數(shù)的過程包括以下步驟：

1)對(duì)分離出的孔隙區(qū)域進(jìn)行形態(tài)學(xué)開閉運(yùn)算處理，選擇性地間斷或連接孔隙；

2)對(duì)經(jīng)過形態(tài)學(xué)處理后的孔隙區(qū)域作骨架細(xì)化處理，獲得孔隙區(qū)域的拓?fù)浣Y(jié)構(gòu)；

3)基于孔隙區(qū)域的拓?fù)浣Y(jié)構(gòu)計(jì)算孔隙區(qū)域的孔隙參數(shù)。

根據(jù)本發(fā)明的實(shí)施例，上述孔隙參數(shù)至少包括孔隙區(qū)域的面積、凸面積、比表面、有效長度、迂曲長度、等效寬度、最大內(nèi)切圓半徑、等效直徑中的一種。

根據(jù)本發(fā)明的實(shí)施例，上述步驟s300包括以下小步驟：

1)通過對(duì)ct掃描圖中似圓狀或者橢圓狀的溶孔孔隙、長條狀的裂縫型孔隙進(jìn)行初步篩選來建立特征向量，為后續(xù)的支持向量機(jī)訓(xùn)練提供訓(xùn)練樣本；

2)采用支持向量機(jī)進(jìn)行樣本分類訓(xùn)練，進(jìn)而自動(dòng)化地完成對(duì)孔隙區(qū)域中的裂縫和溶孔的分類識(shí)別。

根據(jù)本發(fā)明的實(shí)施例，efl為縱截面中有效長度至少大于5mm的裂縫長度總和。

進(jìn)一步地，若tc1>10，碳酸鹽巖孔隙為裂縫型；若1<tc1<10，碳酸鹽巖孔隙為溶孔型；若tc1<1，碳酸鹽巖孔隙為基質(zhì)型。

與現(xiàn)有技術(shù)相比，本發(fā)明的一個(gè)或多個(gè)實(shí)施例可以具有如下優(yōu)點(diǎn)：

本發(fā)明利用孔隙分類系數(shù)的概念來判定碳酸鹽巖儲(chǔ)層的孔隙類型，擺脫了現(xiàn)有技術(shù)中依靠專業(yè)人員經(jīng)驗(yàn)定性判斷碳酸鹽巖儲(chǔ)層類型的傳統(tǒng)模式，對(duì)儲(chǔ)層裂縫、溶孔發(fā)育程度做出科學(xué)的量化評(píng)價(jià)，為批量判定儲(chǔ)層巖心孔隙類型提供了有力的技術(shù)支持。此外，本發(fā)明在獲取孔隙分類系數(shù)的過程中，同時(shí)也能夠批量地獲得碳酸鹽巖巖心孔隙的分布情況、發(fā)育程度以及形狀特征等重要信息，對(duì)烴類運(yùn)移與油氣滲流規(guī)律研究提供了基礎(chǔ)參數(shù)。

本發(fā)明的其它特征和優(yōu)點(diǎn)將在隨后的說明書中闡述，并且，部分地從說明書中變得顯而易見，或者通過實(shí)施本發(fā)明而了解。本發(fā)明的目的和其他優(yōu)點(diǎn)可通過在說明書、權(quán)利要求書以及附圖中所特別指出的結(jié)構(gòu)來實(shí)現(xiàn)和獲得。

附圖說明

附圖用來提供對(duì)本發(fā)明的進(jìn)一步理解，并且構(gòu)成說明書的一部分，與本發(fā)明的實(shí)施例共同用于解釋本發(fā)明，并不構(gòu)成對(duì)本發(fā)明的限制。在附圖中：

圖1是本發(fā)明實(shí)施例中判斷碳酸鹽巖孔隙類型的方法流程圖；

圖2示例性地顯示了本發(fā)明實(shí)施例中對(duì)巖心進(jìn)行ct掃描獲取截圖的示意圖；

圖3(a)和3(b)示例性地顯示了本發(fā)明實(shí)施例中對(duì)ct掃描圖像進(jìn)行預(yù)處理的前后對(duì)比圖；

圖4示例性地顯示了本發(fā)明實(shí)施例中采用自適應(yīng)中值濾波法提高圖像信噪比的工作原理；

圖5示例性地顯示了本發(fā)明實(shí)施例中利用正態(tài)分布曲線對(duì)ct掃描圖像灰度概率密度進(jìn)行擬合的示意圖；

圖6(a)～6(c)示例性地顯示了本發(fā)明實(shí)施例中通過孔隙度約束算法分離出的孔隙區(qū)域與基質(zhì)區(qū)域；

圖7示例性地顯示了本發(fā)明實(shí)施例中種子生長算法中生長點(diǎn)的生長過程；

圖8(a)和8(b)示例性地顯示了本發(fā)明實(shí)施例中通過種子生長算法分離出的孔隙區(qū)域與基質(zhì)區(qū)域；

圖9(a)和9(b)示例性地顯示了本發(fā)明實(shí)施例中通過孔隙區(qū)域邊緣提取算法分離出的孔隙區(qū)域與基質(zhì)區(qū)域；

圖10(a)和10(b)示例性地顯示了本發(fā)明實(shí)施例中獲得的形態(tài)學(xué)圖像、骨架圖像以及方位角玫瑰圖；

圖11示例性地顯示了本發(fā)明實(shí)施例中獲得的形態(tài)學(xué)圖像和歐式距離圖；

圖12(a)和12(b)示例性地顯示了本發(fā)明實(shí)施例中兩種無量綱參數(shù)團(tuán)組合下的svm訓(xùn)練的分類測(cè)試結(jié)果；

圖13示例性地顯示了不同大小的孔隙分類系數(shù)所對(duì)應(yīng)的不同孔隙類型的巖心的參考圖像。

具體實(shí)施方式

為使本發(fā)明的目的、技術(shù)方案和優(yōu)點(diǎn)更加清楚，以下結(jié)合附圖和實(shí)施例對(duì)本發(fā)明作進(jìn)一步地詳細(xì)說明。應(yīng)當(dāng)說明的是，本發(fā)明中的示意性實(shí)施方式及其說明僅用于解釋本發(fā)明，并不作為對(duì)本發(fā)明的限定。

本發(fā)明提供的基于掃描圖像判斷碳酸鹽巖孔隙類型的方法，主要包括以下四大步驟。

s100，對(duì)碳酸鹽巖的巖心掃描圖像進(jìn)行包括彩色轉(zhuǎn)灰度和提高信噪比的預(yù)處理。

對(duì)巖心進(jìn)行ct圖像掃描，掃描方式可參見附圖2，然后對(duì)掃描后的圖像進(jìn)行預(yù)處理。預(yù)處理包括將彩色圖像轉(zhuǎn)為灰度圖像以及提高圖像的信噪比。具體地，使用photoshop軟件打開ct掃描圖像，利用磁性套索工具對(duì)巖心區(qū)域進(jìn)行勾勒，導(dǎo)出到新的空白文檔中保存為通用格式文件(如圖3所示)。將預(yù)處理后的圖像轉(zhuǎn)為灰度圖像，并對(duì)灰度圖像進(jìn)行有效的降噪處理，以便為后續(xù)的圖像分割提供基礎(chǔ)。具體實(shí)施方式可以如下：在matlab軟件中利用imread函數(shù)導(dǎo)入路徑下預(yù)處理后的ct圖像，將導(dǎo)入的8位無符號(hào)rgb圖像轉(zhuǎn)化為灰度圖像，數(shù)據(jù)體由三維變?yōu)槎S。在本實(shí)例中，優(yōu)選采用自適應(yīng)中值濾波的方法對(duì)圖像進(jìn)行降噪處理。自適應(yīng)中值濾波的工作原理是在中值濾波基礎(chǔ)上進(jìn)行選擇性優(yōu)選，對(duì)濾波鄰域的中心像素點(diǎn)進(jìn)行判斷：若被判斷的像素點(diǎn)為鄰域極值，則對(duì)中值濾波的統(tǒng)計(jì)排序進(jìn)行降噪處理；否則跳過此被判斷的像素點(diǎn)，進(jìn)行下一像素點(diǎn)的領(lǐng)域判斷(如圖4所示)。當(dāng)然，也可以采用例如中值濾波、高斯濾波、低通濾波、高通濾波或維納濾波等方式來提高圖像的信噪比，本發(fā)明對(duì)此不做限定。

s200，對(duì)步驟s100處理后的圖像進(jìn)行分割，達(dá)到區(qū)分圖中的孔隙區(qū)域與基質(zhì)區(qū)域的效果。

分析圖像中的孔隙區(qū)域，直到達(dá)到孔隙區(qū)域與基質(zhì)區(qū)域完全區(qū)分開來的效果。在本發(fā)明中，可以采用三種不同的圖像分割算法來區(qū)分圖像中的孔隙區(qū)域和基質(zhì)區(qū)域：

i)基于正態(tài)分布的孔隙度約束算法；

ii)以基質(zhì)區(qū)域?yàn)樯L點(diǎn)的種子生長算法；

iii)孔隙區(qū)域邊緣提取算法。

這三種方法既可以單獨(dú)實(shí)施，也可以組合實(shí)施。

基于正態(tài)分布的孔隙度約束算法是一種能夠針對(duì)同一測(cè)試巖心多張ct掃描圖像進(jìn)行批量閥值約束的算法。該算法的核心在于尋找距離灰度峰值μ的固定距離λ·σ位置處的灰度值作為圖像分割的閥值，因此確定恒定系數(shù)λ(通常λ的數(shù)值范圍在0～3之間)便成為了算法的關(guān)鍵，確定恒定系數(shù)λ的具體步驟如下：

1)統(tǒng)計(jì)同一巖心不同ct截面處的灰度概率密度曲線，利用正太分布曲線對(duì)灰度概率密度曲線進(jìn)行擬合，求出正態(tài)分布曲線中的均值μ、方差σ。

式中，a為峰值伸展系數(shù)，b為基線偏移量，μ為正態(tài)分布均值，σ為正態(tài)分布方差。

圖5示例性地顯示了一條擬合曲線，從中可以看出正態(tài)分布f(x)對(duì)ct掃描圖像具有良好的適應(yīng)性。在不同掃描圖像中，基于正態(tài)分布曲線擬合得到的峰值μ和方差σ均是不同，為了充分考慮每張圖像拍攝背景環(huán)境差異而造成的對(duì)比度、亮度等因素影響，每張圖片的閥值均取在相對(duì)于峰值μ的相同位置處的灰度值作為圖像分割的閥值，進(jìn)而區(qū)分孔隙空間和基質(zhì)區(qū)域。

2)在給定的方差距離范圍內(nèi)[λmin，λmax]，令λ為初始方差距離λ0。此處，λ0的取值范圍在0～3之間。

3)分別計(jì)算出該方差距離下的不同ct截面處所對(duì)應(yīng)的圖像分割閥值μi-λσi(i＝1,2,……n，n為掃描圖像個(gè)數(shù))，然后分別計(jì)算獲得二值圖像中黑色區(qū)域面積，統(tǒng)計(jì)出不同截面下的面孔率spi。

4)計(jì)算不同截面下的面孔率的平均值spave。

5)判斷面孔率的平均值spave與孔隙度φ的差值的絕對(duì)值是否小于預(yù)先給定的條件參數(shù)ε：

若|spave-φ|＜ε，條件中止，且

如果spave＜φ，則λ＝λmax，結(jié)束；

如果spave＞φ，則λ＝λmin，結(jié)束；

若|spave-φ|≥ε，則重復(fù)執(zhí)行步驟3)～5)，直到滿足條件中止。

通過以上步驟便可獲得能夠匹配巖心孔隙度的恒定系數(shù)λ，進(jìn)而獲取不同圖像的分割閥值，閥值大小具體為μi-λσi(i＝1,2,……n，n為掃描圖像個(gè)數(shù))，最終根據(jù)同一測(cè)試巖心不同截面位置處的ct圖像便可以高效地區(qū)分出孔隙區(qū)域與基質(zhì)區(qū)域。

圖6示例性地顯示了通過孔隙度約束算法分離出的孔隙區(qū)域與基質(zhì)區(qū)域。其中，圖6(a)是基于裂縫型ct掃描圖像分離出的孔隙區(qū)域與基質(zhì)區(qū)域；圖6(b)是基于孔洞型ct掃描圖像分離出的孔隙區(qū)域與基質(zhì)區(qū)域；圖6(c)是基于基質(zhì)型ct掃描圖像分離出的孔隙區(qū)域與基質(zhì)區(qū)域。

當(dāng)然，也可以采用另外一種分割方法：以基質(zhì)區(qū)域?yàn)樯L點(diǎn)的種子生長算法。以基質(zhì)區(qū)域?yàn)樯L點(diǎn)的種子生長算法是將預(yù)先通過人機(jī)交互的方式定義的生長點(diǎn)作為起始生長點(diǎn)，按照相似準(zhǔn)則要求進(jìn)行生長，預(yù)選生長點(diǎn)的生長過程可參見圖7(圖7中種子生長相似準(zhǔn)則設(shè)為1，在8鄰域中尋找)。在實(shí)施例中，種子生長點(diǎn)設(shè)置在發(fā)育普遍基質(zhì)像素上，相似生長閥值設(shè)置為40(最高不超過50)，然后對(duì)灰度圖像進(jìn)行孔隙區(qū)域分離，獲得的二值圖像中的黑色區(qū)域(即為0區(qū)域)為孔隙區(qū)域，白色區(qū)域則為巖石骨架區(qū)域。

圖8示例性地顯示了通過種子生長算法分離出的孔隙區(qū)域與基質(zhì)區(qū)域。其中，圖8(a)是基于裂縫型ct掃描圖像分離出的孔隙區(qū)域與基質(zhì)區(qū)域；圖8(b)是基于孔洞型ct掃描圖像分離出的孔隙區(qū)域與基質(zhì)區(qū)域。

而孔隙區(qū)域邊緣提取算法在此可以用來驗(yàn)證其他分割算法(例如上述兩種算法)的分割效果。例如，在本實(shí)施例中，可以使用matlab內(nèi)置的edge函數(shù)調(diào)用具有較好抗噪性和較高識(shí)別精度的canny算子，優(yōu)選地將參數(shù)調(diào)取灰度范圍設(shè)置為[40，70]，調(diào)用方差設(shè)置為1.25，對(duì)孔隙區(qū)域邊緣進(jìn)行提取。

圖9示例性地顯示了通過孔隙區(qū)域邊緣提取算法分離出的孔隙區(qū)域與基質(zhì)區(qū)域。其中，圖9(a)是基于裂縫型ct掃描圖像分離出的孔隙區(qū)域與基質(zhì)區(qū)域；圖9(b)是基于孔洞型ct掃描圖像分離出的孔隙區(qū)域與基質(zhì)區(qū)域。

s300，對(duì)步驟s200分離出的孔隙區(qū)域提取孔隙參數(shù)，并識(shí)別孔隙區(qū)域中的裂縫和溶孔。

基于步驟s200分離出的孔隙區(qū)域圖像提取孔隙參數(shù)，以便為后續(xù)的孔隙類型分布、發(fā)育程度等信息建立無量綱參數(shù)團(tuán)進(jìn)而定量表征提供基礎(chǔ)。該步驟主要分為區(qū)域特征提取和裂縫、溶孔識(shí)別兩大階段。其中：

(一)區(qū)域特征提取(也即孔隙參數(shù)提取)主要包括形態(tài)學(xué)處理、骨骼細(xì)化處理和圖像特征計(jì)算三個(gè)小步驟。

1)形態(tài)學(xué)處理：對(duì)分割出的孔隙區(qū)域進(jìn)行形態(tài)學(xué)開閉運(yùn)算處理，選擇性地間斷或連接孔隙，以對(duì)分割出的孔隙區(qū)域完成形態(tài)學(xué)上的校正和平滑處理，以及去除無效噪聲點(diǎn)。在本發(fā)明中，優(yōu)選采用開閉運(yùn)算聯(lián)用的方式對(duì)孔隙區(qū)域進(jìn)行校正。具體實(shí)施方式為matlab中使用strel函數(shù)建立方形結(jié)構(gòu)元素，設(shè)置正方形邊長像素為2，分別使用imopen函數(shù)和imclose函數(shù)對(duì)二值圖像進(jìn)行形態(tài)學(xué)操作，即先使用開運(yùn)算對(duì)分割二值圖像進(jìn)行形態(tài)學(xué)處理，然后再用閉運(yùn)算進(jìn)行處理，選擇性地連接或斷開圖像中孔隙區(qū)域，達(dá)到孔隙區(qū)域的形態(tài)平滑處理的效果。

2)骨架細(xì)化處理：對(duì)分割出的孔隙區(qū)域的中軸骨架、端點(diǎn)、內(nèi)部孔洞等重要拓?fù)浣Y(jié)構(gòu)進(jìn)行剝離。在本實(shí)施例中，具體實(shí)施過程如下：在3×3的鄰域范圍內(nèi)，將前景色(多為黑色)設(shè)為1，背景色(多為白色)設(shè)為0，按照以下規(guī)則進(jìn)行選擇性的刪除：

①2≤nz(p)≤6，nz(p)表示p點(diǎn)周圍8個(gè)像素點(diǎn)(p上、p下、p左、p右、p左下、p左上、p右上、p右下)中1的數(shù)目；

②z0(p)＝1，z0(p)表示p周圍像素點(diǎn)中相鄰像素相異點(diǎn)的數(shù)目；

③p上×p左×p右＝0或z0(p)≠1；

④p上×p左×p下＝0或z0(p左)≠1。

同時(shí)滿足以上規(guī)則時(shí)，中心像素點(diǎn)p會(huì)由黑色變?yōu)榘咨?，也即刪除前景色。

通過剔除連通域的冗余邊緣像素而剝離出連通域中的骨架結(jié)構(gòu)，能夠著重反應(yīng)出骨架的形狀特征。圖10示例性地顯示了由ct掃描圖像獲得的形態(tài)學(xué)圖像和骨架圖像。其中，圖10(a)是基于裂縫型ct掃描圖像獲得的形態(tài)學(xué)圖像和骨架圖像；圖10(b)是基于孔洞型ct掃描圖像獲得的形態(tài)學(xué)圖像和骨架圖像。

3)圖像特征計(jì)算：主要是孔隙參數(shù)的計(jì)算提取。通過該步驟可以獲得孔隙區(qū)域的面積和長度兩大類型的參數(shù)，可以包括面積、凸面積、比表面、有效長度、迂曲長度、等效寬度、最大內(nèi)切圓半徑、等效直徑、離心度和方位等10余項(xiàng)參數(shù)。以下是這些參數(shù)的具體計(jì)算方法：

(1)孔隙區(qū)域的面積和凸面積

區(qū)域參數(shù)的計(jì)算需要先對(duì)孔隙區(qū)域進(jìn)行標(biāo)注。具體地，可以使用bwlabel函數(shù)對(duì)孔隙區(qū)域做序號(hào)標(biāo)注，然后依據(jù)標(biāo)注的序號(hào)使用regionprops函數(shù)調(diào)用面積算子area求得區(qū)域孔隙像素面積；凸面積的計(jì)算可以使用regionprops函數(shù)調(diào)用凸面積算子convexarea進(jìn)行求解，求得的面積和凸面積均為圖像像素面積。

(2)孔隙區(qū)域方向

孔隙區(qū)域方向需要借助骨架細(xì)化算法所獲得孔隙區(qū)域的拓?fù)浣Y(jié)構(gòu)來獲取，即在孔隙拓?fù)浣Y(jié)構(gòu)的基礎(chǔ)上進(jìn)行方向判斷。具體實(shí)施方法為：對(duì)孔隙區(qū)域進(jìn)行bwlabel標(biāo)注，然后調(diào)用bwmorph函數(shù)中的thicken算子對(duì)單一孔隙區(qū)域進(jìn)行侵蝕處理，在獲取侵蝕后的骨架后使用regionprops函數(shù)中的orientation算子獲取孔隙區(qū)域的方向展布。圖10是示例性地顯示了基于骨架圖像拓?fù)浣Y(jié)構(gòu)所獲得的表征孔隙區(qū)域方向的方位角玫瑰圖。其中，圖10(a)是裂縫型的方位角玫瑰圖；圖10(b)是孔洞型的方位角玫瑰圖。

(3)有效長度和迂曲長度

孔隙的有效長度用于評(píng)價(jià)孔隙區(qū)域的直線距離，它需要通過孔隙區(qū)域所包羅的矩形長寬來計(jì)算，具體計(jì)算方法如下：

式中，l、w為區(qū)域外接矩形的長和寬，θ為孔隙方位角。

具體實(shí)施時(shí)，上式中矩形的長和寬可以利用regionprops函數(shù)中的boundingbox算子來獲取，然后結(jié)合已知的方向參數(shù)(孔隙方位角)，利用上式中g(shù)ivens矩陣做坐標(biāo)旋轉(zhuǎn)即可獲得。

迂曲長度可以基于由骨架細(xì)化算法獲得的孔隙區(qū)域中軸骨架，通過下式計(jì)算獲得：

式中，xi、yi分別表示位于中軸骨架中的x方向、y方向的坐標(biāo)，i表示位于中軸骨架上的不同位置的像素點(diǎn)。

根據(jù)有效長度和迂曲長度可以進(jìn)一步獲取迂曲度等參數(shù)，從而為孔隙區(qū)域的方向展布、規(guī)則程度做出有效判斷。

(4)等效寬度和等效直徑

等效寬度根據(jù)孔隙面積和迂曲長度聯(lián)立計(jì)算確定，具體計(jì)算方法如下：

等效寬度＝孔隙面積/迂曲長度

利用等效寬度可以將連通孔隙區(qū)域中隨位置變化的寬度數(shù)值均一化。

進(jìn)一步地，在孔隙面積的基礎(chǔ)上，還可以將不規(guī)則孔隙區(qū)域模擬成等面積圓，其直徑稱為等效直徑，具體計(jì)算方法如下：

利用不規(guī)則形狀孔隙區(qū)域的等效直徑能夠分辨不規(guī)則孔隙凈面積的大小。

(5)最大內(nèi)切圓半徑

最大內(nèi)切圓半徑是孔隙評(píng)價(jià)中的重要參數(shù)，可以通過歐氏距離圖獲得孔隙區(qū)域中不同位置處的內(nèi)切圓半徑數(shù)值，然后取其中的最大值作為最大內(nèi)切圓半徑。具體實(shí)施時(shí)，可以使用bwdist函數(shù)對(duì)二值圖像的背景色區(qū)域進(jìn)行歐式距離圖的轉(zhuǎn)化，然后選取其中的最大值作為最大內(nèi)切圓半徑。圖11示例性地顯示了由ct掃描圖像獲得的形態(tài)學(xué)圖像和歐氏距離圖。

(6)比表面

孔隙區(qū)域的比表面是孔隙區(qū)域非常重要參數(shù)，它可以通過對(duì)孔隙區(qū)域總周長的計(jì)算獲得，具體計(jì)算方法如下：

sv＝4/π·ba

式中，ba為圖像單位面積內(nèi)中邊界長度，sv為單位體積內(nèi)的表面積。

具體實(shí)施時(shí)，可以首先對(duì)形態(tài)學(xué)處理后的二值圖像做sobel算子邊緣提取，對(duì)邊緣提取后的邊緣像素?fù)Q算后的長度進(jìn)行統(tǒng)計(jì)加和運(yùn)算，再乘以2/π，其數(shù)值大小即為圖像中所有孔隙區(qū)域的比表面。

(二)裂縫、溶孔識(shí)別主要包括建立特征向量、svm(支持向量機(jī))分類識(shí)別訓(xùn)練兩個(gè)小步驟。

1)建立特征向量。該步驟通過對(duì)似圓狀或者橢圓狀的溶孔孔隙、長條狀的裂縫型孔隙(即具有典型特征的孔隙)進(jìn)行篩選來建立特征向量，以為后續(xù)的支持向量機(jī)訓(xùn)練提供訓(xùn)練樣本。

首先對(duì)巖心圖像特征參數(shù)進(jìn)行提取，選取經(jīng)過圖像分割、形態(tài)學(xué)處理后的ct掃描圖像中具有典型特征的長條狀裂縫和似圓狀孔洞參數(shù)為特征向量。例如，表1是166h巖心經(jīng)過區(qū)域特征參數(shù)提取后孔隙類型的初步分類結(jié)果，表中各個(gè)孔隙連通域通過初步判別確定其孔隙類型，其中類型1為裂縫型，類型2為基質(zhì)型，類型3為溶孔型。以相同的方法對(duì)其余編號(hào)巖心進(jìn)行區(qū)域特征參數(shù)提取和孔隙類型初步標(biāo)定等步驟的處理后，在大量的孔隙連通域中分別挑選出具有典型特征的長條狀裂縫型和似圓狀孔洞孔洞型連通域特征參數(shù)建立標(biāo)準(zhǔn)特征向量，如表2所示為篩選出的標(biāo)準(zhǔn)特征向量。最后利用篩選出的標(biāo)準(zhǔn)特征參數(shù)計(jì)算無因次變量組成特征向量，無因次參數(shù)團(tuán)計(jì)算方法具體如下：

(1)長寬比和縱橫比

長寬比和縱橫比均為評(píng)價(jià)孔隙空間條狀形態(tài)的發(fā)育程度的重要參數(shù)，區(qū)別在于計(jì)算中分別使用的參數(shù)不同，長寬比采用計(jì)算參數(shù)為迂曲長度，具體計(jì)算方法：

長寬比＝迂曲長度÷等效寬度

與之相對(duì)應(yīng)的為縱橫比，計(jì)算參數(shù)則變?yōu)橛行чL度，具體計(jì)算方法為：

縱橫比＝有效長度÷等效寬度

(2)迂曲度

迂曲長度主要評(píng)價(jià)孔隙空間中中軸骨架的彎折程度，其涉及參數(shù)主要為迂曲長度和有效長度兩大參數(shù)，具體計(jì)算方法為：

迂曲度＝迂曲長度÷有效長度

迂曲度的數(shù)值變化均在1以上，數(shù)值越大反應(yīng)彎折程度越大。

(3)離心度

離心度是評(píng)價(jià)孔隙空間形狀特征的重要參數(shù)，計(jì)算方法為將孔隙空間近似成橢圓狀，分別計(jì)算得到近似橢圓長軸a和短軸b的長度大小，離心度具體計(jì)算：

離心度的數(shù)值范圍介于0～1之間，數(shù)值越接近于1，孔隙空間形狀越扁；相反，數(shù)值越接近于0的時(shí)候，孔隙空間形狀越接近于圓形。通常，在具體實(shí)施過程中可以通過調(diào)用matlab內(nèi)編譯函數(shù)regionprops函數(shù)中的eccentricity算子計(jì)算得到。

(4)形狀因子

形狀因子主要用于描述孔隙邊緣的平滑和規(guī)則性，其數(shù)值大小用來表征孔隙形狀均勻規(guī)則的程度。由于形狀因子可以對(duì)孔隙的形態(tài)作定量化的表征，因此它是后續(xù)判斷裂縫、溶孔類型的重要參數(shù)。形狀因子的計(jì)算方法如下：

形狀因子＝4π·面積/周長²

(5)實(shí)心度

依據(jù)區(qū)域面積和凸面積的大小可以進(jìn)一步確定區(qū)域孔隙實(shí)心度。區(qū)域孔隙實(shí)心度的具體計(jì)算方法如下：

實(shí)心度＝面積/凸面積

實(shí)心度的大小能夠反應(yīng)形狀規(guī)則程度、內(nèi)部孔洞發(fā)育情況，為后續(xù)的裂縫、溶孔識(shí)別提供重要參數(shù)。

通過以上無量綱參數(shù)團(tuán)的計(jì)算后，可生成由無量綱參數(shù)構(gòu)成的標(biāo)準(zhǔn)特征向量，為后續(xù)的svm(支持向量機(jī))的樣本學(xué)習(xí)提供基礎(chǔ)，無量綱標(biāo)準(zhǔn)特征向量具體可見表3。

表1166h巖心不同連通域的特征參數(shù)及孔隙分類結(jié)果

表2裂縫型和孔洞型連通域特征向量

表3無因次特征向量

2)svm(支持向量機(jī))分類識(shí)別訓(xùn)練。該步驟使用支持向量機(jī)對(duì)孔隙類型進(jìn)行自動(dòng)化的分類篩選，從而為后續(xù)孔隙的分割提取、類型判斷提供極大便利。在本發(fā)明的實(shí)施例中，優(yōu)選采用兩兩參數(shù)組合的方式對(duì)裂縫、溶孔類型分別作訓(xùn)練，依次評(píng)價(jià)分類效果。具體地，分別優(yōu)選出形狀因子&離心度、實(shí)心度&離心度的分類方法，具體分類效果可見圖12。

圖12示例性地顯示了在兩種無量綱參數(shù)團(tuán)組合下的svm訓(xùn)練的分類測(cè)試結(jié)果。其中，圖12(a)為形狀因子&離心度的參數(shù)組合下svm針對(duì)溶孔型孔隙和裂縫型孔隙的分類訓(xùn)練效果圖；圖12(b)為實(shí)心度&離心度的參數(shù)組合下svm對(duì)溶孔型孔隙和裂縫型孔隙的分類訓(xùn)練效果圖。從圖12中可以獲知，“+”代表溶孔型孔隙的標(biāo)準(zhǔn)特征向量，“*”代表裂縫型孔隙的標(biāo)準(zhǔn)特征向量，在形狀因子&離心度或者實(shí)心度&離心度這兩種參數(shù)組合下的類型分界線均具有良好的分類效果。此外，“☆”代表著待測(cè)試的溶孔型孔隙類型，“◇”代表著待測(cè)試的裂縫型孔隙，通過svm分類訓(xùn)練可發(fā)現(xiàn)，“☆”和“◇”均落在了與之對(duì)應(yīng)的孔隙類型的范圍內(nèi)，其反映出通過形狀因子&離心度或者實(shí)心度&離心度的參數(shù)組合下的svm分類訓(xùn)練是可靠的。

s400，根據(jù)步驟s300獲得的孔隙參數(shù)計(jì)算孔隙分類系數(shù)，根據(jù)孔隙分類系數(shù)的大小判斷巖心的孔隙類型，所述孔隙類型包括裂縫型、溶孔型和基質(zhì)型。

應(yīng)當(dāng)說明的是，本發(fā)明中著重側(cè)重評(píng)價(jià)的是裂縫型孔隙、溶蝕型孔隙以及孔隙基質(zhì)型碳酸鹽巖之間的分類，對(duì)溶蝕型孔隙類型未做進(jìn)一步的細(xì)分評(píng)價(jià)。在此，以“溶孔”一詞作為溶蝕型孔隙類型的典型。

在獲得孔隙參數(shù)的基礎(chǔ)上，通過計(jì)算孔隙分類系數(shù)來將單一孔隙類型的判斷上升到巖心孔隙類型的判別。在本發(fā)明中，碳酸鹽巖儲(chǔ)層孔隙分類系數(shù)是用以評(píng)價(jià)孔隙區(qū)域中裂縫型孔隙與溶孔型孔隙之間相對(duì)大小的參數(shù)，其值tc1為：

式中，elc為橫切面裂縫連通系數(shù)；ell為縱切面裂縫發(fā)育系數(shù)；pcr為裂縫和孔洞相對(duì)發(fā)育程度；為巖樣孔隙度；且

式中，fl為橫截面中最大裂縫長度；cl為巖心長度；efl為縱截面中有效長度大于指定長度的裂縫長度總和；cr為巖心直徑；fp為裂縫面積占孔隙區(qū)域百分?jǐn)?shù)；hp為孔洞面積占孔隙區(qū)域百分?jǐn)?shù)。

由于計(jì)算過程中參數(shù)較多，在具體實(shí)施時(shí)可以通過編寫函數(shù)typeporeexport文件將先關(guān)參數(shù)導(dǎo)出，導(dǎo)出的參數(shù)直接用于計(jì)算孔隙分類系數(shù)。在本實(shí)施例中，efl設(shè)定為縱截面中有效長度大于5mm的裂縫長度總和，那么根據(jù)上式計(jì)算得到的孔隙分類系數(shù)tc1：

若tc1>10，碳酸鹽巖孔隙為裂縫型；

若1<tc1<10，碳酸鹽巖孔隙為溶孔型；

若tc1<1，碳酸鹽巖孔隙為基質(zhì)型。

圖13示例性地顯示了不同大小的孔隙分類系數(shù)所對(duì)應(yīng)的不同孔隙類型的巖心的參考圖像。

以上所述，僅為本發(fā)明的具體實(shí)施案例，本發(fā)明的保護(hù)范圍并不局限于此，任何熟悉本技術(shù)的技術(shù)人員在本發(fā)明所述的技術(shù)規(guī)范內(nèi)，對(duì)本發(fā)明的修改或替換，都應(yīng)在本發(fā)明的保護(hù)范圍之內(nèi)。