林婷婷，慧 芳，蔣川東，林 君

吉林大學(xué)地球信息探測(cè)儀器教育部重點(diǎn)實(shí)驗(yàn)室／儀器科學(xué)與電氣工程學(xué)院，長(zhǎng)春 130021

1 引 言

磁共振測(cè)深（Magnetic Resonance Sounding，簡(jiǎn)稱MRS）是一種應(yīng)用核磁共振現(xiàn)象探測(cè)地下水的地球物理方法［1－2］，與其他物探方法相比具有非侵入性和能夠直接定量探測(cè)地下水的特點(diǎn).在磁共振測(cè)深方法［3－4］中，地下水中氫質(zhì)子在地磁場(chǎng)的環(huán)境下處于平衡狀態(tài)，在地面鋪設(shè)線圈中通入當(dāng)?shù)乩獱栴l率的交變電流產(chǎn)生激發(fā)磁場(chǎng)，氫質(zhì)子受磁場(chǎng)激發(fā)產(chǎn)生能級(jí)躍遷，激發(fā)停止后在線圈中接收到氫質(zhì)子由激發(fā)態(tài)恢復(fù)到平衡態(tài)過程中產(chǎn)生的MRS信號(hào).

地面MRS多弛豫時(shí)間反演起步于2005年.Mohnke［10］提出了時(shí)間步長(zhǎng)反演方法，利用信號(hào)隨時(shí)間變化的測(cè)深曲線，得到不同時(shí)間下含水量與深度的關(guān)系，實(shí)現(xiàn)了多弛豫時(shí)間反演.但該方法過度依賴于含水量時(shí)間曲線，結(jié)果極不穩(wěn)定.Mike等［11］考慮反演模型垂向平滑問題，研究了一種穩(wěn)定的QT算法，不過由于平滑約束的固有特點(diǎn)，出現(xiàn)含水層邊界面刻畫模糊，反演結(jié)果不清晰、迭代耗時(shí)長(zhǎng)等缺陷.

本文針對(duì)目前國(guó)內(nèi)外磁共振反演方法的不足，結(jié)合MRS信號(hào)整體性與多指數(shù)特性，借鑒分層反演的思想，首次提出并實(shí)現(xiàn)了磁共振全弛豫數(shù)據(jù)分層多指數(shù)反演方法.采用多指數(shù)擬合方式將MRS信號(hào)所有采集數(shù)據(jù)一次性帶入反演，直接建立起地下探測(cè)深度、弛豫時(shí)間、含水量三個(gè)重要參數(shù)的關(guān)系.該方法通過構(gòu)建泛函極小值方程，利用不等式約束的空間信賴域算法進(jìn)行優(yōu)化，可快速、精確的反演獲得各含水層位置、弛豫時(shí)間及其對(duì)應(yīng)含水量等信息.

為全面驗(yàn)證本文提出的反演方法，文中首先建立了多含水層單一孔隙和多含水層多孔隙兩種不同的地下水模型.反演分別采用平滑、分塊和分層多指數(shù)三種方法依次進(jìn)行，并從模型逼近程度，目標(biāo)值擬合誤差等方面展開對(duì)比分析.最后，筆者對(duì)白城四方坨子測(cè)試點(diǎn)的磁共振野外實(shí)測(cè)數(shù)據(jù)進(jìn)行反演解釋，估算了涌水量［12－13］的大小，與鉆井資料進(jìn)行一致性比對(duì).本研究通過模型數(shù)據(jù)與實(shí)測(cè)數(shù)據(jù)，充分證實(shí)了該反演方法的優(yōu)越性.

2 磁共振數(shù)據(jù)反演方法

2.1 磁共振信號(hào)初始振幅反演

磁共振是具有非零磁矩的原子核在交變磁場(chǎng)的作用下發(fā)生共振吸收的現(xiàn)象.計(jì)算地下水中氫原子核的磁共振角頻率ω0公式為

其中，B0是測(cè)試地點(diǎn)的靜地磁場(chǎng)；γp＝0.267518 Hz／nT是氫質(zhì)子的旋磁比；f0是當(dāng)?shù)貧滟|(zhì)子的拉莫爾頻率［14］.為測(cè)量MRS信號(hào)，在測(cè)試地點(diǎn)鋪設(shè)線圈，線圈的直徑大小取決于所要探測(cè)的最大深度.線圈中通入當(dāng)?shù)乩獱栴l率的交變電流，在地下形成激發(fā)磁場(chǎng).公式為

其中，I0為電流強(qiáng)度，激發(fā)持續(xù)時(shí)間0＜t＜τ，由脈沖矩q＝I0τ體現(xiàn)激發(fā)強(qiáng)度.激發(fā)停止后，地面接收線圈感應(yīng)到磁共振響應(yīng)信號(hào)，方程為

