GPS软件接收机基础.pdf

《GPS软件接收机基础.pdf》由会员分享，可在线阅读，更多相关《GPS软件接收机基础.pdf（82页珍藏版）》请在淘文阁 - 分享文档赚钱的网站上搜索。

1、 4 第第 2 章章 坐标系统和时间系统坐标系统和时间系统 GPS 卫星环绕地球运行， 卫星的观测量是卫星位置/速度和接收机位置/速度函数，为了能够描述观测量，有必要定义适宜的坐标系统和时间系统。 2.1 参考坐标系统参考坐标系统 2.1.1 地心惯性坐标系统地心惯性坐标系统(ECI) 为了测量和确定 GPS 卫星的轨道，使用地心惯性坐标系统是很方便的6。在 ECI 系统中，原点位于地球的质心，GPS 卫星遵循牛顿的运动和引力定律。在一个典型的 ECI 坐标系中，xy 平面选为与地球的赤道面重合，x 轴相对于天球永久固定地指向某一方向，z 轴垂直于 xy 平面指向北极点，y 轴的选择使得坐标系

2、构成一右手坐标系。地面站确定与预测 GPS 卫星的轨道就是在 ECI 坐标系中进行的。 地球的形状是一个扁状的椭球，由于日月对于地球赤道的突出部分的引力作用，使得地球的赤道面相对于天球运动。因为 x 轴的定义是相对于天球的，而 z 轴的定义是相对于赤道面的，这样地球运动的不规则性将使得如上定义的ECI 系统并非真正的惯性系统。此问题的解决办法为：定义某一瞬间的坐标轴指向。GPS 的 ECI 坐标系统以 UTC（USNO）时间 2000 年 1 月 1 日凌晨的赤道面方向作为其基准，x 轴选为由地球质心指向春分点方向，y 轴和 z 轴的定义和上面一样。由于坐标轴的指向是固定的，因此以此种方式定义

3、的 ECI 坐标系统对于 GPS 来说可以认为是惯性的。 2.1.2 地心地固坐标系统地心地固坐标系统(ECEF) GPS 接收机的位置计算常在一个随地球旋转的坐标系统中进行，此种坐标系统称为地心地固坐标系统。在这样的坐标系统中，可方便地计算接收机的纬度，经度和高度。如同 ECI 坐标系统一样，GPS 使用的 ECEF 坐标系统的 xy平面也与地球的赤道面重合。不同的是，在 ECEF 坐标系统中，x 轴指向 0 经度，y 轴指向东经 90 度。因而，x 轴和 y 轴将随地球自转而旋转。在 ECEF 系 5 统中，z 轴垂直赤道面指向地理北极点。 在计算 GPS 接收机位置之前， 有必要把根据星

4、历信息计算出来的 ECI 坐标系中的卫星位置和速度转换到 ECEF 中，这种转换可以通过旋转矩阵来实现，具体细节请参考 2.1.4 小节。一般情况下，计算出来的接收机位置都是用 ECEF坐标表示，对于某些应用需要把这种坐标转换成接收机的经度、纬度和高度。为实现这种转换， 需要定义一个描述地球的标准物理模型， 在 GPS 中选用的是世界大地 84 坐标系统（WGS-84） 。下面给出的是该标准下的地球物理模型的参数： 地球长半轴长： 6378.137akm= 地球扁率 ： 1/298.257 223 563f = 各类参数间变换关系为： ()()222221112211bfbafabeffaff

5、aebf= = = 其中，e为地球的偏心率，e为地球的二阶偏心率。 2.1.3 站心坐标系站心坐标系 站心坐标系又称地方坐标系(Local Coordinate System)，属于左手笛卡尔坐标系统。如 图 2-1 所示，(), x y z表示站心坐标系，它的原点位于地方点()1111,P x y z， z轴垂直于该点地球椭球的切平面指向上方， x轴在该切平面内指向正北， y轴在该切平面内指向东。对于该坐标系中的任意一点2P，A表示方位角，Z是天顶距，E是仰角，d是2P点在站心坐标系中的向径d的长度。其中，A从正北方向按顺时针方向开始计起的；Z是 z轴和向径d的夹角；E的取值为：当2P点在切

