李忠亞 胡敏章 郝洪濤 李 輝
1 中國地震局地震研究所,武漢市洪山側路40號,430071 2 中國地震局地震大地測量重點實驗室,武漢市洪山側路40號,430071 3 湖北省地震局,武漢市洪山側路48號,430071
根據中國地震臺網測定,2017-08-08四川省西北部九寨溝縣發生MS7.0地震,震中位于33.2°N、103.82°E,震源深度約20 km。中國地震局地球物理研究所、美國地質調查局和相關文獻計算結果均表明,本次九寨溝地震為走滑型地震[1-2]。九寨溝地震周邊地形和構造復雜(圖1),近10 a來,其南部和北部先后發生2008年汶川MS8.0地震、2013年蘆山MS7.0地震、2013年岷縣-漳縣MS6.6地震和2016年門源MS6.4地震等破壞性強震。本次地震的發震斷層為樹正斷裂[1,3],其位于塔藏斷裂、岷江斷裂和虎牙斷裂西北延長線的交會區域,也是東昆侖斷裂帶由走滑到虎牙斷裂逆沖構造的過渡區域[4]。
地震孕震過程中的質量遷移必然伴隨著重力場的改變,震前重力場變化可以由GRACE重力衛星和地表重力觀測得到[5-6]。通過絕對重力儀和相對重力儀組合的方式可以進行流動重力觀測,不同期觀測得到的重力變化特征可以作為強震發生前的可靠前兆[7]。目前利用流動重力資料判斷地震發生的可能地點和時間主要依據重力場時空動態演化特征和相關經驗信息。識別重力場動態演化規律時,目標區域經常會出現多個“正負重力變化過渡帶”等震前重力場變化特征,但受限于目前對地震孕育機制的認識和觀測點位分布等因素,多數情況下難以對其進行取舍。本文嘗試通過將重力變化數據進行球諧展開并濾波,壓制局部重力場變化特征,突出主要演化規律,以便于目標特征的識別。

F1:東昆侖斷裂帶;F2:阿萬倉斷裂;F3:塔藏斷裂;F4:光蓋山-迭山斷裂;F5:西秦嶺斷裂帶;F6:岷江斷裂;F7:虎牙斷裂;F8:龍門山斷裂帶;F9:鮮水河斷裂帶;F10:龍日壩斷裂圖1 青藏高原東緣地形、主要斷裂帶和城市分布Fig.1 The distribution of topography, main fault zones and cities in eastern Tibetan plateau
空間任意一點(θ,λ)的引力位球諧展開表達式[8]為:
(1)

(2)
重力儀在地球表面觀測到的重力變化δg為:
δg=Δg+βur
(3)
式中,Δg為空間固定點的重力變化,β為布格梯度,ur為徑向位移。采用球模型近似,空間固定點重力變化Δg和引力位變化ΔV有如下關系:
(4)
當Δg已知時,可以直接進行球諧展開:
(5)
式(4)和式(5)的球諧系數滿足:
(6)
流動重力觀測結果可以很好地反映地震孕育過程中產生的重力變化[9-11]。本文重力變化數據來自中國地震局重力學科管理部,時間變化范圍為2015-04~2017-04。每期重力數據處理過程是采用絕對重力結果作為起算基準,然后將絕對重力和相對重力進行聯合平差,計算過程中采用絕對重力基準對相對重力儀格值進行解算以削弱格值誤差。不同期數據之間作差即可得到相應時間段內的原始重力變化。
為了獲取研究區域重力場動態變化圖像,需對重力變化數據進行后處理,流程如下:1)剔除極少數誤差較大的數據,最終使用的數據點位平均精度優于15 μGal。2)由于提供的重力變化數據是地表值,而根據上節方法計算處理時需要空間固定點重力變化值,為此進行布格梯度歸算。布格梯度改正包含垂直形變效應和布格層效應。垂直梯度可取 -3.086 μGal/cm,地殼密度采用平均值2.67 g/cm3,由布格改正公式得到布格改正為1.1 μGal/cm。布格梯度是垂直改正和布格改正之和,大小為 -2.0 μGal/cm。研究區域每年隆升速率變化范圍在1~6 mm/a[12],以最大隆升速率6 mm/a計算,每年因為垂直形變引起的重力變化絕對值不超過1.2 μGal,該數值小于相對重力儀觀測到的重力變化,因此本文不區分重力儀觀測的地表重力變化和空間固定點重力變化。3)實際數據處理時,選擇整個南北地震帶重力觀測數據(重力觀測點位分布見文獻[13])。將區域內重力變化數據進行格網化,區域外重力變化值設為0,根據式(5)進行球諧分析,得到60階球諧系數。利用球諧系數進行求和并作3°范圍平滑處理,得到的重力變化見圖2。
圖2顯示,九寨溝地震前2 a青藏高原東緣重力場動態變化劇烈。下面具體分析重力場變化特征:

圖2 青藏高原東緣重力場動態演化Fig.2 Gravity field changes in eastern Tibetan plateau
1)2015-04~2015-10期間,全區域絕大部分地區為正重力變化。在九寨溝、黑水、汶川、丹巴和雅安一帶是重力增加峰值區域,最大正重力變化達到20 μGal。除莊浪東北地區重力變化接近0外,其余地區重力變化幅度沿著峰值區域向西北和東南兩側遞減。峰值區域走向與龍門山斷裂帶走向基本一致,說明重力場變化與構造運動密切相關。
2)2015-10~2016-04期間,大部分地區重力變化趨勢與2015-04~2015-10相反,表現為負重力變化。巴顏喀拉塊體東部和柴達木塊體的貴南至澤庫一帶是主要的重力減小地區,周邊的川滇塊體東北部和華南地塊西北部重力變化平緩,但在漳縣和西和東北地區出現正重力變化。區域整體重力場變化特征表現為,巴顏喀拉塊體東緣沿東北方向至莊浪附近,負重力變化逐漸過渡為正重力變化,形成負-正重力變化過渡帶。
3) 2016-04~2016-10期間,重力場變化最明顯的特點是巴顏喀拉塊體負重力變化在2015-10~2016-04期間變化基礎上繼續減小,并且負重力變化分布向西南延伸至川滇塊體東北部。阿壩、馬爾康至新龍一帶是負重力變化的峰值區域,最大負重力變化達到-35 μGal。在柴達木塊體和華南塊體內部,重力變化分布特征分別是向北側和東側逐漸由負重力變化向正重力變化過渡,在貴南附近、資陽和遂寧東南區域出現重力正變化特征。區域整體重力場變化特征為川滇塊體東北部負重力變化向東部、東南部和北部逐漸過渡為正重力變化,形成重力變化過渡帶。
4)2016-10~2017-04期間,重力場變化在前2期基礎上出現反轉回調,具體表現為2個特征:一是負重力變化,在阿壩、若爾蓋、岷縣和漳縣西北地區,重力變化由0逐漸減小,在貴南附近數值達到 -20 μGal;二是正重力變化,在汶川、綿陽、平武和廣元一帶出現正重力變化峰值區域,最大重力變化為18 μGal,并且正重力峰值區域走向與龍門山斷裂帶走向一致。在重力變化峰值區域兩側,即西北方向至阿壩、若爾蓋、岷縣和漳縣一帶,東南方向至遂寧和資陽一帶,正重力變化幅度逐漸減小。區域重力場變化特征是貴南附近的負重力變化沿東南方向逐漸過渡至正重力變化,形成重力變化過渡帶。
地震孕育過程中伴隨著質量遷移和能量交換,可通過分析震前重力場動態演化特征來捕捉地震前兆信息。分析圖2可知,九寨溝地震前,震源周邊重力場演化可以分為3個階段。第1階段為重力增大階段,2015-04~2015-10巴顏喀拉塊體東部邊界重力場均表現為正變化,沿著邊界向四周正變化幅度逐漸減小。第2階段為重力減小階段,2015-10~2016-04重力場開始出現負變化,2016-04~2016-10重力場繼續減小,并且變化幅度增大。該階段內,九寨溝地震震源均位于正負重力變化過渡帶上。第3階段為重力場恢復調整階段,巴顏喀拉塊體東部表現出與第2階段截然不同的重力變化,即重力增加。并且在該階段,震源仍然位于正負重力變化過渡帶上。仔細對比第2階段和第3階段的正負過渡帶可以發現,前者由負向正過渡的方向為東北向,后者為東南向,兩者方向發生約90°旋轉。正負重力變化過渡帶旋轉,表明巖石圈內部物質遷移劇烈。祝意青等[13]指出,九寨溝地震之前震源周邊重力場變化原因有2個:1)2013年四川蘆山7.0級地震和甘肅岷縣-漳縣6.6級地震使青藏高原東部巖石圈物質重新分布;2)川西高原及鄰區深部殼幔邊界與上地幔物質和能量強烈交換。九寨溝地震前重力場變化體征顯示,區域正負重力變化過渡帶空間走向發生旋轉后,緊接著發生了九寨溝地震。因此,捕捉正負重力變化過渡帶空間走向旋轉變化信息對于尋找地震可能發生的時間具有重要前兆意義。
本文分析重力場動態演化規律時,將經重力網平差后的數據進行球諧展開至60階,并進行3°范圍平滑濾波,這樣處理后可以在圖2(b)、2(c)、和2(d)所示的3個時間段重力場變化圖像中清晰顯示重力正負變化過渡帶。3期正負重力變化過渡帶的交會區域基本集中在四川西北部,九寨溝地震震源即位于該區域內部。破壞性強震孕育會引起重力場在較大范圍內出現改變,本文中所用的球諧分析方法可以較好地突出大尺度重力場變化特征,抑制局部細節特征,這種特點有利于確定地震可能發生的大概位置。已有研究表明,強震易發生在重力變化正負梯度過渡帶的零值線附近[9,11,14-15],本文給出的結果顯示,九寨溝MS7.0地震震源位于正負過渡帶內,但與零值線存在一定距離,這與本文中采用的球諧分析技術突出主要重力場變化特征、觀測誤差、點位分布等因素有關。實際采用流動重力觀測資料進行地震前兆分析時,建議采用重力網平差后直接得到的結果和球諧分析方法得到的結果同時進行分析。
本文采用球諧分析方法研究九寨溝地震前重力場動態變化特征,取得的主要認識有:
1)2015-04~2017-04期間,震源周邊區域重力場變化經歷增大、減小和恢復調整3個階段,后2個階段在震源附近出現正負重力變化過渡帶,且重力變化過渡帶空間走向出現旋轉現象。重力場演變特征較好地反映了九寨溝地震的前兆現象。
2)本文使用的球諧分析方法,可以有效突出區域重力場演化的主要特征,便于地震前兆特征識別。但同時也會抑制局部細節變化,實際分析重力場變化數據時,建議與傳統方法同時使用。