球谐分析方法是表示全球范围地磁场分布及其长期变化的一种数学方法。该方法由高斯于1838年首先提出。自高斯理论问世以来,地磁场的解析表达方式有了很大进展。高斯理论的目的是把地磁场表示为该点球坐标的函数(它不管地磁场形成的物理原因),虽然这一理论就其实质没有对地球磁场的成因作出解释,但是它可以解决有关地磁场的结构问题,因此得到了广泛的应用。
(一)球谐函数
磁源外部点的磁位U应满足拉普拉斯方程:
∇2U=0
勘探重力学与地磁学
在球坐标下,解的形式为
勘探重力学与地磁学
对(6-14)式来说,当n=0,就意味着单极子的场,而磁性体是不存在单极子的情况,所以,n=0的项不存在。
若以
勘探重力学与地磁学
代入(6-14)式,则U为
勘探重力学与地磁学
这就是常用的形式,式中:a表示地球的半径; 和 是内源场引起的球谐系数; 是 外源场引起的球谐系数。
(二)地磁场的球谐系数
经过分析研究认为,在地磁场观测的现有精度范围内,地球基本磁场几乎全部来源于地球内部。假定地磁场是由地球内部场源所产生的,则(6-15)式中就有 ,所以磁位可写成
勘探重力学与地磁学
由此,地磁要素的表达式可写为
勘探重力学与地磁学
其中:
勘探重力学与地磁学
勘探重力学与地磁学
以上是用球坐标(r,θ,λ)的表示,一个地磁台至少要观测这三个分量。
(6-17)~(6-19)式中:a是地球的半径();r是地心到观测点的径向距离;θ表示地理余纬度;λ从格林尼治起算的地理经度; (cosθ)是施密特归一化(即正交型)n阶m次缔合勒让德函数,它按下面形式定义:
勘探重力学与地磁学
其中:
勘探重力学与地磁学
从上述分析可看出,当知道了球谐系数 和 就可以用上述公式计算地球表面或在它的外部(r﹥a)上的任意点的地磁要素值。下面介绍求解球谐系数的方法。
当在地球表面测出了X,Y,Z三个分量的数据以后,以这些点的各磁场分量的值以及该点经度λ和余纬度θ代入式(6-17),(6-18),(6-19)式中,则在式中只有 和 为未知数。因为以上这些方程是以无穷级数形式表示的,所以在实际计算时还必须取有限项数。如果限制级数展开到n阶项,那么球谐系数 和 的总个数就为N=n(n+2)。由此欲求解出 和 ,就要有不少于N个方程式才行。这样,若是对地磁场的某一分量观测,就要有N个以上的观测数据;若是对三个分量进行观测,就要有N/3个以上的测点数据。一般在实际计算时,为了提高计算精度,取的测点数都远远大于要求的值的数,并且用最小二乘法来解出这些系数值。在求出球谐系数 和 之后,就可以反过来计算出地面上不同λ和θ点的磁场。
(三)球谐级数各项的物理意义
尽管球谐分析是一种纯粹的数学形式,但是球谐级数的每一项都有它的一定的物理意义。
1.中心偶极子磁场
在球谐级数中,当n=1时,磁位方程可写成:
勘探重力学与地磁学
另外,偶极子的磁位或者均匀磁化球体的磁位可表示为
勘探重力学与地磁学
式中:θm表示磁余纬,如图6-16所示。
图6-16 地磁场径向偶极子模型示意图
设地球表面任意观测点P的地理坐标为(φ,λ),地磁北极Nm点的地理坐标为(φ0,λ0)。N为地理北极,则在球面三角形△NmPN中,(λ-λ0)角即两个经度之差。(λ-λ0)角的两个邻边 弧与 弧的大小分别等于(90°-φ)和(90°-φ0)。因此,根据球面三角余弦定理可有
cosθm=sinφsinφ0+cosφ0cosφcos(λ-λ0)
把cosθm代入(6-22)式即得
勘探重力学与地磁学
地球磁矩m与磁化强度M的关系为
勘探重力学与地磁学
式中:a表示地球半径。令
勘探重力学与地磁学
则有
由(6-21)式与(6-24)式对比可知,当cosθ=sinφ 时,两式的形式完全相同,而(6-21)式是由球谐分析得到的,(6-24)式是直接由均匀磁化球体的磁位公式推导的。由此可见,n=1的球谐分析的磁位表达式就相当于中心偶极子的磁位。同时可看出, 相当于轴向中心偶极子的磁矩m0(指向南极),另两个相当于在赤道平面上的两个中心偶极子磁矩,即 指向东经180°, 指向东经90°(图6-17)。
图6-17 各种磁矩之间的关系
当φ0=90°,即磁轴与地理轴重合时,(6-24)式仅剩下一项:
勘探重力学与地磁学
这就是最粗略的轴向中心偶极子的磁位公式,由(6-23)式可得
勘探重力学与地磁学
两边乘以4πa3,即得地球粗略的磁矩:
勘探重力学与地磁学
当n=1时,由(6-24)式可得磁场的各分量为
勘探重力学与地磁学
以上分析可知,球谐级数的第一项,就磁位来说,相当于地球中心偶极子产生的磁位。其磁轴与地理轴一般不重合,只有φ0=90°时才重合,中心偶极子模型在古地磁学中得到广泛应用。
2.偶极子磁场与非偶极子磁场
在地磁场的球谐级数表达式中,只有n=1的项称为中心偶极子场,除去n=1的项,其余的项统称为非偶极子磁场。非偶极子场实际上是多个偶极子(或叫多极子)磁场。例如n=2项代表四极子,n=3项代表八极子等。
(四)国际地磁参考场(IGRF)
国际地磁参考场是地球基本磁场的定量描述。提出国际地磁参考场是为了计算地球基本磁场有一个约定的标准磁场,以便使各方面的研究工作结果得以统一。
国际地磁参考场是某一时期地球基本磁场及其长期变化的数学模型(由于地磁场是随时间变化的,所以有一系列的数学模型)。它由球谐(高斯)系数及相应的公式组成。1968年10月22~25日在华盛顿召开了关于“地球基本磁场的描述”座谈会,在会上正式提出了1965年的国际地磁参考场。国际地磁和高空物理协会(IAGA)的世界地磁测量部(WMS)于1968年10月28日,以及该委员会的执行委员会于1969年2月先后承认了这一国际地磁参考场(IGRF 1965)。由于地磁场不是恒定的,而且在不断地变化着(磁场的这种变化称为长期变化),所以球谐系数也是随时间系统变化的。国际地磁参考场也包括了表征地磁场年变率特点的 和 系数,即是系数 和 对时间的导数。
在国际地磁参考场(IGRF 1965)公布之后的20多年时间里,几经修改、补充,形成了十来个确定的主磁场模型。国际上规定每5年发表一次球谐系数及绘制一套世界地磁图。IGRF表示确定的地磁参考场,其高斯系数今后不再修改;每5年改变一次模型,即通过年变率的调整取得。表6-1给出了1990~1995年间IGRF球谐系数。
表6-1 1990~1995年代IGRF 球谐系数
续表
如上所述国际地磁参考场由一组高斯球谐系数( , )和年变率系数( , )组成。不同的模型,其系数一般也不同。利用各个模型的系数,算出相应的中心倾斜偶极子的磁矩值,可以了解地磁场长期变化情况以及各个模型间的衔接关系。
历代的球谐系数可以通过有关文献查到。建立国际地磁参考场的球谐系数是由准球面平均半径计算获得的,若要考虑地球形状为旋转椭球体时,则需要采用国际天文协会(IAU)的国际天体椭球坐标,取赤道半径为,扁率为1/。利用球谐系数经地心坐标转换可以求得椭球体的参考场。这对大范围的磁测是必要的。
(五)区域地磁场模型
地磁场模型分为全球模型与区域模型。区域地磁场模型是表示或者描述地球表面某一地区(如某一国家或某一大洲)地磁场时空分布的数学方法。计算区域地磁场模型的数学方法是多种多样的,通常用的方法有多项式、球谐分析、偶极源、球冠谐和分析和矩谐分析方法等。不同的国家使用不同的方法建立各自国家或地区的地磁场模型。同一国家在不同年代采用不同的计算方法。例如,加拿大学者在1985年以前使用泰勒多项式方法计算加拿大地磁场模型,1985年以后则使用球冠谐分析方法计算加拿大地磁场模型。这里介绍三种应用比较广泛的计算区域地磁场模型的方法。
1.多项式拟合法
建立区域地磁场模型应用最早,现在仍被广泛应用的分析方法是多项式拟合法。
多项式拟合法是将地磁要素以多项式表示为经、纬度的函数,或平面坐标的函数;表达式中不包括径向距离(或垂向距离)的项。这种计算方法简单易行,利用模型计算地磁场各要素比较快捷。其研究地区可达数百万平方千米,也可研究数十平方千米的尺度。多项式的阶数一般选为3左右。模型所刻画的最小波长与阶数及研究地区的大小有关,可用以下方法估计:一个n阶多项式在任何涉及的区间L内最多只有n个零点,因此用一个n阶多项式近似地表示的最小波长的估计值λ=L/[(n-1)/2]。一般在互相垂直的两个方向(如东西—南北)采用相同的阶数,但若研究地区两个方向的跨距不同,或考虑到磁场分布特点,采用不同的阶数较合理。
多项式拟合法是一种纯数学的方法,没有考虑地磁各要素间的几何约束(各要素间应满足的一定关系)以及物理约束(矢量场的旋度、散度应为零)。此外,也不能用这个模型求得观测平面上部空间的磁场。
中国、美国、日本及其他一些国家都采用或曾经采用泰劳(Tayler)多项式建立各自的地磁场模型。
美国的一些学者用多项式拟合区域地磁场,用的阶数较高,为7~9阶。由于美国本土东西长而南北短,有人在两个方向上用不同的多项式来进行拟合。
中国学者在20世纪80~90年代,利用我国历年来的地磁观测资料,以及部分国外地磁台站资料、确定的国际地磁参考场(DGRF)值,以三阶泰勒多项式拟合方法,建立了中国地区,,和的基本磁场模型,记为CHINAMF。以磁偏角D为例,地磁场磁偏角D的地理分布,用地理经度(λ)和地理纬度(φ)的泰勒多项式来表示:
勘探重力学与地磁学
式中:λ0,φ0代表展开原点的经度和纬度,并取φ0(N)=36°,λ0(E)=106°;a0,a1,…,a9是待定的系数,用最小二乘法求出。
表6-2列出了1980年中国地区基本磁场模型的系数。
表6-2 1980年中国地区基本磁场模型的系数
中国地区基本磁场(主磁场)模型CHINAMF1980与DGRF1980比较结果表明,在我国新疆北部边境和黑龙江北部边境,两种模型相差较大(如总场强度的差值可达400nT),而在其他边界地区,两种模型吻合较好[总场强度的差值为几十纳特]对于中国陆地,其精度是高的,但对于中国海域,其精度就要差一些。我国的一些学者认为,各个年代的主磁场模型可以作为研究中国地区磁异常的正常背景场,但在使用时要注意上述问题。各地的正常磁场值,可用(6-17)式、(6-18)式和(6-19)式及表示I,H,T的类似公式以及表6-2中列出的系数算出。
2.矩谐分析法(RHA)
应用矩谐分析方法建立区域地磁场模型主要是为了能够反映较短波长的地磁场特征。Alldredge于1981年指出了用球谐分析的方法表示区域地磁场的困难之后,提出用矩谐分析的方法表示区域磁异常。
矩谐分析方法采用直角坐标系。为了进行矩谐分析,首先必须将地理坐标系(地心的或大地测量的)中各个测点位置的坐标及磁场观测值转换到直角坐标系统,一般以所要分析区域的中心作为直角坐标系的原点(图6-18)。然后将磁位或场的分量展成正交的谐和函数级数(正弦、余弦及一个幂函数),用最小二乘法求出系数。一般是先将观测值减去IGRF值,得到残值,利用残值进行矩谐分析。
图6-18 矩谐分析法直角坐标系
图中P点是地面上的一点(余纬θ,东经λ),也是平面直角坐标系x,y,z的原点,Q是研究区中的任意点。
在直角坐标系中,磁位的拉普拉斯方程解,可写成
勘探重力学与地磁学
式中:j=q-i+1;v=2π/Lx;w=2π/Ly;u=[(iv)2+(jw)2]1/2;Lx,Ly是研究地区的边长;A,B,C是常数;Dij,Eij,Fij,Gij均为待定常数。磁场的分量可由下式求出:
T=-μ0∇U (6-28)
由于矩谐分析是用一个平面来近似球面,因此建模的面积不能太大,一般在(3000×3000)km2左右。此种方法在研究区域的边缘会发生振荡现象,需采取相应的技术进行处理。
Alldredge用RHA利用12个地磁台的数据分析了欧洲的地磁场,研究区为矩形区域,东西长3600km,南北宽2800km。Barton用RHA建立了澳大利亚地磁参考场(AGRF 1985),他使用了86个台站的资料,研究区东西为4800km,南北为5200km。我国徐文耀、朱岗昆也采用RHA对我国及邻近地区地磁场进行分析,矩形区域东西长8000km,南北宽6000km,使用了86个台站的数据。