一種基于測井電成像圖像的砂礫巖巖性自動識別方法

袁 陽， 侯國強， 馬濤濤， 周 輝， 董漢強， 潘 濤

(1.中國石油集團測井有限公司吐哈分公司，新疆維吾爾自治區 吐魯番 838200；2.延長油田油氣勘探公司 石油勘探開發部，陜西 延安 716000)

砂礫巖油藏埋藏深，儲集性能非均質性強，復雜的巖性導致砂礫巖地層的油氣儲集空間類型變化多樣。復雜的巖石骨架削弱了儲集流體的常規測井響應，造成常規測井信息識別、劃分砂礫巖巖性、評價儲層難度很大。

由于電成像測井資料分辨率高，可得到砂礫巖地層中礫巖顆粒形狀、大小、空間位置及砂巖和泥質條帶的直觀描述，已經應用到砂礫巖儲層的巖性識別研究中。但以往的研究主要是巖性的定性識別，表現為基于電成像圖像的圖像-巖性識別模式，需要解釋人員根據模式人工做出巖性的判斷、劃分，沒有實現計算機巖性自動定量識別[1-4]。目前，文獻[5]采用數字圖像處理技術對砂礫巖圖像進行特征提取，再基于人工神經網絡自動識別巖性；文獻[6]提取了電成像圖像中的均值、方差、分形維數、包絡和非均質5條特征曲線，分別利用相關性和神經網絡兩種方法進行巖性自動識別；文獻[7]通過圖像灰度轉化，構建灰度共生矩陣，計算樣本的對比度、相關度、熵、均勻度和能量5個特征，利用貝葉斯判別分析法對巖性進行自動判別；文獻[8]通過計算圖像灰度值的平均值和方差，以及方差曲線的變差方差根，構建巖性的定量判別標準模板，進行巖性定量識別；文獻[9]和文獻[10]從電成像測井圖像中提取定量表征礫石含量的結構指示參數，結合常規測井資料建立交會圖，定量自動識別巖性。但這類基于圖像模板特征提取的識別方法，未考慮圖像的井眼覆蓋率，圖像中極板間的空白條帶未充填，使得樣本圖像的特征計算結果不準確；并且這類方法依賴建立的圖像標準模板，標準模板的有限性導致識別各種多變的未知圖像時，出現準確率低、泛化能力差等問題。

本文以圖像中礫巖顆粒提取為核心，首先采用MCA(Morphological Component Analysis)算法充填圖像中極板間的空白條帶，生成全井眼圖像,再進行基于邊緣流的多尺度礫巖顆粒邊緣提取和基于曲線演化的礫巖顆粒區域分割，得到初步分割后的礫巖顆粒區域。但由于電成像測量技術、成像技術的限制和復雜井況，初步分割的顆粒區域與人工視覺判別結果有很大區別，因此在處理井段引入區塊地質沉積狀況時，特別是在將相應沉積條件下的礫巖顆粒磨圓度這一先驗知識作為約束條件時，自頂向下地驅動砂礫巖顆粒區域的合并和分裂。最后，構建砂礫巖粒度譜和巖性精細剖面，得到砂礫巖巖性的精細評價成果。本文方法應用于YJ油田深層砂礫巖難動用油藏的測井評價中，為優選射孔井段和壓裂層段提供了準確依據。

1 電成像圖像空白條帶充填

由于電成像測井在大多數井眼條件下都存在空白條帶，達不到全井眼覆蓋，為保證礫巖顆粒提取的準確度，采用MCA算法[11]對空白條帶進行充填。MCA將圖像分解為幾何結構un和紋理ut兩種形態成分，并假設兩個形態成分在某個特定的基庫或過完備字典下是稀疏的，且基庫或過完備字典之間具有不相干性。MCA模型表示如下[12]：

(1)

圖1 MCA算法具體實現步驟

圖2 原始成像圖及處理結果示例

MCA算法具體實現結果如圖2(b)所示，可以看到圖2(a)的原始圖像經算法充填后，較好地保持了礫巖顆粒邊緣的連續性。本文以每次迭代得到的充填圖像與原始圖像在掩膜模板上的平均像素灰度誤差作為迭代結束的條件。表1給出了圖2(a)原始圖像分別在迭代5，10，20，30，40，50和100次時的平均像素灰度誤差值。

