智利西部沿海处于纳斯卡板块和南美板块之间，由于纳斯卡板块以约8 cm·a^-1的收敛速度向南美板块俯冲(Angermann et al.，1999)，该地区大地震频发.2014年4月1日23∶46(UTC)在智利西北部Iquique近海地区发生了震级为M_W8.1的地震，震中位为19.642°S，70.817°W(图 1，红色五角星位置)，震源深度25 km(USGS National Earthquake Information Center-NEIC: http://earthquake.usgs.gov/regional/neic/).地震发生后，Iquique地震的震源性质和发震构造吸引了许多学者的关注和研究(An et al.，2014；Lay et al.，2014；Ruiz et al.，2014；Schurr et al.，2014；Kato et al.，2014；Yagi et al.，2014；Bürgmann et al.，2014；Hayes et al.，2014；Liu et al.，2015).国际上一些著名的地震科研机构和学者采用不同的数据和方法计算得到Iquique地震的震源机制解，见表 1.其中震源机制解1是GCMT(The Global Centroid-Moment-Tensor Project)通过计算单位矩张量元的格林函数，将格林函数与震源时间函数做卷积，运用最小二乘求得最优矩张量解.GCMT采用体波(截断周期为50 s)，长周期面波(截断周期为50 s)和地幔波(截断周期为200 s)计算了Iquique地震的震源机制解.震源机制解2、3是美国地质调查局(USGS)的结果，USGSCMT(USGS Centroid Moment Tensor Solution)采用体波和长周期面波反演震源机制，USGSWPHASE(USGS WPhase Moment Tensor Solution)采用W-phase(100~1000 s)的波形资料反演震源机制.对于中级地震，USGSBODY(USGS Body-Wave Moment Tensor Solution)采用远场台站记录的体波数据反演震源机制，比如2013年的芦山地震.单纯用体波数据反演的标量地震矩偏小(Jiang et al.，2014a).通常情况下USGSBODY不提供大地震的震源机制解结果，如2011年日本Tohoku M_W9.1地震等，包括此次Iquique M_W8.1地震.Yagi等(2014)利用48个台站的远震P波数据，考虑了格林函数的不确定性反演了Iquique地震的震源机制解.Lay等(2014)利用远程P波、SH波及海啸数据进行有限断层反演，得到的标量地震矩为1.7×10²¹ N·m.An等(2014)利用3个深海评估和海啸台站的海啸数据，反演得到了Iquique地震的滑动分布，并得到该地震的标量地震矩为1.22×10²¹ N·m.从上述研究结果来看，不同的研究机构和学者给出的Iquique地震的震源机制解相差较大，特别是标量地震矩M₀之间的差异显著，最大值约为最小值的2倍，而这种差异不容忽视.

图 1

Fig. 1

图 1 震中位置和16个连续重力观测台站的位置 Fig. 1 The locations of the hypocenter and 16 continuous gravity stations

表 1(Table 1)

表 1 智利Iquique地震震源机制解一览表 Table 1 Focal mechanism solutions of Iquique earthquake 序号 走向(°) 倾角(°) 滑动角(°) 深度/km MW M0/N·m 1 355 15 106 21.6 8.1 1.9×1021 2 340 14 74 33 8.1 2.1×1021 3 358 12 107 25.5 8.2 2.35×1021 4 350 15 106 22 8.1 1.5×1021 5 357 18 109 21.6 8.1 1.7×1021 6 350 15 107 20.1 8.0 1.22×1021 注：文献来源 震源机制解1为Global-CMT Project的结果；震源机制解2和3分别为USGS提供2个不同的结果：USGSCMT、USGSWPHASE；4为Yagi et al.(2014)的结果；5为Lay et al.(2014)的结果；6为An et al.(2014)的结果.

表 1 智利Iquique地震震源机制解一览表 Table 1 Focal mechanism solutions of Iquique earthquake

地震在瞬时释放巨大能量、产生大面积地表破裂的同时，激发大量的地震波，并且导致地球整体的震动.地球的自由振荡一旦被激发，将以驻波形式穿过地球的内部介质.由于自由振荡能量的大小与震源的破裂方式和破裂程度密切相关(Geller and Stein，1979；Stein et al.，2005)，因此地震后的自由振荡信号可用于约束地震震级和检验地震的震源机制解，它也是目前可对震源机制解进行总体检验和约束的非常有效的方法.利用自由振荡信号约束M_W9.0级以上大地震的震源机制解已经得到了不错的效果，比如2004-12-26苏门答腊大地震后，利用应变仪、长周期地震仪及超导重力仪记录到的球型低频模态₀S₀、₀S₂、₀S₃、₀S₄等，径向模₀S₀、₁S₀等与理论模拟结果对比，发现相应的标量地震矩是GCMT结果的2.5~2.67倍，并且得出断层的破裂长度、持续时间及破裂速度(Stein et al.，2005；Park et al.，2005；Lambotte et al.，2006).2011-03-11Tohoku大地震后，薛秀秀等(2012)利用₀S₀振型对Tohoku大地震的多个震源机制解进行分析和约束.对于中级地震来讲低频态的信号不易提取，Jiang等(2014)根据2013-04-20芦山地震的4种不同震源机制解计算了芦山地震频段为2.3~5 mHz的自由振荡信号，发现利用远场体波反演的地震矩偏小.Hu等(2014)为了评估芦山地震标量地震矩的正确性，对比了17个台站和5个矩张量解的观测值与模拟值，也得到了同样的结论.目前，利用自由振荡信号约束震源机制解所用到的观测数据均来自超导重力仪(SG)和宽频带地震仪，当震级较大时，SG会发生限幅情况，因此在计算过程中必须去掉同震时的观测段，对于M_W8.1 Iquique地震来说，中国大陆的几台超导重力仪均出现了限幅情况，而宽频带地震仪主要用于记录高频信息.

通过数字地震观测网络、中国大陆构造环境监测网络、中国地震背景场等国家重大科学工程项目，在中国大陆建成了一个集数据的自动汇集与数据库式管理，仪器的远程监控和数据的图形化分析处理为一体的连续重力观测台网，其平均空间分辨率约400 km.观测仪器包括PET/gPhone、DZW、GS15、TRG-1和GWR超导等，除了SG以外，这些弹簧重力仪也可以观测到微小的重力变化，分布广泛且连续的弹簧重力仪观测数据为检测地球自由振荡信号及利用自由振荡信号约束震源机制解提供了很好的数据基础.本文选取连续重力观测台网中分布广泛的16个弹簧重力仪观测数据(如图 1)，利用自由振荡信号观测约束Iquique地震的震源机制解，并有效的评估哪一个震源机制解最符合地震实际释放的能量，基于符合较好的重力站观测资料，定量研究和分析影响自由振荡振幅的因素.

地球自由振荡中的球型振荡不仅有横向分量，还有径向分量，因此会引起地球的体积变化，从而导致重力发生变化.中国连续重力观测台网提供采样率为1 Hz的重力观测数据用于地球动力学研究.Iquique地震发生后，连续重力观测台网的重力仪清晰的记录了重力变化.提取自由振荡信号之前，首先对各台站的观测资料进行预处理，利用Tsoft(Vauterin，1998)重力潮汐数据预处理程序，去掉尖峰、突跳、间断等.采用Eterna调和分析软件，计算大气导纳因子，并在重力观测中扣除大气的影响(Wenzel，1996).

r₀处的点源地震激发台站的简正模位移可以利用双力偶地震矩心矩张量理论和各项同性的地球模型模拟计算.在r处的位移u(r，t)为(Gilbert，1973)：

地球自由振荡计算值的振幅主要取决于地震发生的位置和地震的震源机制解.模拟计算台站自由振荡需要的参数有：地震发生的时间、位置、走向、倾角、滑动角、标量地震矩、台站的位置等.在参考的震源机制解中，地震的震源机制模型将震源等效为双力偶点源，根据双力偶地震矩心矩张量理论，由地球格林函数对矩张量加权求和，可计算矩张量源激发产生的地球自由振荡位移(江颖等，2014b，2015).

将地震激发的地球自由振荡信号的模拟值与相应的观测值相比对，可评估震源机制解是否能正确反映地震的实际破裂情况.由于不同的地球模型对自由振荡振幅的影响小于0.7%(薛秀秀等，2012)，几乎可以忽略.基于PREM地球模型(Dziewonski et al.，1981)，结合震源机制解，将格林函数对矩张量加权求和，可以得到接收台站的合成地震图(Geller and Stein，1979).

本文选取分布广泛的16个弹簧重力仪观测数据进行分析，台站包括：北京(BJ)、福州(FZ)、高台(GT)、格尔木(GEM)、海拉尔(HLE)、拉萨(LS)、兰州(LZ)、牡丹江(MDJ)、凭祥(PX)、琼中(QZ)、狮泉河(SQH)、炭山(TS)、乌加河(WJH)、乌什(WS)、西安(XA)、漳州(ZZ)，见图 1.为了避免环境噪声在简正模振幅上的影响，选取16个台站采样间隔为1 s，时间长度为16h的重力观测数据，观测数据和模拟数据在经过汉宁窗处理后进行离散傅里叶变换获取线性振幅谱.连续观测的弹簧重力仪会存在较大的零漂，选取的数据长度为16 h，零漂的影响较小，在数据处理过程中已经将该影响去掉.图 2显示了16个重力台站重力仪的自由振荡观测值和基于6个震源机制解计算的模拟值在1.5~5.3 mHz频段的结果.在选取的频段内，弹簧重力仪测得的自由振荡信号有较高的信噪比.从图 2可以明显看出基于震源机制解6的模拟值明显偏小，而基于GCMT和USGS的震源机制解模拟计算的自由振荡信号和观测值符合的较好.

图 2

Fig. 2

图 2 在1.5~5.3 mHz频段内16个弹簧重力仪观测数据与基于6个震源机制解算的模拟值 Fig. 2 Simulations and observations in 1.5~5.3 mHz frequency band made by 16 continuous gravity stations

为了更详细的反映比较结果，基于16个重力仪台站选用1.5~5.3 mHz频段中所有的球型自由振荡的模拟值和观测值进行比较，利用残余方差F(不符合度)和比例因子S(符合度)作为判断标准：

这里O_i和S_i是观测模的振幅和模拟值的振幅，n是所有振型的数量.当比例因子接近1的时候，其震源机制解对应的M₀是可靠的，可以描述震源的大小和强度.

另外，计算了基于各震源机制解结果的F和S的平均值及标准偏差，如表 2.由表 2可知，基于震源机制解3的自由振荡模拟结果与实测结果最接近，其不符合度F的平均值为0.090，比例因子S为0.992，最接近1，F和S的标准偏差也最小，分别为0.055和0.065，最能反映地震释放能量的大小.基于震源机制解1和2的模拟值与观测值的残余方差也较小，分别为0.098和0.091，显示了较高的符合度.但是，基于震源机制解6的模拟值与观测值相差较大，其不符合度F的平均值为0.214，比例因子S达到1.610.震源机制解中各个参数对地球自由振荡信号振幅的影响需要进行进一步的计算与分析.

表 2(Table 2)

表 2 基于6个震源机制解的F与S的平均值和标准差 Table 2 The average and the standard deviation of misfit factors and scaling factors 序号 F σF S σS 1 0.098 0.068 1.078 0.116 2 0.091 0.061 1.024 0.070 3 0.090 0.055 0.992 0.065 4 0.145 0.085 1.337 0.210 5 0.099 0.067 1.054 0.121 6 0.214 0.088 1.610 0.285

表 2 基于6个震源机制解的F与S的平均值和标准差 Table 2 The average and the standard deviation of misfit factors and scaling factors

分别计算震源机制解(标量地震矩、走向、倾角、滑动角和深度)对自由振荡振幅的影响.基于USGSWPHASE的震源机制解计算的模拟值与观测值最符合，由图 3的结果发现，16个台站中格尔木台站观测到的自由振荡信号较接近模拟值，其对应的F为0.047，S为0.996.因此，基于格尔木台站接收到的1.5~5.3 mHz频段的自由振荡信号平均振幅做定量分析.

图 3

Fig. 3

图 3 16个台站模拟值和观测值的残余方差F(A)和比例因子S(B) Fig. 3 Misfit factors (A) and scaling factors (B) at 16 continuous gravity stations

基于USGSWPHASE的震源机制解结果，分别给定标量地震矩、走向、倾角、滑动角和深度一个变化范围，计算不同的震源机制解因素对自由振荡振幅的影响.根据表 1，标量地震矩的变化范围为1.22×10²¹~2.35×10²¹ N·m，保证走向、倾角、滑动角和深度不变，分别计算不同的标量地震矩对应的自由振荡信号平均振幅，结果如图 4a，随着标量地震矩的不断增加，自由振荡的平均振幅也在明显增大，增量达到3.33×10³ nms^-2/Hz；走向的变化范围为240°~358°，保证M₀、倾角、滑动角和深度不变，分别计算不同的走向对应的自由振荡平均振幅，结果如图 4b，随着走向增大，自由振荡平均振幅也在慢慢增大，增加了2.25×10³ nms^-2/Hz；倾角的变化范围为12°~18°，保证M₀、走向、滑动角和深度不变，分别计算不同的倾角对应的自由振荡平均振幅，结果如图 4c，随着倾角的不断增加，振幅在慢慢减小，减小了3.08×10³ nms^-2/Hz；滑动角变化范围为74°~109°，保证M₀、走向、倾角和深度不变，分别计算不同的滑动角对自由振荡平均振幅的影响，结果如图 4d，随着滑动角的不断增加，振幅也在慢慢增大，增量达到2.13×10³ nms^-2/Hz；深度变化范围为20~33 km，保证M₀、走向、倾角和滑动角不变，分别计算不同的深度对应的内核平动振荡振幅的最大值，结果如图 4e，在20~33 km的范围内，振幅有小范围波动，差异仅为0.51×10³ nms^-2/Hz.通过以上数值实验说明M₀对自由振荡振幅的影响最大，而走向、倾角、滑动角和深度对内核平动振荡振幅也有一定影响，但是影响相对较小.因此，将6个震源机制解中的M₀均设置为相同的数值2.35×10²¹ N·m，模拟计算格尔木台站接收到的自由振荡信号，结果如图 5，当M₀为2.35×10²¹ N·m时，基于6个震源机制解得到的自由振荡信号相差不大，且和观测值均符合的较好.进一步证明了，M₀是影响自由振荡振幅的最大因素.An等(2014)仅利用海啸数据反演得到的震源机制解中M₀偏小是造成模拟值和观测值偏差较大的主要原因.

图 4

Fig. 4

图 4 震源机制解对自由振荡振幅的影响 Fig. 4 The influences of free oscillation′s amplitudes caused by focal mechanism solutions

图 5

Fig. 5

图 5 格尔木台重力观测数据与基于改进的6个震源机制解解算的模拟值(其中，M0均为2.35×1021 N·m) Fig. 5 Observations and modified simulations in 1.5~5.3 mHz frequency band made by GEM gravity stations (M0 is supposed 2.35×1021 N·m)

通过比较分析Iquique地震激发的自由振荡信号的实测值与模拟值，结果表明An等(2014)等仅利用海啸数据反演得到的震源机制解中M₀偏小是造成模拟值和观测值偏差较大的主要原因，说明对于本次的Iquique地震仅用海啸数据做反演，有一定的局限性，这是由于海啸数据对浅层地震的分辨率较高，而对于相对较深的逆冲型地震做数据分析的时候，并不一定能得到最佳解.利用海啸数据反演震源机制主要是快速的得到破裂面如何破裂，在哪个方向上会有更大的破坏，而对于标量地震矩的估算会有一定偏差.另外，Lay等(2014)从112个全球地震台站选取200~1000 s的周期范围的160个地面运动记录的W-phase数据反演了Iquique地震，得到的标量地震矩为1.99×10²¹N·m，而另外的一组结合海啸数据得到的震源机制解M₀的为1.7×10²¹ N·m，海啸数据的加入造成了反演的标量地震矩结果偏小.Yagi等(2014)仅采用了远震P波数据，反演结果同样偏小，可能是由于远场体波记录不到一些缓慢破裂的过程，且远场体波容易受到其他震相的污染，特别是在主震后发生大量余震的情况下.对于大地震，USGS的地震档案中一般不提供利用体波反演的震源机制解.USGSWPHASE选取长周期(周期为100~1000 s)方位角为5°~20°的近场数据和30°~90°的远场数据联合反演得到的震源机制解最佳.因此，利用长周期面波信号可以较好的约束Iquique地震的标量地震矩.

本文基于地球自由振荡分析了Iquique地震的震源机制解.结果表明利用弹簧重力仪记录的1.5~5.3 mHz球型简正模观测可以有效地对Iquique地震的震源机制解进行总体评估.从研究结果来看，地震的断层走向、倾角、滑动方向角和震源深度对自由振荡的振幅有一定影响，但是影响较小，M₀对自由振荡振幅的影响最大.另外，基于USGSWPHASE反演得到的震源机制解的自由振荡模拟值与观测值符合最好，其反演的标量地震矩可以较好地反映Iquique地震释放的能量，而仅利用海啸数据反演的标量地震矩偏小，联合近远场观测数据反演可显著改善震源机制解.因此，利用自由振荡的方法可以较好地判断地震震源机制和有效约束地震的震级.

另外，地球自由振荡观测只能对震源机制解进行总体综合评估，而不能单独对各个参数的影响单独评估.如果能利用其他大地测量手段，如近场GPS观测约束震源倾角、滑动方向角，结合地球自由振荡的观测和模拟可更细致地评估震源机制解.