郭文斌 嘉世旭 林吉焱 邱 勇

1 中國地震局地球物理勘探中心，鄭州市文化路75號(hào)，450002

重力－地震聯(lián)合反演方法可分為同步反演和順序反演兩類。常用的順序反演是一種人機(jī)交互反演［1－3］，主要依據(jù)地震速度結(jié)構(gòu)模型確定地下空間各地質(zhì)單元的形態(tài)，根據(jù)密度－速度轉(zhuǎn)換公式計(jì)算各地質(zhì)單元初始密度值，然后反演各地質(zhì)單元真實(shí)密度值。若反演無法收斂，或結(jié)果明顯不合理，則根據(jù)工作經(jīng)驗(yàn)及對(duì)區(qū)域地質(zhì)構(gòu)造的認(rèn)識(shí)調(diào)整各地質(zhì)單元形態(tài)，重新反演，直到結(jié)果符合要求為止。但若欠缺相關(guān)工作經(jīng)驗(yàn)，則會(huì)將對(duì)區(qū)域構(gòu)造的錯(cuò)誤認(rèn)識(shí)帶到反演結(jié)果中。針對(duì)這一問題，本文采用擬BP 神經(jīng)網(wǎng)絡(luò)物性反演算法，并輔以小波多尺度分析手段提升垂向分辨率，實(shí)行重力－地震的自動(dòng)聯(lián)合反演，以減少研究者的主觀影響，并將該方法用于遂寧－阿壩深地震測(cè)深剖面，驗(yàn)證其效果。

1 重力地震聯(lián)合反演

反演目標(biāo)遂寧－阿壩剖面［4］位于四川盆地中部，長約500km，其位置如圖1中虛線所示。該剖面跨越多個(gè)地質(zhì)構(gòu)造單元，其地殼速度、界面非均勻構(gòu)造模型及解釋研究可參見文獻(xiàn)［4］。

重力－地震聯(lián)合反演的第一步是建立速度－密度的轉(zhuǎn)換關(guān)系，統(tǒng)一物性參數(shù)，將地震獲得的速度結(jié)構(gòu)模型轉(zhuǎn)換為密度結(jié)構(gòu)模型。參照物性反演原理，將反演區(qū)域以長方形單元均勻剖分，再將速度結(jié)構(gòu)轉(zhuǎn)換為密度結(jié)構(gòu)，映射到反演區(qū)域。為準(zhǔn)確描述遂寧－阿壩速度模型的結(jié)構(gòu)特征，用大小為10km×1km 的長方形單元格將反演區(qū)域剖分為50×60個(gè)單元格（圖2）。

目前反演區(qū)域沒有很好的波速－密度關(guān)系經(jīng)驗(yàn)公式［5］，本文使用的反演算法受模型初值影響較?。?－7］，可以選用馮銳等提出的波速－密度關(guān)系計(jì)算初始密度模型［1－2］：

轉(zhuǎn)換結(jié)果如圖3所示，其基本特征與地震速度結(jié)構(gòu)特征一致。

聯(lián)合反演中，地震反演結(jié)果的另一個(gè)重要作用就是提供約束，盡可能降低重力反演的多解性。本文通過控制各單元格密度的變化范圍來約束反演方向，即根據(jù)速度結(jié)構(gòu)模型確定各圈層包含單元格的變化范圍，每迭代m次后，對(duì)各單元密度值作如下調(diào)整：

圖1 遂寧－阿壩剖面（虛線）位置圖Fig.1 The position of Suining－Aba profile（the dotted line）

圖2 遂寧－阿壩剖面速度結(jié)構(gòu)模型及網(wǎng)格剖分Fig.2 The velocity structure and mesh generation of Suining－Aba profile

圖3 遂寧－阿壩剖面密度結(jié)構(gòu)模型Fig.3 The density structure of Suining－Aba profile

2 反演算法

