巖心孔隙結(jié)構(gòu)的空間分辨核磁表征

張 巖, 黃曉瑩

(中國石油大學(北京)地球物理學院,北京 102249)

2000年以來,拉普拉斯核磁共振技術(shù)(縱向弛豫時間T1、橫向弛豫時間T2以及擴散)迅速發(fā)展,在石油勘探領(lǐng)域得到廣泛應(yīng)用,有效獲取了巖石內(nèi)部結(jié)構(gòu)及流體信息[1-3]。但該方法無法精確定位復雜結(jié)構(gòu)巖石局部孔隙結(jié)構(gòu)及流體種類。核磁共振成像技術(shù)能夠快速、無損地對巖心進行可視化表征,并提取巖心內(nèi)部空間信息[4-8]。將以上兩種方法結(jié)合,可有效表征復雜孔隙介質(zhì)局域結(jié)構(gòu)及流體信息,該方法也被稱為空間定位核磁。Liu[9]通過小角度扳轉(zhuǎn)脈沖編碼快速獲得了T1信息,并將其與成像技術(shù)相結(jié)合,提取了樣品局部弛豫信息。Li等[10]將T2弛豫測量與一維相位編碼成像技術(shù)相結(jié)合,實驗獲取了巖心樣品的空間分辨的T2剖面。Xiao等[11]通過采取快速成像的方法,用相位編碼成像技術(shù)獲得二維空間分辨T2分布譜。Zhang等[12]通過采用適中分辨率的成像圖以及將成像過程和弛豫測量過程同時進行,大大縮短了實驗時間,獲得了二維空間分辨D-T2關(guān)聯(lián)譜。通過同樣的方法,Zhang等[13]獲得了二維空間分辨Gint2D-T2關(guān)聯(lián)譜,研究了樣品局部結(jié)構(gòu)與內(nèi)部梯度場、擴散和弛豫的關(guān)聯(lián)以及外加場強對實驗結(jié)果的影響。另外,Zhang等[14]在低場核磁儀器上測量了多孔材料的空間分辨的孔徑(a)-弛豫(T2)關(guān)聯(lián)譜,獲取了多孔樣品的空間分辨信息,提取了局部區(qū)域孔徑及表面弛豫弛豫率信息。之前的空間定位核磁多是在高場強核磁條件下進行的,但巖心中含有鐵磁性物質(zhì),由于鐵磁的存在,會影響磁場的場強分布,進而影響實驗結(jié)果,如在成像中出現(xiàn)重影、偽影等現(xiàn)象。所以在場強較高的情況下,對于天然巖心,不能很好地觀測巖心內(nèi)部結(jié)構(gòu)和流體的分布情況。在場強較低的情況下,又有信噪比低的問題。筆者在低場核磁條件下系統(tǒng)地進行空間定位核磁研究,并通過將多個像素點相加以獲取更高的局部區(qū)域信號強度的方法,以提升實驗的信噪比。同時,對非均質(zhì)性進行系統(tǒng)研究。

1 基本原理與數(shù)據(jù)處理方法

1.1 核磁成像及與拉普拉斯核磁聯(lián)合的基本原理

在巖石物理研究中,磁共振成像(MRI)具有重要意義,可以無損表征多孔介質(zhì)內(nèi)部結(jié)構(gòu)以及流體的分布情況。Lauterbur和Mansfield完成首個MRI實驗[15-16]。通過將磁場梯度G疊加到主磁場B0上,產(chǎn)生了空間相關(guān)磁場,其拉莫爾頻率和空間位置的關(guān)系可以描述為ω0(r)=γ(B0+G·r),其中γ為旋磁比。在樣品的體積元dV之內(nèi),位置r的NMR信號元dS可以描述[17]為

dS(G,t)=C(r)ρ(r)dVexp[i(γB+γG·r)t].

(1)

式中,ρ(r)為局部自旋密度;C(r)為由其他核磁共振參數(shù),如T1、T2和擴散系數(shù)等引起的對比因子。

如果梯度足夠強,橫向磁化矢量的散相將由梯度主導,可忽略T2效應(yīng)。當NMR信號“共振”時,上述積分可以表述為

