周俊杰，陳聰，喻翔，高玲舉，陳濤

(核工業(yè)北京地質(zhì)研究院，北京 100029)

0 引言

地球物理三維反演是獲取深部地質(zhì)結(jié)構(gòu)信息的重要技術(shù)之一，是近年來礦產(chǎn)勘查工作攻深找盲的熱點研究問題之一[1-2]。重力異常對橫向密度變化敏感度高，在礦產(chǎn)勘查中常用于推斷地層起伏以及局部巖體分布等，但垂向分辨能力較差，通常需要引入先驗信息來約束，以期改善重力反演結(jié)果的質(zhì)量。因此，約束的構(gòu)建和引入是重力三維反演中的關(guān)鍵問題。Li 和Oldenburg[3]、Williams[4]、姚長利等[5]較早的介紹了重力三維約束反演的方法原理，通過理論模型和典型實測案例說明了重力三維反演中約束的重要性。隨后，關(guān)于重力三維約束反演的相關(guān)研究不斷涌現(xiàn)。祁光等[6]、蘭學(xué)毅等[7]使用交互方式研究了位場反演中的先驗信息約束問題，并在礦集區(qū)得到成功應(yīng)用。郭冬等[8]對三維密度填圖進行了研究，認為地表地質(zhì)信息和鉆孔信息對提升重力反演的垂向分辨能力具有幫助。羅凡等[9]在深部磁鐵礦勘查中利用先驗信息約束的重磁反演圈定了礦體范圍。陳長敬等[10]采用電法資料構(gòu)建約束進行重力三維反演，利用電法垂向分辨率高的特點改善了重力反演結(jié)果。

重力反演的約束方式有多種，如參考模型約束[11]、物性范圍約束[12]、方向約束[13]、聚焦約束[14]以及結(jié)構(gòu)約束[15]等，可適用不同種類、不同性質(zhì)的先驗信息的引入。相較于其他方式，參考模型具有理論簡單，穩(wěn)定可靠等特點，但其構(gòu)建流程較為繁瑣。使用參考模型進行約束的主要步驟為：廣泛收集已知資料并轉(zhuǎn)化為地下空間物性信息，構(gòu)建參考模型并衡量相應(yīng)權(quán)重，最后代入反演計算得到結(jié)果。在實際應(yīng)用中，常見的可利用先驗信息包括地表地質(zhì)圖、物性統(tǒng)計數(shù)據(jù)、鉆孔資料以及其他地球物理成果等。筆者針對參考模型的構(gòu)建及約束的引入問題進行分析，通過理論算例和實測資料分析探討重力三維約束反演的效果，結(jié)果顯示參考模型法可有效提升反演質(zhì)量，值得在實踐中推廣應(yīng)用。

1 重力三維約束反演方法

重力三維約束反演利用重力數(shù)據(jù)及相關(guān)地質(zhì)信息約束，計算地下密度三維分布特征。Li 和Oldenburg[16]提出的重力三維約束反演方法在重力勘探領(lǐng)域得到廣泛應(yīng)用，其基本思想是將地下三維空間剖分為包含若干直立長方體單元的規(guī)則網(wǎng)格，通過求解每個單元物性值來獲取地下密度分布情況。本文采用以此理論為基礎(chǔ)的UBC-GIF GRAV3D程序進行反演計算。為便于先驗信息的引入，反演算法在滿足數(shù)據(jù)擬合的同時也引入了模型擬合，其目標函數(shù)可簡寫為：

φ(m)=‖Wd[d0-g(m)]‖2+

β[α1‖Wm(m-mref)‖2+α2‖D(m-mref)‖2]，

(1)

