对《转换函数与汶川大地震关系的初步研究》一文的分析*

龚绍京 马 骥 刘双庆 栗连弟

（天津市地震局，天津300201）

对《转换函数与汶川大地震关系的初步研究》一文的分析*

龚绍京 马 骥 刘双庆 栗连弟

（天津市地震局，天津300201）

我们不能苟同该文作者的一些概念，认为该文用的计算公式以及所得的结果存在较多问题，其所得转换函数的量级以及与地震对应关系的分析叫人难以理解。

复转换函数；帕金森矢量；电性结构；复数最小二乘法

1 关于计算方法

关于复转换函数的计算公式，我们发现一些作者的描述概念不清。例如该文［1］及引文［2-4］的作者仅求出了功率谱，利用功率谱求出的振幅谱和（1）式求出了复转换函数或其模|A|、|B|或帕金森矢量系数a，b。仅用功率谱是无法算出复转换函数的，因此，他们求出的不是真正的复转换函数及其模，也不是真正的帕金森矢量系数。

1.1 复转换函数的计算方法

求复转换函数可有两种计算方法：（a）对某一事件三个分量的时间序列，求它们的自谱和互谱，然后利用下面公式求复转换函数：

式中，Phh、Pdd是自谱，即功率谱。Pzh、Pzd、Phd等是互谱。A＝Au＋i Av，B＝Bu＋iBv。下标u代表实部，v代表虚部。功率谱互谱方法由于只用一个事件即可求转换函数，所以要求样本长度较长。一般为24小时，也有用12小时的。即使如此，由于连续的长时间扰动多发生在磁暴期间，源场准均匀的假设很难满足，所以计算结果涨落较大。有人用5个事件的转换函数求平均。我们从对青光台短周期事件的时间序列分析［6］中发现，源场的扰动有一定的“惯性”，因此，认为5个事件太少，结果很难理想。（b）对一系列事件三个分量的时间序列进行离散傅里叶变换，求出各周期成分的余谱、求积谱。然后，用复数最小二乘法推导的Everett-Hyndman公式求A，B。实数的最小二乘法是要使∑e2达到最小，复数最小二乘法是要使∑e¯e达到最小，e为残差，∑为求和符号。公式如下：

1.2 用振幅谱和公式（1）可否作为求Parkinson矢量公式或（3）式的近似？

丁鉴海在《地震地磁学》［2］一书中提到：“帕金森矢量及转换函数的求解公式如下：对给定的各周期对应的振幅谱ΔZ、ΔH和ΔD生成数据矩阵，用最小二乘法解矩阵方程（5-45）（注：即（1）式），求转换函数A、B及相应的其他变量：A＝Au＋i Av，……”。我们认为这个提法不妥，我们不知道光有振幅谱没有相位谱怎么能求出复转换函数A＝Au＋i Av，……？因为由振幅谱ΔZ、ΔH和ΔD生成的数据矩阵中根本就没有虚数。

还要说明帕金森矢量和转换函数概念的来源完全不同，前者来源于电磁感应的基本原理，并由事实证明地磁变化矢量有限定在一个平面上的趋势（而非该书［2］中所提的扰动矢量，因扰动矢量另有含义）。后者来源于数学和工程控制领域，表达在复频率域内输出和输入的函数关系。在作法上，帕金森和威斯都要求量取相同时间间隔内三个分量的变化幅度，尤其帕金森矢量不一定量取三个分量到极值点的变化幅度——振幅。因三个分量的起止点都相同，才能用三个分量来描述地磁变化矢量。因此，从概念上说，振幅谱的ΔZ、ΔH、ΔD与地磁变化矢量的三个分量ΔZ、ΔH、ΔD完全是两码事。只有当三个分量间的相位差可忽略时，才可量取三个分量的振幅（而非振幅谱）并近似求出帕金森矢量系数，也就是所谓的实转换函数a、b。因此，从概念上说，用振幅谱和（1）式求出的不是帕金森矢量系数。

现在我们来分析可否用振幅谱和（1）式求出复转换函数或其模？

按二元回归分析原理［5］，（1）式的系数a、b应由下式求得：

按前述作者的定义，ΔZ、ΔH、ΔD是振幅谱，都是实数。我们可以证明，当（3）式的Z、H、D都是实数时，即三个分量间不存在相位差时，（3）式即退化成（4）式（证明从简）。但实际情形中，Z、H、D三个分量的变化从来都不是同相位的。

