基于分形理論探究碳酸鹽巖CT圖像二值化最佳閾值

吳國銘 李熙喆 高樹生 劉華勛 吳國銓 安為國

(①中國科學院大學滲流流體力學研究所，河北廊坊 065007； ②中國石油勘探開發研究院廊坊分院，河北廊坊 065007； ③蘭州工業學院電氣工程學院，甘肅蘭州 730050)

基于分形理論探究碳酸鹽巖CT圖像二值化最佳閾值

吳國銘*①李熙喆②高樹生②劉華勛②吳國銓③安為國②

(①中國科學院大學滲流流體力學研究所，河北廊坊 065007； ②中國石油勘探開發研究院廊坊分院，河北廊坊 065007； ③蘭州工業學院電氣工程學院，甘肅蘭州 730050)

吳國銘，李熙喆，高樹生，劉華勛，吳國銓，安為國.基于分形理論探究碳酸鹽巖CT圖像二值化最佳閾值.石油地球物理勘探，2017,52(5):1025-1032.

CT掃描構建數字巖心的一個重要環節是灰度圖像的二值化。本文基于圖像處理分析和分形理論計算孔滲參數，為CT圖像二值化提供了更準確、更適用的約束條件。孔洞型碳酸鹽巖不同形態和尺度的微觀孔隙累積數量與半徑分布遵循冪律關系，具有統計意義上的分形特征，且分維值與孔隙度之間存在非線性定量關系。通過自定義四連通像素值梯度突變算法與分形理論相結合的方法，統計分析了孔隙形態、數量等參數并計算了孔隙度，其平均誤差小于7.2%；改進“先二值化后邊緣識別”的常規方法，運用濾波降噪和邊緣識別算子識別標定CT灰度圖像微觀孔隙，既保證了精度又提高了效率。根據這兩種圖像處理技術可確定最佳灰度閾值，實現CT圖像二值化。該方法處理的二值圖像較好地保留了不同尺度微觀孔隙結構的形態和分布特征，因此可廣泛應用于常規巖心重構及后續的三維數字巖心的構建。

碳酸鹽巖 微觀孔隙 CT圖像 分形 孔隙度 最佳閾值

1 引言

分形理論的基本思想是客觀事物具有自相似的層次結構，局部與整體在形態、功能、信息、時間和空間等方面完全相同或具有統計意義上的相似性[1-3]。定量刻畫這種自相似結構復雜程度的參數是分數維[4]。分形理論在石油工業領域取得了一系列的成果及應用[5-10]。碳酸鹽巖油氣藏儲層縫、洞發育[11-13]，非均質性強，研究表明裂縫和孔洞等孔隙結構具有自相似性[14-16]，如何準確獲取分維值成為關鍵。通過巖心CT掃描圖像分析微觀孔隙結構的分維值得到了業界的關注，文獻[17-20]基于計盒維數法分析計算CT灰度圖像和二值化圖像不同孔隙結構的分維值，但未見針對所得結果進一步的推廣應用。將分維值與圖像處理技術相結合，確定最佳二值閾值，實現圖像孔隙與固相最優分割，目前尚未看到相關文獻。

數字巖心可用于研究多孔介質內部微觀孔隙結構及流動規律，能從微觀角度揭示油氣儲集和滲流機理，獲取儲層宏觀參數，為油氣藏高效開發提供有力的技術指導[21-25]。CT成像技術是一種構建三維數字巖心的重要方法，通過對巖心樣品不同斷層X射線衰減程度的測量生成連續的反映微觀孔隙信息的斷面圖像，并將處理后的二維圖像疊加形成三維數字巖心。CT掃描圖像的合理分割可對數字巖心準確性產生重要影響，在圖像二值化處理中如何準確分割固相與孔隙，取決于分割閾值的合理確定。最大間距法、分水嶺算法和神經網絡等均可實現圖像二值分割，區分固相與孔隙，但上述方法僅針對灰度信息進行處理，未考慮圖像所含的物理信息，如孔隙度等，導致分割結果與實際存在偏差。雖然以實測孔隙度為約束條件確定的分割閾值能更準確地反映儲層微觀孔隙特征[26-27]，但對于形狀不規則的巖心樣品，往往無法通過實驗獲取其孔隙度，給二值閾值的確定帶來不便，本文基于分形理論，通過對部分圖像的預處理，獲取微觀孔隙的分維數，最終得到理論孔隙度作為分割閾值的約束條件，為CT圖像的二值化處理提供新的方法和思路。

