范建柯丁志峰郭慧麗蘇道磊張斌，

1） 中國(guó)北京 100081 中國(guó)地震局地球物理研究所

2） 中國(guó)山東青島 266071 中國(guó)科學(xué)院海洋研究所

3） 中國(guó)山東青島 266071 中國(guó)科學(xué)院海洋地質(zhì)與環(huán)境重點(diǎn)實(shí)驗(yàn)室

4） 中國(guó)濟(jì)南 250001 濟(jì)南市地震監(jiān)測(cè)中心

引言

山東地區(qū)位于中國(guó)大陸東部，自中生代以來經(jīng)歷了復(fù)雜的地質(zhì)構(gòu)造運(yùn)動(dòng)，如中生代早期華北地塊與揚(yáng)子地塊的碰撞、中生代中晚期的華北克拉通破壞以及新生代以來的拉張變形等.陸地區(qū)域內(nèi)發(fā)育沂沭斷裂帶（郯廬斷裂帶山東段）和聊考斷裂帶兩大斷裂系統(tǒng)，北部海域發(fā)育張家口——蓬萊斷裂帶，導(dǎo)致地震活動(dòng)頻繁，并伴有大地震的發(fā)生，如1668年郯城M8.5地震，造成了重大的人員傷亡和經(jīng)濟(jì)損失.

2003年6月5日和2020年2月18日，山東省最大的兩個(gè)城市青島和濟(jì)南分別發(fā)生了一次M4.1地震（圖1），雖然未導(dǎo)致較大的人員傷亡和經(jīng)濟(jì)損失，但由于其發(fā)生于人口均近千萬的兩大城市，且有感范圍較大，余震次數(shù)較多，而上述兩個(gè)地區(qū)通常被認(rèn)為屬于“弱震”區(qū)，因此，這兩次地震仍然引發(fā)了非常大的關(guān)注.前人綜合余震重定位結(jié)果、震源機(jī)制解等數(shù)據(jù)，推測(cè)這兩次地震可能屬于相對(duì)完整巖體條件下的一次新破裂活動(dòng)（潘元生等，2004）或區(qū)域構(gòu)造應(yīng)力作用下附近斷裂或次級(jí)派生斷裂活動(dòng)的結(jié)果（張斌等，2020）.另外，自2000年以來，山東半島地區(qū)相繼發(fā)生了嶗山震群、乳山震群和長(zhǎng)島震群，但目前上述發(fā)震區(qū)尚未有高精度的層析成像數(shù)據(jù)，速度結(jié)構(gòu)與發(fā)震機(jī)理仍不清楚.因此，有必要查明上述地區(qū)的淺部及深部速度結(jié)構(gòu)形態(tài)，探討速度結(jié)構(gòu)與地震發(fā)震機(jī)理之間的關(guān)系，為后期的防震減災(zāi)提供理論支持.

圖1 山東及周邊地區(qū)構(gòu)造分區(qū)圖（修改自蘇道磊等，2016）Fig.1 Tectonic settings in and around Shandong area（modified from Su et al，2016）

1 數(shù)據(jù)與方法

本文所用數(shù)據(jù)來源于山東及鄰區(qū)96個(gè)臺(tái)站（圖2a）記錄的區(qū)域內(nèi)地震事件的P波和S波到時(shí)數(shù)據(jù)，主要分為兩部分：第一部分為1975年1月至2014年1月發(fā)生的天然地震事件（蘇道磊等，2016）；第二部分為2016年1月至2019年12月發(fā)生的天然地震事件.原始數(shù)據(jù)共包含7 271 個(gè)地震事件.

圖2 臺(tái)站（a）及地震（b）分布圖Fig.2 Distribution of seismic stations （a） and local earthquakes （b） used in this study

為保證反演結(jié)果的準(zhǔn)確性，需要對(duì)地震事件進(jìn)行嚴(yán)格挑選，設(shè)定標(biāo)準(zhǔn)如下：① 每個(gè)地震事件至少被4個(gè)臺(tái)站接收；② 震相走時(shí)殘差絕對(duì)值小于3.0 s；③ 重定位前后，發(fā)震時(shí)刻偏差小于2.0 s，水平向偏差小于6 km，震源深度偏差小于8 km.最終篩選出4 652個(gè)地震事件，包括3萬6 482個(gè)P波震相到時(shí)和3萬2 600個(gè)S波震相到時(shí)用于成像反演（圖2b）.