6、平面以上时为正值，在切平面以下时为负值。 6 y(East)x(North)z(Up)AZd()2, P x y zE()1111,P x y z 图 2-1 站心坐标系统 2.1.4 不同坐标系间坐标变换不同坐标系间坐标变换 在GPS数据处理中，常常需要把一个向量的坐标在不同坐标系统间进行转换。下面首先介绍坐标变换的基础，然后给出几种常见的坐标变换公式。 2.1.4.1 坐标变换与旋转矩阵坐标变换与旋转矩阵 如果两个笛卡尔坐标系统共原点且都是左手或都是右手坐标系，那么二者之中的任何一个坐标系都可通过三次连续的旋转转换到另一坐标系统中。这三个旋转矩阵为： Equation Chapter (N

7、ext) Section 1Equation Section (Next) ( )( )( )123100R0cossin,0sincoscos0sinR010,sin0coscossin0Rsincos0001= (2.1) 其中，代表旋转角， 从旋转轴的正方向看如果逆时针旋转则该角取正值。1R， 7 2R和3R分别称作x，y和z轴的旋转矩阵。对于任何旋转矩阵R，有( )( )()1TRRR=；也就是说，旋转矩阵是正交矩阵，其中1R和TR分别表示矩阵的R的逆和转置。 对于两个具有不同原点和不同单位长度的笛卡尔坐标系统，一般的转换公式表示为 0Rnold=+xxx (2.2) 或者 000Rn

8、oldnoldnoldxxxyyyzzz=+ 其中，是比例因子，R是转换矩阵。nx 和oldx分别表示同一向量新的和旧的坐标；0 x 表示平移向量，它是旧的坐标系的原点在新的坐标系中的坐标。 2.1.4.2 ECEF 与大地坐标坐标系间的坐标变换与大地坐标坐标系间的坐标变换 ECEF坐标与大地坐标间转换有不止一种实现方法，比如文献4、6和7中分别介绍了三种不同的实现，其中文献7中的推导最为详细。为方便计，这里仅给出文献6中的变换公式。如 图2-2所示。 xyz(),uuuUxyzhu 图 2-2 ECEF 与大地坐标系 (),uuuU xyz表示接收机在ECEF坐标系中的位置，表示接收机所处位

9、置的经度，表示接收机所处位置的纬度，h表示接收机所处位置相对于参考地 8 球球面的高度。 已知接收机在ECEF中的坐标， 求接收机在大地坐标系中的坐标按 表2-1进行计算6。 已知接收机在大地坐标系中的坐标参数，和h，按照下式可以计算得接收机在ECEF坐标系中的坐标参数。 ()()()2222222coscoscos11tansinsincos11tan1sinsin1sinuuuahexayhezaehe+=+u (2.3) 其中，, a e分别表示地球的长半轴长和偏心率。 表 2-1 根据 ECEF 参数计算大地坐标系下的经度，纬度和高度 (1) 22uurxy=+ (2) 222Eab=

10、 (3) 2254uFb z= (4) ()222221uGreze E=+ (5) 423e FrcG= (6) 3212sccc=+ (7) 22131FPsGs=+ (8) 412Qe P=+ (9) ()()22220111111212u2PezPe rrarQQQQ= +P 9 (10) ()2220uUre rz=+ (11) ()()222201uVre rez=+ (12) 20ub zzaV= (13) 21bhUaV= (14) 20arctanuzezr+= (15) arctan, 0180arctan, 00180arctan, 00uuuouuuouuuuyxxyx

11、xyxyx=+= (3.5) 其中，对于C/A码，61 023 10cf.Hz=。 对于接收时刻ut，由接收机的本地时钟提供，回顾式（2-8） ，伪距表示为 ()()usudratmscttrctttt=+ + (3.6) 假设通道延迟drt在形成伪距观测量时，我们已经提前消除掉了，则上式变为 17 ()()()ususiontropcttrcttctt= =+ + + (3.7) 这里定义的伪距还有一些其它误差因素没有考虑，这些误差将在第4章予以说明。 3.2 载波相位及其测量方法载波相位及其测量方法 3.2.1 载波相位的定义载波相位的定义 载波相位，又称载波差拍相位，是接收机接收的带有多

12、普勒的卫星载波信号与接收机产生的恒定频率信号的差拍信号的相位。这种观测量作为相关通道的副产品而得到或产生自平方通道。平方通道用接收到信号去乘以自己从而得到载波的二次谐波，并且平方后的信号不再含有码的调制11。 由于载波的波长比P码和C/A码的波长都短，因而载波相位的测量量的精度要比用伪码测到的伪距的精度高。对于GPS的L1载波，波长大约是20厘米。根据经验法则，相位测量可精确到1%的波长，这意味着可获得2毫米的精度。 载波相位测量的主要缺点与整周模糊度有关。获得初始的卫星与接收机间的整周模糊度是比较困难的。 一个高质量的GPS接收机大多数时候， 在卫星与接收机间相对位置变化时仍能维护整周计数。