我们来看看L11、L22、L12、L10、L20与N、W、X、P、O的对应关系。设三个分量某周期成分的余谱和求积谱为：Zuj，Zvj，Huj，Hvj，Duj，Dvj，则按他们的定义，ΔZj。此时：

同理，P与L10，O与L20亦是不可比的。

最好的情况是Z和H的相位差可忽略，但D与Z、H存在明显的相位差，例如昌黎的情况。此时可设经谱分析后某周期的谱成份为：

即Zj和Hj都只有实部而无虚部。j为事件的序数。式中k也表达事件的序数，第一项是j＝k时的和。第二项是j≠k时的和，只有当第二项中的虚部可忽略时，（3）式的分母才能相当于（4）式的分母。

但这是不可能的。至于分子部分，更是不能比。为此我们找到几组经FFT运算后求出的Zuj、Zvj、Huj、Hvj、Duj、Dvj值。用由复数最小二乘法推导的（3）式求出Au、Av、Bu、Bv、|A|、|B|，又用振幅谱按（4）式求出a、b。用实际资料进行对比，以验证可否用振幅谱和（1）式求复转换函数及其模。结果列于表1。差距很大，事实说明上述引文中的方法是有问题的。

表1 按（3）式求出的结果与按（4）式求出的结果的对比

表1的文件名中，第一、二个字母代表台站，He代表菏泽，L代表崙坪，Q代表泉州。因此，LQ文件对应的数值不是A、B，而是水平场转换函数，且最后一行是做错动试验的结果。

从表1看出，不仅数值差别大，符号甚至相反。由振幅谱和（1）式求出的a、b也有正负，并不是转换函数的模，也不是复数。另要说明，人为的忽略正负而取a、b的模（即文［1，7］中的|A|、|B|），这样的做法是不正确的。模|A|、|B|没有明确的物理含义。一些文章画出了它们的逐月变化［7］，由于忽略了正负，图中的曲线会显得平坦一些。但转换函数是要构成帕金森矢量的，并要由此矢量指出导电率高的一方并由此分析震源区电性的变化。如果由‘模'来组成帕金森矢量，那么矢量只能永远指向西南方了［8］。可我们引用的3个图中，矢量可以指向四面八方。

2 关于转换函数的量级问题

我们知道帕金森（威斯）矢量表达地下电性的横向不均匀性，并指向（背向）导电率高的一方。它的典型表现有所谓的海岸效应（图1和图2），内陆异常（图3）和电流通道。这种横向不均匀性愈大，矢量的值也愈大。也就意味着构成矢量的a、b（或Au、Bu）值中的一个或两个会比较大。既然它们是表达地下电性结构的参量，那么对一个固定的地点，在正常情况下应该有大体一定的取值。

表2列出了我们收集到的各台站的取值。从表2和图1～图3可以看出：在内陆，只有在高导层的隆起区或断裂带附近矢量才比较大，如昌黎；在海岸，靠近深海的地方其矢量大于浅海附近。在日本东部的深海沟，在深1 000多米的J1和2 000多米的S1的斜坡上，C＝（A2u＋B2u）1/2达到1.7和1.9。这也是我们查到的最大C值。在陆地上，柿岗的a值达0.6～0.7，这与该地处于日本海岸线的拐弯处，朝南面向太平洋有关。

再看作者作出的结果。其纵坐标是取转换函数值的对数，|A|的最大坐标值分别是1.14，0.9，0.82，0.6，5.24，6.6。那么转换函数的最大值应是纵坐标值的反对数。作者没有特别说明是自然对数，我们只能按常用对数算。算出来转换函数达到13～3×106，这样的结果令人难以理解（表3）。

图1 帕金森矢量的海洋效应［9，15］（矢量长度为：L＝sin I，I为优势面的倾角）

图2 日本海沟处的帕金森矢量［10］（矢量方向按帕金森的定义，但长度为：C＝（A2u＋B2u）1/2）

3 与地震的对应

表3列出了各台至汶川地震的大致距离及异常的最大量级。按作者的说法都有反应。这里有几个问题：（1）异常的可靠性？因转换函数的值不可靠，因而异常也难以叫人置信。（2）对应关系非常牵强。从袁宝珠等［1］的图上看，地震前和地震后都有异常，如果作者再多处理几年资料，会不会也是这样跳上跳下，不定时地出现异常？（3）作者认为所有列出的台都有异常，这么远的距离，这么大的量级，这种异常真叫人难以想象。

