1. 引言
形变学科是地震前兆监测和地震预报的主要手段之一，具有明确的物理含义以及预测机理，因构造地震孕育过程的应变积累和震时应变释放都伴随有地壳形变发生，通过形变观测能够获取震间长期应变积累、震前及临震应力应变状态等直接有效的地壳形变动态信息，在地震综合预测中，形观测手段发挥了重要的作用 [1]。同时，形变异常落实复杂性和异常甄别差异性造成很多异常无法判断。本文通过对泾源卧龙山垂直摆倾斜仪异常变化跟踪落实，排查仪器状态和环境干扰等影响因素，多角度、多方法分析研究异常资料，客观评价此次异常变化 [2]。
2. 台站背景资料
台站地理、地质构造情况
固原台泾源卧龙山观测点位于宁夏最南端的固原市泾源县城东南角卧龙山上，地理位置为35.5˚N、106.34˚E，海拔高度1890 m。泾源县地处六盘山东麓腹地，地势西北高，东南低，属低山丘陵区，海拔高度范围1608~2942米(图1(a))。位于六盘山下，泾水源头的泾源县，因泾河发源于此而得名，素有“秦风咽喉、关陇要地”之称。其东与甘肃省平凉市相连，南与甘肃省华亭县、庄浪县接壤，西与宁夏隆德县毗邻，北与宁夏原州区、彭阳县交界。气候区域属温带半湿润区，为森林草原类型气候，呈现出“春寒、夏凉、秋短、冬长”的特点。年均气温6.9℃，日照时数2370小时，年均降水量为641.5毫米。泾源县地形、地貌奇特。地势自西向东倾斜，地貌可分为侵蚀构造石山区，剥蚀构造丘陵区和侵蚀堆积河谷平川区三种类型。境内水资源以泾河为主，年均径流量2.4亿立方米。
地震地质构造上，泾源卧龙山测点位于宁夏南部弧形构造区的清水河断陷盆地内，其西侧为六盘山东麓断裂，东侧为云雾山断裂，两条断裂皆为NNW走向，相距不超过20 km，两条断裂分别为清水河盆地的东、西边界，测点以北40 km左右为已探明的牛首山–固原断裂。距离测点最近的断裂为六盘山东麓断裂，约4 km，该断裂沿六盘山东麓延伸，总体走向北北西，倾向南西西，向南止于陕西宝鸡东，长度近200 km (图1(b))。六盘山断裂带处在海原断裂带的南东端，其走向与海原断裂带有较大的夹角，晚新生代的构造变形以北北西至近北南走向分布的逆冲断裂及其伴随的褶皱为主要特征(丁国瑜等，1991)。由于处在青藏块体东北缘和鄂尔多斯西南缘交汇处，受青藏块体的NE-NEE向推挤，及鄂尔多斯块体的阻挡之下，区域内主要发育着北西西、北西、北北西以及近南北向的全新世活动断裂(祁连山东段–海原断裂、六盘山汇聚构造断裂、西秦岭北缘断裂等) [3]，构成了3个构造变形组合系统：海原–六盘山走滑–挤压缩短转换系统、陇山–宝鸡局部拉张系统和西秦岭构造系挤压旋转系统(张培震等，2003；郗欽文等，1994)。区域地质构造复杂，新构造和剥蚀作用下，基岩裸露，地貌上断层崖面和水系有规则拐弯，水平扭动表现明显。地层岩性虽然为第四系洪积沙土，但基岩完整质密 [4]。
(a) 
(b)
Figure 1. Topographic distribution map (a) and regional geological map (according to Deng Qidong et al.) (b) of Wolong Mountain in Jingyuan
图1. 泾源卧龙山地形分布图 (a)和区域地质图(据邓启东等) (b)
3. 异常概述
异常形态特征描述
4月28日起，固原台(泾源卧龙山)垂直摆倾斜(VP) NS向、EW向分别出现转折加速上升、转折加速下降变化。5月14日6时，NS向转折下降，自4月28日起累积上升幅度为1271.729 × 10−3''，上升了52.2%；5月15日20时，EW向由转折加速下降转为上升，自4月28日起累计下降幅度为1030.01 × 10−3''，上升了18.3%。5月22日22时，NS向再次转折上升，自5月15日累积下降幅度为316.093 × 10−3''，下降了27.2%；5月21日22时，EW向再次转折下降，自5月15日起累积下降幅度为732.415 × 10−3''，下降了15.9%。截至5月26日，NS向经历了上升–下降–上升的过程，累积上升幅度达1221.205 × 10−3''，上升了50.2%，EW向经历了下降–上升–下降的过程，累积下降幅度达705.64 × 10−3''，下降了15.5% (图2、表1)。
Figure 2. Long-term (Figure 1(a)) and short-term (Figure 1(b)) hourly value curves of vertical pendulum tilt in Guyuan (Wolong Mountain, Jingyuan)
图2. 固原(泾源卧龙山)垂直摆倾斜长期(图1(a))及短期(图1(b))整点值曲线
Table 1. Statistics of vertical pendulum inclination in Guyuan (Wolong Mountain, Jingyuan)
表1. 固原(泾源卧龙山)竖值摆倾斜变化统计
4. 异常调查
4.1. 异常产生原因的初步判断与调查步骤
固原(泾源卧龙山) VP垂直摆倾斜NS向、EW向自2021年4月28日起出现异常，经仪器厂家鉴定，仪器观测系统工作正常；同时，此次异常变化与VS垂直摆在2013年庐山7.0级地震前、2015年尼泊尔8.1级地震前、2017年九寨沟7.0级地震前出现的异常类似，且异常出现前周边环境无明显变化，5月22日发生青海玛多7.4级地震，初步认为此次异常为地震前兆异常，并应考虑其与青海玛多7.4级地震的关系 [5]。针对此次前兆异常，台站成立了异常核实小组，对现场进行核实调查，并制定了异常核实工作方案：① 观测系统工作状态检查；② 仪器标定记录和工作日志信息核实；③ 气象因素调查；④ 异常变化分析和性质判断。
4.2. 仪器工作状态检查
经过现场仪器工作状态检查，仪器供电电压稳定，指示灯和面板数据均正常，仪器数据采集及传输系统正常。观测仪器底脚螺丝不存在生锈。仪器标定正常，5月11日NS向标定失败可能与当时天气状况有关。武汉地震科学仪器研究院工程师5月23日对仪器调零装置进行检测，电极正常，打开调零装置进行检测重新安装后，远程调零功能恢复。仪器厂家工程师判定仪器观测系统工作正常，与此次VP摆异常变化无关 [6]。
查阅相关日志，无影响此次异常的干扰发生。由于格值会导致潮汐因子和数据趋势等出现相应的变化(闫伟等，2015)。通过日志整理及标定结果检查，工作人员每年6月和12月均会按时对仪器进行标定和巡检。图3为VP型垂直摆2018年12月~2020年12月NS向和EW向格值标定结果，5次标定格值变化不大，格值误差及校准幅度均符合《细则》要求。异常发生前未进行标定，使用格值为2020年12月28日的格值，NS向为0.06764，EW向为0.06625，与之前格值误差分别为−0.368%、−0.27%，均小于2%，符合《细则》要求。分析认为格值变化与此次异常无关 [7]。
Figure 3. The change of grid value in Guyuan (Jingyuan Wolong Mountain) from 2018 to 2020
图3. 固原(泾源卧龙山) 2018年~2020年格值变化情况
Figure 4. Long-term change curve of temperature, air pressure, cave temperature and rainfall in Wolong Mountain, Jingyuan
图4. 泾源卧龙山气温、气压、洞温、降雨量长期变化曲线图
Figure 5. Comparison of temperature, air pressure and VP vertical pendulum changes at the Wolong Mountain measuring point in Jingyuan
图5. 泾源卧龙山测点气温、气压与VP垂直摆变化对比图
4.3. 气象因素的调查
泾源卧龙山测点的气象三要素观测仪为ZKGD3000-M型气象三要素综合观测仪，其中降雨量仅2019年正常观测，气温、气压观测正常。图4(a)~(c)分别为2019年以来泾源卧龙山气温、气压、洞温曲线图，可以看出，从长期来看，气温、气压、洞温均呈现规律年变，2021年4月28日以来无显著异常变化；图4(d)给出了2019年VP垂直摆与降雨量对比图，可以看出，降雨量显著集中的3月、4月，VP摆仅呈现正常年变，无明显短期变化，且经调查，此次异常期间，固原并无显著降雨发生；分析认为，此次固原(泾源卧龙山) VP垂直摆异常变化与气温、气压变化无关。图5给出了2021年3月1日~3月31日气温、气压与VP摆NS向、EW向对比图，可以看出，从短期来看，VP垂直摆与气温正相关，与气压反相关；此次固原(泾源卧龙山) VP垂直摆NS向经历了上升–下降–上升的过程，EW向经历了下降–上升–下降的过程，气温、气压虽然与此次异常的大幅度变化无关，但可能与异常持续过程中的小幅度波动有关，下文将对气温、气压变化对VP摆的影响予以剔除 [8]。
5. 形变异常变化分析与性质判定
5.1. 观测资料质量监控
中国地震台网中心地球物理台网学科组每月对全国各台站连续率、完整率及精度进行统计，2021年1月~4月泾源VP垂直摆完整率、连续率均在99.9%以上，两个测项精度集中在0.05~0.06附近，精度不高，但1~4月无显著变化(表2)。
Table 2. Statistics on the continuity rate, completeness rate and accuracy of Jingyuan VP vertical pendulum tilt from January to April in 2021
表2. 泾源VP垂直摆倾斜2021年1~4月连续率、完整率、精度统计
5.2. 干扰因素排除
5.2.1. 卧龙山水库对VP垂直摆影响的数值模拟
基于均匀的、各向同性的半无限弹性体的点状载荷模型引起周边垂向位移场变化的解析解，闫伟等(2015)推导了水平倾斜场的解析解计算方法，给出了二维、三维不规则形状载荷对地面或地下某点的垂直位移和水平倾斜的解析解计算方法 [9]。本文使用三维不规则形状载荷模型，计算了卧龙山水库对泾源VP垂直摆的影响。通过91卫图助手获取水库边界经纬度、转换成成直角坐标系下坐标，绘制水库与台站的相对位置图。台站周边的岩性为第三系粉砂岩，取杨氏模量E = 4 × 107、μ = 0.25。通过数值模拟，计算了水库蓄水变化最大的时段——2020年8月18日至2021年3月15日，蓄水量由54万m3减至11万m3，高差变化为7 m时引起的泾源VP垂直摆变化分别为70.8132 × 10−3''、−101.2294 × 10−3''，远小于此次异常的变化幅度，认为卧龙山水库蓄水与此次异常无关 [10]。
5.2.2. 气象因素对VP垂直摆的多元线性回归分析
如前所述，气温、气压不是引起此次固原(泾源卧龙山) VP垂直摆异常的原因，但从短期来看，可能影响异常持续进程，为了剔除气温、气压变化对VP垂直摆短期变化的影响，本文选取无异常、数据较为稳定的2019年1月1日~2021年3月31日进行多元回归拟合，得到气温、气压对VP垂直摆观测的影响系数a、b和常数c，公式如下：
气温、气压与VP垂直摆的相关系数及多元线性回归拟合结果表3和图6所示，由表3可知，NS向、EW向与气温、气压的相关系数均较小，多元线性回归拟合的R2亦较小，表明VP垂直摆受气温、气压的影响并不显著 [11]。
Table 3. Relationship among VP vertical pendulum, temperature, air pressure and cave temperature
表3. VP垂直摆与温度、气压和洞温的关系
Figure 6. Results of multiple linear regression fitting
图6. 多元线性回归拟合结果图
5.3. 固体潮汐分析
固体潮汐能反应地壳结构变化，大量实验结果表明，潮汐因子等参数的时间序列变化能在一定程度上反映周边介质物性的非线性变化过程(郗欽文等，1994；唐久安等，1995；张晶等，2003)。正常情况下一个地方的潮汐因子会在一定幅度内上下波动(刘冬英等，1998；侯晓真等，2018)。本文使用Nakai拟合检验、潮汐变化分析、去潮汐和零漂分析对VP垂直摆的固体潮进行分析 [12]。
5.3.1. Nakai拟合检验
Nakai拟合检验是固体潮观测检验的常用方法，应用Nakai拟合模型，可以对固体潮的整点值观测数据进行拟合检验，求取固体潮观测振幅和理论振幅之比(振幅因子)、时间滞后等特征参量。如图7给出了固原(泾源卧龙山) VP垂直摆的Nakai拟合检验结果，从图中可以看出，此次异常发生期间NS向非潮汐拟合常数，EW向拟合方差、振幅因子、非潮汐拟合常数存在异常 [13]。
5.3.2. 潮汐变化分析
1966年，维尼狄可夫提出了基于计算机数字滤波器进行滤波的方法，使用该法对固体潮资料进行调和分析，即被称为Venedikov调和分析法。潮汐变化分析为日滑动的月调和分析：计算采用观测整时值序列，通常以48小时为滑动步长，30天为窗长，采用月长度调和分析，解算各主要潮波的潮汐因子、相位滞后等参数。图8为VP垂直摆的潮汐变化分析结果，从图中可以看出，此次异常期间，NS向、EW向M2波、半日波潮汐因子变化存在低值异常 [14]。
5.3.3. 去潮汐、零漂分析
去除潮汐、零漂后，NS向、EW向日差值、消除潮汐改正和零漂后的高频信息两项存在异常(图9)。
Figure 7. Nakai fitting test of VP vertical pendulum
图7. VP垂直摆Nakai拟合检验
Figure 8. Analysis of tidal variation of VP vertical pendulum
图8. VP垂直摆潮汐变化分析
Figure 9. Analysis of VP vertical swing to tide and zero drift
图9. VP垂直摆去潮汐、零漂分析
5.3.4. 固体潮汐分析总结
由5.1~5.3可知，固原台(泾源卧龙山) VP垂直摆NS向存在5项、EW向存在7项固体潮汐异常，表4为对VP垂直摆固体潮汐异常的归纳总结。
Table 4. Summary of Guyuan tidal anomalies for the VP vertical pendulum at Guyuan Station (Wolong Mountain, Jingyuan)
表4. 固原台(泾源卧龙山) VP垂直摆固原潮汐异常总结
5.4. 同台及周边台站资料分析
5.4.1. 同台资料分析
泾源卧龙山测点架设有水管仪、伸缩仪、VP垂直摆三套形变仪器，目前水管仪变化正常，伸缩仪EW向为自2018年起被列为年度异常。泾源伸缩仪EW向2018年6月以来持续上升，呈张性变化，2020年趋势有所转平，张性减弱，同时2021年年初以来出转折下降–上升变化，呈现压缩–拉张的短期波动变化。此次异常与2017年九寨沟地震前异常类似，2017九寨沟地震前3年泾源卧龙山洞体应变EW向由拉张转为压缩变化，至震前半年压缩速率减小，随后出现加速压缩、拉张等短期波动变化 [15]。2021年5月22日发生青海玛多7.4级地震，认为此次异常为青海玛多7.4级地震的前兆异常的可能性较大(图10)。
Figure 10. Variation curve of daily mean value of strain in Jingyuan cave
图10. 泾源洞体应变日均值变化曲线
5.4.2. 周边台站资料分析
海原小山钻孔应变第一、第二、第三分量2021年4月11日出现明显转折变化，4月15日后逐渐恢复。海原小山钻孔应变以往多次出现类似变化，其中1次为显著强降雨造成，3次与周边中强地震对应，分别为2013年漳县–岷县6.6级地震、2016年中宁Ml4.4级地震、2017年固原4.6级地震。此次异常经现场核实后排除了降雨影响，认为该转折变化为前兆异常，2021年5月22日发生青海玛多7.4级地震，海原小山台距离震中，681 km，分析认为该异常变化与玛多7.4级地震的关系较大(图11)。
Figure 11. Variation curve of the whole point value of borehole strain in Haiyuan Xiaoshan
图11. 海原小山钻孔应变整点值变化曲线
5.5. GNSS基线变化分析
利用中国地震局第一监测中心提供的GNSS测站时间序列数据，计算了泾源卧龙山测点附近GNSS基线变化，分析了单站和周边GNSS基线时间序列异常基线分布图见图12。鉴于GSPL (甘肃平凉)站在2020年起出现异常，与GSPL站有关的基线必然出现异常，故未选择GSPL (甘肃平凉)站进行基线时间序列分析。
图13为泾源卧龙山附近异常GNSS测站时间序列图，从图中可以看出GSPL (甘肃平凉) 2020年2月以来转折加速，由缓慢南向运动转为加速北向运动；由缓慢西向运动转为加速东向运动；整体由北东向运动转为南西向运动。NXHY (宁夏海原)站2020年2月以来转折加速，由北向运动转为南向运动；由缓慢西向运动转为加速东向运动，整体由北西向运动转为东南向运动。GSPL和NXHY站测站坐标时间序列呈准同步变化。
图14为泾源卧龙山附近异常GNSS基线和方位角时间序列图，从图中可以看出，GSDX (甘肃定西)-GSQS (甘肃清水)、GSJN (甘肃静宁)-GSQS (甘肃清水)、GSJN (甘肃静宁)-GSQS (甘肃清水)基线时间序列和方位角均在2020年年初出现转折下降异常，由拉张转为挤压，并发生方位角偏转；GSLX (甘肃陇西)-GSTS (甘肃清水)基线时间序列在2020年初转折上升，由压缩转为拉张。
GNSS测站和基线时间序列分析的结果表明，2020年初以来，泾源卧龙山附近地壳形变发生了显著变化，主要以由张性运动转为压性运动为主 [16]。
Figure 12. Distribution map of stations and GNSS baselines near the Wolong Mountain measuring point in Jingyuan
图12. 泾源卧龙山测点附近测站和GNSS基线分布图
Figure 13. Time series diagram of GSPL station (a) and NXHY station (b) near Wolong Mountain in Jingyuan
图13. 泾源卧龙山附近GSPL站 (a)和NXHY站 (b)时间序列图
Figure 14. Anomaly map of GNSS baseline at Wolongshan station in Jingyuan
图14. 泾源卧龙山测站GNSS基线异常图
5.6. 震例分析
VP垂直摆架设于2018年11月，此前，同一观测墩架设VS垂直摆，VP垂直摆是VS垂直摆的升级，观测原理相同，故VS垂直摆和VP垂直摆在异常变化上可进行对比分析。自2007年3月VS垂直摆架设以来，共出现3次异常与此次异常类似的变化，仅为EW项，分别为2012年9月12日~2013年9月2日、2015年4月17日~11月14日、2017年6月28日至8月27日，异常持续过程中，分别发生2013年4月20日芦山7.0级地震、2015年4月25日尼泊尔8.1级地震和2017年8月8日九寨沟7.0级地震。下面对三次地震前的异常形态进行具体阐述。
1) 2013年4月20日芦山7.0级地震
2013年4月20日芦山7.0级地震位于青藏块体东缘，距离泾源卧龙山测点655 km。2012年9月12日起，固原(泾源卧龙山) VS垂直摆EW向转折加速下降，至2012年12月26日恢复上升，2013年1月1日再次转折下降，2013年3月27日再次恢复上升，2013年4月7日转平，转平13天后，2013年4月20日发生芦山7.0级地震(距异常发生220天)，地震发生后，异常并未立刻恢复，87天后，2013年7月16日恢复上升，2013年9月2日恢复到异常变化前水平。此次异常经历了转折下降–恢复上升–转折下降–恢复上升–平稳变化–恢复上升的复杂过程，自异常发生至发生芦山7.0级地震，经历了220天，并于震后135天恢复到异常发生前，异常最大幅度为1857.38 × 10−3''，上升了15倍(图15(a))。
2) 2015年4月25日尼泊尔8.1级地震
2015年4月25日尼泊尔8.1级地震位于青藏块体东缘，距离泾源卧龙山2192 km。2015年4月17日起，固原(泾源卧龙山) VS垂直摆转折加速下降，8天后，4月25日发生尼泊尔8.1级地震，5天后，4月30日开始恢复上升，至11月14日，211天恢复至异常前状态。此次异常经历了转折下降—恢复上升的过程，异常最大幅度为306.369 × 10−3''，上升了2倍(图15(b))。
3) 2017年8月8日九寨沟7.0级地震
2017年8月8日九寨沟7.0级地震位于青藏块体东缘，距离鲸鱼卧龙山测点345 km。2017年6月28日起，固原(泾源卧龙山) VS垂直摆转折加速下降，40天后，2017年8月7日恢复上升，8月8日发生九寨沟7.0级地震(图15(c))。此次异常经历了下降–上升的过程，异常最大幅度为236.701 × 10−3''，自异常开始至恢复至异常前状态历时60天(图15(c))。
综合对比三次震例，均仅发生在EW向，NS向无异常，异常变化的显著特征均为转折加速下降，对应的地震均在站点的西南方向，青藏高原东缘；但异常持续时间、异常进程、异常幅度、差异较大。此次异常为NS向、EW向同步转折加速异常，EW向异常为转折加速下降，与以往异常类似，但在异常进程、异常幅度与以往异常不同，2021年5月22日发生青海玛多7.4级地震，在泾源卧龙山测点西南方向，距离测点747 km，位于青藏高原东缘，认为此次异常为青海玛多7.4级地震的前兆异常的可能性较大，但震后异常持续，目前没有恢复迹象，需继续关注异常发展变化 [17]。
Figure 15. Anomaly change of VS-type vertical pendulum in Guyuan (Wolong Mountain, Jingyuan)
图15. 固原(泾源卧龙山) VS型垂直摆异常变化图
6. 异常变化分析与判断
本文通过对仪器状态、环境干扰、数据分析、地壳变形等资料的多方面分析，结果如下：
① 泾源VP垂直摆2021年5月以来多次出现无法远程调零现象，经仪器厂家检测重新安装后恢复功能，并判定与此次异常无关。格值无明显变化，无人为干扰，数据观测质量无显著变化。
② 对周边环境进行调查，可能对此次VP摆造成影响的干扰因素为气温、气压和卧龙山水库，通过对气温、气压与VP垂直摆进行多元线性回归分析、卧龙山水库对VP垂直摆的影响进行数值模拟分析，排除了环境干扰。
③ 利用Nakai拟合检验、潮汐变化分析、去潮汐和零漂分析对VP垂直摆进行固体潮汐分析，发现存在多项固体潮汐异常，表明此次VP摆异常不仅为趋势转折变化，也同时出现固体潮汐变化。
④ 同台泾源伸缩仪器EW项自2018年6月出现张性变化异常，附近测点海原小山钻孔应变第一分量、第二分量、第三分量自2021年4月11日出现明显转折变化。
⑤ GNSS测站时间序列和基线时间序列的结果表明，自2020年年初开始，泾源卧龙山测点周边地壳形变发生显著变化。
⑥ 震例分析的结果表明：此次VP垂直摆异常与之前三次VS垂直摆异常在异常变化形态上存在相似之处，但在异常进程、异常幅度等方面有较大不同。
综合分析认为，固原台(泾源卧龙山)VP垂直摆为前兆异常，但是否是青海玛多7.0级地震的地震前兆异常，尚需看后续异常是否恢复 [18]。
7. 结论
综合现场检查与分析结果，目前已基本排除仪器故障、观测环境干扰、气象因素与人为干扰等干扰情况；VP垂直摆NS向、EW向在趋势变化和固体潮变化两个方面均存在异常。
2018年6月以来，泾源伸缩仪EW项存在张性异常，反映测点附近应变发生变化；2020年以来出现多项GNSS异常，反映测点附近地壳形变发生变化；2021年4月11日海原小山钻孔应变出现趋势转折异常，反映测点附近应变出现短临变化；此次异常与以往三次地震前兆异常有雷同之处 [19]。因此，认为此次出现的短期变化具有较高的可信度，判定为A类异常；2021年5月22日发生青海玛多7.4级地震，此次异常为青海玛多地震的前兆异常可能性较大，但目前异常扔持续，无恢复迹象，且2021年5月30日内蒙古鄂托克前旗发生ML3.6级地震，可能预示着未来3个月及稍长时间，宁夏南部及邻区存在发生中强地震的可能，还应继续密切关注异常发展，注意宁夏南部及邻区发生5.0级以上地震及青藏块体东缘发生7.0级以上地震的可能 [20]。
基金项目
基金资助：2022年度宁夏地震科研基金课题(NX202205)。