13、然而，由于各种原因，像信号噪声、天线阻挡，都可以导致整周跳变。在多数时候，艰苦的后处理允许检测周跳和对周跳进行修正。可是，整周跳变限制了载波相位测量量的实时应用。 根据定义，载波相位的表达式如下 ( )( )jjuuustt= (3.8) 其中，( )jst表示卫星在卫星时间st发射的载波信号的相位，( )uut表示接收机在接收机时间ut的载波信号的相位，二者的单位都是周(cycles)。在有些的文献中，载波相位定义成( )jst减去( )uut，不过这并不重要。 高稳振荡器在短时间间隔内，相位和频率满足如下关系 ()( )0tttft+=+ (3.9) 令usttt=，有 ( )( )()0

14、juususttftt=+ (3.10) 根据上式我们可以得到 18 ( )( )()0juususttftt= (3.11) 根据伪距测量一节的说明，我们知道接收时刻和发射时刻满足如下关系 ssiontropuurttttttc+= (3.12) 需要注意是，载波相位测量中的电离层对于相位传播的影响与对伪码传播的影响符号相反，相对于光在真空中的传播速度一个超前一个滞后，关于这一点的详细推导可以参考文献4、5和6等。 根据上式，有 ususiontroprttttttc=+ (3.13) 这样，带有误差的载波相位的数学模型为 ()()000usiontropfrfttfttc= + + (3.

15、14) 在实际应用中，载波相位在某一历元的测量基于接收机再生载波与卫星信号载波的对齐，而并不知道哪一周期表示理想的周期同步12。因此，总的相位total包含测量得到的小数相位分量( )Fr和从初始锁定历元0t到时刻t的相位的整数周期计数( )Int，以及在初始历元0t未知的整数周期N ( )()( )00totalFrInt;t ,tN t=+ (3.15) 未知的周期计数N通常被称作整周模糊度(integer cycle ambiguity)。只要接收机在观测阶段保持对载波的连续跟踪，那么在接收机与每颗卫星间只有一个整周模糊度。然而，如果有失锁，将会有周跳引入。 真 正 为 接 收 机 观

16、测 的 是( )()0measuredFrInt;t ,t=+。 因 此 ，( )0totalmeasuredN t=+。在某一瞬间，某颗卫星与接收机间的观测方程可以写作 ()()( )0000totalusiontropmeasuredfrfttfttN tc= + +=+ 进而 ()()( )0000measuredusiontropfrfttfttN tc= + + (3.16) 等式两边乘以载波波长0Lc f=，并定义 Lmeasured = (3.17) 则可得到以米为单位的相应的载波相位方程 ()()usiontropLrcttcttN =+ + + (3.18) 19 与式(3.

17、16)相比，方程增加了一个整周模糊度，而且电离层的误差项符号相反。 3.2.2 载波相位的测量方法载波相位的测量方法 这一小节在上面的基础上给出一种载波相位的测量方法。假设接收到的第j颗卫星信号表示为 ( )()12jRjstA(t )cosf t= (3.19) 其中，RjTjdjfff=+表示接收到的信号频率， 其中Tf是发射频率， 对于L1载波，不考虑钟差，该值为0f =1575.42MHz，对于L2载波，该值为1227 600f =.MHz， df是多普勒频率，j表示第j颗卫星，1A(t )是信号幅度、伪码和数据位合成的信号的简写，而且这里只代表信号的某一支路。 假设接收机产生的本地载

18、波信号为 ( )()202rgstA cosf t=+ (3.20) 其中0是初始相位，gf表示本地载波的频率。 若0Tgfff，测量的载波相位根据前面的定义表示 ()()00000tttmeasured00djdjtttfdtffdtf dt=+= + (3.21) 其中，0对应0，单位为cycle。 根据多普勒的定义，参见3.3小节。measured实际上表示的是接收机与卫星j间的距离自0t开始的变化量，正值表示距离在增加，负值表示距离在减小。从这里可以知道，接收机只要对载波的多普勒连续积分就得到了载波相位观测量。 参照图3-3，如果输入的信号是零中频信号，即( )rst与( )jst混频

19、后的差频信号。实际上，具体实现时接收机并不需要真的产生一个频率为0f的信号，而是可以通过多级混频和采样等方式实现。 下面给出载波相位累积器的实现原理框图。 20 CARRIER_NCO(32bit)dF加法器INTE_CYCLE_CNT(32bit)CLOCK频率字更新输入 图 3-5 载波相位累积器 寄存器dF存放多普勒频率字，用载波跟踪环的处理结果来更新。载波相位累积器包含两部分，CARRIER_NCO作为载波NCO，用以产生本地的再生载波来剥离仅有多普勒频率的载波，同时表示小数周；INTE_CYCLE_CNT表示整周期数。之所以两个计数器合在一起与载波频率字累加，是因为载波频率字表示多普

