基於三維FIB‑SEM圖像的頁巖孔隙定量表徵方法與流程

2023-10-22 05:16:41 2

本發明涉及石油天然氣地質勘探開發領域，特別涉及一種在頁巖油氣勘探中，採用fib-sem(聚焦離子束-掃描電鏡)識別頁巖內部三維孔隙的類型和定量表徵其結構信息。

背景技術：

頁巖的儲油和儲氣能力與頁巖中的孔隙度和有效孔隙度等因素息息相關，而孔隙類型和大小的定量描述更是直接決定油和氣在孔隙內的賦存狀態、運移機理和能力。同時，頁巖在沉積演化過程中，孔隙的演變始終是關注的焦點之一，是我們研究頁巖油氣沉積演化史的重點。因此對頁巖儲層的孔隙研究與刻畫有意義重大。

頁巖內的孔隙按形成的機理和位置，可以簡單地分為有機孔隙和無機孔隙。有機孔隙是有機質成熟過程中在其內部形成的孔隙，分布在有機質周圍和內部，而無機孔隙是由原始的孔隙演變而來或後期在無機礦物內部或礦物顆粒之間形成的，分布在無機礦物顆粒間或礦物內部。還有部分孔隙介於有機質和無機礦物顆粒之間。

目前對於頁巖儲層孔隙的分析已開展了大量研究，其中主要的測試技術手段有壓汞法和氣體吸附法，但這兩種測試方法的結果展示的都是頁巖整體的孔隙特徵，不能將頁巖中有機質孔隙和無機礦物孔隙進行區分。

而近幾年普遍採用的三維掃描技術有：fib-sem技術、納米ct技術。通過這些技術，可以清晰看到孔隙的三維空間分布。如果想要達到區分有機孔隙和無機孔隙並分別定量描述，往往是通過分割出孔隙後，人工根據經驗區分有機孔隙和無機孔隙。這種方法操作簡單，但是十分耗時，難以準確，推廣困難。

另外，通過掃描電鏡技術，也可以先獲取真實頁巖巖心在不同平面上的二維孔隙掃描圖像，然後將圖像分割轉換得到頁巖巖心二值圖像；並利用模擬的方法構建有機孔隙數字巖心，並通過疊加，建立有機、無機孔隙的數字巖心。這種通過二維數據去模擬三維空間的方法，雖然使用的機時較少，但是不能真正地反應三維空間有機孔隙和無機孔隙的分布，從而對評價頁巖儲油儲氣能力、孔隙內的油和氣的賦存狀態、運移機理和能力帶來局限性。

技術實現要素：

針對上述問題，本發明提供一種基於fib-sem圖像的頁巖有機孔隙和無機孔隙的識別和分析新方法。

本發明是通過以下技術方案實現的：

基於三維fib-sem圖像的頁巖孔隙定量表徵方法，所述方法包括：

圖像採集：利用fib-sem對頁巖中的無機質、有機質和孔隙進行微米和納米級尺度的高精度的三維成像，得到三維fib-sem圖像；

圖像預處理：將三維fib-sem圖像的各二維圖像進行對齊、對圖像亮度不均勻性進行校正、濾波去除圖像的噪聲；

圖像分割：將預處理後的圖像進行分割，把圖像劃分孔隙區、有機質區、無機質區；

孔隙分析：對圖像孔隙區中各孔隙進行標註，根據各孔隙周圍有機質區、無機質區的分布情況將孔隙分為有機孔隙、無機孔隙、邊界孔隙；

定量表徵參數：分別對有機孔隙、無機孔隙、邊界孔隙進行分析，得到三維孔隙的各種表徵參數值。

進一步地，所述孔隙分析具體為，

對三維圖像的各xy方向二維圖像進行孔隙連通性分析、提取各孔隙邊緣；

一個區域中任意兩個像素能夠連通，則認為這是一個孔隙；

二維連通性分析依據為，像素p和q，q在像素p同一層周圍的4或8個像素，則像素p和q是連通的；

對各二維圖像中孔隙提取邊緣；

對邊緣進行z方向連通性分析：三維連通性分析依據為，相鄰兩層面的孔隙邊緣像素p1和q1，如果q1在像素p1上下相鄰層周圍的2或10或18個像素，則像素p1和q1在z方向是連通的；

對各孔隙進行標註：對每個三維孔隙進行標註；

計算孔隙影響因子n：

其中，num_o為每個孔隙周圍有機質像素的數目；num_i為每個孔隙周圍無機質像素的數目；其中0≤n≤1；

孔隙類型判斷：

