虞良斌 高 飛 屈小川 安家春 綦子民
1 合肥工業(yè)大學(xué)土木與水利工程學(xué)院,合肥市屯溪路193號(hào),230009 2 武漢大學(xué)中國(guó)南極測(cè)繪研究中心,武漢市珞喻路129號(hào),430079
地震是地殼快速釋放能量過(guò)程中造成的地面振動(dòng),因此,可以利用震前溫度或電離層異常擾動(dòng)的動(dòng)態(tài)變化來(lái)監(jiān)測(cè)、預(yù)測(cè)地震[1-3]。但很多研究在分析震前異常征兆時(shí)忽略了溫度和電離層異常擾動(dòng)在時(shí)間和空間上的內(nèi)在聯(lián)系和發(fā)展。
本文以2016-11-13新西蘭南島凱庫(kù)拉發(fā)生的MW7.8大地震(42.53°S,173.05°E)為例,利用數(shù)值天氣模式資料和GNSS觀測(cè)數(shù)據(jù)對(duì)震前溫度及電離層異常發(fā)生、發(fā)展過(guò)程及內(nèi)在關(guān)聯(lián)性進(jìn)行初步探討,以期為震前異常征兆監(jiān)測(cè)提供參考與借鑒。
顧及接收機(jī)端和衛(wèi)星端的硬件延遲,建立如式(1)的區(qū)域VTEC多項(xiàng)式模型,利用載波相位平滑偽距[4]和最小二乘法[5]求解模型系數(shù):
9.524 37(ΔBS+ΔBR)cosz′
(1)
式中,n、m為階數(shù),Eik為待求解的多項(xiàng)式模型系數(shù),φ0為測(cè)區(qū)中心點(diǎn)的地理緯度,S0為測(cè)區(qū)中心點(diǎn)(φ0,λ0)在該時(shí)段中央時(shí)刻t0的太陽(yáng)時(shí)角,z′是穿刺點(diǎn)處的天頂距,ΔBS為待求的衛(wèi)星硬件延遲,ΔBR為待求的接收機(jī)硬件延遲。
從GeoNet網(wǎng)站(http:∥www.geonet.org.nz)下載2016年新西蘭南島地區(qū)14個(gè)GPS連續(xù)跟蹤站的觀測(cè)數(shù)據(jù),跟蹤站分布見(jiàn)圖1(因NLSN站和WITH站位置接近,故只顯示了WITH站)。為驗(yàn)證解算VTEC值的準(zhǔn)確性,基于式(1)求得各測(cè)站的VTEC,并與IGS提供的GIM結(jié)果進(jìn)行對(duì)比。圖2給出2016-11-01~07 CMBL站的VTEC值與 IGS提供的GIM數(shù)據(jù)??梢钥闯?,兩者吻合度較高,最小差值和最大差值分別為0.02 TECu和3.1 TECu,7 d的RMS僅為0.713 TECu。表明本文解算的VTEC值精度較高,能夠用于震前電離層異常檢測(cè)研究。
圖1 GPS連續(xù)跟蹤站與震中位置Fig.1 GPS continuous tracking stations and epicenter location
圖2 CMBL站建模的VTEC與IGS的GIM對(duì)比Fig.2 Comparison between modeled VTEC and GIM of IGS at CMBL station
滑動(dòng)四分位距法[6]是當(dāng)前VTEC異常提取中常用的有效方法之一。依據(jù)式(2)、(3)確定VTEC值上、下閾值U上、U下:
U上=Q2+k(Q3-Q1)
(2)
U下=Q2-k(Q3-Q1)
(3)
式中,Q1、Q2和Q3分別為下、中和上四分位,k為倍常數(shù),決定閾值的大小。當(dāng)k取1.5時(shí),上下閾值差約為標(biāo)準(zhǔn)差的2倍,此時(shí)異常探測(cè)的置信度為95%[7],且k越大,模型對(duì)異常值的敏感度越低。當(dāng)某時(shí)刻檢測(cè)數(shù)據(jù)大于或者小于閾值時(shí),判定該時(shí)刻VTEC出現(xiàn)異常,大于上閾值為正異常,小于下閾值為負(fù)異常。
ECMWF最新發(fā)布的ERA5再分析資料具有高時(shí)間分辨率(1 h)和高空間分辨率(約31 km)[8]。本文下載震前15 d以及震后7 d的溫度數(shù)據(jù),數(shù)據(jù)范圍以震中為中心,覆蓋整個(gè)地震發(fā)生區(qū)域(169°~177.25°E,37.75~47°S)。
圖像差值法直觀易懂,且能在一定程度上排除地形、地面和氣象等客觀因素帶來(lái)的干擾,是檢測(cè)溫度異常的常用方法之一,但該方法必須選取合適的參考背景值[9]。天體引潮力作為一種能夠在震前計(jì)算的地殼形變物理參數(shù),其反映的地質(zhì)構(gòu)造作用與地震熱異常本質(zhì)是一致的,因此可以作為選取參考背景溫度值的依據(jù)[10]。
綜上,本文使用震前引潮力及ERA5的溫度數(shù)據(jù),利用圖像差值法提取震前溫度異常。
基于滑動(dòng)四分位距法提取14個(gè)測(cè)站震前20 d及震后3 d(2016-10-25~11-16)的VTEC異常值。圖3(e)給出了CMBL站的VTEC異常擾動(dòng)值。從圖中可以看出,該站于10-28、11-03、06、09分別出現(xiàn)3.5 TECu、1.3 TECu、0.8 TECu、1.8 TECu的正異常。考慮到VTEC同樣會(huì)受到日地空間活動(dòng)的影響而產(chǎn)生異常,因此需排除與日地空間活動(dòng)有關(guān)的VTEC異?,F(xiàn)象。從美國(guó)NGDC和日本京都地磁中心下載同期的4個(gè)反映日地空間活動(dòng)的指數(shù),并依據(jù)表1中的常用設(shè)置[11-12]將其分成3個(gè)級(jí)別,分別表示低、中和高3種日地活動(dòng)水平。
圖3 2016-10-25~11-16日地空間活動(dòng)指數(shù)Fig.3 The indices of terrestrial space activityfrom October 25 to November 16, 2016
根據(jù)表1及圖3可以看出,F(xiàn)10.7指數(shù)在23 d內(nèi)都未發(fā)生異常變化,其余3種活動(dòng)指數(shù)在2016-11-04~09、15、16這8 d內(nèi)都處于正常水平,而在其余時(shí)段內(nèi)3者都發(fā)生了或強(qiáng)或弱的異常變化,因此可以認(rèn)為這8 d內(nèi)的VTEC異常與此次地震相關(guān)。CMBL站在2016-11-06、09的正異常與此次地震有關(guān)。
表1 各日地空間活動(dòng)指數(shù)的活動(dòng)水平
為了進(jìn)一步研究各測(cè)站VTEC異常與地震的關(guān)聯(lián)性,給出各測(cè)站VTEC異常擾動(dòng)時(shí)間及量級(jí)(圖4)。結(jié)合圖3日地活動(dòng)水平可知,除了2016-11-06、09,其余時(shí)間段的VTEC異常都與日地空間活動(dòng)有關(guān)。在2016-11-06,除了震中西南方向的HOKI和MRBL站外,其余12個(gè)測(cè)站都出現(xiàn)VTEC正異常,最大正異常發(fā)生在CMBL站,約1.8 TECu;而在2016-11-09,除了震中東北方向的WGTN和TORY站外,其余12個(gè)測(cè)站都出現(xiàn)VTEC正異常,最大正異常發(fā)生在WITH站,約2.5 TECu。對(duì)比各測(cè)站的VTEC異常值可知,本次地震最大VTEC異常值出現(xiàn)在震中偏北的位置,與本次地震能量釋放方向一致。此外,相同方位測(cè)站VTEC異常值大小隨著震中距的增大而略微減小。
圖4 各測(cè)站出現(xiàn)VTEC擾動(dòng)時(shí)間及量級(jí)Fig.4 Time and magnitude of VTEC disturbances at each station
引潮力是周期性變化的,由其誘發(fā)的地震多發(fā)生在周期振幅相位最大值或最小值處。從圖5可以看出,新西蘭地震發(fā)生時(shí)引潮力處于高峰,反映了引潮力對(duì)地應(yīng)力達(dá)到臨界點(diǎn)的活動(dòng)構(gòu)造具有誘發(fā)作用。2016-11-07引潮力處于低谷,因此可作為本次地震溫度異常提取的參考背景值。將2016-11-08~17的溫度值與2016-11-07同時(shí)刻的溫度值相減,就可得到同區(qū)域不同氣壓層的溫度變化圖像。為避免日間太陽(yáng)輻射的影響,選取夜間(世界時(shí)12:00)的溫度數(shù)據(jù)進(jìn)行震前溫度異常檢測(cè)[13]。
圖5 凱庫(kù)拉MW7.8地震天體引潮力時(shí)序變化曲線Fig.5 Time series curve of astro-tidal-triggering of Kaikula MW7.8 earthquake
圖6為研究區(qū)2016-11-08~17的溫度異常時(shí)空變化??梢钥闯觯鹬斜毕蚝蜄|向位置在震前及震后都出現(xiàn)異常增溫。由于地震熱異常在自下而上的傳導(dǎo)過(guò)程中逐漸減弱,因此隨高度上升,溫度異常增幅呈現(xiàn)出明顯的降低,越接近地表面,增溫現(xiàn)象越明顯。2016-11-08~10所有氣壓層的異常增溫區(qū)域和幅度都逐漸增大,同時(shí)向震中移動(dòng),并在震前3 d(2016-11-10)達(dá)到最大,可達(dá)6 K。隨后2016-11-11~12異常增溫出現(xiàn)衰減;2016-11-13地震發(fā)生前,溫度異常僅增加3 K。地震發(fā)生后,由于地表大面積破裂,形成熱通道,釋放大量熱量,故地震當(dāng)天(2016-11-13)和震后1 d(2016-11-14)增溫面積大幅度增加,增溫異常區(qū)覆蓋整個(gè)震中區(qū)域,且增溫幅度達(dá)到12 K。震后2 d(2016-11-15)增溫面積較前1 d(2016-11-14)稍有減少,但局部地區(qū)增幅仍然達(dá)到10 K。2016-11-16~17增溫快速衰減,且逐漸恢復(fù)平靜。綜上,此次地震前后溫度異常過(guò)程為:增溫(2016-11-08)-增溫加強(qiáng)(2016-11-09~10)-增溫衰減(2016-11-13)-再次增溫-達(dá)到增溫高峰(2016-11-14~15)-快速衰減(2016-11-16~17)。這與2001-11-14昆侖山MS8.1地震[13]、2013-04-20蘆山MS7.0地震[14]震前溫度出現(xiàn)的起始增溫-增溫加強(qiáng)-高峰-快速衰減異常特征相似。
基于天體引潮力提取的震前增溫異常具有增加-削弱-高峰-衰減的變化特點(diǎn),這與震前巖石破裂遵循的應(yīng)力微破裂-破裂-閉鎖規(guī)律有著很強(qiáng)的相似性[15]。在不同的氣壓層,溫度異常分布具有相似性,但隨著氣壓的不斷抬升,溫度異常逐漸衰減。本文同樣處理了與引潮力相位變化相似的2個(gè)時(shí)段的溫度數(shù)據(jù),并未發(fā)現(xiàn)溫度異常現(xiàn)象,因此認(rèn)為,此次溫度異常與地震活動(dòng)有關(guān)。2016-11-08出現(xiàn)溫度異常,隨后2016-11-09出現(xiàn)電離層異常,且電離層異常區(qū)域與溫度異常區(qū)域具有高度一致性(圖7),這進(jìn)一步表明震前溫度異常與電離層異常存在一定關(guān)聯(lián)。對(duì)于震前電離層異常以及溫度異常的成因至今沒(méi)有確切的解釋。電離層異常主要有聲重力波說(shuō)和電磁輻射說(shuō),其中,電磁輻射說(shuō)[16]認(rèn)為地球巖石的破裂會(huì)導(dǎo)致內(nèi)部的帶電粒子溢出進(jìn)而引起電磁場(chǎng)的異常。對(duì)于臨震低空大氣增溫的成因[17]主要有2種學(xué)說(shuō):一是大氣靜電場(chǎng)的極化作用;二是異常大氣的增溫效應(yīng),即在地表微破裂后會(huì)溢出N2、CO2等異常增溫氣體。震前溫度異常并不是單一因素引起的,而是多種因素共同作用的復(fù)雜過(guò)程。在此次地震中,震前5 d(2016-11-08)出現(xiàn)的溫度異常很大程度上是由于地表微破裂溢出的帶電粒子在大氣靜電場(chǎng)的極化作用下導(dǎo)致的,這種作用在自下而上的傳遞過(guò)程中導(dǎo)致了電離層異常;而震后1~2 d(2016-11-14~15)的大面積溫度異??赡苁怯捎诘乇淼耐耆屏褜?dǎo)致溫室氣體溢出,在自下而上的傳遞過(guò)程中并不會(huì)達(dá)到電離層,因此震后未監(jiān)測(cè)出明顯的電離層異常。這種震前同區(qū)域出現(xiàn)的溫度異常和電離層異??勺鳛榈卣甬惓z測(cè)的征兆,以提高地震預(yù)測(cè)的準(zhǔn)確性。
本文依據(jù)ERA5再分析溫度數(shù)據(jù)和電離層VTEC數(shù)據(jù)捕捉2016-11-13新西蘭MW7.8地震發(fā)生前的空間大氣擾動(dòng)現(xiàn)象,再次驗(yàn)證了強(qiáng)震來(lái)臨前會(huì)發(fā)生電離層異常和大氣溫度異常,且兩者并不是孤立存在的。文中大氣溫度異?,F(xiàn)象始于2016-11-08,隨著時(shí)間推移具有增溫-增溫加強(qiáng)-增溫衰減-增溫高峰-快速衰減的變化特征。由于地震熱異常在自下而上的傳導(dǎo)過(guò)程中逐漸減弱,因此增溫幅度隨高度呈衰減趨勢(shì)。此后,同方位區(qū)域于2016-11-09出現(xiàn)電離層異?,F(xiàn)象,表明震前大氣擾動(dòng)存在一定的時(shí)空關(guān)聯(lián)性,這一特征有助于強(qiáng)震預(yù)測(cè)。
致謝:感謝中國(guó)地震臺(tái)網(wǎng)中心馬未宇老師提供天體引潮力數(shù)據(jù)。