附录E. 模型与算法
E.6 单点定位
RTKLIB采用迭代加权最小二乘法(LSE)用于“Single”(单点定位)模式(无论是否包含SBAS校正)。
E.6.1 线性最小二乘
假设给定测量向量,它可以被建模为一个包含未知参数向量和一个随机观测误差向量的线性方程:
最小二乘代价函数定义为观测误差的平方和:
通过使用 (E.6.1) 和 (E.6.2) ,代价函数可以重写为(计算时注意代价值为标量,而标量的转置为自身):
为了最小化代价函数,的梯度应该为零。然后
这就是最小二乘法的“正规方程”(Normal Equation):
为了解正规方程,我们可以通过最小二乘法得到估计的未知参数向量:
如果给定了每个测量的权重,可以使用权重矩阵重写代价函数 (E.6.3) 。
为了最小化代价函数,我们可以通过加权最小二乘法(WLS)以类似最小二乘法原本的方式获得估计的未知参数向量:
加权最小二乘法(WLS)的权重矩阵通常定义为:
其中是第个观测误差的先验标准差。
E.6.2 高斯-牛顿迭代
如果观测值不是以线性模型给出的,测量方程可以由一个通用的非线性向量函数表示为:
其中是参数向量的测量向量函数。该方程可以通过在初始参数向量处使用泰勒(级数)展开:
其中是关于在处的偏导数矩阵:
假设初始参数足够接近真实值,并且忽略泰勒级数的高阶项(第二项及其后续项)。我们可以将 (E.6.9) 近似为:
然后可以得到以下线性方程:
通过将线性加权最小二乘法 (E.6.8) 应用于 (E.6.13) ,我们可以得到非线性加权最小二乘法的正规方程:
因此,我们可以通过以下方式获得估计的未知参数向量:
如果初始参数不够接近真实值,我们可以迭代改进估计参数如下:
如果迭代收敛,我们可以得到最终的估计参数为:
迭代最小二乘法通常被称为高斯-牛顿法。请注意,对于具有较大非线性的病态(ill-conditioned)测量方程,简单的高斯-牛顿方法并不总是能够收敛。在这种情况下,我们应该采用另一种策略来处理这种非线性最小二乘问题。对于非线性最小二乘法,最流行的方法是非线性最小二乘的LM(Levenberg-Marquardt)方法。
E.6.3. 接收机位置和钟差估计
如果定位模式 (Positioning Mode) 设置为 “Single” ,将会通过以下单点定位过程,通过逐历元的方式获得最终解。对于一个历元时间,未知参数向量 定义为:
伪距测量向量 可以表示为:
其中 是伪距测量。如果配置选项"Ionosphere Correction"设置为"Iono-Free LC",则使用附录 E.5(E.5.7) 中定义的无电离层组合(线性组合)伪距。其他情况下,则仅使用 伪距。
图E.6-1 单点定位的卫星几何结构
单点定位的测量方程及其偏导数矩阵构成如下:
其中几何距离 和视距(LOS)向量 由 E.3 (3.4) 和 E.3 (3.5) 给出,结合卫星和接收机的位置。卫星位置 和钟差 则根据配置选项“Satellite Ephemeris/Clock”从 E.4 中描述的GNSS卫星星历和时钟推导而来。
为了求解测量方程以获得接收机位置和接收机钟差,RTKLIB采用了迭代加权最小二乘法(LSE),如下所示:
对于迭代加权最小二乘法的初始参数向量 ,在单点定位的第一个历元中可以全部设置为0。一旦获得解,该位置将用于下一个历元接收机位置的迭代初值。对于权重矩阵 ,RTKLIB使用以下公式:
其中:
:卫星系统误差因子
(1: GPS, Galileo, QZSS和Beidou, 1.5: GLONASS, 3.0: SBAS)
:码/载波相位误差比率,伪距定位和载波相位定位使用的是同样的加权公式,这里以载波相位为基础,然后通过系数转换到了伪距上
:载波相位误差因子 和 (m)
:星历和钟差标准差(m)
:电离层校正模型误差标准差(m)
:对流层校正模型误差标准差(m)
:码偏差误差标准差(m)
手册与代码中权计算的差异
RTKLIB代码中是这样的:
// var = fact^2*R^2*(a^2 + (b^2/sin(el) + c^2*(10^(0.1*(snr_max-snr_rover)))) + (d*rcv_std)^2)
注意,这里的和都需要带上平方,手册中并没有很好的说明这一点,另外代码中还包含SNR的部分。
对于星历和钟差的标准差,RTKLIB中使用了URA(用户测距精度)或类似的指标。通过几次迭代,通常情况下解会收敛,并获得估计的接收机位置 和接收机钟差 。
估计的接收机钟差 没有明确输出到结果文件中。它被合并在定位解的时间标签中。这意味着定位解的时间标签表示的并不是接收机的时间标签,而是GPST中测量到的真实信号接收时间(包含钟差修正的观测时间)。
E.6.4 接收机速度和钟漂估计
如果给出了多普勒频率观测值,可以按照以下过程估计接收机速度和钟漂。对于一个历元时间,速度估计的未知参数向量 定义为:
其中 和 分别是接收机在ECEF(地心地固)坐标系中的速度(m/s)和接收机钟漂(s/s)。相应的速率测量向量 可以表示为:
其中 是卫星 的 多普勒频率观测值。
Q1:多普勒观测值的正负
RINEX 中的多普勒观测值 可以表示为:
根据多普勒效应的物理定义:
- 当接收机和卫星靠近时, , ,多普勒记录为正值。
- 当接收机和卫星远离时, , ,多普勒记录为负值。
绝大部分接收机满足上述规律,可以通过查看观测文件中前后两个历元的伪距或载波相位是增加还是减小来判断多普勒的符号定义。
RTKLIB总是使用 频段的多普勒频率观测值。这些测量方程及其偏导数矩阵构成如下:
这些方程中卫星相对于接收机的速度 从以下公式推导:
其中 , 。通过使用与位置估计类似的迭代最小二乘法,我们可以获得接收机的速度和钟漂:
其中权重矩阵 设置为 (非加权最小二乘法)。
E.6.5 解的验证与RAIM FDE
6.3中描述的接收机位置估计可能由于未建模的观测误差而产生无效解。为了检验结果是否有效并防止无效解,RTKLIB在完成位置估计后可以采用以下验证方式。如果验证失败,将发出警告信息并拒绝该解。
其中 是待估参数的数目, 是观测量的数目。 是自由度 和 的卡方分布。GDOP是几何精度因子(dilution of precision)。 ,可以配置选项“Reject Threshold of GDOP”中进行相关设置。
除了上述解的检验之外,RTKLIB在版本2.4.2中还增加了RAIM(接收机自主完好性监控)FDE(故障检测与排除)功能。如果启用了配置选项“RAIM FDE”,那么在式(E.6.33)中的卡方检验失败后,RTKLIB将通过逐个排除可见卫星来进行重新估计。在重新进行了所有尝试后,选择具有最小归一化平方残差 的估计位置作为最终解。在这种方案中,由于卫星故障、接收机故障或较严重的多路径所导致的无效观测值将被视为异常值并被剔除掉。请注意,此功能需要2颗冗余的可见卫星,这意味着至少需要6颗可见卫星才能获得最终解。