其中，ωr和φr分別代表接收信號(hào)的角頻率和相位；平均弛豫時(shí)間T＊2反映地下平均有效孔隙大??；信號(hào)的初始振幅E0表示在激發(fā)停止時(shí)刻磁共振信號(hào)的幅值大小，它與地下含水量直接相關(guān).磁共振信號(hào)的包絡(luò)稱為自由感應(yīng)衰減（Free Induction Decay，簡(jiǎn)稱FID）信號(hào)，呈指數(shù)規(guī)律衰減.Weichman［15］等人提出通過方程（4）描述FID信號(hào)，建立信號(hào)初始振幅與含水量的關(guān)系：

其中，L代表了地下水探測(cè)的深度范圍；K1D（q，z）是一維核函數(shù)［16］，描述不同脈沖矩下不同深度處核磁共振信號(hào)的響應(yīng)，反映了激發(fā)磁場(chǎng)的穿透能力及分辨率，體現(xiàn)地下磁共振響應(yīng)的靈敏度，如圖1所示，橫軸代表激發(fā)脈沖矩，縱軸代表地下深度，顏色代表了響應(yīng)強(qiáng)度；M0是氫原子核的磁化強(qiáng)度，代表平衡條件下單位體積內(nèi)的磁矩；B⊥是激發(fā)磁場(chǎng)垂直于地磁場(chǎng)的分量；0≤m（z）≤1代表地下含水量情況.按深度z將方程（4）離散化，得到了初始振幅與含水量之間的矩陣關(guān)系為

平滑反演方法是將地下從地表至最大探測(cè)深度通過等間距或遞增間距的方式使地下含水層呈層狀分布，公式為

采用吉洪諾夫正則化［17－18］方法計(jì)算每個(gè)含水層對(duì)應(yīng)的含水量m（zi）及平均弛豫時(shí)間（zi），i＝1，…，b.其吉洪諾夫泛函極小值方程為

分塊反演的整個(gè)過程與平滑反演類似，不同之處在平滑反演預(yù)設(shè)地下含水層按深度z呈層狀分布，而分塊反演是已知地下含水層個(gè)數(shù)，各含水層深度未知，需要反演獲得，其結(jié)果使得各含水層之間分界明顯，易于判斷含水層位置.

圖1 傳統(tǒng)一維核函數(shù)圖Fig.1 Traditional one dimensional kernel matrix

2.2 磁共振弛豫信號(hào)的多指數(shù)描述

目前，磁共振平滑反演與分塊反演方法通過平均弛豫時(shí)間僅將地下環(huán)境假想為單一孔隙.然而，核磁共振測(cè)井及實(shí)驗(yàn)室核磁共振［19］技術(shù)指出，每一地層中探測(cè)到的MRS信號(hào)可視為不同孔隙產(chǎn)生的弛豫信號(hào)疊加.根據(jù)磁共振表面弛豫擴(kuò)散理論，孔隙內(nèi)流體的弛豫時(shí)間與孔隙大小以及形狀密切相關(guān)，較小孔隙對(duì)應(yīng)信號(hào)的弛豫時(shí)間短，反之，較大孔隙對(duì)應(yīng)信號(hào)的弛豫時(shí)間長(zhǎng)［20］.由上述分析可知，地面磁共振信號(hào)應(yīng)具有多指數(shù)特性，而多指數(shù)FID信號(hào)可通過方程（10）表述：

其中，T2i（i＝1，2，…n）表示地下環(huán)境中第i類孔隙物質(zhì)對(duì)應(yīng)的弛豫時(shí)間，Ei表示第i類弛豫單元對(duì)總MRS信號(hào)的貢獻(xiàn).

分析多指數(shù)方式描述磁共振信號(hào)對(duì)獲取信號(hào)初始振幅的意義.圖2直觀給出單指數(shù)方式和多指數(shù)方式針對(duì)含有噪聲的FID仿真信號(hào)的擬合情況，黑色點(diǎn)即以信噪比以10dB為例的FID信號(hào)；空心圓圈為初始振幅值；黑色及灰色實(shí)線分別為采用單指數(shù)擬合及多指數(shù)擬合FID信號(hào)的結(jié)果.由圖可見，單指數(shù)擬合后得到的信號(hào)初始振幅低于模型設(shè)定值；而多指數(shù)方式獲得的初始振幅與模型設(shè)定值吻合.從趨勢(shì)而言，采用多指數(shù)擬合能夠更好的逼近FID信號(hào)整個(gè)衰減過程.

圖2 磁共振FID仿真信號(hào)的單指數(shù)和多指數(shù)擬合對(duì)比Fig.2 The comparison of mono－exponential and multi－exponential fitting of FID model data