2 分維值與孔隙度相關性分析

碳酸鹽巖固體顆粒級配良好，直徑λ跨度范圍約為三個量級，且連續分布，如從若干微米到若干毫米，具有統計自相似性和標度不變性[28]，伴生的孔洞滿足累積數量與直徑分布的冪律關系，為分形多孔介質。

根據分形理論，在d維歐式空間中的分形體的某個測量物理量滿足[21]

(1)

式中：F(γ)為待測物理量；L0為參考線性常數；γ為測量線性尺度；d為歐式空間維數；D為分形維數。

在無標度區間內，多孔介質單位體積中半徑大于λ的孔隙累積數與半徑分布滿足

Nc(λ≤L)=Aλ-D

(2)

式中：Nc為單位體積中半徑大于λ的孔隙累積數；λ為孔隙半徑；L為測量尺度；A為比例系數。對式(2)中的λ求微分，可得單位體積分形多孔介質中孔隙半徑的概率密度函數

(3)

推導出孔隙度[21,22]與分形維數的關系式為

(4)

式中：φ為孔隙度；λmax為孔隙最大半徑；λmin為孔隙最小半徑。

綜合式(3)、式(4)可得[29]

(5)

若已知多孔介質單位體積孔隙數量N、孔隙最大半徑λmax和最小半徑λmin、分維值D，即可計算出巖心孔隙度φ。

3 CT圖像處理分析

3.1 CT成像技術簡述

X射線CT是一種借助于計算機將X射線斷層掃描的斷面重現出來的成像技術，可實現非透明物體組成及結構的無損檢測，當強度均勻的射線束透照射物體時，由于內部不同組分的密度、結構等的差異，使不同部位透射射線強度不同，能量衰減程度滿足比爾定律

(6)

式中：I0、I分別為X射線穿過物體前、后的能量；μi、Li分別為第i種組分的衰減系數和透射路徑長度。

X-CT機的基本結構包括：X射線源、樣品夾持器和X射線檢測器(圖1)。X射線源發射X射線，經由載物臺上的樣品(保持旋轉狀態)吸收后，由X射線檢測器接收衰減后的信號并傳送給計算機，最后輸出樣品掃描圖像。

圖1 X射線CT原理圖

選取4塊碳酸鹽巖儲層巖心(基本物性參數見表1)進行CT掃描(縱向掃描間隔為0.076mm)，每塊巖心得到1300余張圖像(圖2)，每塊巖心等間隔選取400張圖像作為預處理對象，計算其分維值和孔隙度，從而確定最佳分割閾值。

表1 巖心物性參數表

圖2 巖心原始CT圖像

3.2 盒維數法確定分維值

本次實驗中巖心掃描后得到的原始CT圖像為RGB模式，為便于運用計盒維數法確定微觀孔隙結構分維值，將RGB圖像轉換為像素值為0～255的灰度圖像(圖像大小不變，像素點總數仍為990×990)，可通過灰度值表征掃描區域的不同組分、結構等信息，如固相顏色較淺、灰度值較大，微觀孔隙顏色較深、灰度值較小。

CT灰度圖像中，灰度值可視為沿Z軸方向高度為z的微小長方柱體，這樣，掃描區域的灰度值具有高低起伏的空間形態(圖3)，每個像素點的自由度為3(指平面坐標x、y和灰度值z)。

圖3 CT圖像中各像素點的灰度值柱狀分布圖(局部)

針對由一系列微小長方柱體構成的灰度空間，采用計盒維數法確定微觀孔隙結構分維值[17]，其基本原理如下： 用高度為hi的微小長方柱體(盒子)度量起伏不定的灰度空間，如某點灰度值為zi，所需盒子的數量為

(7)

(8)

圖關系曲線

基于1號巖心400張預處理灰度圖像得到的400個分維值D分布情況如圖6所示，可知分維值D分布在(1.6411，1.6507)范圍內，沿巖心軸向分布的波動性較小。

圖曲線線性部分