綜上所述，該次重力－地震聯(lián)合反演計(jì)算量巨大，并且模型初值并不十分準(zhǔn)確。本文采用的擬BP神經(jīng)網(wǎng)絡(luò)非線性反演算法，使用一個(gè)固定值代替梯度值作為每次迭代的調(diào)整量，能夠有效節(jié)省存儲(chǔ)空間，加快收斂速度，并且具有跳出局部最小的能力［6－8］，其反演過程詳見文獻(xiàn)［6－8］。針對(duì)文獻(xiàn)［8］提到的深部垂向分辨率略有不足的問題，本文引入小波多尺度變換進(jìn)行改善。

利用小波多尺度分析及Mallat算法，重力異常可分解為：

其中，AJΔg（x）為重力異常的j階接近，即重力異常的低頻部分，通常是由深部構(gòu)造或區(qū)域構(gòu)造引起的，DiΔg（x）是j階分解后的各階細(xì)節(jié)異常，對(duì)應(yīng)重力異常的高頻部分，多為淺部或局部異常，DiΔg（x）不隨j的增大而變化。多尺度反演中多選用近似部分為反演對(duì)象，首先在j尺度下反演得到包含主要深部信息和區(qū)域信息的模型，將其作為j－1尺度下反演的初始模型。隨著尺度的減小，反演對(duì)象的淺部信息逐漸增加，模型的局部及淺部特征逐步顯現(xiàn)，當(dāng)j＝0時(shí)，AJΔg（x）＝Δg（x），此時(shí)反演結(jié)果即可包含所有信息。因此，小波多尺度變換可以在反演過程中實(shí)現(xiàn)自動(dòng)的異常分離，從而提高重力的垂向分辨率［9－11］。

3 反演結(jié)果分析

為防止目標(biāo)函數(shù)取值過小導(dǎo)致過度擬合或者無法收斂，反演時(shí)設(shè)目標(biāo)函數(shù)Φ＜0.01或者目標(biāo)函數(shù)值無明顯變化時(shí)為迭代停止條件。目標(biāo)函數(shù)初始值為0.08，經(jīng)過1 500次迭代降至0.01，目標(biāo)函數(shù)收斂，反演過程耗時(shí)343s（酷睿i3 處理器，2G 內(nèi)存）。重力異常曲線的擬合效果如圖4所示。可以看出，反演結(jié)果與實(shí)測(cè)異常曲線有較好的擬合效果。

圖4 重力異常擬合曲線Fig.4 The fitting curve of gravity

聯(lián)合反演結(jié)果如圖5所示，其密度結(jié)構(gòu)特征與剖面的地質(zhì)解釋［4］（黑實(shí)線）有很好的一致性，基本保留了地震數(shù)據(jù)的界面分層信息。與初始模型相比，反演結(jié)果有較明顯的變化：四川盆地各圈層內(nèi)部密度較為均一，而川西北高原中下地殼密度相對(duì)混雜，相鄰各圈層互有滲透；龍門山褶皺帶下方15～30km 處的低密度填充變?yōu)橄鄬?duì)高密度填充，并且其密度與下地殼相近；西北高原地下介質(zhì)整體密度降低，四川盆地中下地殼密度變大，上地殼與川西北高原上地殼密度接近，使得龍門山褶皺帶與四川盆地交界處的斷裂特征更加明顯。總體而言，剖面地下空間的界面分層信息與地震速度結(jié)構(gòu)模型基本吻合，但橫向變化特征更為明顯。

圖5 聯(lián)合反演結(jié)果Fig.5 The result of joint inversion

結(jié)合剖面所處的地質(zhì)環(huán)境特征及相關(guān)資料［3］可知，本文反演結(jié)果有更貼近實(shí)際的地質(zhì)學(xué)解釋：四川盆地中下地殼密度較大，地殼分層明顯，說明該區(qū)域巖石完整，整體性較好，川北高原中下地殼向東運(yùn)動(dòng)受其阻礙，使得各圈層內(nèi)巖石擠壓，沿龍門山斷裂帶上涌，導(dǎo)致川西北高原地下介質(zhì)整體密度降低，相鄰各圈層界線模糊（擠壓作用使巖石破碎），龍門山褶皺帶中上地殼有高密度物質(zhì)填充（下地殼物質(zhì)上涌）。