本文采用TOMOG3D三維層析成像反演方法（Zhaoet al，1992）.該方法采用三維網(wǎng)格節(jié)點(diǎn)表示空間的速度分布，空間內(nèi)每一點(diǎn)的速度值由周圍八個(gè)節(jié)點(diǎn)的速度值進(jìn)行線性插值獲得，允許三維空間內(nèi)存在間斷面，通過偽彎曲算法快速準(zhǔn)確地計(jì)算射線路徑和走時(shí)，并能夠同時(shí)處理近遠(yuǎn)震及后續(xù)震相到時(shí)數(shù)據(jù).

合適的一維初始速度模型對(duì)最終反演結(jié)果的準(zhǔn)確性至關(guān)重要.蘇道磊等（2016）分別測(cè)試了三個(gè)初始速度模型：① 魯西地區(qū)地殼速度模型，同時(shí)考慮莫霍面起伏（嘉世旭，張先康，2005）；② 山東地震臺(tái)網(wǎng)定位用地殼速度模型；③ 華北地區(qū)地殼速度模型（陳立華，宋仲和，1990）.測(cè)試結(jié)果表明，魯西地區(qū)地殼速度模型的走時(shí)殘差均方根最小，也更符合真實(shí)地層情況，而華北地區(qū)地殼速度模型走時(shí)殘差遠(yuǎn)大于其它兩個(gè)速度模型.因此，本文僅對(duì)比前兩個(gè)初始速度模型的走時(shí)殘差，并同時(shí)考慮莫霍面起伏的影響.其中，莫霍面埋深數(shù)據(jù)來自于CRUST1.0模型（Laskeet al，2013）和鄭宏等（2021）利用接收函數(shù)反演得到的山東地區(qū)莫霍面埋深等（圖3）.初始S波速度模型由P波速度除以1.732得到.

圖3 基于不同模型反演所得的莫霍面埋深（a，b）及本研究所用莫霍面埋深（c）Fig.3 Moho depths inverted from different models （a，b） and used in this study （c）

經(jīng)過計(jì)算，上述兩個(gè)初始速度模型的殘差均方根分別為0.841 s和0.860 s，走時(shí)殘差分布如圖4所示.對(duì)比發(fā)現(xiàn)，魯西地區(qū)速度模型要優(yōu)于山東地震臺(tái)網(wǎng)定位用地殼速度模型，與蘇道磊等（2016）的結(jié)論相似.因此，后續(xù)的層析成像反演采用魯西地區(qū)速度模型，并考慮莫霍面起伏的影響.

圖4 使用魯西地區(qū)地殼速度模型（a）和山東地震臺(tái)網(wǎng)定位用地殼速度模型（b）在地震重定位前后所得的走時(shí)殘差分布圖Fig.4 Travel time residual distribution based on crustal velocity model in Luxi area （a） and Shandong seismic network （b） before and after relocation

本研究中，初始模型三維網(wǎng)格節(jié)點(diǎn)橫向間隔為0.2°×0.2°，縱向上分別在1，10，20和30 km深度處設(shè)置節(jié)點(diǎn)層.反演前，基于初始速度模型和原始到時(shí)數(shù)據(jù)對(duì)所有地震事件進(jìn)行了重新定位.重定位前后的走時(shí)殘差分布如圖4所示，重定位后的總體走時(shí)殘差均方根由0.841 s降低到 0.649 s，其中 P 波走時(shí)殘差均方根由 0.792 s降低到 0.620 s，S 波走時(shí)殘差均方根由0.893 s降低到0.680 s，說明地震定位精度有了較大幅度的提高.走時(shí)殘差絕大多數(shù)位于±2 s以內(nèi)，因此選擇走時(shí)殘差在±2 s以內(nèi)的震相參與最終的成像反演.反演采用帶阻尼和平滑因子的最小二乘（least squares QR-factorization，縮寫為 LSQR）方法（Paige，Saunders，1982）得到最終的三維速度結(jié)構(gòu).通過大量的測(cè)試，P波和S波速度反演采用的最優(yōu)阻尼和平滑因子均為 5.0 和 50.0 （圖5）.

