郭文斌 嘉世旭 林吉焱 邱 勇
1 中國地震局地球物理勘探中心,鄭州市文化路75號,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)及對區(qū)域地質(zhì)構(gòu)造的認(rèn)識調(diào)整各地質(zhì)單元形態(tài),重新反演,直到結(jié)果符合要求為止。但若欠缺相關(guān)工作經(jīng)驗(yàn),則會將對區(qū)域構(gòu)造的錯誤認(rèn)識帶到反演結(jié)果中。針對這一問題,本文采用擬BP 神經(jīng)網(wǎng)絡(luò)物性反演算法,并輔以小波多尺度分析手段提升垂向分辨率,實(shí)行重力-地震的自動聯(lián)合反演,以減少研究者的主觀影響,并將該方法用于遂寧-阿壩深地震測深剖面,驗(yàn)證其效果。
反演目標(biāo)遂寧-阿壩剖面[4]位于四川盆地中部,長約500km,其位置如圖1中虛線所示。該剖面跨越多個地質(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個單元格(圖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é)果的另一個重要作用就是提供約束,盡可能降低重力反演的多解性。本文通過控制各單元格密度的變化范圍來約束反演方向,即根據(jù)速度結(jié)構(gòu)模型確定各圈層包含單元格的變化范圍,每迭代m次后,對各單元密度值作如下調(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
綜上所述,該次重力-地震聯(lián)合反演計(jì)算量巨大,并且模型初值并不十分準(zhǔn)確。本文采用的擬BP神經(jīng)網(wǎng)絡(luò)非線性反演算法,使用一個固定值代替梯度值作為每次迭代的調(diào)整量,能夠有效節(jié)省存儲空間,加快收斂速度,并且具有跳出局部最小的能力[6-8],其反演過程詳見文獻(xiàn)[6-8]。針對文獻(xiàn)[8]提到的深部垂向分辨率略有不足的問題,本文引入小波多尺度變換進(jìn)行改善。
利用小波多尺度分析及Mallat算法,重力異??煞纸鉃椋?/p>
其中,AJΔg(x)為重力異常的j階接近,即重力異常的低頻部分,通常是由深部構(gòu)造或區(qū)域構(gòu)造引起的,DiΔg(x)是j階分解后的各階細(xì)節(jié)異常,對應(yīng)重力異常的高頻部分,多為淺部或局部異常,DiΔg(x)不隨j的增大而變化。多尺度反演中多選用近似部分為反演對象,首先在j尺度下反演得到包含主要深部信息和區(qū)域信息的模型,將其作為j-1尺度下反演的初始模型。隨著尺度的減小,反演對象的淺部信息逐漸增加,模型的局部及淺部特征逐步顯現(xiàn),當(dāng)j=0時,AJΔg(x)=Δg(x),此時反演結(jié)果即可包含所有信息。因此,小波多尺度變換可以在反演過程中實(shí)現(xiàn)自動的異常分離,從而提高重力的垂向分辨率[9-11]。
為防止目標(biāo)函數(shù)取值過小導(dǎo)致過度擬合或者無法收斂,反演時設(shè)目標(biāo)函數(shù)Φ<0.01或者目標(biāo)函數(shù)值無明顯變化時為迭代停止條件。目標(biāo)函數(shù)初始值為0.08,經(jīng)過1 500次迭代降至0.01,目標(biāo)函數(shù)收斂,反演過程耗時343s(酷睿i3 處理器,2G 內(nèi)存)。重力異常曲線的擬合效果如圖4所示。可以看出,反演結(jié)果與實(shí)測異常曲線有較好的擬合效果。
圖4 重力異常擬合曲線Fig.4 The fitting curve of gravity
聯(lián)合反演結(jié)果如圖5所示,其密度結(jié)構(gòu)特征與剖面的地質(zhì)解釋[4](黑實(shí)線)有很好的一致性,基本保留了地震數(shù)據(jù)的界面分層信息。與初始模型相比,反演結(jié)果有較明顯的變化:四川盆地各圈層內(nèi)部密度較為均一,而川西北高原中下地殼密度相對混雜,相鄰各圈層互有滲透;龍門山褶皺帶下方15~30km 處的低密度填充變?yōu)橄鄬Ω呙芏忍畛?,并且其密度與下地殼相近;西北高原地下介質(zhì)整體密度降低,四川盆地中下地殼密度變大,上地殼與川西北高原上地殼密度接近,使得龍門山褶皺帶與四川盆地交界處的斷裂特征更加明顯??傮w而言,剖面地下空間的界面分層信息與地震速度結(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)動受其阻礙,使得各圈層內(nèi)巖石擠壓,沿龍門山斷裂帶上涌,導(dǎo)致川西北高原地下介質(zhì)整體密度降低,相鄰各圈層界線模糊(擠壓作用使巖石破碎),龍門山褶皺帶中上地殼有高密度物質(zhì)填充(下地殼物質(zhì)上涌)。
圖5中存在個別較為明顯、但顯然不可信的局部細(xì)小構(gòu)造特征,這主要是由于反演過程中完全拋棄了人工干涉,并且未對結(jié)果在分辨率尺度上進(jìn)行合適的平滑濾波,導(dǎo)致未能將低于分辨能力的假異常完全抹除。
本文反演所得密度結(jié)構(gòu)圖(圖5)較好地保留了速度結(jié)構(gòu)圖(圖2)劃分的圈層構(gòu)造,且與該剖面的地質(zhì)學(xué)解釋具有很好的一致性。兩者的差異主要體現(xiàn)在少見的“低速-高密度”對應(yīng)關(guān)系:速度結(jié)構(gòu)圖中下地殼的低速通道在密度結(jié)構(gòu)圖中消失不見,龍門山褶皺帶中上地殼的低速異常在密度結(jié)構(gòu)圖中表現(xiàn)為高密度異常。這是由于地震波速的變化對地質(zhì)界面、巖石破碎程度相當(dāng)敏感,而重力對巖性的變化更為敏感。下地殼低密度通道的消失說明了該區(qū)域并無明顯的巖性變化,速度結(jié)構(gòu)圖中低速通道主要由巖石破碎而非巖性變化引起;中上地殼相應(yīng)位置的高密度填充,說明該區(qū)域包含下地殼物質(zhì),由于巖石破碎,其在速度結(jié)構(gòu)圖中體現(xiàn)為低速。重力-地震聯(lián)合反演結(jié)果的這種“低速-高密度”對應(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ì)多個參數(shù)特征對剖面作出更合理、可靠的地質(zhì)解釋。與傳統(tǒng)的人機(jī)交互聯(lián)合反演相比,該方法簡單、快捷,反演過程無需人工參與,受研究者主觀影響較小。遂寧-阿壩剖面的反演結(jié)果也顯示了該方法的不足之處,即無法完全抹除個別不可分辨的假異常。針對這一問題,可以考慮借鑒最小構(gòu)造反演[12]相關(guān)方法進(jìn)行改進(jìn)。
[1]劉皓,方盛明,嘉世旭.地震力聯(lián)合反演及效果——以天津-北京-赤城地震探測剖面為例[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)惠芬,張若水.三維位場的快速反演方法及程序設(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].大地測量與地球動力學(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