2017, Vol. 37
2. 中国地质大学(武汉)地球物理与空间信息学院，武汉市鲁磨路388号，430074;
3. 云南省地震局，昆明市知春街249号，650041
川滇菱形块体由西南边界的红河断裂、北东边界的鲜水河断裂、安宁河-则木河断裂及小江断裂所围成，北西向的曲江断裂(即通海-峨山断裂)位于小江断裂与红河断裂交汇部的西北侧。自有记载以来，曲江断裂及其附近地区发生5.0级以上地震10多次，其中1913年峨山7.0级地震、1970年通海7.7级强震相隔57a，而震中相距仅38 km(表 1)。
|
表 1 1913年以来曲江断裂5级以上地震目录 Tab. 1 M≥5 earthquakes catalog in Qujiang fault since 1913 |
调查发现，通海地震破裂带基本是连续的，走向北西60°，倾向以北东为主，倾角50°~85°，右旋错动为主兼有逆冲，较大的水平错动为2~2.7 m，逆冲量一般在0.5 m以下[1]。用测震资料获得的震源机制倾向基本为南西，阚荣举等[2]给出走向306°，倾向南西，倾角87°；张之立等[3]给出的倾向是212°；而陈颙[4]认为走向291°，倾向北东，倾角80°；王椿镛等[5]根据地形变资料反演的倾向是南西，倾角84°，平均走向滑距为-2.24 m，破裂长度为52 km；Zhou[6]认为，通海地震至少由相距较远的两次子事件组成，破裂长度90 km。
通海地震前震区布设有较密集的三角测量点，震后利用GNSS对这些三角点进行复测，可以监测到更多的近场形变资料。如果同震形变很大，震前三角测量误差引起的同震形变误差可以忽略，则可以对地震破裂滑动分布进行强有力的约束。该方法在2008年汶川地震破裂研究中发挥了很大的作用[7]。本文以1970年云南通海地震为例，一方面检验评估用GNSS对震前三角点进行复测在研究7级地震破裂中的作用及需要注意的问题，另一方面也希望能够通过对通海地震的研究，对1970年通海地震破裂获得一些新的认识。
1 观测资料和数据处理根据文献[8]，通海震区震前只有1958年为测图而施测的一期Ⅱ等补充网三角观测成果。1970年地震后复测了36条边，边长精度为10~20 cm。假设北盘远处两相邻点为不动点，按经典平差方法计算断裂两侧相对位移矢量，反映了曲江断裂两侧存在顺时针的水平错扭。错动范围沿曲江断裂、长度约50 km，两侧三角点相对位移为2m左右，离开断裂20 km仍有比较大的相对位移。但是，该形变结果含有不动点选取误差的影响。
王椿镛等[5]假设南北两侧远处的两个参考点具有一个确定的位移值，该位移值是在求解断层破裂模式的迭代过程中给出的。每次用破裂模型计算该理论位移值，并加到参考点的原坐标中，得到参考点的新坐标；然后按经典平差方法求解其他点的位移值。因而，最后的平差结果依赖于地震断层参数解。
刘玉权等[9]在经典自由网平差的基础上采用相似变换获得形变场，参考点不稳定性带来的误差能够得到比较好的消除和削弱，计算结果与王椿镛等[5]的结果极为相似。可以看出，曲江断裂发生了明显的相对右行错动，错距达2 m左右。断裂两盘所显示的应变特征相似、量级相当，说明两盘点位移应该是对称的。
图 1是早期三角测量结果，精度低，复测困难。我们于2007-01用GPS对震区20多个三角点进行复测, 每点观测不少于6~8 h, 数据采样率30 s。用BERNESE GNSS软件对GPS数据进行处理，获得在ITRF2005参考框架下的单日解坐标，属于WGS84坐标系统。而震前的三角点坐标属于北京54坐标系统，两个坐标系统需要进行七参数转换。从表 1地震目录可见，两期观测间没有较大地震；据文献[8]，通海地震震后形变在1 a多后逐渐平缓，曲江断裂两侧震后水平变形只有17 mm，远远小于其同震形变影响。因此，可以忽略两期观测间小地震同震形变以及震后形变影响。
|
图 1 根据三角测量计算的1970年通海地震同震形变与模拟值[5] Fig. 1 Measurement and modeling coseismic displacements derived from triangulation caused by the 1970 Tonghai earthquake[5] |
根据中国地壳运动速度场[10]扣除两期观测间的长期构造运动后，可以获得通海地震同震形变场。由于没有震前震后公共点控制，所以坐标参数转换需要作一些假定。经过实验发现，如果假定断裂同一侧(例如北侧)离开断裂相同(或相近)距离的两个相邻点具有某个理论变形，即以断裂同一侧的相邻两点作为参数转换控制点，则可以获得比较合理的同震形变场。该结果与文献[5, 9]的形变分布特征大致相同，但在断裂东南端的变形较大(图 2)。
|
图 2 根据三角测量和GPS测量计算的1970年通海地震同震形变与模拟值 Fig. 2 Measurement and modeling coseismic displacements derived from triangulation and GPS caused by the 1970 Tonghai earthquake |
一般情况下假设震源区为半无限弹性空间，则矩形断层位错引起的地表变形可用半无限空间的弹性位错模型[11]计算，主要与断层几何特征七参数(长、宽、深、走向、倾向、水平坐标)和滑动参数(走滑量、倾滑量、张性分量)有关，以形变观测值为约束可以反演断层破裂参数。假设破裂由少数几个面积比较大的断裂组成，则反演的目的是求得破裂的大体几何形状和平均意义的滑动量。如果将破裂细分为更多的子断层，则反演的目的是为了获得精细的破裂滑动分布。
两种反演方法都应该满足观测数据拟合度和滑动分布粗糙度最小：
|
(1) |
式中，d是形变观测值，主要是水平位移观测值；W是观测值的权矩阵，是观测值方差-协方差D的逆矩阵，即D = W－1W；G是形变格林函数，用半无限空间弹性位错模型[11]计算；s是待求的未知数矢量，一般指子断层滑动矢量。
假设地震断层只有一个矩形断层，我们首先使用非线性优化的模拟退火算法研究该矩形断层的平均长、宽、倾角、位置和大致平均滑动量。根据张四昌等[1]的研究，通海地震断层与曲江断裂基本吻合，完全产生于曲江断裂的破碎带范围之内。因此，地震断层的走向、地表迹线、埋深可以固定或者约束在较小的活动范围之内，主要研究地震断层的长、宽、倾角、位置和滑动量。根据地质调查的先验信息设定模拟退火算法的搜索区间，断层长度在20~100 km内搜索，宽度给予较大的活动范围如5~40 km，上边界埋深0~10 km，倾角30~89°，均衡走滑量-3~0 m，倾滑量-1~1 m。由于地质调查获得的破裂倾向有南西方向，也有北东方向，所以我们也设置南西和北东两个倾向。
分别以文献[5]和本文的同震形变结果为约束，采用模拟退火算法反演单断裂矩形断层破裂参数。两个反演结果都获得超过80°甚至接近89°的倾角，因而推测破裂倾向南西和北东都是有可能的，目前的约束数据可能还不足以分辨其倾向差异，也可能有些破裂倾向南西，有些破裂倾向北东。为了简便，后面将采纳大部分震源机制解给出的倾向南西的结果。破裂上边界埋深一般在0~2 km以内，可以认为破裂出露地表。两个反演结果获得的破裂长度都接近100 km。均衡倾滑量在0.5 m以内，走滑量一般大于1.5 m，本文的同震形变反演的均衡走滑量比文献[5]稍大，反映了破裂高倾角、走滑为主的特征。
3 地震破裂滑动分布基于地震破裂几何形状的研究结果，我们将地震破裂迹线固定为与地质调查的破裂重合，破裂倾角取84°，倾向南西。断层长度有意向两端延长变为140 km，下边界深度延长为26 km。将断层划分为3 km×3 km的315个子断裂，如此则式(1)中各子断裂的滑动量与地表形变之间是线性关系，可以用约束变量的最小二乘法来求解子断层滑动量。
在地震破裂滑动分布反演中，滑动分布平滑因子的选取会对最终破裂滑动分布产生显著影响。根据不符值和粗糙度的折衷(trade-off)曲线可以选取滑动分布平滑因子，选取的原则是不能不恰当地增加观测值不符值，一般可以选取折衷曲线拐点或转折点的平滑因子。但实际上很难准确可靠地选取平滑因子，所以可以采用第三方观测数据、其他先验信息或指标，如地表错动量，对反演结果进行有益的检验。本文观测数据较少，所以采用地质调查位错量，将破裂滑动约束在3 m以内，以增加平滑因子选取的合理性。
以文献[5]的同震形变为约束，反演获得的破裂滑动分布如图 3所示，以本文同震形变为约束反演获得的破裂滑动分布如图 4所示。两个结果的最大破裂都位于测震资料给出的震中附近，破裂长度、深度也大致相同。假定剪切模量为30 GPa，两个模型计算的矩震级Mw=7.4，破裂长度都接近100km，但是第一个结果显示的破裂幅度较小，而第二个结果显示的破裂幅度较大。两个结果的西北段破裂分布比较接近，而东南段差异较大。图 3中东南段破裂较小，而图 4中东南段破裂较大，与两个同震形变场的形变分布特征相对应。
|
图 3 以文献[5]同震形变为约束反演的破裂滑动分布 Fig. 3 Slip distribution model constrained by the coseismic deformations in document[5] |
|
图 4 根据GPS与三角测量同震形变为约束反演的破裂滑动分布 Fig. 4 Slip distribution model constrained by the coseismic deformations derived from triangulation and GPS |
两个反演结果差异的可能原因如下：
1) 尽管GNSS没有误差积累，三角测量误差较大，但是本文联合两者进行同震变形计算缺乏严密的坐标转换公共点约束。虽然以断裂北侧两个相邻点作为坐标转换控制点，并考虑了该两个控制点的形变理论值，但是可能还是含有一些系统误差，并在断裂东南段引起较大的误差，导致反演结果的误差。
2) 文献[5]可能低估了断裂东南端的变形，从而引起反演结果的差异。
3) 也许两个形变结果都含有一定的误差。
本文分析同震形变时根据文献[8]的形变量级，假设可以忽略通海震后变形的影响。但是通海震后变形还是存在的，并且其空间分布不均匀，靠近断层的测点影响明显。因此，如果要可靠评估通海地震的破裂滑动分布、震级大小以及曲江断裂的地震复发周期，还需要进行如下工作：
1) 收集更多的观测数据，例如采用GNSS对震区三角点进行更多的复测；
2) 应该在断裂北侧或南侧，甚至南北两侧几乎没有通海地震变形的地区，观测若干个三角点作为坐标转换控制点，并对GPS联合三角测量计算同震形变进行精度评估；
3) 可靠地评估震后变形对同震形变及破裂滑动分布的影响。
4 结语本文分别以三角测量获取的通海地震同震形变和联合GPS与三角测量获取的同震形变为约束数据，先采用模拟退火算法反演单断裂矩形断层破裂模型参数，然后用约束变量的最小二乘法反演地震滑动分布模型。结果表明，模拟退火算法获得的单矩形断裂模型的倾角都超过80°，甚至接近89°，破裂出露地表，破裂长度都接近100 km，破裂以走滑为主，走滑量一般大于1.5 m。以约束变量的最小二乘法反演获得的两个破裂滑动分布模型中，矩震级Mw=7.4，最大破裂都位于测震资料给出的震中附近，破裂深度也大致相同，两个结果的西北段破裂分布比较接近。但是，三角测量结果确定的东南段破裂幅度较小，而根据GPS和三角测量确定的东南段破裂幅度较大。其原因可能是，本文联合GPS和三角测量进行同震变形计算时，缺乏严格的坐标转换公共点约束，引起一些较大的误差；也可能文献[5]低估了断裂东南端的变形，从而引起反演结果的差异；也有可能是两个形变结果都含有一定的误差，甚至包含震后变形的影响。如果要可靠地评估通海地震的破裂滑动分布、震级大小以及曲江断裂的地震复发周期，还需要收集更多的观测数据，并在几乎没有通海地震变形的地区观测若干三角点作为坐标转换控制点，并对GPS联合三角测量计算同震形变进行精度评估，同时需要评估震后变形的影响。
| [1] |
张四昌, 刘百箎. 1970年通海地震的地震地质特征[J]. 地质科学, 1978(4): 323-335 (Zhang Sichang, Liu Baichi. Seismic Geological Characteristics of Tonghai Earthquake in 1970[J]. Scientia Geologica Sinica, 1978(4): 323-335)
(  0) |
| [2] |
阚荣举, 张四昌, 晏凤桐, 等. 我国西南地区现代构造应力场与现代构造活动特征的探讨[J]. 地球物理学报, 1977, 20(2): 96-109 (Kan Rongju, Zhang Sichang, Yan Fengtong, et al. Present Tectonic Stress Field and Its Relation to the Characteristics of Recent Tectonic Activity in Southwestern China[J]. Acta Geophysica Sinica, 1977, 20(2): 96-109)
( 0) |
| [3] |
张之立, 刘新美. 1970年1月5日云南通海地震的震源特征[J]. 地球物理学报, 1982, 25(5): 441-447 (Zhang Zhili, Liu Xinmei. The Focal Characteristics of the Tonghai Earthquake, Yunnan Province, on January 5, 1970[J]. Acta Geophysica Sinica, 1982, 25(5): 441-447)
( 0) |
| [4] |
陈颙. 测定浅震震源参数的宏观方法[J]. 地球物理学报, 1975, 18(4): 246-255 (Chen Yong. On the Determination of Focal Parameters of Shallow Earthquakes by Means of Macroseismic Method[J]. Acta Geophysica Sinica, 1975, 18(4): 246-255)
( 0) |
| [5] |
王椿镛, 朱成男, 刘玉权. 用地形变资料测定通海地震的地震断层参数[J]. 地球物理学报, 1978, 21(3): 191-198 (Wang Chunyong, Zhu Chengnan, Liu Yuquan. Determination of Earthquake Fault Parameter for the Tonghai Earthquake from Ground Deformation Data[J]. Acta Geophysica Sinica, 1978, 21(3): 191-198)
( 0) |
| [6] |
Zhou H, Allen C R, Kanamori H. Rupture Complexity of the 1970 Tonghai and 1973 Luhuo Earthquakes, China, from P-Wave Inversion, and Relationship to Surface Faulting[J]. Bulletin of the Seismological Society of America, 1983, 73(6): 1585-1597
( 0) |
| [7] |
Wang Q, Qiao X J, Lan Q G, et al. Rupture of Deep Faults in the 2008 Wenchuan Earthquake and Uplift of the Longmenshan[J]. Nature Geosci, 2011(4): 634-640
( 0) |
| [8] |
国家地震局地震测量队.1970年云南通海地震的地形变特征[J].地球物理学报, 1975, 18(4): 240-245 (The Geodetic Survey Brigade for Earthquake Research, National Seismological Bureau.The Characteristics of the Crustal Deformation Associated with the Tonghai Earthquake, Yunnan, in January 1970[J]. Acta Geophysica Sinica, 1975, 18(4): 240-245)
( 0) |
| [9] |
刘玉权, 黄震民, 杨来宝. 1970年通海地震的水平位移研究[J]. 地震研究, 1984, 7(2): 187-195 (Liu Yuquan, Huang Zhenmin, Yang Laibao. Study on the Horizontal Displacement Caused by the Tonghai Earthquake in 1970[J]. Journal of Seismological Research, 1984, 7(2): 187-195)
( 0) |
| [10] |
Zhao B, Huang Y, Zhang C H, et al. Crustal Deformation on the Chinese Mainland during 1998-2014 Based on GPS Data[J]. Geodesy and Geodynamics, 2015, 6(1): 7-15
( 0) |
| [11] |
Okada Y. Surface Deformation due to Shear and Tensile Faults in a Half-Space[J]. Bull Seismol Soc Amer, 1985, 75(4): 1135-1154
( 0) |
2. Institute of Geophysics & Geomatics, China University of Geosciences, 40 Hongshance Road, Wuhan 430074, China;
3. Yunnan Earthquake Administration, 249 Zhichun Street, Kunming 650041, China