利用InSAR技術(shù)研究新疆阿什庫勒火山群現(xiàn)今活動性

季靈運 許建東 趙 波 萬 園 申歡歡

1)中國地震局地質(zhì)研究所，活動構(gòu)造與火山重點實驗室，北京 100029

2)中國地震局第二監(jiān)測中心，西安 710054

3)中國地震局地質(zhì)研究所，北京 100029

0 引言

阿什庫勒火山群位于新疆于田縣以南250km處的阿什庫勒盆地內(nèi)部(圖1)。阿什庫勒盆地位于康西瓦斷裂和阿爾金斷裂的交會處，構(gòu)造運動活躍，中強地震頻繁發(fā)生(圖1)。阿什庫勒火山群由10余座主火山和數(shù)十個子火山組成，區(qū)域海拔超過4 600m，是世界上海拔最高的火山區(qū)之一，熔巖面積達200km2(劉若新，2000)。新疆日報曾報道阿什庫勒火山群于1951年5月27日爆發(fā)過(新疆日報，1951)，但存在爭議(鄧萬明，1989;劉嘉麒等，1990;劉若新，2000)。2011年5月中國地震局地質(zhì)研究所和新疆維吾爾自治區(qū)地震局對阿什庫勒火山群開展了科學(xué)考察，在阿什火山北側(cè)火口沿發(fā)現(xiàn)了疑似1951年噴發(fā)遺留的噴氣孔(許建東等，2011)。如果1951年的噴發(fā)事件真的存在，那么阿什庫勒火山群現(xiàn)今是否還在活動?是否具有潛在的噴發(fā)危險性?另外，2008年3月21日，阿什庫勒火山群南部山區(qū)發(fā)生了MS7.3強烈地震(陳學(xué)忠等，2008;圖1)，本文稱之為2008年于田地震，此次強震事件對火山活動是否有影響?

圖1 阿什庫勒火山群區(qū)域構(gòu)造與地震震中分布圖Fig.1 Tectonic and epicentral distribution map around Ashikule volcanic group.

阿什庫勒火山群地區(qū)自然條件極其惡劣，研究程度較低，目前的研究主要集中在火山地質(zhì)與年代學(xué)方面(Liu et al.，1989;劉嘉麒等，1990;許建東等，2011)，迄今為止，還沒有地表形變方面的報道。惡劣的氣候環(huán)境使得人類難以接近，因而以地面控制點為基礎(chǔ)的常規(guī)地表形變監(jiān)測手段如GPS、水準(zhǔn)測量很難實施。近年來，作為一種全新的空間對地觀測手段，InSAR(Interferometric Synthetic Aperture Radar)技術(shù)因其具有不需要地面控制點，且能夠以高精度、低成本、快速度的優(yōu)勢完成大范圍、高空間分辨率的“面”的地表形變信息，而廣泛應(yīng)用于地學(xué)研究的各個領(lǐng)域(Massonnet et al.，1993，1995;單新建等，2006;孫建寶等，2008)。

本文基于Envisat ASAR和ALOS PALSAR衛(wèi)星影像資料，利用干涉圖堆疊技術(shù)(Stacking)和短基線集合成孔徑雷達干涉測量技術(shù)(SBAS-InSAR，Small BAseline Subsets-Interferometric Synthetic Aperture Radar)提取阿什庫勒火山群地區(qū)的地表形變，了解并分析其在2008年于田地震前后的活動特征，填補該火山群現(xiàn)今活動性未知的空白，進而為該火山群的噴發(fā)預(yù)測、災(zāi)害評估等工作提供基礎(chǔ)資料。

1 SAR影像資料與InSAR數(shù)據(jù)處理

1.1 SAR影像簡介

為了檢測阿什庫勒火山群地區(qū)在2008年3月于田地震前后可能存在的形變信息，我們購買和收集到了Envisat ASAR和ALOS PALSAR兩種衛(wèi)星影像資料。其中Envisat ASAR影像共41景(1個航跡，降軌)，其中于田地震之前20景，震后21景，時間跨越2003—2010年;ALOS PALSAR影像共8景(1個航跡，升軌)，其中震前3景，震后5景，時間跨越2007—2009年。

1.2 InSAR數(shù)據(jù)處理方法

首先，對于所有影像，不指定公共主影像，進行自由組合干涉?？紤]到阿什庫勒火山群地區(qū)地處高原無人區(qū)，地表無植被覆蓋，長時間間隔仍能保持較高的相干性，但山區(qū)地形起伏較大，宜采用較短垂直基線的干涉圖，所以給定干涉限制條件為:時間間隔短于5a，垂直基線＜300m(Envisat ASAR)與＜3 000m(ALOS PALSAR)。