20、勒频率，而多普勒频率是有正有负的数量，所以这里数据都用补码表示，只要位数足够长，该整周计数器不存在溢出问题，因为它的最大变化范围就是卫星信号传播时间，取70毫秒的话，对应整周计数器（L1载波）的比特位数为 ()( )101070 15754202loglog26.727=位。 所以测量的载波相位可以表示 322measuredCARRIER_ NCOINTE _CYCLE _CNT ,cycle=+单位： (3.22) 需要说明，上式中计数器的内容都按二进制补码解释。 3.3 伪距变化量及其测量方法伪距变化量及其测量方法 3.3.1 伪距变化量的定义伪距变化量的定义 伪距变化量是指两个相邻时刻

21、伪距的差值，即 ( )()1kktt= (3.23) 21 3.3.2 伪距变化量的测量方法伪距变化量的测量方法 根据载波相位的定义（参见3.3小节） ，相邻时刻的载波相位的差对应的就是伪距变化量。当然，如果仅为得到伪距变换量，则没有必要对载波相位连续积分，只需在观测间隔内积分即可。 3.4 多普勒频率及其测量方法多普勒频率及其测量方法 3.4.1 多普勒频率的定义多普勒频率的定义 在接收机天线处，接收的频率，根据经典的多普勒方程表示为 ()1rRTffc=v ai (3.24) 其中，Tf是发射的卫星信号频率，rv是卫星相对于接收机的速度向量，a是由接收机指向卫星的方向的单位向量，c是光速。

22、点积rv ai表示卫星沿a方向速度，相对速度由下式给出 rs=vxu?- - (3.25) 其中，sx?表示卫星的速度，u?表示接收机的速度，二者都是在相同的ECEF坐标系中。多普勒频率为 ()sloslosdRTTTLvvfffffcc= = = xu a? i (3.26) 其中， ()lossv=xu a? i=，称为伪距变化率或者视线速度，正值表示卫星远离接收机， 负值表示卫星靠近接收机， 在这种定义下其符号与多普勒频率相反；LTc f=为载波的波长。 因为losv=?，所以伪距变化量 ( )()1111kkkkkktttlosLdLmeasuredkmeasuredktttdtv d

23、tf dttt= ? (3.27) 22 3.4.2 多普勒频率的测量方法多普勒频率的测量方法 多普勒频率观测量是用来计算接收机速度的观测量13。在测速精度要求不高的情况下，多普勒频率可以直接从多普勒频率字得出（参见 图3-5） ，在这种方法中多普勒频率的估计误差，主要受接收机动态和接收机噪声影响。 第二种方法，就是通过对载波相位微分运算获得。对于实时动态应用，理想的微分器应该具有宽带的频率响应以覆盖所有的动态，另外它还要有尽可能短的群时延。已经提出的载波相位的微分器设计方法，可以分为如下几类： （1）曲线拟合 （2）卡尔曼滤波 （3）泰勒级数近似 （4）使用傅立叶级数和加窗技术的FIR滤波器

24、 （5）使用Remez交换算法的FIR最优滤波器 在设计微分器时需要考虑的问题主要有三个。首先，由于系统的动态性，多普勒的动态范围很宽，因此需要用一个宽带或全带微分器。第二个问题与信号的相关性有关。当载波相位的采样率低于1Hz时，载波相位信号可以看作白噪声，但是当采样率变高，信号间的时间相关性必须予以考虑。最后，对于面向实时的应用，微分运算可能由于缺乏未来信号的信息而受到影响。因此，微分器的设计标准可总结如下： （1）在低频带幅频特性要足够精确，依据系统的动态在宽带范围尽可能地接近理想的微分器； （2）相位响应是线性或接近线性的； （3）滤波器系数的平方和能最小化（保证对噪声的放大效应最小，证

25、明参见文献13） ； （4）容易实时实现，即要求滤波器具有因果性和低阶次。 根据上述标准，在文献13中提出了一种一阶IIR滤波器，即Al-Alaoui一阶微分器。设计IIR滤波器不能像FIR滤波器设计那样能从已知的频率响应计算冲激响应。许多IIR滤波器都能够从相应的模拟滤波器设计导出，然后再把模拟滤波器转换到采样的z-平面。A1-Alaoui设计数字微分器的新颖方法在于通过使用积分器来设计模拟微分器。其依据是：在模拟信号处理中，微分器经 23 常通过颠倒模拟积分器的传输函数得到。 设计A1-Alaoui数字微分器的一般步骤如下 （1）设计一个积分器，使其具有与想要设计的微分器相同的范围和精度；

