2018年阿拉斯加Mw7.9級(jí)地震同震電離層擾動(dòng)研究

林清瑩, 于忠海, 張永明, 陳 鵬, 張 碩

(1. 濟(jì)南市勘察測(cè)繪研究院, 濟(jì)南 250101；2. 西安科技大學(xué)測(cè)繪科學(xué)與技術(shù)學(xué)院, 西安 710054；3. 山東農(nóng)業(yè)大學(xué)資源與環(huán)境學(xué)院, 山東 泰安 271018)

0 引言

地震是由于地殼快速釋放能量所造成的震動(dòng)，對(duì)地表具有很強(qiáng)的破壞性，給人類造成嚴(yán)重的危害，因此探測(cè)地震引起的異常現(xiàn)象具有重要的意義。地震發(fā)生時(shí)釋放的能量不僅會(huì)破壞地球表面，還會(huì)傳播至電離層并引起電子密度(Ne)和總電子含量(Total Electron Content，TEC)的變化等異常現(xiàn)象。自1964年美國(guó)阿拉斯加大地震中首次探測(cè)到地震電離層異常擾動(dòng)以來，電離層異常與地震之間的耦合關(guān)系已經(jīng)引起了廣泛關(guān)注[1]。隨著大量全球和區(qū)域衛(wèi)星導(dǎo)航系統(tǒng)(Global Navigation Satellite System，GNSS)網(wǎng)的建設(shè)和使用，GNSS技術(shù)已廣泛地應(yīng)用于地球動(dòng)力學(xué)研究中，突破了傳統(tǒng)地震監(jiān)測(cè)技術(shù)的時(shí)空限制。

通過全球定位系統(tǒng)(Global Positioning System，GPS)觀測(cè)數(shù)據(jù)探測(cè)發(fā)現(xiàn)，地震發(fā)生時(shí)激起的地震波向上傳播至電離層高度，會(huì)引起電離層總電子含量的變化，產(chǎn)生同震電離層擾動(dòng)(Coseismic Ionospheric Disturbances,CID)現(xiàn)象。E.Astafyeva等研究了1994年Kurile Mw8.1級(jí)地震同震電離層異常現(xiàn)象，發(fā)現(xiàn)CID出現(xiàn)的時(shí)間與距震中的距離存在相關(guān)性[2]。通過CID的傳播速度，可以研究電離層擾動(dòng)的激發(fā)條件，如陳鵬和陳家君發(fā)現(xiàn)2004年蘇門答臘Mw8.6級(jí)和Mw8.2級(jí)地震CID傳播速度分別為0.73km/s與0.69km/s，與聲波在0～450km內(nèi)的平均傳播速度一致[3]。Liu J.Y.等發(fā)現(xiàn)1999年Chi-Chi Mw7.6級(jí)地震引起的兩種CID傳播速度分別為2.5km/s和370m/s，認(rèn)為是由瑞利波和重力波引發(fā)的擾動(dòng)[4]。CID的傳播方向也是一個(gè)重要特征，K.Heki等發(fā)現(xiàn)因受地磁場(chǎng)的影響,北半球的CID表現(xiàn)出向南傳播的趨勢(shì)，擾動(dòng)幅度一般不超過2個(gè) TECU，擾動(dòng)峰值點(diǎn)出現(xiàn)在震后10～20min，持續(xù)幾分鐘后逐漸消失[5-6]，M.N.Cahyadi等發(fā)現(xiàn)破裂帶走向也會(huì)對(duì)CID傳播方向產(chǎn)生影響[7-8]。同震地表垂直位移被認(rèn)為是引發(fā)CID現(xiàn)象的主要原因，如Jin S.等發(fā)現(xiàn)瑞利波引起的CID主要受到垂直地面運(yùn)動(dòng)的影響，并會(huì)直接影響CID波形[9-10]。E. Astafyeva等發(fā)現(xiàn)CID的幅度和持續(xù)時(shí)間與地震的震級(jí)存在明顯的正相關(guān)關(guān)系[11]。此外，低軌衛(wèi)星等空基數(shù)據(jù)也被用于CID研究。Ma X.等利用cosmic數(shù)據(jù)，發(fā)現(xiàn)汶川地震期間TEC和NmF2時(shí)空分布變化存在同震響應(yīng)[12]。

