煤巖CT圖像的孔隙度和比表面積測量方法

金智敏 周宏偉 薛東杰

摘?要：為實(shí)現(xiàn)基于CT圖像的煤巖孔隙度和比表面積快速測量，提出了一種基于MATLAB圖像處理的測量方法。以平煤十二礦深部采集的煤巖為研究對象，應(yīng)用Nano CT技術(shù)對樣品進(jìn)行了掃描，并分別借助Avizo和MATLAB軟件計(jì)算孔隙度和比表面積，研究了Avizo三維建模和MATLAB圖像處理2種方法的優(yōu)缺點(diǎn)。發(fā)現(xiàn)2種方法計(jì)算出的結(jié)果誤差僅在10%左右，且以復(fù)合Cotes公式的誤差估計(jì)最小，孔隙度和比表面積的誤差分別為9.392%和8.247%.結(jié)果表明：傳統(tǒng)的基于Avizo三維建模的測量方法精確度較高，同時還能獲取樣品的三維孔隙網(wǎng)絡(luò)模型和孔隙的形狀、走向及孔徑分布，但計(jì)算過程復(fù)雜，對于復(fù)數(shù)個試樣必須進(jìn)行多次計(jì)算，且對CT圖像有較高的質(zhì)量要求。基于MATLAB圖像處理的測量方法則通過簡單的數(shù)值積分實(shí)現(xiàn)了快速計(jì)算煤巖的孔隙度和比表面積，避免了三維建模，極大地降低了計(jì)算復(fù)雜度，且能夠一次性處理復(fù)數(shù)個試樣，較好地滿足了實(shí)際工程需要。且對于CT圖像具有更高的容錯率，每次計(jì)算僅需其中數(shù)張掃描質(zhì)量效果較好的CT圖像即可。

關(guān)鍵詞：煤巖;CT圖像;孔隙度;比表面積;MATLAB;Avizo;

中圖分類號：TP 391

文獻(xiàn)標(biāo)志碼：A

文章編號：1672-9315（2020）01-0133-08

DOI：10.13800/j.cnki.xakjdxxb.2020.0118開放科學(xué)（資源服務(wù)）標(biāo)識碼（OSID）：

Measurement of porosity and specific surface

area of coal rock using CT images

JIN Zhi-min，ZHOU Hong-wei，XUE Dong-jie

（School of Mechanics and Civil Engineering，China University of Mining and Technology（Beijing），Beijing 100083，China）

Abstract：In order to realize the rapid measurement of porosity and specific surface area of coal using CT images，a measurement method based on MATLAB image processing is proposed.Taking the coal rock obtained from the deep part of Pingdingshan Coal Mine Group 12 as the research object，the sample was scanned by Nano-CT，and the porosity and specific surface area were calculated by Avizo and MATLAB software respectively.The advantages and disadvantages of Avizo 3D modeling and MATLAB image processing were studied.It is pointed that the errors calculated by the two methods are only about 10%，with the error done by the composite Cotes formulathe smallest.The errors of porosity and specific surface area are 9.392% and 8.247%，respectively.The results indicate that the traditional measurement method based on Avizo 3D modeling is more accurate，and the 3D pore network model and the shape，orientation and size distribution of pores can be obtained.However，the calculating process is comparatively complicated，and it is necessary to perform multiple calculations for a plurality of samples，with high quality requirements for CT images.The measurement method based on MATLAB image processing proposed in this paper realizes the rapid calculation of porosity and specific surface area of coal by simple numerical integration，which avoids 3D modeling，reduces the computational complexity greatly，and processes multiples at a time.The sample satisfies the engineering needs well.Moreover，for CT images，there is a higher fault tolerance rate，and only a few CT images with better scanning quality effects are needed for each calculation.Key words：coal rock;CT image;porosity;specific surface area;MATLAB;Avizo

0?引?言

煤巖是典型的多孔介質(zhì)，力學(xué)上常表現(xiàn)出非均勻性、非連續(xù)性和各向異性。其內(nèi)部結(jié)構(gòu)可視為由固體基質(zhì)和孔（裂）隙2部分組成，其中孔（裂）隙的微觀結(jié)構(gòu)特征表現(xiàn)為孔隙與喉道的幾何形狀、體積、分布及連通關(guān)系，與煤巖的力學(xué)性質(zhì)、滲流能力密切相關(guān)[1]。在煤礦瓦斯治理等實(shí)際工程問題中，孔隙度和滲透率等參數(shù)都具有極其重要的意義。Kozeny-Carman（KC）方程是多孔介質(zhì)滲流領(lǐng)域著名的半經(jīng)驗(yàn)公式之一，即