0≤n≤k1，孔隙為有機質孔隙；k1＜n＜k2，孔隙為邊界孔隙；k2≤n≤1，孔隙為無機質孔隙；

其中0≤k1＜k2≤1；

進一步地，所述對圖像亮度不均勻性進行校正具體為：

採用徑向基核函數各xy二維圖像進行灰度校正，然後將fib-sem三維圖像進行旋轉，對zx和zy進行二維圖像灰度校正，從而實現三維圖像的灰度校正。

進一步地，所述將三維fib-sem圖像的各二維圖像進行對齊具體為：通過變化相鄰兩張二維圖像的相對位置，計算其灰度差的平方和，得到相鄰圖像的差異程度，然後選取差異程度最小相對位置為對齊位置；依次循環，對所有二維圖像進行配準，得到位置對齊後的三維圖像。

進一步地，所述濾波去除圖像的噪聲具體為：

利用三維濾波方法濾除灰度圖像的噪聲；

所述三維濾波方法包括均值濾波、中值濾波、維納濾波、高斯濾波。

進一步地，所述圖像分割具體為：

其中，f(x,y,z)為原圖像的灰度值；g(x,y,z)為分割後圖像的灰度值；

t1、t2為分割域值；0≤t1＜t2≤255；式中a、b和c是不同的介於[0,255]灰度值。

進一步地，所述圖像分割還包括異常點剔除，將分割後孔隙中小於3像素的點剔除。

進一步地，所述表徵參數包括孔隙的體積、表面積、中心位置、長度、寬度、高度、傾角、直徑。

本發明的有益技術效果：

(1)本發明在頁巖油氣勘探中，採用fib-sem識別孔隙的類型，判斷孔隙是無機孔隙或有機孔隙，並對其結構進行分析表徵，解決現有技術中沒有方法實現對頁巖孔隙在三維空間有效分類和定量表徵的問題，為後續評價頁巖儲油儲氣能力、孔隙內的油和氣的賦存狀態、運移機理和能力提供基礎。

(2)該方法具有普適性，可用於各種頁巖的分析。

(3)本發明所用方法簡單、可用於同時包含有機質、無機礦物和孔隙的fib-sem圖像中，準確率高達100％。

附圖說明

圖1為本發明基本流程圖；

圖2為本發明實施例頁巖的fib-sem三維圖像；

圖3a為本發明實施例頁巖的配準三維圖像；

圖3b為本發明實施例頁巖的灰度校正三維圖像；

圖3c為本發明實施例頁巖的濾波平滑三維圖像；

圖4a為本發明實施例頁巖的分割三維圖像；

圖4b為本發明實施例頁巖的孔隙三維分布圖；

圖4c為本發明實施例頁巖的有機質三維分布圖；

圖5a為本發明實施例頁巖的有機孔隙三維分布圖；

圖5b為本發明實施例頁巖的無機孔隙三維分布圖；

圖5c為本發明實施例頁巖的邊界孔隙三維分布圖；

圖6a為本發明實施例頁巖的孔隙的直徑——數量柱狀分布圖；

圖6b為本發明實施例頁巖的孔隙的直徑——體積柱狀分布圖；

圖6c為本發明實施例頁巖的有機孔隙、無機孔隙和邊界孔隙的直徑——數量柱狀分布圖；

圖6d為本發明實施例頁巖的有機孔隙、無機孔隙和邊界孔隙的直徑——體積柱狀分布圖。

具體實施方式

為了使本發明的目的、技術方案及優點更加清楚明白，以下結合附圖及實施例，對本發明進行進一步詳細描述。應當理解，此處所描述的具體實施例僅僅用於解釋本發明，並不用於限定本發明。

相反，本發明涵蓋任何由權利要求定義的在本發明的精髓和範圍上做的替代、修改、等效方法以及方案。進一步，為了使公眾對本發明有更好的了解，在下文對本發明的細節描述中，詳盡描述了一些特定的細節部分。對本領域技術人員來說沒有這些細節部分的描述也可以完全理解本發明。

實施例1

本實施例針對頁巖樣品進行孔隙分析，包括如下步驟：

1.fib-sem獲取圖像

fib-sem一種集掃描電鏡(sem)的成像功能和聚焦離子束(fib)的切割功能為一體的雙束系統，主要用於對巖石中的礦物、有機質和孔隙進行微米和納米級尺度的高精度的三維成像。其成像原理為：樣品表面與離子束保持垂直；使用fib離子束對樣品表面進行刻蝕，露出觀察面，使用sem電子束對觀察面進行掃描，收集二次電子或背散射電子得到高解析度成像；然後根據剝蝕厚度的要求，設置離子束能量參數對觀察面進行剝蝕，剝蝕之後使用電子束對新的觀察面進行成像，重複這個步驟直到成像完成。最終得到一系列圖像。

