利用地壳位移场观测资料,运用弹性力学理论可以得到地壳应变场。连续清晰的水平应变场对揭示地震孕育过程,应变场的中短期前兆异常演化特征,中强地震的发震地点、时间和震级的确定均有积极作用(阿提内·约克塞尔,2001)。目前研究地壳应变场的方法主要分为2类:第一类方法是对岩石的压力测试和数学模拟方法。陈颙等(1984)在研究受压岩石断层面形成过程中变形场的时空规律时发现岩石表面变形场压力增加、应力应变由分散到逐步集中的变化特征; 马瑾等(1979)在利用数字模拟方法研究断层闭锁区附近应变场演化中,分析总结了震源体周围震前应变场的演化特征,并对震中位置、发震时间和震级分别给出了判断方法。第二类方法是运用GPS观测资料开展地壳应变场前兆异常特征研究。自20世纪80年代以来,随着中国大陆GPS观测网的建立,运用GPS观测资料开展地壳应变场前兆异常特征研究在实践和理论方面都有了很大发展,逐步形成了对我国大陆内部较大构造块体形变和应变的监测能力。在大量研究结果中,不同研究者给出了一些GPS地壳运动与变形图像。张希和江在森(1998,1999)通过最小二乘配置法建立中国大陆水平运动速度场模型,获得了基于连续介质假设的中国大陆水平应变场初步结果。王琪等(1998)运用红河断裂的GPS监测资料计算了该区域构造应力场; 江在森等(2008)以华北为例研究了区域应变时空特征及其与强震关系,发现中强地震通常发生在剪切应变率的高值区及其边缘,尤其是与构造变形背景相一致的剪应变率高值区。早期主要是使用GPS水准观测资料,在中强地震地壳应变前兆异常特征研究领域取得了许多成果(江在森,2003a,b,c,2000,2001,2006; 王琪等,1998,戴黎明等,2010),这些成果对中强地震危险区的判定起到了不可或缺的作用。虽然运用GPS观测资料可以得到连续空间应变场结果,但由于GPS连续观测站台网稀疏,区域站复测周期为2年,因此在时间上存在不连续的缺陷。本文试图运用在时间连续性和空间分布密度都优于GPS基准站的青海及其周边数字地震台网记录的地脉动资料,借鉴GPS观测资料计算水平应变场的方法开展青海及周边区域应变场的计算,给出2016年1月21日门源6.4级地震前20 d到震前1 h的水平应变场,对门源6.4级地震的应变场时空前兆异常特征进行初步研究。

本文选取了青海区域台网中的67个地震台(图1)2016年1月1～21日每天0时所记录的地脉动资料。通过粗选去掉了断记较多和标定有问题的数据,并去除了仪器响应,对干扰信号进行了滤波处理。采用中离散DB(Daubechies)小波分析软件对观测数据进行了分析处理,小波阶数为10,小波系数为4(图2)。预处理之后,进行功率谱计算分析。以信号强度大小确定参与区域应变场的计算的数据,最终选定0.01～0.15 Hz频段的数据计算区域应变场。

以2016年1月4日00～01时AXX台地脉动资料为例,通过小波分析,从分出的10阶地脉动信号,分别计算功率谱,通过功率谱值大小,确定地脉动信号强度最大的频率范围,从而确定带通滤波的上下限频率。AXX台计算结果功率谱最大频率范围是0.01～0.15 Hz(图2i-2)。

通过小波分析分解观测数据,运用经验协方差函数,得到网格化(0.3°×0.3°)地脉动位移空间相关的水平位移场,再利用面膨胀计算公式得到水平应变空间分布。按时间先后为序,给出门源6.4级地震前20 d(每天00时)(图4)和震前14 h(图5)青海及周边区域应变场时间、空间和强度的演化进程。震前20 d时间演化进程总体分为2个阶段:第一个阶段为2016年1月1～14日; 第二个阶段为2016年1月15～21日。第一阶段应变场空间分布特征主要表现为张性(红色)与压性(蓝色)区域四象限或相间分布,且主要集中分布在青海北部和东南部,青海西南部应变仅有一些微弱的变化。高应变梯度带和四象限张压相间应变场出现频率较高的区域主要集中在青海及周边区域二级活动块体(张培震等,2003)边界带上的阿尔金断裂带、祁连断裂带、东昆仑断裂带以及玉树甘孜断裂带在青川两省交界的东段(图3a,图4)。