表1分析了不同信噪比下分別采用單一指數(shù)擬合和多指數(shù)擬合的初始振幅擬合誤差.通過對(duì)比可以發(fā)現(xiàn)，在不同信噪比下，多指數(shù)擬合信號(hào)初始振幅的準(zhǔn)確度均高于單一指數(shù)擬合，尤其在低信噪比（1～2dB）的情況下仍能比較理想的進(jìn)行初始振幅擬合.初始振幅的大小直接決定了反演的準(zhǔn)確性，因此，采用多指數(shù)方式描述磁共振信號(hào)對(duì)反演解釋將具有重要意義.

2.3 分層多指數(shù)反演方法

由前所述，平滑反演和分塊反演方法均通過初始振幅－脈沖矩曲線實(shí)現(xiàn)地下水垂向分布反演，忽略了大部分弛豫信號(hào)，導(dǎo)致各脈沖矩下信息獨(dú)立，造成較大反演偏差.考慮到核函數(shù)的積分特性，即信號(hào)在相鄰時(shí)間相鄰脈沖矩下的信息并非完全獨(dú)立，本文由此提出了分層模式與多指數(shù)結(jié)合的磁共振全弛豫信號(hào)反演方法，能夠直接建立地下含水層位置、各層含水量大小及其映射弛豫時(shí)間的關(guān)系.

表1 不同信噪比下初始振幅擬合誤差Table 1 Initial amplitude fitting error in different SNR

2.3.1 弛豫信號(hào)多指數(shù)的正演推導(dǎo)

為實(shí)現(xiàn)弛豫信號(hào)的反演方法，需要重新推導(dǎo)磁共振弛豫信號(hào)的正演方程，建立弛豫信號(hào)與地下含水量之間的關(guān)系.公式為

其中，mi（T2i，z）＞0，i＝1，2，…，n代表地層z處不同弛豫時(shí)間 T2i對(duì)應(yīng)的含水量大小.K1D（q，z，t，T2i）是本文建立的弛豫信號(hào)的一維核函數(shù)矩陣，宏觀上說，新核函數(shù)通過不同弛豫時(shí)間T2的e指數(shù)因子刻畫了在多指數(shù)前提下脈沖矩隨著時(shí)間變化的磁共振響應(yīng).根據(jù)弛豫時(shí)間與含水地層對(duì)應(yīng)關(guān)系［21］，確定30ms＜T2＜1000ms，在此范圍內(nèi)選取合適的T2譜分布作為先驗(yàn)信息.

為使該正演過程更加清晰，筆者在圖3繪制了新的核函數(shù)，橫軸代表弛豫時(shí)間T2，縱軸代表信號(hào)時(shí)間t.m行n列的子核函數(shù)矩陣能夠精細(xì)描述出不同信號(hào)采集時(shí)間下、不同弛豫時(shí)間T2對(duì)應(yīng)的隨深度及脈沖矩變化的MRS響應(yīng).圖3中可見，隨著時(shí)間增加，同一T2下對(duì)MRS的響應(yīng)逐漸衰減，小T2對(duì)應(yīng)的響應(yīng)衰減快，大T2對(duì)應(yīng)的響應(yīng)衰減慢，這與本文2.2中的磁共振表面弛豫擴(kuò)散理論吻合.

將方程（11）離散化得到多指數(shù)下磁共振FID信號(hào)與含水量的關(guān)系，用矩陣形式表示為

2.3.2 分層控制下的目標(biāo)函數(shù)

磁共振分層多指數(shù)反演之目的是要獲得從地表至最大探測(cè)深度的地下水分布情況及賦存狀態(tài)，即求解方程（12）中的mT2，z，對(duì)該矩陣的求解過程也是決定反演收斂速度的關(guān)鍵.由于采用全弛豫信號(hào)進(jìn)行反演，涉及的數(shù)據(jù)量龐大，同時(shí)為避免光滑連續(xù)模型中出現(xiàn)的迭代緩慢及水層分界模糊現(xiàn)象，文中采取了分層控制方式突出含水邊界，經(jīng)有限步的計(jì)算就可以得到精確解.

將待反演的含水層個(gè)數(shù)設(shè)為控制量，分層數(shù)為k.若地下實(shí)際存在的含水層個(gè)數(shù)為θ，則地層分層數(shù)至少為2θ＋1.僅當(dāng)含水層個(gè)數(shù)小于激發(fā)脈沖矩個(gè)數(shù)時(shí)能獲取較好的反演效果，固而確定k的取值范圍為2θ＋1≤k＜g，g代表激發(fā)脈沖矩個(gè)數(shù).對(duì)含水層層頂深度zt和層底深度zb求解，就能夠確定含水層的位置及厚度，它們是構(gòu)成本文分層反演的核心，方程（13）給出了分層多指數(shù)反演中含水量與深度及弛豫時(shí)間的關(guān)系：