(2)

該方程式體現(xiàn)了一種傅里葉變換。為了清楚說明這一點,Mansfield引入了波矢量k=γGt/2π。式(2)可以改寫為

(3)

(4)

因此,用于成像的信號采集是k空間采樣[18]。通過將成像信號進行傅里葉變換,即可得到成像圖。

在此基礎(chǔ)上,增加測量拉普拉斯核磁的脈沖序列,通過系統(tǒng)性改變回波間隔、恢復時間及梯度強度等,完成兩種方法的聯(lián)合,其信號可由如下公式表示:

S(tE,td)=S(0)?f(T1,T2)exp(-td/12)×

exp(-ntE/T2)dT1dT2ρ(r)exp(-kr)dr.

(5)

式中,f(T1,T2)為空間分辨拉普拉斯譜;tE為回波間隔;td為恢復時間;S(0)為初始信號;S(tE,td)為衰減后的信號。

1.2 非均質(zhì)性定量表征數(shù)據(jù)處理方法

首先,引入了變異函數(shù)的概念[19-20],設(shè)區(qū)域化變量Z(x)滿足二階平穩(wěn)條件或本征假設(shè),則變異函數(shù)的計算表達式為

(6)

式中,h為兩個樣本點的空間分隔距離;Z(xi)和Z(xi+h)分別為空間位置xi和xi+h處的觀測值;N(h)為分隔距離為h時的樣本對數(shù)。計算的變異函數(shù)值實質(zhì)上反映的是分隔距離為h的兩個位置測得的參數(shù)之間的平均平方差,有了采樣數(shù)據(jù)及變異函數(shù)計算公式就可以獲知距離h的區(qū)域化變量變異性。

變異函數(shù)γ(h)是一個單調(diào)遞增的函數(shù),當h超過某一數(shù)值后,γ(h)不再繼續(xù)單調(diào)地增大,而往往穩(wěn)定在一個極限值γ(∞)附近。設(shè)Z(x)具有各向同性的變異函數(shù)γ(h),選用指數(shù)模型對變異函數(shù)γ(h)隨h變化的比例關(guān)系進行擬合,擬合公式為

γ(h)=C0+C1(1-exp(-h/a)).

(7)

當h=0時,變異函數(shù)γ(h)≠0,而等于一個常數(shù)C0,這種現(xiàn)象稱為“塊金效應(yīng)”。其中C0稱為塊金常數(shù),反映區(qū)域化變量在小于抽樣尺度h時所具有的變異性。(C0+C1)為基臺值,C1為拱高。

2 實 驗

實驗中共測量3塊人造巖心:巖心N孔隙度為10%,滲透率為10×10-3μm2;巖心P孔隙度為18%,滲透率為30×10-3μm2;巖心B是一塊人造砂巖,孔隙度為25%,滲透率為3 μm2。實驗前用飽和儀將巖心抽真空,然后加壓飽和3 d,這里所用的飽和流體為NaCl溶液。

實驗脈沖序列如圖1所示。選擇a序列對巖心N進行T2-MRI實驗,設(shè)置好相應(yīng)的圖像參數(shù)和序列參數(shù),每次只改變回波間隔tE參數(shù),第一個tE值設(shè)置為儀器允許的最小值,tE=6 ms,之后按6 ms的1～10倍逐次變化:tE=6, 12, 18, 24, 30, 36, 42, 48, 54, 60 ms,共10組實驗,進而得到T2衰減曲線。設(shè)置采樣的累加次數(shù)為32。

圖1 空間分辨脈沖序列

在T1-MRI實驗中使用反轉(zhuǎn)恢復成像脈沖序列(b)對飽和鹽水的巖心進行成像實驗。設(shè)置采樣的累加次數(shù)為32,回波間隔tE=6 ms,通過改變恢復時間(td),對樣品進行多組成像實驗,進而得到T1恢復曲線。

如圖2所示,左圖為巖心N的矢狀面定位像,右圖為巖心N的軸向面定位像,5條紅色長條表示選取該巖心的5個圓截面層面進行成像。