圖5 P 波（a，b）和 S 波（c，d）成像中不同阻尼和平滑因子對(duì)應(yīng)的三維速度模型標(biāo)準(zhǔn)差與走時(shí)殘差均方根關(guān)系曲線Fig.5 Trade-off curves between the standard deviation of the 3-D velocity model and the root-mean-square of travel time residual according to the damping （a and c） and smoothing （b and d） parameters for P- （a，b） and S-wave （c，d） tomographies

2 分辨率測(cè)試

在分析成像結(jié)果之前，需要對(duì)成像結(jié)果的可靠性進(jìn)行評(píng)估，一般采用檢測(cè)板（checkerboard）方法（Zhaoet al，1992）.在該方法中，首先將三維空間內(nèi)相鄰網(wǎng)格節(jié)點(diǎn)分別設(shè)置±3%的速度擾動(dòng)，然后利用相同的地震和臺(tái)站分布計(jì)算理論走時(shí)，并在計(jì)算理論走時(shí)過程中加入標(biāo)準(zhǔn)差為0.1 s的隨機(jī)誤差來檢測(cè)計(jì)算穩(wěn)定性.隨后對(duì)得到的理論走時(shí)在初始一維速度模型的基礎(chǔ)上進(jìn)行反演，通過對(duì)比三維網(wǎng)格節(jié)點(diǎn)處擾動(dòng)值的恢復(fù)情況對(duì)成像結(jié)果的可靠性進(jìn)行評(píng)估.如果反演后擾動(dòng)值與反演前擾動(dòng)值的分布相似，則說明分辨率較好.

圖6顯示了橫向間隔分別為 0.33°×0.33°，0.4°×0.4°和 0.5°×0.5°的分辨率測(cè)試結(jié)果.結(jié)果顯示，對(duì)于P波和S波成像結(jié)果，研究區(qū)內(nèi)大部分地區(qū)的分辨率能達(dá)到0.4°×0.4°，部分地區(qū)如魯西南地區(qū)和山東半島地區(qū)可達(dá)到0.33°×0.33°的分辨率.其中：網(wǎng)格間距為0.33°×0.33°時(shí)，P波和S波檢測(cè)板擾動(dòng)值恢復(fù)準(zhǔn)確性分別達(dá)到74%和81%，振幅恢復(fù)達(dá)到70%的比例約為47%和59%，振幅恢復(fù)達(dá)到100%的比例約為32%和41%；網(wǎng)格間距為0.4°×0.4°時(shí)，P波和S波檢測(cè)板擾動(dòng)值恢復(fù)準(zhǔn)確性分別達(dá)到76%和82%，振幅恢復(fù)達(dá)到70%的比例約為52%和62%，振幅恢復(fù)達(dá)到100%的比例約為36%和41%；網(wǎng)格間距為0.5°×0.5°時(shí)，P波和S波檢測(cè)板擾動(dòng)值恢復(fù)準(zhǔn)確性分別達(dá)到77%和82%，振幅恢復(fù)達(dá)到70%的比例約為55%和64%，振幅恢復(fù)達(dá)到100%的比例約為37%和42%.

圖6 檢測(cè)板分辨率測(cè)試結(jié)果Fig.6 Results of checkerboard resolution tests

3 結(jié)果與討論

