(1)清水温泉井

清水温泉位于清水县城北东8 km的汤浴河谷内(34.45°N,106.13°E),海拔高程1 445 m。地势北东高,南西低,属低山丘陵。其构造位置属于陇西系南段、秦岭北缘大断裂的北侧,六盘山断裂带西南侧(图1)。台基岩性为通渭—清水断裂晚古生代花岗岩,该温泉是深部热水沿岩石裂隙运动,在一定地质构造条件下出露于地表而形成的。清水温泉井深443.3 m,水温观测仪器为SZW-1A型数字式温度计,水温探头位于井下185 m处。2007年6月仪器入网,产出正常数据产品。

(2)临夏井

临夏井位于临夏市东北方向约7 km的折桥乡(35.63°N,104.26°E),海拔1 830 m。地形开阔平坦,周围没有明显干扰。第四纪覆盖层厚10～20 m,为亚粘土卵石,下伏第三系基底花岗岩。井深200.1 m,2006年11月开始水温数字化观测,2007年入网,2007年6月产出正常数据。观测仪器为SZW-1A型数字化水温仪,水温探头位置在井下185 m处的第三系花岗岩中。

图1 研究区域构造背景及井点分布图
Fig.1 Regional tectonic background and distribution of well points in study area

赵刚等(2009)对全国近277个观测井的地热正常动态进行了系统研究,将其划分为降温型、升温型、波动型、稳定型、跳变型、长周期型6种基本类型。据统计,在整个地下流体台网中,降温型、升温型和波动型占70%左右。在了解其正常动态特征基础上,才能够提取地震前兆异常信息。井水温度的背景值,一般取决于井孔所在地的大地热流(q)或地温梯度(dT/dz)。一般来说,热流密度大的地区地温高,地温高的地方水温也高。根据地温成因,一个地区不同深度的地温状态是不同的。一般由地表到一定深度段上,地温受太阳辐射热的影响而随时间有规律起伏,如随昼夜、季节变化等,该变化带称为变温带。变温带的厚度一般为十几米至几十米。变温带以下有一个面,该面上的深度不受太阳辐射的影响,一般常年不变,称为恒温带。恒温带以下水温,主要受地温梯度的控制,即随着深度(Z)、水温(T)升高。

清水温泉井水温自架台产出数据以来基本呈现趋势上升变化(图2a),年变幅度达到0.05 ℃,长期动态特征属于漂移上升型,探头位置处于恒温带下,所测为深部水温变化,不受地表温度影响。日变幅在0.02 ℃左右,月变幅约为0.03 ℃。

临夏井水温自2007年6月产出数据以来基本呈现上升趋势变化(图2d),年变幅度达到0.03 ℃,长期动态特征属于漂移上升型,探头位置很深,处于恒温带下,所测温度为深部水温变化,不受地表温度影响。日变幅在0.01 ℃,水温观测无潮汐效应,月变幅约为0.02 ℃,该井日变幅较大,几乎为年变幅的30%,和清水温泉井水温相似。

(1)直线回归分析

变量之间的关系可分为两种:一种是函数关系,另一种叫相关关系,又叫统计相关。水温观测所蕴含的信息成分比较复杂,非单一因素引起,是多种成分的混合、叠加,因而无法找到一个确定的函数关系,但数据的变化还是有一定的规律可循。假定水温是一个随时间变化的量,并与时间变化有相关关系,且按照:

y=f(x).(1)

的函数形式相关。通过实测的观测数就可以计算出几个数值,用直线来大致表示两个变量的关系:

y=a+bx,(2)

其中y为水温,x为时间,直线斜率b为回归系数,a为常数,a,b的值计算公式为

b=(∑^N_t=1x_ty_t-1/N∑^N_t=1x_t∑^N_t=1y_t)/(∑^N_t=1x²_t-1/N∑^N_t=1x²_t),(3)

a=y^--bx^-.(4)

在直线回归分析公式中,自变量x为时间序列,因变量y为水温测值。

(2)斜率K值法

斜率K值法原理上就是趋势速率法,在本文中求取的K值就是水温日变化速率。如果水温数据以一个恒定的速率上升,则计算得出的斜率K值为一个相对固定的值,那么K值曲线则是一条直线; 若高于或低于这个值,K值曲线偏离直线的时段都属于异常。

清水温泉井和临夏井长期动态类型都属于长期漂移上升型,但对于一个一直随时间呈趋势上升变化的量,单从原始曲线上无法较直观的观察到这种变化,因此计算两个井水温动态变化的斜率,可以提取地震异常。由于清水温泉、临夏井日变幅较大,对两台日均值数据先进行11日滑动处理,再以11日滑动均值计算斜率。通过反复计算,笔者采用步长为30日的斜率K值,以抑制短周期跳变信息,突出中长周期成分信息:

K_i=(T_i+m-T_i)/m, i=1,2,3…….(5)

式中,K为斜率,T为水温,i为时间序列,m为步长,取值为30。

(3)从属函数法

前兆观测资料的异常变化形态虽然是多种多样的,如突跳、转折、加速等,但其本质都是观测曲线随时间变化的斜率变化的量。从属函数就是表征一条曲线变化的斜率和跳动程度(叶秀薇,2004)。通过对资料的分析比较,用11点滑动平均值作为数据序列计算的效果较好,明显优于直接使用月均值提取异常的方法(杨兴悦等,2006)。其从属函数为

μ_i=[1+α/(|k_i||r_i|)]^-1.(6)

式中,k_i为观测值变化斜率; r_i 为滑动平均值M(t)与时间t的相关系数,它反应了观测序列内在质量的好坏; α为经验常数。

通过多次程序计算,清水温泉井α取值为0.000 3,拟合窗长m值取9,计算出的μ值较佳; 临夏台α取值为0.000 4,拟合窗长m值取9,计算出的μ值较佳。

图2 清水、临夏两井水温数据经过处理后的曲线图
(a)清水井11日平滑日均值;(b)清水井斜率K值;(c)清水井从属函数;(d)临夏井11日
平滑日均值;(e)临夏井斜率K值;(f)临夏井从属函数
Fig.2 The curve of the processed water temperature data recorded by Qingshui and Linxia Well(a)smooth day mean for 11 days in Qingshui Well;(b)slope K-value in Qingshui Well;(c)subordinate function in Qingshui Well; (d)smooth day mean for 11 days in Linxia Well;(e)slope K-value in Linxia Well;(f)subordinate function in Linxia Well

2007年6月2日至2008年12月31日,清水温泉井水温一直处于快速上升变化状态,在汶川8.0级地震前后,没有明显的突变、转折等异常形态出现,但水温上升变化很显著,呈前期加速、后期减缓的异常变化(图2a),且与对数函数曲线形态很相似,这种在汶川地震前后水温曲线呈对数变化形态的水温井还有在甘肃东南部的成县井、武都樊坝井、ZK801井等,由于篇幅有限,此文不再详述。2009年9月1日至2011年3月11日水温变化相对以前速度放慢,日本9.0级地震发生后,水温测值明显回落,持续了约70天后,于2011年5月20日后恢复正常,数据以比较缓慢而相对固定速度上升,曲线形态呈现直线斜形上升变化。

按年对两井整点值数据进行线性回归分析,计算出回归方程(表1),2007年直线斜率b值最大,为0.000 2,即水温数据基本以每小时 0.000 2℃的速率在上升; 2008年降为0.000 1℃,在b值变化大的2007～2008年间,发生了汶川8.0级地震,2009年、2010年b值相同且较小,2011年b值变为很小的负值,且在2011年3月11日的发生了本9.0级大地震。2012年b值为正,但b值很小,仅为0.000 01。

表1 井水温直线回归方程a值与b值
Tab.1 The a-value and b-value in linear regression equation of well water temperature

2007年6月2日至2008年12月31日,K值以稳定的速率由高值逐渐下降(图2b),但相对其后的K值,此时段K值仍处于高值区,说明在此时段水温总体上升较快,但上升的速率是逐渐减小的; 2009年9月1日至2011年3月11日K值下降速度放缓; 3月11日日本9.0级地震发生后,K值出现一个明显的大幅度下降—回升的变化形态; 2011年5月20日后K值没有出现明显的上升、下降等变化形态,基本稳定在0.000 01附近。

2007年6月7日至2009年1月23日,从属函数μ值呈现明显的高值变化(图2c),汶川地震前出现的高值变化更多、更剧烈,μ最高值出现在2007年6月初,为0.62,震后有所衰减,但也很明显,异常变化一直持续到2009年1月24日,μ值接近为0。2009年11月26日至2010年4月5日,μ值又出现高值,最大幅度出现在2009年12月30日,μ为0.58,随后于2010年4月14日发生了青海玉树7.1级地震,μ值异常变化显著。在2011年3月11日日本9.0级地震前6天,μ值在2011年3月4日、3月5日出现明显大于背景值的高值变化,约为0.35; 2012年8月12日新疆和田6.2级地震前8天μ值也出现了0.34的较高数值变化。

临夏水温自2007年6月15日至2007年12月20日基本为缓慢下降变化(图2d),21日后开始由降转升,2008年1月9日开始加速上升,直到6月17日,水温数据累积变化达0.083℃,约为正常年变幅的6倍,在此时段发生了汶川8.0级大地震,异常变化形态非常显著,出现了明显的非对称“V”字型变化(何案华等,2012)。2008年6月18日后水温以指数形式回落,7月7日后开始以缓慢而稳定的速度呈上升变化,其后的2008年12月29日、2009年3月15日有两次比较明显的“V”字型下降变化,这可能是汶川地震的震后效应,但持续时间不长,变化幅度不大。临夏水温自2009年3月15日出现最后一个“V”字型变化后,水温以稳定的速率呈直线斜形上升变化。

