基于InSAR同震形变观测反演2010年新西兰南岛Mw7.1Darfield地震同震破裂分布

罗旭巍，孙建宝，沈正康，＊，Yo Fukushima

1 北京大学地球与空间科学学院，北京 100871

2 中国地震局地质研究所地震动力学国家重点实验室，北京 100029

3 Disaster Prevention Research Institute，Kyoto University，Uji，Kyoto，Japan

1 引 言

2010年9月4日新西兰南岛Canterbury平原发生了Mw7.1级地震，震源深度约为10km，震中位于Darfield东南10km处，地震造成2人受伤［1］.但它的余震，2011年2月21日Christchurch Mw6.3级地震造成了180余人死亡［2］.

新西兰岛屿横跨太平洋板块和澳洲板块（图1），地壳运动十分活跃.新西兰南岛和北岛之间，太平洋板块以约43mm／a的速度向澳洲板块俯冲；而在南岛的西南端，澳洲板块以35mm／a的速度插入太平洋板块下方［3］（图1）.由于上述两大板块汇聚在南岛附近且存在斜向运动，新西兰南岛南端形成了一条东西走向的阿尔卑斯走滑断裂带，GPS测量数据显示其存在约40mm／a的右旋走滑速率［4－5］.除了阿尔卑斯走滑断裂带吸收的震间滑动外，还有约10mm／a的跨边界带缩短分布在阿尔卑斯山走滑断裂带东南部的山地内［5］.因而新西兰成为世界上地震频繁的国家之一（每年超过15000个ML≥2的地震）［7］，历史上发生的大震包括1855年的Mw8.2 Wairarapa地震和1931年的 Mw7.8Hawke′s Bay地震.

2010年的Mw7.1Darfield地震发生在地震活动相对寂静的地区，发震断层震前并未为人所知.地震之后一些地质学家和地球物理学家迅速进行了相关研究，调查了断裂的地表破裂轨迹［8］，给出了地震的震源机制解［7，9］.同时，利用地震波形以及GPS和InSAR位移场反演了该地震的同震破裂分布［1，4，10］.本文中，我们利用一套时空基线更好的InSAR卫星观测数据，结合地质调查，通过最大后验概率反演方法同时确定断层的几何参数，光滑参数，不同InSAR数据集的权重，以及这次地震的同震滑动分布.最后，将形变反演结果与野外观测结果进行系统比较，并研究该地震对周围地区应力场的加载作用.

图1 新西兰地区构造形变主要特征（箭头表示板块相对运动方向）Fig.1 Main features of plate tectonics in New Zealand（Arrows show relative plate motion directions）

2 InSAR数据

本研究采用了日本空间探测局（JAXA）的ALOS卫星上PALSAR传感器获得的雷达数据，在覆盖地震发震时间范围内，选取了三对垂直基线较短且完整覆盖地震形变区的数据源（表1），包括两对升轨数据和一对降轨数据.

表1 卫星数据参数Table 1 Parameters of Raw Data

选择的两条升轨数据都满足垂直基线较短，震后时间较短，主辅图间隔时间较短三个要素，有助于保持震前震后雷达信号的干涉相干性，且同震形变场中混入的震后形变较小，从而获得高信噪比和高精度的同震形变场.与前人所采用的雷达数据不同的是，本研究所用的336轨道数据垂直基线为302m，成像时间间隔为6个月，而前人采用的数据雷达数据成像时间间隔达到了2年6个月［10］.另外，在数据处理过程中，我们采用了二阶雷达图像配准算法，使得雷达干涉信号的相干性有所增强.特别是本研究所用的336轨道数据，只有在使用二阶配准算法的情况下才能得到完整的形变场，而前人使用的一阶配准模型无法得到有效的形变信号.因此，本研究中的数据选择和处理方面比前人有明显改善，而这些显著改善对于更好约束断层的几何特征和滑动分布非常关键［11－12］.

