基于背景噪聲面波的吉林省東南部S波速度結(jié)構(gòu)成像

董向欣，安百州，畢鳳儀，于平

吉林大學(xué) 地球探測科學(xué)與技術(shù)學(xué)院，長春 130026

0 引言

吉林省地處東北亞地理中心位置，地貌形態(tài)差異明顯，可分為東南部山地和西北部平原兩大地貌。其地質(zhì)構(gòu)造較為復(fù)雜，分布有郯廬斷裂帶北段主干的依蘭—伊通斷裂和郯廬斷裂帶北段分支的敦化—密山斷裂，松遼盆地以及長白山火山、龍崗火山等全新世活動火山。太平洋板塊向西俯沖于歐亞板塊下，俯沖帶在上地幔的水平運移影響了長白山火山、鏡泊湖火山的火山活動[1]。位于吉林省東南部、中朝邊境的長白山火山屬于巨型的板內(nèi)層狀復(fù)式火山[2]，由數(shù)以百計的火山錐和廣袤的熔巖臺地組成了總面積達12 000 km2的龐大火山群，已有多次大規(guī)模噴發(fā)的地質(zhì)歷史，今后可能有中等規(guī)模的噴發(fā)危險[3]。位于吉林省靖宇縣與輝南縣交界處、長白山脈西北側(cè)的龍崗山脈，是長白山脈的支脈，由160余座火山組成了龍崗火山群，是中國最典型的單成因火山作用發(fā)育特征的地區(qū)之一[4]。龍崗火山群的火山活動始于新近紀，鼎盛于第四紀，休眠于距今約1 600 a±。龍崗火山群和長白山火山群分別展現(xiàn)了單成因和多成因的不同火山成因類型。20世紀90年代以來，前人在吉林省東南部地區(qū)利用不同的地球物理學(xué)方法技術(shù)開展了探測研究，包括人工地震[5--6]、天然地震接收函數(shù)[7--10]、大地電磁測深[11--15]、背景噪聲成像等[16--21]，顯示在長白山火山下方存在低速異常、高泊松比異常、低阻異常和重力負異常等特征，但對于長白山火山地殼巖漿房的空間分布位置和部分熔融程度等地質(zhì)學(xué)問題的認識仍存在較大分歧。近年來，王武等[17]和范興利等[18]在長白山地區(qū)的背景噪聲成像研究，探討了長白山火山巖漿房的位置和地殼熔融程度等爭議問題，但所使用的成像方法均為需要生成群速度或相速度圖的中間步驟的兩步法，具有一定的局限性。

利用多個地震臺站的背景噪聲記錄計算互相關(guān)函數(shù)并進行長時間疊加，獲得面波的經(jīng)驗格林函數(shù)并提取群速度或相速度的頻散曲線，最終進行成像獲得地下速度結(jié)構(gòu)，是近年來發(fā)展快速且應(yīng)用廣泛的地震成像方法。背景噪聲的理論可以追溯到Aki[22]的文獻。背景噪聲成像方法同基于地震面波層析成像的傳統(tǒng)方法相比，具有射線分布均勻、不需震源信息、不受短周期面波衰減影響和分辨率高等優(yōu)點，更適用于地震活動較少的地區(qū)。Shapiro[23]最早使用背景噪聲面波信號反演得到美國加州的高分辨率瑞利面波群速度分布圖，隨后該方法逐漸應(yīng)用于全球多個地區(qū)不同尺度的地殼速度的研究。Fang et al.[24]提出的基于射線追蹤的面波走時直接反演橫波速度結(jié)構(gòu)的方法，考慮了復(fù)雜介質(zhì)下的射線路徑彎曲的影響，減小了計算誤差，有利于較小尺度的地質(zhì)構(gòu)造復(fù)雜地區(qū)的面波反演。為研究吉林省東南部地區(qū)S波速度結(jié)構(gòu)，驗證直接反演方法在多火山發(fā)育區(qū)的適用性，本研究使用國家測震臺網(wǎng)的吉林、黑龍江和遼寧臺網(wǎng)的連續(xù)波形數(shù)據(jù)，利用臺站對互相關(guān)函數(shù)獲得瑞利面波的經(jīng)驗格林函數(shù)，提取群速度和相速度頻散曲線，用面波直接成像方法反演得到吉林省東南部深度為0～30 km的S波速度結(jié)構(gòu)圖像。

1 數(shù)據(jù)預(yù)處理及頻散曲線提取