經(jīng)過反演，P波和S波走時(shí)殘差均方根由反演前的0.620 s和0.680 s分別降低到0.378 s和0.417 s，減少近40%.利用P波、S波反演結(jié)果計(jì)算得到了泊松比成像結(jié)果，如圖7所示.蘇道磊等（2016）的P波層析成像結(jié)果揭示出：在1——10 km深度切片上，沂沭斷裂帶沿構(gòu)造走向表現(xiàn)出強(qiáng)烈的橫向不均一性，高低速異常交替出現(xiàn)，膠東半島北部海域、膠萊盆地和濟(jì)陽坳陷主要表現(xiàn)為低速異常；20——30 km深度，魯西地區(qū)存在較大范圍的低速異常（蘇道磊等，2016）.上述結(jié)果與本文的P波成像結(jié)果非常一致（圖7a）.另外，本文的S波成像結(jié)果與Li等（2018）基本一致，如在10 km深度上，沂沭斷裂帶北部以低速異常為主（圖7b），20——30 km山東半島表現(xiàn)為大范圍的低速異常（圖7a）.但由于本文應(yīng)用了更多的P波和S波到時(shí)數(shù)據(jù)，因此本文的成像結(jié)果具有更高的分辨率，如在1 km和10 km深度，沂沭斷裂帶南部西側(cè)的兩個(gè)低速異常區(qū)相互分離，斷裂帶東側(cè)表現(xiàn)為高速異常（圖7a），這在蘇道磊等（2016）的研究中揭示得并不明顯.同時(shí)，本文的檢測(cè)板測(cè)試結(jié)果也要優(yōu)于蘇道磊等（2016），魯西南和山東半島地區(qū)分辨率能夠達(dá)到0.33°×0.33°.山東半島地區(qū)中下地殼存在大范圍低速異常，說明該地區(qū)可能存在較強(qiáng)烈的地幔上涌（李志偉等，2006）.本文同時(shí)利用Liu和Zhao （2018）提出的方法計(jì)算了1——30 km的P波、S波和泊松比的平均值.平均泊松比異常與鄭宏等（2021）利用接收函數(shù)反演得到的山東地區(qū)的泊松比分布基本一致，進(jìn)一步說明了本文成像結(jié)果的可靠性.

圖7 山東地區(qū)不同深度 P 波速度擾動(dòng)（a）、S 波速度擾動(dòng)（b）及泊松比相對(duì)擾動(dòng)（c）結(jié)果底部3幅圖分別為1——30 km的P波、S波速度擾動(dòng)及泊松比相對(duì)擾動(dòng)的平均值Fig.7 Map views of P-wave （a），S-wave （b） and Poisson’s ratio （c） perturbations in Shandong area The averages of P-wave，S-wave and Poisson’s ratio anomalies from 1 to 30 km are shown at the bottom

泰山作為山東地區(qū)的最高點(diǎn)，自新生代以來經(jīng)歷了多期快速抬升（李理，鐘大賚，2006）.接收函數(shù)研究結(jié)果顯示，泰山地區(qū)具有較薄的地殼厚度（32 km左右）和較大的泊松比（0.27），且沒有明顯的方位各向異性（鄭宏等，2021）.本文的成像結(jié)果顯示，泰山地區(qū)下方存在明顯的低速異常（圖8左下四幅圖），這說明該地區(qū)仍存在較強(qiáng)的地幔上涌，導(dǎo)致泰山新生代以來出現(xiàn)顯著抬升，并且現(xiàn)在仍處于抬升階段（鄭宏等，2021）.另外，在淺部1——10 km深度（圖8左上四幅圖），泰山北部為高速異常，南部為低速異常，這與泰山的巖性和地質(zhì)構(gòu)造是一致的：以泰山山前斷裂為界，斷裂北側(cè)為泰山主體，主要出露前寒武紀(jì)泰山變質(zhì)雜巖，以隆升為主；山前斷裂南側(cè)以沉降為主，主要為泰安——萊蕪盆地巨厚的第三系碎屑巖及第四系沉積物（李理，鐘大賚，2006）.

圖8 長(zhǎng)清地震及泰山周邊不同深度處 P 波速度擾動(dòng)（a）、S 波速度擾動(dòng)（b）及泊松比相對(duì)擾動(dòng)（c）結(jié)果Fig.8 Map views of P-wave （a），S-wave （b） and Poisson’s ratio （c） perturbations around Changqing earthquake and Mount Tai