2016年门源6.4级地震震中位于GPS给出的应变梯度带附近,表明中长期应变场在地震前有一定的前兆反映(图3a)。对比图3a,b可以看出,研究区内大部分6级以上地震与主要断裂空间分布(张培震等,2003)和应变场梯度变化分布基本一致。地脉动解算的应变场空间分布格局反映研究区域内短期地壳应变场的活动特征,其反映出高应变值在空间上的变化特征与使用中国地震局第一形变测量中心给出的GPS速度场解算的中长期应变场(图3a)、应变梯度带空间较为一致,

图3 2013～2015年GPS资料得到应变场(a)和二级活动块体及主要断层空间分布(b)
Fig.3 Strain field obtained by GPS data between 2013 and 2015(a)and spatial distribution of two level active blocks and main faults

图4 2016年1月1～20日每天00时应变场空间分布(a)1月1日;(b)1月2日;(c)1月3日;(d)1月4日;(e)1月5日;(f)1月6日;(g)1月7日;(h)1月8日; (i)1月9日;(j)1月10日;(k)1月11日;(l)1月12日;(m)1月13日;(n)1月14日;(o)1月15日; (p)1月16日;(q)1月17日;(r)1月18日;(s)1月19日;(t)1月20日
Fig.4 Spatial distribution of strain field at 00 a.m.from Jan.1 to Jan.20,2016

图5 2016年1月20日11时～21日00时整点应变场空间分布(a)1月20日11时;(b)1月20日12时;(c)1月20日13时;(d)1月20日14时;(e)1月20日15时; (f)1月20日16时;(g)1月20日17时;(h)1月20日18时;(i)1月20日19时;(j)1月20日20时; (k)1月20日21时;(l)1月20日22时;(m)1月20日23时;(n)1月21日00时
Fig.5 Spatial distribution of hour strain field from 11 a.m.on Jan.20 to 00 a.m.on Jan.21,2016

即地脉动所得出的应变场异常空间分布与中长期(活断裂和GPS资料)所指示的地震孕震区在空间上有一定的相关性。在第一阶段,门源6.4级地震震中附近应变场在空间上分布有张压相间的四象限分布的特征,且应变最大值的区域始终在震中附近。从图4a～m可以看出,1月1～13日,研究区域内的应变强度虽有起伏波动,但变化不大; 1月14日研究区域内应变强度升高达到最大值(图4n),升高幅度达到背景值的4倍(图6),此时为门源地震前6天。第二阶段,1月15～20日,研究区内高应变梯度带进一步向震中附近收缩,且四象限分布更为明显(图4o-t)。从图5可以看出,门源6.4级地震前14 h的应变场空分布,结果显示应变强度略有回升,高应变梯度带主要集中在震中附近。门源地震区域应变场强度在时间进程上特征,表现为2015年1月1～13日研究区域内应变场强度变化较为平缓且应变场强度较低,1月14日应变场强度突然升高,达到背景值4倍,1月15～21日应变场强度逐步降低至震前水平(图6)。

图6 2016年1月1～21日应变强度时间序列
Fig.6 Time series of strain intensity from Jan.1 to Jan.21,2016

本文通过小波分析提取青海区域数字地震台网记录的地脉动观测资料低频信号,计算门源6.4级地震震前20 d地壳应变场。从区域应变场空间分布、应变场强度随时间变化特征,分析门源6.4级地震前前兆异常特点。从区域应变场的角度,总结对该地震的发震地点和发震时间的短期预测方法。研究结果显示,对门源6.4地震前发震地点的判定,在应变场空间上表现为在震中附近出现张压相间应变梯度带,而且张性与压性区域的应变强度大致相当。这一结果与陈颙等(1984)利用脉冲激光全息干涉术研究岩石受压破裂过程所观测到的岩石变形特征一致; 马瑾等(1995)在研究演示受压过程中,应变场也有类似的张压相间四象限分布特征。发震时间的判定在应变强度随时间进程上表现为,随着发震时刻的临近,区域应变场在震前6 d左右升高,大约是背景值的4倍,随后恢复到高值前的变化水平。

利用区域数字地震台网地脉动观测资料计算得到的区域应变场,其高应变梯度带在空间分布上与GPS资料给出的中长期应变场基本一致,其张压性有时也不一致,这可能与两种观测资料所反映的地壳运动周不同有关,这方面的问题还有待于做进一步研究。

本文利用相对于其它形变观测资料,在时间上连续性好、采样率高、台网分布密相对度较高的区域数字地震台网地脉动观测资料获取应变场,是对GPS观测资料获取应变场的有益补充,对跟踪判定中强地震发震地点和发震时间具有一定意义。