文章快速检索  
  高级检索
地震同震滑动分布反演平滑因子的确定
王乐洋1,2,3 , 赵雄1,2     
1. 东华理工大学测绘工程学院, 江西 南昌 330013;
2. 流域生态与地理环境监测国家测绘地理信息局重点实验室, 江西 南昌 330013;
3. 江西省数字国土重点实验室, 江西 南昌 330013
摘要:针对大地测量地震同震滑动分布反演中平滑因子(又称正则化参数)的确定问题,一般采用模型粗糙度与数据拟合残差之间的折中曲线来确定(为便于区分,本文将该方法称为"L曲线")。本文在L曲线的基础上提出一种确定平滑因子的新方法——折中相交曲线法。模拟试验反演结果表明,利用折中相交曲线法确定的平滑因子反演地震滑动分布各参数精度要优于L曲线法。利用折中相交曲线法确定平滑因子反演拉奎拉与台湾美浓实际地震,并与L曲线法反演结果进行对比分析。分析结果表明,利用折中相交曲线法确定的平滑因子反演拉奎拉与台湾美浓实际同震滑动分布各参数结果均在国内外其他学者研究的范围内,并且利用折中相交曲线法确定平滑因子较L曲线法具有计算效率高、无需依赖数据拟合度、确定平滑因子大小更为合适等优点。
关键词:地震同震滑动分布    平滑因子    折中相交曲线法    L曲线法    拉奎拉地震    美浓地震    
Determination of Smoothing Factor for the Co-seismic Slip Distribution Inversion
WANG Leyang1,2,3 , ZHAO Xiong1,2     
1. Faculty of Geomatics, East China University of Technology, Nanchang 330013, China;
2. Key Laboratory of Watershed Ecology and Geographical Environment Monitoring, NASG, Nanchang 330013, China;
3. Key Laboratory for Digital Land and Resources of Jiangxi Province, Nanchang 330013, China
Foundation support: The National Natural Science Foundation of China (Nos. 41664001; 41204003; 41874001); The Support Program for Outstanding Youth Talents in Jiangxi Province (No. 20162BCB23050); The National Key Research and Development Program(No. 2016YFB0501405)
First author: WANG Leyang (1983—), male, PhD, associate professor, majors in geodetic inversion and geodetic data processing.E-mail:wleyang@163.com
Abstract: For the determination of the smoothing factor (also known as the regularization parameter) in the co-seismic slip distribution inversion, the compromise curve between the model roughness and the data fitting residual is generally used to determine (in order to distinguish the method proposed in this paper, the method is called "L curve" according to its shape). Based on the L-curve, the eclectic intersection curve as a new method are proposed to determine the smoothing factor in this paper. The results of the simulated experiment show that the inversion accuracy of the parameters of the seismic slip distribution with the smoothing factor determined by the eclectic intersection curve method is better than that of the L curve method. Moreover, the eclectic intersection curve method and the L curve method are used to determine the smoothing factor of L'Aquila and Taiwan Meinong earthquake slip distribution inversion respectively, and the inversion results are compared and analyzed. The analysis results show that the L'Aquila and Meinong of Taiwan actual earthquake slip distribution results are in the range of other scholars at home and abroad, and compared with the L curve method, the eclectic intersection curve method has advantages of high computation efficiency, no need to depend on data fitting degree and more appropriate of smoothing factor and so on.
Key words: co-seismic slip distribution     smoothing factor     eclectic intersection curve method     L curve method     the L'Aquila earthquake     the Meinong earthquake    

在地震机制研究领域中,利用大地测量资料研究地震机制是地学领域的前沿研究热点[1]。大地测量数据(如GPS、InSAR数据等)因观测范围广、精度高,反演出的滑动分布结果更为精细、全面,近年来常被国内外学者用于进行地震同震滑动分布反演研究[2-8]。在震源机制研究中,位错模型是反演震源参数的有效模型[9-11],一般采用弹性半空间矩形位错模型进行地震同震滑动分布反演[12-13]

在地震同震滑动分布反演过程中,地表形变位移与滑动分布参数之间为线性关系。滑动分布参数与地表形变位移的线性模型依赖于格林函数矩阵的构造,而格林函数矩阵元素通常是确定断层的几何参数后,将断层破裂面进一步离散成若干断层子块,由每个子块的单位走滑与倾滑引起的地表形变位移构成。另外,格林函数矩阵在地震滑动分布反演过程中通常是病态的,一般采用Tikhonov正则化方法[14]解决。

在正则化过程中,平滑因子的确定是关键问题之一,不同的平滑因子会使反演结果大不相同。目前主要方法有广义交叉核实法(generalized cross validation, GCV)[15]、岭迹法[16-17]、L曲线法[18-20]。GCV法理论严密,但有时GCV方法函数绘制的曲线过于平缓,很难定位到最优平滑因子[21];岭迹法计算简单,但是在确定平滑因子的时候具有一定的主观性;在地震同震滑动分布反演中,目前最常用的平滑因子确定方法为L曲线法,但是L曲线法在确定平滑因子过程中具有以下几个缺点:①无法从L曲线图中根据横纵坐标轴直观读取平滑因子的大小;②L曲线法过于依赖数据拟合精度[22],求解过程可能不收敛[23];③计算时间过长。