對于每個符合條件的干涉對，進行常規(guī)兩軌法D－InSAR處理，使用美國NASA公布的SRTM DEM(分辨率約為90m)去除地形相位部分。目標(biāo)是提取火山區(qū)可能存在的大范圍地表形變異常，所以干涉圖的多視因子選定4×20(距離向×方位向)。則對于Envisat ASAR影像生成的干涉對，每個像元約為80×80m2，對于ALOS PALSAR，每個像元約為60×60m2。干涉圖相位解纏采用最小費用流算法(Minimum Cost Flow)，僅解纏相干系數(shù)值＞0.5的像元。對含有基線誤差的干涉對，基于高程信息和解纏相位，采用非線性最小二乘平差方法去除基線誤差(Rosen et al.，1996;Lu，2007)。對明顯含有與地形相關(guān)的大氣相位延遲(Atmospheric Phase Screen，APS)的干涉對，根據(jù)無形變區(qū)域的解纏相位和高程信息，建立大氣相位延遲模型，去除了與地形相關(guān)的大氣相位延遲干擾。

削弱各種噪聲，提取緩慢變形信號的一種常用的方法是干涉圖堆疊技術(shù)(Zebker et al.，1994，1997;Peltzer et al.，2001)。干涉圖堆疊技術(shù)通過對覆蓋同一地區(qū)的多個干涉圖進行加權(quán)平均(一般以每幅干涉圖的時間間隔為權(quán)重)，實現(xiàn)削弱每個干涉對所含的大氣延遲相位以及軌道殘余誤差，達到提高最終形變結(jié)果可靠性的目的。

為了揭示阿什庫勒火山群地區(qū)在于田地震后地表形變隨時間的演化特征，本文采用了短基線集算法求解地表形變時間序列。短基線集算法是由Berardino等(2002)提出來的，后來由Jung等(2008)進行了精化改進。短基線集算法以分布式散射體(DS，Distribute Scatterers)為研究對象，允許以相對短的空間基線和時間間隔對多景影像進行自由組合干涉，增加了時間采樣率并提高了干涉圖的相干性，削弱了空間失相干的影響，降低了外部DEM引入的誤差，從而達到提高形變結(jié)果精度的目的。由于大氣延遲相位在空間域上高度相關(guān)，而在時間域則是隨機的(Ferretti et al.，2000;Berardino et al.，2002;Jung et al.，2008)，短基線集算法通過時間域的高通濾波和空間域的低通濾波來去除大氣相位延遲。精化的短基線集算法在以下4個方面提高了時間序列形變結(jié)果的可靠性:1)選取高質(zhì)量的干涉圖(不含相位解纏誤差，且相干性好)估計外部DEM誤差和線性形變相位;2)對形變時間序列和誤差進行迭代計算;3)利用有限差分平滑方法削弱時間域的噪聲;4)對參考點本身存在的誤差進行了校正(Jung et al.，2008)。

2 InSAR形變結(jié)果分析

2.1 2008年于田地震前阿什庫勒火山群地區(qū)地表形變特征

通過逐一分析2008年于田地震前的所有干涉圖，沒有發(fā)現(xiàn)明顯的形變信息，選取18個相干性較好的干涉圖，采用干涉圖堆疊技術(shù)，計算得到了阿什庫勒火山群地區(qū)2003—2007年的平均形變速率圖(圖2)，其中黑色三角形表示各個小火山口的位置。可以看出，阿什庫勒火山群地區(qū)相干性較好，南部山區(qū)由于終年積雪導(dǎo)致失相干現(xiàn)象較嚴(yán)重，但溝谷地帶仍能保持較好的相干性。形變?yōu)檠乩走_視線方向(LOS，line of sight)，可以看出在2008年3月于田地震發(fā)生之前，阿什庫勒火山群地區(qū)沒有顯示出明顯的地表形變，表明該區(qū)在2003—2007年期間相對比較穩(wěn)定。

圖2 阿什庫勒火山群2003—2007年平均形變速率圖(LOS)Fig.2 Average deformation rate of Ashikule volcanic group during 2003 －2007(Line of Sight).黑色三角形與圖1相同;紅色圓點表示參考點

2.2 2008年于田地震后阿什庫勒火山群地區(qū)地表形變特征