圖6 分維值沿巖心軸向分布

綜合400張預處理圖像，得到平均分維值，見表2。分維值大小可表征微觀孔隙發育程度及非均質性強弱，分維值分布與巖心特征相吻合。如1、3號巖心孔洞、裂縫較發育，非均質性較強，其平均分維值較大； 2號巖心部分發育裂縫，其平均分維值次之； 4號巖心孔隙相對不發育，其平均分維值最小。

3.3 微觀孔隙結構邊緣提取

為獲取單位體積多孔介質中的微觀孔隙總數N，需準確識別灰度圖中微觀孔隙結構邊緣，統計每張圖像中單連通孔隙數量。

在邊緣檢測前，運用鄰域大小為3×3的中值濾波器[31]對灰度圖像進行降噪處理，可在不改變微觀孔隙與固相對比度的前提下，有效地提高邊緣檢測精度。

在灰度圖像中，邊緣是指研究對象的邊界，即反映灰度值劇烈變化的曲線，邊緣檢測的目的就是尋找這些灰度值劇烈變化的像素點，可通過MATLAB的edge函數實現灰度圖中微觀孔隙結構邊緣的有效檢測。edge函數邊緣檢測方法的優點在于使用了兩種不同的閾值分別檢測強邊緣和弱邊緣，只有在強、弱邊緣相連時，輸出圖像中才包含弱邊緣，真正的弱邊緣易于檢測，且有效避免了噪聲填充。edge函數求取邊界的算子較多，如sobel、prewitt、log、roberts和zeroscross等，本次實驗運用canny算子完成邊界檢測。

BS-240VET全自動生化分析儀（深圳邁瑞生物醫療電子股份有限公司）；Bio-Rad D-10 HbA1C分析儀（美國Bio-Rad公司）；α1-微球蛋白（α1-MG）試劑盒（德靈公司）；血清尿清蛋白/尿肌酐（UALB/Cr）、視黃醇結合蛋白（RBP）、血清胱抑素-C（Cys-C）（南京威科特曼公司產品）；同型半胱氨酸（Hcy）、尿β2微球蛋白（U-β2-MG）試劑盒（寧波瑞源生物科技有限公司）；糖化血紅蛋白（HbA1C）試劑盒（美國Bio-Rad公司）。C-反應蛋白（CRP）（寧波美康生物科技股份有限公司）。

canny算子是一種較新的邊緣檢測算子，具有較高的邊緣檢測性能，得到了越來越廣泛的應用。canny算子邊緣檢測的基本原理如下： ①求取高斯濾波器f與圖像函數g(x，y)的卷積(式9)； ②計算圖像梯度的幅值和方向； ③對梯度圖像應用非極大值抑制，確定梯度局部最大值

(9)

對于canny算子更深入的理論解析不做贅述，詳細推導過程可參考文獻[33]。

圖7展示了運用edge函數的canny算子提取得到的灰度圖像中微觀孔隙結構的邊緣信息，將原始灰度圖像中絕大部分孔隙結構包含在內(其余更微小的孔隙結構限于圖像分辨率無法表征)，實現了單連通孔隙的有效識別與準確劃分。

圖8所示為1號巖心每張預處理圖像(共計400張)中微觀孔隙數量沿巖心軸向分布特征，溶蝕段孔隙密集發育，從圖像中識別出的微觀孔隙數量較多； 致密段孔隙稀疏分布，從圖像中識別出的微觀孔隙數量較少，單張圖像最大的孔隙數量可達最小值的3倍之多，體現了碳酸鹽巖較強的非均質性。獲取孔隙數量后，可進一步計算得到單位體積多孔介質中的孔隙總數N(見表2)。

3.4 微觀孔隙結構參數計算

由圖7b可知，canny算子輸出圖像為二值圖像，孔隙邊界點為1，其余點均為0。針對該二值圖像，采用自定義4連通像素值梯度突變算法，自主編程實現單連通孔隙結構參數(如面積)的計算，自定義像素值梯度計算公式為

=f(i,j+1)-f(i,j)

(10)

式中：i、j為像素點坐標；f為像素值。在二值圖中，固相和孔隙相中的像素值梯度均為0，而在兩相邊界處像素值梯度出現突變，為±1。