针对以上问题,本文提出一种确定地震同震滑动分布反演平滑因子的新方法——折中相交曲线法(eclectic intersection curve method,EI)(简称EI曲线法)。通过模拟试验以及拉奎拉与台湾美浓实际震例反演验证EI曲线法确定平滑因子的优势与可行性。

1 EI曲线法确定平滑因子的基本原理 1.1 同震滑动分布基本方程

在同震滑动分布反演中,地表同震位移与断层滑动量之间为线性关系,如式(1)所示

(1)

式中,d表示地表形变数据;G表示格林函数矩阵;m为滑动参数;ε为随机观测误差。

在同震滑动分布反演中,通常将断层面划分成n个矩形单元,基于矩形位错模型计算矩形单元的单位走滑量和单位倾滑量引起地表k个观测点的形变值,构成了k×2n个格林函数矩阵G。将断层面离散化处理,直接用式(1)求解会出现病态问题,从而使滑动分布解不稳定,因为此时格林矩阵G是严重秩亏的。为了避免这种情况,通常对断层面施加一定的平滑约束条件,一方面增加虚拟观测数据解决系数矩阵G的秩亏问题,另一方面约束了相邻断层单元间出现大的梯度变化。

本文采用目前较成功用于断层面平滑约束的拉普拉斯二阶平滑矩阵对断层单元进行约束[24]。根据Tikhonov正则化原理,相应的准则为

(2)

由式(1)、式(2)可得地震滑动分布解为

(3)

式中,α为平滑因子;T为拉普拉斯二阶平滑矩阵(R=TTT为正则化矩阵);Ω(m)为稳定泛函;||·||2为·的二范数;P1P2为实际观测数据与虚拟观测数据内部观测权阵。其中平滑因子α,一般由拟合残差||dGm||P12和平滑度||Tm||P22之间的折中曲线来确定,即本文中L曲线法。

1.2 L曲线法的基本原理

由式(2)可知,||dGm||P12和||Tm||P22均是关于α的函数,L曲线法的基本原理则是取一系列α值, 以||dGm||P12为纵坐标,||Tm||P22为横坐标绘图,得到多组(||Tm||P22,||dGm||P12)坐标。

在地震同震滑动分布反演中,平滑因子的选择是根据拟合残差||dGm||P12与平滑度||Tm||P22之间的折中曲线来确定。在L曲线中, 折中点一般在L曲线图中曲率最大点附近,通过试错,调试选取出一个较为合适的点。如图 1所示,在折中点α左侧目标函数(式(2))倾向于拟合残差||dGm||P12的大小,而在折中点α右侧目标函数(式(2))倾向于平滑度||Tm||P22大小。如何确定平滑因子α的大小是地震同震滑动分布反演的主要问题之一。由图 1及L曲线法定义可知,L曲线法在确定平滑因子的过程中存在以下几点不足:①无法根据横纵坐标直接从图 1中获知平滑因子α的大小;②L曲线法绘制过程中过度依赖于||dGm||P12与||Tm||P22之间的数据拟合精度;③若要获取合适的平滑因子,需要多组(||Tm||P22,||dGm||P12)坐标,使计算效率大大降低。针对以上几点,本文提出一种确定平滑因子的新方法——EI曲线法。

图 1 L曲线图 Fig. 1 L curve

1.3 EI曲线法基本原理

平滑因子一般是由拟合残差||dGm||P12与平滑度||Tm||P22之间的折中曲线来确定。上文提到在L曲线中,该折中点位于L曲线上曲率最大点附近,而本文提出的EI曲线法中折中点则是定义为:存在某点α使||dGm||P12与||Tm||P22的值相等,则称该点为折中点。即在该点左侧所有||dGm||P12的值要大于||Tm||P22的值,在该点右侧所有||Tm||P22的值要大于||dGm||P12。但是在地震滑动分布反演中,实际观测数据与虚拟观测数据内部观测权均是根据先验信息所定,例如一般实际观测数据如GPS数据则是根据观测中误差定权。通常此时实际观测数据和虚拟观测数据的内部观测权阵P1P2可能是不准确的[25],因为它们的单位权方差大小没有统—,即σ01σ02。故不能直接以初始||dGm||P12与||Tm||P22相等时的α作为所求平滑因子大小。

针对以上问题,EI曲线法首先利用赫尔默特方差分量估计法[1, 25]统一反演过程中各类数据单位权方差,在统一单位权方差过程中需要加入权修正因子[26]ab进行修正,使σ01=σ02,此时相应的目标函数可以表达如下

(4)

式中,ab为实际观测数据与虚拟观测数据内部权阵P1P2的权修正因子。

