引言

磁异常场是地壳和上地幔顶部所有低于居里点温度的磁性岩石所产生的磁场（也称为岩石圈磁场或地壳场）（Purucker，Whaler，2007）。磁异常场相对主磁场可以被认为是静态场，在地球上不同位置的强度不一。在近地表高度，磁异常场的幅值一般为几十到几百nT，个别地区如俄罗斯库尔斯克磁异常区和非洲中部的班吉磁异常区可达数千nT。对全球和区域磁异常场的系统研究在板块假说、地质构造、航海、采矿勘探等领域的理论和应用方面具有重要价值，从而对理解地球的组成和演化具有重要意义（冯彦，2011）。

科研中通常使用固定地磁台、流动地磁、海洋、航空和卫星磁测数据建立模型并研究磁异常场。最近几十年除了地磁台站、流磁复测点及航空磁测外，使用较多的是以高密度、高精度、全天候为优点的卫星磁测数据。1999年Ørsted卫星（Neubert et al，2001）的发射拉开了现代卫星矢量地磁场测量的帷幕，从2000年起陆续发射了SAC-C，CHAMP （Reigber et al，2002）卫星和Swarm星座（Friis-Christensen et al，2006）用于测量地磁场，随之一系列重要的全球模型不断出现，如综合模型（comprehensive model，缩写为CM）系列（Sabaka et al，2004，2015，2020）、卫星岩石圈磁场模型（lithospheric magnetic field model，缩写为MF）系列（Maus et al，2007，2008）、融合了卫星及近地表磁测的NGDC-720模型（Maus，2010）、CHAOS模型系列（Finlay et al，2015，2016，2020）、LSC-1 （Olsen et al，2017）以及以近地面磁测为主、长波长部分由卫星模型代替的EMAG2v3模型（Meyer et al，2017）等。最新由Thébault等（2021）建立的地磁场模型，截断阶数高达1050阶，能计算得到约38 km波长的磁异常场。这些地磁模型具有的截断阶数更高、空间分辨率更高、内外源场分离更合理等优点。然而从对区域磁异常场的研究角度而言，则需更加注重地表或者岩石圈浅层的磁异常信息，该信息主要由地表处的地磁短波长部分组成，这些信息会随着高度呈r⁻³ （r为场点到源点的距离）衰减，卫星数据都是在磁层电离层之间（高度约300—800 km）获得，因此在卫星高度所测得的仅为大尺度、场源较深部分的信息。尽管一些全球模型如NGDC-720等是结合了地表、卫星数据等进行建模，但是地面数据分布严重不均，一些小尺度的磁异常信息并不能够均匀地覆盖测量，因此结合研究区域内的高密度近地面数据和区域模型进行中小尺度的磁异常研究是一种可行的研究方法。

针对卫星数据及所建模型的局限，利用区域分布稠密的地面测点，结合各种经验或者半经验的区域模型是研究区域磁异常的常用方法。常用的区域模型有球冠谐模型（Haines，1985）、泰勒多项式模型（Le Mouel，1969）、矩谐模型（Alldredge，1981）等。通过上述模型，安振昌（2003）和高金田等（2006）对二十世纪九十年代以来中国大陆的磁异常场进行了深入研究。然而这些工作仅考虑了地面数据，若联合卫星数据建模则可更好地获取磁异常场的三维空间特性，例如：Ou等（2013）基于CHAMP卫星数据和改进的球冠谐模型研究了中国地区的岩石圈磁场模型；Feng等（2016a）基于CHAMP卫星数据及球冠谐模型反演了卫星高度的岩石圈磁场，还利用三维泰勒多项式（three-dimensional Taylor polynomial，缩写为3DTP）（Feng et al，2016b）和三维曲面样条（three-dimensional surface spline，缩写为3DSS）模型（冯彦等，2018）研究了中国地区的岩石圈磁场。这些工作为区域三维磁异常场的研究奠定了基础。

一般认为磁异常场为静态场，基本不随时间发生变化，然而我国1960—2000年期间每隔10年左右进行一次的全国性地磁测量显示，由于每次所使用的测量仪器、测量技术、测量点的数量和位置分布不均匀，或建模技术不同，不同时期的区域磁异常存在差异，而这些差异也是我们需要关注的方面。由于泰勒多项式模型采用最小二乘法进行拟合，适用于观察一些中等尺度磁异常的分布及变化情况，本研究拟通过中国大陆各年代的磁测数据，利用泰勒多项式模型，结合CM4模型（Sabaka et al，2004），计算得到各年度磁异常的分布，以分析磁异常差异，并研究造成差异的原因。为了进一步分析结果，将利用最新的3DSS模型结合地面数据来比较分析。

