基于Sentinel-1A的2021年雙湖MW5.7地震同震破裂模型

夏 濤 劉小利 余鵬飛 鄧德貝爾 樂子揚 高天琪

1 中國地震局地震大地測量重點實驗室，武漢市洪山側路40號，430071 2 防災科技學院，河北省三河市學院街465號，065201

2021-03-30西藏雙湖縣黃水湖西岸發生MW5.7地震，中國地震臺網中心(CENC)給出震中位置為34.38°N、87.68°E，震源深度為10 km。雙湖縣位于青藏高原腹地的羌塘塊體中部，自1970年有地震記錄以來，羌塘塊體中部共發生過14次5級以上地震，本次地震是周邊100 km范圍內的首個5級以上地震(圖1(a))。由于自然環境惡劣、交通限制等原因，羌塘塊體中部測震臺和GNSS連續站非常稀疏，地質學資料有限，距離震中200 km范圍內除震中西側約1 km處有一條近SN向、傾角不明的無名正斷層外(本文稱黃水湖斷層)，尚未有其他已知斷層[1-2]。受限于震中區域地質構造、余震和同震GPS等資料較少，發震斷層的確定較難，CENC未給出此次地震的震源機制解，而GCMT和USGS給出的震源機制解差異較大，發震斷層面選擇不統一。

斷裂數據來源于文獻[1-2]；震源機制解源于哈佛大學GCMT；歷史地震源于CENC；GPS速度場源于文獻[3]圖1 羌塘塊體中部歷史地震及斷裂分布Fig.1 The distribution of the faults and historical earthquakes of the central Qiangtang block

InSAR可大范圍、高精度地獲取地震同震形變場，是確定中強地震發震構造的重要手段之一。本文利用Sentinel-1A SAR-C地震前后升軌和降軌4景數據獲取2021年雙湖地震同震形變場，基于半空間位錯模型確定發震構造幾何學參數，計算不同深度上的同震靜態庫侖應力變化，并結合羌塘塊體構造背景討論本次地震的發震機制及其對區域未來地震危險性的影響，為深入理解羌塘塊體中部構造活動性和地殼變形機制提供參考依據。

1 區域構造背景

GPS速度場(圖1(b))顯示，青藏高原內部各次級塊體呈現不同的地殼形變機制[4-5]。班公湖-怒江縫合帶以北的羌塘塊體從西向東滑移速率逐漸減小，羌塘塊體中部至東部滑移速率相差可達10.3 mm/a[6]，表現出明顯的空間差異性；班公湖-怒江縫合帶以西在上地殼物質持續向東擠出的構造作用下，沿羌塘內部小型正斷層和共軛走滑斷裂發生彌散型變形特征，構造活動以EW向延伸擴展為主；班公湖-怒江縫合帶以東則直接向東流出，使得塊體東邊界不斷向東擴張，發生剛性塊體變形，呈現明顯的物質快速轉移運動特征[4]。

目前，對羌塘塊體的斷裂活動性仍缺乏定量研究，特別是地處青藏高原腹地的塊體中部，僅Wang等[7]利用時序InSAR技術發現羌塘塊體中部存在持續向東的滑移運動，速率約為4 mm/a。2021年雙湖發震發生于羌塘塊體中部的黃水湖西岸，其發震構造和變形機制不明，變形特征是呈彌散型還是剛性整體向東擠出均有待進一步研究。

2 InSAR數據處理

本文利用Sentinel-1A影像(表1)，基于兩軌差分干涉獲取本次雙湖地震的干涉影像，提取相應的同震形變場?；陂_源雷達數據處理平臺GMTSAR[8]進行以下數據處理：1)采用20∶4的多視比例；2)通過GACOS移除影像中的大氣誤差[9]；3)根據精密軌道和SRTM提供的30 m DEM去除軌道誤差和平地相位；4)采用Goldstein濾波方法和最小費用流解纏方法分別對干涉圖進行濾波和相位解纏處理；5)地理編碼。

表1 干涉對信息Tab.1 Interferometry pairs information