式中?K為滲透率，μm2;為多孔介質(zhì)的孔隙度;C為KC常數(shù)（一般取5）[2];S為比表面積，μm?-1.由式（1）可知，多孔介質(zhì)的滲透率可由其孔隙度和比表面積估算出。

孔隙度和比表面積有多種測定方法，常用的有低壓氣體吸附法（LPGA）[3-4]、壓汞法（MIP）[5]、核磁共振技術(shù)（NMR）[5]、聚焦離子束-掃描電鏡技術(shù)（FIB-SEM）[6]、計(jì)算機(jī)斷層掃描技術(shù)（CT）[6-7]、掃描電鏡技術(shù)（SEM）[8]等。其中計(jì)算機(jī)斷層掃描技術(shù)（Computed Tomography，CT）由于具備無損檢測物體內(nèi)部結(jié)構(gòu)、組分及損傷情況的能力，近年來廣泛應(yīng)用于煤炭行業(yè)。目前有許多學(xué)者在這方面做了大量工作。Zhang等基于CT掃描圖像利用Avizo軟件對某無煙煤進(jìn)行了三維重構(gòu)并研究其孔隙度和滲透率的演化[9]。

Mayo等利用CT技術(shù)研究了不同煤樣中Xe，Kr和CO2氣體吸附量隨時間的變化，并對氣體擴(kuò)散曲線和擴(kuò)散模型進(jìn)行了分析[10]。王剛等通過CT技術(shù)和三維重構(gòu)技術(shù)建立了6種煤樣的真實(shí)模型并進(jìn)行模擬試驗(yàn)，對非達(dá)西滲流中各參數(shù)對滲流的影響進(jìn)行了研究[11]。李偉等結(jié)合顯微CT技術(shù)和多孔介質(zhì)逾滲理論，對煤中孔隙進(jìn)行了三維可視化表征，定量研究了不同煤體結(jié)構(gòu)煤的孔隙連通性和滲透能力的差異[12]。

目前對于實(shí)驗(yàn)室尺度的煤樣，主要采用顯微CT技術(shù)表征其孔隙-裂隙結(jié)構(gòu)。煤巖中瓦斯氣體吸附、擴(kuò)散和緩慢層流滲透的主要空間是微米級和納米級孔隙，然而顯微CT掃描精度一般最高只有幾十微米，難以滿足對這些孔隙的測量。Nano CT精度較顯微CT高得多，但掃描樣本尺度只有微米級，使用常規(guī)方法難以快速測定其孔隙度和比表面積。鑒于此，提出了一種基于MATLAB數(shù)字圖像處理技術(shù)的孔隙度和比表面積測量方法。該方法基于少量的二維CT圖像信息，通過簡單的數(shù)值積分快速計(jì)算煤巖的孔隙度和比表面積，無須三維重構(gòu)，為定量研究煤巖的孔隙結(jié)構(gòu)提供了新的思路。并與基于Avizo三維建模的傳統(tǒng)方法對比，結(jié)果表明，借助MATLAB處理二維CT圖像來測量煤巖的孔隙度和比表面積，能較好地滿足工程需要。

1

基于Avizo三維建模的孔隙度及比表面積測定

試驗(yàn)煤巖取自平煤十二礦深部，埋深約1 000 m.煤巖作為一種低滲巖體，其孔隙尺度一般在1 μm左右。試驗(yàn)前已用ACTIS300-320/225工業(yè)CT檢測系統(tǒng)對煤巖進(jìn)行成像，分辨率達(dá)30 μm時仍難見孔隙，因此改用Nano CT掃描煤巖的孔隙分布。

Nano CT是在Micro CT基礎(chǔ)上發(fā)展出來的精度更高的計(jì)算機(jī)斷層掃描技術(shù)。其原理與Micro CT基本一致，均是根據(jù)樣品中不同部分對X射線的吸收和透過率差別，測量X射線在不同角度穿過樣品時的衰減，并采用濾波反投影法（FBP）[13]、有序子集-最大期望法（OSEM）[14-15]、代數(shù)重建法（ART）[16]等算法進(jìn)行圖像重建，將光信號轉(zhuǎn)變成數(shù)字信號，從而獲得樣品的斷面CT圖像[17]。CT圖像的分辨率一般取決于射線等效束寬BW，即

