0 引言

地球电场简称地电场，主要由固体地球外部电流体系在地球内部感应产生的分布于地表的快变化电场和地下介质中物理、化学作用引起的慢变化电场2部分组成（崔腾发等，2013），地球电场既受场源控制又受地下局部介质电性结构的影响，进而表现出局部性差异，其中地电场日变化是每日出现的变化成分（克拉耶夫，1954；孙正江等，1984；傅承义等，1985；李金铭，2005；杜学彬等，2007）。关于地电场随时间变化的研究取得了有益进展，如：黄清华等（2006）定量分析日本新岛台地电场固体潮汐的频谱特征及震前异常变化；张学民等（2007）探讨了昌黎台地电场固体潮汐波远震前变异现象及机理；田山等（2009）利用大地电场同一测向长、短极距比值法，排除来自大地电场及自然电场远场的变化信息，突出了与台站局部电性不均匀性有关的地震前异常信息；谭大诚等（2010）分析了与固体潮和空间S_q电流系有关的地电场变化。但是，目前关于地电场的研究工作仍主要集中在认识地电场变化波形、周期特征等方面，较少使用模型来描述其变化规律。本文在王丽等（2015）对区域地电日变化模型研究基础上，利用泰勒多项式和曲面样条函数建立江苏及邻区的地电场日变化模型，以研究江苏及邻区地电场日变化规律和局部干扰特征。

本研究欲建立江苏省及邻区地电场日变化（即随时间变化）模型，期望能基于模型进一步认识各台站地电场日变化规律和差异性，并识别夹杂在地电场正常变化中的干扰。

1 台站概况与数据选择

江苏省及相邻的山东省、安徽省现有13个地电台站长期定点观测，分别为江苏省的南京、新沂、高邮、海安地电台；山东省的郯城、安丘、菏泽、陵阳、无棣大山、邹城、乳山地电台；安徽省的蒙城、嘉山地电台。台站分布见图 1。

由于地电场观测环境和装置系统（电极、外线路、仪器等）的长期稳定性会影响观测数据质量，为此采用以下原则筛选地电观测数据：①地电观测数据质量合格；②考虑到地电日变化对潮汐响应的逐日变化，选取各台站2013-2017年相同农历日或前后相差一天的NS测道和EW测道观测数据，作为拟合模型的样本数据及用模型描述的实测数据，磁情指数分别按K_p = 1、2，K_p = 3、4和K_p = 5取值。

2 理论模型

当前我国地电场仅观测水平分量，无垂直分量，不能使用地磁学中的球谐函数等方法建立三维模型，故采用泰勒多项式（Taylor Polynomial，文中缩写为TP）和曲面样条函数（Surface Spline，缩写为SS），选定区域若干台，按地电场单分量E_x（北分量）、E_y（东分量）随时间变化的二维空间（即平面）建立模型，并描述地电日变化。地电日变化明显响应太阳、月球对地球的潮汐作用，按幅度大小，主要日变化成分依次为半日波（主太阴半日波M₂波周期和主太阳S₂波周期）、全日波（太阳日周期）（黄清华等，2006；叶青等，2007；谭大诚等，2010；崔腾发等，2013）；地电日变化存在逐日相位差和半月潮现象，且受到磁情天气影响。因此，文中以一年中农历日期的地电场原始观测曲线作为样本曲线建立模型，并采用拟合模型曲线来描述不同年份或同年份中同等农历日期、同等磁情指数的地电日变化。地电场变化对地下介质电导率结构响应敏感，因此，采用台站某分量（即某测道的测量数据）建立的日变化模型，只适用于描述该台站该分量在同农历日期、同等磁情指数下的地电日变化。

2.1 泰勒多项式方法

应用泰勒多项式方法，安振昌等（1994）建立中国地区地球主磁场模型和东亚地磁场模型，柳士俊等（2011）建立区域地磁场模型。本文将此方法用于地电场研究中，建立地电日变化模型。泰勒多项式基本公式为

$ F = \sum\nolimits_{n = 0}^N {\sum\nolimits_{m = 0}^n {{A_{nm}}{{\left({\varphi - {\varphi _0}} \right)}^{n - m}}{{\left({\lambda - {\lambda _0}} \right)}^m}} } $

(1)