2.3.3 空間信賴域反演

人們求解含水量反演這類非線性規(guī)劃問題的方法主要有擬牛頓法和共軛梯度法等.但這些方法過于依賴初值選取，有時(shí)還會(huì)限于局部極值，無法找到全局極值.鑒于傳統(tǒng)方法的缺點(diǎn)，我們引入空間信賴域算法［21－22］，這是現(xiàn)代應(yīng)用較廣泛的地球物理反演方法，基本思想在于：通過一系列信賴域子問題的最優(yōu)值逼近最優(yōu)化問題的解.該方法穩(wěn)定的數(shù)值性能，適合于求解病態(tài)最優(yōu)化問題，并具備全局收斂性.對(duì)于無約束或等式約束的優(yōu)化問題，該方法已被證實(shí)是可靠而強(qiáng)適的.筆者通過構(gòu)造信賴域子問題，得到不等式約束優(yōu)化問題的信賴域算法.

建立目標(biāo)函數(shù)的最小二乘形式，方程為

其中，grad（J）d和 Hess（J）d分別代表泛函J 在迭代點(diǎn)～md的梯度和Hessian陣，分別記為

（1）不等式約束的施加

將不等式約束施加于反演方程能夠有效的約束信賴域空間，提高反演效率.

（2）正則化參數(shù)的確定

圖3 分層多指數(shù)一維核函數(shù)圖Fig.3 Layered multi－exponential one dimensional kernel matrix

其中，I為單位矩陣.為保證算法順利進(jìn)行，調(diào)節(jié)正則化參數(shù)α使Hess（J）d＋αI是正定的，則ξ唯一的由方程（21）決定，有

然后判斷ξd是否落在信賴域半徑Δd內(nèi)，滿足‖ξd‖≤Δd則接受該點(diǎn)，并轉(zhuǎn)入下一步迭代，否則按αd＋1＝αd·γ（γ＞1）調(diào)節(jié)參數(shù)，使每次迭代的目標(biāo)函數(shù)為一單調(diào)下降序列，直至 （15）式迭代收斂到一個(gè)穩(wěn)定數(shù)值解.

圖4為分層多指數(shù)反演流程圖.

圖5是本文歸納的平滑反演、分塊反演和分層多指數(shù)反演的總體流程框圖.

3 模型對(duì)比

3.1 模型實(shí)驗(yàn)一

構(gòu)造單一弛豫時(shí)間T2多含水層模型，分別設(shè)置含水層弛豫時(shí)間及含水量，比對(duì)平滑、分塊、分層多指數(shù)三種方法反演效果.

模型數(shù)據(jù)是在如下條件下通過正演計(jì)算獲得：

（1）邊長(zhǎng)100m的單匝方形線圈，采用收發(fā)共線方式；

（2）假設(shè)地下為均勻半空間，電阻率為500μm.地磁場(chǎng)強(qiáng)度B0＝46970nT，地磁傾角60°；

圖6 模型一Fig.6 Model one

（3）拉莫爾頻率為2000Hz，激發(fā)脈沖矩從0.4As到10As分為20次激發(fā)地下0m到100m；

（4）假設(shè)地下存在三個(gè)含水層（圖6），第一層設(shè)置于地下淺層5m至15m，含水量為10%，弛豫時(shí)間T2為100ms；第二層設(shè)置于地下中層25m到35m，含水量為10%，弛豫時(shí)間T2為300ms；第三層設(shè)置于地下深層55m到75m，含水量為10%，弛豫時(shí)間T2為200ms；

（5）由于磁共振方法采集信號(hào)水平在納伏級(jí)，易受到噪聲的干擾，導(dǎo)致反演精度下降.因此，為使仿真結(jié)果接近實(shí)際并驗(yàn)證該方法的適用性，在正演模型中加入5%Guass白噪聲.設(shè)置反演地層個(gè)數(shù)k＝10；初始正則化參數(shù)α＝100.

圖7給出了模型一的反演結(jié)果.（a）為含噪磁共振數(shù)據(jù)，藍(lán)色實(shí)線為不同脈沖矩下的信號(hào)衰減曲線，紅色圈為不同脈沖矩下接收時(shí)刻信號(hào)的幅值E1（td）.圖（b）、（c）中左圖分別為平滑反演和分塊反演得到的含水量垂向分布，右圖為平均弛豫時(shí)間垂向分布，紅色虛線為模型值，藍(lán)色實(shí)線為反演值.（d）中左圖為分層多指數(shù)反演得到的地下含水量分布，橫軸代表弛豫時(shí)間T2，縱軸代表深度，顏色代表含水量；右圖為總含水量結(jié)果，紅色虛線為模型值，藍(lán)色實(shí)線為反演值.