GPS數(shù)據(jù)為探測(cè)電離層活動(dòng)提供了高時(shí)空分辨率數(shù)據(jù)支持，這有助于深入探測(cè)地震電離層異常效應(yīng)，促進(jìn)地震-電離層耦合機(jī)制的研究。本文對(duì)2018年1月23日在阿拉斯加灣科迪亞克島東南部地區(qū)發(fā)生的Mw7.9級(jí)地震進(jìn)行了研究，據(jù)美國(guó)地質(zhì)調(diào)查局(United States Geological Survey，USGS)公布，此次地震是由于太平洋板塊淺層巖石圈斷層走滑造成的，沿西南偏西或西北偏北方向發(fā)生破裂。本文采用美國(guó)連續(xù)運(yùn)行參考站(Continuously Operating Reference Station，CORS)發(fā)布的GPS觀測(cè)數(shù)據(jù)，解算孕震區(qū)周圍傾斜總電子含量(Slant Total Electron Content，STEC)數(shù)據(jù),利用奇異譜分析(Singular Spectrum Analysis，SSA)的方法提取電離層擾動(dòng)信號(hào)，對(duì)CID現(xiàn)象和較高緯度地區(qū)衛(wèi)星選取進(jìn)行進(jìn)一步分析，并利用SWARM衛(wèi)星垂直總電子含量(Vertical Total Electron Content，VTEC)數(shù)據(jù)對(duì)震中上空的總電子含量變化進(jìn)行研究。

1 數(shù)據(jù)和處理方法

1.1 電離層數(shù)據(jù)