1.   建模数据和模型

1.1   建模数据

本研究使用中国大陆地区（18°N—54°N，73°E—136°E）地面实测数据进行建模。中国科学院地质与地球物理研究所和中国地震局分别在1960，1970，1980，1990和2000年左右进行全国范围的地磁数据测量，并将数据分别归算至1960.0 （即1960年1月1日0点，下同），1970.0，1980.0，1990.0和2000.0年（冯彦，2011）。各年代的测点分布有所不同，主要以磁偏角D、磁倾角I、磁场水平分量H和总强度F为主。此外，为了考察磁异常的整体分布，我们还将1960.0—1990.0年的测点归算至2000.0年，进行整体研究，具体测点分布见图1。

图  1  中国大陆地区在1960—2000年的测点分布

Figure  1.  Distribution of surveying points over Chinese mainland during 1960−2000

下载: 全尺寸图片 幻灯片

由图1可见，1960.0，1970.0，1980.0，1990.0和2000.0年的测点分别为445，1887，255，139和156个，共2 882个地面磁测点。这些数据点主要分为地磁台站和复测点两类。对于各年代的重复测点的选取，考虑到CM4模型相比于IGRF采用了更多的高精度数据，同时对于主磁场有更高的截断阶数，能更准确地表示主磁场（Maus et al，2009），故利用每个实测点的磁场强度F与该点CM4主磁场强度的绝对差值进行选取，保留绝对差值较小的测点；然后考察北向分量X、东向分量Y、垂直分量Z以及磁偏角D和磁倾角I的残差（实测值减去对应CM4模型的主磁场值），若残差大于500 nT或0.5°，则将该异常点剔除（冯彦，2011），最终分别得到405，1761，199，128和138个无重复的磁测数据。另外，为了考察磁异常的整体分布，假设岩石圈磁场不随时间变化，将1960.0，1970.0，1980.0，1990.0年代的测点基于IGRF12归算至2000.0年，剔除重复点后，共得到2494个有效数据点。

1.2   模型

地磁场在三维空间中是关于时间和位置的连续分布函数，可通过泰勒函数近似计算（Le Mouel，1969），其表达式如下：

式中：$ W $为任一地磁分量；$ N $为模型最大截断阶数；$ \varphi，\lambda $分别为各磁测点的地理纬度和经度；$ {\varphi _0} $和$ {\lambda _0} $为多项式模型的展开原点（其它各点相对此点进行计算），一般为研究区域中心点，也可为研究区域内的任意一点，将所研究地区的原点设为经纬度方向的中点，即$ {\varphi _0} ＝ {36^ \circ }{\text{N}}，{\lambda _0} ＝ {104.5^ \circ }{\text{E}} $；$ {A_{nm}} $为模型系数，每个模型共有$ （ {N ＋ 1} ） （ {N ＋ 2} ） /2 $个系数，考虑到测点数远多于待定系数，故利用最小二乘法进行拟合。

泰勒多项式模型具有计算简单、使用方便的特点，能够较好地表示区域地磁场中等尺度的空间变化特征（冯彦，2011）。由于模型基于最小二乘法拟合计算，会对一些小尺度的磁异常分布进行平滑拟合，因此该模型并不适用于反映小尺度空间的磁异常变化，而更适用于反映地核场以及部分中长波长岩石圈磁场（岩石圈中层及深层部分）分布。泰勒模型通过最小二乘原理进行拟合建模，理论上，测点越多建模效果越好，因此适合在测点稠密区进行区域或者局部建模。

泰勒模型存在的问题是截断阶数的选取：阶数太高会导致计算缓慢，产生龙格现象；阶数太低则会产生较大的拟合误差。本研究结合以往的经验（Feng et al，2015），分别计算了模型1—10阶的模拟值及其与对应CM4模型主磁场的均方偏差，并观察其变化趋势，在数值稳定的前提下，还考虑到所建模型需要反映区域磁异常而不能太过平滑，因此综合考虑了偏差的稳定性和数值变化，最终截断阶数选为7阶。

区域模型存在的另一个问题是边界效应，如果不进行控制则会产生边界畸变。因此在不同年代的境外均匀选取了约30个CM4模型的磁异常值作为补充点来控制边界异常。对于数据归算后的整体分布，由于数据较多，选取150个境外测点进行控制。