表2列出了采用平滑、分塊以及分層多指數(shù)三種反演方法處理模型一數(shù)據(jù)的反演結(jié)果，對(duì)比了反演含水層的層頂深度、層底深度、含水量及弛豫時(shí)間，并對(duì)總含水量的擬合情況進(jìn)行了誤差分析.

圖7 模型一數(shù)據(jù)反演結(jié)果（a）模型一數(shù)據(jù)；（b）平滑反演結(jié)果；（c）分塊反演結(jié)果；（d）分層多指數(shù)反演結(jié)果.Fig.7 Inversion results of model one data（a）Data set of model one；（b）Result of smooth inversion；（c）Result of block inversion；（d）Result of layered multi－exponential inversion.

通過圖7及表2可以得出結(jié)論：平滑反演和分塊反演能夠得到地下含水量及平均弛豫時(shí)間的垂向分布，但由于方法的固有缺陷，無法定位深層含水量的準(zhǔn)確位置，總含水量擬合數(shù)值偏小，誤差較大.分層多指數(shù)反演能夠得到地下隨深度及弛豫時(shí)間T2的含水量分布情況.同時(shí)它精確區(qū)分了淺層、中層、深層的含水構(gòu)造，與模型設(shè)置基本一致，總含水量反演結(jié)果準(zhǔn)確，擬合誤差僅為2.59%.

3.2 模型實(shí)驗(yàn)二

構(gòu)造多弛豫時(shí)間T2含水層模型.模型數(shù)據(jù)是在如下條件下通過正演計(jì)算獲得：

表2 模型一數(shù)據(jù)反演結(jié)果Table 2 Inversion results of model one

（1）邊長(zhǎng)為100m的單匝方形線圈，采用收發(fā)共線方式；

（2）假設(shè)地下為均勻半空間，電阻率為500μm.地磁場(chǎng)強(qiáng)度B0＝46970nT，地磁傾角為60°；

（3）拉莫爾頻率為2000Hz，激發(fā)脈沖矩從0.4As到10As分為20次激發(fā)地下0m到100m；

（4）假設(shè)地下存在兩層含水層（圖8），第一層設(shè)置于地下淺層20m至40m，含有兩種弛豫時(shí)間成分，T21為100ms，對(duì)應(yīng)含水量為5%，T22為400ms，對(duì)應(yīng)含水量為5%；第二層設(shè)置于地下深層60m至80m，含有兩種弛豫時(shí)間成分，T21為200ms，對(duì)應(yīng)含水量為5%，T22為600ms，對(duì)應(yīng)含水量為5%；

（5）模型中加入5%Guass白噪聲.設(shè)置反演地層個(gè)數(shù)k＝8，初始正則化參數(shù)α＝100.

圖8 模型二Fig.8 Model two

圖9 模型二數(shù)據(jù)反演結(jié)果（a）模型二數(shù)據(jù)；（b）平滑反演結(jié)果；（c）分塊反演結(jié)果；（d）分層多指數(shù)反演結(jié)果.Fig.9 Inversion results of model two data（a）Data set of model two；（b）Result of smooth inversion；（c）Result of block inversion；（d）Result of layered multi－exponential inversion.

表3 模型二數(shù)據(jù)反演結(jié)果Table 3 Inversion results of model two

圖9和表3給出了三種方法處理模型二數(shù)據(jù)的反演結(jié)果.結(jié)論如下：平滑反演和分塊反演無法準(zhǔn)確獲得深部含水層位置信息，并且僅能獲得平均弛豫時(shí)間，無法區(qū)分同一含水層的不同弛豫時(shí)間分布，含水量擬合誤差分別達(dá)到4.79%及3.59%.而分層多指數(shù)反演能精確實(shí)現(xiàn)淺層、深層含水定位，與模型基本相符.同時(shí)，在淺層和深層區(qū)分出給出不同弛豫時(shí)間T2及其映射的含水量分布情況，含水量擬合誤差僅為2.86%.表3中給出最大局部含水量所對(duì)應(yīng)的弛豫時(shí)間T2.此外，圖9（d）中右圖所示的總含水量則代表每一含水層多個(gè)弛豫時(shí)間所映射的含水量之和，與模型設(shè)置吻合.

需要指出的是，隨著探測(cè)深度的增加，磁共振響應(yīng)的分辨率逐漸降低，T2反演精度下降，向相鄰的T2值均勻擴(kuò)散，這是導(dǎo)致深層兩個(gè)弛豫時(shí)間映射含水量偏小的本質(zhì)原因.但從總含水量的數(shù)值而言，分層多指數(shù)反演更接近于模型，顯著優(yōu)于傳統(tǒng)反演方法.