選取國家測震臺網(wǎng)數(shù)據(jù)備份中心[25]共39個固定地震臺站從2019年1月1日至2020年12月31日的連續(xù)地震背景噪聲數(shù)據(jù)。其中吉林測震臺網(wǎng)26個臺站、黑龍江測震臺網(wǎng)7個臺站、遼寧測震臺網(wǎng)6個臺站。研究區(qū)內(nèi)的主要地質(zhì)構(gòu)造和臺站位置分布如圖1所示。儀器類型為長周期儀器，波形采樣率為100 Hz，實際數(shù)據(jù)處理使用了垂直分量的地震記錄。

圖1 研究區(qū)主要地質(zhì)構(gòu)造和臺站分布圖Fig.1 Distribution map of main geological structures and stations in study area

1.1 數(shù)據(jù)預(yù)處理

數(shù)據(jù)預(yù)處理采用了Bensen et al.[26]介紹的由背景噪聲數(shù)據(jù)測量瑞利面波群速度和相速度頻散曲線的方法和處理流程。首先對所有地震臺站的單臺數(shù)據(jù)按天截取為相同長度，重采樣為1 Hz。為降低環(huán)境噪聲和儀器噪聲及人為干擾造成的波形記錄傾斜和零飄現(xiàn)象的影響，對重采樣后的數(shù)據(jù)做了去平均值和去線性趨勢處理。然后去除儀器響應(yīng)，目的是消除不同地震儀器固有工作頻率的干擾。帶通濾波周期范圍為3～40 s。最后做了時域歸一化和頻譜白化處理，時域歸一化的目的是除去波形記錄中能量較強的地震和爆破等干擾，頻譜白化的目的是抑制單頻信號和微震的干擾，拓寬背景噪聲的頻帶。

對任意一個臺站每天的垂直分量波形數(shù)據(jù)進行互相關(guān)運算，選取時間長度為-200～200 s。將該臺站對24個月的單日互相關(guān)運算結(jié)果線性疊加，生成最終的互相關(guān)函數(shù)。 假定有臺站A和臺站B，波場由臺站A激發(fā)，被臺站B接收，臺站A和臺站B記錄的背景噪聲的互相關(guān)函數(shù)對時間求導(dǎo)可獲得臺站B的經(jīng)驗格林函數(shù)。

如果在空間上的噪聲源是均勻分布的，那么背景噪聲記錄的互相關(guān)函數(shù)的正向和負向分支應(yīng)該是嚴格對稱的，這樣求導(dǎo)得到的經(jīng)驗格林函數(shù)與精確的經(jīng)驗格林函數(shù)具有相同的相位信息。由于實際地球介質(zhì)中噪聲源的不均勻空間分布和隨機性，不能保證由實際背景噪聲記錄的互相關(guān)函數(shù)的正向和負向分支的對稱性。因此，用背景噪聲記錄的互相關(guān)函數(shù)的正向和負向分支反序平均的結(jié)果求導(dǎo)計算經(jīng)驗格林函數(shù)。這同時可以提高計算結(jié)果的信噪比。圖2是所有臺站從2019年1月疊加至2020年12月(共24個月)得到的經(jīng)驗格林函數(shù)。圖中的結(jié)果顯示了清晰的面波信號，具有較高的對稱性和信噪比。

紅色虛線之間為群速度為2～5 km的瑞利面波信號范圍，紫色虛線外側(cè)為噪聲范圍。圖2 所有臺站經(jīng)驗格林函數(shù)Fig.2 Empirical Green function of all stations

1.2 群速度和相速度頻散曲線提取

本研究采用姚華建等[27]提出的一種基于圖像分析的雙臺面波頻散曲線快速提取方法和由此編寫的Matlab交互式處理程序(EGFAnalysisTimeFreq_version_2015)，該程序設(shè)計了一個圖形用戶界面(GUI)，可快速準確地從經(jīng)驗格林函數(shù)或互相關(guān)函數(shù)中提取出雙臺瑞利面波群速度和相速度頻散曲線。筆者以遼寧地震臺網(wǎng)的西豐臺(LN.XFN)和吉林地震臺網(wǎng)的榆樹臺(JL.YST)為例，介紹了頻散曲線提取的主要過程。

如圖3所示，生成西豐臺和榆樹臺的瑞利面波群速度頻譜能量圖。為滿足遠場假設(shè)，獲得理論上真實可靠的頻散曲線，應(yīng)用了臺站間沿大圓路徑的距離至少為2倍瑞利面波波長的質(zhì)量控制標準。圖3中的紅色直線表示臺站間沿大圓路徑的距離等于2倍波長時對應(yīng)的周期，紅色的實心圓表示滿足臺站間距不小于2倍波長的的質(zhì)量控制標準的頻散點。