通過逐一分析2008年于田地震后的所有干涉對，發(fā)現(xiàn)阿什庫勒火山群和于田地震發(fā)震地區(qū)均存在不同程度的形變信息。通過進一步對比，發(fā)現(xiàn)形變與時間呈明顯的線性關(guān)系。我們首先選取了Envisat 43個、ALOS PALSAR 6個相干性較好的干涉圖，采用干涉圖堆疊技術(shù)獲得了震后的平均形變速率圖(圖3)。可以看出，阿什庫勒火山群地區(qū)整體有沿雷達視線方向隆升的態(tài)勢，但量級不大，＜1cm/a。南部山區(qū)沒有明顯的形變趨勢，但是于田地震同震地表破裂帶(徐錫偉等，2011)北端地區(qū)呈現(xiàn)了明顯的沿雷達視線方向下沉的趨勢，且量級較大，最大達到約4cm/a(雷達視線方向)。

考慮到大氣相位延遲、外部DEM誤差，選擇短基線集算法求解形變時間序列，便于直觀分析和更好地揭示研究區(qū)域震后形變隨時間的演化特征。選取Envisat相干性較好的43幅干涉圖，ALOS PALSAR的6幅干涉圖，分別獲取了研究區(qū)域基于2種衛(wèi)星影像資料的形變時間序列特征。圖4給出了所使用的干涉對的時間與空間垂直基線，實線代表Envisat ASAR影像干涉對，虛線代表ALOS PALSAR影像干涉對。可以看出，2種衛(wèi)星影像資料各自生成了一個干涉圖子集。而ALOS PALSAR的時間被包含在Envisat ASAR干涉對的時間內(nèi)。圖5列出了根據(jù)Envisat影像獲取的阿什庫勒火山群地區(qū)震后形變時間序列，即所有影像成像時間相對于參考影像(Envisat為20080420，ALOS PALSAR為20080711)的累積形變。從整體上看，阿什庫勒火山群地區(qū)震后有沿雷達視線方向隆升的態(tài)勢。值得注意的是，于田地震同震地表破裂帶北端存在明顯的震后形變，表現(xiàn)為明顯的沿雷達視線方向下沉的形變趨勢，形變與時間呈較強的線性關(guān)系。

圖3 阿什庫勒火山群于田地震震后平均形變速率圖(LOS)Fig.3 Average deformation rate of Ashikule volcanic group after the Yutian earthquake(Line of Sight).

為了更加清晰直觀地表現(xiàn)相應(yīng)研究目標(biāo)的形變時間序列演化，以阿什山火山以及于田地震同震地表破裂帶北端為例，研究它們的形變演化特征。圖6a，b分別列出了相應(yīng)目標(biāo)的形變時間序列，結(jié)合圖3，5，可以得到如下結(jié)論:

(1)鄰近參考點地區(qū)在時間上比較穩(wěn)定，表明參考點選擇合理。

(2)阿什山火山在Envisat(降軌，即衛(wèi)星由N向S飛行)干涉結(jié)果上顯示沿雷達視線方向隆升(即朝向衛(wèi)星方向)，而在 ALOS PALSAR(升軌，即衛(wèi)星由S向N飛行)干涉結(jié)果上則顯示沿雷達視線方向下沉(即遠(yuǎn)離衛(wèi)星方向)。然而，InSAR僅能觀測到一維形變信息，即雷達視線方向形變是地表3方向(E、N、上)形變的共同貢獻。根據(jù)雷達入射角以及衛(wèi)星飛行方位角，可以計算得到雷達視線方向位移與地表3方向位移的幾何關(guān)系，本文中，Envisat為［0.3719，－0.080，0.922］，ALOS PALSAR 為［－0.616，－0.110，0.780］?？梢钥闯觯?種衛(wèi)星影像均對垂向位移最敏感，而對SN方向位移的敏感性最差。另外，注意到2種衛(wèi)星影像對EW方向位移的表現(xiàn)是相反的，因而可以推斷阿什山火山的形變以水平方向分量為主。通過對比圖5，6可以得出，以阿什山火山為代表的阿什庫勒火山群地區(qū)在震后有整體向E運動的趨勢，量級較小，累積最大約1cm。這種形變分布不屬于火山活動的形變特征，應(yīng)為區(qū)域塊體水平構(gòu)造運動的反映。

圖4 干涉對基線組合圖Fig.4 Temporal and perpendicular baselines of InSAR pairs.

圖5 阿什庫勒火山群于田地震震后形變時間序列(LOS)Fig.5 Deformation time series in Ashikule volcanic group after Yutian earthquake.