4 野外數(shù)據(jù)實(shí)例

經(jīng)模型數(shù)據(jù)驗(yàn)證分層多指數(shù)方法后，本文對(duì)野外實(shí)測(cè)MRS數(shù)據(jù)進(jìn)行反演處理，證明其適用性.

研究組利用吉林大學(xué)自主研發(fā)的JL－MRS磁共振地下水探測(cè)儀于白城市鎮(zhèn)賚縣四方坨子嫩江岸邊進(jìn)行了150m地下水探測(cè)試驗(yàn).測(cè)試點(diǎn)地磁場(chǎng)55936nT，采用邊長(zhǎng)150m方形線圈，激發(fā)頻率為2383Hz.鉆探結(jié)果表明含水層見于井深3m到21m、63m到106m和126m到148m處，含水層厚度共83m，為孔隙承壓水，靜止水位4m，當(dāng)水位下降至7m時(shí)的涌水量為320t／h.

圖10a為實(shí)測(cè)數(shù)據(jù)，平均信噪比達(dá)34.68dB（相當(dāng)于信號(hào)中僅含有1.85%的噪聲，在核磁共振反演中，信噪比高于6.02dB數(shù)據(jù)即為真實(shí)可信［23］）.圖10（b—d）是分別采用平滑、分塊和分層多指數(shù)方法對(duì)數(shù)據(jù)進(jìn)行反演的結(jié)果，圖10e為鉆井成圖.表4給出了白城測(cè)點(diǎn)鉆井地層結(jié)構(gòu)與三種反演方法的結(jié)果對(duì)比情況.平滑反演得到地下存在兩層含水層，淺層位于地下5.25m至38.25m，深層存在含水層但無法確定其具體位置，其平均弛豫時(shí)間T＊2分布平緩過渡，無法與含水層準(zhǔn)確對(duì)應(yīng).分塊反演顯示地下主要存在兩層含水層，淺層位于地下6.75m至29.25m，同樣無法準(zhǔn)確體現(xiàn)深層含水層的位置，與平滑反演相較，T＊2有明確劃分，然而深層含水層僅有較大弛豫時(shí)間247.2ms的反演解釋結(jié)果與野外實(shí)際地層情況有一定偏差.上述兩種方法反演的含水層均為連續(xù)變化，未見相對(duì)隔水層出現(xiàn).分層多指數(shù)方法的初始正則化參數(shù)α＝100，反演設(shè)置地層個(gè)數(shù)k＝10（參數(shù)優(yōu)化后，E0－q曲線擬合誤差最?。?結(jié)果顯示地下主要存在兩個(gè)含水層，淺層位于地下3.75m至33.75m，從弛豫時(shí)間上看地表土對(duì)應(yīng)的T2為53ms；隔水構(gòu)造出現(xiàn)在33.75m與41.25m之間，該隔水構(gòu)造對(duì)應(yīng)右圖所示的較小的含水量分布；中間層亞粘土地層位于地下41.25m至65.25m，對(duì)應(yīng)的弛豫時(shí)間T2為85ms比理論值略大，因?qū)嶋H探測(cè)中，亞粘土屬于束縛水構(gòu)造（一般認(rèn)為磁共振探測(cè)為自由水），反演擬合易造成偏差，而成為上下兩個(gè)含水層弛豫時(shí)間值的過渡.深層含水層層頂位置反演準(zhǔn)確，從弛豫時(shí)間上看，地下63.75m至150m主要存在兩個(gè)弛豫時(shí)間成分，分別是35ms和196ms，能夠?qū)?yīng)于鉆井結(jié)果63m到148m的巖礫石結(jié)構(gòu)與亞粘土.隨著深度的增加，激發(fā)場(chǎng)逐漸減弱，導(dǎo)致了磁共振響應(yīng)分辨率隨著深度的增加而逐漸降低，使得分層多指數(shù)方法亦不能明確的給出第三層出水層，但其弛豫時(shí)間的分布差異提示在主要含水層中或有隔水介質(zhì)出現(xiàn).

圖10 白城測(cè)點(diǎn)數(shù)據(jù)反演結(jié)果（a）白城測(cè)點(diǎn)數(shù)據(jù)；（b）平滑反演結(jié)果；（c）分塊反演結(jié)果；（d）分層多指數(shù)反演結(jié)果；（e）測(cè)點(diǎn)鉆探成井圖.Fig.10 Inversion results of Baicheng test（a）Dataset of Baicheng test；（b）Smooth inversion result；（c）Block inversion result；（d）Layered multi－exponential inversion result；（e）Drilling result.