同時，為提高測量精度，獲得盡可能準確的頻散曲線，應(yīng)用了滿足測量信噪比(SNR)的質(zhì)量控制標準。信噪比的定義參考Bensen et al.[26]，為一定周期或頻率的面波信號時窗內(nèi)包絡(luò)的最大振幅與一段時間后的噪聲時窗內(nèi)包絡(luò)的振幅均方根之比。本研究的噪聲時窗在信號時窗之后的100 s選取。圖3中的藍色圓環(huán)表示信噪比至少為5的瑞利面波群速度和相速度頻散點。

a.瑞利面波群速度頻散曲線測量圖，黑色曲線為各周期的信號最大振幅的連線，紅色圓點表示各周期的信號最大振幅所在的位置，藍色圓環(huán)表示提取的滿足信噪比要求的頻散點，紅色直線表示臺站間沿大圓路徑的距離等于2倍波長時對應(yīng)的周期；b.瑞利面波相速度頻散曲線測量圖，黑色曲線、紅色圓點、藍色圓環(huán)和紅色直線的意義與圖3a中的相同。紫色曲線為圖3a中提取的群速度頻散曲線，綠色曲線為根據(jù)保留的相速度頻散點計算得到的群速度頻散曲線。圖3 瑞利面波群速度和相速度頻散曲線測量示例Fig.3 Example of Rayleigh wave group velocity and phase velocity dispersion curve measurement

根據(jù)以上兩個質(zhì)量控制標準，舍棄了小于2倍波長的臺站間距對應(yīng)的周期和信噪比小于5的周期對應(yīng)的頻散點，只保留了圖3a中紅色實心圓和藍色圓環(huán)同時所在的群速度頻散點，以保證群速度測量的可靠性。

然后，基于已保留的群速度頻散點，提取相速度頻散曲線。在相速度頻譜能量圖(圖3b)中存在多個相速度頻散的分支，只選擇離群速度頻散曲線最接近且能量最強的分支用來提取相速度頻散曲線。同樣地，只保留圖3b中紅色實心圓和藍色圓環(huán)同時所在的相速度頻散點。紫色曲線表示上一步提取的頻散曲線，綠色曲線表示根據(jù)保留的相速度頻散點計算得到的群速度頻散曲線，兩者較為接近，說明根據(jù)已保留的群速度頻散點提取的相速度頻散曲線可靠。

1.3 射線路徑分布

本研究共選取研究區(qū)域內(nèi)39個臺站，按照n(n-1)/2 公式計算，理論上可以得到741條頻散曲線。經(jīng)過1.2節(jié)的數(shù)據(jù)篩選，實際反演中拾取674條瑞利面波相速度的混合路徑頻散曲線進行面波直接成像，射線路徑分布如圖4所示。

圖4 射線分布圖Fig.4 Ray distribution map

2 S波速度結(jié)構(gòu)反演

2.1 反演方法

采用Fang et al.[24]提出的基于射線追蹤的三維淺層地殼結(jié)構(gòu)的面波頻散直接反演方法，反演吉林省東南部地區(qū)三維地下S波速度結(jié)構(gòu)。相比于傳統(tǒng)的兩步法，如Barmin et al.[28]提出的面波成像方法，該直接反演方法不需要生成群速度圖或相速度圖的步驟，使用基于頻率相關(guān)的射線追蹤和基于小波的稀疏約束反演。射線追蹤使用Rawlinson et al.[29]提出的快速行進法，計算每個周期的面波旅行時和源與接收器之間的射線路徑，避免了在大多數(shù)面波層析研究中使用的大圓傳播的假設(shè)?；谛〔ǖ南∈杓s束反演[30]考慮了具有平滑變化地震特性的準層狀介質(zhì)，通過網(wǎng)格點下的一維剖面來表示三維剪切波速度模型，這些網(wǎng)格點是使用基于小波的稀疏約束成像方法從所有頻散數(shù)據(jù)中同時確定的。波速模型的小波系數(shù)用迭代加權(quán)最小二乘算法來估計，并且在迭代時，使用新獲得的波速模型來更新面波射線路徑和數(shù)據(jù)靈敏度矩陣。

2.2 檢測板測試