为了验证结果，我们将泰勒模型计算所得的磁异常值与CM4模型和3DSS模型的磁异常场相比较。CM4模型最大的特点是能够进行不同地磁场源如地核场、岩石圈磁异常、电离层场、磁层场及感应场的估算并分离，因此可直接用来移除场源并比较分析。3DSS是区域三维磁场建模，主要在研究空间经过的所有测点进行插值，并取得了较好效果。上述两种模型的建模方法可参阅相关文献（Sabaka et al，2004；冯彦等，2018）。

2.   结果

现通过7阶泰勒模型分析1960.0，1970.0，1980.0，1990.0以及2000.0年的中国大陆地区岩石圈磁场分布情况。所有年代复测点的日变化及其它扰动均基于临近台站进行了校正。将各年代实测值减去对应年代的CM4主磁场值后得到磁异常场值，利用泰勒模型对磁异常场进行反演得到模型系数，然后正演得到研究区域内1°为间隔的网格值。

2.1   1960—2000年的X分量磁异常分布

各年代X分量磁异常分布如图2a所示，可见：1960.0年X分量磁异常除西藏西南部和新疆西部外，中国大陆地区皆为负值，磁异常强度随纬度增加，由−120 nT增加到−20 nT左右，极值点出现在贵州地区，强度为−170 nT；1970.0年负值区继续扩展，中国大陆皆为负值，极值点出现在四川北部地区，强度为−90 nT左右，强度随经度增加从−40 nT减至−80 nT左右；到1980.0年负值区继续向西部移动，极值点出现在西藏西南部，强度约为−240 nT，东部的正值区开始出现，以105°E附近的零值线为界几乎东西各半；1990.0年负值区域进一步向北及东北部扩展，而正值区则向西南移动，磁异常强度随纬度增加而降低，从20 nT降至−100 nT左右；及至2000.0年，中国大陆的负磁异常向西南移动，强度有所减弱，而在东海附近的正磁异常则向东北及南部省份移动，极值点出现在湖南、广东交界处，强度为100 nT左右。

图  2  基于泰勒多项式模型的X分量（a）和Y分量（b）的磁异常分布图

Figure  2.  Geomagnetic anomaly distributions of component X （a） and Y （b） based on Taylor polynomial model

下载: 全尺寸图片 幻灯片

2.2   1960—2000年的Y分量磁异常分布

各年代磁异常Y分量分布如图2b所示，由该图可见：1960.0年的Y分量磁异常在中国大陆地区除了东北（黑龙江、吉林和内蒙古）、新疆以及南方部分区域外，都以正磁异常值为主，强度随经度由0 nT增加到40 nT左右；从1970.0年起南方的负磁异常开始与东北方的异常相联结，正值区有向西移动的趋势，强度随经度增加由60 nT降为0 nT左右；异常经过十年（1980.0年）后变化不大，只是黑龙江地区的负极值点的强度有所增强，中国大陆负值区有所减小，极值点位于四川与陕西交界处；磁异常随着时间逐步扩展，及至1990.0和2000.0年，负值区大幅扩展，而位于新疆北部的负值区变化不大，原先位于东北和中部的负值区联结起来并逐步向东北移动，呈现东西向的条带分布。

2.3   1960—2000年的Z分量磁异常分布

各年代磁异常Z分量磁异常分布如图3所示，由该图可见：Z分量在整个研究期间都以负值为主，强度变化不大；两个较为显著的磁异常极值点主要出现在东北黑吉辽与内蒙古交界处以及1960.0年的新疆与西藏交界处，而前者的位置几乎不变，后者随着时间向北移动，并固定在新疆与俄罗斯交界处；各年代的Z分量磁异常基本都随经度增加而降低，幅度在±200 nT之内。

图  3  基于泰勒多项式模型的Z分量磁异常强度分布图

Figure  3.  Geomagnetic anomaly distributions of component Z based on Taylor polynomial model

下载: 全尺寸图片 幻灯片

2.4   基于区域模型与全球模型和区域三维模型的分析比较