式中?BW為射線等效束寬，μm;d為探測器孔徑，μm;a為射線源焦點(diǎn)尺寸，μm，M=L/λ;M為放大倍數(shù);L為射線源到探測器距離，m;λ為射線源到掃描中心距離[18]，m.射線等效束寬越小，CT圖像的分辨率越高。

CT圖像的銳利度又與半影H有關(guān)，即

式中?H為半影，μm.如圖1所示，CT圖像由本影和半影組成，其中本影指投影時沒有任何光線到達(dá)的區(qū)域，半影指只有部分光線到達(dá)的區(qū)域。半影越小，CT圖像的銳利度越高。由式（2）和式（3）可知，CT圖像的分辨率受到射線源焦點(diǎn)尺寸、CCD分辨率（即探測器孔徑d）及放大倍數(shù)的約束，而在射線源焦點(diǎn)尺寸和CCD分辨率一定的情況下，任意增大放大倍數(shù)以期提高CT圖像的分辨率，可能會導(dǎo)致嚴(yán)重的半影模糊。因此，提高CT圖像分辨率的主要途徑就是縮小射線源焦點(diǎn)和提高CCD分辨率。

試驗(yàn)儀器選用Xradia Ultra-XRM L200立體顯微鏡，該CT檢測系統(tǒng)擁有65和15 μm 2種視場模式?；竟ぷ髟砣鐖D2所示。同步輻射源發(fā)出的X射線經(jīng)過單色器分離出單色X光后，再通過聚焦波帶片聚焦，并用小孔過濾掉零級光和高級衍射光，這樣在理論上就產(chǎn)生了一個焦點(diǎn)無窮小的射線源。通過小孔的X光照射到焦點(diǎn)附近的樣本上，再通過物鏡波帶片和顯微物鏡放大成像，同時使用相襯環(huán)提高透明物體的清晰度，以提高CCD分辨率。

試驗(yàn)選用掃描直徑為65 μm的大視場模式，測試電壓8 kV.最終所測得的圖像大小為1 024×1 024像素，每個像素點(diǎn)的邊長約為0.063 μm.試驗(yàn)最終測得1014張二維CT層析圖像，依次編號分別為1～1 014.為了使結(jié)果更加精確，截取了樣品內(nèi)部較清晰處的內(nèi)接正方體以保證樣品的代表性。選用正方體不僅便于統(tǒng)計(jì)孔隙度與孔隙結(jié)構(gòu)參數(shù)，而且立體展示效果更佳。最終所截取的正方體邊長為601像素，即只選取了其中601張601×601像素的圖像進(jìn)行分析。如圖3所示，以編號702的圖像為例，圖中黃色邊框內(nèi)即為分析范圍。

試驗(yàn)所得的圖片均保存為tiff格式。所以在將圖片導(dǎo)入三維圖像數(shù)據(jù)分析軟件Avizo前，需將得到的tiff圖片進(jìn)行格式變化。載入數(shù)據(jù)后，還需對二維CT層析圖像進(jìn)行圖像處理。一個標(biāo)準(zhǔn)的巖心圖像處理流程包括：灰度圖像濾波處理、灰度圖像的二值化、基于分水嶺算法的圖像分割、結(jié)果的分析以及滲流方向上的孔隙連通性判斷等等[20]。在Avizo軟件中，濾波降噪常用的命令模塊包括Median filter，Non-local Means filter，Anisotropic Diffusion filter等。選用各向異性擴(kuò)散濾波函數(shù)（Anisotropic Diffusion filter）對圖像進(jìn)行處理，并添加Interactive Thresholding命令模塊進(jìn)行圖像分割，其中閾值參數(shù)為59.仍以編號702的圖像為例，處理前后的圖像如圖4所示。圖4（a）即為原始切片圖像，圖4（b）則為圖4（a）經(jīng)過Avizo濾波降噪、閾值分割和二值化等圖像處理后得到的圖像。其中藍(lán)色部分為孔隙或裂縫，黑色部分為基質(zhì)礦物。