图3 波兰的威斯矢量分布和海西期基底等深线图［11］

表2 各台短周期变化参量的量级（包括存在异常的时段）

我们认为寻找地震前兆是件很严肃的事。怎样才算是异常？怎样才算是地震前兆？应该有严格而统一的标准和有效的检验方法。例如，该文中异常非常多，震前震后都有。我们要问：“在没有地震的正常年份，是否也常常有这种跳上跳下的‘异常'？作者是否认真地思考过，并尽可能多处理一些年份的资料以验证这一点？”。同时，转换函数的大小也是一件很重要且包含重要物理意义的事。多年前王锜［12］的文章中也是出现了转换函数在当地不可能出现的量级，Au在0.2～0.8之间，Bu在0.6～1.2之间。如果王锜的结果正确，则该地应该是一个比昌黎——渤海地区更大的导电率异常区，这样的异常区在内陆还没有发现过。陈伯舫［13］和龚绍京等［14］已经对其进行了分析。

表3 各台的震中距及该文中异常的最大量级

由于转换函数有明确的物理含义，许多人都想做这方面的工作。我们认为应该有组织地交流与学习这方面的知识，并对各种不同的作法进行分析比较，从而选择一种较好的作法。

（作者电子信箱，龚绍京：caogong2003＠hotmail.com）

［1］袁宝珠，陈化然，张素琴，等.地磁转换函数与汶川大地震关系的初步研究.国际地震动态，2009（7）：69-75

［2］丁鉴海，卢振业，黄雪香.地震地磁学.北京：地震出版社，1994：270-273

［3］曾小苹，林云芳.地磁短周期变化异常对中国中强地震的响应.地震，1995，1：29-36

［4］林云芳，曾小苹，续春荣，等.地磁方法在地震预报中的应用.地震地磁观测与研究，1999，20（6）：35-44

［5］中国科学院数学研究所统计组.常用数理统计方法.北京：科学出版社，1973：100-106

［6］龚绍京.青光台地磁短周期事件的时间序列分析.地震，1983，1：6-10

［7］郑在壮，沈瑞童.地磁短周期转换函数在地震预报中的应用.地震地磁观测与研究，2010，31（3）：13-17

［8］龚绍京.广东省地磁台的帕金森矢量及广州台的系数在河源地震前后的时间变化.地震研究，1987，10（5）：575-582

［9］Parkinson W D.The influence of continents and oceans on geomagnetic variations.Geophys.J.Int.，1962，6（4）：441-449

［10］Yukutake T.太平洋西北部地磁台阵研究的初步报告.地磁短周期变化译文集.国外地震科技情报，1987，增刊1：27-38

［11］Untiedt J.欧洲中部和南部的电导率异常.地磁短周期变化译文集.国外地震科技情报，1987，增刊1：13-20

［12］王锜.与1983年11月7日荷泽5.9级地震可能有关的荷泽台地磁转换函数异常变化.地震学报，1988，10：49-57

［13］陈伯舫.评《与1983年荷泽5.9级地震有关的荷泽地磁转换函数异常变化》.地震学报，1989，11：211-212

［1 4］龚绍京，杨桂君，田山，等.荷泽5.9级地震前后荷泽台转换函数随时间变化的研究——兼与王锜同志商榷.地震学报，1991，13：113-120

［15］陈伯舫.日本鹿屋台地磁转换函数的变化.华南地震，2003，23（1）：8-12

Comment on“Study on the relationship between the Wenchuan strong earthquake and the geomagnetic transfer function”by Baozhu Yuan et al.

Gong Shaojing，Ma Ji，Liu Shuangqing，Li Liandi
（Earthquake Administration of Tianjin Municipality，Tianjin，300201，China）

Baozhu Yuan et al.published their paper on the title by“Study on the relationship between the Wenchuan strong earthquake and the geomagnetic transfer function”in 2009.However，we do not agree on the authors about some concepts presented in their paper.The formula and the corresponding result calculated by them are questionable.Therefore，the relation of the magnitude of the transfer function to the earthquake analyzed in that paper is unimaginably queer.

complex transfer functions；Parkinson vector；electrical structure；complex least sq

P315.72＋1；

A；

10.3969/j.issn.0235-4975.2012.08.006

2011-07-25；

2012-01-13。

uare method