一种基于蒙特卡洛采样的GNSS单频单历元姿态确定方法与流程

技术特征：

1.一种基于蒙特卡洛采样的GNSS单频单历元姿态确定方法，其特征在于：其步骤如下：

步骤一：准备工作

首先，给出m+1个接收机对n+1个GNSS卫星的载波观测量给出站际-星际双差观测线性化方程如下所示：

式(1)中，Φ为m条基线向量的载波双差观测量组成的n×m阶矩阵，每一列向量表示一条基线对应的n个双差观测量；B为m条未知基线向量在GNSS参考坐标系下的坐标值组成的3×m阶矩阵；G是从卫星到接收机方向的双差单位向量在GNSS参考坐标系下的表示组成的n×3阶矩阵；λ是载波波长；Z是双差整周模糊度组成的n×m阶矩阵，每一列向量表示一条基线对应的n个双差整周模糊度；V观测噪声组成的n×m阶矩阵；符号vec(·)表示把矩阵的列向量按列序号从小到大的顺序依次从上到下重新排列成一列，组成新的列向量；Q是vec(V)的协方差矩阵；取Q为：

其中σ是载波噪声的标准差；符号表示克罗内克积；

蒙特卡洛取样方法：表示概率分布函数p_x(x)的随机测量，是取样点，是各点对应的权重，并且N_s是取样次数，那么概率分布函数p_x(x)近似表示为

$p_x\left(x\right)\≈\underset{i=1}{\overset{N_s}{\mathrm{\Σ}}}{w}_i\mathrm{\δ}(x-x_i)...\left(3\right)$

其中，δ(·)表示狄拉克δ函数；

步骤二：建立多变量的GNSS姿态模型

定义本地坐标系：原点在主接收机位置，x轴沿着第一条基线方向，y轴在第一条基线和第二条基线决定的平面内并且垂直于x轴，z轴方向由右手定则决定；m条基线向量在本地坐标系中的表示已知，则由m条基线向量组成的3×m坐标矩阵F如下所示：

$\begin{array}{c}m\≥3,F=\left[\begin{array}{ccccc}{f}_{11}& f_{21}& f_{31}& ...& f_{m1}\\ 0& f_{22}& f_{32}& ...& f_{m2}\\ 0& 0& f_{33}& ...& f_{m3}\end{array}\right]\\ m=2,F=\left[\begin{array}{cc}{f}_{11}& f_{21}\\ 0& f_{22}\\ 0& 0\end{array}\right]\end{array}...\left(4\right)$

其中，F中每一列向量表示一条基线向量在本地坐标系下的表示，即：在本地坐标系下，副接收机1的坐标为(f₁₁,0,0)，副接收机2的坐标为(f₂₁,f₂₂,0)，副接收机3的坐标为(f₃₁,f₃₂,f₃₃)，副接收机m的坐标为(f_m1,f_m2,f_m3)；

令R表示从本地坐标系到GNSS参考坐标系的坐标转移矩阵，则

B＝RF···············(5)

将其代入式(1)即得到多变量的GNSS姿态模型：

其中，3阶标准正交方阵R和双差模糊度矩阵为未知量；将R中的参数用四元数表示，则(6)式写为：

其中，q＝(q₀,q₄)^T，q₀＝(q₁,q₂,q₃)^T，并且满足q^Tq＝1，

$R\left(q\right)=\left[\begin{array}{ccc}{q}_1^2-q_2^2-q_3^2+q_4^2& 2\left(q_1{q}_2+q_3{q}_4\right)& 2\left(q_1{q}_3-q_2{q}_4\right)\\ 2\left(q_1{q}_2-q_3{q}_4\right)& -q_1^2+q_2^2-q_3^2+q_4^2& 2\left(q_2{q}_3+q_1{q}_4\right)\\ 2\left(q_1{q}_3+q_2{q}_4\right)& 2\left(q_2{q}_3-q_1{q}_4\right)& -q_1^2-q_2^2+q_3^2+q_4^2\end{array}\right]...\left(8\right)$

q为未知量；

步骤三：蒙特卡洛取样方法构造模糊度的概率分布函数

根据(7)式，双差模糊度表示为：

$Z=\frac{1}{\mathrm{\λ}}(\mathrm{\Φ}-GR(q)F-V)...\left(9\right)$

观测噪声遵从均值为0协方差为Q的正态分布；vec(V)的概率分布函数为：

如果已知四元数q有上界q^u和下界q^l，根据先验姿态信息，q遵从均匀分布：