EI曲线法基本原理是,在同一坐标系下以α为横坐标,分别以a||dGm||P12b||Tm||P22为纵坐标绘制关于α的曲线图,从而取使a||dGm||P12b||Tm||P22值相等时的α为EI曲线法所确定平滑因子,一般EI曲线如图 2(a)所示。为了减少计算时间本文设置阈值,当a||dGm||P12b||Tm||P22之间的绝对值差值小于阈值(一般大小为10-3~10-6)时输出曲线(如图 2(b)所示),这样便减少了计算时间,提高了计算效率。

图 2 EI曲线图 Fig. 2 EI curve

综上,利用EI曲线法确定地震滑动分布平滑因子步骤可总结如下:

(1) 利用方差分量估计法统一实际观测数据与虚拟观测数据单位权方差:即使σ01=σ02,并求得ab

(2) 绘制EI曲线图,给定平滑因子α初值、步长、终值,以α为横坐标,以a||dGm||P12b||Tm||P22为纵坐标,分别绘制关于α的曲线图,当a||dGm||P12b||Tm||P22之间的绝对值差值小于阈值时输出曲线。

图 2的EI曲线图可以看出,本文提出的EI曲线法在确定平滑因子时较L曲线法具有以下几点优势:①在EI曲线图中可以直观地根据横坐标读取平滑因子α值的大小;②由于直接绘制||dGm||P12、||Tm||P22关于α的曲线,故无需依赖||dGm||P12与||Tm||P22之间的数据拟合精度;③当α取相同初始值、步长、终值时,由于EI曲线法在||dGm||P12与||Tm||P22值近似相等处停止计算,直接输出,从而本文方法较L曲线法大大地减少了计算时间。

2 模拟试验

为了验证本文方法的可行性,本文作了以下模拟试验。模拟试验中,模拟断层面几何参数如下:断层面几何中心位置X=0km,Y=0km,断层顶深0km,长度30km,宽度30km,走向角60°,倾角45°,滑动角45°。模拟GPS 3方向观测数据如图 3所示,根据点位与震中距离设置形变点密度,并给形变点施加观测误差N(0,32mm2)。

图 3 GPS 3方向模拟观测点 Fig. 3 Distribution of three directions simulated points of GPS

在模拟试验中,本文分别利用L曲线法、EI曲线法确定平滑因子,模拟试验中平滑因子α均以0.01为初值、0.001为步长、1为终值绘制L曲线图、EI曲线图(其中,L曲线法确定平滑因子大小为0.179,EI曲线法确定的平滑因子大小为0.095。在利用EI曲线法确定平滑因子过程中,利用方差分量估计法确定实际观测数据与虚拟观测数据内部权阵修正因子大小分别为a=1,b=0.0037),并将两种方法确定的平滑因子进行地震滑动分布反演。利用L曲线法、EI曲线法确定平滑因子如图 4(a)(d)所示;两种方法反演地震滑动分布结果以及两种方法反演地震滑动分布与模拟值的差值分布如图 4(b)(e)(c)(f)所示。

图 4 模拟试验两种方法反演地震滑动分布结果及解的残差 Fig. 4 Slip distribution results of simulation experiment and the residuals of two methods

模拟试验预设参数以及利用L曲线法、EI曲线法反演滑动分布各参数结果见表 1。两种方法确定平滑因子反演地震滑动分布部分参数结果与模拟值的差异,以及EI曲线法反演地震滑动分布部分参数结果,较L曲线法精度改进程度见表 2

表 1 模拟试验反演地震滑动分布结果 Tab. 1 The slip distribution inversion results in simulation experiments
方法 平滑因子
α
计算时间
/s
最大滑动量
/m
平均滑动量
/m
平均滑动角
/(°)
矩震级
/Mw
地震矩
/(1018N·m)
均方根误差
/mm
预设参数 1.4721 0.3104 45 6.5822 8.3806
L曲线法 0.179 126.6230 1.3049 0.3168 45.9012 6.5881 8.5546 3.0
EI曲线法 0.095 18.8387 1.3722 0.3138 45.8826 6.5853 8.4725 2.9

表 2 模拟试验两种方法反演各参数结果与模拟值的差异 Tab. 2 Comparison the results of each parameters and simulation values of two methods in simulation experiments
与模拟值的差值 最大滑动量
/m
平均滑动量
/m
平均滑动角
/(°)
L曲线法 0.1672 0.0064 0.9012
EI曲线法 0.0999 0.0034 0.8826
L-EI* 5.23% 0.95% 0.04%
注:L-EI*表示利用EI曲线法确定平滑因子反演地震滑动分布各参数结果较L曲线法改进程度。