本研究使用美国JPL实验室开发的Repeat Orbit Interferometry Package［13］（ROI－PAC）软 件包以及Gamma InSAR处理软件包，处理获得同震形变干涉条纹（图2），同时用雷达幅值信息互相关的方法得到雷达斜距向和方位向的水平位移数据（图3（c，d））.InSAR相位通过全局最优算法解缠，获得了覆盖地震形变区的包含三分量位移信息的雷达视线向图2（Line－of－Sight，LOS）位移［14］.

差分干涉雷达相位信息中包括若干分量：φ＝φdef＋Δφε＋φatm＋Δφorb＋φn.其中φdef是形变相位，为地表形变在LOS向引起的相位变化；φatm是大气噪音，由大气延迟作用引起的相位误差；Δφorb是由轨道误差引起的剩余相位；Δφε是高程剩余相位，由高程模型DEM中的误差引起的相位，φn是系统噪音，为处理过程中引起的误差相位.φdef形变相位是有用相位，其它的相位的贡献都被视为形变测量中的误差源.高程剩余相位Δφε与垂直基线和高程模型有关.本文采用精度为3″的SRTM高程模型［15］，并将其插值到1″用于消除干涉相位中的地形贡献，其引入的误差相对较小，对同震形变场的影响可以忽略.系统噪音φn是信号噪音、插值误差和配准误差等引入的相位误差［13］，在相干性较好的数据中，系统噪音对干涉相位的贡献非常小，一般可以忽略不计.因此本文处理的数据中，大气误差φatm和 轨道误差φorb是主要的误差来源.

观测信号中的形变相位φdef，大气延迟相位φatm，轨道误差相位Δφorb都具有空间相关性，但是各自的空间相关尺度不一样，可以利用其空间相关尺度的不同将这几种信号区分开［16－17］.轨道误差信号具有长波特征，其尺度大于区域性形变特征尺度.同震形变主要集中在近场，选择一个足够远的区域使得同震形变近似为零，利用最小二乘方法将此远场区域误差用一个线性或者二阶的剩余轨道面进行拟合，得到剩余轨道相位，并从干涉形变场中消除，从而消除轨道误差的影响.另外，大气噪声中有一部分与地形起伏线性相关［13］，因此可以将上述改正的远场信号再用最小二乘方法拟合其与地形的线性相关性，将得到的大气分量从总的形变场中进一步消除掉.

3 同震破裂模型反演

对于线性、非线性混合反演的问题，目前有两种比较常见的概率反演算法，分别是完全贝叶斯概率反演方法（Fully Bayesian Inversion，以下简称FBI）［18］和线性－非线性混合蒙特卡洛反演方法（Mixed linear－nonlinear Monte－Carlo method，以下简称MCI）［19］.前者将所有参数，包括断层几何参数和滑动参数，数据相对权重和光滑因子等都作为非线性参数反演；而后者将滑动参数作为线性参数（基于最小二乘方法）、其它参数作为非线性参数（基于蒙特卡洛方法）反演.FBI方法对所有的未知量都采用蒙特卡洛随机采样进行试探，以期获得最终的收敛解，其求解大量未知量的地球物理问题时，存在耗时长、收敛慢的问题.也正是因为这一不足，Fukuda和Johnson又将FBI方法改良为了 MCI.在 MCI中，未知量被分成线性未知量和非线性未知量两类.MCI算法大大降低了未知力量的个数，需要通过蒙特卡洛方法反演的未知参量的个数，因而能够有效减少迭代次数.两种方法使用相同的参数采样方法，即马尔科夫链蒙特卡洛（Markov Chain Monte Carlo，MCMC）方法.它们有一个共同的问题，即参数采样步长选择不当可能使得反演不收敛.如果采样步长设定过大，反演参数采样可能趋向于小概率区域，得到一个滑动分布过于粗糙但数据拟合很好的解.如果采样步长过小，每一步的采样缺乏效率，需要大量的尝试才能得到收敛解.实际操作中，往往需要多次尝试才能最终确定合理步长.对于存在较大误差的InSAR观测数据，很难在较短的时间内获得合理的反演结果.鉴于此，我们最终采用了最大后验概率反演方法（three－step Maximum－A－Posteriori Probability inversion）.