圖5中存在個(gè)別較為明顯、但顯然不可信的局部細(xì)小構(gòu)造特征，這主要是由于反演過程中完全拋棄了人工干涉，并且未對(duì)結(jié)果在分辨率尺度上進(jìn)行合適的平滑濾波，導(dǎo)致未能將低于分辨能力的假異常完全抹除。

4 結(jié) 語

本文反演所得密度結(jié)構(gòu)圖（圖5）較好地保留了速度結(jié)構(gòu)圖（圖2）劃分的圈層構(gòu)造，且與該剖面的地質(zhì)學(xué)解釋具有很好的一致性。兩者的差異主要體現(xiàn)在少見的“低速－高密度”對(duì)應(yīng)關(guān)系：速度結(jié)構(gòu)圖中下地殼的低速通道在密度結(jié)構(gòu)圖中消失不見，龍門山褶皺帶中上地殼的低速異常在密度結(jié)構(gòu)圖中表現(xiàn)為高密度異常。這是由于地震波速的變化對(duì)地質(zhì)界面、巖石破碎程度相當(dāng)敏感，而重力對(duì)巖性的變化更為敏感。下地殼低密度通道的消失說明了該區(qū)域并無明顯的巖性變化，速度結(jié)構(gòu)圖中低速通道主要由巖石破碎而非巖性變化引起；中上地殼相應(yīng)位置的高密度填充，說明該區(qū)域包含下地殼物質(zhì)，由于巖石破碎，其在速度結(jié)構(gòu)圖中體現(xiàn)為低速。重力－地震聯(lián)合反演結(jié)果的這種“低速－高密度”對(duì)應(yīng)關(guān)系更完美地詮釋了該地區(qū)“褶皺造山帶下地殼塑性流變介質(zhì)在四川盆地下地殼剛性結(jié)構(gòu)強(qiáng)力阻擋下被迫轉(zhuǎn)向向上，在松潘－甘孜褶皺帶東側(cè)邊緣形成由下地殼發(fā)起、有緩而陡向上逆沖貫穿地殼的巨型上升流”［4］的地質(zhì)構(gòu)造特征。

綜上所述，本文所使用的重力－地震聯(lián)合反演方法合理可信，能夠結(jié)合不同方法的優(yōu)點(diǎn)，使用地殼內(nèi)部介質(zhì)多個(gè)參數(shù)特征對(duì)剖面作出更合理、可靠的地質(zhì)解釋。與傳統(tǒng)的人機(jī)交互聯(lián)合反演相比，該方法簡單、快捷，反演過程無需人工參與，受研究者主觀影響較小。遂寧－阿壩剖面的反演結(jié)果也顯示了該方法的不足之處，即無法完全抹除個(gè)別不可分辨的假異常。針對(duì)這一問題，可以考慮借鑒最小構(gòu)造反演［12］相關(guān)方法進(jìn)行改進(jìn)。

［1］劉皓，方盛明，嘉世旭.地震力聯(lián)合反演及效果——以天津－北京－赤城地震探測(cè)剖面為例［J］.地震學(xué)報(bào)，2011，33（4）：443－450（Liu Hao，F(xiàn)ang Shengming，Jia Shixu.Jiont Inversion of Seismic and Gravity Data and Its Effect：Taking the Tianjin－Beijing－Chicheng Seimic Profile as Example［J］.Acta Seismologica Sinica，2011，33（4）：443－450）

［2］馮銳，嚴(yán)惠芬，張若水.三維位場(chǎng)的快速反演方法及程序設(shè)計(jì)［J］.地質(zhì)學(xué)報(bào)，1986，60（4）：390－403（Feng Rui，Yan Huifen，Zhang Ruoshui.The Rapid Inversion of 3－D Potential Field and Program Design［J］.Acta Geologica Sinica，1986，60（4）：390－403）