表 1给出了利用L曲线法和EI曲线法确定平滑因子反演地震同震滑动分布各参数的结果及两种方法确定平滑因子所用时间;表 2给出了两种方法反演地震同震滑动分布各参数结果与模拟值的差异。在各参数反演精度上,从表 1表 2可以看出,利用EI曲线法确定的平滑因子反演地震滑动分布结果在最大滑动量、平均滑动量、矩震级等参数反演结果均要优于L曲线法反演结果;在计算时间上,在相同的初值、步长、终值情况下,EI曲线法比L曲线法在计算时间上减少了107.7843s,大大地提高了计算效率;在均方根误差上,EI曲线法较L曲线法减小0.1mm,进一步说明EI曲线法确定平滑因子反演地震滑动分布结果精度要高于L曲线法。

3 实际震例反演 3.1 拉奎拉地震滑动分布反演

2009年4月6日,意大利中部拉奎拉地区发生Mw6.3级地震。据美国地质调查局(USGS)网站发布,地震震中位置为(13.334°E,42.334°N),震源深度为8.8km。地震造成1000余人的伤亡,大量的房屋被损毁或破坏[27]。本文从文献[6, 28]中获取拉奎拉地震InSAR降轨数据位移1254个,升轨位移数据1282个,并利用文献[29]中通过非线性反演方法反演得到的拉奎拉地震断层的走向角、倾角及断层中心位置等参数来固定断层破裂面。在此基础上将断层破裂面沿断层走向、倾向把断层长度、宽度扩展至30km、27km,并将破裂面延伸至地表,将断层面均匀剖分成1.5×1.5km大小的矩形单元,共得到了360个矩形断层单元。分别采用EI曲线法与L曲线法确定平滑因子进行地震滑动分布反演。

在拉奎拉实际震例中,平滑因子α均以0.05为初值、0.001为步长、1为终值绘制L曲线图、EI曲线图(其中,L曲线法确定平滑因子大小为0.588与文献[28]中用L曲线法确定的平滑因子大小0.6,两者大小相近;EI曲线法确定的平滑因子大小为0.159。在利用EI曲线法确定平滑因子过程中,利用方差分量估计法确定实际观测数据与虚拟观测数据内部权阵修正因子,大小分别为a=1,b=0.0024)。利用EI曲线法和L曲线法确定平滑因子如图 5(a)(c)所示;两种方法反演拉奎拉地震成果及两种方法反演滑动分布结果差值分布如图 5(b)(d)(e)所示。

图 5 两种方法反演拉奎拉地震滑动分布结果及差值分布 Fig. 5 The slip distribution inversion results of L'Aquila earthquake and the residuals of two methods

利用L曲线法、EI曲线法确定的平滑因子反演拉奎拉地震滑动分布各参数结果见表 3;近年来国内外部分其他学者反演拉奎拉地震滑动分布各参数结果见表 4

表 3 用L曲线法、EI曲线法确定平滑因子反演拉奎拉地震滑动分布结果 Tab. 3 The slip distribution inversion results of L'Aquila earthquake by EI and L curve methods
方法 平滑因子
α
计算时间
/s
最大滑动量
/m
平均滑动量
/m
平均滑动角
/(°)
矩震级
/Mw
地震矩
/(1018N·m)
均方根误差
/mm
L曲线法 0.588 145.5253 1.005 0.1535 -101.3043 6.3478 3.7294 14.6
EI曲线法 0.159 19.8612 1.02 0.1468 -103.9264 6.3349 3.5676 13.8

表 4 拉奎拉地震断层滑动参数 Tab. 4 The slip parameters for L'Aquila earthquake
来源 L曲线法 EI曲线法 文献[28] 文献[6] 文献[27] 文献[30] 文献[31] USGS
最大滑动量/m 1.01 1.02 0.95 1.07 0.66 0.49 1.1
平均滑动角/(°) -101.3 -103.9 -96.4 -102.8 -105 -98 -98.5
矩震级/Mw 6.35 6.33 6.34 6.32 6.23 6.32 6.36 6.29
地震矩/(1018Nm) 3.73 3.57 3.63 3.43 2.8 3.2 3.9 3.4

图 5(b)(d)可知,利用L曲线法、EI曲线法确定的平滑因子进行拉奎拉地震反演时,其反演结果均显示拉奎拉同震滑动分布区域主要分布在地下4~15km,该结果同文献[28]研究成果一致;对比表 3表 4可以看出,本文两种方法确定的平滑因子反演拉奎拉地震滑动分布结果均在其他学者研究范围内。对于最大滑动量,利用L曲线法和EI曲线法反演最大滑动量大小分别为1.01、1.02m,文献[6]联合了不同波长、不同入射角的Enviast和ALOS卫星差分干涉数据,并联立了震区部分GPS三维形变数据,通过弹性三角位错模型,反演得到滑动分布的最大滑动量为1.07m,对于其最大滑动量略大于本文反演结果,这可能与约束地表的数据源有关。文献[28]利用最小二乘和总体最小二乘反演最大滑动量均为0.95m,略小于本文两种方法反演结果;对于平均滑动角,本文两种方法反演的结果分别为-101.3°、-103.9°,该研究结果均在其他学者研究范围内,并且与文献[6]反演结果相近;对于矩震级参数,本文两种方法反演结果相近,且在其他学者研究范围内。另外,文献[27]联合ASAR数据反演得到的断层滑动量为0.66m,给出的矩张量为2.80×1018N·m;文献[31]利用拉奎拉震中附近的GPS观测站数据,研究了拉奎拉地震的同震滑动分布和震后的余震滑动分布,得到的最大滑动量为1.1m,矩张量为3.90×1018N·m。