本文采用CORS系統(tǒng)發(fā)布的GPS觀測(cè)數(shù)據(jù)，反演阿拉斯加地震孕震區(qū)周圍的震發(fā)1h內(nèi)STEC數(shù)據(jù)，處理后精度達(dá)到0.01個(gè)TECU[10]，可用于阿拉斯加震后的電離層擾動(dòng)探測(cè)。CORS系統(tǒng)在美國(guó)有1900多個(gè)觀測(cè)站，并且在阿拉斯加地區(qū)分布非常密集，為研究電離層異常擾動(dòng)提供了充足的數(shù)據(jù)，CORS系統(tǒng)在阿拉斯加地區(qū)的站點(diǎn)分布如圖1所示。利用采樣間隔為30s和15s的GPS數(shù)據(jù)反演的TEC中均發(fā)現(xiàn)了異常現(xiàn)象，為了提高精度，本文選取了采樣間隔15s的GPS觀測(cè)數(shù)據(jù)(ftp://data-out.unavco.org/pub/rinex/obs/)作為實(shí)驗(yàn)數(shù)據(jù)，進(jìn)行CID探測(cè)。

圖1 CORS站點(diǎn)分布、震源球表示震中位置

利用GPS相位觀測(cè)數(shù)據(jù)計(jì)算相對(duì)STEC的關(guān)系式

STECΔ=

(1)

其中，f1、f2為載波頻率；L1、L2為載波相位觀測(cè)值；λ1、λ2為載波的波長(zhǎng)；N1、N2為相位觀測(cè)值的整周模糊度；ε為觀測(cè)噪聲。由于本文采用高精度GPS數(shù)據(jù)，求取STEC的相對(duì)精度很高，為了減少電離層水平梯度的影響，衛(wèi)星的截止高度角選取為15°。

采用歐洲航天局(European Space Agency，ESA)發(fā)布的SWARM衛(wèi)星Level2級(jí)TEC數(shù)據(jù)(https://swarm-diss.eo.esa.int/)，VTEC數(shù)據(jù)采樣間隔為1Hz，將相同經(jīng)緯度格網(wǎng)點(diǎn)的VTEC數(shù)據(jù)求均值，作為相應(yīng)高度的電離層總電子含量。

1.2 太陽(yáng)-地磁活動(dòng)

太陽(yáng)和地磁環(huán)境的變化會(huì)造成電離層活動(dòng)的異常干擾[14-16]，在分析阿拉斯加地震同震電離層異常效應(yīng)前，應(yīng)首先對(duì)太陽(yáng)-地磁活動(dòng)狀況進(jìn)行判斷。Dst指數(shù)可以反映地磁活動(dòng)變化情況，正常情況下在-30nT～30nT范圍內(nèi)，而太陽(yáng)輻射指數(shù)F10.7可以反映太陽(yáng)活動(dòng)的變化情況，通常太陽(yáng)活動(dòng)平靜時(shí)F10.7指數(shù)小于90SFU。本文利用日本京都地磁數(shù)據(jù)中心提供的磁暴時(shí)變化指數(shù)Dst(http://wdc.kugi.kyoto-u.ac.jp/dstae/index.html)數(shù)據(jù)和空間環(huán)境預(yù)報(bào)中心提供的太陽(yáng)輻射指數(shù)F10.7(http://www.sepc.ac.cn/F107Index_chn.php)的觀測(cè)數(shù)據(jù)對(duì)太陽(yáng)-地磁活動(dòng)狀況進(jìn)行分析。

圖2給出了包括震發(fā)當(dāng)天在內(nèi)10天的地磁Dst指數(shù)和太陽(yáng)射電量F10.7的變化。由圖2可知，在包括震發(fā)在內(nèi)的10天時(shí)間里，|Dst|<30，地磁場(chǎng)比較平靜；同時(shí)太陽(yáng)射電量F10.7指數(shù)<90SFU，太陽(yáng)活動(dòng)在這個(gè)時(shí)期很穩(wěn)定，說明震發(fā)當(dāng)天太陽(yáng)-地磁活動(dòng)對(duì)電離層造成的影響很小。

(a) Dst指數(shù)

1.3 方法

從STEC中提取地震電離層擾動(dòng)信號(hào)，需要去除由GPS衛(wèi)星的沿軌運(yùn)動(dòng)、電離層的空間變化和時(shí)間變化等引起的電離層TEC背景的趨勢(shì)項(xiàng)。因此，本文采用SSA的方法進(jìn)行提取處理，該方法通過對(duì)一維時(shí)間序列的分解與重構(gòu)，識(shí)別出原始序列中的不同信號(hào)，提取出趨勢(shì)分量、周期分量與噪聲分量等。

(1)構(gòu)建遲滯矩陣

對(duì)反演出的震發(fā)時(shí)段UTC 9-10h STEC數(shù)據(jù)，按照15s的采樣頻率，將1h的STEC數(shù)據(jù)長(zhǎng)度N定義為240，將其組成x1,x2,…,xN的時(shí)間序列{x}，建立時(shí)滯矩陣。

(2)奇異值分解(Singular Value Decomposition，SVD)

定義矩陣S=XXT，XT為X的轉(zhuǎn)置矩陣，設(shè)λ1,…,λM為矩陣S的特征值，U1,…,UM分別為λ1,…,λM對(duì)應(yīng)的特征向量，其中λ1≥…≥λM≥0。

(3)分組

將初等矩陣{1,…,d}劃分為m個(gè)不相交的子集I1,I2,…,Im，其中I={i1,…,ip}。軌跡矩陣X的SVD可以表達(dá)為：XI=XI1+…+XIm。分組即指對(duì)I1,I2,…,Im的確定過程。

(4)對(duì)角平均化

將分組得到的矩陣轉(zhuǎn)換為一系列新的長(zhǎng)度為N的時(shí)間序列。長(zhǎng)度為N的時(shí)間序列定義為重建成分(Reconstruction Components，RC)，原始序列即所有重建成分之和。圖3給出了前15重構(gòu)成分的w-correlation圖。由圖3可知，分解的前5項(xiàng)重構(gòu)向量與其他重構(gòu)向量獨(dú)立性較好。

圖3 STEC數(shù)據(jù)SSA前15重構(gòu)成分的w-correlation

(5)獲取擾動(dòng)信號(hào)

利用SVD得到l個(gè)重構(gòu)分量和對(duì)應(yīng)的特征向量λi(1 ≤i≤l)，將重構(gòu)分量的特征值λK按照大小排列λ1≥λ2≥ …≥λl≥ 0。特征值的大小對(duì)應(yīng)的特征向量代表了TEC信號(hào)的變化趨勢(shì)的大小。考慮到前5個(gè)重構(gòu)分量的特征值累積貢獻(xiàn)率超過95%，截取前5個(gè)較大的特征值，找到對(duì)應(yīng)的xk之和，即可以重建出充分反映TEC原序列的整體特征，即

(2)

定義包含電離層短周期信號(hào)的重構(gòu)STEC時(shí)間序列主成分為STECmain，為提取擾動(dòng)信號(hào)提供高精度的背景參考值。

用原始STEC時(shí)間序列減去STECmain，得到電離層STEC擾動(dòng)信號(hào)

ΔSTEC=STEC-STECmain

(3)

2 結(jié)果與分析

2.1 CID現(xiàn)象

選取震中附近測(cè)站和GPS PRN05號(hào)衛(wèi)星的STEC數(shù)據(jù)，基于SSA方法提取震發(fā)時(shí)刻UTC 9:00-10:00 ΔSTEC時(shí)間序列。圖4給出了各測(cè)站的ΔSTEC時(shí)間序列和衛(wèi)星軌跡分布圖。由圖4(a)可知，在震后0.5h內(nèi)，測(cè)站ab07、ac42、av06和av10的ΔSTEC序列出現(xiàn)了明顯的N形波擾動(dòng)，且先后出現(xiàn)2個(gè)峰值點(diǎn)。在測(cè)站ac42的ΔSTEC序列中，UTC 9:39開始出現(xiàn)擾動(dòng)，經(jīng)過了2次由負(fù)到正的變化，第1個(gè)峰值點(diǎn)出現(xiàn)在UTC 9:43，擾動(dòng)幅度達(dá)到了0.15個(gè)TECU;經(jīng)過幾分鐘后，在UTC 9:49出現(xiàn)了第2個(gè)峰值點(diǎn)，擾動(dòng)幅度為0.06個(gè)TECU，在UTC 9:54后恢復(fù)平靜，整個(gè)擾動(dòng)過程持續(xù)了15min左右。第1個(gè)峰值點(diǎn)的SIP(Sub-Ionospheric Point)距離震中177.98km，第2個(gè)峰值點(diǎn)的SIP距離震中186.7km，且SIP位于震中西南方向區(qū)域,如圖4(b)所示。

(a) ΔSTEC時(shí)間序列

2.2 觀測(cè)衛(wèi)星選取

由于地震產(chǎn)生的聲波會(huì)向各個(gè)方向傳播，能量在傳播過程中不斷衰減，因此在基于STEC擾動(dòng)分析同震電離層異常現(xiàn)象過程中，不僅要保證測(cè)站對(duì)于觀測(cè)衛(wèi)星有連續(xù)的觀測(cè)信號(hào)，還要選取合適的GPS觀測(cè)衛(wèi)星，才能確保阿拉斯加地震同震電離層異常探測(cè)的有效性。

本文選取測(cè)站av40的STEC擾動(dòng)結(jié)果進(jìn)行分析，如圖5(a)所示，在UTC 9:00-10:00經(jīng)過測(cè)站上空的6顆具有連續(xù)觀測(cè)信號(hào)的衛(wèi)星中，只有PRN05號(hào)衛(wèi)星在UTC 9:40探測(cè)出了明顯的擾動(dòng)現(xiàn)象，擾動(dòng)幅度達(dá)到0.12個(gè)TECU，峰值點(diǎn)出現(xiàn)在UTC 9:44和UTC 9:49，擾動(dòng)過程持續(xù)了15min，其他觀測(cè)衛(wèi)星未發(fā)現(xiàn)明顯的擾動(dòng)現(xiàn)象。圖5(b)給出了6顆衛(wèi)星的SIP軌跡，可以看出只有PRN05號(hào)衛(wèi)星SIP軌跡在地震附近區(qū)域。

(a) ΔSTEC時(shí)間序列

由于聲波高度角低時(shí)，其在對(duì)流層和中氣層中接近水平方向傳播；高度角高時(shí)，聲波沿著垂直方向到達(dá)電離層，再沿水平方向在電離層中傳播[5]，只有高度角較大的聲波才能傳播到電離層高度引起電離層擾動(dòng)。聲波只有通過衛(wèi)星和測(cè)站的連線路徑，才能探測(cè)到震中上空附近電離層擾動(dòng)情況，而高緯度地區(qū)測(cè)站觀測(cè)到的衛(wèi)星數(shù)目較少且SIP大多距離測(cè)站超過2000km，因此選擇合適的GPS衛(wèi)星對(duì)該地區(qū)同震電離層研究特別重要。

2.3 CID傳播速度和擾動(dòng)幅度

為了進(jìn)一步分析CID傳播特征，利用高度為350km的電離層單層模型獲得CID的位置信息。在分析STEC擾動(dòng)過程中，由衛(wèi)星接收機(jī)連線與電離層薄層的交點(diǎn)可以獲得電離層穿刺點(diǎn)(Ionosph-eric Pierce Point, IPP)位置。IPP在地面上的投影點(diǎn)為層下點(diǎn)SIP，計(jì)算CID峰值點(diǎn)處的SIP經(jīng)緯度，根據(jù)CID出現(xiàn)的時(shí)間和對(duì)應(yīng)的SIP點(diǎn)震中距即可計(jì)算CID的水平傳播速度。

圖6給出了2018.1.23 UTC 9:00-10:00 震中附近ΔSTEC擾動(dòng)序列的時(shí)間-距離圖。根據(jù)CID峰值點(diǎn)出現(xiàn)的時(shí)間和相應(yīng)的位置，利用最小二乘進(jìn)行線性擬合，斜率近似等于震后UTC 9:00-10:00出現(xiàn)的CID傳播速度。

圖6 阿拉斯加地震ΔSTEC變化走勢(shì)(UTC 9：00-10：00),藍(lán)色虛線表示震發(fā)時(shí)刻，黑線代表根據(jù)CID峰值點(diǎn)擬合的趨勢(shì)線

由圖6可以看出，CID按照特定的速度進(jìn)行傳播，并且發(fā)現(xiàn)CID振幅隨著距離增加而減小。CID在震后8min距離震中42.91km處最先被探測(cè)到，出現(xiàn)在震后10～20min內(nèi)，和2011年Mw9.1級(jí)日本地震、2008年Mw7.8級(jí)汶川地震的CID出現(xiàn)時(shí)間一致[17]。CID的傳播速度為(2.62±0.04)km/s，符合地震產(chǎn)生的瑞利波向上傳播到電離層高度引起的CID水平傳播速度[5]。

2.4 CID方向性

利用GPS跟蹤站的觀測(cè)信號(hào)，可以得到GPS跟蹤站與PRN05、23、31衛(wèi)星的SIP軌跡，同時(shí)得到跟蹤站與衛(wèi)星連線的STEC擾動(dòng)信息。本文基于阿拉斯加地震震中附近CORS網(wǎng)測(cè)站，并分析這些站點(diǎn)上空的電離層擾動(dòng)情況，結(jié)果如圖7所示。

圖7 UTC 09:00-10:00阿拉斯加地區(qū)SIP點(diǎn)分布圖，震中上方和西南方向擾動(dòng)情況及SIP分布，顏色深淺代表電離層擾動(dòng)幅度

圖7給出了震發(fā)當(dāng)天UTC 9:00-10:00在震中上方和西南方向探測(cè)到的電離層異常和SIP軌跡，可以看到明顯的CID，擾動(dòng)幅度達(dá)到0.16個(gè)TECU，出現(xiàn)在距離震中42.9～488km范圍內(nèi)。圖4(b)給出了震中東北方向的電離層擾動(dòng)和SIP軌跡，并沒有探測(cè)到明顯的CID信號(hào)，CID傳播表現(xiàn)出了明顯的方向差異性。

由圖7可以看出，CID從震中沿著西南方向傳播，一方面由于瑞利波的傳播受到地震破裂帶的影響，引起CID傳播方向和破裂方向保持一致；另一方面，由于CID現(xiàn)象的存在主要是通過帶電粒子的運(yùn)動(dòng)形式反映的，地磁場(chǎng)施加的庫(kù)侖力使其更容易沿著磁力線運(yùn)動(dòng)，因此在北半球的CID表現(xiàn)出向南傳播的趨勢(shì)。

為了進(jìn)一步分析瑞利波傳播至電離層高度引起的擾動(dòng)情況，圖8給出了震中上方和西南方向CID峰值點(diǎn)分布位置，可以看出地震上空的擾動(dòng)值最大，達(dá)到0.16個(gè)TECU，CID幅度沿著西南方向隨著距離增加而減小，然后逐漸消失，這進(jìn)一步證明了瑞利波引起的CID沿著西南方向傳播，并且CID能量以較快速度衰減，符合瑞利波特性[17]。

圖8 震中上方和西南方向CID峰值點(diǎn)示意圖，震源球代表震發(fā)位置，顏色深淺代表峰值點(diǎn)幅度

2.5 SWARM衛(wèi)星數(shù)據(jù)

利用軌道高度為460km的SWARM-A衛(wèi)星觀測(cè)得到的VTEC數(shù)據(jù)，選取了包括震發(fā)當(dāng)天在內(nèi)前后3天的觀測(cè)數(shù)據(jù),對(duì)震發(fā)上空電離層異常進(jìn)行研究。SWARM衛(wèi)星在固定軌道高度運(yùn)行，在經(jīng)過同一區(qū)域時(shí)相對(duì)于地面高度一致，且穿過電離層高度。

圖9(a)給出了衛(wèi)星經(jīng)過的SIP軌跡的位置及震中位置的分布圖，可以看到,SWARM-A衛(wèi)星經(jīng)過震中和西南方向區(qū)域上空，持續(xù)時(shí)間為3min左右。震發(fā)當(dāng)天從UTC 20:08(黑色軌跡)開始經(jīng)過震中區(qū)域，與相同區(qū)域前一天UTC 20:43(紅色軌跡)和后一天UTC 19:17(藍(lán)色軌跡)總電子含量進(jìn)行對(duì)比。圖9(b)給出了前后3天在相同高度444～446km上的VTEC，可以看出震發(fā)當(dāng)天與前一天和后一天相比，總電子含量明顯變小。在排除了太陽(yáng)和地磁活動(dòng)影響的情況下，電子含量明顯降低可能是由于地震造成的。

(a) 衛(wèi)星軌跡及震中位置

3 結(jié)論

本文基于CORS提供的GPS觀測(cè)數(shù)據(jù)，解算震中上空高精度的STEC數(shù)據(jù)，利用SSA方法提取電離層擾動(dòng)信息，分析阿拉斯加Mw7.9級(jí)地震CID現(xiàn)象。

分析表明，震后8min出現(xiàn)了明顯的電離層擾動(dòng)，先后出現(xiàn)2個(gè)峰值點(diǎn)，擾動(dòng)幅度達(dá)到0.16個(gè)TECU，擾動(dòng)時(shí)間持續(xù)了15min。CID在距震中600km的范圍內(nèi)傳播，傳播速度為2.62km/s,與瑞利波水平方向傳播速度一致，利用頻譜分析發(fā)現(xiàn)異常時(shí)刻中心頻率為4.87mHz,符合瑞利波頻率。在阿拉斯加地震的CID現(xiàn)象中，發(fā)現(xiàn)擾動(dòng)具有明顯的方向性，CID傳播方向受到地震破裂帶和地磁場(chǎng)的影響。CID在西南方向傳播的過程中，能量以較快速度衰減，符合瑞利波特性。結(jié)合聲波傳播高度角和衛(wèi)星SIP軌跡，發(fā)現(xiàn)緯度較高地區(qū)的CID探測(cè)需要先選取合適的觀測(cè)衛(wèi)星。最后,利用SWARM衛(wèi)星的VTEC數(shù)據(jù)，發(fā)現(xiàn)震發(fā)當(dāng)天震中上空區(qū)域總電子含量明顯減小，進(jìn)一步證明了阿拉斯加地震存在電離層擾動(dòng)現(xiàn)象。

上述研究表明，利用SSA方法處理由GPS數(shù)據(jù)反演得到高精度的STEC數(shù)據(jù)，可以探測(cè)電離層異常效應(yīng)。此外，需要充足的高時(shí)空分辨率的電離層數(shù)據(jù)進(jìn)行更詳細(xì)的探測(cè)，并結(jié)合大量震例研究地震釋放能量方式對(duì)電離層的影響，進(jìn)一步了解地震電離層耦合機(jī)制。