由表1看到，算法迭代到20次時，誤差急劇下降，隨著迭代次數的增加，誤差變化不大，說明算法迭代20次以后，充填結果基本相同，實驗結果也證明了這一點。因此，通過適當地選擇誤差閾值(或兩次迭代的差異值)，可以自適應地結束迭代，得到可以接受的充填結果，從而避免固定較大迭代次數的方法帶來的計算冗余。

表1 不同迭代次數下的平均像素灰度誤差值

2 基于邊緣流的多尺度礫巖顆粒邊緣提取

邊緣流算法[13]在初始尺度計算邊緣流能量，并將各個像素點的邊緣流能量投影到水平、垂直方向，進行累加，得到初始邊緣流向量場；隨著尺度提升，當初始向量場幅度小于設定幅度閾值時，當前尺度的向量場等效于初始向量場；當初始向量場與當前向量場角度小于設定的角度閾值時，則將當前向量場疊加上初始向量場。在多個尺度重復上述過程，直至到達最大設定尺度，得到最終的多尺度邊緣流向量場。最后，根據Possion方程，由最終的多尺度邊緣流向量場的散度計算得到圖像的邊緣流函數，結果如圖2(c)所示，邊緣流函數通過歸一化轉換為灰度圖像顯示，可以看到灰度值較小的礫巖顆粒邊緣，但有些邊緣并不連續。

3 基于曲線演化的礫巖顆粒區域分割

基于曲線演化的砂礫巖顆粒區域分割算法[14]，是根據上節計算得到的邊緣流函數，通過求解水平集方程，計算區域分割單像素曲線，其特點是可以將不閉合的斷裂邊緣進行封閉。水平集方程如下所示[14]：

(2)

式中，φ為迭代計算曲線演化結果圖像；t為時間；S為邊緣能量場；▽為取散度；w為平滑權重系數；C為邊緣流函數；K為一個表示收縮或擴張的常數。最后對曲線演化結果圖像分別以高斯函數水平、垂直方向的一階導數為卷積核，進行卷積運算，將幅度值大于設定閾值的像素點選作邊緣點，結果如圖2(d)所示，可以看到藍色的單像素寬度曲線封閉了各個礫巖顆粒區域。

4 區塊地質沉積約束條件下的砂礫巖顆粒區域合并與分裂

從圖2(d)中可以看到，各個區域被閉合的藍色單像素曲線區分開來，但有的礫巖顆粒區域明顯與人眼視覺判別結果不同，出現明顯的陡峭凹陷。圖3為砂礫巖顆粒區域分裂示例。從圖3(a)的區域分割結果可以看出，同樣各個區域被閉合的藍色單像素曲線區分開來，但有的礫巖顆粒區域存在明顯的粘連情況，兩個鄰近的礫巖顆粒被劃分為一個區域。這些礫巖顆粒的形狀特征出現的概率在當前地質沉積條件下很低。

分析當前區塊地質沉積狀況，由此得到礫巖顆粒磨圓度的統計規律，作為約束條件知識驅動顆粒區域的合并與分裂過程。由顆粒磨圓度判定凹陷的顆粒分割區域是否進行相鄰區域的合并，以及明顯粘連的顆粒分割區域是否進行區域的分裂，由一個區域分裂為多個區域，并將每個分裂后的區域作為單獨的一個礫巖顆粒。具體方法是采用自適應曲率計算技術，定位礫巖顆粒區域中邊界上的凹點、凹點對和鄰近凸點。由鄰近的非礫巖區域的相對位置，合并鄰近區域與該礫巖區域。同時，根據定位的凹點對之間的弦長和計算的對應弧長，與弧弦比閾值比較，將該礫巖區域分裂為多個礫巖顆粒，這樣就解決了初步分割區域的不完整和粘連情況。