fib-sem得到的圖像就是一系列的二維圖像，三維圖像是通過各個二維圖像在空間不同高度下的堆砌。

通過fib-sem獲取的圖像大小為：1188*1192，每個像素代表10nm。即掃描區域大小長為11.88um，寬為11.92um。切割的厚度為10nm，共切割了842張，即切割的物體高為8.42um。圖2為本發明實施例頁巖的fib-sem三維圖像。對儀器本身型號沒有要求。在此本發明使用的是蔡司生產的crossbeam540聚焦離子束-掃描電子顯微鏡。本發明使用二次電子信號成像，也可以是背散射電子成像。

2.圖像預處理

由於fib-sem獲得的三維圖像是分別採集的，因此每張二維圖像與下一張圖像之間有時會有一定的偏移。因此在對fib-sem的三維圖像分析之前，需要對三維圖像的位置進行對齊。最簡單的辦法是：根據圖像，人工拖拽將相鄰圖像對齊，然而fib-sem圖像一般有幾百或上千張圖像，人工工作量很大。另一種方法是自動對齊方法，由於fib-sem切片厚度較薄(一般只有幾納米至幾十納米)，因此相鄰兩張二維圖像具有一定的近似性，可對相鄰兩張二維圖像進行圖像配準。由於相鄰二維圖像只存在上下左右的偏移，因此可以通過變化相鄰兩張圖像的相對位置，計算其差異程度，然後選取差異程度最小的狀態，則認為該相對位置是校正好的位置。依次循環，採用上述方法對所有圖像進行配準，得到位置校正對齊後的三維圖像。採用兩種方法相結合的方式，可以將圖像很好地對齊，最後對校正後的圖像進行剪裁，去掉校準過程中的多餘部分，使得三維圖像保持規範的長方體。

其次，由於fib-sem掃描樣品觀察面前方空間有時不夠開闊，有些電子信號被前方遮擋，有時會出現觀察面底部信號比上部弱，最終掃描圖像的亮度由上至下逐步變暗。這種亮度的不均勻性，對圖像分割會產生很大幹擾。因此需要對圖像的灰度進行校正。對此，本發明採用的方法是徑向基核函數(radialbasisfunction，簡稱rbf)，它可以看作是一個高維空間中的曲面擬合(逼近)問題，學習是為了在多維空間中尋找一個能夠最佳匹配訓練數據的曲面，用這個訓練的曲面來處理(比如分類、回歸)後續新的數據。rbf的本質思想是應用遞歸技術反向傳播學習算法，這種技術在統計學中被稱為隨機逼近。rbf裡的basisfunction(徑向基函數裡的基函數)是在神經網絡的隱單元裡提供了一個函數集，該函數集在輸入模式(向量)擴展至隱空間時，為其構建了一個任意的「基」。這個函數集中的函數就被稱為徑向基函數。但是，rbf灰度校正針對的是二維圖像，若z方向需要進行灰度校正，只需將fib-sem三維圖像進行旋轉，然後進行灰度校正，再旋轉回來即可。

最後，在保證原始圖像所有重要細節信息特徵情況下，將fib-sem得到的帶有一定噪聲的灰度圖像進行濾波處理。一般三維濾波方法有均值濾波、中值濾波、維納濾波、高斯濾波等，本發明採用meanshift(均值偏移)濾波函數，meanshift濾波是一種有效的統計迭代算法，使得每一個點漂移到密度函數的局部極大值點。它能夠在實現圖像的自適應濾波過程中，有效地濾除複雜背景下噪音的幹擾，同時不損失目標信息。它的一個重要特徵是，在迭代的過程中，圖像結構沒有發生改變。因此，meanshift濾波具有非常好的邊緣保持性，在圖像降噪的同時能夠最大程度地保留圖像的細節特徵。

本實施例首先進行圖像配準，先採用人工方法，粗略地將圖像進行配準；然後再進行自動配準；圖3a為本發明實施例頁巖的配準三維圖像。配準裁剪多餘部分後物體大小為：11.00*11.00*8.42um3。其次，進行灰度校正，圖3b為本發明實施例頁巖的灰度校正三維圖像，可以看出圖像有了很好的均一性。最後採用meanshift濾波函數進行濾波，圖像子塊為3*3*3，圖3c為本發明實施例頁巖的濾波平滑三維圖像。可以看出通過以上三個步驟，圖像的錯位，灰度都得到很好的校正，同時三維圖像在去除噪聲的同時很好地保留圖像的細節特徵。