3.2 台湾美浓地震滑动分布反演

2016年2月6日3时57分(北京时间)我国台湾省高雄市美浓区发生了Mw6.4(global centroid moment tensor, GCMT)地震(下称美浓地震)。此次地震给台湾南部地区带来巨大的人员伤亡及财产损失。本文采用文献[32]获取美浓GPS站位移数据225个参与计算美浓地震滑动分布反演。利用文献[32]中均匀滑动分布反演的断层走向、倾向、断层深度等参数,作为本文通过多峰值颗粒群算法[29]进行非线性断层参数反演的初始参数大小。通过非线性反演最终确定适合本文数据的美浓地震断层参数走向角为strike=303.6°,倾角为dip=16.2°,震源中心深度depth=12.3km。将断层长和宽扩展至44、54km,并将其按照2km×2km矩形划分为594个单元。为了避免少量滑动单元出现方向矛盾问题,采用非负约束最小二乘法进行线性滑动分布求解[33]。分别采用EI曲线法与L曲线法确定平滑因子进行地震滑动分布反演。

在美浓实际震例反演中,平滑因子α均以0.05为初值、0.002为步长、1为终值绘制L曲线图、EI曲线图(其中,L曲线法确定的平滑因子大小为0.445,EI曲线法确定的平滑因子大小为0.328。在利用EI曲线法确定平滑因子过程中,利用方差分量估计确定实际观测数据与虚拟观测数据内部权阵修正因子,大小分别为a=1,b=0.8311)。利用EI曲线法和L曲线法确定平滑因子如图 6(a)(c)所示;两种方法反演美浓地震结果及两种方法反演结果差值分布如图 6(b)(d)(e)所示。

图 6 两种方法反演美浓地震滑动分布结果及差值分布 Fig. 6 The slip distribution inversion results of Meinong earthquake and the residuals of two methods

利用L曲线法、EI曲线法确定的平滑因子反演美浓地震滑动分布各参数结果见表 5;近年来国内外部分其他学者反演美浓地震滑动分布各参数结果见表 6

表 5 用L曲线法、EI曲线法确定平滑因子反演美浓地震滑动分布结果 Tab. 5 The slip distribution inversion results of Meinong earthquake by EI and L curve methods
方法 平滑因子
α
计算时间
/s
最大滑动量
/m
平均滑动量
/m
平均滑动角
/(°)
矩震级
/Mw
地震矩
/(1018N·m)
均方根误差
/mm
L曲线法 0.445 118.4569 0.4529 0.0517 45.3159 6.3441 3.6828 12.9
EI曲线法 0.328 43.7275 0.4793 0.0521 43.4510 6.3466 3.7140 12.8

表 6 美浓地震断层滑动参数 Tab. 6 The slip parameters of Meinong earthquake
来源 L曲线法 EI曲线法 文献[32] 均匀滑动 GCMT USGC 文献[34] 文献[35]
最大滑动量/m 0.453 0.479 0.517/0.553
平均滑动角/(°) 45.32 43.45 51.5 -51.5 21 38 45
矩震级/Mw 6.34 6.35 6.38 6.33 6.40 6.40 6.39 6.52

图 6(b)(d)可以看出,利用L曲线法、EI曲线法确定的平滑因子进行美浓地震反演时,其反演结果均显示美浓地震滑动分布区域主要分布在地下9~14km,该结果同文献[32]反演结果一致。表 6列举了美浓地震的部分研究成果,表 5表 6显示对于最大滑动量,本文两种方法反演结果分别为0.453、0.479m。文献[32]联合InSAR和GPS数据反演美浓地震结果显示最大倾滑量为0.517m、最大走滑量为0.553m,略大于本文两种方法反演结果;对于平均滑动角,本文两种方法反演的结果分别为45.32°、43.45°,其中,利用L曲线法反演结果与文献[34]反演结果相近。而EI曲线法反演结果在文献[34]与文献[32]反演结果之间;对于矩震级参数,本文两种方法反演结果相近,且在其他学者研究范围内;文献[35]反演的矩震级结果为Mw6.52,对于其反演结果较其他学者反演结果偏大的原因可能与没有较好地拟合GPS数据有关[32]

对于本文所做的模拟试验以及拉奎拉与台湾美浓实际震例试验,作出以下几点分析:

(1) 对于两种方法确定平滑因子反演地震滑动分布各参数精度,由模拟试验可知,利用EI曲线法确定平滑因子反演滑动分布的各参数结果与模拟值的差值均小于L曲线法。模拟试验及两个实际震例反演的均方根误差结果表明EI曲线法反演的均方根误差均小于L曲线法,从而进一步说明EI曲线法确定的平滑因子反演地震滑动分布精度要高于L曲线法。

(2) 对于计算效率,在模拟试验、拉奎拉及美浓实际震例中EI曲线法确定平滑因子时间较L曲线法分别减少107.7843、125.6641、74.7294s。总体上,EI曲线法较L曲线法效率提高2~7倍,故利用EI曲线法确定平滑因子较L曲线法可以有效地提高计算效率。

(3) 对于L曲线图与EI曲线图,由于EI曲线图横坐标即为平滑因子α的大小,故从EI曲线图上可以直观地读取平滑因子的大小,而从L曲线图中却不能直观读取α大小;L曲线图绘制时需要依赖||dGm||2和||Tm||2之间的数据拟合精度,EI曲线法则不需要。

4 结论

本文提出的EI曲线法为地震滑动分布反演中平滑因子的确定增加了一条新的途径。本文通过模拟试验及拉奎拉与美浓实际地震反演表明,EI曲线法较L曲线法确定平滑因子具有定位准确、反演精度高、计算效率高,以及EI曲线图横坐标直观反映平滑因子大小、无需依赖于数据拟合度等优势。本文提出的EI曲线法可以应用在地震滑动分布反演中,在其他领域的应用有待于进一步的研究与探索。

致谢: 感谢武汉大学温扬茂博士提供意大利拉奎拉地震的InSAR数据。


