在已实现spp程序基础上改rtk，附C++源码（事后和实时流）

原文地址：在已实现spp程序基础上改rtk，附C++源码（事后和实时流）

上面这个链接现在可以下载啦～ 换了地址～ 因为报告我已经搞丢了 所以只分享我的事后处理以及实时流的计算程序给大家。

事先说明😚，本程序实现基于老师linux程序修改，提供代码不能直接运行，只是里面的内容可能会帮助读者思考rtk流程原理，附报告，对于想要学习rtk、理解rtk或是只是想简单实现rtk的同学都适用。看完报告后，你应该对rtk与spp的关系有更加深刻的理解，也应该更清楚最小二乘和模糊度固定在rtk中的作用。
由于程序比较复杂，不同学习能力的人吸收起来可能差异较大，毕竟这是我们花了一个多月才学会吸收的程序（现在再写的花估计几天就能撸出来了😂），所以不提供额外指导~
首先，应该进行时间同步处理，将基准站和流动站数据的时间统一。
该模块主要是修改了示例Rx3ObsData.hpp/.cpp文件，在其中重写了readRecord函数，修改如下
//read obs data with lastepochtime of base
bool Rx3ObsData::readRecord(std::fstream &strm, CommonTime &lastBaseEpoch)
noexcept(false)
{
readRecord(strm);
if(currEpoch>lastBaseEpoch){
return true;
}else{
return false;
}
}

该函数主要是应用于流动站数据的读取，相比原函数添加了新参数lastBaseEpoch，该参数存储了上一次读取的基准站的时间信息，以此来判断本历元是否需要读取新的基准站数据，即当当前读取的流动站时间大于上一基准站时间时，返回true，小于时返回false，在主函数中新增调用如下
// time synchronize
CommonTime lastEpochOfBase = rxDataBase.currEpoch;
bool getBase = false;
while(!rxStreamRover.eof()){
getBase = rxDataRover.readRecord(rxStreamRover,lastEpochOfBase);
if(getBase){
rxDataBase.readRecord(rxStreamBase);
lastEpochOfBase = rxDataBase.currEpoch;
// cout<<“now we get the new base station epoch:”<<currEpochBase<<endl;
}
}

以流动站数据流为外层分辨率，每一步都读取流动站信息，返回是否读取新的基准站数据，若是，读取；若否，循环继续。

然后，对观测方程进行线性化处理。
观测方程线性化主要是利用基站观测卫星的卫地距与伪距或载波相位距离之差，流动站的该差值与基准站的该差值对应作差即得到了线性化的残差结果，利用TypeID存储在satmap中，供之后最小二乘时计算L矩阵使用，单独作C11卫星的残差结果图如下

图中以卫地距和伪距之差为例，作了C2和C6频率的残差图，由于采用的流动站数据和基准站数据一致（零基线数据），两曲线重合。残差绝对值的最大值达到了14m左右，但大部分历元还是处在5m以内。

再然后，最小二乘计算浮点解。
将线性化结果按照map存储顺序进行组合，得到最小二乘残差矩阵L，计算参数矩阵A，赋权阵P，对于载波相关量赋值1/(0.0030.003)，对于伪距相关量赋值1/(0.30.3)，由最小二乘公式计算得到dx矩阵，这里得到的dx矩阵由xyz方向的基线向量改正数、接收机钟差、模糊度组成，维度为5+2*satnum，步骤具体公式见后文算法公式推导流程。
计算得到的基线向量改正数结果可视化如下

为进行对比，作原spp结果与参考坐标残差图如下

由上图可以看到，spp结果精度大部分时刻都在10m以内，一部分时刻达到了20~30m的量级，基线改正数与spp结果差在数值上基本对应，且符号相反，由此可以知道基线改正已经在很大程度上修正了spp的误差。
作浮点解精度结果图如下

可以看到零基线的结果已经在很高的精度了（1cm以内）。
至此，rtk浮点解计算完成，以上结果都是由零基线数据计算得到的，所以精度在较高的水平，而对于短基线的结果，精度则在0.5m以内，为达到更高的精度，接下来需要对计算得到的模糊度进行固定。

然后，进行模糊度固定（短基线情况）。
模糊度固定直接利用老师提供的Lambda方法，floatAmb矩阵使用的是我们dx中第5位之后的单差模糊度浮点解与基准线单差模糊度浮点解之差，而floatAmbCov即双差模糊度的协方差矩阵，该矩阵我们通过单差模糊度协方差矩阵，利用协方差传播率得到双差模糊度协方差矩阵，传播率中的系数矩阵即M阵为双差系数矩阵，基准星项为-1，对角线元素为1，当对角线对应卫星即为基准星时，设置对角线元素为0，通过协方差传播率计算得到的双差模糊度协方差矩阵，需要进一步处理，对基准星项赋值一个较小值，我设置的是0.0001（参照老师提供的模糊度解算例子程序）。
模糊度解算函数如下
ARLambda AR;
VectorXd floatAmb;
floatAmb=lsqRtk.getDifBlur();
MatrixXd floatAmbCov=lsqRtk.getQ();
VectorXd intAmb(floatAmb.size());
intAmb.setZero();
intAmb = AR.resolve(floatAmb, floatAmbCov);
bool isFixed = AR.isFixed();

以起始历元数据为例，浮点解双差模糊度与固定解双差模糊度部分结果截图如下

可以看到，模糊度固定前后差异基本上在1以内。
我对模糊度解算中得到的ratio值可视化如下。

  由图可知，ratio值一般在50~250的区间内波动，而ratio值反映的是模糊度固定的准确程度，在一般应用中，都设置当ratio值大于3时，模糊度解是可信的，并且在短基线解算中，ratio值一般在50以上，由此可知，我们的模糊度固定结果还是在较高的可信水平。

最后，利用固定好的模糊度得到固定解。

  由以上公式计算得到的dx_值前三项为固定后的基线改正，后两项为固定后的接收机钟差，利用dx_前三项对基线改正后，利用改正后的基线加上基准站位置，就得到了我们的RTK解算结果——固定解。绘制固定解相关曲线如下。

上图即固定解结果与参考坐标作差得到的残差示意图，由图可以看到，固定解的精度一般在1cm以内，在某些特殊时刻，固定解的精度达到了2cm左右，整体来说，我的精度也达到了rtk精度的一般要求。

后面就是实时流相关的报告展述了。
报告中有公式推导和程序流程等，就不在这里罗列了，放付费了😭。