川滇地區深部結構的三維地震層析成像研究

王海燕

摘 要：本文利用二維介質中的射線追蹤完成在地球內部地震波的三維的射線追蹤，根據走時殘差，建立地震波旅行時層析成像系統方程組。通過川滇地區寬頻帶地震儀所記錄的遠震P波走時殘差，獲得該區域的三維P波速度結構，同時結合當地地質構造和已有的層析成像結果，對本論文的成像結果進行分析，并對青藏高原東緣地區深部速度結構和殼慢動力學特點進行討論與研究。

關鍵詞：三維地震層析成像；射線追蹤；地殼結構

1 前言

地球內部蘊含著豐富的資源的同時其內部的構造活動又給人類帶來巨大的災難，而利用天然地震所產生的地震波信息，已經獲得了許多對地球科學產生深遠影響的發現，比如地球內部成層結構的建立、深源地震的發現等[1-5]。由于地震波的速度取決于介質的密度、物質組成、彈性模量、溫度等因素，因而地震波場所攜帶信息的差異反應了地球內部介質特性的變化，地震層析成像技術就是利用地震波的到時、波形或頻散特征等，建立地球內部的三維速度結構影像。但是由于研究對象的特殊性，如震源定位不是很準確，地震波射線的透射角度有限，在非均勻介質中地震波傳播的復雜性等種種原因，致使地震層析成像技術很難達到醫學CT的效果。

2 IASP91地球模型中的射線追蹤

如果我們把橢圓的地球近似看做球形，且內部成層狀構造，可以認為地球是由具有對稱性的無限多個均勻的同心圈層構成，地震波在地球介質內部傳播，因而，可以利用二維介質中的射線追蹤完成地震波在地球內部的射線追蹤，包括地震波在圈層介質分界面上的反射波等震相的追蹤。

根據IASP91地球模型的速度參數表可知，其內部速度是地球半徑的函數，并且每個地層的速度梯度不等，為保證地層速度變化的連續性，網格結點之間的速度值利用三次樣條插值獲得，進而對介質的縱向非均勻特征進行描述。

針對上述地球模型，我們可以將地震波在地球內部傳播的縱波及橫波速度分別表示為

（1）

在二維介質的情況下，地震波傳播的射線軌跡滿足程函方程

（2）

可以推導出二維情形下射線軌跡滿足的一階微分方程組。

在數值計算中，利用數值微分可以寫為

（3）

其中（x0，z0）為射線起始點E的坐標（圖1），為射線的離源角，及為速度場在x、z兩個坐標軸方向上的速度梯度，為給定的時間步長。射線經過時間后，射線將到達A點，其空間坐標為（x，z），離源角為。當A點位于某一界面時，需要根據斯奈爾定律來計算經過折射后的離源角。

圖1 射線路徑示意圖

（1）當射線由震源E以離源角α0出發，按規定的步長dx逐步前進，到達A點，計算A點的坐標，以及走時TA，并判斷點A是否到達或者超過第一層下界面。

（2）若A點沒有達到下一層界面，則根據點所在地層的層序數，點到地心的距離R1，獲取射線當前所在位置的速度值以及震相。

（3）根據層序數判斷是射線是否到達地層的最底層，并對地層或塊體內部進行追蹤，確定射線傳播的路徑及旅行時間，并返回射線終點的坐標及離源角α等參數。

（4）從A點出發，繼續以步長dx前進，到達B點，計算B點的坐標以及走時TB，同樣進行判斷B點所處位置，若依然沒有到達或者超過界面，則重復以上計算，直到射線到達和超過界面為止；若射線超過界面，則采用“二分法”求射線與界面的交點，以及交點處界面傾角，并判斷是否進行廻折波射線追蹤。

（5）如果到達的不是目的界面，則要計算出入射線與法線的夾角、折射線與橫軸的夾角等相關參數；射線到達目的界面后，則進行反射波上行射線追蹤，直到追蹤到地面為止。

（6）改變震源處的離源角，重復上述步驟就得到了一條條經過地下介質到達地面的射線。

3 層析成像的基本原理

地震波旅行時層析成像是典型的地球物理反演問題，依據地震波穿過介質內部所需的時間來重建介質的速度結構。

通常情況下，射線的走時t 可以寫成：

（4）

其中，i是射線條數，l是射線元，v、s為分別介質的速度、慢度函數。

如果假設介質的慢度分布函數為，其中為介質的慢度初始模型，為介質的慢度擾動量，那么，根據（4）式有

（5）

從而有

（6）

所以，由上式可知，如果知道地震波的走時擾動量（走時殘差），就可以研究介質中慢度（速度）的擾動分布特征。將積分離散化則層析成像系統的方程組具有如下形式

（7）

寫成矩陣，有

D=CU （8）

