利用 SBAS-D InSAR技术提取腾冲火山区形变时间序列＊

季灵运 王庆良 崔笃信 胡亚轩 郝 明 李煜航 秦姗兰

(中国地震局第二监测中心,西安 710054)

利用 SBAS-D InSAR技术提取腾冲火山区形变时间序列＊

季灵运 王庆良 崔笃信 胡亚轩 郝 明 李煜航 秦姗兰

(中国地震局第二监测中心,西安 710054)

基于 6景 JERS-1 L波段 SAR影像,利用小基线集-合成孔径雷达干涉测量技术,通过线性形变相位、非线性形变相位、大气延迟相位以及地形残差相位的分离,提取了腾冲火山地区 1995—1997年间地表形变时间序列(雷达视线向),与 2003—2004年的 GPS观测结果对比表明,SBAS-D InSAR技术提取地壳形变的精度可达亚厘米级。时间序列形变显示胆札-高田断裂两侧形变差异性显著,可能与其下方存在的地壳岩浆囊的活动有关。打鹰山地区地表形变揭示其下方可能存在隐伏断裂。

SBAS-D InSAR;JERS-1;形变时间序列;大气延迟相位;腾冲火山区

1 引言

腾冲火山区南北长约 90 km、东西宽约 50 km,共有 68座新生代火山[1],其中全新世有过活动的有3座(黑空山、马鞍山、打鹰山)[2]。前人研究结果表明,腾冲火山区现今仍有较强的活动性,具有再次喷发的危险[3-10]。对腾冲火山的形变监测,云南省地震局分别于 1997年和 2003年共埋设了 94个水准点,并进行了多期观测,并于 2002、2003、2004年进行了 GPS组网观测,李成波等[6]基于 GPS水平形变资料计算了地壳面膨胀率和面应变,认为压力源位于腾冲西南方向。李春光等[7]计算了 1998—1999年的垂直形变,发现打鹰山附近形成了一个明显的隆起区,而其南面则呈下降趋势。胡亚轩[8,9]、施行觉[10]等对 1998、1999、2000 3期水准数据进行了分析,结果显示该区岩浆活动比较活跃,且压力源中心位置不定,在破火山口周围迁移,他们还利用Mogi模型对岩浆源参数进行了反演,结果表明火山活动受断层控制。由前人的形变研究成果不难得出,腾冲火山区近年来存在形变异常。所以,跟踪腾冲火山区的形变演化,对了解岩浆源的活动特征具有重要的意义。然而,常规的 GPS、水准测量由于存在作业周期长、仅能获取点的形变信息等缺陷,很难全局把握火山的形变活动规律,近年发展起来的 InSAR技术可以实现全天候、高分辨率、大空间尺度的连续观测,D InSAR技术已在地表形变监测方面得到了广泛应用[11-13],SBAS-D InSAR技术[14]克服了常规D InSAR监测中的时空失相干因素影响,在城市缓慢地表形变[15]、断层蠕动[16]、火山区域缓慢地表形变[17-20]监测等方面表现出了极大的潜力和优势。所以本文尝试利用 SBAS-D InSAR技术提取腾冲火山区形变时间序列,探讨腾冲火山活动的状态。

2 SBAS-D InSAR技术原理与数据处理流程

SBAS算法中,差分干涉相位定义为[14]：

其中 x和 r为像元坐标,λ为雷达波长,Δd为雷达视线方向地表形变,B⊥为垂直基线,θ为 SAR视角, Δz为地形残差,Δφatm为大气延迟相位,Δφn为其他噪声相位。显然,SBAS算法将差分干涉相位分为地表形变相位、地形残差相位、大气延迟相位以及其他噪声相位 4部分。其数据处理流程如图 1所示。

图1 SBAS-D InSAR数据处理流程Fig.1 Flow chart of SBAS-D InSAR processing