3.1 最大后验概率反演简介

形变场的最大后验概率反演［21］过程分为三步.

第一步，采用自适应模拟退火法（Adaptive Simulated Annealing，以下简称 ASA）［22］搜索形变模型的最大后验概率，将满足最大后验概率的断层滑动模型作为形变场反演结果，反演中模型的后验概率公式沿用了Fukuda和Johnson［18］在MCI方法中的定义.每一步采样中，断层倾角、数据权重和平滑因子作为非线性量随机采样，而断层错动量则通过最小二乘法求得.不断重复这一过程，直到得到稳定的收敛结果.反演第一步是全局最优化反演，其效率远远高于FBI和MCI中的MCMC采样方法.ASA反演中自动确定参数步长和各个参数在反演中的敏感性，因此该步反演能够很快收敛到最大后验概率模型附近，同时局部后验概率极大值和孤立的参数采样区域都被逐一测试和跳过.

第二步，将第一步得到的模型作为初始值，采用类似于FBI的方法对各个未知量进行了带约束的非线性反演，将断层几何参数和滑动参数均作为非线性参数反演.第一步中为了反演快速收敛到合理的区域，断层滑动参数在反演中用最小二乘方法一次反演得到，但最小二乘法无法对滑动参数施加约束，可能会在一些数据约束不够强的区域（如断层较深部位），或者观测数据局部噪声过大时，出现不合理的解，如局部的滑动方向与周围区域的滑动方向相反，需要进一步的修正；因此第二步中将断层上每一个矩形元的滑动参数均作为非线性参数在给定的范围内逐一单独采样反演，便于对滑动参数的值进行进一步的修正，反演进一步趋于系统解的最大后验概率.

第三步，将第二步中得到的断层几何模型固定下来，只反演错动量及数据权重、平滑因子等非线性参数.该步反演可以避免断层几何参数和滑动参数之间的相互影响（trade－off），同时进一步使得反演趋于最大后验概率模型.

3.2 断层分段依据和断层反演模型设置

对于主要发震断层Greendale断层，前人做了很多卓有成效的调查工作，Elliot等［10，23］已经准确地描述出了本次地震在地表的破裂轨迹（图3a）.本研究中的两条升轨干涉图信号具有很高的信噪比，形变干涉条纹连续清晰.337轨道的干涉条纹图中（图2b）可以看到，形变条纹中存在一条近东西走向的分界线，分界线两侧具有不同的条纹色彩变化模式，显示它们具有不同的LOS移动方向.该分界线与野外观测发现的Greendale断层一致.Greendale断层以南的条纹呈现以断层为界的椭圆形，且条纹非常连续，表明这个区域的形变主要由Greendale断层的滑动产生.断层北侧，可以明显看到形变干涉条纹的不规则扭曲，表明这个区域内可能存在多条断层，并产生了破裂.

从337轨道的信号相干性（图3b）图上，可以看到6个低相干性条带，其中1—4是贯通的，与地质调查的Greendale断层基本一致（图3a），5号低相干带位于Greendale断层的西南端，可以解释在337、631（图2b，图2c）两条轨道干涉图中在（171.9°E，43.57°S）附近出现的单独闭合条纹带.6号低相干带接近仪器震中，从336、337两幅干涉图（图2a，图2b）中可以明显看到几条与周围颜色反向的条纹，说明此处存在一条较小断层.以上几条断层在337轨道斜距向位移场中也能够看到（图3c）.从337轨道方位方向上（图3d）还能看到一条位移突变分界线，相同位置上336、337的干涉图条纹向南弯曲，显示此处存在一条微弱滑动的断层，在对应位置添加了一条断层单元7.根据余震分布［24］，将断层两端稍加延伸，使得整个断层模型东西总长约45km.破裂主体的Greendale断层，根据走向变化由西向东分成四段，分别编号为1至4，三个隐伏断层由西向东分别编号为5至7.经过多次测试之后，发现断层上的错动量基本集中在10km深度以内，同时，考虑到InSAR数据对断层错动的分辨率［11］，最终将断层深度固定为12km，相关参数如表2所示.断层面划分成约1km×1km的矩形位错元，倾角作为未知量参与非线性反演.