26、 （2）颠倒积分器的传输函数； （3）为了使传输函数稳定，把落在单位圆外的极点反射到圆内； （4）使用落在单位圆外的极点的倒数补偿幅度。 Al-Alaoui基于非最小相位数字积分器设计了一个一阶IIR微分器，它的有效范围是0.78倍的Nyquist频率。 这个积分器由矩形积分器和梯形积分器合成。假设加权因子分别为34和14，理想的积分器近似为 ( )( )( )()()000314413174142181IRTHzHzHzTzTTzzzz=+=+= (3.28) 其中，0T是载波相位的采样间隔。 反射零点7z = ，使用它的倒数1 7，并通过乘以7r=来补偿幅度，得到一个最小相位数字积分器，它

27、的传输函数为 ( )017781IzTHzz+= (3.29) 颠倒上面的传输函数，得到稳定的Al-Alaoui IIR一阶微分器 ( )8171 7D0zHzTz/=+ (3.30) 这个微分器的频率特性如 图3-6所示。它在0至0.78的全部频带内接近理想的微分器，并且具有很好的线性相位响应。由于这个滤波器的延迟只有半个采样周期，因此它适合实时应用。 根据上述传输函数，可以得到提取多普勒频率微分器为 ( )( )()()01811177dmeasuredmeasureddfnnnfnT= (3.31) 采样间隔0T如何选呢？0T的选择应满足下面的约束 24 0012maxmaxmaxLaT

28、vvT (3.32) 其中，maxa表示接收机与卫星相对运动的最大的加速度；maxv表示允许的速度误差；maxv表示接收机与卫星相对运动的最大速度，因此maxLv表示最大的多普勒频率，该式反映的是采样定理的要求，即采样率不小于信号中的最高频率的2倍。 图 3-6 Al-Alaoui 一阶 IIR 滤波器的频率响应 25 第第 4 章章 误差源误差源 码相位和载波相位测量过程中的误差依据与卫星、信号传播和接收机的关系可以分为三类，即与卫星有关的误差、与信号传播有关的误差和与接收机有关的误差。在参考文献9中有一章专门讨论误差源。在这里，仅简单说明主要误差源的特性、大小和消除方法。 4.1 与卫星有

29、关的误差与卫星有关的误差 4.1.1 卫星时钟误差卫星时钟误差 GPS卫星配有原子钟用以控制星上的所有定时操作， 像广播信号的产生等。虽然这些钟是高度稳定的，但是卫星钟差st仍可与GPS系统时间偏差接近1毫秒。卫星钟差引入的测距误差在3.0米(1)量级。主控站确定出钟差修正参数后上传给卫星，卫星在通过导航电文转发给用户6。接收机通过如下多项式使用这些参数来估计卫星的钟差Equation Section (Next) ()()2012sffocfocrtaattattt=+ (4.1) 其中， ()()()()()()0ff1f2ocraa/attt=2时钟偏差 sec时钟漂移率 sec sec

30、频率漂移率 sec/sec钟差数据的参考时间 sec当前时间历元 sec相对论效应修正 sec 由于上述参数是真正的卫星钟差的拟合估计，所有仍有误差残留。 4.1.2 星历预测误差星历预测误差 卫星的星历估计由地面控制段计算，然后上传给卫星，由卫星广播给用户。就像卫星钟差参数一样，这些估计的星历估计也还有残留的误差，如果把这些误差的影响投影到接收机与卫星的视线方向，这个投影的误差被称作有效伪距 26 误差，它的大小在4.2米(1)量级。 在GPS定位中，根据不同的要求，处理卫星轨道误差的方法主要三种5。 （1）忽略轨道误差。这一方法，主要用于单点定位，这时卫星轨道实际存在的误差，将成为影响定位

31、精度的主要因素之一。 （2）采用轨道改进法处理观测数据。这一方法的基本思想是，在数据处理中，引入表征轨道偏差的改正参数，并假设在短时间内这些参数为常量，将其作为待估量与其他未知参数一起求解。这种方法主要用于大地测量等精密定位，由于未知数多，计算量大，因此不适合实时性要求较高的导航应用。 （3）同步观测值求差。这一方法，是利用在两个或多个观测点，对同一卫星的同步观测值求差，以减弱卫星轨道误差的影响。这种方法，适用于精密相对定位，也不适合单点定位。 4.1.3 选择可用性选择可用性(SA) 对于SPS的用户来说，最大的误差源是SA。SA是美国国防部为了降低用户的导航解精度有意引入的误差源。SA于1