通過Avizo重構(gòu)的煤巖的三維數(shù)字模型如圖5所示，該模型選用編號200～800的二維CT圖像沿y軸排列并重構(gòu)。煤巖的孔隙結(jié)構(gòu)模型則如圖6所示。Avizo軟件主要基于分水嶺算法識別孔隙邊界并對各孔隙進(jìn)行標(biāo)記，其基本原理為將CT圖像視為具有拓?fù)浣Y(jié)構(gòu)的地形圖，用各點(diǎn)的灰度值表示該點(diǎn)的海拔高度，則每一個局部極小灰度值及其影響范圍構(gòu)成了一個盆地，而各盆地的交匯處則形成分水嶺。圖5和圖6的三維重構(gòu)模型中總計(jì)有6 360個孔隙。在Avizo中添加Label Analysis命令對各孔隙進(jìn)行單獨(dú)分析，通過統(tǒng)計(jì)各孔隙中的體素個數(shù)及邊界上的像素個數(shù)來表征孔隙的體積及表面積，并按孔隙大?。吹刃е睆剑┓謩e對孔隙數(shù)量、表面積、體積的分布進(jìn)行統(tǒng)計(jì)，結(jié)果如圖7所示。最終計(jì)算得到孔隙的體積為3 600.286 μm?3，表面積為26 334.324 μm?2，所以該煤巖的孔隙度為6.633%，比表面積為0.485 μm?-1.

2?基于MATLAB圖像處理的孔隙度及比表面積測定

由二維CT圖像，可以測出該圖像所對應(yīng)截面上的孔隙周長與面積。從圖4可以看出，將二維CT圖像沿y軸堆積即可得到煤巖的三維數(shù)字模型，各圖像上的孔隙周長與孔隙面積均可視為軸向坐標(biāo)y的函數(shù)。設(shè)試樣整體體積為V，孔隙體積為Vp，孔隙內(nèi)表面積為δ，垂直于截面方向的長度為l.再設(shè)截面上孔隙面積為Sp，孔隙周長為Pp.由此可得孔隙率為

式中?V，Vp分別為試樣整體體積和孔隙體積，μm3;δ為孔隙內(nèi)表面積，μm2;l為垂直于截面方向的長度，μm;Sp為再設(shè)截面上孔隙面積，μm2;Pp為孔隙周長，μm.

對于式（4）和式（5）中的積分部分，可以利用機(jī)械求積方法近似求解，將離散的CT圖像數(shù)據(jù)連續(xù)化，即

式中?yk為求積節(jié)點(diǎn);Ak為節(jié)點(diǎn)yk的權(quán)，僅與節(jié)點(diǎn)yk的選擇相關(guān)。因此，只需數(shù)張二維CT層析圖像，測出它們的孔隙周長與孔隙面積，即可通過數(shù)值積分的方法計(jì)算出試樣整體的孔隙率與比表面積。

為了方便與前一節(jié)的對比，同樣從編號200～800共601張圖像中選取了編號為200，275，350，425，500，575，650，725，800的二維CT圖像。由于原始圖像位深度為32位，MATLAB中的imcrop函數(shù)無法識別，故需將其導(dǎo)入畫圖、ACDSee等圖像編輯工具中轉(zhuǎn)存為位深度24位的bmp格式后再導(dǎo)入MATLAB。圖像在該過程中僅丟失了代表透明度的Alpha通道信息，對原始數(shù)據(jù)并無影響。再參照圖3剪裁其中邊長為601像素的正方形以便于分析，其中該正方形左上角像素在原圖中的坐標(biāo)為（240，200）。

CT圖像中的噪聲普遍存在。為了有效地抑制噪聲的干擾，必須對原始圖像進(jìn)行濾波以突出圖像中的有效信息。濾波主要可分為3類：線性濾波[21]、中值濾波[22-23]和自適應(yīng)濾波[24-25]。由于消除噪音的同時還需保持圖像細(xì)節(jié)清晰，且圖像中需統(tǒng)計(jì)的細(xì)節(jié)如點(diǎn)、線、尖頂部等較多，因此選擇自適應(yīng)濾波對原始圖像進(jìn)行處理。自適應(yīng)濾波是一種新型的信號處理方法，它是一種基于最小均方誤差準(zhǔn)則的最優(yōu)估計(jì)，對高斯白噪聲的去除效果尤其明顯。它既能一定程度上克服線性濾波后圖形細(xì)節(jié)模糊的問題，也能避免中值濾波導(dǎo)致的圖像細(xì)節(jié)缺失對統(tǒng)計(jì)結(jié)果的影響。