其中：第一項為重力數(shù)據(jù)擬合目標函數(shù)；第二項為模型目標函數(shù)，內(nèi)含最小模型項和光滑模型項，β為正則化參數(shù)。兩項目標函數(shù)均用L2范數(shù)定義。Wd為觀測數(shù)據(jù)標準差對角陣，d0為重力數(shù)據(jù)向量，m為待求模型向量，g為正演算子。模型項中，mref為參考模型向量，D為三維光滑度矩陣，α1,α2分別為最小項和光滑項的權(quán)重因子。在數(shù)據(jù)項和模型項目標函數(shù)的共同作用下，反演迭代所得結(jié)果將在滿足數(shù)據(jù)擬合的基礎(chǔ)上，同時符合參考模型的特征，并具有光滑特性。參考模型作為已知向量，已經(jīng)包含在目標函數(shù)中，需要在反演前預(yù)先構(gòu)建。

2 參考模型構(gòu)建方法

參考模型是指和待反演模型具有相同網(wǎng)格架構(gòu)，且內(nèi)含地下空間已知物性信息的三維數(shù)據(jù)體。根據(jù)式(1)，反演迭代解將盡可能的貼近參考模型特征，以此來達到約束的效果。為表明參考模型約束的作用，使用典型理論模型進行分析測試。圖1所示為測試采用的理論模型，其東向、北向和垂向網(wǎng)格規(guī)模為100×1×30，其中剖面方向單元長度為10 m，走向方向單元長度1 000 m，為降低邊界效應(yīng)影響，在模型體邊緣部位均進行了擴邊處理。模型由淺地表低密度體和起伏基底構(gòu)成，各地質(zhì)單元密度參數(shù)如表1所示。取平均密度為2.6×103kg/m3，對剩余密度進行三維反演。首次反演采用無約束方式，即不考慮先驗信息，參考模型為0的情形，如圖2a所示；假設(shè)已知基底的大致深度信息，將參考模型設(shè)置為層狀模型，代入反演所得結(jié)果如圖2b所示；假設(shè)除基底信息外，還已知淺地表密度信息，將其納入?yún)⒖寄Ｐ?，代入反演所得結(jié)果如圖2c所示。3種情況均使用GRAV3D程序的默認反演參數(shù)，僅對比使用不同參考模型時反演效果。

圖2 不同參考模型下反演結(jié)果對比Fig.2 Comparison of inverted results under different reference model constraints

表1 理論密度模型參數(shù)設(shè)置

圖1 理論密度模型Fig.1 Synthetic density model

從理論模型反演試驗可以看出，參考模型約束對于改善反演結(jié)果質(zhì)量具有成效。當不使用參考模型(即參考模型為0)時，反演僅在正則化條件下對觀測數(shù)據(jù)進行擬合，受制于重力反演體積效應(yīng)的影響，其結(jié)果僅在橫向上較好地分辨出基底隆起位置，但其層位難以進行解釋。當已知測區(qū)內(nèi)基底平均深度時，可根據(jù)蓋層和基底密度特點構(gòu)建一個層狀模型，反演結(jié)果則會在滿足數(shù)據(jù)擬合的同時盡可能符合參考模型特征，相較于無約束反演，在分辨出基底隆起橫向位置的同時還能反演出基底和蓋層的分界面，但基底隆起處密度界面模糊，不利于隆起的進一步解釋。當已知淺地表信息時，可在層狀參考模型基礎(chǔ)上再引入局部低密度體。值得注意的是，兩者的可靠程度不同，通常淺地表信息更加準確，可賦予高權(quán)重值，反演結(jié)果會更加貼近高權(quán)重位置的參考模型密度特征。在本例中，引入地表信息后，模型權(quán)重取地表30 m內(nèi)從100衰減至10，其余部位權(quán)重均為1。反演結(jié)果可看出，地表密度分布更符合設(shè)定信息，而層狀基底的隆起部位更加明顯，埋深也與已知數(shù)值相近?？梢姡?將不同的地質(zhì)信息以參考模型方式引入到重力反演中，可有效提高反演質(zhì)量，得到更加可靠的反演結(jié)果。