为了考察整个研究期间磁异常场的分布趋势，现将1960.0，1970.0，1980.0，1990.0年基于实测值的岩石圈磁场值按照IGRF12模型的长期变化统一归算至2000.0年，对于异常点和重复点按照以下准则进行筛选：利用每个测点的总强度与对应于该点的CM4地核场强度的绝对差值进行重复测点筛选，保留绝对差值较小的测点，接着剔除残差相对于CM4模型X，Y，Z分量大于500 nT以及磁偏角D及磁倾角I大于0.5°的异常点。最后得到的有效数据点共2494个（冯彦，2011）。建模时需考虑截断阶数的选取和边界效应的控制两个问题，综合考虑均方偏差和磁异常形态，截断阶数依然选为7，通过均匀添加150个境外补充点降低了边界效应。最后利用7阶泰勒多项式模型模拟并绘制了1960.0—2000.0中国大陆地区磁异常场7个量的整体分布图（图4），并将其与CM4模型绘制的7个量的磁异常图相比较。

图  4  基于泰勒模型 （左）和CM4模型（右）的磁异常分布

（a） X分量；（b） Y分量

Figure  4.  Geomagnetic anomaly distribution based on Taylor polynomial model （left） and CM4 （right）

（a） Component X；（b） Component Y

下载: 全尺寸图片 幻灯片

图  4  基于泰勒模型 （左）和CM4模型（右）的磁异常分布

（c） Z分量；（d） 总强度F；（e） 水平分量H；（f） 磁偏角 D；（g） 磁倾角I

Figure  4.  Geomagnetic anomaly distribution based on Taylor polynomial model （left） and CM4 （right）

（c） Component Z；（d） Total intensity F；（e） Horizontal component H；（f） Declination D；（g） Inclination I

下载: 全尺寸图片 幻灯片

图4显示基于7阶泰勒模型的磁异常分布相较于CM4模型略显简单，磁异常较少。首先考察泰勒模型分布情况（左列），X分量整体为负值分布，并在青海、四川、西藏等地区的交会处有强度约为−100 nT的负极值点，强度随经度增加由−20 nT降低为−40 nT左右；Y分量以105°E附近的零值线为界，正负异常基本东西对半，强度随经度增加由40 nT降低至大约0 nT；Z分量和总强度F不管强度还是分布都高度一致，为负值，极值点在青海附近出现，为−150 nT左右，强度随纬度增加由−50 nT降低为−100 nT左右；H分量的分布变化与X分量类似；磁偏角D则与Y分量类似，强度随经度增加而降低；磁倾角I在大陆地区也以负值为主，在东北、西北和西南少数地区存在正值，强度随经度增加由0.5°降低至−0.5°左右。

基于两模型绘制的7个量分布显示，X，Z，H分量和总强度F有较大差异，主要体现在分布和强度上，这是由于利用泰勒多项式模型建立的是低阶区域模型，不满足正交性。图4显示中国大陆地区主要以负磁异常为主，尤其是中部和东部地区都为负磁异常。两模型的Y，D和I的分布具有一致性，比较可以看出泰勒模型的拟合尺度更大，而CM4模型能反映约600 km波长的磁异常，且卫星数据很密集，因此中大尺度的磁异常分布较为均匀。仔细观察也可以看出两种模型的相似之处，如泰勒模型的Y分量在中国大陆地区南部的分布形态为正负相间，在福建东南部和西藏大部分地区为正值，而在南部大部分地区为负异常，北部的新疆北部和黑龙江大部仍为负值，只是分布形态略为简单；磁偏角D与Y分量变化相似；泰勒模型的磁倾角I在新疆南部有正异常分布，且与北部的负异常相邻，这个特征从CM4模型中也可看出，新疆东北部的正异常也具有一致性。从磁异常东北方向的分布看，泰勒模型以负值为主，而CM4模型则表现为正负相间。整体而言，由于两者的建模数据和分辨率存在差异，相关性不高，但从Y，D，I的分布来看，两者在局部磁异常方面存在一定的相关性。

为了进一步考察所绘制磁异常的真实性，我们选用冯彦等（2018）最新开发的三维曲面样条3DSS模型，结合地面及航磁数据（Xiong et al，2016）进行总强度F的计算，具体为1936—2000年的地面数据共3143个，航磁数据3064个，一共6207个数据点进行三维建模，绘制了距地面0.85 km高度处（所有地面测点的平均高度）的磁异常分布，并与泰勒模型和CM4模型比较，结果如图5所示。

图  5  总强度F的磁异常分布

（a） 3DSS模型；（b） 泰勒模型；（c） CM4模型；（d） 3DSS模型的实测点分布

Figure  5.  Geomagnetic anomaly distribution of total intensity F

（a） 3DSS model；（b） Taylor model；（c） CM4 model；（d） All modeling points of 3DSS

下载: 全尺寸图片 幻灯片