對原始圖像進(jìn)行濾波后，可以獲得更平滑的灰度直方圖信息。以編號200的CT圖像為例，其自適應(yīng)濾波前后的灰度直方圖變化如圖8所示。從圖8（a）可以看出，濾波前該灰度直方圖有3個波峰，但其中2個分別在左端（暗部）和右端（亮部）產(chǎn)生了溢出。而濾波后該灰度直方圖兩端溢出消失，波形雖有一定的偏移，但波峰和波谷的偏差不大，如圖8（b）所示。

當(dāng)基于雙峰直方圖的谷底閾值法推廣到三峰直方圖時，可以得到：三值化時若灰度直方圖具有3峰，3峰間的2個波谷則為2個閾值點(diǎn)[26]。因此選取了圖8（b）中2個波谷作為閾值點(diǎn)。這樣對CT圖像進(jìn)行三值化處理，不僅可以區(qū)分孔（裂）隙和固體基質(zhì)，還能區(qū)分固體基質(zhì)中的煤基質(zhì)和煤雜質(zhì)。第1個波谷為[30，80]區(qū)間（圖8（b）中紅色區(qū)間），在此區(qū)間中圖像可近似為水平線，其誤差不超過1%;第2個波谷則在225附近。因此，可將2個閾值點(diǎn)分別近似取為55和225，即灰度值在[0，55]的像素點(diǎn)代表孔（裂）隙，灰度值在[56，225]的像素點(diǎn)代表煤基質(zhì)，灰度值在[226，255]的像素點(diǎn)代表煤雜質(zhì)。

圖9（a）和（b）分別給出了編號200的CT原始切片圖像及其三值化結(jié)果。圖9（b）中黑色區(qū)域?yàn)榭祝眩┫?，灰色區(qū)域?yàn)槊夯|(zhì)，白色區(qū)域?yàn)槊弘s質(zhì)。其三值化后的灰度直方圖則如圖9（c）所示。但由于文中主要測量的是CT圖像的孔隙周長與孔隙面積，因此只需對原始圖像進(jìn)行二值化處理即可。按照灰度特性將圖像劃分為背景和目標(biāo)2部分，其中將孔（裂）隙視為目標(biāo)，煤基質(zhì)和煤雜質(zhì)等固體基質(zhì)視為背景。閾值分割點(diǎn)的灰度值取55，二值化后的圖像則如圖10所示。其中黑色像素點(diǎn)代表孔（裂）隙，白色像素點(diǎn)代表固體基質(zhì)。圖10黑白顏色對調(diào)后所有像素的像素值求和的結(jié)果即為編號200圖像的孔隙面積，單位為像素。

利用for循環(huán)語句對圖10中所有黑色像素點(diǎn)的4鄰域像素值進(jìn)行求和并映射到一個新的601×601矩陣中，該矩陣中的任一元素表示的是所對應(yīng)像素點(diǎn)在孔隙邊界上的長度，單位為像素邊長。該矩陣所有元素之和即為圖10中的孔隙周長。編號200，275，350，425，500，575，650，725，800的維CT圖像的孔隙周長與孔隙面積見表1，其中每個像素點(diǎn)邊長約為0.063 μm.

得到了各CT圖像的孔隙周長與孔隙面積，則可通過式（4）和式（5）計(jì)算煤樣的孔隙度和比表面積。式（6）中常用的機(jī)械求積公式有復(fù)合梯形公式、復(fù)合Simpson公式、復(fù)合Cotes公式等，分別記為

最終各數(shù)值積分方法計(jì)算出的孔隙度、比表面積見表2.

若以Avizo三維建模測定的孔隙度及比表面積為精確值，復(fù)合梯形公式、復(fù)合Simpson公式、復(fù)合Cotes公式計(jì)算出的孔隙度、比表面積誤差則見表3.從表3可以看出，各數(shù)值積分方法計(jì)算出的孔隙度和比表面積誤差均是復(fù)化Cotes公式最小，復(fù)化Simpson公式次之，復(fù)化梯形公式最大。各誤差均在10%左右，表明提出的測量方法能較好地估計(jì)煤樣的孔隙率和比表面積。