参考文献
[1] 许才军, 邓长勇, 周力璇. 利用方差分量估计的地震同震滑动分布反演[J]. 武汉大学学报(信息科学版), 2016, 41(1): 37–44.
XU Caijun, DENG Changyong, ZHOU Lixuan. Coseismic Slip Distribution Inversion Method Based on the Variance Component Estimation[J]. Geomatics and Information Science of Wuhan University, 2016, 41(1): 37–44.
[2] FUNNING G J, PARSONS B, WRIGHT T J, et al. Surface Displacements and Source Parameters of the 2003 Bam (Iran) Earthquake from Envisat Advanced Synthetic Aperture Radar Imagery[J]. Journal of Geophysical Research:Solid Earth, 2005, 110(B9): B09406.
[3] JIANG Zaisen, WANG Min, WANG Yanzhao, et al. GPS Constrained Coseismic Source and Slip Distribution of the 2013Mw6.6 Lushan, China, Earthquake and Its Tectonic Implications[J]. Geophysical Research Letters, 2014, 41(2): 407–413. DOI:10.1002/2013GL058812
[4] 许才军, 刘洋, 温扬茂. 利用GPS资料反演汶川Mw7.9级地震滑动分布[J]. 测绘学报, 2009, 38(2): 195–201, 215.
XU Caijun, LIU Yang, WEN Yangmao. Mw7.9 Wenchuan Earthquake Slip Distribution Inversion from GPS Measurements[J]. Acta Geodaetica et Cartographica Sinica, 2009, 38(2): 195–201, 215.
[5] 张国宏, 屈春燕, 宋小刚, 等. 基于InSAR同震形变场反演汶川Mw7.9地震断层滑动分布[J]. 地球物理学报, 2010, 53(2): 269–279.
ZHANG Guohong, QU Chunyan, SONG Xiaogang, et al. Slip Distribution and Source Parameters Inverted from Co-seismic Deformation Derived by InSAR Technology of Wenchuan Mw7.9 Earthquake[J]. Chinese Journal of Geophysics, 2010, 53(2): 269–279. DOI:10.3969/j.issn.0001-5733.2010.02.005
[6] 温扬茂, 何平, 许才军, 等. 联合Envisat和ALOS卫星影像确定L'Aquila地震震源机制[J]. 地球物理学报, 2012, 55(1): 53–65.
WEN Yangmao, HE Ping, XU Caijun, et al. Source Parameters of the 2009 L'Aquila Earthquake, Italy from Envisat and ALOS Satellite SAR Images[J]. Chinese Journal of Geophysics, 2012, 55(1): 53–65. DOI:10.6038/j.issn.0001-5733.2012.01.006
[7] 李志才, 张鹏, 金双根, 等. 基于GPS观测数据的汶川地震断层形变反演分析[J]. 测绘学报, 2009, 38(2): 108–113, 119.
LI Zhicai, ZHANG Peng, JIN Shuanggen, et al. Wenchuan Earthquake Deformation Fault Inversion and Analysis Based on GPS Observations[J]. Acta Geodaetica et Cartographica Sinica, 2009, 38(2): 108–113, 119. DOI:10.3321/j.issn:1001-1595.2009.02.003
[8] WANG Leyang, ZHAO Xiong, GAO Hua. A Method for Determining the Regularization Parameter and the Relative Weight Ratio of the Seismic Slip Distribution with Multi-source Data[J]. Journal of Geodynamics, 2018(118): 1–10.
[9] OKADA Y. Surface Deformation due to Shear and Tensile Faults in a Half-space[J]. Bulletin of the Seismological Society of America, 1985, 75(4): 1135–1154.
[10] OKADA Y. Internal Deformation due to Shear and Tensile Faults in a Half-space[J]. Bulletin of the Seismological Society of America, 1992, 92(2): 1018–1040.
[11] SUN Wenke, OKUBO S, VANÍČEK P. Global Displacements Caused by Point Dislocations in a Realistic Earth Model[J]. Journal of Geophysical Research:Solid Earth, 1996, 101(B4): 8561–8577. DOI:10.1029/95JB03536
[12] 孙建宝, 徐锡伟, 沈正康, 等. 基于线弹性位错模型及干涉雷达同震形变场反演1997年玛尼Mw7.5级地震参数-Ⅰ.均匀滑动反演[J]. 地球物理学报, 2007, 50(4): 1097–1110.
SUN Jianbao, XU Xiwei, SHEN Zhengkang, et al. Parameter Inversion of the 1997 Mani Earthquake from InSAR Co-seismic Deformation Field Based on Linear Elastic Dislocation Model-I. Uniform Slip Inversion[J]. Chinese Journal of Geophysics, 2007, 50(4): 1097–1110. DOI:10.3321/j.issn:0001-5733.2007.04.017
[13] 刘洋, 许才军, 温扬茂, 等. 2008年大柴旦Mw6.3级地震的InSAR同震形变观测及断层参数反演[J]. 测绘学报, 2015, 44(11): 1202–1209.
LIU Yang, XU Caijun, WEN Yangmao, et al. The InSAR Coseismic Deformation Observation and Fault Parameter Inversion of the 2008 Dachaidan Mw6.3 Earthquake[J]. Acta Geodaetica et Cartographica Sinica, 2015, 44(11): 1202–1209. DOI:10.11947/j.AGCS.2015.20140628
[14] TIKHONOV A N. Regularization of Incorrectly Posed Problems[J]. Soviet Mathematics Doklady, 1963, 4(1): 1624–1627.
[15] GOLUB G H, HEATH M, WAHBA G. Generalized Cross-Validation as a Method for Choosing a Good Ridge Parameter[J]. Technometrics, 1979, 21(2): 215–223. DOI:10.1080/00401706.1979.10489751
[16] 王振杰.大地测量中不适定问题的正则化解法研究[D].武汉: 中国科学院测量与地球物理研究所, 2003.
WANG Zhenjie. Research on the Regularization Solutions of Ill-posed Problems in Geodesy[D]. Wuhan: Institute of measurement and Geophysics of the Chinese Academy of Sciences, 2003. http://cdmd.cnki.com.cn/article/cdmd-80057-2004041260.htm
[17] 王振杰, 欧吉坤. 用L-曲线法确定岭估计中的岭参数[J]. 武汉大学学报(信息科学版), 2004, 29(3): 235–238.
WANG Zhenjie, OU Jikun. Determining the Ridge Parameter in a Ridge Estimation Using L-curve Method[J]. Geomatics and Information Science of Wuhan University, 2004, 29(3): 235–238.
[18] HANSEN P C. Analysis of Discrete Ill-posed Problems by Means of the L-curve[J]. SIAM Review, 1992, 34(4): 561–580. DOI:10.1137/1034115
[19] HANSEN P C, O'LEARY D P. The Use of the L-curve in the Regularization of Discrete Ill-posed Problems[J]. SIAM Journal on Scientific Computing, 1993, 14(6): 1487–1503. DOI:10.1137/0914086
[20] FAN Qingbiao, XU Caijun, YI Lei, et al. Implication of Adaptive Smoothness Constraint and Helmert Variance Component Estimation in Seismic Slip Inversion[J]. Journal of Geodesy, 2017, 91(10): 1163–1177. DOI:10.1007/s00190-017-1015-0
[21] 鲁洋为, 王振杰. 用U曲线法确定岭估计中的岭参数[J]. 导航定位学报, 2015(3): 132–134, 138.
LU Yangwei, WANG Zhenjie. Determining the Ridge Parameter in Ridge Estimation Using U-curve Method[J]. Journal of Navigation and Positioning, 2015(3): 132–134, 138.
[22] 袁强强, 沈焕锋, 李平湘, 等. 自适应正则化多幅影像超分辨率重建[J]. 中国图象图形学报, 2010, 15(12): 1720–1727.
YUAN Qiangqiang, SHEN Huanfeng, LI Pingxiang, et al. Adaptively Regularized Muti-frame Image Super-resolution Reconstruction[J]. Journal of Image and Graphics, 2010, 15(12): 1720–1727.
[23] XU Peiliang. Truncated SVD Methods for Discrete Linear Ill-posed Problems[J]. Geophysical Journal International, 1998, 135(2): 505–514. DOI:10.1046/j.1365-246X.1998.00652.x
[24] JONSSON, ZEBKER H, SEGALL P, et al. Fault Slip Distribution of the 1999Mw 7.1 Hector Mine, California, Earthquake, Estimated from Satellite Radar and GPS Measurements[J]. Bulletin of the Seismological Society of America, 2002, 92(4): 1377–1389. DOI:10.1785/0120000922
[25] 王乐洋, 许才军, 张朝玉. 一种确定联合反演中相对权比的两步法[J]. 测绘学报, 2012, 41(1): 19–24.
WANG Leyang, XU Caijun, ZHANG Chaoyu. A Two-step Method to Determine Relative Weight Ratio Factors in Joint Inversion[J]. Acta Geodaetica et Cartographica Sinica, 2012, 41(1): 19–24.
[26] WANG Leyang, XU Guangyu. Variance Component Estimation for Partial Errors-in-variables Models[J]. Studia Geophysica et Geodaetica, 2016, 60(1): 35–55. DOI:10.1007/s11200-014-0975-2
[27] WALTERS R J, ELLIOTT J R, D'AGOSTINO N, et al. The 2009 L'Aquila Earthquake (Central Italy):A Source Mechanism and Implications for Seismic Hazard[J]. Geophysical Research Letters, 2009, 36(17): L17312. DOI:10.1029/2009GL039337
[28] 王乐洋, 李海燕, 温扬茂, 等. 地震同震滑动分布反演的总体最小二乘方法[J]. 测绘学报, 2017, 46(3): 307–315.
WANG Leyang, LI Haiyan, WEN Yangmao, et al. Total Least Squares Method Inversion for Coseismic Slip Distribution[J]. Acta Geodaetica et Cartographica Sinica, 2017, 46(3): 307–315. DOI:10.11947/j.AGCS.2017.20160212
[29] 李海燕.震源位错模型参数反演方法研究[D].南昌: 东华理工大学, 2016.
LI Haiyan. Research on Inversion Method for Hypocenter Parameter Base on Dislocation Model[D]. Nanchang: East China University of Technology, 2016. http://cdmd.cnki.com.cn/Article/CDMD-10405-1016758485.htm
[30] ANZIDEI M, BOSCHI E, CANNELLI V, et al. Coseismic Deformation of the Destructive April 6, 2009 L'Aquila Earthquake (Central Italy) from GPS Data[J]. Geophysical Research Letters, 2009, 36(17): L17307. DOI:10.1029/2009GL039145
[31] CHELONI D, D'AGOSTINO N, D'ANASTASIO E, et al. Coseismic and Initial Post-seismic Slip of the 2009 Mw6.3 L'Aquila Earthquake, Italy, from GPS Measurements[J]. Geophysical Journal International, 2010, 181(3): 1539–1546.
[32] 王乐洋, 高华, 冯光财. 2016年台湾美浓Mw 6.4地震震源参数的InSAR和GPS反演[J]. 地球物理学报, 2017, 60(7): 2578–2588.
WANG Leyang, GAO Hua, FENG Guangcai. InSAR and GPS Inversion for Source Parameters of the 2016 Mw6.4 Meinong, Taiwan Earthquake[J]. Chinese Journal of Geophysics, 2017, 60(7): 2578–2588.
[33] FENG Guangcai, LI Zhiwei, SHAN Xinjian, et al. Geodetic Model of the 2015 April 25 Mw 7.8 Gorkha Nepal Earthquake and Mw7.3 Aftershock Estimated from InSAR and GPS Data[J]. Geophysical Journal International, 2015, 203(2): 896–900. DOI:10.1093/gji/ggv335
[34] HUANG M H, TUNG H, FIELDING E J, et al. Multiple Fault Slip Triggered Above the 2016 Mw 6.4 Meinong Earthquake in Taiwan[J]. Geophysical Research Letters, 2016, 43(14): 7459–7467. DOI:10.1002/2016GL069351
[35] LEE S J, YEH T Y, LIN Y Y. Anomalously Large Ground Motion in the 2016 ML 6.6 Meinong, Taiwan, Earthquake:A Synergy Effect of Source Rupture and Site Amplification[J]. Seismological Research Letters, 2016, 87(6): 1319–1326. DOI:10.1785/0220160082
http://dx.doi.org/10.11947/j.AGCS.2018.20170724
中国科学技术协会主管、中国测绘地理信息学会主办。
0

文章信息

王乐洋,赵雄
WANG Leyang, ZHAO Xiong
地震同震滑动分布反演平滑因子的确定
Determination of Smoothing Factor for the Co-seismic Slip Distribution Inversion
测绘学报,2018,47(12):1571-1580
Acta Geodaetica et Cartographica Sinica, 2018, 47(12): 1571-1580
http://dx.doi.org/10.11947/j.AGCS.2018.20170724

文章历史

收稿日期:2017-12-20
修回日期:2018-08-26

相关文章

工作空间