文献[14-17]的研究表明,SBAS算法的精度可达毫米级,与 GPS观测结果也吻合得较好。SBAS算法由于在影像自由组合干涉时限制了时间基线和空间基线,保证了每幅干涉图的高相干性,又利用奇异值分解的方法将多个干涉图子集联合进行最小二乘求解,增加了时间采样。SBAS算法将外部DEM误差产生的相位分离出来,降低了使用外部DEM所引入的误差。充分考虑到大气延迟相位对形变结果的影响,SBAS算法通过时空域滤波的方法削弱了大气延迟相位。相对于永久散射体技术,小基线集技术获取到的形变序列在空间上更为连续,从而可以应用于监测地壳长时间缓慢变形。

3 SBAS-D InSAR技术在腾冲火山区的应用

3.1 SAR影像的选取

腾冲火山区地貌复杂,所处纬度偏低,植被覆盖茂密,很容易造成 SAR影像失相干,所以本文选用L波段 JERS-1 SAR数据,L波段数据能够较好地克服植被覆盖所造成的失相干问题。数据处理结果表明,在腾冲火山区,时间基线仍不超过一年才具有较好的相干性。6景 JERS-1 SAR影像,获取时间分别为 1995年 4月 10日、1995年 8月 20日、1995年 11月 16日、1996年 3月 27日、1996年 8月 6日、1997年 4月 27日,给定时间基线为 365天,空间垂直基线为 2 km作为限制条件,进行自由干涉组合,得到6个可用的干涉对,时间跨越了 1995—1997年。可用干涉对的时间基线和空间基线关系如图 2所示。

图2 干涉对基线组合Fig.2 Baseline combination of all the interferomeric pairs

3.2 常规 D InSAR技术处理

对每个符合时间、空间基线限制条件的干涉对均利用两轨法 D InSAR技术处理,处理平台为瑞士GAMMA遥感公司开发的GAMMA软件。处理步骤主要包括：原始 SAR信号成像、图像配准与干涉、外部DEM模拟、基线估计与精化、相干系数计算、干涉图自适应滤波、相位解缠。外部地形数据使用美国NASA(美国宇航局)公布的 SRT M DEM。另外,为了削弱噪声,增加相干性,同时为了检测相对大尺度的形变信息,对干涉图进行了多视处理,多视因子为9：18。相位解缠采用最小费用网络流算法,为了与GPS观测资料进行比对,采用腾冲县城附近为解缠起始区域。

3.3 SBAS技术处理

对所有符合条件的 6个差分干涉对的解缠相位进行小基线集技术处理。首先,通过差分干涉相位估计得到研究区域的线性形变相位,并转换为 LOS向形变速率,其中 3座全新世火山的线性形变速率分别为黑空山 -2.8 mm/a,打鹰山 -6.6 mm/a,马鞍山 -0.82 mm/a;同时得到了残差地形相位,发现其值均为 -0.1～0.1 rad,换算为 LOS向的形变约为 -1.4～1.4 mm(对于 JERS-1)。然后从原始差分干涉相位中减去线性形变相位和残差地形相位,再对残余相位进行解缠,利用奇异值分解方法求解非线性形变相位,最后联合线性形变相位和非线性形变相位得到时间序列累积形变相位。对于大气延迟相位的估计,首先在时间域上进行均值滤波处理,然后在空间域上采用小波分解的方法,分离出高频成分,从而实现了大气延迟相位的分离。图 3是1995年 8月 20日的大气延迟相位图,从图 3可以看出,最大达 1.5 rad。

图3 大气延迟相位(1995-08-20)Fig.3 I mage on the atmospheric phase delay(1995-08-20)

3.4 SBAS形变结果的质量评价

图4为LOS向形变时间序列。根据雷达成像的几何关系,地面的上、东、北 3方向形变对 LOS向形变的贡献可表示为：

其中,θ为雷达脉冲入射角,φ为卫星轨道方位角。