3?結(jié)?論

1）傳統(tǒng)的基于Avizo三維建模的測量方法精確度較高，同時還能獲取樣品的三維孔隙網(wǎng)絡(luò)模型和孔隙的形狀、走向及孔徑分布，但計(jì)算過程復(fù)雜，對于復(fù)數(shù)個試樣必須進(jìn)行多次計(jì)算，且對CT圖像有較高的質(zhì)量要求。

2）提出的基于MATLAB圖像處理的測量方法則通過簡單的數(shù)值積分實(shí)現(xiàn)了快速計(jì)算煤巖的孔隙度和比表面積，避免了三維建模，極大地降低了計(jì)算復(fù)雜度，且能夠一次性處理復(fù)數(shù)個試樣，更加簡便經(jīng)濟(jì)，較好地滿足了工程需要。且對于CT圖像具有更高的容錯率，每次計(jì)算僅需其中數(shù)張掃描質(zhì)量效果較好的CT圖像即可。

參考文獻(xiàn)（References）：

[1]Song L，Ning Z F，Duan L.Research on reservoir characteristics of Chang 7 tight oil based on nano-CT[J].Arabian Journal of Geosciences，2018，11（16）：472.

[2]徐?鵬，邱淑霞，姜舟婷，等.各向同性多孔介質(zhì)中Kozeny-Carman常數(shù)的分形分析[J].重慶大學(xué)學(xué)報(bào)，2011，34（4）：78-82.XV Peng，QIU Shu-xia，JIANG Zhou-ting，et al.Fractal analysis of Kozeny-Carman constant in the homogenous porous media[J].Journal of Chongqing University，2011，34（4）：78-82.

[3]Zhang R，Liu S，Bahadur J，et al.Changes in pore structure of coal caused by coal-to-gas bioconversion[J].Scientific Reports，2017，7（1）：3840.

[4]Rodrigues C F，Sousa M J L D.The measurement of coal porosity with different gases[J].International Journal of Coal Geology，2002，48（3）：245-251.

[5]Li X，Kang Y，Haghighi M.Investigation of pore size distributions of coals with different structures by nuclear magnetic resonance（NMR）and mercury intrusion porosimetry（MIP）[J].Measurement，2018，116：122-128.

[6]Liu S，Sang S，Wang G，et al.FIB-SEM and X-ray CT characterization of interconnected pores in high-rank coal formed from regional metamorphism[J].Journal of Petroleum Science and Engineering，2017，148：21-31.

[7]Karacan C O，Okandan E.Adsorption and gas transport in coal microstructure：investigation and evaluation by quantitative X-ray CT imaging[J].Fuel，2001，80（4）：509-520.

[8]宮偉力，李?晨.煤巖結(jié)構(gòu)多尺度各向異性特征的SEM圖像分析[J].巖石力學(xué)與工程學(xué)報(bào)，2010，29（4）：2681-2689.GONG Wei-li，LI Chen.Multi-scale and anisotropic characterization of coal structure based on SEM image analysis[J].Chinese Journal of Rock Mechanics & Engineering，2010，29（4）：2681-2689.

[9]Zhang G，Ranjith P G，Perera M S A，et al.Characterization of coal porosity and permeability evolution by demineralisation using image processing techniques：A micro-computed tomography study[J].Journal of Natural Gas Science & Engineering，2018，56（8）：384-396.

[10]Mayo S，Josh M，Kasperczyk D，et al.Dynamic micro-CT study of gas uptake in coal using Xe，Kr and CO2[J].Fuel，2018，212：140-150.

[11]王?剛，楊鑫祥，張孝強(qiáng)，等.基于CT三維重建的煤層氣非達(dá)西滲流數(shù)值模擬[J].煤炭學(xué)報(bào)，2016，41（4）：931-940.WANG Gang，YANG Xin-xiang，ZHANG Xiao-qiang，et al.Numerical simulation on non-Darcy seepage of CBM by means of 3D reconstruction based on computed tomography[J].Journal of China Coal Society，2016，41（4）：931-940.