式中，F表示测点（即台站）的地电场某分量测值；φ和λ代表地电场台站的地理纬度和经度，φ₀和λ₀为泰勒多项式模型展开原点的纬度和经度，展开原点一般选在区域中心；N为多项式截断阶数，共有(N + 1)(N + 2)/2个模型系数；A_nm为多项式模型系数，通过测量点所选样本测值及φ、λ，采用最小二乘法求得。

2.2 曲面样条函数方法

安振昌等（1982）应用曲面样条函数建立局部地区地磁场长期变化模型。曲面样条函数方法基本公式如下

$ w(x, y) = {a_0} + {a_1}x{\rm{ + }}{a_2}y\sum\nolimits_{i = 1}^n {{F_i}r_{^i}^2\ln \left({{r^2} + \varepsilon } \right)} $

(2)

式中，w(x，y)是地电场某分量测值；a₀、a₁、a₂与F_i为待定系数；x、y为台站纬度和经度；$r_{^i}^2{({x_i} - {x_0})^2} + {({y_i} - {y_0})^2}$，其中x₀、y₀为所选区域中心点的纬度和经度；ε为控制曲面曲率变化的经验性调节参数，视实际情况适当选取。当曲面曲率变化较大时，ε取值较小，反之，取值较大；对平坦曲面，一般取ε = 1-10⁻²。设N个测点的地电场同一分量观测数据为w_j，对应测点经纬度坐标为(x_i，y_i)，则建立曲面样条函数模型的过程即为采用最小二乘法确定待定系数的过程。

3 计算结果分析

以一年中农历日期的地电场原始观测曲线作为样本曲线建立模型，并采用拟合模型曲线来描述不同年份或同年份中同等农历日期、同等磁情指数的地电日变化，以期能够建立江苏及邻区的地电日变化模型。用泰勒多项式和曲面样条函数得到的地电场变化模型曲线几乎一致，故文中仅给出泰勒多项式计算结果。

3.1 模型拟合

对于所选的地电台站数据，分别对比分析不同K_p指数下（K_p = 1、2，K_p= 3、4，K_p= 5）泰勒多项式模型拟合精度。

在K_p = 1、2时，选取蒙城、南京、新沂台NS向2016年1月26日-27日及2017年1月14日-15日地电数据；在K_p = 3、4时，选取高邮、海安、郯城和新沂台EW向2013年9月25日及2015年10月2日-3日地电数据；在K_p= 5时，选取安丘、高邮、无棣和新沂台NE向2016年3月14日-15日地电数据。计算以上所选数据每日模型曲线与样本曲线误差（差值绝对均值），结果见表 1和图 2、图 3、图 4。

由表 1可见：当K_p = 1、2时，3个台站平均误差为0.03-0.08 mV/km，小于1.7%（见表 1平均值一栏数据，代表每日模型曲线与样本曲线误差，即表中的差值绝对均值比日均值）；当K_p= 3、4时，4个台站平均误差为0.09-0.16 mV/km，小于2.5%；当K_p= 5时，4个台站平均误差为0.15-0.43 mV/km，小于3.3%。因此，随着K_p指数的增加，任一时段不同方向泰勒多项式拟合曲线精度越来越差。

由图 2、图 3、图 4可见：①在K_p = 1、2时，地电日变化较平稳，拟合模型曲线相对平滑，图 2所示，蒙城、新沂台2峰1谷地电日变化形态与南京台1峰1谷的变化形态，拟合曲线更加清晰；②当K_p = 3、5时，如图 3中4个台站样本曲线在不同日期存在不同程度的地电脉动，使用泰勒多项式拟合，可消除部分局部干扰，但拟合曲线没有图 2光滑；③当K_p = 5时，拟合曲线可大致显示地电日变化动态，但局部干扰并未被完全消除。

由图 2-图 4及表 1可见，在不同磁情指数下，各台站泰勒多项式模型曲线与样本曲线在变化波形和幅值大小上均较吻合，变化动态一致，但模型曲线更加光滑，可见泰勒多项式拟合模型曲线可代表台站实际观测值，且随着K_p指数的降低，模型拟合精度依次增高，说明在空间电磁环境平静时，选用K_p指数较低的样本观测值来拟合某地电台某测项分量的模型曲线，可作为其正常背景曲线。