2020年2月18日，濟(jì)南長(zhǎng)清發(fā)生了一次MS4.1地震，震源深度約為2.7 km，并引發(fā)近40次的余震活動(dòng)（張斌等，2020）.本文的成像結(jié)果顯示，該地震震中位于P波、S波高低速異常和高低泊松比異常過渡帶，可能與該地震有關(guān)的長(zhǎng)清斷裂也位于P波、S波和泊松比高低異常過渡帶（圖8上半部和圖9）.震中東部的高速異常與該地區(qū)的背景噪聲成像結(jié)果一致，可能與濟(jì)南侵入巖體有關(guān)（雷霆，2020）.震源機(jī)制解顯示，該地震具有正斷兼走滑性質(zhì)，其余震多為走滑性質(zhì)（張斌等，2020）.另外，余震重定位結(jié)果顯示，主震兩側(cè)的破裂呈明顯不對(duì)稱分布，以西北側(cè)破裂為主（張斌等，2020）.因此，濟(jì)南長(zhǎng)清MS4.1地震可能是區(qū)域構(gòu)造應(yīng)力下長(zhǎng)清斷裂發(fā)生左旋走滑運(yùn)動(dòng)的結(jié)果（張斌等，2020），這與本文的層析成像結(jié)果相一致.另外，主震的西北向不對(duì)稱破裂的形成可能是由于震中東側(cè)高速異常體的存在阻礙了地震的東向破裂.

圖9 過長(zhǎng)清地震的剖面 AA′（a）的 P 波速度擾動(dòng)（b）、S 波速度擾動(dòng)（c）和泊松比擾動(dòng)（d）Fig.9 Vertical cross section of P-wave （b），S-wave （c） and Poisson’s ratio （d） perturbations along a profile AA′ （a） across the Changqing earthquake shown in the inset map

2003年6月5日青島嶗山地區(qū)發(fā)生了ML4.1地震（圖10a左側(cè)），2004年11月1日附近地區(qū)又發(fā)生一次ML3.6地震，這兩次地震都形成了震群序列.根據(jù)定位結(jié)果，這兩次震群序列大致呈北西走向，與附近的主要斷裂近似垂直、與次級(jí)斷裂展布方向基本一致但并不重合，可能屬于相對(duì)完整巖體條件下的一次新的破裂（潘元生等，2003，2005）.根據(jù)本文的成像結(jié)果，這兩次地震及震群序列發(fā)生在P波高低速異常過渡帶、S波低速異常區(qū)和高低泊松比過渡區(qū)且偏向于高泊松比區(qū)域（圖10a——c 1km），同時(shí)深部存在明顯的P波和S波低速異常（圖10a，b 10——30 km），這與1995年日本神戶地震（ML7.2）震源區(qū)的速度和泊松比結(jié)構(gòu)非常相似（Zhaoet al，1996）.另外，震中附近水資源豐富，同時(shí)為花崗巖侵入形成地區(qū)，且構(gòu)造裂隙發(fā)育（趙廣濤等，1996）.因此本文認(rèn)為，可能是深部的流體填充相對(duì)完整巖體內(nèi)的裂隙并引發(fā)破裂，從而導(dǎo)致地震的發(fā)生.需要注意的是，雖然這兩次地震及震群序列與本地區(qū)北東走向的主斷層沒有直接關(guān)系，但本文的成像結(jié)果顯示，該地區(qū)北東走向的主斷層基本位于低速區(qū)內(nèi)以及高低泊松比異常的過渡帶（圖10），因此仍需警惕該地區(qū)發(fā)生震級(jí)更大的中強(qiáng)地震的可能性.