主震升降軌干涉圖和同震形變場如圖2所示。可以看出，形變場以近SN向分界呈三瓣式(L1、L2、L3)分布。L1部分在升降軌形變場中最大視線向(line of sight，LOS)形變量分別為0.53 cm、-0.46 cm，二者符號相反、量值接近；升軌時遠離、降軌時靠近L2部分，指示L1存在NWW向水平運動。L2部分在升降軌形變場中皆呈現負值，LOS向形變量分別為-6.78 cm、-4.73 cm，表明該區域整體為垂直沉降運動。L3部分在升降軌形變場中LOS形變量分別為-1.65 cm、1.45 cm，符號相反、量值接近；升軌時接近、降軌時遠離L2部分，指示L3存在近東向水平運動。圖2(g)顯示，沿AB測線得到的剖面T114A在從L1向L2過渡時(L1與L2之間的虛線)，其梯度明顯大于剖面T121D，可能是由于T121D的成像視角與干涉圖長軸方向近乎平行，導致真實地表形變在形變場中的投影分量較少，限制了LOS形變的觀測精度。表2中最佳均勻滑動模型的傾角為本次雙湖事件形變場長軸為近SN走向，主要以垂直向形變為主。因InSAR對SN方向形變分量不敏感，因此在忽略SN向形變分量的情況下，利用不同視角的升降軌同震形變場提取震區EW水平位移場和垂直向位移場[10](圖2(c)、2(f))更適用于本次地震。計算得出，EW向和垂直向形變場的殘差(RMS)分別為0.4 cm和1.3 cm，在InSAR精度范圍內，說明2.5D形變場具有一定的可靠性。從圖2(c)水平向形變場來看，L1表現為西向水平位移，最大位移值為2.5 cm；L2、L3則為東向水平位移，最大位移值為5.8 cm，且出現在L2區間；從圖2(f)垂直向形變場來看，最大沉降量達11.8 cm，出現在L2區間；L1、L3部分呈微弱抬升，最大抬升量為2 cm。此外，對比3組形變場來看，在L1與L2的過渡區域，與T121D形變場相比，T114A表現出更為明顯的近SN向長軸走向。由此推斷，造成T114A和T121D沿測線形變量變化梯度差距較大(圖2(g))的可能原因是T121D成像角度與形變場近SN向長軸幾乎垂直，對形變場水平向的位移并不敏感，且L1部分本身水平位移量較小。

暖色代表接近衛星方向或抬升運動，冷色代表遠離衛星方向或沉降運動圖2 雙湖地震同震形變場Fig.2 Coseismic deformation field of the Shuanghu earthquake

3 斷層滑動分布反演

為精細地確定發震斷層位錯滑移特征，采用以下步驟進行計算：1)以同震形變場為輸入數據，利用蒙特卡洛算法搜索最佳斷層參數，獲得發震斷層的均勻滑動模型；2)以最佳均勻滑動模型參數為先驗條件約束，利用Okada位錯模型反演同震滑動分布。

3.1 均勻滑動反演

為提高最佳斷層參數的搜索效率，確保有效收斂，使用四叉樹采樣算法[11]對升降軌同震形變場進行降采樣處理。處理后得到494個升軌、422個降軌形變點位，以此為蒙特卡洛搜索約束[12]，獲得收斂性較好的馬爾科夫鏈(圖3)，并統計特征最佳的均勻滑動模型斷層參數(表2)。

紅線表示最佳斷層模型參數圖3 蒙特卡洛-馬爾科夫鏈聯合概率密度分布Fig.3 Distribution of joint probability density of Monte Carlo-Markov chain

表2 均勻滑動模型Tab.2 Uniform slip model

表2中最佳均勻滑動模型的傾角為51°，走向為36.6°，走向與主震干涉圖NNE向長軸走向一致；斷層走滑量為0.39 m，傾滑量為-0.61 m，指示雙湖地震可能是一次正斷為主、兼右旋走滑的地震，但已有研究并未表明該區域存在右旋走滑斷層[13]。因此，有必要構建滑動分布模型，進一步正演精細的三維形變特征以確定發震斷層幾何學參數。

3.2 斷層滑動分布反演

以均勻滑動模型得到的最佳斷層模型參數為先驗條件，基于Okada位錯模型[14]進行同震斷層滑動分布反演[15]，使用Crust1.0模型，反演約束方程為：

式中，G為聯系觀測值與彈性位錯模型的格林函數矩陣，d為觀測值，?2為拉普拉斯算子，γ2為平滑因子。

根據均勻滑動模型，設定以下初始參數：斷層長18 km、寬20 km、傾角50°、走向30°、滑動角范圍-130°～-50°。將斷層按1 km×1 km劃分為360個子斷層，并且設定平滑因子γ2為0.07。最終得到的斷層滑動模型為：震中34.37°N、87.71°E，震源深度6.51 km，走向33.0°、傾角50°、平均滑動角-74°、最大滑動量0.26 m(圖4)。從圖4(a)～4(c)可見，破裂集中在沿走向8～16 km(滑動量大于10 cm)、地下深度5～15 km處，斷層滑動以傾滑為主，兼有微弱左旋分量。假定剪切模量為30 GPa，計算得到地震矩張量為4.5×1017Nm，相當于震級MW5.74。