在實際應(yīng)用中，重力三維反演的關(guān)鍵在于參考模型的構(gòu)建。一般可先在不引入先驗信息的情況下開展無約束反演，獲得地下密度分布初步認識；之后結(jié)合測區(qū)地質(zhì)圖、物性資料和相關(guān)勘探成果，根據(jù)巖體和地層的大致分布特點和密度統(tǒng)計規(guī)律構(gòu)建參考模型。對于較可靠的地表信息可賦予高權(quán)重，而無法確認是否準確的深部推測資料可賦予低權(quán)重。此外，鉆孔信息也可賦予較高權(quán)重加以利用。實驗表明，通常情況下高權(quán)重設(shè)置為低權(quán)重的100倍左右較為適宜，反演結(jié)果既能具有參考模型特征，又能接受迭代擾動，滿足正常的觀測數(shù)據(jù)擬合要求。

3 張麻井礦區(qū)應(yīng)用實例

3.1 礦區(qū)概況

筆者結(jié)合張麻井礦區(qū)實例來分析參考模型在重力三維反演中的構(gòu)建及應(yīng)用。張麻井礦區(qū)(圖3a)位于沽源火山斷陷盆地西南緣，區(qū)域面積約64 km2，盆地基底為太古宇紅旗營子群中深變質(zhì)巖系，研究區(qū)中未出露；蓋層為侏羅系上統(tǒng)及白堊系下統(tǒng)火山熔巖、火山碎屑巖系。其中，區(qū)內(nèi)西北部大面積出露張家口組二段(J3z2)粗面巖，中部及南部多為張家口組三段(J3z3)流紋巖、角粒熔巖，西南部出露張家口組三段(J3z3)凝灰?guī)r，東北部出露花吉營組(K1h1)安山巖等。區(qū)內(nèi)地表第四系覆蓋廣闊，為下賦地層界線的劃定帶來一定困難。礦床主要受構(gòu)造、巖漿活動、熱液蝕變等因素綜合控制，摸清成礦地質(zhì)結(jié)構(gòu)對深部找礦具有較大價值[17-18]。

1—第四系；2—花吉營組安山巖；3—張家口組二段粗面巖；4—張家口組三段第8層流紋巖；5—張家口組三段第7層角粒熔巖；6—張家口組三段第6層凝灰?guī)r；7—張家口組三段第5層板狀流紋巖；8—張家口組三段第1～4層凝灰?guī)r；9—次流紋斑巖；10—花崗斑巖；11—硅化帶；12—地質(zhì)界線；13—實測斷層；14—推測斷層；15—剖面；16—鉆孔

在測區(qū)內(nèi)開展了1∶25 000大比例尺重力測量，對布格重力異常進行分離后得到剩余布格重力異常，其形態(tài)可基本反映2 km深度范圍內(nèi)地層走勢和巖體分布[19]。研究區(qū)西北部呈現(xiàn)高值異常，表征地下高密度地層較厚，而中部低值異常則表征低密度地層較厚。在研究區(qū)開展了巖石密度測量統(tǒng)計，并結(jié)合前人研究成果，獲得本區(qū)巖石物性主要特征(表2)。蓋層侏羅系上統(tǒng)及白堊系下統(tǒng)巖石密度較低，其中張家口組流紋巖(J3z3)密度為2.34×103kg/m3，凝灰?guī)r(J3z3)密度為2.41×103kg/m3，花吉營組安山巖(K1h)密度為2.41×103kg/m3。張家口組角粒熔巖(J3z3-7)和粗面巖(J3z2)密度較高，分別為2.51×103kg/m3和2.52×103kg/m3。由于區(qū)內(nèi)基底資料較為缺乏，據(jù)以往資料太古宇混合變質(zhì)基底(Arh)密度為2.69×103kg/m3，遠高于火山巖密度。此外，次流紋斑巖(Cλπ)也顯示出密度高的特征，為2.59×103kg/m3。

3.2 無參考模型約束反演