表4 白城四方駝子鉆井結(jié)果與反演結(jié)果對(duì)比Table 4 The comparison of drilling result and inversion results at BaichengSifangtuozi

在實(shí)測(cè)數(shù)據(jù)處理中，分層多指數(shù)反演在含水層的位置上與鉆井結(jié)果一致性最好，且能夠通過弛豫時(shí)間T2分布來明確地下巖石結(jié)構(gòu).分別利用三種反演方法進(jìn)行涌水量評(píng)估，發(fā)現(xiàn)平滑反演和分塊反演計(jì)算的涌水量值偏小，而分層多指數(shù)反演得到該測(cè)點(diǎn)的涌水量約為327t／h，與鉆井結(jié)果320t／h基本相符，說明了分層多指數(shù)方式能夠準(zhǔn)確評(píng)價(jià)地下水賦存情況.

5 結(jié) 論

通過對(duì)比平滑反演、分塊反演和分層多指數(shù)反演方法處理模型數(shù)據(jù)及野外實(shí)測(cè)數(shù)據(jù)算例，本文得出如下結(jié)論：

（1）分層多指數(shù)反演采用多指數(shù)方式來擬合FID信號(hào)，摒棄了平滑與分塊反演中參數(shù)提取環(huán)節(jié)，降低了傳統(tǒng)單一指數(shù)擬合引入的誤差對(duì)反演結(jié)果的影響；將弛豫信號(hào)所有數(shù)據(jù)一次性帶入反演計(jì)算，建立起相鄰脈沖矩之間及相鄰時(shí)間下信號(hào)的相互關(guān)系，提高了反演結(jié)果的準(zhǔn)確性與穩(wěn)定性.

（2）與平滑與分塊反演得到的含水量垂向分布相比較，分層多指數(shù)反演能夠得到隨深度及弛豫時(shí)間T2變化的含水量分布情況，更加直觀的反映出地下水的賦存狀態(tài)，適用于多巖性解釋的復(fù)雜地質(zhì)環(huán)境.

（3）分層多指數(shù)反演采用設(shè)定地下分層個(gè)數(shù)的方式，以及不等式約束下的空間信賴域算法，縮短了反演所需時(shí)間，使得各含水層之間分界明顯，不僅有助于實(shí)際應(yīng)用中準(zhǔn)確判斷含水層所處的位置，還避免了預(yù)先設(shè)置復(fù)雜反演控制參數(shù)導(dǎo)致的人為結(jié)果偏差.

本文突破了傳統(tǒng)地面MRS技術(shù)僅局限于平均孔隙度解釋的瓶頸，提供了MRS在復(fù)雜地質(zhì)構(gòu)造下進(jìn)行全數(shù)據(jù)反演解釋的新思路，而隨著多通道［24］MRS儀器的問世，該反演方法可進(jìn)行深入研究，從而拓展應(yīng)用到二維全數(shù)據(jù)的反演解釋之中.

（References）

［1］Legchenko A，Valla P.A review of the basic principles for proton magnetic resonance sounding measurements.Journal of Applied Geophysics，2002，50（1－2）：3－19.

［2］Legchenko A，Baltassat J M，Beauce A，et al.Nuclear magnetic resonance as a geophysical tool for hydrogeologists.Journal of Applied Geophysics，2002，50（1－2）：21－46.

［3］潘玉玲，張昌達(dá).地面核磁共振找水儀理論與方法.北京：中國(guó)地質(zhì)大學(xué)出版社，2000.Pan Y L，Zhang C D.The Theory and Methods of Water Detecting（in Chinese）.Beijing：Publishing House of China University of Geosciences，2000.

［4］林君，段清明，王應(yīng)吉等.核磁共振找水儀原理與應(yīng)用.北京：科學(xué)出版社，2011.Lin J，Duan Q M，Wang Y J，et al.Theory and Design of Magnetic Resonance Sounding Instrument for Groundwater Defection and Its Applications（in Chinese）.Beijing：Science Press，2011.

［5］張昌達(dá)，潘玉玲.關(guān)于地面核磁共振方法資料巖石物理學(xué)解釋的一些見解.工程地球物理學(xué)報(bào)，2006，3（1）：1－8.Zhang C D，Pan Y L.Some views on petrophysical interpretation of SNMR data.Chinese Journal of Engineering Geophysics（in Chinese），2006，3（1）：1－8.

［6］Legchenko A V，Shushakov O A.Inversion of surface NMR data.Geophysics，1998，63（1）：75－84.

［7］Mohnke O，Yaramanci U.Smooth and block inversion of surface NMR amplitudes and decay times using simulated annealing.Journal of Applied Geophysics，2002，50（1－2）：163－177.