3.2 模型验证

为验证不同K_p指数下泰勒多项式模型拟合曲线的可靠性及实用性，选用不同年份相同农历日期同测道地电场原始观测值，与图 2-图 4中相同条件下的泰勒多项式拟合曲线进行对比，结果见图 5、图 6、图 7，计算与上述拟合数据差值（除去建立模型时所用观测值），得到差值绝对均值，结果见表 2。差值计算公式如下

$ \Delta F = {F_E} - {F_{\rm{0}}} $

式中，F_E为地电场实际观测值，F₀为模型拟合值。

当K_p=1、2时，蒙城、南京和新沂台NS测道短极距分钟值2014年1月17日-18日（农历腊月十七至十八）原始曲线与图 2中给出的泰勒多项式拟合曲线进行对比，结果见图 5，可知3个台站实测曲线与上述模拟地电日变化波形较吻合，均呈现2峰1谷的变化形态。由表 2可见，差值绝对均值除2014年1月18日新沂台为2.1外，其他均在1.2上下变化，说明3个台站观测环境、观测系统等稳定性在2-3年内未发生明显改变，但应略有差异，表现在其差值绝对平均值比表 1略大。

当K_p= 3、4时，高邮、海安、郯城和新沂台EW测道短极距地电场日变化分钟值2013年9月24日（农历八月二十，K_p= 3）及2016年9月21日（农历八月二十一，K_p= 4）与原始曲线，与图 3中给出的泰勒多项式拟合曲线进行对比，结果见图 6，可见4个台站实测地电场日变化与上述拟合结果存在一定偏差，差值绝对值在3.0 mV/km上下波动（表 2）。由图 6(a)可见，在农历八月二十，1 100-1 440 min时段拟合与原始数据同步差值较大，查询观测日志和工作日志后，确认为外空地电脉动现象所致，其中新沂台因当日暴雨及电极老化干扰叠加，导致全天差值相对较大。由图 6(b)可见，4个台站实测曲线与上述拟合地电场日变化波形偏差更大，除高邮台在3.0 mV/km上下波动外，其他3台均大于4.0 mV/km，其中新沂台达6.34 mV/km；在400-600 min和800-1 100 min时段内，三峡至常州高压直流输电线路因调试过程中电流不平衡对地漏电，对4个台站造成干扰，导致同步差值较大。

当K_p= 5时，利用安丘、无棣、高邮和新沂台NE测道长极距地电场分钟值2012年2月27日（农历二月初六）实测地电场日变化曲线，与图 4中给出的泰勒多项式拟合曲线进行对比，结果见图 7。由图 7可见，4个台站实测地电场日变化波形与上述拟合曲线基本吻合，但表 2显示，其差值绝对值均大于3.13 mV/km，与K_p= 1-4时的差值相比，明显偏大，与地电暴时段地电场存在高频突跳现象具有一致性。

4 结论及讨论

通过上述分析，可以得到以下结论。

（1）使用不同年份、同农历日期、同等磁情天气的样本地电场数据建立拟合模型，其中泰勒多项式所建立的拟合模型较为可靠，在磁静日等同等条件下建立的地电日变化模型是可行的。随着磁情指数增大，误差越来越大，因此，在建立模型时应尽量采用磁情指数较小的平静日地电数据作为样本数据。

（2）使用分钟值曲线建立的地电场日变化模型可反映地电场日变化特征，2峰1谷和1峰1谷的变化形态清晰（叶青等，2007），可作为正常背景变化模型。使用模型曲线描述不同年份、同农历日期、同等磁情天气的实测地电场日变化，模型曲线与实测曲线误差较小，且能反映一天当中台站观测环境、测量系统等引起的地电变化。

（3）使用同一个区域地电台网中不同台站同方向测道的样本曲线建立各测道地电场变化模型，可代表各测道地电场变化自身特征，也能反映同方向测道样本曲线的共性特征。

各测道模型曲线与样本曲线的一致性表明了各测道地电场变化的自身特征；在计算模型系数解线性方程组的求逆过程中，各测道地电场观测数据对模型系数计算是相互关联的，此关联性可反映区域地电台网各测道地电场变化的共性特征，即区域地电台网同方向测道空间上的共有特征。