圖7 canny算子邊緣識別效果圖

圖8 微觀孔隙數量沿巖心軸向分布特征

遍歷二值圖像矩陣，若第i行第j列的像素值為0，則自第j列始逐列計算像素值梯度; 若出現突變，確定第i行的左右邊界，第i行也是該單聯通孔隙相的上邊界。

判斷第i行第j列的左邊界點下方(第i+1行第j列)像素值： ①若為0，則以第i+1行第j列像素點為始，分別向左、向右計算不同列的像素值梯度，通過突變確定第i+1行的左右邊界； ②若為1，則判斷第i行第j+1列像素點下方(第i+1行第j+1列)像素值。重復上述過程，確定第i+1行的左右邊界。

如果循環至第n行，0值下方均為1，可確定第n行為該單聯通孔隙相的下邊界，循環結束。這樣即可找到一個完整的單連通孔隙區域，并對其所含全部像素點及邊界點分別標記。

運用上述算法，可準確計算單張二值圖像相互獨立的微觀孔隙結構面積(即所含像素點總數量)。但由于真實孔隙的位置、形態、大小等復雜多樣，難以獲其特征參數，可將各孔隙視為與其等面積的圓，進而計算其等效半徑λ。表2為統計不同巖心預處理圖像得到的微觀孔隙等效半徑最大值λmax和最小值λmin。圖9為孔隙等效半徑分布特征，最大等效半徑約為80像素點，即5.976mm，最小等效半徑約為2像素點，即0.156mm，表明孔隙形態、位置和大小差異性較大，具有較強的非均質性； 同時，較好地實現了微小孔隙的識別和標定，程序穩定性和準確性較好。

圖9 孔隙等效半徑沿巖心軸向分布特征

由表2可知，對于縫洞型碳酸鹽巖儲層，普遍發育橫截面積較大的孔隙，如1號巖心的孔隙最大等效半徑為5.976mm，具有較好的儲集特性，分析CT圖像得到的認識與常規物性實驗結果相一致，充分體現了可視化方法的直觀性和準確性。

上述所有過程均由自編的MATLAB程序自動批量處理并輸出計算結果。

表2 計算孔隙度所需各參數值

3.5 確定孔隙度和最佳分割閾值

經過前述分析計算，已經獲取了微觀孔隙平均分維值D、單位體積多孔介質中的孔隙總數N和孔隙等效半徑λmax、λmin。現在，將各參數代入式(6)，即可得到巖心孔隙度(表3)，計算所得孔隙度和實驗孔隙度相近，最大相對誤差為8.7%，表明基于理論公式計算的孔隙度精度較高，符合儲層實際，該孔隙度可作為CT灰度圖像二值分割的約束條件，特別是針對未開展常規孔隙度測試的巖心，其CT圖像二值分割最佳閾值滿足[32]

f(k*)

(11)

式中：Imin、Imax分別為CT圖像最小和最大灰度值；P(i)為像素點i的灰度值；φ為計算孔隙度；k*為最佳閾值。

將得到的計算孔隙度代入式(11)，確定的最佳分割閾值見表3，分割后的二值圖像如圖10所示，較好地保留了不同尺度微觀孔隙結構的形態和分布特征，為后續構建三維數字巖心奠定堅實基礎。

從理論上來說，同一儲層段的巖心樣品CT圖像應取相同的最佳閾值，但由表3可知，4塊巖心對應的最佳閾值并不相同，而是存在微小差異。逐次選取表3中的一個最佳閾值，對4塊巖心CT圖像開展相同閾值下的二值分割，統計孔隙度相對誤差，結果如圖11所示。

表3 孔隙度及相對誤差和最佳閾值

圖10 最佳閾值確定的二值圖像

圖11 不同閾值對應的孔隙度相對誤差

由圖11可知，采用相同閾值對CT圖像分割后，4塊巖心所得孔隙度相對誤差差異較為明顯，不適宜對同一儲層段的巖心樣品的CT掃描圖像采用完全相同的閾值開展二值分割。上述結果可能是由于不同巖心樣品中礦物成分及顆粒膠結程度的差異與CT掃描中的噪聲干擾等所導致。

4 結論與認識