圖4為砂礫巖顆粒區域合并示例。圖4(a)是圖2(d)使用并查集提取的各個分割區域，得到各區域的包含像素點，并以不同顏色表示，根據礫巖灰度閾值，判定是否為礫巖顆粒。提取礫巖顆粒邊界像素，并以順時針排序。自適應計算邊界曲率，定位凹點及其鄰近凸點，凹點和其鄰近凸點形成一黃色矩形區域，如圖4(b)所示。如果鄰近區域的質心在該矩形區域內，則合并鄰近區域；如圖4(c)中所示的左下方淺綠色表征的礫巖顆粒區域，經過區域合并處理后，形成了完整的礫巖顆粒區域，符合當前礫巖顆粒磨圓度的設定。圖3(b)為圖3(a)使用并查集提取的各個分割區域，計算礫巖區域各個凹點對之間的弦長以及對應的弧長，并根據凹點的次數權重，使用權重大于1的凹點對對應弦，切分該礫巖區域。圖3(b)中大面積的深藍色區域表明存在明顯的粘連情況；如圖3(c)中所示，該區域由紅色的凹點對之間相連的弦，分裂為兩個不同的區域(分別用深紅色與深藍色表征)，這樣，每個區域的形狀特征符合當前顆粒磨圓度的知識約束。

圖4 砂礫巖顆粒區域合并示例

5 粒度譜構建及精細巖性剖面計算

砂礫巖粒度譜的構建和巖性精細剖面的計算采用滑動窗口方法，如圖5所示，統計某深度點對應滑動窗口內，礫巖顆粒區域面積的大小，并由此劃分區域面積的BIN值(面積大小的統計分組數目，每個BIN(Binary In Net)代表了各自的面積大小范圍)，統計屬于每個BIN值范圍內的礫巖顆粒的數目，構成粒度譜，作為解釋成果直觀地展現砂礫巖顆粒的大小分布。如圖6所示，精細巖性剖面計算首先根據礫巖顆粒大小閾值，將礫巖顆粒細分為細礫巖、中礫巖和粗礫巖，并根據滑動窗口內各巖性顆粒面積除以滑動窗口面積，得到精細巖性剖面中各細分巖性的百分比，而砂巖巖性百分比由礫巖顆粒占比的殘余值計算，并形成曲線，由解釋軟件進行充填繪制，形成直觀的砂礫巖精細巖性剖面。

圖5 砂礫巖粒度譜構建原理

圖6 砂礫巖精細巖性剖面計算原理

6 應用效果分析

砂礫巖巖性精細評價成果圖及與NRA巖性識別方法對比如圖7所示。

圖7 砂礫巖巖性精細評價成果圖及與NRA巖性識別方法對比

圖7(b)給出了Y920井3555.517～3557.808 m井段的精細評價綜合處理成果圖，可以看到構建的粒度譜和精細巖性剖面，與圖7(a)所示的NRA(Neutron Resistivity Acoustic)巖性識別方法[15]成果圖進行對比可以看到，在解釋井段，兩者的對應性較好；解釋井段上部礫巖總組分偏少，對應NRA方法的礫狀砂巖；下部礫巖總組分逐漸增加，對應NRA方法的礫巖巖性。粒度譜直觀地展現了相應井段礫巖顆粒的大小分布，巖性剖面相對于傳統的礫巖、礫狀砂巖和含礫砂巖的巖性劃分，更加精細地區分了礫巖顆粒的大小，而不是籠統地根據井段中所有礫巖部分所占面積劃分砂礫巖巖性，不考慮礫巖顆粒的獨立性和特殊性。

圖8為Y936井3676～3677.6 m井段處理成果圖，巖性剖面中以細礫巖和中礫巖為主，對比常規解釋結論(油層)，表明細礫巖、中礫巖的物性、含油性較好，符合YJ油田砂礫巖體區塊的儲層“四性”關系，據此優選壓裂設計井段，壓裂成功實施后，獲得初產近20 t的產能，目前日穩產油10 t。

圖8 Y936井3676～3677.6 m井段處理成果圖

7 結束語

本文的巖性自動識別方法通過圖像預處理以及在地質沉積約束條件知識驅動下的精細圖像處理技術，有效地提取了圖像上的礫巖顆粒，并據此構建粒度譜和精細的巖性剖面，為用戶提供了直觀、準確的砂礫巖評價成果。

本文方法可在今后的老井復查和新井解釋中發揮重要作用，為開拓服務市場和油田的勘探開發提供重要的技術支撐。