其中為走時殘差，由實際觀測值得到；為介質分塊中的慢度分布，為待求向量；C為不同單元格內所有線段的長度所形成的一個大型稀疏系數矩陣，由射線的傳播路徑確定。

（9）

4 地震層析成像的應用

青藏高原的形成原因一直備受國內外地學家的關注，尤其是青藏高原巖石圈的深部結構。2003-2004年期間，美國麻省理工大學與國土資源部成都地質礦產研究所合作，在三江地區進行了寬頻帶地震觀測，本文利用遠震P波走時殘差，獲得該區域的三維P波速度結構。

我們研究的區域范圍為24.0°N～31.0°N，99.5°E～104.0°E，本文使用了25臺寬頻帶地震儀所記錄的地震數據，從中挑選出了信噪比較高的714個地震事件，震級大于4.5級，震中距為25.0°～180.0°，震源深度為0～641km，并且每一個遠震事件至少有10個臺站具有清晰的地震記錄，每個記錄的走時殘差絕對值小于5s。本文按照上文闡述的射線追蹤以及層析成像方法對研究區域寬頻帶地震數據記錄進行地震波旅行時層析成像。

根據眾多學者對三江地區的研究資料 [5-10] 可知，該區域地殼比較厚，莫霍面的深度達到了70km左右，并且該地殼有上、中、下三層結構，其厚度分別為10km、25km、35km，所以本文在IASP91地球模型基礎上稍微做了修改，修改后的IASP91地球模型參數我們稱為初始模型。網格在經度方向和緯度方向上均以1.0°間距作水平方向的劃分，在深度方向上以10、20、30、40、50、65、80、100、120、150、200、250km進行劃分。

我們對三江地區觀測數據進行層析成像，則由0～125 km深度處的速度擾動水平剖面圖中可以看出，川西高原地殼和上地幔的速度存在明顯的橫向不均勻性，川中塊體的速度相對較高，而川青塊體以及川滇菱形塊體的速度相對較低，并且川滇菱形塊體相對較低的速度區域呈南北向分布；大型斷裂帶兩側存在明顯的速度差異，如龍門山斷裂帶、安寧河斷裂帶、金沙江斷裂帶等為研究區域中幾大塊體的速度分界線。

縱觀圖2的4個深度剖面，我們發現，在金沙江斷裂帶和安寧河斷裂帶之間的研究區域，隨著深度的增加，特別是深度150～250km區間，低速度區域的面積逐漸增大。在川滇菱形塊體東側的低速區域，已穿過其速度分界線龍門山斷裂、安寧河斷裂帶延伸至四川盆地內部。根據曾融生、Blackman、Sliver 等[11-14]提出青藏高原物質有橫向流動的可能性，推斷出現這一現象的原因可能為：印度板塊向青藏高原的下部俯沖，使青藏高原的地殼抬升，在巨大的擠壓下，增多的地殼物質向青藏高原東部地區的上地幔流動。其中東流的物質在龍門山斷裂帶附近遇到了四川盆地的阻擋，部分改向東南方向流動，但是仍有部分向下侵入到四川盆地的西南區域，致使區域 P 波速度較低。并且在深部斷裂帶呈現負的速度異常，更有助于地殼塊體沿斷裂的側向擠出。

圖3 P 波震相速度擾動水平剖面

從102°E經度P波震相的速度擾動垂直剖面（圖4）可以看出，該剖面自南往北依次穿過滇西塊體、滇中塊體、川西塊體、川東南塊體、川青塊體。在0～150 km深度上，23°N～31°N范圍內的地殼和上地幔速度存在明顯的不均勻性，低速和高速相間。產生這一現象的原因可能是印度板塊在與歐亞板塊碰撞、俯沖、擠壓過程中，由于地殼應力場的非均勻性和鮮水河斷裂帶的存在，在高速的四川盆地邊緣發生了斜上方的逆沖運動造成的。在深度150km以下，該剖面的速度值逐漸變小，且低速度區域的面積逐漸增大。這一結果也在一定程度上驗證了川西高原在下地殼至上地幔的范圍內存在較軟物質的可能性。

圖4 102°E經度P 波震相速度擾動垂直剖面

圖5 29°N緯度P 波震相速度擾動垂直剖面

從29°N緯度P波震相的速度擾動垂直剖面（圖5）可以看出，在0～100 km深度范圍內，地殼和上地幔的速度存在橫向的不均勻性。在川滇菱形塊體內101°E～102°E范圍的80km深度以上，該處存在一個高速區域。100km以下，川滇菱形塊體低速度區域的面積逐漸增大。并且低速體在西部以金沙江斷裂為界，東部以小江斷裂為界，據此我們推測，在青藏高原受到印度板塊NNE方向的擠壓后，川滇地區受到連帶影響，由于斷裂帶的巖石張力的破碎，使四川盆地沿著小江斷裂逆向俯沖（圖5箭頭所示）。