32、990年3月25日正式开始实施。这种人为干扰，通过()epsilon和()delta两种技术实现。技术干扰卫星星历数据，通过降低GPS广播的轨道参数的精度来降低利用C/A码进行实时单点定位的精度；技术是通过抖动卫星时钟来实现的。 不过当前的用户可以不必在考虑如何降低这种影响了，因为SA已于2000年在美国前总统克林顿的授权下停止了干扰。 4.1.4 相对论效应相对论效应 爱因斯坦的广义和特殊相对论理论都是伪距测量过程要考虑的因素。卫星时钟同时受到广义和特殊相对论的影响。为了补偿这两种效应的影响，卫星时钟在发射前已经调整到10.22999999545MHz。这样，位于海平面的用户观测到的频率将会

33、是10.23MHz；因此，用户不需修正这一效应。用户确实需要对另一种相对论的周期性效应予以修正，这种效应的产生是由卫星轨道的微小的离心性导致的。 这种周期性的效应一半由卫星相对于ECI坐标系的速度的周期性变化导致，一半由卫星在它引力场中势能的周期性变化导致。 27 当卫星位于近地点，卫星的速度加快势能变小，二者都导致卫星的时钟变慢。 当卫星位于远地点， 卫星速度减慢势能变大， 二者都导致卫星的时钟变快。这种效应可以由下面的修正项予以补偿 rktFe a sinE= (4.2) 其中， 104 442807633 10F.sec/ mea= =卫星的轨道偏心率卫星轨道的主半轴长卫星轨道的偏近点角

34、kE 4.2 与信号传播有关的误差与信号传播有关的误差 与卫星信号传播有关的误差，主要包括大气层折射误差、多径效应和地球自转效应。 4.2.1 电离层折射的影响电离层折射的影响 GPS卫星信号和其他电磁波信号一样，当其通过电离层时，将受到这一介质弥散特性的影响，使信号的传播路径发生变化。假设，由于电离层引起电磁波信号传播路径的变化为()iono,kSkgp,=或者 分别表示群速度和相速度，则电离层对伪码相位和载波相位的影响的分别为（一阶近似） 2240 340 3iono,piono,g.TECSf.TECSf (4.3) 其中，UsereSVTECn dl=表示信号传播路径上的电子总量，en

35、是电子密度。可见电离层对信号的影响取决于电子总量TEC和信号的频率。 对于GPS卫星信号来说，在夜间当卫星处于天顶方向时，电离层折射对信号传播路径的影响最小，将小于5米；而在日间正午前后，当卫星接近地平线时，对其影响大于150米，为了减弱其影响，在GPS定位中通常采用的措施包括： 28 （1）利用双频观测 由于电离层的影响是信号频率的函数，所以，利用不同频率的电磁信号进行观测，便可以确定其大小，从而对观测量予以修正。 （2）利用电离层模型加以修正 对于单频的GPS接收机用户来说，为了减弱电离层的影响，一般采用由导航电文提供的电离层模型，即Klobuchar模型。Klobuchar模型在中纬度地

36、区，在平均意义上能消除掉50%的电离层延迟。 （3）利用同步观测求差 这一方法，是利用两台或多台接收机，对同一组卫星的同步观测值求差，以减弱电离层折射的影响。尤其当观测点距离较近时，由于卫星到达不同接收机的路径相近，所以，大气层对它们影响接近，所以不同接收机对相同卫星的同步观测值求差，可显著地减弱电离层折射的影响。 上面提到的三种方法，对于单频单点定位的接收机来说，只有第二种方法可取。 因此这里给出利用GPS导航电文里的电离层模型参数修正伪距观测量的算法14。 表 3-1 Klobuchar 电离层模型修正算法 电离层修正模型为 ()()()24995 0 1011 572245 0 101

37、57ionoxxF.AMP, x.TsecF. , x.+= 其中， ionoT是对于L1频率的；如果用户使用的是L2频率，修正项必须乘以()()()222121575 42 1227 670 66LLff.=， ()30000nnmn,AMPAMPsecif AMP,AMP = ()()2504800tx radiansPER= 29 30; 72,000 72,000, 72,000nnmnPERPERif PERPER = 21 0 16 0 0 53F.E=+ n和n都是卫星发送电离层模型参数，n=0,1,2,和3。 其他需要计算的方程为： ()()()()()()40 0641 61