圖2 巖心N的矢狀面與冠狀面定位像

3 實驗結(jié)果與討論

3.1 巖心N非均質(zhì)性的定性分析

對成像數(shù)據(jù)進行傅里葉變換以得到每個層面的圖像,圖3為tE=6 ms的情況下巖心N的5個層面的成像結(jié)果(分辨率為1.6 mm,像素點為32×32)。根據(jù)成像圖,可以從宏觀尺度上得到巖心含水量分布、孔隙分布及非均質(zhì)信息。

圖3 tE=6 ms時巖心N的5個層面的成像

從圖3中可以看出,巖心N整體的非均質(zhì)性較強,第1層信號較弱,孔隙最少,含水量最少;第2層信號較強,孔隙最多,含水量最多;第1、2、3層的右半部分比左半部分的孔隙多,含水量相對也多;第4、5層的孔隙集中分布在層面的中部區(qū)域,且從宏觀上看這兩個層面的孔隙分布和含水量分布情況最為相似。

在此基礎(chǔ)上,將拉普拉斯核磁共振方法和核磁共振成像方法結(jié)合起來,對巖心N每個層面進行空間分辨的定位分析,進而得到樣品內(nèi)部的詳細結(jié)構(gòu)及流體分布信息。具體程序如下所述:

在巖心N成像實驗測量的5個層面上各選取12個點,如圖3所示,坐標分別為

x=[17 20 23 17 20 23 17 20 23 17 20 23];

y=[13 13 13 16 16 16 19 19 19 22 22 22].

然后將每個點及周圍8個點的信號幅度相加,得到12個區(qū)域的信號幅度,由此有效增加了信噪比。將每個區(qū)域不同tE時間的信號幅度連成一條衰減曲線,共可得到的12條衰減曲線。經(jīng)過T2反演之后,得到每個區(qū)域的T2分布。同樣,通過反轉(zhuǎn)恢復成像脈沖序列對巖心N進行T1-MRI實驗,共可得到12條恢復曲線,得到每個區(qū)域的T1分布。

如圖4所示:

圖4 巖心N的5層小區(qū)域與整塊巖心的T1和T2分布

(1)第1、2層都是同一列的3個區(qū)域的孔隙結(jié)構(gòu)較為相似,最左邊的3個區(qū)域孔隙相對較少,中部區(qū)域孔隙較多,前3列區(qū)域孔徑分布較為一致,而最后一列的3個區(qū)域較其他區(qū)域來說,孔徑分布范圍更廣,含有一些孔徑稍大的孔隙。

(2)第3層中前兩列區(qū)域孔徑分布較為一致,最左邊的3個區(qū)域孔隙相對較少,區(qū)域4～6孔隙較多;區(qū)域7～9的孔徑相對其他區(qū)域要小一些,而最后一列的3個區(qū)域較其他區(qū)域來說,孔徑分布范圍更廣,含有一些孔徑稍大的孔隙。

(3)第4層和第5層的孔隙結(jié)構(gòu)和流體分布的特點較為相似,但與第1、2、3層的差別較大,孔隙主要集中分布在中部區(qū)域,左右兩邊的區(qū)域孔隙相對較少,前兩列區(qū)域孔徑分布較為一致,區(qū)域2、3的孔徑較前兩列的其他區(qū)域來說略大。區(qū)域7～9的孔徑相對其他區(qū)域要小一些,而最后一列的3個區(qū)域較其他區(qū)域來說,孔徑分布范圍更廣,含有一些孔徑稍大的孔隙。

可以發(fā)現(xiàn),同一層面各局部區(qū)域的T1、T2分布有差異,層面與層面之間同坐標區(qū)域的T1、T2分布也有較大差異,局部區(qū)域的T1、T2分布與整塊巖心的T1、T2分布也有較大差異,由此可以定性地表征巖心N的非均質(zhì)性。

3.2 巖心P非均質(zhì)性的定性分析