本文基于縫洞型碳酸鹽巖CT圖像分析處理，探究了微觀孔隙結構分形特征和圖像二值分割最佳閾值的確定，主要結論如下。

(1)CT掃描技術能實現不同尺度微觀孔隙結構的可視化表征，為研究多孔介質內部孔隙網絡和數字巖心的構建提供了條件。

(2)碳酸鹽巖儲層屬于固相分形多孔介質，其中伴生的縫洞等微觀孔隙結構滿足累積數量與半徑分布的冪律關系，即具有統計意義上的分形特征；分維值大小可表征微觀孔隙發育程度及非均質性強弱，微觀孔隙等效半徑范圍越大，結構越復雜，分維值越大，且分維值與孔隙度之間存在非線性定量關系。

(3)對灰度圖像的適當降噪處理及邊緣識別算子的合理運用，可實現對灰度圖像微觀孔隙的有效標識，有別于先二值化后邊緣識別的常規方法，在保證精度的同時提高了效率；自定義4連通像素值梯度突變算法可有效統計微觀孔隙結構參數，準確獲取孔隙最大和最小等效半徑。

(4)分維值等參數計算所得孔隙度與實驗孔隙度的相對誤差較小，表明通過理論公式所得孔隙度精度較高，可作為確定二值分割最佳閾值的約束條件，以此所得的二值圖像中較好地保留了不同尺度微觀孔隙結構的形態和分布特征，為后續構建三維數字巖心奠定了堅實基礎。

(5)本文為探究實驗孔隙度未知的巖心CT圖像如何實現準確的二值化處理提供了新的角度和方法。

[1] Mandelbrot B B.The Fractal Geometry of Nature.Freeman San Francisco，1983.

[2] 翟云芳.滲流力學.北京：石油工業出版社，2009，152-156. Zhai Yunfang.Seepage Flow Mechanics.Petroleum Industry Press，Beijing，2009，152-156.

[3] 聶篤憲.隨機分形裂縫模型的研究[學位論文].廣東廣州：廣東工業大學，2005. Nie Duxian.Study of Random Fractal Fracture Model [D].Guangzhou University of Technology，Guangzhou，Guangdong,2005.

[4] 鄭紅軍.裂縫性儲集層的分形分維應用研究[學位論文].四川成都：西南石油大學，2003，2-12. Zheng Hongjun.Application of fractal dimension of fractured reservoirs [D].Southwest Petroleum University，Chengdu，Sichuan,2003，2-12.

[5] 施澤進，羅蟄潭，彭大鈞等.油氣藏體積分布的分形特征及其意義.大地構造與成礦學，1995，19(1)：71-72. Shi Zejin，Luo Zhetan，Peng Dajun et al.Fractal pro-perty of size distributions of oil-gas deposites.Geotectonica et Metallogenia，1995，19(1)：71-72.

[6] 潘葆芝，薛林福.分數維及其在測井地質解釋中的應用.測井技術，1992，16(3)：215-217. Pan Baozhi，Xue Linfu.Fractal dimension and its application to geologic interpretation of log data.Well Logging Technology，1992，16(3)：215-217.

[7] 鄭繼明.分形和小波變換在油氣預測中的應用.預測，1996，(6)：51-52. Zheng Jiming.Fractal and wavelet transform in application of oil-gas prediction.Forecasting，1996，(6)：51-52.

[8] 隋少強，宋麗紅，李淑敏.分形理論在儲集層橫向預測和油氣層判斷中的應用.新疆石油地質，2001，22(4)：342-343. Sui Shaoqiang，Song Lihong，Li Shumin.The application of fractal theory in lateral prediction of reservoir and hydrocarbon diagenesis.Xinjiang Petroleum Geo-logy，2001，22(4)：342-343.

[9] 蔣加鈺，段玉順.分形理論在油氣檢測中的應用.石油物探，1995，34(1)：49-50. Jiang Jiayu，Duan Yushun.Application of fractal theory to hydrocarbon detection.GPP，1995，34(1)：49-50.

[10] 李春華，魯港，毛俊莉等.分形理論在井間薄儲集層描述中的應用.石油學報，1998，19(4)：27-28. Li Chunhua，Lu Gang，Mao Junli et al.Application of fractal theory in the inter-well thin reservoir description.Acta Petrolei Sinica，1998，19(4)：27-28.