38、70 01370 0220 114 32 10miiiuii.cos. semi-circlessin A semi-circlescos. semi-circlesE.t.GPS time sec=+=+=+=+ 其中 要求086400t ，因此，如果86400t,减去86400秒；如果0t，则604800kktt=； 如果302400kt ，则604800kktt=+。 计算平近点角： ( )0kkMMn t=+ (5.40) 根据平近点角和偏近点角的关系式 kkkMEesinE= (5.41) 采用迭代的办法计算偏近点角。 根据 42 2111kkkkkke sinEsinvecosEc

39、osEecosvecosE= (5.42) 计算真近点角kv。 （2）计算升交角距和轨道摄动改正项 计算升交角距 kkv=+ (5.43) 计算升交角距、轨道半径和轨道倾角的摄动修正项 ()()()()()()222222kuskuckkrskrckkiskickC sinC cosrC sinC cosiC sinC cos=+=+=+ (5.44) （3）计算观测时刻的升交角距、卫星的轨道半径及轨道倾角 ()01kkkkkkkkkuraecosEriii ti=+=+=+ +修正的升交角距 修正的轨道半径 修正的轨道倾角 ? (5.45) （4）计算卫星在轨道直角坐标系中的坐标 0pkkp

40、kkpxr cosuyr sinuz= (5.46) 其中x轴指向升交点，z轴垂直于轨道面，y轴与x,z 形成右手坐标系。 （5）计算升交点的经度 ()00kekeett = + ? (5.47) 其中，e?是指地球的自转角速度，约为-57.2921151467 10 rad / s。 （6）计算卫星在 ECEF 坐标系中的位置 卫星在ECEF坐标系中的位置，只需两个坐标旋转即可得到，图 5-1 中的()TTTX ,Y ,Z这里指的就是ECEF坐标系， 回顾2.1.4小节， 先绕x轴顺时针旋转ki ，然后在绕z轴瞬时针旋转k，从而 ()()310sppkpkkskkppkpkkspkxxx c

41、osy cosi sinyRRiyx siny cosi coszy sini = + (5.48) 43 5.2.2 GPS 卫星速度的计算卫星速度的计算 从多普勒观测方程中知道，为了计算接收机的速度需要知道卫星的速度，然而 GPS 的接口控制文件 ICD-GPS-200 中并没有像计算卫星位置那样给出计算卫星速度的算例。因为，速度是位置坐标的对时间的导数，因此可以仿照前面计算卫星位置的方法， 计算卫星的速度。 在下面给出了 GPS 卫星速度的计算步骤18。 由式(5.40)，有 kMn=? (5.49) 根据式(5.41)，等式两边同时求导得 ()1kkkkkkMEe cosEEe cos

42、EE= = ? (5.50) 所以 1kkkMEe cosE= ? (5.51) 根据式(5.42)，有 ()1kkkcosEe cosvcosve+ =+ (5.52) 对上式，两边同时对时间求导，有 ()1kkkkkkkksinEEe cosvcosEe sinvvsinvv+ = ? (5.53) 所以，真近点角的导数为 ()()11kkkkkksinEEe cosvve cosEsinv+ = ? (5.54) 根据式(5.43)，有升交角距的导数 kkv=? (5.55) 根据式(5.44)，有摄动修正项的导数 ()()()()()()222222222kuskuckkkrskrck

43、kkiskickkuC cosC sinrC cosC siniC cosC sin=? (5.56) 根据式 (5.45)，有 44 kkkkkkkkkuura e sinEEriii=+= += +? (5.57) 根据式 (5.46)，有 pkkkkkkkpkpkkkkkkkpkxr cosur sinuur cosuyuyr sinur cosuur sinuxu=+=+? (5.58) 根据式 (5.47)，有 ke = ? (5.59) 根据式 (5.48)，有卫星在ECEF中的速度为 spkpkkpkkkskspkpkkpkkkskspkpkkxx cosy cosi siny

44、sini siniyyx siny cosi cosy sini cosixzy siniy cosii + = + +? (5.60) 5.3 最小二乘导航解算最小二乘导航解算 5.3.1 算法描述算法描述 在 5.1 小节已经给出了三种观测量分别对应的观测方程组，即 H =x+ (5.61) Hdx?=+ (5.62) H=+gx (5.63) 为了减小观测误差的影响， 可以使用比模型中未知参数数目多些的观测量。因此， 当使用多于 4 颗卫星的观测量时， 问题就变成了求解超定线性方程组20。 为了描述更一般化的问题，而不只限于用户位置和速度的估计，假设n维观测量nyR与p维未知参量pxR满