本文选用的 GPS站点的形变速率值基本含盖了 InSAR数据所覆盖的范围。将 GPS观测的 3方向速率分量带入式 (2),得到其在 LOS向的形变速率,与 InSAR形变速率值的比较如表 1所示。

表1 GPS与 I nSAR获得的LOS向形变速率值的比较Tab.1 Comparison between deformation velocities(LOS) measured by GPS and I nSAR

由表1可以看出,虽然 GPS观测时间与 SAR成像时间不同,得到的速率不同,但是两者的一致性很好,最大差值约 1 cm/a,差值的标准差为 0.65 cm/ a。对比结果表明 SBAS技术对地表形变的监测精度可达亚厘米级。

3.5 腾冲火山区形变时间序列分析

图4直观地表现了腾冲火山区 1995—1997年的地表形变 (LOS向)随时间的演化过程。从全局来看,腾冲火山区发生了较大范围的区域性地表变形,大部分地区以下降为主要形变特征,其中下降相对较大的地方出现在黑空山西南以及腾冲县城以南地区,累积最大约 8 cm;研究区的西南部在 1995—1996年发生了隆升,最大达 4 cm,1996年后又回落;从局部来看,个别断裂在形变场上有明显显示,反映出断裂两侧有明显的差异活动,例如西南部的北西向胆札-高田断裂以及中部的南北向固东-腾冲断裂。图 5给出了垂直于胆札-高田断裂走向的时间序列形变剖面曲线 (图 4A-A’剖面),可以看出,断裂东北侧远离断层逐渐下沉,而西南侧远离断层则表现为隆升,断裂两侧存在明显的形变梯度,是断裂活动的反映。这一结果与李春光等[7]计算得到的 1998—1999年形变特征一致,与野外近场数字化台站记录结果所显示的火山地震活动主要发生在打鹰山-马鞍山一带以及腾冲县城西南部地区[1]也是一致的,表明该断裂现今仍在活动。值得注意的是,胆札-高田断裂的西南侧在 1996年后发生了约 4 cm的不规则形变变化,体现了该断裂的较强活动性[4]。

从图 4还可以定性地得出：1995—1997年,相对于腾冲县城,打鹰山的东北和西南两侧表现为明显的差异性形变特征,马鞍山地区形变不明显,黑空山一直处于下降状态。图 6给出了穿过打鹰山的形变时间序列剖面曲线(图 4B-B’剖面),从图 6可以看出,打鹰山的东北侧表现为上升,而西南侧表现为下降,累积最大差异达 10 cm,这种巨大的差异性运动特征预示着打鹰山下方可能存在北西向隐伏断裂,且活动性较强。这一结果与 1998—2000年的水准测量结果[8]一致,1998—1999年与 1999—2000年打鹰山的南北两侧的垂直形变差异性非常明显,推测压力源的位置不是固定不变的,沿破火口周围迁移。表 2列出了马鞍山和黑空山两座火山的时间序列形变量值,可以看出,两地区均表现为缓慢下降,但马鞍山地区下降幅度相对较小,累积不到 1 cm,表明其在 1995—1997年活动性不显著。黑空山附近表现为连续性下沉,累积最大超过 3 cm,表明其在 1995—1997年活动性不强。

图4 地表形变(LOS向)时间序列(1995—1997年)Fig.4 Ti me seriers of defor mation on LOS(1995-1997)

图5 胆札-高田断裂形变时间序列剖面Fig.5 Cross section on time series of defor mation of Danzha-Gaotian fault

图6 打鹰山形变时间序列剖面Fig.6 Cross section on time series of deformation of Dayingshan volcano area

表2 两座全新世活动火山的形变时间序列Tab.2 Ti me series of deformation of two Holocene volcanoes

综合火山区的整体形变时间序列和 3座全新世活动火山的形变时间序列演化特征分析,可以得出：腾冲火山区现今地壳运动比较活跃,形变特征复杂多样;打鹰山地区活动性比较显著,下方可能存在隐伏断裂;胆札-高田断裂两侧形变相对较大,差异性显著,地球物理资料揭示其下方存在活动性较强的岩浆囊可能是导致其差异性形变的主要因素。