2007年临夏井b值为负值(表1),2008年转为正值,且b值较大,为0.002,在此期间发生汶川8.0级地震。自2009年至今,b值变化不大。

临夏水温自2007年6月15日至2008年6月17日,K值基本为加速上升变化,2008年6月18日后急速回落,下降到一个低值后恢复到正常变化状态,在此时段发生了汶川8.0级大地震,K值异常变化形态非常显著(图2e)。2008年7月7日后,K值基本稳定在0.000 2附近,没有趋势上升和下降变化形态,这种变化持续至今。

临夏井μ值自2007年6月起基本上呈现上升型高值异常变化,汶川8.0级地震后不久,于2008年6月21日达到最高值0.788,然后从高值快速回落到背景值0(图2f),此后至2009年3月22日期间还出现了3次μ值高值变化。青海玉树地震前没有出现μ值高值变化,2011年3日11日日本9.0级地震前11天,μ值出现了持续4天的高达0.32的变化。

两水温观测井在汶川地震前都有明显的前兆异常,水温测值震前都出现快速上升变化; 但异常持续时间不同,临夏井震后1个多月基本恢复正常,而清水温泉井则持续到2008年底。

正常情况下,一个上升漂移型的井水温会以稳定的速率上升,因而斜率K值几乎保持不变,曲线形态呈现为直线形态; 但当有干扰因素或有地震前兆异常出现时,K值就会有大的变化。两井在汶川地震前后斜率K值都处于高值状态,但清水温泉是由高值快速降低,而临夏井则是加速上升变化,震后K值变化的时间长度也不同,清水温泉井持续到2008年底,临夏井震后一个多月恢复正常。

从属函数是对测值斜率变化的另一种反映方式,虽然和K值变化有一定的相似性,但只要参数取值合理,异常变化则会更显著。两水温观测井μ值在汶川地震前出现了多次大于0.5的异常变化,异常变化都非常明显,震后都出现了μ值的高值变化; 但清水温泉井震后μ值逐渐减小,而临夏井在震后μ值则呈现增大趋势,一个多月后达到最高值然后快速减小。

(1)对于长期漂移上升(下降)型且日动态变化很大的水温观测井,采用一定的数学方法进行异常的提取较合理,斜率K值法、从属函数法是两种非常好的异常提取方法。

(2)汶川8.0级地震前清水温泉井和临夏井都出现比较明显的异常变化。但两井的异常形态具有差异性,清水温泉井地震异常形态呈现为“对数变化形态”,而临夏井则呈现为非对称“V”字型变化形态。震后两口井都有继承性异常变化,但结束时间不同步。

(3)清水温泉井在青海玉树7.1级地震前出现了测值跳动、K值波动,μ值高值变化等异常变化形态,临夏井却没有出现异常。日本9.0级地震前两口井都出现幅度较小的异常变化,但震后形态不同,清水温泉井震后测值明显下降,K值降低,而临夏井在震后立即恢复正常。

(1)清水温泉井和临夏井在汶川地震前都出现了明显的地震前兆异常。虽然两井都距汶川约500 km,但都在青藏块体的东北缘,因而印度洋板块对欧亚板块的挤压所引起的青藏地块应力场的增加,也同样在这两个地方能反映出来。在强震活动的高潮阶段,构造块体内各地震带(包括未参与强震活动的地震带)都有应力增强的过程(耿鲁明,石耀霖,1993)。如果把青藏地块应力场的增加看做是“因”,则汶川地震的发生,和清水、临夏应力场的增加则是3个不同形态的“果”。两口井在汶川地震后都继承性的延续了一定时段震前异常形态,这可能是汶川地震后累积的应力并未完全释放完毕,后续地应力释放引起了区域应力场的变化。

(2)清水温泉井和临夏井相距约290 km,而青海玉树7.1地震前只有清水温泉井出现了地震异常。这可能与两井所处构造部位不同有关,清水温泉井位于NNW、NWW向断层的交汇部位附近,处在鄂尔多斯地块、华南地块、青藏地块的交接部位。因为其构造部位的特殊性,按照“多点场应力集中模型”(马宗晋,1980),清水温泉井所处位置可能是应力集中点,对区域应力场的变化反应更灵敏。

(3)日本9.0级地震前后,清水温泉井原始曲线、从属函数异常都比较明显,震后不久水温数据恢复正常,变化速率也即斜率K值基本保持在一个基值附近。临夏井只有从属函数出现异常,但不显著。日本地震距这两口井3 000 km,这样的远兆异常在地下流体观测中不断被发现(车用太等,1999)。