转换函数方法在磁台磁场异常变化分析中的应用

转换函数方法在磁台磁场异常变化分析中的应用

“转换函数”方法地磁场是地球上最重要、最基本的物理场所之一。一百多年以前,人们就已注意到地磁异常变化与地震的关系,国外做了许多这方面的研究工作。1966年以来,我国也逐步建成了比较完善的地磁观测系统,在震磁效应理论研究和地震预测方面取得了很大进展(丁鉴海等,1994)。经过多年的实践与研究,地磁已成为地震预测不可缺少的手段之一。“转换函数”方法是国内用于预测地震的地磁方法之一,林云芳等(1999),曾小苹等(1995)和龚绍京等(1991,1997)基于地磁台站磁场连续记录仪器磁变仪产出的磁照图经数据采集后的资料,均使用“转换函数”方法在地震前提取到了地磁场的异常变化。近几年来,随着越来越多地磁台站的数字化,可直接利用地磁台站记录的数字化数据来计算“转换函数”。李琪等(2005)对模拟数据和数字化数据从时间序列和转换函数计算结果两方面进行了对比,发现二者基本一致。2006年7月4日11时56分河北省文安县发生5.1级地震,震中距离天津静海地磁台仅52km,成为研究地震前后转换函数变化过程的极好震例。本文利用静海地震台的地磁数据,分析了文安地震前后静海地磁台转换函数参数的变化,并与距离文安地震震中298km的泰安地磁台的计算结果进行对比,试图找到距离震中距离远近不同的台站在转换函数参数的变化规律上有何相同与不同之处。1地下电性环境分析转换函数转换函数的基本原理在许多文献中都有详细说明,在此参照曾小苹等(1998)撰写的《地磁方法预报灾害》中的转换函数部分内容进行简单的描述。由于电磁感应效应,磁感应强度F的分量ΔZ、ΔH和ΔD满足以下关系ΔΖ=AΔΗ+BΔD(1)ΔZ=AΔH+BΔD(1)式中的A和B为转换函数。转换函数A和B的时空变化特征为地下电性构造的函数,因而可用于研究局部和部分地区地壳和上地幔电性结构及其变化。据(1)式对地面各测点所记录到的天然磁场信号的分析得到的转换函数A和B,是位置(λ,φ)、深度(h)、频率(ω)或周期(T)和时间(t)的函数,即A=A(λ,φ,h,ω,t)B=B(λ,φ,h,ω,t)(2)A=A(λ,φ,h,ω,t)B=B(λ,φ,h,ω,t)(2)对某一测点(λ0,φ0),某一频率(ωj)的地磁场F,由电磁波的趋肤效应即可得到该测点对应ω0的穿透深度hj=(2/ωjμσ)12(3)hj=(2/ωjμσ)12(3)处的电性信息。其中μ和σ分别为介质的磁导率和电导率。由(2)式、(3)式可知,分析研究转换函数随时间t的变化,则可了解和监视地震活动区地下电性环境的变化。由(1)式已知,某测点处、某一频率(或周期)j的地磁短周期变化量满足关系ΔΖj=AjΔΗj+BjΔDj(4)ΔZj=AjΔHj+BjΔDj(4)其中复转换函数(暂去掉角标j)A=Ar+iAiB=Br+iBi(5)A=Ar+iAiB=Br+iBi(5)角标r和i分别表示实部和虚部。A、B之模|A|=(A2r+A2i)12|B|=(B2r+B2i)12(6)|A|=(A2r+A2i)12|B|=(B2r+B2i)12(6)A、B的相角为φA=tg-1(Ai/Ar)φB=tg-1(Bi/Br)(7)φA=tg−1(Ai/Ar)φB=tg−1(Bi/Br)(7)帕金森(Parkinson)矢量的磁方位角αp、倾角Ip和长度Cp为αp=tg-1(Br/Ar)Ιp=tg-1(A2r+B2r)12Cp=sinΙ(8)αp=tg−1(Br/Ar)Ip=tg−1(A2r+B2r)12Cp=sinI(8)解矩阵方程(4)的总方差σz=[n∑j=1(¯ΔΖ-ΔΖj)2n-1]12(9)其中¯ΔΖ=n∑j=1ΔΖj/n2各周期对应的谱值数据采集、分析处理流程按如下顺序进行:数字化数据→高通滤波→谱分析→按周期重组事件→计算转换函数参量。(1)对三要素分钟值进行Fourier谱分析,得出各周期对应的谱值。(2)选出从6～164min之间以2min为间隔递增的80个周期的转换函数参数。(3)对给定的各周期对应的谱值ΔZj、ΔHj和ΔDj生成数据矩阵,用最小二乘法,按(4)～(9)式求解,得出转换函数A,B及其参量:转换函数的模|A|和|B|,转换函数实部参量Ar、Br和虚部参量Ai、Bi,实部参量和虚部参量的差Ac、Bc,实部参量和虚部参量的相角Ao、Bo,方位角αp以及总方差σz等。36年6月10日第1组总方差z异常点的出现距离文安地震震中最近的天津静海地磁台是国家基准地磁台站,观测环境和数据质量都很好。本文中用到的两个地磁台的数据都是GM3磁通门磁力仪记录的分钟值采样数字化数据。图1给出了静海地磁台2006年1月1日～2006年9月30日总方差σz日均值时间序列曲线,图中虚线的值为均值加上2倍标准偏差。如果把超出2倍标准偏差的数据点视为异常点,可看出异常点首先在短周期出现,并逐步过渡到长周期。从时间上看,周期为6min的总方差σz于2006年2月中旬,文安地震前4个半月出现了第一组异常点。然后在3月中旬和4月上旬仍是周期为6min的总方差σz出现了两组刚好超出2倍标准偏差的异常点。在6月上旬,距离发震时间不到一个月时,周期为6min、36min、115min和156min的总方差σz同时出现了异常点。2006年6月底到地震发生这段时间内,周期为36min、66min、115min和156min的总方差σz都出现了异常点,而周期为6min的总方差σz无异常出现。文安地震之后到7月底这段时间,周期为66min、115min和156min的总方差σz都有异常点出现,而周期为6min和36min的总方差σz都无异常出现。如果仅关注各周期异常最大的点,则可以看到异常极大点出现的时间随着周期的增长而向后推移。震后短周期的异常点消失,而长周期的异常点继续存在了一段时间,到7月底所有周期的异常点都消失,恢复平静状态。图2给出了静海地磁台2006年1月1日～2006年9月30日转换函数A的模|A|的日均值时间序列曲线,图中虚线的值为均值加上3倍标准偏差。从图中可看出,2月下旬首先在短周期出现(T=6min)了第一组异常点,其后依次在周期为116min和66min出现异常点,而在文安地震发震前后,除了周期为6min以外的数据曲线,都出现了异常点。图3给出了静海地磁台2006年1月1日～2006年9月30日转换函数B的模|B|的日均值时间序列曲线,图中虚线的值为均值加上3倍标准偏差。除了2006年初出现在周期为116min的异常点之外,从文安地震发震前后异常点的出现时间来看,仍是异常点出现的时间随着周期的增长而向后推移。用同样的方法分析了山东泰安地磁台的数据,计算结果见图4和图5。2006年距离泰安地磁台300km范围内发生了两个地震,分别是2006年4月9日的ML4.7地震(震中位置为115.4°E,35.7°N,震中距168km)和文安地震。从图4和图5中可看到,这两个地震前转换函数参数都出现了异常点,而ML4.7地震前出现的异常(周期为36min和66min)在时间上是同步的,而在文安地震前出现的异常尽管时间上并不完全同步,但也相差不远。此外,在文安地震震后也有比较大的异常出现,尤其在图4中表现尤为明显。从泰安地磁台转换函数参数的变化来看,尽管不像静海地磁台那样震前有异常点出现时间随周期的增长向后推移的规律,但在震前震后都有明显的异常出现。4异常点的出现时间本文分析了河北文安地震前后静海和泰安地磁台的转换函数参数的日均值时间序列曲线,发现①静海地磁台总方差σz数据异常点首先在短周期出现,随着时间向后推移逐渐在比较长的周期出现,尤其是异常的极大点表现得更为明显。在静海地磁台转换函数的模|B|中也有类似的规律出现;②尽管静海地磁台转换函数的模|A|数据异常点出现规律不明显,但异常点仍是最先出现在短周期(T=6min),而最后出现于长周期(T=156min);③泰安地磁台的转换函数参数(总方差σz和转换函数的模|A|)异常点出现时间随周期的增长向后推移的规律不是特别明显,但震前震后都有很明显的异常出现;④震后一个半月内,静海和泰安地磁台的所有转换函数参数都仍然有异常点出现。曾小苹等(1998)发现,如果震前所有周期的总方差σz的数值都增大,预报的可信度会增加。龚绍京等(1998)也提出了判断异常的“成片”原则,即几个相邻的周期同时出现异常。在本文分析的震例中也同样观察到了类似异常“成片”出现现象,这与过去的震例基本一致。但是,从细节上看,距离震中最近的静海台的转换函数参数每个周期的异常极大点出现的时间不一致。总的来说,出现时间从短周期向长周期迁移,地震发生之后的异常更多出现在长周期。在过去的转换函数计算震例中并没有发现类似的现象。从转换函数的原理来分析,不同周期的转换函数参数对应不同深度的地下电性介质的变化,周期越长,对应的地下介质越深。数据异常点的出现时间从短周期逐渐向长周期过渡,有可能是地下介质由浅至深发生了变化。也许由于静海地磁台距离震中比泰