4 结论

1)相对于常规 D InSAR技术,SBAS-D InSAR技术分离了大气延迟相位以及地形残差相位等干扰因素,从而能够获得高精度的地表形变时间序列 (LOS向),腾冲火山区的试验结果表明 SBAS-D InSAR技术能够以亚厘米级的精度检测区域地表形变。

2)胆札-高田断裂两侧形变差异性显著,可能与其下方存在的地壳岩浆囊的活动有关。

3)打鹰山地区地表形变揭示其下方可能存在隐伏断裂。

致谢 感谢路中博士提供 JERS-1 SAR影像!

1 皇甫岗,姜朝松.腾冲火山研究[M].昆明：云南科技出版社,2000.(Huangfu Gang and Jiang Chaosong.Study on Tengchong volcanic activity[M].Kunming：Yunnan Press of Science and Technology,2000)

2 李大明,李齐,陈文寄.腾冲火山区上新世以来的火山活动[J].岩石学报,2000,16(3)：362-370.(Li Daming, LiQi and ChenWenji.Volcanic activities in the Tengchong volcano area since Pliocene[J].Acta Petrologica Sinica, 2000,16(3)：362-370)

3 叶建庆,等.腾冲火山区地震群的活动特征[J].地震地质,2003,25(增)：128-137. (Ye Jian qing,et al. Characteristics of earthquake cluster activity in Tengchong volcanic area[J].Seismology and Geology,2003,25(supp l.)：128-137)

4 赵慈平,冉华,陈坤华.由相对地热梯度推断的腾冲火山区现存岩浆囊 [J].岩石学报,2006,22(6)：1 517-1 528.(Zhao Ciping,Ran Hua and Chen Kunhua.Presentday magma chambers in Tengchong volcano area inferred from relative geother mal gradient[J].Acta Petrologica Sinica,2006,22(6)：1 517-1 528)

5 楼海,等.云南腾冲火山区上部地壳三维地震速度层析成像[J].地震学报,2002,24(3)：243-251.(Lou Hai, et al.Three dimensional seis mic velocity tomography of the upper crust in Tengchong volcanic area,Yunnan province [J].Acta Seismologica Sinica,2002,24(3)：243-251)

6 李成波,等.腾冲火山区的 GPS形变特征[J].地球物理学进展,2007,22(3)：765-770.(Li Chengbo,et al.Research on the character of the GPS defor mation in the Tengchong volcano area[J].Progress in Geophysics,2007,22 (3)：765-770)

7 李春光,等.腾冲火山区的形变特征[J].地震研究,2000, 23(2)：165-171.(Li Chunguang,et al.Deformation feature in Tengchong volcano areas[J].Journal of Seis mological Research,2000,23(2)：165-171)

8 胡亚轩,等.腾冲火山区地表垂直形变分析[J].大地测量与地球动力学,2003,(2)：37-41.(Hu Yaxuan,et al. Analysis of vertical deformation in Tengchong volcanic area [J].Journal of Geodesy and Geodynamics,2003,(2)：37 -41)

9 胡亚轩,王庆良,赵慈平.腾冲火山区形变分析[J].国际地震动态,2008,352(4)：42-47.(Hu Yaxuan,Wang Qingliang and Zhao Ciping.Analysis on the deformation in Tengchong volcanic region[J]. Recent Developments in World Seismology,2008,352(4)：42-47)

10 施行觉,等.以垂直形变资料反演腾冲火山区岩浆活动性的初步研究 [J].地震研究,2005,28(3)：256-261. (Shi Xingjue,et al.Vertical deformation inversion ofmagma activity in Tengchong volcanic area[J].Journalof Seismological Research,2005,28(3)：256-261)