乳山地區(qū)位于大別——蘇魯超高壓變質(zhì)帶的東部，可能存在大量的中生代巖漿侵入體（郭敬輝等，2005）.自2013年10月開始，記錄到的乳山震群地震數(shù)已超過1萬余次，最大震級(jí)為M5.0，目前震群活動(dòng)可能已結(jié)束.震群周邊的斷裂大多以NE-NNE向?yàn)橹鳎诘闹囟ㄎ唤Y(jié)果和震源機(jī)制解顯示，乳山地震序列主要呈NW向展布，絕大多數(shù)地震發(fā)生于中上地殼，發(fā)震斷層可能為傾角近直立的左旋走滑斷層，與區(qū)域內(nèi)的主要斷裂并不一致（張斌等，2017）.本文的成像結(jié)果顯示，乳山震群周邊1 km深度以高速異常和高泊松比異常為主，10 km深度主要為高低速異常過渡區(qū)，20——30 km深度以低速異常和低泊松比異常為主（圖10）.另外，震群附近地?zé)豳Y源較為豐富（田粟，2012）.因此，本文認(rèn)為，乳山地區(qū)深部存在熱地幔物質(zhì)上涌，內(nèi)部含有的流體注入淺部相對(duì)較完整的侵入體內(nèi)或侵入體之間，在區(qū)域應(yīng)力場(chǎng)的作用下導(dǎo)致侵入體的破裂或侵入體之間隱伏斷裂的活動(dòng)，從而引發(fā)乳山震群.這與前人對(duì)乳山震群進(jìn)行的震源譜參數(shù)反演、震中空間分布、震中隨時(shí)間的演化規(guī)律以及震源區(qū)應(yīng)力狀態(tài)的研究結(jié)果是一致的（鄭建常等，2016；王鵬，2019）.

圖10 山東半島地區(qū)不同深度 P 波速度擾動(dòng)（a）、S 波速度擾動(dòng)（b）及泊松比相對(duì)擾動(dòng)（c）結(jié)果Fig.10 Map views of P-wave （a），S-wave （b） and Poisson’s ratio （c） perturbations in and around Shandong Peninsula

長(zhǎng)島震群的發(fā)生可能主要受控于NWW向的張家口——蓬萊斷裂，這是一條深大斷裂，可能已切穿莫霍面甚至巖石圈，成為地幔熱物質(zhì)或基性物質(zhì)上涌的通道（張嶺等，2007），斷裂主要以正斷兼走滑運(yùn)動(dòng)為主（索艷慧等，2013）.歷史上，附近海域曾發(fā)生過M6.0和M7.0左右的大地震（王志才等，2006）.本文的成像結(jié)果顯示（圖10），張家口——蓬萊斷裂帶的地殼速度結(jié)構(gòu)特征與沂沭斷裂帶具有非常大的相似性（蘇道磊等，2016），地殼速度和泊松比結(jié)構(gòu)在1——10 km深度處具有非常強(qiáng)烈的橫向不均一性，斷裂一側(cè)高低速異常交替分布，斷裂帶整體位于高低速異常與高低泊松比異常的過渡帶，20 km以下以低速和低泊松比異常為主，但部分地區(qū)仍有顯著的高泊松比異常，可能反映了地幔熱物質(zhì)的上涌.速度和泊松比異常在10 km左右發(fā)生明顯變化，這與前人得到的b值在8.5 km左右發(fā)生轉(zhuǎn)折是一致的（申金超等，2019）.本文認(rèn)為，深部地幔熱物質(zhì)沿?cái)嗔褞嫌克a(chǎn)生的構(gòu)造應(yīng)力，和/或地幔熱物質(zhì)內(nèi)部含有的流體沿?cái)嗔焉嫌炕蚯秩肓严秾?dǎo)致了長(zhǎng)島震群、甚至周邊地區(qū)的強(qiáng)震活動(dòng).

4 結(jié)論

利用山東及鄰區(qū)的地震臺(tái)站記錄的P波和S波到時(shí)數(shù)據(jù)反演獲得了研究區(qū)內(nèi)高精度的縱橫波速度結(jié)構(gòu)和泊松比異常分布形態(tài).研究結(jié)果揭示了研究區(qū)內(nèi)的地殼結(jié)構(gòu)具有強(qiáng)烈的橫向不均一性.2020年濟(jì)南長(zhǎng)清MS4.1地震可能是區(qū)域構(gòu)造應(yīng)力下長(zhǎng)清斷裂發(fā)生左旋走滑運(yùn)動(dòng)的結(jié)果，震中東側(cè)的高速異常體可能阻礙了地震的東向破裂.2003年青島嶗山ML4.1地震可能是由深部的流體填充相對(duì)完整巖體內(nèi)的裂隙并引發(fā)破裂所致.嶗山震群、乳山震群和長(zhǎng)島震群可能與深部流體有非常強(qiáng)的相關(guān)性.