表2 断层单元相关参数Table 2 Parameters of fault segments

3.3 位移场离散采样

干涉雷达获得的形变场是连续的，一幅干涉图像中包含了大约106个数据点，像元之间具有较强的相关性，若利用全部的观测数据点进行反演，会造成计算困难.本文采用四叉树离散采样法［25］对观测数据离散采样.该方法利用形变量变化梯度值为图像采样依据，采样后的观测数据，既能有效提高计算效率，同时不会在近场剧烈形变区域丢失信息，可以代表整个形变场的特征，降低了噪音对形变场的影响［26］.经过测试，选用形变梯度0.06为阈值，即当形变区域梯度大于0.06时将区域分割为4个小块，直至形变梯度小于0.06为止，同时忽略了干涉图中相关系数小于0.25区域内的数据点.最终从336、337和631三个轨道的干涉图里分别提取出了1152，2513，733个数据点加入反演.

4 反演结果及讨论

4.1 断层几何形态及同震破裂分布

通过三步反演，最终得到了令最大后验概率趋于最大值的解.每一步反演的相关结果如表3所示.

表3 三步反演的结果Table 3 Results of the three－step method

本次事件总共释放的地震矩达到了5.0×1019N·m，对应的矩震级约为7.1.模型的滑动分布如图4a、表4所示.

表4 各断层整体错动特征Table 4 Overall Fault Slip Parameters

图4 （a）各断层面上位错分布图；（b）模型中Greendale断层地表错动与其他结果对比图Fig.4 （a）Slip distribution on each fault segment；（b）Comparison of surface rupture calculated in our model with other data

对于Darfield地震的同震破裂模型，已经有两篇研究成果发表［4，10］.其中 Beavan等采用了 GPS数据（尚未公开）以及 ALOS 337轨道获得的InSAR数据进行联合反演，Elliot等采用了ALOS三个轨道的InSAR数据以及P波和SH波形数据分别反演.本文同样采用了三个轨道的InSAR数据，但部分数据不同于Elliot等的模型中所用数据［10］，本研究中的336轨道震前、震后数据离发震时间更近，因此获得的数据具有更高的性噪比，对于研究断层东段的破裂有很大帮助.

本研究的反演结果显示，Greendale断层总长度为45km，1号子断层位于它的最西端，倾角为东北倾70°，最东边有一个最大值约3.5m的破裂集中区，往西逐渐减少.最西端（第1段）的3km范围滑动量接近于0，与地表观测吻合较好.Greendale断层中央的2号和3号子断层是主要破裂段，倾角分别为南倾87°和85.5°.2号子断层面上的破裂较集中，滑动发生在深度5km以上区域，最大滑动约为7.0m.3号子断层面上的滑动最为显著，出现了整个断层的最大滑动，约8.5m，该段部分区域的破裂在深度9.0km处仍达到了3.1m.Greendale断层东段，即4号子断层，是一个南倾69.5°的断层，最大破裂约为4.2m，发生在与3号断层相连的地方.在它的末端，出现了一部分的逆冲分量，最大约为1.0m.

断层段6长度为6km，倾角为59°（倾向东南）.反演显示该段以逆冲为主，最大逆冲分量约为5.3m，整个断层面上都有比较大的滑动，并且一直延伸到10km深度，该断层段是此次地震中初始破裂所在位置，因此区域地震波形反演得到了以逆冲为主的初始破裂震源机制解［1，9］.该断层的西南端与Greendale断层相连.