其中，

$\underset{q\∈\[q^l,q^u\]}{\∫}cdq=1...\left(12\right)$

根据q^Tq＝1，得

$p_q\left(q\right)=\left\{\begin{array}{cc}{c}^{\′}& q_0\∈\[q_0^l,q_0^u\],q_4=\±\sqrt{1-\left|\right|q_0|{|}^2}\\ 0& q_0\∉\[q_0^l,q_0^u\],q_4=\±\sqrt{1-\left|\right|q_0|{|}^2}\end{array}\right....\left(13\right)$

并且，

$\underset{\left|\right|q_0\left|\right|\≤1}{\underset{q_0\∈\[q_0^l,q_0^u\]}{\∫}}{c}^{\′}{\mathrm{dq}}_0=1...\left(14\right)$

其中，是q^l的前三个分量，是q^u的前三个分量；在缺少先验姿态信息的情况下，取q^l＝[1,1,1]^T，q^u＝[1,1,1]^T；根据概率分布函数p_q(q)和p_vec(V)(v)对q和v进行N_s次取样

$\begin{array}{c}\begin{array}{cc}{q}_i~p_q\left(q\right);& v_i~p_{vec\left(V\right)}\left(v\right)\end{array}\\ i=1,2,...,N_s\end{array}...\left(15\right)$

每次取样点对应的权重为

$w_i=\frac{1}{N_s},i=1,2,...,N_s...\left(16\right)$

根据(9)式得出模糊度的样本点则

$p_{vec\left(Z\right)}\left(z\right)\≈\underset{i=1}{\overset{N_s}{\mathrm{\Σ}}}{w}_i\mathrm{\δ}(z-z_i)...\left(17\right)$

模糊度的期望和协方差为

$\stackrel{\‾}{z}=\underset{i=1}{\overset{N_s}{\Σ}}{w}_i{z}_i...\left(18\right)$

$P=\underset{i=1}{\overset{N_s}{\Σ}}{w}_i(z_i-\stackrel{\‾}{z}){(z_i-\stackrel{\‾}{z})}^T...\left(19\right)$

这里得出的期望和方差将用于步骤四的LAMBDA算法中；

步骤四：LAMBDA算法搜索出模糊度的候选值

由于LAMBDA算法是已有的现成方法，这里不做详细说明；通过LAMBDA算法能搜索出N_c个模糊度的候选值

步骤五：计算模糊度最优整数解并确定姿态

由模糊度的候选值根据(9)式算得各候选值对应的姿态矩阵

${\hat{R}}_k=G^{+}(\mathrm{\Φ}-\mathrm{\λ}{\hat{Z}}_k)F^{+},k=1,2,...,N_c...\left(20\right)$

X⁺表示矩阵X的伪逆；其中四元数的候选值由下式得到

$\begin{array}{c}{q}_0=\frac{1}{4q_4}\left[\begin{array}{c}{R}_{23}-R_{32}\\ R_{31}-R_{13}\\ R_{12}-R_{21}\end{array}\right]\\ q_4=\frac{1}{2}{(trR+1)}^{\frac{1}{2}}\end{array}...\left(21\right)$

其中{R_ij,i,j＝1,2,3}是矩阵R的元素；则模糊度的最终整数解和姿态矩阵为

$\[\hat{q},\hat{Z}\]=\underset{k=1,...,N_c}{min}\left|\right|vec(\mathrm{\Φ}-GR({\hat{q}}_k)F-\mathrm{\λ}{\hat{Z}}_k)\left|\right|...\left(22\right)$

通过以上步骤，用蒙特卡洛方法求得了模糊度的期望和协方差，然后将其用于求解模糊度的LAMBDA方法中进行模糊度的解算，减少了计算模糊度的计算量，摆脱了对伪距测量量的依赖，从而能简便快捷的确定姿态，减弱了对环境条件的要求。

2.根据权利要求1所述的一种基于蒙特卡洛采样的GNSS单频单历元姿态确定方法，其特征在于：在步骤四中所述的“LAMBDA算法搜索出模糊度的候选值”，其搜索的过程如下：

先不考虑模糊度的整数限制，直接用最小二乘法解出模糊度的浮点解，把得出的协方差矩阵构造转换矩阵Z，采用Z整数转换矩阵将模糊度和协方差矩阵进行转换，进行整周模糊度搜索，求得模糊度的整数解及其协方差矩阵，再对求得的模糊度和协方差矩阵进行Z逆变换，得到模糊度的整数解及其协方差矩阵。