根據(jù)測區(qū)基本情況，將地下介質(zhì)剖分為138×152×35的規(guī)則網(wǎng)格，單元長度為50 m。為削弱邊界效應(yīng)，各方位向外各擴展10個網(wǎng)格。使用2.57×103kg/m3作為背景密度值，結(jié)合物性統(tǒng)計數(shù)據(jù)，將剩余密度限定在(-0.5 ～0.5)×103kg/m3之間。反演使用GRAV3D程序，除物性范圍和參考模型外均使用默認參數(shù)。圖4為無參考模型約束的反演結(jié)果。從密度切片可見，地下密度分布在西北較大，東南較小，與剩余布格重力異常特征相符，但密度分層不明顯，深部密度較低，與地質(zhì)認識不符。此外，淺地表密度形態(tài)較為雜亂，和巖體露頭、出露地層特征差異很大。這說明無約束反演結(jié)果分辨能力較弱，可靠性低，不利于開展進一步的地質(zhì)解釋工作。

圖4 無參考模型約束反演結(jié)果外觀(a)和切片(b)Fig.4 Perspective(a) and sliced view(b) of inverted result without reference model constraint

3.3 參考模型構(gòu)建

3.3.1 物性填充法構(gòu)建地表約束

根據(jù)巖石密度特征和地質(zhì)圖相關(guān)層位幾何信息可構(gòu)建地表物性約束。利用自主研發(fā)的地表信息數(shù)字化軟件將巖層密度值填充到對應(yīng)位置，其網(wǎng)格大小和前述三維網(wǎng)格的平面設(shè)置相同。首先根據(jù)地質(zhì)圖將第一層地表物性進行數(shù)字化填充，之后順層推斷。推斷完畢后，將平面網(wǎng)格按照深度關(guān)系合為一體，組裝為三維網(wǎng)格模型。圖5所示為根據(jù)測區(qū)地質(zhì)圖對地表巖性進行數(shù)字化編碼，進而整合為三維地表參考模型結(jié)果，結(jié)合表2中統(tǒng)計的各類巖性密度值進行物性填充，整合后即可用作參考模型。

圖5 基于地質(zhì)圖的物性數(shù)字化和三維模型構(gòu)建Fig.5 Geological-map-based density digitization and 3D model construction

表2 研究區(qū)巖石密度特征

3.3.2 分層賦值法構(gòu)建地層約束

在測區(qū)內(nèi)收集了密度界面反演結(jié)果和二維人機交互正反演結(jié)果[19]，可將其作為深部地層信息引入到反演中。首先將密度界面按照三維網(wǎng)格參數(shù)進行重采樣，隨后轉(zhuǎn)化為網(wǎng)格模型，在對應(yīng)地層單元上賦予相應(yīng)的物性代號值(圖6a、b)，之后結(jié)合表2中物性統(tǒng)計值進行填充。對于多層地層，需注意地層整合及穿插關(guān)系，斷層處的地層錯斷也需注意，以免發(fā)生物性填充錯誤。由于深部地層信息通常為推斷結(jié)果，在反演時應(yīng)賦予不高于地表信息的權(quán)重。如有鉆孔信息，則應(yīng)在鉆孔處賦予較高的權(quán)重。

3.3.3 模型融合

地表信息和地層信息構(gòu)建為三維模型后，將兩者進行拓撲融合(圖6c、d)，主要思路為權(quán)重大的物性單元將權(quán)重小的代替。此外還需注意地表和露頭的連續(xù)性，可采用剖面掃描方式逐個進行檢查。之后還需將相應(yīng)的權(quán)重模型也進行融合，確保準確度高的物性單元配以高權(quán)重。

圖6 基于地層信息的模型構(gòu)建及其與地表模型的融合Fig.6 Strata-based model construction and fusion with surface model

3.4 參考模型約束反演結(jié)果

圖7為對張麻井鈾礦研究區(qū)的重力三維約束反演結(jié)果， GRAV3D程序參數(shù)設(shè)置除上述參考模型外，其他均與圖4相同。反演初始均方差為1.3×107，經(jīng)12次迭代反演，最終均方差為2.8×103，收斂曲線呈對數(shù)下降并趨于穩(wěn)定，數(shù)據(jù)擬合非常理想。此外，該模型也具有參考模型提供的地表地質(zhì)和地層起伏基本特征，同時密度的空間分布信息更為豐富，對中淺層密度不均勻體也有很好的體現(xiàn)。圖8a、b分別為使用和不使用參考模型得到的反演結(jié)果對比。盡管無約束結(jié)果也能描述基底的隆起，但在深部無法反映出火山蓋層，淺部物性分布雜亂；和無約束情況相比較，約束反演結(jié)果則有了很大改善：反演結(jié)果反映的地層露頭形狀和地質(zhì)圖相匹配，密度幅值符合統(tǒng)計規(guī)律，呈現(xiàn)出隨深度增加，密度逐漸增大的特點，且地層具有明顯的分層特征，呈現(xiàn)西北隆起，東南凹陷的特點，與重力異常特征相符。此外，中淺層隱伏巖體在反演剖面上也反映較為明顯，呈現(xiàn)出次高密度異常。