［3］袁惟正，徐新忠，雷江鎖，等.大別山地震波速剖面的重力擬合及花崗巖帶［J］.中國地質(zhì)，2003，30（3）：235－239（Yuan Weizheng，Xu Xinzhong，Lei Jiangsuo.Simulation of Gravity Anomaly of the Velocity Profile from Wild－Angle Reflection Seimic Profiling and the Granite Belt in the Dabie Mountatins［J］.Geology in China，2003，30（3）：235－239）

［4］嘉世旭，劉保金，徐朝繁，等.龍門山中段及兩側(cè)地殼結(jié)構(gòu)與汶川地震構(gòu)造［J］.中國科學(xué)：地球科學(xué)，2014，44（3）：497－509（Jia Shixu，Liu Baojin，Xu Zhaofan.The Crustal Structures of the Central Longmenshan along and Its Margins as Related to the Seismotectonics of the 2008Wenchuan Earthquake［J］.Science China Earth Sciences，2014，44（3）：497－509）

［5］Barton P J.The Relationship between Seismic Velocity and Density in the Continental Crust：A Useful Constrain［J］.Geophys J R Astr Soc，1986，87（1）：195－208

［6］郭文斌，朱自強(qiáng)，魯光銀.重力梯度張量的擬BP神經(jīng)網(wǎng)絡(luò)反演［J］.中南大學(xué)學(xué)報(bào)：自然科學(xué)版，2011，42（12）：3 798－3 803（Guo Wenbin，Zhu Ziqiang，Lu Guangyin.Quasi－BP Neural Network Inversion of Gravity Gradient Tensor［J］.Journal of Central South University：Science and Technology 2011，42（12）：3 798－3 803）

［7］郭文斌，朱自強(qiáng)，魯光銀.重力異常的BP神經(jīng)網(wǎng)絡(luò)三維物性反演［J］.地球物理學(xué)進(jìn)展，2012，27（2）：409－416（Guo Wenbin，Zhu Ziqiang，Lu Guangyin.3－D Gravity Inversion for Physical Properties Using BP Network［J］.Progress of Geophysics，2012，27（2）：409－416）

［8］管志寧，侯俊生，黃臨平，等.重磁異常反演的擬BP神經(jīng)網(wǎng)絡(luò)方法及其應(yīng)用［J］.地球物理學(xué)報(bào)，1998，41（2）：242－251（Guan Zhining，Hou Junsheng，Huang Linping，et al.Inversion of Gravity and Magnetic Anomalies Using Pseduo－BP Neural Network Method and Its Application［J］.Chinese Journal of Geophysics，1998，41（2）：242－251）

［9］邱寧，何展翔，昌言君.分析研究基于小波分析與譜分析提高重力異常的分辨能力［J］.地球物理學(xué)進(jìn)展，2007，22（1）：112－120（Qiu Ning，He Zhanxiang，Chang Yanjun.Ability of Improving Gravity Anomaly Resolution Based on Multiresolution Wavelet Analysis and Power Spectrum Analysis［J］.Progress of Geophysics，2007，22（1）：112－120）

［10］侯遵澤，楊文采，劉家琦.中國大陸地殼密度差異多尺度反演［J］.地球物理學(xué)報(bào)，1998，41（5）：642－651（Hou Zunze，Yang Wencai，Liu Jiaqi.Multiscale Inversion of the Density Contrast within the Crust of China［J］.Chinese Journal of Geophysics，1998，41（5）：642－651）

［11］吳立新，楊明芝，趙衛(wèi)明，等.利用重力多尺度分解資料反演青藏高原東北緣地殼厚度［J］.大地測(cè)量與地球動(dòng)力學(xué)，2011，31（1）：19－23（Wu Lixin，Yang Mingzhi，Zhao Weiming.Crust Thickness Inversed from Multi－Scale Decomposition of Bouguer Gravity Anomalies in Northeastern of Qinghai－Tibet Plateau［J］.Geodesy and Geodynamic，2011，31（1）：19－23）

［12］Zelt C A，Barton P J.Three－Dimensional Seismic Refraction Tomography：A Comparison of Two Methods Applied to Data from the Faeroe Basin［J］.Journal of Geophysical Research Atmospheres，1998，103（B4）：7 187－7 210