5号断层段也是一个逆冲断层，倾角为50°（倾向西北）.它的破裂主要出现在中部5km以上的较小区域内，最大的滑动为3.8m，相当于一个Mw＝6.2级地震，最大滑移量比Beavan等［4］的模型值略大，与Elliott等［10］的模型值比较一致，能够很好地解释干涉条纹中5号断层段附近的独立封闭形变条纹.

7号断层段在InSAR干涉图上没有明显的显示，因此并没有加入到Beaven等［4］的模型中，计算显示的破裂也很小，最大破裂不超过2.0m，以走滑破裂为主，对整个形变场贡献较小.

4.2 LOS位移场残差

图5是正演得到 的三条轨道上预测的InSAR形变场以及残差.336和337轨道，信号较好，获得了大范围清晰形变场，模型正演产生的条纹整体形态与观测比较吻合.与前人的结果［4，10］比较，本文的模型能够更好地解释形变观测.特别是，前人的破裂模型无法有效消除6号断层段东侧的剩余残差，而从本文的模型残差图中可以看到，无论是336轨道还是337轨道，该区域近场残差最大仅为20.0cm左右.

631轨道基线较长，且震前震后雷达数据获取时间间隔较长，因而相干性不高，只得到了不连续的形变条纹，但能看到大体形变趋势.从模型计算预测的形变条纹看，其基本趋势与观测数据是一致的，特别是5、6号逆冲断层段引起的单独条纹带比较明显.Greendale断层东西两端存在一部分残差，说明形变场还存在一部分形变无法用现有模型有效解释，该残差可能是由于形变场中包含的非弹性形变信息，如沙土液化等因素所致.沙土液化现象在本次地震中大量存在，也造成了一定的地面形变［4］.

4.3 与前人结果的比较

断层几何参数方面，本研究所用的4条断层位置与Elliott等［10］所用的断层模型基本一致，而与Beavan等［4］所用模型仅相差一条整体破裂较小的7号断层段.从断层几何形态来看，本文的模型与前人研究给出的Greendale断层倾角变化一致，6号断层的倾角在本文的模型和Elliott等的反演结果［10］中均为约50°，而Beavan等［4］计算给出的倾角为75°.考虑到该段以逆冲断层活动为主，较缓的倾角可能更加合理.

从错动分布来看，我们计算得到靠近震源的断层，即6号断层，整个断层面上都有较大的错动，在约3km深度的地方存在3.0m左右的滑动.从残差上也能看到，这种滑动分布能比较合理地解释6号断层东侧的位移.Elliott等［10］的结果显示在这条断层上只在4km深度以下的部分区域才存在较明显的错动，可能低估了这条断层段的滑动量.对于Greendale断层，Beavan等［4］，Elliott等［10］以及我们的结果都显示Greendale断层的最大破裂出现在中部偏东的区域.Beavan等［4］的模型最大滑动量（约5.0m）远小于Elliott等模型的约8.2m以及我们的结果约8.5m.这可能是Beavan等［4］所采用的数据缺乏足够的近场信息造成的.在4号子断层的中部3～4km深度，本文模型显示存在超过4m的破裂，超过了其它两篇文献中模型的滑动量，该差异可能来源于本文采用的336轨道数据的形变信号优于Elliott等［10］采用的形变数据，因而本文模型能提供4号子断层附近更多的破裂信息.

我们选取Greendale断层0～1km深度断层面上的位移近似表征它的表面滑动，与地质调查结果［8］、航拍数据等遥感数据［10］以及 Elliott等Greendale断层模型中0～1km深度上的位移［10］进行对比（图4b）.结果表明，观测数据之间匹配较好，地表的最大破裂出现在2号和3号子断层交界的区域，最大值在5.5～6m之间，之后向两边减小.在我们的模型中，地表的最大错动出现在3号断层的最西端，滑动量为6.2m，与观测结果基本吻合，而在向两边减小的趋势也与观测结果大致吻合.Elliott等［10］的模型中可能由于过强的滑动约束使断层滑动过于光滑，该模型中的地表错动在1号、2号、3号子断层都小于地表野外观测数据，而在4号子断层上又高于地表观测，峰值约4.0m出现在3、4号子断层交界的地方，比观测结果偏东7km左右.因此，本文的模型结果比Elliott等［10］的反演结果更合理，能更好解释地表位移数据.在2号子断层中部偏东部分，本文模型中给出的地表位移高于地表观测结果2.0m左右，这可能与2号断层和6号断层在此处相邻，它们的滑移存在较高相关造成模型误差增大有关.