對巖心P進行定位分析。如圖5所示(分辨率為1.6 mm,像素點為32×32),巖心P整體的非均質(zhì)性較強,孔隙主要集中分布在各層面的中部區(qū)域和左半邊區(qū)域,這些區(qū)域流體含量較多,信號較強,而各層面的右半邊區(qū)域孔隙相對較少,流體含量較少,信號較弱。

在巖心P成像實驗測量的5個層面上各選取9個點,如圖5所示,坐標分別為

圖5 tE=6 ms時巖心P的5個層面的成像

x=[14 14 14 17 17 17 20 20 20];

y=[17 20 23 17 20 23 17 20 23].

然后將每個點及周圍8個點的信號幅度相加,得到了9個區(qū)域的信號幅度,由此有效增加了信噪比。將每個區(qū)域不同tE時間的信號幅度連成一條衰減曲線,共可得到的9條衰減曲線。經(jīng)過T2反演之后,得到每個區(qū)域的T2分布。同樣,通過反轉(zhuǎn)恢復成像脈沖序列對巖心P進行T1-MRI實驗,共得到9條恢復曲線及每個區(qū)域的T1分布(圖6)。

從圖6中可以看出:

圖6 巖心P的5個層面的小區(qū)域與整塊巖心的T1和T2分布

(1)第1層。在該層面上選取的3列區(qū)域中,同一列的小區(qū)域的孔隙結(jié)構(gòu)和流體含量較為相似,不同列區(qū)域間存在較大差異,中部3個小區(qū)域孔隙最多,含有較多大孔隙,而左右兩邊區(qū)域的大孔隙較少,右邊區(qū)域孔隙最少。但該層面各區(qū)域的孔徑分布范圍基本一致。

(2)第2層。在該層面上同一列區(qū)域的孔隙結(jié)構(gòu)和流體含量較為相似,不同列區(qū)域間存在較大差異,中部3個區(qū)域孔隙最多,含有較多大孔隙;而左右兩邊區(qū)域的大孔隙較少,但與第1層同坐標區(qū)域相比,含有的大孔隙略多一些。但該層面各區(qū)域的孔徑分布范圍基本一致。

(3)第3層。該層面的中部區(qū)域孔隙最多,左邊區(qū)域次之,右邊區(qū)域孔隙最少,中部3個小區(qū)域含有的大孔隙最多,區(qū)域3、6、9含有一些其他區(qū)域沒有的孔徑很小的孔隙。與第1、2層相比,第3層的中部區(qū)域含有的大孔隙要少一些。

(4)第4層。該層面各區(qū)域的孔徑分布范圍有較大差異,左邊區(qū)域的孔徑分布范圍最窄,右邊區(qū)域的孔徑分布范圍稍微大一些,而中部區(qū)域的孔徑分布范圍最廣,含有一些其他區(qū)域沒有的孔徑較大的孔隙,且總的來說中部區(qū)域含有的大孔隙也最多。

(5)第5層。相對于其他層面來說,該層面各小區(qū)域間的T1和T2分布的差異較小,各區(qū)域的孔徑分布范圍基本一致,但區(qū)域8孔徑分布范圍較其他區(qū)域要小一些。相對于左邊區(qū)域和中部區(qū)域來說,右邊區(qū)域的孔隙相對少一些,特別是孔徑較小的孔隙含量明顯要少一些。

總的來說,巖心P的非均質(zhì)性較強,同一層面各局部區(qū)域的T1、T2分布有較大差異,層面與層面之間同坐標區(qū)域的T1、T2分布也有所差異,局部區(qū)域的T1、T2分布與整塊巖心的T1、T2分布也有較大差異,由此可以定性地表征巖心P的非均質(zhì)性。

3.3 巖心B非均質(zhì)性的定性分析

對成像原始數(shù)據(jù)進行傅里葉變換以得到每個層面的圖像。圖7為tE=6 ms的情況下巖心B的5個層面的成像結(jié)果(分辨率為1.6 mm,像素點為32×32)。從宏觀尺度上可看出,巖心B整體較為均質(zhì),各層面的孔隙結(jié)構(gòu)和流體分布很相似,各層面內(nèi)的孔隙數(shù)量和含水量都很多;而相對于其他層面來說,第3層的孔隙數(shù)量和含水量略少一些。