[11] 楊柳，沈亞，管俊亞等.多維數據裂縫檢測技術探索及應用.石油地球物理勘探，2016，51(增刊1)：58-61. Yang Liu，Shen Ya，Guan Junya et al.Fracture detection based on multi-dimension data.OGP，2016，51(S1)：58-61.

[12] 吉學林.碳酸鹽巖裂縫型儲層預測.石油地球物理勘探，2013，48(增刊1)：95-98. Ji Xuelin.The characterization of fractures in carbo-nate reservoirs by full azimuth seismic data.OGP，2013，48(S1)：95-98.

[13] 凌云，高軍，孫德勝等.裂縫儲層的地震預測技術與應用實例.石油地球物理勘探，2015，50(1)：92-96. Lin Yun，Gao Jun，Sun Desheng et al.A case study of fractured reservoir prediction.OGP，2015，50(1)：92-96.

[14] 文慧儉，閆林，姜福聰等.低孔低滲儲層孔隙結構分形特征.大慶石油學院學報，2007，31(1)：16-18. Wen Huijian，Yan Lin，Jiang Fucong et al.Fractal property of low porosity and permeability reservoirs’ pore structure.Journal of Daqing Petroleum Institute，2007，31(1)：16-18.

[15] 王自明，宋文杰，戴勇等.利用分形模擬建立裂縫型碳酸鹽巖儲層滲透率變異函.天然氣工業，2007，27(11)：73-74. Wang Ziming，Song Wenjie，Dai Yong et al.Application of fractal simulation technique to modeling of fractured carbonate reservoir.Natural Gas Industry，2007，27(11)：73-74.

[16] 鄧攀，陳孟晉，楊泳.分形方法對裂縫性儲集層的定量預測研究和評價.大慶石油地質與開發，2006，25(2)：19-20. Deng Pan，Chen Mengjin，Yang Yong.The application of fractal approach to the quantitative estimation research and evaluation of fractured reservoir.Petroleum Geology and Oilfield Development in Daqing，2006，25(2)：19-20.

[17] 彭瑞東，楊彥從，鞠楊等.基于灰度CT圖像的巖石孔隙分形維數計算.中國科學，2011，56(26)：2258-2264. Peng Ruidong，Yang Yancong，Ju Yang et al.Fractal dimension of pore structure based on gray CT images.Science China Press，2011，56(26)：2258-2264.

[18] 何娟，劉建立，呂菲.基于CT數字圖像的土壤孔隙分形特征研究.土壤，2008，40(4)：663-665. He Juan，Liu Jianli，Lv Fei.Characterizing fractal characteristics of soil pores based on high-resolution digital CT images.Soils，2008，40(4)：663-665.

[19] 馮杰，于紀玉.利用CT掃描技術確定土壤大孔隙分形維數.灌溉排水學報，2005，24(4)： 27-28. Feng Jie，Yu Jiyu.Determination fractal dimension of soil macropore using computed tomography.Journal of Irrigation and Drainage， 2005， 24(4)： 27-28.

[20] 陳學松，吳愛祥，楊保華等.基于CT圖像法的砂巖散體孔隙結構分形特征研究.金屬礦山，2007，37(10)：75-76. Cheng Xuesong，Wu Aixiang，Yang Baohua et al.Study on fractal characteristics of pore structure of bulk rock-ore in dump learching based on CT image.Metal Mine，2007，37(10)：75-76.

[21] 王金波.巖石孔隙結構三維重構及微細觀滲流的數值模擬研究[學位論文].北京：中國礦業大學(北京)，2014，2-6. Wang Jinbo.3D Reconstruction of Porous Rock and Numerical Simulations of Fluid Flow at Mesoscale Levels [D].China University of Mining and Technology (Beijing)，Beijing,2014，2-6.

[22] 聶昕.頁巖氣儲層巖石數字巖心建模及導電性數值模擬研究[學位論文].北京：中國地質大學(北京)，2014，19-24. Nie Xin.Digital Core Modeling and Numerical Study of Electrical Conductivity of Shale Gas Reservoir Rock [D].China University of Geosciences(Beijing)，Beijing,2014，19-24.