4.4 Darfield地震对周围区域应力场的影响

Darfield地震是该地区应力场高度集聚后的调整，同时会对其附近应力场中余震的发生产生一定影响［27］.Darfield地震的余震既有逆冲型地震，也有走滑型地震，反映了该地区比较复杂的区域应力状态.我们选取Greendale断层的震源机制（走滑为主）作为接收断层上地震的震源机制，计算Darfield地震对于周围区域应力场的影响［28－29］.接收断层的参数选如表5.

图5 位移条纹和模型残差（由缠绕条纹表示，一个周期：－11.81－11.81cm）Fig.5 Deformation fringes and residues of each path（Showed by wrapped fringes，value in one cycle is from －11.81cm to 11.81cm）

表5 选定的接收断层参数Table 5 Parameter of receiving fault

图6a中显示了这一系列断裂对6.0km深度范围应力场的加卸载作用.与余震的分布对照可以看到，约有65%的余震发生距在发震断层3km距离范围内，这是由于主震发生之后断层附近积累的应力没有完全释放，而发震后的断层面容易发生破裂，因此在断层面上发生了后续破裂.如Greendale断层西侧5号断层段附近，虽然主震减小了该区域的应力，但仍能明显地看到断层面上发生了密集余震.远离断层的余震带主要是以43.55°S，172.05°E为中心的一个矩形，这恰好与库仑应力的一个增加区域（红区）相吻合，因此我们推测，这个区域内可能存在一些隐伏断层，震后释放了Darfield地震施加的部分应力.Greendale断层东南端的蓝区也出现了较多的余震，这个余震区域呈带状分布，向东北延伸到2011年Christchurch地震的发震断层，我们认为这部分小震预示着区域应力释放的方向从Greendale断层转变到一条西南—东北走向的断层上，Christchurch地震的发生可能与此有关.

2011年Christchurch地震区域正好处于一个库仑应力增大的区域.选定截面为Christchurch地震的断层面，可以看到在断层面上库仑应力增加值在0.5～1.2bar之间，在震源处［7］增大量约为1.0 bar（图6b）.可见仅仅从库仑应力的变化来考虑，Darfield地震的发生对于Christchurch地震起到了一定的触发作用.

4.5 区域应力场分配与震源机制

图6 （a）断层破裂造成周围区域库仑应力场变化（圆点为震后一个月内余震分布）；（b）Darfield地震引起的Christchurch地震断层面上的库伦应力改变Fig.6 （a）The change of the Coulomb stress around the faults（The dots are the after－shocks in one month）；（b）The change of the Coulomb stress on the fault of the Christchurch earthquake caused by the Darfield earthquake

澳洲板块和太平洋板块的边界带贯穿新西兰南岛.根据GPS数据计算发现（图1），在此地区太平洋板块相对澳洲板块存在约42mm／a方向为南偏西76°的运动［3，29］.边界带内的阿尔卑斯断裂带的走向为北偏东60°［7］，存在约40mm／a的右旋走滑［3］.因此边界带内阿尔卑斯断裂带以外仍有约10mm／a垂直于边界带的相对缩短，即北偏西约60°的压应变在新西兰南岛内积累.Darfield地震破裂主要发生在走向33°的6号逆冲断层段，以及近东西走向的Greendale右旋走滑断层上，它是板块边界带区域应力的一种表现.区域震源机制解显示此次地震的初始破裂是发生在6段断层段上的逆冲破裂［1，9］，其最大主压应力的方向垂直于断层的走向，因此6号断层段的震源机制对应于北偏西57°的压应力.其后发生的破裂主要为近东西走向的Greendale断层上的右旋滑移，其震源机制对应于北偏西约45°的压应力.因此这二个主破裂段释放的应力场就是板块边界带应变分配在该地区残留的构造应力场，Darfield地震的复杂性体现了区域构造应力场的特点.