3.圖像分割

所謂圖像分割指的是根據灰度、顏色等特徵把圖像劃分成若干個互不交迭的區域，並使這些特徵在同一區域內呈現出相似性，而在不同區域間呈現出明顯的差異性。常用的分割方法有：閾值分割、區域生長法、分水嶺分割等。其中閾值分割以其實現簡單、計算量小、便於並行、速度快等優勢，成為三維圖像分割中最常用的一種方法。閾值分割方法的本質是輸入圖像f到輸出圖像g的如下變換：

其中，t1、t2為閾值。閾值確定後，將閾值與像素點的灰度值逐個進行比較，然後給出每個閾值區間內對應的圖像區域。通過閾值分割，我們可以初步將fib-sem掃描電鏡下頁巖樣品的孔隙、有機質和無機礦物區分開。對於裂縫，在此我們認為是一種特殊的孔隙，不進行單獨區分。

但是由於閾值分割不考慮孔隙真實存在的可能性，因此，需要將異常點剔除，例如：分割出小於3像素的點是否為孔隙，在此我們認為此情況下不是孔隙，需要將其剔除。然後重新得到分割結果。

對圖3c採用閾值分割的方法，分割得到孔隙、有機質和無機礦物的分布圖。在此，本發明實施例在公式(1)中的閾值t1、t2分別為：66、98。即孔隙的灰度值範圍：0～66，有機質灰度值範圍：67～98，無機質灰度值範圍：99～255。然後將異常點小於3個像素的孔隙剔除。圖4a為本發明實施例頁巖的分割三維圖像，其中黑色部分代表孔隙，灰色部分代表有機質，白色部分代表無機礦物。圖4b為本發明實施例頁巖的孔隙三維分布圖，圖4c為本發明實施例頁巖的有機質三維分布圖。其中灰色部分為分別為孔隙和有機質，白色部分為其它物質。為了凸顯三維可視化立體效果，在圖中加入了光照和陰影設計，其中黑色部分為陰影。

4.二維孔隙分析

分割後得到孔隙，首先需要定義什麼是一個孔隙，以及什麼是有機孔隙，什麼是無機孔隙。

在此，我們認為不是每一個單獨的像素都定義為一個孔隙，而是一個區域中任意兩個像素能夠連通，則認為這是一個孔隙，這與數學上的連通定義相一致，簡單來說：像素p和q，如果q在像素p同一層周圍的8個像素，或者上下層周圍各9個像素，則稱這兩個像素p和q是26連通的。接下來需要判斷每個孔隙是有機孔隙還是無機孔隙，在此根據孔隙周圍是有機質還是無機礦物佔優(含量所佔的比重較大)，如果是有機質，則該孔隙為有機孔隙；如果是礦物，則該孔隙為無機孔隙；如果孔隙周圍的有機質和無機含量相差不大，則稱該孔隙為邊界孔隙。

通過上述定義，首先需要對孔隙的連通性進行分析，從而得到有多少個孔隙。但是，由於三維fib-sem圖像的像素十分大(一般為1000*1000個像素或更大)，層數也多(一般為500～1000層)，直接採用三維連通性分析將十分耗時，因此該方法可行性不大。對此，我們先對每個層面的二維孔隙進行連通性分析，然後進行邊緣提取，再對邊緣進行z方向的連通性分析，通過這種方式的處理，需要分析的數據量將減少很多。

對二維孔隙進行分析，首先進行連通性分析，像素p和q，如果q在像素p周圍的8個像素，則稱這兩個像素是8連通的。從而對每個二維層面的每個孔隙進行標註。標註後的孔隙的像素值為1～m1(m1為孔隙總個數)。標註的目的是為了區分每個孔隙，若像素彼此為一個孔隙，則標註的數值為相同，若像素屬於不同孔隙，則標註的數值不同。

對於一個孔隙，由於內部孔隙周圍仍是孔隙，對判斷孔隙類型沒有影響。因此，只需要對孔隙的邊緣進行判斷，在此採用提取邊緣，常用的邊緣提取模板有laplacian算子、roberts算子、sobel算子、laplacian-gauss算子、kirsch算子和prewitt算子等。

然後對二維孔隙邊緣進行接觸面的類型判斷，統計周圍有機質和無機礦物的數量，分別得到每個孔隙周圍有機質和無機礦物像素的數目num_o和num_i。為了減少運算量，本發明將這兩個影響因子轉化為一個，如採用下面公式(2)的變換，將兩個影響因子轉化為在區間[0,1]的一個影響因子n，以便於下一步三維孔隙的分析。若n值越接近1，則說明孔隙周圍的無機礦物相對含量越多；n值越接近0，則說明孔隙周圍的有機質含量越多；n值越接近0.5，則說明孔隙周圍的有機質和無機礦物的含量相差不大。