綜合分析反演結(jié)果可知，地表第四系密度較低，在研究區(qū)廣泛分布，但其埋深有限，多在50 m以內(nèi)。北側(cè)出露張家口組二段粗面巖(J3z2)，呈現(xiàn)出密度中等偏高的特點，且密度隨深度加大而逐漸變大。南側(cè)出露多為流紋巖(J3z3)，僅西南側(cè)為凝灰?guī)r(J3z3)，其密度比粗面巖較低，角粒熔巖(J3z3-7)在西大山以西部分出露，和圍巖相比密度值略高。東北部出露安山巖(K1h)，其密度相對較低。從圖7b剖面可知，基底和火山蓋層均呈現(xiàn)北側(cè)隆起、中南部凹陷，南側(cè)隆起的特征，其中南側(cè)隆起比參考值更高，但低于北側(cè)隆起。三維密度模型對巖體分布也有反映，在張麻井礦床處顯示局部高密度特征，表征了次流紋斑巖的賦存；西大山西側(cè)角粒熔巖也顯示出局部次高密度特征。西大山隆起處也出現(xiàn)次高密度特征。根據(jù)SN向密度剖面AA′結(jié)果(圖8)，測區(qū)北部為地層隆起，地表有張家口組二段粗面巖出露；而中部和南部主要為張家口組三段流紋巖。結(jié)合巖層密度特征和地質(zhì)特征，可推斷F45下盤張家口組二段粗面巖頂界面埋深約為海拔600～1 000 m之間；太古界混合變質(zhì)基底頂界面埋深在海拔100 m以下。在F45上盤，變質(zhì)基底頂界面埋深在海拔100～1 000 m之間。據(jù)鉆孔資料，張麻井ZKB-1在海拔810 m處見張家口組二段粗面巖，密度剖面解釋結(jié)果與鉆孔揭露位置相符。綜上所述，參考模型約束反演由于先驗信息的引入而變得更加可靠。相比無約束反演的分辨能力更好、質(zhì)量更高；相較于收集的界面反演資料，反演結(jié)果信息更為豐富，在體現(xiàn)地層錯斷和刻畫中淺層巖體賦存方面優(yōu)勢更為明顯。

圖7 參考模型約束反演結(jié)果外觀(a)和切片(b)Fig.7 perspective and sliced view of inverted result with reference model constraint

圖8 AA′剖面反演結(jié)果對比Fig.8 Comparison of inversion results on AA′ profile

4 結(jié)論

1) 將地質(zhì)圖、物性統(tǒng)計資料以及相關(guān)推斷成果作為參考模型約束引入到重力三維反演，能顯著提升反演結(jié)果質(zhì)量。借助物性填充、分層構(gòu)建等方法可有效構(gòu)建參考模型，并根據(jù)準確度賦予相應(yīng)的權(quán)重，是一種實用的先驗信息引入方法，值得在實踐中推廣應(yīng)用。

2) 筆者將地表約束信息、地層約束信息整合為參考模型，實現(xiàn)對反演結(jié)果的改善，并用鉆孔信息驗證結(jié)果的正確性。在實際應(yīng)用中，鉆孔信息也可作為高準確度的約束信息引入到反演中。如何快捷的構(gòu)建鉆孔參考模型并使之實用化，值得今后進一步研究。

致謝：衷心感謝王萬銀教授、兩位審稿專家和責任編輯為本文提出的寶貴意見!