5 结 论

通过完整覆盖发震区域且信噪比较高的SAR数据对，获得了高空间分辨率的同震形变场.根据形变信息，可以对Darfield地震的断层形态和同震破裂分布进行研究.本文中断层错动模型在地表的滑动，相比前人结果，与地质调查和遥感观测数据更为吻合.根据错动模型计算得到的库伦应力变化场跟远离发震断层的余震有较好的对应.2011年2月21号的Christchurch地震震源处于Darfield地震同震应力增加区域内，库仑应力增加约1.0bar，证明Darfield地震对于余震Christchurch地震存在应力加载效应.

（References）

［1］Gledhill K，Ristau J，Reyners M，et al.The darfield（Canterbury，New Zealand）Mw7.1earthquake of september 2010：apreliminary seismological report.Seismological Research Letters，2011，82（3）：378－386，doi：10.1785／gssrl.82.6.378.

［2］Reyners M. Lessons from the destructive Mw6.3 Christchurch，New Zealand，earthquake.Seismol.Res.Lett.，2011，82（3）：371－372.

［3］Norris R J，Cooper A F.Late quaternary slip rates and slip partitioning on the alpine fault，New Zealand.J.Structural Geol.，2001，23（2－3）：507－520.

［4］Beavan J，Samsonov S，Motagh M，et al.The MW7.1 Darfield（CANTERBURY）earthquake：Geodetic observations and preliminary source model.Bull.N.Z.Soc.Earthq.Eng，2010，43：228－235.

［5］Wallace L M，Beavan J.A large slow slip event on the central Hikurangi subduction interface beneath the Manawatu region，North Island，New Zealand.Geophysical Research Letters，2006，33（11），doi：10.1029／2006GL026009.

［6］Wallace L M，Beavan J，McCaffrey R，et al.Balancing the plate motion budget in the South Island，New Zealand using GPS，geological and seismological data.Geophysical Journal International，2007，168（1）：332－352.

［7］USGS http：／／earthquake. usgs. gov／earthquakes／eqinthenews／2010／us2010atbj／

［8］Quigley M，Dissen R V，Villamor P，et al.Surface rupture of the greendale fault during the darfield （canterbury）earthquake，New Zealand：initial findings.Bulletin of the New Zealand Society for Earthquake Engineering，2010，43（4）：236－242.

［9］Holden C，Beavan J，Fry B，et al.Preliminary source model of the Mw7.1Darfield earthquake from geological，geodetic and seismic data ／／ Proceedings of the Ninth Pacific Conference on Earthquake Engineering Building an Earthquake－Resilient Society，Auckland，New Zealand，2011.

［10］Elliott J R，Nissen E K，England P C，et al.Slip in the 2010—2011Canterbury earthquakes，New Zealand.Journal of Geophysical Research，2012，117，doi：10.1029／2011JB008868.

［11］Delouis B，Giardini D，Lundgren P，et al.Joint inversion of InSAR，GPS，Teleseismic，and Strong－Motion data for the spatial and temporal distribution of earthquake slip：Application to the 1999izmit mainshock.Bulletin of the Seismological Society of America，2002，92（1）：278－299.

［12］Delouis B，Nocquet J M，Vallée M，et al.Slip distribution of the February 27，2010 Mw＝8.8maule earthquake，central chile，from static and high－rate GPS，InSAR，and broadband teleseismic data.Geophysical Research Letters，2010，37（17），L17305，doi：10.1029／2010GL043899.