圖(c)中，視線方位角為150°，高度角為20°，加粗實線為斷層上盤；圖(d)中，箭頭表示水平運動方向圖4 斷層滑動分布模型及模擬三維形變場Fig.4 Fault slip distribution model and the simulated 3D deformation field

圖4(d)是基于滑動分布模型重建的矢量形變場。水平向位移顯示，斷層西側向西移動，東側則微弱東移，呈左旋走滑特征；垂直向形變場存在明顯的NNE向分界線，以沉降運動為主。顯然，模擬的矢量形變場整體形變特征與2.5D形變場基本一致，且更加精細地展示了斷層在不同維度上的滑動量值。

圖5是基于均勻滑動模型、滑動分布模型得到的正演模型及其擬合殘差?？梢钥吹剑琓114A的2個模型與觀測值擬合度更好，RMS均為0.3 cm，殘差主要出現在形變場沉降區，主要是因黃水湖干擾而缺乏部分近場觀測信號導致的，滑動分布模型以清晰的NNE向分界線指示了發震構造行跡。T121D的均勻滑動模型和滑動分布模型擬合的RMS分別為0.3 cm和0.6 cm，但均勻滑動模型的指示意義相對不佳。聯合T114A得到的斷層走向參數約束后，可以改善T121D滑動分布模型在斷層走向的擬合效果(圖5(g)、5(i))，進一步表明T121D近乎平行于斷層走向的飛行方向導致該方向上的形變獲取受限，進而使得在該方向上的形變量低于T114A。

黑色矩形框表示滑動分布模型斷層地表投影，加粗實線表示斷層上盤，(a)～(e)和(f)～(j)分別對應T114A和T121D圖5 同震形變場觀測值、模型值與殘差Fig.5 The observation values, model values and residual of coseismic deformation field

4 討 論

均勻滑動模型和滑動分布模型獲得的震源機制參數見表3?？梢钥闯觯琁nSAR反演的模型反映的是地震破裂面的總體特征，與USGS、GCMT給出的遠震體波模型結果基本一致，僅震源深度相差較大，這是由于遠震體波模型地震矩張量反演中波形擬合誤差隨深度變化，深度結果誤差相對較大[16-17]。

表3 2021年雙湖MW5.7地震震源參數Tab.3 Focal parameters of the 2021 Shuanghu MW5.7 earthquake

測震學和大地測量學結果都顯示，2021年雙湖MW5.7地震為一次正斷為主、兼具走滑分量的事件。但受特定方位的遠震地震資料信噪比低等影響，USGS和GCMT給出的節面I的滑動角以及InSAR均勻滑動模型和滑動分布模型的滑動角分別表現出右旋或左旋走滑，干擾了走滑方向的確定。已有地質學資料顯示，羌塘塊體中部已知斷層未見右旋走滑性質[13]。由于震中區域缺少GPS連續觀測站，測震臺也非常稀疏(圖1)，主震后2個月內僅發生8次余震(表4，數據來源于CENC地震目錄)，無法提供更多資料約束發震節面。

表4 2021年雙湖MW5.7主震及余震信息Tab.4 2021 Shuanghu MW5.7 main shock and aftershocks information

為了解決上述問題，采用靜態庫侖應力轉移與同震形變場、余震的對應關系來討論不同性質發震節面的差異，進而判定發震構造性質，評估本次地震對區域未來地震危險性的影響?？紤]到InSAR滑動分布模型給定的左旋節面解與GCMT節面I參數更接近，為對比分析，采用Crust1.0模型、摩擦系數0.4、楊氏模量30 GPa、泊松比0.25，分別以InSAR左旋節面解(走向33°、傾角50°，滑動角-74°)和GCMT右旋節面解(走向179°、傾角54°，滑動角-112°)為接收面計算不同深度上的靜態庫侖應力，結果分別見圖6(a)～6(g)、6(h)～6(n)。