[23] 劉學峰.基于數字巖心的巖石聲電特性微觀數值模擬研究[學位論文].山東青島：中國石油大學(華東)，2010，2-20. Liu Xuefeng.Numerical Simulation of Elastic and Electrical Properties of Rock Based on Digital Cores [D].China University of Petroleum (East China)，Qingdao,Shandong,2010，2-20.

[24] 張挺.基于多點地質統計的多孔介質重構方法及實現[學位論文].安徽合肥：中國科學技術大學，2009，3-6. Zhang Ting.Research on the Reconstruction of Po-rous Media Based on Multiple-point Geostatistics [D].University of Science and Technology of China，Hefei,Anhui,2009，3-6.

[25] 趙秀才.數字巖心及孔隙網絡模型重構方法研究[學位論文].山東青島：中國石油大學(華東)，2009，18-22. Zhao Xiucai.Numerical Rock Construction and Pore Network Extraction [D].China University of Petroleum (East China)，Qingdao,Shandong,2009，18-22.

[26] 耿沖，楊永飛，高瑩.不同分辨率巖石CT圖像的優化處理方法.科學技術與工程，2014，14(2)：2-3. Geng Chong，Yang Yongfei，Gao Ying.Optimization of image processing method based on rock CT images of different resolutions.Science Technology and Engineering，2014，14(2)：2-3.

[27] 姜蓉蓉，肖沙里，高家利.改進分水嶺算法的彩色圖像分割技術.計算機工程與應用，2009，45(19)：172-173. Jiang Rongrong，Xiao Shali，Gao Jiali.Color image segmentation based on improved watershed transform.Computer Engineering and Applications，2009，45(19)：172-173.

[28] 孔祥言.高等滲流力學.安徽合肥：中國科學技術大學出版社，2010. Kong Xiangyan.Higher Seepage Mechanics.University of Science and Technology of China Press，Hefei,Anhui,2010.

[29] 劉俊亮，田長安，曾燕偉等.分形多孔介質孔隙微結構參數與滲透率的分維關系.水科學進展，2006，17(6)：813-815. Liu Junliang，Tian Chang’an，Zeng Yanwei et al.Fractal dimensionality dependence of micro structural parameters and permeability in fractal porous media.Advances in Water Science，2006，17(6)：813-815.

[30] 姜虹，范宜仁，陳華.應用分形理論計算碳酸鹽巖地層膠結指數m值.測井技術， 2012， 36(3)： 251-252. Jiang Hong，Fan Yiren，Chen Hua.On calculation of cementation indexmin carbonate formation using fractal theory.Well Logging Technology，2012，36(3)：251-252.

[31] 周麗，莊眾，成景旺等.利用自適應迭代多級中值濾波法分離海上多震源混合波場.石油地球物理勘探，2016，51(3)：434-436. Zhou Li，Zhuang Zhong，Cheng Jingwang et al.Multi-source blended wavefield separation for marine seismic based on an adaptive iterative multi-level median filtering.OGP，2016，51(3)：434-436.

[32] 趙秀才，姚軍，房克榮.合理分割巖心微觀結構圖像的新方法.中國石油大學學報(自然科學版)，2009，33(1)：65-66. Zhao Xiucai，Yao Jun，Fang Kerong.A new reasonable segmentation method for microstructure image of re-servoir rock.Journal of China University of Petroleum(Edition of Natural Sciences),2009，33(1)：65-66.

[33] 張強，王正林.精通MATLAB圖像處理.北京：電子工業出版社，2009. Zhang Qiang，Wang Zhenglin.Proficient MATLAB in Image Processing.Publishing House of Electronics Industry，Beijing,2009.

(本文編輯：朱漢東)

吳國銘 博士研究生，1989年生；2012年本科畢業于中國石油大學(華東)應用物理學專業，隨后在中國科學院大學滲流流體力學研究所流體力學專業取得碩博連讀資格，主要從事油氣田開發方面的學習和研究。

1000-7210(2017)05-1025-08

P631

A

10.13810/j.cnki.issn.1000-7210.2017.05.016

*河北省廊坊市44號信箱，065007。 Email：guomwu@163.com

本文于2016年10月18日收到，最終修改稿于2017年8月16日收到。

本項研究受國家重大科技專項(2011ZX05013-002)資助。