从图5可以看出，三种模型绘制的总强度F分布各有差异。3DSS模型在国内以负值区居多，图中的白色部分为超过1000 nT的异常。3DSS模型由于结合了地面数据和中国地区的航磁数据，能反映出研究地区的中小尺度磁异常分布，具有较高的可信度。3DSS模型在国内以小尺度磁异常为主，除了在新疆南部和西藏北部有大尺度的正异常，其它地区以负值居多，其中云南及其北部和南部都为负磁异常区，这与泰勒模型一致，而在东北三省和湖南、湖北、广西、浙江、江西等地区也都是负异常，但尺度都不大，这是样条函数反映的小尺度磁场信号的特点。三维模型与CM4模型的对应较差，说明地面测点确实能在二维及三维模型中体现磁异常信息，如果结合更多近地面的磁场数据，则可以反映更多地表浅层的磁异常信息。

3.   讨论与结论

通过1960—2000年期间以十年为间隔的地磁测量数据，结合经典的泰勒多项式模型方法分析了各年代的磁异常场的分布，可以看出测点的不同相应导致了磁异常场的不同，之后将泰勒多项式模型结果与CM4模型相比较，并与结合了地面及航磁数据的3DSS模型进行验证，得到以下结论。

通过各年代比较，分析了泰勒多项式模型X，Y，Z分量在不同年代的磁异常变化趋势，表现为：X分量随时间呈现强度降低，分布向西南收缩的趋势；Y分量磁异常在中国大陆呈现正值区向西部移动、负值区逐步占据优势的趋势，强度由正变负，逐步减弱；Z分量在中国大陆呈现负值为主，形态较为稳定，两个位于东北和新疆西藏交界处的极值点几乎固定不变。综合来看，X和Y分量分别平行于地表，受到测点、设备及测量因素的影响而产生了磁场变化，而径向分量Z若是测量方法和手段差异不大，则其形态分布不会产生较大变化，这也解释了为何大多数卫星磁异常都用Z分量来描述。X和Y分量随时间变化可能有两点原因，一是每次测量的数据点位数量和位置不同，数量相差近十倍，这会导致不同年代反演时相同或相近测点的权重不同，从而直接影响建模分布；二是由于板块运动，如地震、火山活动导致短期内的地质变化，在局部地区产生新的磁性岩石，从而改变区域磁性岩石的形成和分布，导致了磁异常场的变化。

将所有数据归算至2000.0年进行建模绘图（图4），分析表明建模测点的数量和高度不同导致了模型反演和正演的不同。三种模型由于分辨率差异导致了一致性很低。泰勒模型的X，Z，H分量和总强度F几乎均为负值，与CM4模型值存在一定差异，而Y分量、磁偏角D、磁倾角I与CM4模型具有一定的相关性。基于磁异常结果，考虑到3DSS模型、泰勒多项式模型和CM4模型的各自特点，我们认为磁异常的差异主要由数据点数量、位置和建模方法所致。

从本研究所得到的不同年代的磁异常分布比较结果来看确实存在问题，如磁异常的一致性问题。从地质构造角度看，短时期（50年左右）内研究区域不会发生明显的变化，也不会产生明显的磁异常变化，但是局部地区的剧烈地质运动，如汶川、玉树地区地震等，使得对应地区的板块挤压而改变由居里等温面变化导致的岩石圈磁性，从而产生小范围磁异常强度和范围的变化；另外假设建模数据是干净的，故可认为不同时期的磁异常，尤其是X和Y分量的变化，主要是由测点的数量和位置的差异所导致，而不同的建模方法也是导致差异的原因，不同的区域模型对于磁异常拟合的差异并不明显。以上说明除了区域地质活动的剧烈变化，局部地区磁异常研究结果的准确性和精度主要取决于测点的质量和位置。区域磁异常存在连续的变化既有外部原因，如测量仪器、测量手段、环境干扰等，也有一定的内部原因，如地震等。

本研究主要考察的区域磁异常与局部地区的磁性岩石密切相关，而磁性岩石又与地质构造、居里等温面的深度相关，本研究可作为研究震磁关系的一种参考。考虑到我们使用近40年的中国地区稠密的实测点进行区域建模，在研究区域的地面测点密度远超各种全球模型，泰勒模型本身的特点可以基本反映出中等尺度的磁异常信息，因此具有一定的研究意义和实用价值，而对于一些小规模且磁场强度较大的磁异常，可以通过航磁数据和台站等验证。