11 路旭,等.用 I NSAR作地面沉降监测的试验研究[J].大地测量与地球动力学,2002,(4)：66-70.(Lu Xu,et al.Experiment on land subsidence monitoring by InSAR [J].Journal of Geodesy and Geodynamics,2002,22(4)：66-70)

12 乔学军,郭利民.新疆伽师强震区的 InSAR观测研究[J].大地测量与地球动力学,2007,(1)：7-13.(Qiao Xuejun and Guo Limin.Study on Jiashi strong earthquake s warm area by InSAR[J].Journal of Geodesy and Geodynamics,2007,(1)：7-13)

13 乔学军,等.当雄Ms6.6地震的 InSAR观测及断层位错反演[J].大地测量与地球动力学,2009,(6)：1-7. (Qiao Xuejun,et al. Study on dislocation inversion of Ms6.6 Dangxiong earthquake as constrained by InSAR measurement[J].Journal of Geodesy and Geodynamics, 2009,(6)：1-7)

14 Berardino P,et al.A new algorithm for surface deformation monitoring based on small baseline differential SAR interferograms[J].IEEE Transact. Geoscience and Remote Sensing,2002,40(11)：2 375-2 383.

15 Antonio Pepe,et al.On the generation of ERS/ENV ISAT D InSAR time-series via the SBAS technique[J].IEEE Geoscience and Remote sensing letters,2005,2(3)：265-269.

16 Lanaria R,et al.Application of the SBAS-D InSAR technique to fault creep：A case study of the Hayward fault, California[J].Remote Sensing of Environment,2007,109 (1)：20-28.

17 Pepe A,et al.Surface deformation of active volcanic areas retrieved with the SBAS-D InSAR technique：an overview [J].Annals of Geophysics,2008,51(1)：247-263.

18 Chang-Wook Lee,et al.Surface defor mation of Augustine Volcano,1992-2005,from multiple-interferogram processing using a refined SmallBaseline Subset(SBAS)Interferometric Synthetic Aperture Radar(InSAR)approach [A].The 2006 eruption of Augustine Volcano,Alaska：U.S.Geological Survey Professional Paper 1769[C].453 -465.

19 Lundgren P,et al.Gravity and magma induced spreading ofMount Etna volcano revealed by satellite radar interferometry[J].Geophysical Research Letters,2004,31(4)：L04602.

20 Casu F,et al.Surface defor mation dynamics ofMauna Loa and Kilauea volcanoes,Hawaii,revealed by InSAR time series analysis[R].American Geophysical Union,Fall Meeting,2008,V11B-2023.

TI M E SERIES OF DEFORMATI ON IN TENGCHONG VOLCANIC AREA EXTRACTED BY SBAS-D InSAR

JiLingyun,WangQingliang,CuiDuxin,Hu Yaxuan,HaoMing,Li Yuhang and Qin Shanlan
(Second CrustM onitoring and Application Center,CEA,X i’an 710054)

On the basisof the 6 JERS-1 SAR images,the time seriesof deformation in Tengchong volcanic area were obtained by SBAS-D InSAR technique through separating the linear deformation phase,the nonlinear deformation phase,the at mospheric phase screen and theDEM residual errorphase.Comparedwith the GPSmeasurement results,SBAS-D InSAR technique can detect the crustal defor mation at the precision of sub-cm level.The defor mation of the Danzha-Gaotian fault is obviously different bet ween the two sides,which may relate the chamber below. The defor mation pattern of the Dayingshan area indicates that there is probably an buried fault underground.

SBAS-D InSAR;JERS-1;time seriesof deformation;atmospheric phase delay;Tengchong volcanic area

1671-5942(2011)04-0149-06

2011-03-23

国家自然科学基金(40974062);地震行业科研重点专项(200908029)

季灵运,男,1982年生,博士生,主要研究方向为 InSAR技术在地震、火山形变中的应用.E-mail：dinsar010@163.com

P225.1

A