圖7 tE=6 ms時巖心B的5個層面的成像

圖8所示為巖心B的5個層面上12個小區(qū)域與整塊巖心T1、T2的分布?？梢钥闯?巖心B很均質(zhì),各層面不同區(qū)域的孔徑分布范圍基本一致,同一層面各局部區(qū)域的T1、T2分布的差異很小,各層面間同坐標區(qū)域的T1、T2分布的差異也很小,各局部區(qū)域與整塊巖心的T1、T2分布也較為相似,由此可以定性地看出巖心B的非均質(zhì)性很小。

圖8 巖心B的5個層面的小區(qū)域與整塊巖心的T1和T2分布

3.4 砂巖樣品非均質(zhì)性的定量表征

引入地質(zhì)統(tǒng)計學中的理論方法對巖心非均質(zhì)性進行定量分析。將變異函數(shù)的概念引入到對成像數(shù)據(jù)的處理分析中,可以實現(xiàn)對巖心樣品非均質(zhì)性的定量表征。選用多層自旋回波成像實驗得到的成像數(shù)據(jù)進行處理,h為兩個樣本點之間的距離,Z(x)為各樣本點的信號強度,可以得到各巖心的各層面對應(yīng)的變異函數(shù)γ(h)隨h變化的比例關(guān)系,然后對其進行指數(shù)擬合,得到各層面非均質(zhì)特征值a,用a表征巖心各個層面的非均質(zhì)性。圖9所示為巖心N、P和B的5個不同層面分別對應(yīng)的變異函數(shù)γ(h)隨h變化的比例關(guān)系及指數(shù)擬合曲線。

圖9 巖心B、N和P的變異函數(shù)與h的關(guān)系及指數(shù)擬合曲線

表1所示為巖心B、N和P的5個不同層面分別對應(yīng)的a值。可以看出,巖心B各層面的非均質(zhì)性都較小,第2層和第4層相對于其他層面非均質(zhì)性略大。巖心N各層面的非均質(zhì)性均比巖心B的大,其中第2層相對于其他層面非均質(zhì)性更顯著(a值增大了一倍);而相對于巖心B和巖心N,巖心P非均質(zhì)性很強,其中第1、2、3層的非均質(zhì)性很大,而第4層和第5層的非均質(zhì)性略小一些。該方法可能的誤差主要來源于兩個方面:

表1 巖心B、N和P的5個不同層面的a值

①實驗誤差。包括信號較弱(如低孔低滲樣品),成像選層重疊等。通過將多個像素點疊加,增加信噪比,通過在成像選層之間有所間隔,避免重疊,盡量降低實驗誤差。

②擬合誤差。通過選取合適的擬合方法,將擬合誤差降到最低。該方法的優(yōu)勢在于能夠定量評價樣品的非均質(zhì)性,缺點在于對樣品信噪比有一定要求,對于信號太弱的樣品,可能因成像效果不佳,無法進行有效表征。

4 結(jié) 論

(1) 非均質(zhì)性定性表征發(fā)現(xiàn),巖心B較均質(zhì),孔隙尺寸分布較窄,巖心P、 N非均質(zhì)性較強,孔隙尺寸分布較大。由成像圖可知,巖心N左邊含水量較小,孔隙較小,而右邊含水量較大,孔隙較大。巖心P總體與巖心N相反,巖心B信號和孔隙分布較均勻。由局域T1、T2圖可知,巖心N與P在同一層及不同層之間顯示出明顯的非均質(zhì)性,表明其孔隙結(jié)構(gòu)在不同區(qū)域差異較大。

(2) 非均質(zhì)性定量表征發(fā)現(xiàn),巖心P非均質(zhì)性最強,巖心B最均質(zhì)。在不同層面上,巖心N第2層非均質(zhì)較強,其他層均質(zhì)性較弱;巖心P第1～3層非均質(zhì)較強,4～5層均質(zhì)性較弱;巖心B各層之間非均質(zhì)性差別較小。非均質(zhì)性定量與定性表征結(jié)果總體一致性較好。