為檢查反演參數(shù)設(shè)置的合理性，評價反演模型的分辨能力，通常使用模型的檢測板恢復(fù)性測試，對臺站分布與射線路徑覆蓋程度的影響進行檢測。設(shè)置的檢測板棋盤模型深度范圍為地表至62.5 km，沿深度方向劃分為26層，每層厚度為2.5 km，每層的S波速度為1.4 km+深度×0.05 km，根據(jù)檢測板分辨率測試的原理，使用正弦函數(shù)在每層的S波相同速度值上添加了一定的速度擾動量，速度擾動量的幅度范圍為每層的S波相同速度值±0.4 km。如圖5所示，棋盤模型的平面沿經(jīng)度和緯度方向劃分成15×14共210個網(wǎng)格點，網(wǎng)格大小為0.5°×0.5°。在正演時，先按照實際射線路徑計算理論走時，由于受噪聲源分布不均勻因素造成的相速度頻散誤差約為1%[31]，在理論數(shù)據(jù)中加入2%的隨機誤差生成合成數(shù)據(jù)，最后對合成數(shù)據(jù)反演，獲得不同深度的檢測板測試結(jié)果。如圖5的紅色矩形所示，反演結(jié)果的范圍比模型的范圍分別在經(jīng)緯度和深度方向上少1°和1層深度。

圖5給出了初始棋盤模型及不同深度的測試結(jié)果。檢測結(jié)果顯示深度為5～30 km時，反演的棋盤圖案基本重現(xiàn)了給定的初始模型和速度擾動，特別對于研究區(qū)域的中心位置，成像結(jié)果的可信度較高。但模型邊緣重建結(jié)果相對較差，這與缺乏射線路徑覆蓋有關(guān)。隨著深度的增加，反演重建效果逐漸變差，在深度為35～40 km的圖像中表現(xiàn)得較為明顯，這是由于長周期的瑞利面波的相速度對相應(yīng)深度的S波速度的敏感性降低導(dǎo)致。

2.3 反演和討論

一般來說，瑞利面波的相速度對其三分之一個波長左右對應(yīng)的地下深度處的S波速度結(jié)構(gòu)最為敏感。筆者的反演初始速度模型根據(jù)拾取的瑞利面波相速度頻散點一維平均S波速度模型設(shè)置，深度范圍為地表至62.5 km，沿深度方向劃分為26層，每層厚度為2.5 km。模型的平面略大于研究區(qū)平面，沿經(jīng)度和緯度方向劃分成15×14共210個網(wǎng)格點，網(wǎng)格大小為0.5°×0.5°。使用拾取的符合質(zhì)量控制標準的相速度頻散作為輸入數(shù)據(jù)，經(jīng)過反復(fù)試驗，選取了合適的平滑參數(shù)和阻尼參數(shù)，迭代次數(shù)為10次，反演得到了吉林省東南部地區(qū)三維地下S波速度結(jié)構(gòu)。

根據(jù)檢測板可靠性分析結(jié)果，選取研究區(qū)深度為5 km、10 km、15 km、20 km、25 km、30 km的S波速度分布圖像進行討論(圖6)。

由圖6可以看出，研究區(qū)存在造山帶、沉積盆地、斷裂帶及火山等多種地質(zhì)構(gòu)造，同一深度的S波速度值不均勻分布，且速度異常變化范圍較大，表明吉林省東南部地區(qū)中上地殼S波速度結(jié)構(gòu)存在較強的橫向不均勻性。5 km的成像結(jié)果與地表地質(zhì)構(gòu)造相對應(yīng)，速度異常呈現(xiàn)出受郯廬斷裂帶北延的依蘭—伊通斷裂和敦化—密山斷裂等構(gòu)造控制的區(qū)塊化特征。依蘭—伊通斷裂西北方的松遼盆地在上地殼深度表現(xiàn)為整體的低速異常，可能是由于其具有厚度較大的松散沉積層；東南方的廣大山脈地區(qū)表現(xiàn)為大面積的高速異常。

隨著深度的增加(10～30 km)，松遼盆地在中下地殼的低速異常逐漸減弱并有向高速異常轉(zhuǎn)變的趨勢?？赡苁桥璧刂醒刖哂休^大的沉積厚度，邊緣具有較淺的結(jié)晶基底的原因。

長白山火山的低速異常逐漸顯現(xiàn)，在15～20 km表現(xiàn)明顯(朝鮮一側(cè)可能也具有低速異常)并有向下延伸的趨勢，可能是火山下方中部地殼發(fā)育的巖漿房或者地殼深部巖漿向淺部運移的通道，這和范興利等[18]的研究認識相一致。