(3)于田地震同震地表破裂帶北端附近在2種衛(wèi)星影像資料上均表現(xiàn)出了較強的下沉形變態(tài)勢(沿雷達視線方向)，且與時間呈較強的相關(guān)性。考慮到圖5中于田地震發(fā)震斷層兩盤在震后的位移，沿垂直于斷層的走向做了一條剖面(位置見圖3)，以詳細(xì)分析其形變的空間特征。由升降軌道成像幾何關(guān)系可知，SN方向形變在雷達視線方向上最不敏感，所以假定SN方向位移為0，聯(lián)合2種衛(wèi)星的成像幾何關(guān)系(Wright et al.，2004)，求解得到于田地震同震地表破裂帶北端的E向、垂向的形變速率分量，如圖7所示。

圖6 幾個局部區(qū)域的形變時間序列(位置見圖3a)Fig.6 Deformation time series at certain areas(see Fig.3a for locations).a Envisat;b PALSAR

由圖7可以看出，于田地震地表破裂帶北端的震后形變特征為水平方向整體向E運動，最大約5cm/a;垂直方向下沉，且距離斷層越近下沉速率越大，最大約6cm/a。同時，形變速率分量與高程也無明顯的相關(guān)關(guān)系。我們認(rèn)為，此局部地表下沉不屬于火山活動的影響范圍:1)該區(qū)域距離火山群較遠(yuǎn);2)如果確屬巖漿運移引起的地表下沉，那么在未噴發(fā)的情況下，巖漿到達新的存儲地點后，會引起局部地表隆升。因而，阿什庫勒火山群在2008年于田地震后沒有表現(xiàn)出明顯的火山活動引起的地表形變異常。同震地表破裂帶北端的形變異常應(yīng)該是震后塊體應(yīng)力平衡調(diào)整的結(jié)果。

我們收集了阿什庫勒火山群地區(qū)1970—2011年的地震目錄(圖1)，雖然區(qū)域內(nèi)地震活動比較頻繁，但沒有明顯的空間分布規(guī)律?；鹕饺簝?nèi)部地震分布較少，與各個子火山也沒有空間的相關(guān)關(guān)系。這表明，研究區(qū)域的地震活動可能與阿什庫勒火山群無直接的聯(lián)系，應(yīng)該是區(qū)域內(nèi)部發(fā)育的斷層活動引發(fā)的構(gòu)造地震。因此，在過去的40余年間，阿什庫勒火山群地區(qū)沒有發(fā)生與火山活動有關(guān)的地震事件。

圖7 于田地震地表破裂帶北端剖面(圖3)的E向與垂向速度關(guān)系曲線Fig.7 Cross section of average deformation rates along the rupture caused by Yutian earthquake.

3 結(jié)論

(1)2008年于田地震之前，阿什庫勒火山群地區(qū)總體形變不明顯，表明阿什庫勒火山群在2003—2007年期間相對較穩(wěn)定。

(2)2008年于田地震之后，阿什庫勒火山群地區(qū)整體有E向運動的趨勢，到2010年7月，累積最大位移約1cm，為震后塊體應(yīng)力調(diào)整的結(jié)果，非火山活動引起。地表同震破裂帶北端在震后顯示出了較強的E向水平運動和垂向下沉運動的形變態(tài)勢，累積最大約6cm，與火山活動沒有必然的聯(lián)系。

(3)1970—2011年的地震目錄表明，阿什庫勒火山群在過去的40多年以來沒有發(fā)生火山地震活動事件。

陳學(xué)忠，蔣長勝，李燕娥.2008.2008年3月21日新疆于田7.3級地震［J］.國際地震動態(tài)，325(4):18—28.

CHEN Xue-zhong，JIANG Chang-sheng，LI Yan-e.2008.The March 21，2008 Yutian，Xinjiang MS7.3 earthquake［J］.Recent Development in World Seismology，325(4):18—28(in Chinese).

鄧萬明.1989.新疆烏魯克庫勒火山群地質(zhì)考察:“1951年1號火山噴發(fā)”質(zhì)疑［J］.自然資源學(xué)報，4(4):349—354.

DENG Wan-ming.1989.Geological investigation on the volcanic group in the Wulukule(lake)region，Xinjiang:Query the validity of“the eruption of No.1 volcano in 1951”［J］.Acta Natural Resources，4(4):349—354(in Chinese).

劉嘉麒，買買提依明.1990.西昆侖山第四紀(jì)火山的分布與K－Ar年齡［J］.中國科學(xué)(D輯)，(2):180—187.

LIU Jia-qi，YIMING Maimaiti.1990.The distribution of Quaternary volcanoes and its K － Ar isotopic age in western Kunlun Mts［J］.Science in China(Ser D)，(2):180—187(in Chinese).

劉若新.2000.中國的活火山［M］.北京:地震出版社.

LIU Ruo-xin.2000.Active Volcanoes in China［M］.Seismological Press，Beijing(in Chinese).