［13］Rosen P A，Hensley S，Peltzer G，et al.Updated repeat orbit interferometry package released.Eos，Transactions American Geophysical Union，2004，85（5）：47.

［14］Fialko Y，Simons M，Agnew D.The complete（3－D）surface displacement field in the epicentral area of the 1999 Mw7.1 Hector Mine Earthquake，California，from space geodetic observations.Geophysical Research Letters，2001，28（16）：3063－3066.

［15］Farr T G，Rosen P A，Caro E，et al.The shuttle radar topography mission.Rev.Geophys.，2007，45（2），doi：10.1029／2005RG000183.

［16］Hanssen R F.Radar Interferometry－Data Interpretation and Error Analysis（V.2）.Dordrecht：Kluwer Academic Pub.，2001.

［17］孙建宝，徐锡伟，石耀霖等.东昆仑断裂玛尼段震间形变场的InSAR观测及断层滑动率初步估计.自然科学进展，2007，17（10）：1361－1370.Sun J B，Xu X W，Shi Y L，et al.InSAR observation of the inter－seismic deformation of Mani segment in the East Kunlun fault and estimation of the slip rate.Progress in Natural Science （in Chinese），2007，17（10）：1361－1370.

［18］Fukuda J，Johnson K M.A fully Bayesian inversion for spatial distribution of fault slip with objective smoothing.Bulletin of the Seismological Society of America，2008，98（3）：1128－1146.

［19］Fukuda J，Johnson K M.Mixed linear－non－linear inversion of crustal deformation data：Bayesian inference of model，weighting and regularization parameters. Geophysical Journal International，2010，181（3）：1441－1458.

［20］Okada Y.Internal deformation due to shear and tensile faults in a half－space，Bull.Seism.Soc.Am.，1992，82（2）：1018－1040.

［21］Sun J B，Shen Z K，Bürgmann R，et al.A three－step Maximum－A－Posteriori probability method for INSAR data inversion of coseismic rupture with application to the April 14，2010 Mw6.9Yushu，China earthquake.Journal of Geophysical Research：Solid Earth，doi：10.1002／jgrb.50244.

［22］Ingber L.Simulated annealing：Practice versus theory.Math.Comput.Model.，1993，18（11）：29－57.

［23］Quigley M，van Dissen R，Litchfield N，et al.Surface rupture during the 2010 Mw7.1Darfield （Canterbury）earthquake：Implications for fault rupture dynamics and seismic－hazard analysis.Geology，2012，40（1）：55－58.

［24］GeoNet http：／／www.geonet.org.nz／

［25］Jónsson S，Zebker H，Segall P，et al.Fault slip distribution of the 1999 Mw7.1Hector Mine，California，earthquake，estimated from satellite radar and GPS measurements.Bulletin of the Seismological Society of America，2002，92（4）：1377－1389.

［26］孙建宝，徐锡伟，沈正康等.基于线弹性位错模型及干涉雷达同震形变场反演1997年玛尼Mw7.5级地震参数—I.均匀滑动反演.地球物理学报，2007，50（4）：1097－1110.

Sun J B，Xu X W，Shen Z K，et al.Parameter inversion of the 1997Mani earthquake from INSAR co－seismic deformation field based on linear elastic dislocation model—Ⅰ.uniform slip inversion.Chinese J.Geophys.（in Chinese），2007，50（4）：1097－1110.

［27］King G C P，Stein R S，Lin J.Static stress changes and the triggering of earthquakes.Bulletin of the Seismological Society of America，1994，84（3）：935－953.

［28］Lin J，Stein R S.Stress triggering in thrust and subduction earthquakes and stress interaction between the southern San Andreas and nearby thrust and strike－slip faults.Journal of Geophysical Research，2004，109（B2），doi：10.1029／2003JB002607.

［29］Sella G F，Dixon T H，Mao A L.REVEL：A model for recent plate velocities from space geodesy.Journal of Geophysical Research，2002，107（B4），doi：10.1029／2000JB000033.