龍崗火山在上地殼深度(5 km)不具有低速異常，在中下地殼深度(20～30 km)顯現(xiàn)出微弱的低速異常，且位于低速異常范圍的邊緣。從15 km、25 km和30 km的圖像上可以看出，龍崗火山和長白山火山的低速異常處于同一個較大的異常范圍中，但與后者相比程度微弱，這可能是由于龍崗火山巖漿噴發(fā)后的殘留物導(dǎo)致的。

在5～20 km深度，鏡泊湖火山在中上地殼深度也不具有低速異常。但從25 km、30 km的圖像上可以看出，鏡泊湖火山及其東北方在下地殼深度顯示出一定的低速異常，原因可能是火山下方存在程度較低的部分熔融。鏡泊湖火山的速度結(jié)構(gòu)，與龍崗火山、長白山火山的速度結(jié)構(gòu)沒有直接聯(lián)系，而與敦化—密山斷裂有關(guān)。從巖漿類型來看，鏡泊湖火山的噴出物無地殼混染，與長白山火山的噴出物類型不同。

為更好地對比和說明長白山火山和龍崗火山的S波速度結(jié)構(gòu)，筆者在長白山火山附近沿緯度42°N和經(jīng)度128°E分別作了2條剖面，對應(yīng)圖6a中的AA’和BB’線段(圖7)。

黑色矩形區(qū)域為棋盤模型的范圍，紅色矩形區(qū)域為反演結(jié)果的范圍，黑色多邊形區(qū)域為射線路徑覆蓋的范圍。a.5 km深度；b.5 km深度；c.10 km深度；d.15 km深度；e.20 km深度；f.25 km深度；g.30 km深度；h.35 km深度；i.35 km深度。圖5 檢測板測試結(jié)果Fig.5 Test results of check board

黑色三角代表研究區(qū)內(nèi)的新生代火山，紅色線段代表圖7的剖面位置。a.5 km深度；b.10 km深度；c.15 km深度；d.20 km深度；e.25 km深度；f.30 km深度。圖6 S波速度層析成像結(jié)果Fig.6 S-wave velocity tomography results

圖7 S波速度結(jié)構(gòu)剖面圖Fig.7 S-wave velocity profiles

由圖7可以看出，長白山下方的低速異常在地下6～23 km深度范圍，分為一個主要的低速異常體和一個較小的低速異常體，具有向下延伸的趨勢。在12 km深度附近含有少量的相對高速異常。龍崗火山下方7～15 km深度范圍也存在低速異常，與相鄰的長白山火山下方的低速異常是連通的，但沒有向下延伸的趨勢，說明龍崗火山下方可能不存在地殼巖漿囊。

由于廣泛且快速的熔融物分離發(fā)生于深部地殼的熱流條件下，深部陸殼因此是巖漿演化和存儲最理想的位置[32]。結(jié)合圖6c、圖6d和圖7，筆者認為長白山火山的地殼巖漿房主要位于火山的正下方，深度為16～21 km。這支持了Kyong--Song et al.[33]提出的地殼的很大一部分已被巖漿作用改造，與部分熔融相關(guān)的地殼巖漿體可能直接位于火山下方的說法。另外，Zhang et al.[34]觀察到巖漿從基性噴發(fā)到中基性噴發(fā)的成分變化。巖漿房中可能發(fā)生進一步的巖漿分離和化學(xué)分異，導(dǎo)致中酸性熔融物的產(chǎn)生?；瘜W(xué)分異后，來自下地殼巖漿房的密度較低、長英質(zhì)的巖漿可能會上升至淺部，形成橫向和縱向延伸的低速體。因此，從下地殼巖漿房運移的差異巖漿可能是目前在長白山火山觀察到豐富的近地表熱液活動的原因。

3 結(jié)論

(1)吉林省東南部地區(qū)地殼深度的S波速度結(jié)構(gòu)存在較強的橫向不均勻性。上地殼的速度結(jié)構(gòu)與造山帶、沉積盆地、斷裂帶及火山等多種地表地質(zhì)構(gòu)造相對應(yīng)，呈現(xiàn)出區(qū)塊化特征。

(2)長白山火山的地殼巖漿房的可能深度在6～23 km，且位于火山正下方，與地幔熱物質(zhì)上涌有關(guān)。

(3)龍崗火山下方可能不存在活動巖漿囊。鏡泊湖火山下方可能存在與長白山火山相比程度較低的部分熔融。

(4)成像結(jié)果能較好地反映出研究區(qū)地殼深度內(nèi)S波速度結(jié)構(gòu)的特征，直接反演方法在多火山發(fā)育的復(fù)雜構(gòu)造區(qū)域具有良好的適用性。