單新建，張國宏.2006.孕震區(qū)震前D－InSAR干涉形變場動態(tài)演化圖像分析［J］.地震地質(zhì)，28(3):441—446.

SHAN Xin-jian，ZHANG Guo-hong.2006.An analysis of dynamic evolution of preseismic interferometric deformation fields in seismic area［J］.Seismology and Geology，28(3):441—446(in Chinese).

孫建寶，梁芳，沈正康，等.2008.汶川MS8.0地震InSAR形變觀測及初步分析［J］.地震地質(zhì)，30(3):789—795.

SUN Jian-bao，LIANG Fang，SHEN Zheng-kang，et al.2011.InSAR deformation observation and preliminary analysis of the MS8 Wenchuan earthquake［J］.Seismology and Geology，30(3):789—795(in Chinese).

新疆日報.1951.新疆昆侖大阪西溝發(fā)現(xiàn)火山噴發(fā)［Z］.7月5日，第二版.

Xinjiang Daily.1951.Volcano Eruption Occurred in the Western Valley of Kunlun Mountains，Southern Xinjiang［Z］.May 7，page 2(in Chinese).

許建東，趙波，張柳毅，等.2011.新疆阿什庫勒火山群野外地質(zhì)科學(xué)考察［J］.地震地質(zhì)，33(3):499—505.doi:10.3969/j.issn.0253 －4967.2011.03.001.

XU Jian-dong，ZHAO Bo，ZHANG Liu-yi，et al.2011.Field geological exploration of Ashikule volcano group in western Kunlun Mountains［J］.Seismology and Geology，33(3):499—505(in Chinese).

徐錫偉，譚錫斌，吳國棟，等.2011.2008年于田MS7.3地震地表破裂帶特征及其構(gòu)造屬性討論［J］.地震地質(zhì)，33(2):462—471.doi:10.3969/j.issn.0253 －4967.2011.02.019.

XU Xi-wei，TAN Xi-bin，WU Guo-dong，et al.2011.Surface rupture features of the 2008 Yutian MS7.3 earthquake and its tectonic nature［J］.Seismology and Geology，33(2):462—471(in Chinese).

Berardino P，F(xiàn)ornaro G，Lanari R，et al.2002.A new algorithm for surface deformation monitoring based on small baseline differential SAR interferograms［J］.Institute of Electrical and Electronics Engineers Transactions on Geosciences and Remote Sensing，40(11):2375—2383.

Ferretti A，Prati C，Rocca F.2000.Nonlinear subsidence rate estimation using permanent scatterers in differential SAR interferometry［J］.IEEE Transactions on Geoscience and Remote Sensing，38(5):2202—2212.doi:10.1109/36.868878.

Jung H S，Lee C W，Park J W，et al.2008.Improvement of small baseline subset(SBAS)algorithm for measuring time-series surface deformations from differential SAR interferograms［J］.Korean Journal of Remote Sensing，24:165—177.

Liu J Q，Maimaiti Y M.1989.Distribution and ages of Ashikule volcanoes on the West Kunlun Mountains，west China［J］.Bulletin of Glacier Research，7:187—190.

Lu Z.2007.InSAR imaging of volcanic deformation over cloud-prone areas:Aleutian Islands［J］.Photogrammetric Engineering ＆ Remote Sensing，73(3):245—257.

Massonnet D，Rossi M，Carmona C，et al.1993.The displacement field of the Landers earthquake mapped by radar interferometry［J］.Nature，364:138—142.

Massonnet D，Briole P，Arnaud A.1995.Deflation of Mount Etna monitored by spaceborne radar interferometry［J］.Nature，375:567—570.

Peltzer G，CrampéF，Hensley S，et al.2001.Transient strain accumulation and fault interaction in the eastern California shear zone［J］.Geology，29:975—978.

Rosen P，Hensley S，Zebker H，et al.1996.Surface deformation and coherence measurements of Kilauea Volcano，Hawaii，from SIR － C radar interferometry［J］.Journal of Geophysical Rearch，101:23109—23125.

Wright T J，Parsons B E，Lu Z.2004.Toward mapping surface deformation in three dimensions using InSAR［J］.Geophysical Research Letter，31(1):L01607.

Zebker H A，Rosen P A，Goldstein R M，et al.1994.On the derivation of coseismic displacement fields using differential radar interferometry—the Landers earthquake［J］.Journal of Geophysical Research，99(B10):19617—19634.

Zebker H A，Rosen P A，Hensley S.1997.Atmospheric effects in interferometric synthetic aperture radar surface deformation and topographic maps［J］.Journal of Geophysical Research，102:7547—7563.