[12]李?偉，要惠芳，劉鴻福，等.基于顯微CT的不同煤體結(jié)構(gòu)煤三維孔隙精細(xì)表征[J].煤炭學(xué)報(bào)，2014，39（6）：1127-1132.LI Wei，YAO Hui-fang，LIU Hong-fu，et al.Advanced characterization of three-dimensional pores in coals with different coal-body structure by Micro-CT[J].Journal of China Coal Society，2014，39（6）：1127-1132.

[13]Pan X C，Xia D，Zou Y，et al.A unified analysis of FBP-based algorithms in helical cone-beam and circular cone-and fan-beam scans[J].Physics and Medicine &Biology，2004，49（18）：4349-4369.[14]Liu X，Comtat C，Michel C，et al.Comparison of 3-D reconstruction with 3D-OSEM and with FORE+OSEM for PET[J].IEEE Transactions on Medical Imaging，2001，20（8）：804-814.

[15]David S，Burion S，Tepe A，et al.Experimental validation of an OSEM-type iterative reconstruction algorithm for inverse geometry computed tomography[C]//Medical Imaging 2012：Physics of Medical Imaging. International Society for Optics and Photonics，2012，8313（6）：125.[16]Herman G T，Meyer L B.Algebraic reconstruction techniques can be made computationally efficient[positron emission tomography application][J].IEEE Transactions on Medical Imaging，1993，12（3）：600.

[17]Peyrin F，Dong P，Pacureanu A，et al.Micro-and nano-CT for the study of bone ultrastructure[J].Current Osteoporosis Reports，2014，12（4）：465-474.

[18]劉郁紀(jì).X射線工業(yè)CT物理設(shè)計(jì)及圖像重建[D].蘭州：蘭州大學(xué)，2010.LIU Yu-ji.The physical design and image reconstruction of X-ray CT[D].Lanzhou：Lanzhou University，2010.[19]李?光，羅守華，顧?寧.Nano CT成像進(jìn)展[J].科學(xué)通報(bào)，2013，58（7）：501-509.LI Guang，LUO Shou-hua，GU Ning.Research progress of Nano CT imaging[J].Chinese Science Bulletin，2013，58（7）：501-509.

[20]陶?鵬.基于數(shù)字巖心的低滲儲層微觀滲流機(jī)理研究[D].成都：西南石油大學(xué)，2017.TAO Peng.Mechanism of micro seepage in low-permeability reservoirs based on digital core[D].Chengdu：Southwest Petroleum University，2017.

[21]Flamant J，Chainais P，Bihan N L.A complete framework for linear filtering of bivariate signals[J].IEEE Transactions on Signal Processing，2018，66（17）：4541-4552.

[22]孫宏琦，施維穎，巨永鋒.利用中值濾波進(jìn)行圖像處理[J].長安大學(xué)學(xué)報(bào)（自然科學(xué)版），2003，23（2）：104-106.SUN Hong-qi，SHI Wei-ying，JU Yong-feng.Image processing with medium value filter[J].Journal of Changan University（Natural Science Edition），2003，23（2）：104-106.

[23]Huang T，Yang G，Tang G.A fast two-dimensional median filtering algorithm[J].IEEE Transactions on Acoustics，Speech and Signal Processing，1979，27（1）：13-18.

[24]張旭明，徐濱士，董世運(yùn).用于圖像處理的自適應(yīng)中值濾波[J].計(jì)算機(jī)輔助設(shè)計(jì)與圖形學(xué)學(xué)報(bào)，2005，17（2）：295-299.ZHANG Xu-ming，XV Bin-shi，DONG Shi-yun.Adaptive median filtering for image processing[J].Journal of Computer-Aided Design & Computer Graphics，2005，17（2）：295-299.

[25]Althahab，Jumaah A Q.A new robust adaptive algorithm based adaptive filtering for noise cancellation[J].Analog Integrated Circuits and Signal Processing，2018，94（2）：217-231.

[26]鐘江城，周宏偉，任偉光，等.基于CT圖像灰度分布的含雜質(zhì)煤體三值化方法[J].力學(xué)與實(shí)踐，2018，40（2）：140-147.ZHONG Jiang-cheng，ZHOU Hong-wei，REN Wei-guang，et al.A three-value-segmentation method of coal containing inclusion based on gray distribution of computed tomography image[J].Mechanics in Engineering，2018，40（2）：140-147.