1、1利用 SBASInSAR 技术研究印尼 Agung 火山的形变特征与岩浆房参数摘要：基于 ALOS PALSAR 影像，利用小基线集合成孔径雷达干涉测量技术，提取了位于印度尼西亚巴厘岛的 Agung 火山 20072009 年的地表形变时间序列，并基于 Mogi 点源模型和竖直椭球体模型反演了岩浆房参数。结果表明：Agung 火山地区大气延迟相位干扰较严重，Agung 火山在 20072009 年发生了较明显的隆升形变，且与时间呈正相关。竖直椭球体模型能够更好地拟合 InSAR 形变场，岩浆房位于火山体下方约 5 km处。SBASInSAR 结果表明，应加强跟踪监测 Agung 火山的潜在喷
2、发危险性。关键词：SBAS-InSAR；形变时间序列；Agung 火山；岩浆房参数；大气延迟相位 中图分类号：P315727 文献标识码：A 文章编号：1000-0666（2013）03-0313-06 0 引言 Agung 火山位于印度尼西亚的巴厘岛东部，是该岛的最高峰，最高处海拔 3 142 m，构造上位于爪哇海沟以北约 300 km 处，属于印度澳大利亚板块向东南亚大陆下方俯冲的弧后火山。Agung 火山体呈锥状，顶部存有 520 m375 m 的漏斗状破火山口，无植被生长，火山体的侧面陡峭，植被较茂盛（Dilmy，1965） 。历史上，Agung 火山曾多次喷发，在过去的2300 年里
3、，有记录的喷发事件分别发生在 1808、1821（待考证） 、1843、19631964 年（Seach，2012） 。其中，19631964 年的布里尼式喷发夺去了至少 1 000 人的生命，爆炸性喷发所产生的巨大的火山灰云随后笼罩整个地球，造成了巨大的灾难（Frederic，1970） 。喷发历史表明，Agung 火山属于具有巨大喷发危险的当代活动火山，一旦喷发，将会给周围居民乃至全球气候带来巨大灾难。然而，由于当地社会经济条件落后等方面的原因，Agung 火山目前缺乏地表形变监测资料。 近年发展起来的 InSAR 技术因其自身的优势已在火山区地表形变监测方面得到了广泛应用（Massonn
4、et et al，1995；Hooper et al，2004，2007；Peltier et al，2010） 。常规 InSAR 技术容易受到时间、空间失相干的影响，Berardino 等（2002）提出了小基线算法（small baseline subsets，简称 SBAS） ，该算法仅选取短基线干涉对进行时间序列形变求解，削弱了空间失相干的影响，同时大大降低了地形因素引入的误差和大气延迟相位。Jung 等（2008）对小基线算法在 4 个方面进行了精化，提高了时间序列形变结果的可靠性：（1）选取高质量的干涉图（不含相位解缠误差，且相干性好）估计外部 DEM 误差和线性形变相位；（2）
5、对形变时间序列和误差进行迭代计算；（3）利用有限差分平滑方法削弱时间域的噪声；（4）对参考点本身存在的误差进行了校正。目前小基线算法已在火山区缓慢地表形变监测方面表现出了极大的潜力和优势（Lundgren et al，2004；Pepe et al，2008；Casu et al，2008；季灵运等，2011） 。 1InSAR 数据处理与结果 3Agung 火山位于热带地区，植被覆盖茂盛（图 1） ，本文选用了对植被具有较强穿透能力的日本 ALOS PALSAR 影像（波长 236 cm） ，获取时间从 2007 年 7 月至 2009 年 1 月，选取的影像提供了 15 a 的地表形变观测
6、资料。 11InSAR 数据处理 常规的 DInSAR 技术容易受到大气延迟相位的干扰。Zebker 等（1997）的研究表明，云的相对湿度在空间或时间域内变化 20%将会导致1014 cm 的形变误差。对于 Agung 火山区，大气相位延迟在某些干涉图中表现的比较明显。图 2 列出了 3 幅干涉图（已解缠） ，可以看出图 2a和图 2b 中的颜色分布相反，若差分干涉相位不包含误差，则两者反映了相反的形变态势，即在图 2a 干涉对的时间间隔内，火山区域发生了明显的下沉，最大下沉达 75 cm；而在图 2b 干涉对的时间间隔内，火山区域则发生了明显的隆升，最大隆升达 38 cm；图 2c 却未反
7、映出较大的相位差信息。2007 年 710 月火山并无喷发，如此短的时间内（92 d）形变态势发生转换，不符合火山区形变的物理解释。图 2a，b 中的相位分布与地形具有较强的相关性，即靠近火山体顶部，相位的绝对值较大，顶部和侧面的相位差达到 34 弧度，符合对流层水汽垂直分层效应对雷达信号延迟干扰的特征。由此可以得出，图 2a，b 中的较大相位差均由2007 年 8 月 24 日这景影像携带的大气延迟相位所导致。1 地震研究 136卷第 3 期季灵运等：利用 SBASInSAR 技术研究印尼 Agung 火山的形变特征与岩浆房参数图 1Agung 火山周围 SRTM DEM（右上角插图 显示了
8、巴厘岛及其周围岛屿的位置， 4虚线框表示干涉图的范围） Fig1Digital Elevation Map from Shuttle Radar Topography Mission surrounding Agung volcano （the illustration in the right top corner shows the location of Bali island（small rectangle） and its surrounding islands，the dashedline frame shows the range of interferogram）图 2 带有显
9、著大气延迟相位的解缠干涉图 （a）2007-07-092007-08-24；（b）2007-08-242007-10-09；（c）2007-07-092007-10-09 Fig2Unwrapping interferograms with atmospheric delayed phase 考虑到大气延迟相位的干扰，本文选择小基线算法求解时间序列形变。给定垂直基线小于 1 000 m，时间间隔小于 6 个月，自由组合干涉生成 33 幅干涉图。对于每个干涉图，多视因子给定 49，即每个像元约为（3030）m2。对含有基线误差的干涉图，基于已有 DEM，利用最小二乘方法进行了基线精化（Rosen
10、 et al，1996；Lu，2007） 。选取 12 幅高信噪比的干涉图，组成了一个子集（图 3） ，进行小基线算法解算。 12 基于 SBASInSAR 技术的 Agung 火山区形变时间序列利用前述精化的SBASInSAR 技术提取了 Agung 火山的形变时间序列形变场（图 4） 。从整体上看，从 2007 年 7 月 9 日到 2009 年 1 月 11 日，Agung 火山地区表现出了较强的地表隆升趋势，与时间呈明显的线性正相关，且形变从火山体侧部到顶部逐渐增大，552 d 内累积最大形变超过 12 cm。为了表现火5山不同部分的时间序列形变特征，图 5 给出了火山顶部（CO） 、
11、火山体侧部（FL）以及参考点附近地区（NR）的形变时间序列曲线。从图中可以看出，火山顶部和火山体侧部表现出了相似的隆升形变特征，其中火山顶部的累积形变约为 13 cm，火山体侧部的累积形变约为 5 cm。图3InSAR 干涉组合和基线图 Fig3InSAR interference combinations and baselines 2 时间序列累积形变场模拟与分析 火山区的地表形变一般由于其下方岩浆房或者岩浆系统的扰动等所引起的（Dvorak，Dzurisin，图 5 火山体不同部位的形变时间序列（火山体 不同部位见图 4a 中黑色三角） Fig5InSAR LOS deformation
12、 time series for different parts of volcano（the different parts of volcano marked by the black triangles in Fig4a）图 4 基于 SBAS-InSAR技术的 LOS 方向时间序列形变场（2007 年 7 月至 2009 年 1 月） （无颜色的地区表示雷达信号失相干） Fig4Time series deformation fields derived by SBAS-InSAR teohnology in LOS direction from Jul，2007 to Jan，200
13、9（decorrelation of radar signal are uncolored）1997） 。为了解释 20072009 年 Agung 火山区的隆升形变态势，研究其下方可能存在的岩浆房或岩浆系统的动力学机制，本文选取552 d 的累积形变场，分别利用埋置于均匀弹性半空间中的点源 Mogi 模6型（Mogi，1958）和竖直椭球体模型（Yang et al，1988）进行了模拟。考虑到可能残留的轨道误差，模拟过程中引入了两个方向的线性趋势项（Massonnet，Feigl，1998） 。采用下山单纯形算法和蒙特卡洛搜索算法估计最优参数及其误差（Press et al，1992） ，
14、以残差干涉图（原始干涉图与模拟干涉图之差）的均方根误差作为迭代结束的阈值。 图 6 列出了模拟结果，其中图 6a 为 552 d 的累积形变场，为了削弱噪声，将观测形变场（图 4h）进行了均值滤波，图 6b 和 d 分别为 Mogi点源模型和竖直椭球体模型模拟的结果，图 6c 和 e 为相应的残差图。从整体上看，两模型均可以很好地模拟原始观测形变场。表 1 列出了两模型模拟的 Agung 火山区岩浆房的几何参数，两模型模拟的岩浆房的中心位置较一致，均位于火山口的下方，但是竖直的椭球体模型模拟的结果相对较深。为了定量直观地比较两模型的模拟结果，沿火山锥体做了一条剖面（图 6a 中黑线，图 6f
15、为模拟结果） ，从图 6f 中可以看出，Mogi点源模型能够很好地拟合火山体侧面的形变场，而对火山体顶部的形变场拟合较差，而竖直的椭球体模型对火山体各个部位的形变场均能很好地模拟。由此表明，竖直的椭球体模型更加接近于真实的岩浆房几何形状。图 6 基于 Mogi、竖直椭球体模型的形变场模拟结果 （a）时间序列累积形变场（图 4h）的均值滤波图（黑线表示剖面位置） ；（b）Mogi 模型模拟的形变场；（c）Mogi 模型模拟的残差图；（d）竖直椭球体模型模拟的形变场；（e）竖直椭球体模型模拟的残差图；（f）模拟结果 （黑线表示地形，红线表示观测值，绿线和蓝线分别表示 Mogi、竖7直椭球体模型模拟
16、结果） Fig6Deformation field simulation results based on Mogi and vertical spheroid models （a）mean filter of time series accumulation deformation field in Fig4（black line represents the position of the profile） ； （b）deformation field simulated by the Mogi model；（c）residual simulated by Mogi model；（d）de
17、formation field simulated by the vertical spheroid models；（e）residual simulated by vertical spheroid models；（f）simulation results（the black line represents the terrain，the red line represents the observed value，the green and blue lines show the simulation results based on Mogi and vertical spheroid
18、models respectively）综合火山区的形变时间序列与基于点源 Mogi 模型、竖直的椭球体模型模拟实验，笔者认为 Agung 火山区在 20072009 年发生了比较明显的隆升形变，岩浆房中心位于火口下方约 5 km 处（表 1 列出了竖直的椭球体模型得到的深度约为平均海平面以下 2 km，火山体最高处约为 3 km） ，具备再次喷发的危险，需加强跟踪监测。表 1Agung 火山岩浆房几何参数（误差为 95%置信区间） 3 结论 本文基于 ALOS PALSAR 影像，利用 SBASInSAR 技术，提取了位于印度尼西亚巴厘岛的 Agung 火山 20072009 年的地表形变时
19、间序列，弥补8了该火山区一直以来无形变监测结果的现状。结果表明，Agung 火山在20072009 年发生了较明显的隆升形变，且与时间呈正相关，具有再次喷发的危险。基于 Mogi 点源模型以及竖直椭球体模型，利用 552 d 的累积形变场模拟得到了岩浆房的几何参数，结果表明，Agung 火山的岩浆房位于火山体下方约 5 km 处。 Agung 火山的 SBASInSAR 技术试验再次表明了 InSAR 技术在提取区域地表缓慢形变方面的巨大潜力和优势。InSAR 技术的应用，为 Agung 火山区的现今岩浆活动状况和地表形变态势分析提供了重要依据。我国境内分布有多座活动火山，受气候、当地社会经济
20、等条件的限制，还有一部分处于无监测或少监测资料的现状。本文基于 SBAS-InSAR 技术对 Agung火山开展的研究工作为查明和研究我国火山的现今活动特征起到重要的借鉴作用。 参考文献： 季灵运，王庆良，崔笃信，等 2011 利用 SBASDInSAR 技术提取腾冲火山区形变时间序列J.大地测量与地球动力学，31（4）：149-154 Berardino P，Fornaro G，Lanari R，et al2002A new algorithm for surface deformation monitoring based on small baseline differential SA
21、R interferogramsJ.Institute of Electrical and Electronics Engineers Transactions on Geosciences and Remote Sensing，40（11）：2 375-2 383 Casu F，Poland M P，Solaro G，et al2008Surface deformation dynamics of Mauna Loa and Kilauea volcanoes，Hawaii，revealed 9by InSAR time series analysisR.Sanfrancisco：Ameri
22、can Geophysical Union Dilmy A1965Pioneer plants found one year after the 1963 eruption of Agung in BaliJ.Pac Sci，19（4）：498-501 Dvorak J J，Dzurisin D1997Volcano geodesy：the search for magma reservoirs and the formation of eruptive ventsJ.Rev Geophys，35（3）：343-384 Frederic E V1970Atmospheric turbidity
23、 after the Agung eruption of 1963 and size distribution of the volcanic aerosolJ.JGR，75（27）：5 185-5 194 Hooper A，Segall P，Zebker H2007Persistent scatterer InSAR for crustal deformation analysis，with application to Volcano Alcedo，GalapagosJ.JGR，112，B07407，doi：101029/2006JB004763，2007 Hooper A，Zebker
24、H，Segall P，et al2004A new method for measuring deformation on volcanoes and other natural terrains using InSAR persistent scatterersJ.Geophys Res Lett，31（23） ，L23611，doi：101029/2004GL021737，2004 Jung H S，Lee C W，Park J W，et al2008Improvement of small baseline subset（SBAS）algorithm for measuring time
25、series surface deformations from differential SAR interferogramsJ.Korean Journal of Remote Sensing，24（2） ，10165-177 Lu Z2007InSAR imaging of volcanic deformation over cloudprone areas-Aleutian islandsJ.Photogrammetric Engineering & Remote Sensing，73（3）：245-257 Lundgren P，Casu F，Manzo M，et al2004Grav
26、ity and magma induced spreading of Mount Etna volcano revealed by satellite radar interferometryJ.Geophys Res Lett，31（4）：L04602 Massonnet D，Briole P，Arnaud A1995Deflation of Mount Etna monitored by spaceborne radar interferometryJ.Nature，375：567-570 Massonnet D，Feigl K1998Radar interferometry and it
27、s application to changes in the Earths surfaceJ.Rev Geophys，36（4）：441-500 Mogi K1958Relations between the eruptions of various volcanoes and the deformations of the ground surface around themJ.Bull Earthquake Res Inst Univ Tokyo，36：99-134 Peltier A，Bianchi M，Kaminski E，et al2010PSInSAR as a new tool to monitor preeruptive volcano ground deformation：Validation using GPS measurements on Piton de la FournaiseJ.Geophys Res Lett，37（12） ，L12301，doi：101029/2010GL043846 Pepe A，Manzo M，Casu F，et al2008Surface deformation of active volcanic areas retrieved with the SBASDInSAR