45、足如下形式的线性方程组 H=+yx (5.64) 其中，np维矩阵Hn pR是一个已知的线性连接矩阵， 是一个正态分布的随机误差向量，即(), ,00。除此之外，在线性模型中，还假设 ( )HE=yx (5.65) ( )V= y (5.66) 其中，( )Ei是数学期望函数，( )Vi是方差协方差函数。 45 求最小二乘解即寻找使残差平方和最小的参数，也即 2minHxx-y (5.67) 当应用于导航问题时，最小化问题表达为 minH xx位置参数估计 (5.68) minHxx-?速度参数估计 (5.69) minHx位置变化量估计 x x (5.70) 其中，i代表欧几里德向量范数。

46、如果不考虑观测误差的影响，上述问题的解就是简单的最小二乘问题，向文献6所示范的那样，上述问题的解统一地表示为 ()1TTLSH HH=xy (5.71) 这时估计的参数虽然是无偏的，但并非最佳的。 考虑观测误差的影响， 如果以观测量的方差协方差矩阵(VCM)的逆矩阵1作为权矩阵的话，则这时的加权最小二乘（WLS）估计等于一个最佳线性无偏估计(BLUE)。 加权最小二乘法通过最小化观测量间的加权偏差和来获得最佳的解。假设矩阵1HHT是非奇异（如果是非奇异的，则假设正确）的且有至少和未知参数一样多的独立的观测量，那么未知参数x的 BLUE 为 ()1T1T1WLSH HH=xy (5.72) 类似

47、的，用户位置的增量x（相对于线性化点的修正量）的加权最小二乘解是 ()1T1T1WLSHHH= x (5.73) 用户速度x ?的 BLUE 是 ()1T1T1WLSHHH=xd?- (5.74) 用户位置的变化量x的 BLUE 是 ()1T11WLSHHHT=xg (5.75) 用户速度x ?利用的观测量是在观测时刻的多普勒估计，因此它表示是对观测时刻的瞬时速度的估计，而用户位置的变化量x的估计实际上表示的是观测间隔内的平均速度。 加权矩阵，即观测量的协方差矩阵可以通过一般的假设或模型得到，像信号仰角或信号强度决定的方差模型。详细的说明见 5.3.4 节。 46 在 用 户 位 置 计 算

48、中 ， 未 知 的 用 户 位 置 坐 标 和 时 间 偏 差TT uuuuuc tx y zc t=xu 估计可通过把增量x加到线性化点得到 0=+xxx (5.76) 其中，00000 uuuuxyzc t=x，线性化方程通过在点()0000 uuuuxyzt做 Taylor级数展开得到。用户位置的 BLUE 的最优性要求线性化误差是可以忽略的，为了保证这一点，估计过程以最新的估计结果作为新的估计，利用新的估计参数再对方程进行线性化，一直到x足够小为止。 使用参数估计的方法来计算未知参数并把观测量看作随机变量，而并非使用与未知参数一样多的观测量来求得唯一的解，这样做的好处在于能够利用冗余的

49、观测量来提高解算精度和进行质量控制。假设一个正确的模型，观测残差定义为估计的观测值与对应的测量值的差，这个残差值反映了测量值间的一致性程度。异常值(Outliers)，即粗差(gross errors) 严重威胁最小二乘估计的结果。因此，残差对于监视被估参数的质量是非常有用的。在 5.6 小节将详细说明基于残差实现接收机的可靠性监视的算法。 5.3.2 卫星分布的几何特征卫星分布的几何特征 用户位置估计的质量不仅取决于测距精度，也与用户/卫星的观测几何有关。精度稀释因子（DOP）概念提供了一个简单的评估用户/卫星分布质量的参数。不好的卫星分布可能放大随机误差和偏差，因而产生大的位置误差。 DO

50、P 概念的详尽描述可以参见文献6。基本的 DOP 参数可以从协因子(cofactor)矩阵 ()1QH HT=算得，若Q表示为 Qxxxyxzxtyxyyyzytzxzyzzzttxtytzttqqqqqqqqqqqqqqqq= (5.77) 其中，矩阵H是设计矩阵(design matrix)，即式 (5.16)所定义的线性连接矩阵。经常使用的DOP参数包括GDOP(几何DOP)，PDOP（位置DOP） ，和TDOP(时间DOP)，它们定义为 xxyyzzttGDOPqqqq=+ (5.78) xxyyzzPDOPqqq=+ (5.79) 47 ttTDOPq= (5.80) 相对于线性化点