綜上所述，我們推測：印度板塊以NNE方向運動與歐亞板塊碰撞，并向青藏高原下部俯沖，致使青藏高原地殼增厚，并向東強烈擠壓，由于地殼應力場的非均勻性和斷裂帶的存在，在這種擠壓下，川滇地區的地殼發生變形，區內的次級塊體及邊界產生逆沖或者扭轉，對川滇地區的地質構造造成影響。由研究區域大面積的低速體的存在，推測該現象可能是由地殼和上地幔存在高導層、高熱流值造成的。

縱觀整個研究區域，并結合前人的研究成果，我們得出川滇地區地殼結構的總體特征是：在上地殼速度異常分布中，存在明顯的橫向不均勻性，但是總體上川滇地區地殼、上地幔的速度均速度較低；龍門山斷裂帶、安寧河斷裂帶，金沙江斷裂帶在地殼和上地幔一定深度內的速度異常中顯示出構造分界特征，并且深部斷裂帶呈現負速度異常；地殼厚度變化劇烈，地殼和上地幔存在高導層、高熱流值。

5 存在問題

地震波旅行時層析成像包括數據處理，正演模擬和數據反演等環節，而每個環節都有可能導致最終結果的不確定。該方法在應用中主要存在問題具體如下：

（1）反演的數據量不足。由于臺站的數量有限，參與層析成像反演的數據量不足，從而導致成像結果的分辨率不高。

（2）震源定位不準確。由于臺網的稀疏，使得地震的定位存在很大的誤差，尤其是震源的深度，那么對地震層析成像反演結果有很大影響的的旅行時就會變的不準確

（3）初始參考模型引起的誤差。在進行射線追蹤時，由于不知道真實的射線路徑，我們首先引入了參考模型，然后通過逐步迭代對參數進行修正。如果成像反演時利用了折射或者反射震相，則初始參考模型中速度界面信息的準確性尤為重要，因為速度間斷面對這些震相的射線傳播路徑影響很大，所以可靠的速度間斷面信息對地震層析成像反演是至關重要的。

參考文獻

[1] Wang Y X， Jiang M， Nábelek J，et al. The new insight of the velocity structure beneath the Indian-Eurasian continental collision zone. 2006， AGU West. Pac. Geophys. Meet， Suppl.， Abstract S45A～04.

[2] Zhao Dapeng， Akira Hasegawa. P wave tomographic imaging of the crust and upper mantle

beneath the Japan islands [J]. J. Geophys. Res.， 1993， 98（B3）： 4333～4353

[3] 賀世杰，郭峰. 地幔柱的識別和演化綜合述評 [J]. 地球科學進展，2003，18（3）：433～439

[4] 劉建華，劉福田，吳華，等. 中國南北帶地殼和上地幔的三維速度圖像 [J]. 地球物理學報，1989，32（2）：143～151

[5] 王有學，姜枚，熊盛青，等. 西昆侖巖石圈的拆沉作用及其深部構造含義—地震層

析成像及航磁異常證據 [J]. 中國地質，2006，33（2）：299～308

[6] 姚振興，李白基，梁尚鴻，等. 青藏高原地區瑞利波群速度和地殼構造[J]. 地球物理學報，1981，24（3）：287～29

[7] 丁志峰，何正勤，吳建平，等. 青藏高原地震波三維速度結構的研究[J]. 中國地震， 2001，17（2）：202～209

[8] 錢輝，姜枚， Chen Wangping，等. 青藏高原吉隆—魯谷（Hi-Climb）層析成像與引藏碰撞的消減作用[J]. 地球物理學報，2007，50（5）：1427～1436

[9] 吳慶舉，曹融生，趙文津. 喜馬拉雅-青藏高原的上地幔傾斜構造與陸-陸碰撞過程[J]. 中國科學D輯，2004，34（10）：919～925

[10] 賀日政，趙大鵬，高銳，等. 西昆侖造山帶下巖石圈地幔速度結構[J]. 地球物理學報，2006，49（3）：778～787

[11] 曾融生，朱介壽，周兵，等. 青藏高原及其東部鄰區的三維地震波速度結構與大陸碰撞模型 [J]. 地震學報，1992，11（S）：520～533

[12] 曾融生，丁志峰， 吳慶舉. 青藏高原巖石圈構造及動力學過程研究[J]. 地球物理學報， 1994，37（S）：99～116

[13] Blackman D K. Seicmic anisotrophy in the mantel beneath an oceanic spreading center[J]. Nature，1993， 366： 675～677

[14] Sliver P G， Chan W W. Implication for continental structure and evolution from seismic anisotropy [J]. Nature， 1998， 335（1）： 34～39