首先對每個層面的孔隙進行連通性分析，然後進行邊緣提取。在本發明實施例中，我們採用的是sobel算子。對二維孔隙邊緣進行分析，統計周圍有機質和無機礦物的數量，並利用公式(2)中，得到每一個二維孔隙的影響因子n。

5.三維孔隙分析

將步驟4中每個孔隙的邊緣沿著z方向進行連通性分析，得到三維孔隙的連通性分析，統計得本發明實施例孔隙共有23524個。然後對每個孔隙的數據n進行數據統計，得到均值mean_n，這個數據在區間[0,1]內。而且越接近0，表示該孔隙周圍有機質的比重越高，有機質越佔優；越接近1，表示該孔隙周圍無機礦物的比重越高，無機礦物越佔優；越接近0.5，表示該孔隙周圍無機礦物和有機質的比重越相近。

6.三維孔隙類型判斷

根據三維孔隙的均值分析結果，判斷每個孔隙周圍是有機質還是無機礦物所佔的比重大(佔優)，若有機質明顯佔優，則認為該孔隙為有機孔隙；若無機礦物明顯佔優，則認為該孔隙為無機孔隙；若無二者所佔比重相差不大，則認為該孔隙為邊界孔隙。據此為依據，分別可得到有機孔隙、無機孔隙、邊界孔隙的三維空間分布圖。

在本發明實施例中，我們規定均值mean_n不大於0.4時，有機質明顯佔優，該孔隙為有機孔隙；mean_n不小於0.6時，無機礦物明顯佔優，該孔隙為無機孔隙；否則有機質和無機礦物都無明顯佔優，該孔隙為邊界孔隙。得到有機孔隙共計22479個，無機孔隙共計782個，邊界孔隙共計263個。圖5a為本發明實施例頁巖的有機孔隙三維分布圖，圖5b為本發明實施例頁巖的無機孔隙三維分布圖，圖5c為本發明實施例頁巖的邊界孔隙三維分布圖。其中灰色部分為分別為有機孔隙、無機孔隙和邊界孔隙，白色部分為其它物質，為了凸顯三維可視化立體效果，在圖中加入了光照和陰影設計，其中黑色部分為陰影。

7.分析

進一步地，對孔隙進行分析定量表徵。通過數學方法進行分析計算，根據步驟6得到的結果，分別對有機孔隙、無機孔隙、邊界孔隙進行分析，得到三維孔隙的各種表徵參數：比如孔隙的體積、表面積、中心位置、長度、寬度、高度、傾角、直徑等數據。

在此，我們主要關注於孔隙的直徑——數量分布和孔隙的直徑——體積分布。首先對整體孔隙進行分析，孔隙的體積為6.167um3，對應孔隙度為0.61vol％。圖6a為本發明實施例頁巖的孔隙的直徑——數量柱狀分布圖，圖6b為本發明實施例頁巖的孔隙的直徑——體積柱狀分布圖。對圖6a進行分析，可以看出頁巖中孔隙的直徑多分布於0～50nm、50～100nm之間。對圖6b進行分析，可以看出孔隙的體積主要分布在10～50nm、50～100nm和100～200nm之間。

對有機孔隙、無機孔隙和邊界孔隙分別分析定量表徵，得到有機孔隙、無機孔隙和邊界孔隙的體積分別為3.887um3、2.191um3、0.089um3，對應的有機孔隙、無機孔隙和邊界孔隙的孔隙度分別為0.38％、0.22％、0.01％。圖6c為本發明實施例頁巖的有機孔隙、無機孔隙和邊界孔隙的直徑——數量柱狀分布圖，圖6d為本發明實施例頁巖的有機孔隙、無機孔隙和邊界孔隙的直徑——體積柱狀分布圖。

對圖6c進行分析，可以看出頁巖中的有機孔隙和邊界孔隙直徑主要分布於0～50nm和50～100nm之間，但超過400nm的孔隙只存在於有機孔隙中；而無機孔隙主要分布在10～50nm和50～100nm和100～200nm之間。

對圖6d進行分析，可以看出頁巖中有機孔隙的體積主要分布在10～50nm和50～100nm之間；而邊界孔隙主要的體積主要分布在50～100nm和100～200nm之間；無機孔隙的體積主要分布在100～200nm和200～300nm之間。