［8］翁愛華，王雪秋，劉國(guó)興等.導(dǎo)電性影響的地面核磁共振反演.地球物理學(xué)報(bào)，2007，50（3）：890－896.Weng A H，Wang X Q，Liu G X，et al.Nonlinear inversion of surface nuclear magnetic resonance over electrically conductive medium.Chinese J.Geophys.（in Chinese），2007，50（3）：890－896.

［9］戴苗，胡祥云，吳海波等.地面核磁共振找水反演.地球物理學(xué)報(bào)，2009，52（10）：2676－2682.Dai M，Hu X Y，Wu H B，et al.Inversion of surface nuclear magnetic resonance.Chinese J.Geophys.（in Chinese），2009，52（10）：2676－2682.

［10］Mohnke O，Yaramanci U.Forward modeling and inversion of MRS relaxation signals using multi－exponential decomposition.Near Surface Geophysics，2005，3（3）：165－185.

［11］Yaramanci U，Müller－Petke M.Surface unclear magnetic resonance—a unique tool for hydrogeophysics.The Leading Edge，2009，28（10）：1240－1247.

［12］孫淑琴，趙義平，林君等.用地面核磁共振方法評(píng)估含水層涌水量的實(shí)例.地球物理學(xué)進(jìn)展，2008，23（4）：1317－1321.Sun S Q，Zhao Y P，Lin J，et al.Example of evaluating groundwater－saturated layers yield capacity based on the surface NMR method.Progress in Geophys.（in Chinese），2008，23（4）：1317－1321.

［13］Ryom N M，Hagensen T F，Chalikakis K，et al.Comparison of transmissivities from MRS and pumping tests in Denmark.NearSurfaceGeophysics，2011，9（2）：211－223.

［14］Braun M，Hertrich M，Yaramanci U.Complex inversion of surface－NMR signals－extending the limits of model resolution.17th EEGS symposium on the application of Geophysics to Engineering and environmental problems.Extended Abstracts，2004：856－867.

［15］Weichmann P，Lavely E，Ritzwoller M.Theory of surface nuclear magnetic resonance with applications to geophysical imaging problems.Physical Review E，2000，62（1）：1290－1312.

［16］Müller M，Hertrich M，Yaramanci U.Analysis of magnetic resonance sounding kernels concerning large scale applications using SVD.Proc.SAGEEP，EEGS.2006，19（1）：595－603.

［17］唐隆基，李文，鄧陽生.關(guān)于解地球物理中病態(tài)方程的若干問題.地球物理學(xué)報(bào)，1995，38（1）：105－114.Tang L J，Li W，Deng Y S.Some aspects for solving ill－posed equation in geophysics.Chinese J.Geophys.（in Chinese），1995，38（1）：105－114.

［18］張文權(quán)，翁愛華.地面核磁共振正則化反演方法研究.吉林大學(xué)學(xué)報(bào) （地球科學(xué)版），2007，37（4）：809－813.Zhang W Q，Weng A H.On regularization inversion method of surface nuclear magnetic resonance data.Journal of Jilin University （Earth Science Edition）.（in Chinese），2007，37（4）：809－813.

［19］Mohnke O， Yaramanci U.Pore size distributions and hydraulic conductivities of rocks derived from Magnetic Resonance Sounding relaxation data using multi－exponential decay time inversion.Journal of Applied Geophysics，2008，66（3－4）：73－81.

［20］Schiroc M，Legchenko A，Creer G.A new direct noninvasive groundwater detection technology for Australia.Exploration Geophysics，1991，22（2）：333－338.

［21］曾齊，師學(xué)明，武永勝等.一維大地電磁信賴域反演法研究.地球物理學(xué)進(jìn)展，2011，36（3）：885－893.Zeng Q，Shi X M，Wu Y S，et al.A study of onedimensional magnetotelluric sounding inversion using the trust region algorithm.Progress in Geophys.（in Chinese），2011，26（3）：885－893.

［22］田玥，陳曉非.利用擬牛頓法和信賴域法聯(lián)合反演震中分布與一維速度結(jié)構(gòu).地球物理學(xué)報(bào)，2006，49（3）：845－854.Tian Y，Chen X F.Simultaneous inversion of hypocenters and velocity using the quasi－Newton method and trust region method.Chinese J.Geophys.（in Chinese），2006，49（3）：845－854.

［23］Legchenko A.MRS measurements and inversion in presence of EM noise.Boletín Geológico y Minero，118（3）：489－508.

［24］Walsh D.Multi－channel surface NMR instrumentation and software for 1D／2Dgroundwater investigations.Journal of Applied Geophysics，2008，66（3－4）：140－150.