将式（3.1-36）改写成RN?RecosL1?(1?e)tanL122，再代入上式，得
??Ree2tanLtanL??z?? （3.1-42）
2222x?y??1?(1?e)tanL??如记t?tanL，则由上式可构造出求解纬度正切值的迭代公式，如下
??Ree2ti1?z?? （3.1-43） ti?1?2222x?y?1?(1?e)ti???令迭代初值t0?0，一般经过5～6次迭代便可达到足够的数值计算精度，再求反正切即可获得纬度L。
最后，根据（3.1-39）求解高度，得
x2?y2h??RN （3.1-44）
cosL式（3.1-38）、（3.1-43）和（3.1-44）即为由(x,y,z)求解(?,L,h)的算法。
3.2地球的正常重力场
在地球的大地水准体描述中，水准体表面是地球实际重力场的一个等位面，每一点的重力方向均与该点所在等位面相垂直，实际的重力方向一般称为天文垂线，或称真垂线。由于实际地球内部密度分布不均匀，并且表面凹凸不平，大地水准面不规则，因而实际重力的大小和方向也不规则。与地球的几何形状描述类似，也希望使用一个简单的数学函数来拟合地球重力场，这个简单函数表示的重力场就称为正常重力场。
重力是地球万有引力和离心力共同作用的结果，参见图3.2-1，在P点处，重力矢量G是引力矢量F和离心力矢量F?的合力。
ωieP图3.2-1 重力矢量
1．圆球假设下的地球重力
若将地球视为圆球体并且认为密度均匀分布，那么地球引力指向地心，根据牛顿万有引力定律，地球对其表面或外部单位质点的引力大小为
F?GM??2 （3.2-1） r2r其中，G为万有引力常数，M为地球质量，记?=GM为地心引力常数，r为质点至地心的距离。
由于地球绕其极轴存在自转角速率?ie，使得与地球表面固连的单位质点受到离心力的作用，其大小为
2F???ieRcosL （3.2-2）
其中，R为圆球半径，L是地理纬度（在圆球假设中即为地心纬度）。重力是引力与离心力的合力，引力与离心力之间的夹角为π?L，根据余弦定理，在纬度为L的地方上P点的重力大小为
46
gL?F2?F?2?2FF?cos(π?L)22?F2?(?ieRcosL)2?2F?ieRcos2L?(F??R)?2F?RsinL?(?R)sinL22?(F?1/2?ieR)?ieR2??(F??R)?1?2sinL?22(F??R)ie??2ie22??ieR2??ge1?2sinL?ge?1?sinL?gege??22?ieR2ie22ie22ie22 （3.2-3）
22其中，ge?F??ieR为赤道上的重力大小，?ieR/ge为赤道上的离心力与重力加速度的比值。
实践表明，基于圆球假设的重力公式（3.2-3）与实际椭球地球相比，在高纬度地区偏小将近2mg，部分原因归结于实际椭球地球在高纬度地区半径缩小，实际引力增大。为了更精确地计算正常重力值，需要在椭球条件下进行重力推导。 2．旋转椭球假设下的地球重力
对于地球旋转椭球体，假设在椭球表面上重力矢量处处垂直于表面，也就是说，旋转椭球表面为重力的一个等位面，意大利人索密里安（Somigliana）于1929年经过严密推导（过程比较复杂，从略），获得了旋转椭球体的正常重力公式，如下
Regecos2L?Rpgpsin2L （3.2-4） gL?2222RecosL?RpsinL其中，Re和Rp分别为旋转椭球的赤道长半轴和极轴短半轴，ge和gp分别为赤道重力和极点重力，gL为地理纬度L处椭球表面的重力大小。 对于赤道重力ge和极点重力gp，近似有
?33?g?(1?m?mf)e?RR27?ep （3.2-5） ??g??(1?m?6mf)p?Re27?其中，f为旋转椭球几何形状扁率，约等于1/298；并且
m?2?ieRe?/(ReRp)?2?ieRe （3.2-6）
ge为赤道上的离心力与重力加速度的比值，约等于1/288。
记地球重力扁率为
??gp?gege （3.2-7）
其实际值约等于1/189。将式（3.2-5）代入上式，并展开成关于m和f级数形式，可得
66R(1?m?mf)(1?f)(1?m?mf)pgp?gegp77????1??1??13333gegeRe(1?m?mf)1?m?mf2727163333???(1?m?f?mf?mf2)?1?(m?mf)?(m?mf)2???1 （3.2-8）
772727??51715?m?f?mf?m2?2144式（3.2-8）建立了重力扁率?与椭球形状扁率f之间的关系，若忽略高阶小量，近似有
5??m?f （3.2-9）
2
47
上式是利用重力测量方法确定地球形状参数f的基本公式，它最早在1743年由法国数学家克莱罗（A.C. Clairaut）求得，通常称为克莱罗公式。克莱罗在推导公式（3.2-9）时作了如下前提假设：地球是由密度不同的均匀物质层圈组成的椭球体，各椭球面都是重力等位面，且各层密度由地心向外有规律的减小。需要说明的是，克莱罗公式是近似成立的公式，而索密里安公式却是理论上严格成立的，并且后者的前提条件更加宽松，对椭球体密度分布不做任何限制。
若将gp?(1??)ge和Rp?(1?f)Re代入索密里安公式（3.2-4），展开为关于?和f的级数形式，可得
gL?gecos2L?(1?f)(1??)sin2LcosL?(1?f)sinL22222?1?(2f?f)sinL???23?122222????ge?1?(??f??f)sinL1?(2f?f)sinL?(2f?f)sinL???2???8??ge1?(??f??f)sin2L1/2??? （3.2-10）
??1??ge?1??(??f??f)?(2f?f2)?sin2L2???1?3???(2f?f2)2?(??f??f)(2f?f2)?sin4L?2?8?11??ge?1?(???f?f2)sin2L?(?f?f2?22?考虑到如下三角函数恒等式
将上式移项，有
??????)sin4L?sin22L?(2sinLcosL)2?4sin2L(1?sin2L)?4sin2L?4sin4L （3.2-11）
1sin4L?sin2L?sin22L （3.2-12）
4再将式（3.2-12）代入式（3.2-10），忽略高阶小量，可得
111??gL?ge?1?(???f?f2)sin2L?(?f?f2)(sin2L?sin22L)?224??1?? （3.2-13） ?ge?1??sin2L?(2?f?f2)sin22L?8???ge(1??sin2L??1sin22L)其中
?1?(2?f?f2) （3.2-14）
若将克莱罗公式（3.2-9）代入上式，还可得 1?5?1?1??2(m?f)f?f2??(5mf?f2) （3.2-15）
8?2?8式（3.2-13）为索密里安公式的良好近似（最大误差约为12μg），实用中的正常重力模型常常是以该形式给出的，它比式（3.2-4）的计算量稍小些。历史上曾给出过以下一些重要的正常重力模型。
（1）1901年，德国人赫尔默（Helmert）根据当时波斯坦系统的几千个重力测量结果，给出正常重力公式为：
18gL?9.7803?(1?0.005302?sin2L?0.000007?sin22L)相应的参考椭球的扁率f?1/298.3。
(m/s2) （3.2-16）
上式称为赫尔默正常重力公式，其中重力扁率??0.005302?1/188.6，利用克雷诺定理，可以计算出
（2）1909年，美国人海福特（Hayford）根据美国当时的大地测量结果给出了一个参考椭球，它的赤道半径Re?6378388m和几何扁率f?1/297.0；1928年，芬兰人海斯卡宁（Heiskanen）根据当
48
2时的重力测量结果计算出正常场地球模型在赤道上的重力为ge?9.78049m/s。若取地球自转角速率
?ie?7.2921151?10?5rad/s，根据上述四个独立参数，可计算得
gL?9.78049?(1?0.0052884?sin2L?0.0000059?sin22L)(m/s2) （3.2-17）
1930年，国际大地测量与地球物理联合会（IUGG）将上式定为国际正常重力公式。
（3）利用现代卫星测量技术，IUGG于1979年通过了1980大地参考系，与其对应的正常重力公式为
gL?9.780327?(1?0.00530244?sin2L?0.00000585?sin22L)(m/s2) （3.2-18）
（4）1987年，WGS-84（World Geodetic System 1984）大地坐标系给出的地球参数为：半长轴
Re?6378137m，扁率f?1/298.257223563，地心引力常数（含大气层）??3.986004418?1014m3/s2，地球自转角速率?ie?7.2921151467?10?5rad/s。如不考虑大气层影响，可推导得正常重力公式
gL?9.780325?(1?0.00530240?sin2L?0.00000582?sin22L)3．重力与海拔高度的关系
(m/s2) （3.2-19）
在地球表面附近的重力场中，引力与离心力相比前者占主要成分，重力随海拔高度增加而减小，其变化规律与引力随距离增加而减小的规律近似相同。分析高度影响时，不妨将地球近似成圆球且质量集中在地心，地球对高度为h的单位质点的引力为
F?对上式两边同时微分，得
?(R?h)2 （3.2-20）
dF??2其中
?(R?h)3dh??2?R3dh???2dh （3.2-21）
?2?2R? （3.2-22） 3若设地心引力常数??3.986005?1014m3/s2和地球平均半径R?6371km，则有?2?3.08?10?6s?2，这表示在地球表面附近高度每升高1m，引力值（或重力值）将减小3.08?10?6m/s2（约0.3μg）。
综合式（3.2-13）和式（3.2-21），可求得地球表面附近正常重力随纬度和高度变化的实用计算公式，即在大地坐标og(L,?,h)处的重力值，记为
gLh?ge(1??sin2L??1sin22L)??2h （3.2-23）
值得注意的是，上式只给出了重力的大小，其方向应该是og点处的铅垂向（真垂线），而不是该点的地理法向（地理垂线），如图3.2-2所示。 1967年，Heiskanen给出了真垂线和地理垂线之间夹角的近似公式
???3hsin2L （3.2-24） gLh其中，?3?8.08?10?9s?2。据上式计算，在纬度L?45°处，高度每上升1km，?约增加
?98.08?10?1000/?9.8。显然，除了赤道和极点外，只有在旋转椭球表面上（0.17\h?0）真垂线与地
理垂线才能严格重合。
49
?zPoLx
图3.2-2 正常重力场的垂线偏差
当以地理坐标系（g系）作为惯性导航参考坐标系时，为了扣除重力的影响，需将重力矢量投影到地理坐标系，表示为
00?????????hsin2L? （3.2-25） gg???gsin?Lh???3????gLhcos??????gLh??上式的三个分量依次表示重力矢量在东、北和天向的投影值，同样在纬度L?45°处，高度每上升1km，北向重力分量将变化8.08?10?9?1000?0.8ug，这一影响在高空高精度惯性导航中是需要考虑和补偿的。更高精度的且与实际地球相关的重力场计算和补偿详见3.3节介绍。
3.3地球重力场的球谐函数模型
在3.2节中正常重力场描述的是规则地球产生的重力场，但实际地球并不规则，正常重力场只是实际重力场的一个较好的近似，为了更加细致地刻画实际地球的重力场，需引入球谐函数和重力位等概念，这在高精度惯性导航系统的重力场建模和补偿中十分有用。
3.3.1 球谐函数的基本概念
1. 拉普拉斯方程
如果三元函数u(x,y,z)在空间区域?上满足如下偏微分方程 ?2u?2u?2u???0 （3.3-1） ?x2?y2?z2则称u是区域?上的调和函数，或称谐函数，上述方程称为拉普拉斯方程（Laplace equation）。
方程（3.3-1）可简记为
?u?0 （3.3-2）
其中，算符
?2?2?2??2?2?2 （3.3-3）
?x?y?z称为拉普拉斯算子。
参见图3.3-1，球面上一点u，其直角坐标(x,y,z)与极坐标(r,?,?)之间的坐标转换关系为
50