從圖6(a)～6(g)看出，5 km深度處，庫侖應力呈現明顯的對稱特征，應力卸載區呈EW向，應力加載區呈SN向，東側卸載區與南側加載區分界明顯，呈NNE向，與左旋節面走向一致，左旋節面給定的主震即落于該分界處；最大應力增加為0.6 bar，超過應力觸發閾值0.1 bar，存在較大地震風險。7 km深度處，應力卸載區主要沿左旋節面走向進一步擴張，與左旋節面給定的主震位置和SAR觀測到的主要形變區范圍(圖5)基本一致，同時南側應力加載明顯。考慮到InSAR左旋節面給定的主震深度為6.51 km，且在7 km處釋放了28.9%的力矩，推斷該深度為主震應力釋放的主要深度。10 km深度處，SN向應力加載區范圍增大，EW向應力卸載區范圍減小，且呈東西兩瓣沿左旋節面方向反向對稱；最大應力增加達0.9 bar，遠大于應力觸發閾值0.1 bar，發生于2021-04-02、震源深度為9 km的4級余震即落于東側卸載區與北側加載區的分界處；南側應力加載也較明顯，仍存在較大地震風險。15 km深度處，SN向應力加載區進一步擴張并連通，主震所在區域應力持續卸載。在20～30 km深度處仍以應力加載為主，最大應力在20 km深度時增加到0.4 bar，而其他7次余震都集中發生于該深度區間，且均落于庫侖應力大于0.1 bar的應力加載區或邊界。

(a)～(g)為InSAR結果；(h)～(n)為GCMT結果；綠色、黃色小球分別表示不同震級的余震分布,序號與表4對應圖6 不同接收面對應的不同深度上同震靜態庫侖應力變化與主余震分布Fig.6 Coseismic static Coulomb stress changes at different depths corresponding to different receivers and distribution of main shock and aftershocks

從圖6(h)～6(n)看出，應力變化分布與左旋節面解結果存在較大差異，主要表現在4個方面: 1) 5～10 km深度，應力卸載區呈NW向擴張，與右旋節面走向不符；2) 最大應力卸載區與SAR觀測到的主要形變區范圍(圖5)相差較大，且出現在7 km深度處，與表3中GCMT模型給出的主震深度均不一致；3) GCMT給出的主震深度為16.4 km，但在15 km深度處未出現明顯的應力卸載區，相反在主震東側應力加載明顯，而所有余震均未出現在該區域范圍；4)除2021-04-02的4級余震震源深度為10 km左右，其余7次余震震源深度都集中在25～32 km，但在25 km、30 km 深度處余震所在位置絕大部分對應應力卸載區。

綜上分析，由左旋節面作為接收面獲得的靜態庫侖應力變化與SAR觀測得到的形變場、主余震分布更加契合，也符合羌塘塊體區域已有應力場及地質學研究結果[6,13]。可以認為，2021年雙湖MW5.7地震的發震斷層應為正斷兼左旋走滑性質。從圖5、圖6(a)～6(g)看出，主震和2021-04-02的4級余震發生于10 km深度以上的淺部地殼，且分布于黃水湖斷層南段；其他余震均發生于25 km以下深度，且大都散落于黃水湖斷層兩側，集中在斷層面地表投影范圍內。分析圖6中不同深度上庫侖應力變化及余震分布，認為黃水湖斷層南段短時間內地震風險性較小，北段地震風險性需結合更多資料進一步分析。

本次地震發生于羌塘塊體中部的NNE向地塹區，InSAR形變場特征反映這是一次典型的地塹擴張活動。從區域狀態來看[6]，羌塘塊體整體張、壓、剪應變都較強烈，為擴張型地塊。從西至東，塊體SN向擴張率逐漸增加，最大達24.40×10-9/a；而EW向擴張率逐漸減小，低至3.21×10-9/a；塊體中部的主壓、主張、最大剪切和面應變率分別為14.29×10-9/a、25.13×10-9/a、35.41×10-9/a和6.84×10-9/a，表明羌塘塊體持續向東延展。Yu等[18]指出，羌塘塊體東邊界可能已向東延伸至95°E，而不是Taylor等[2]認為的93°E。因此，本次地震進一步指示，羌塘塊體仍持續向東擴張，因EW向拉伸作用導致近SN向小規模正斷層呈彌散型變形和地塹發育。

5 結 語

1)本次雙湖地震震中位置為34.37°N、87.71°E，震源深度6.51 km，發震斷層走向33°、傾角50°、傾向東、滑動角-74°，以傾滑為主兼少量左旋走滑分量，最大滑動量0.26 m；

2)短時間內，黃水湖斷層南段地震風險性較小，北段則需結合更多資料進一步分析；

3)本次地震是在羌塘塊體持續向東擴張的背景下，受EW向拉伸作用使得黃水湖正斷層發生的一次彌散型變形活動，SN向地塹得到進一步擴張。