后续文章:
论文研读2-2：多GNSS双历元纯相位定位-固定模糊度精度增益
论文研读2-3：多GNSS双历元纯相位定位-定位精度分析
仅相位定位中的模糊度解算问题
在卫星导航定位中，载波相位测量是实现高精度定位的基础，但如果仅使用相位测量值将导致定位解算中方程数目不足，从而无法正确解算出模糊度，位置和钟差，但如果使用伪距测量值，又将会引入精度较低的伪码误差。
一种可行的解决方法便是使用双历元的相位观测值，本文基于A study on multi-GNSS phase-only positioning论文，更详细但又不失简约地介绍仅相位定位模型中的模糊度结算问题。
1 定位模型
对于 m m m颗卫星， f f f个频点的观测数据，对于历元 i ( i = 1 , 2 ) i(i=1,2) i(i=1,2),构造双差观测值double-differenced（DD） Δ ϕ i \Delta \phi_i Δϕi
E ( Δ ϕ i ) = [ Λ ⊗ I m − 1 ] a + [ e ⊗ G i ] Δ b i \begin{equation} E(\Delta \phi_i)=[\Lambda \otimes I_{m-1}]a+[e \otimes G_i]\Delta b_i \end{equation} E(Δϕi)=[Λ⊗Im−1]a+[e⊗Gi]Δbi
其中 Δ ϕ i ∈ R f ( m − 1 ) \Delta\phi_{i}\in R^{f(m-1)} Δϕi∈Rf(m−1)
- ⊗ \otimes ⊗为kronecker积，后续的处理中会用到其性质，详细运算规则参考克罗内克积kronecker的定义与常见性质
- a ∈ Z f ( m − 1 ) a\in\mathbb{Z}^{f(m-1)} a∈Zf(m−1)是待求解的整周模糊度， Δ b i ∈ R 3 \Delta b_i\in R^3 Δbi∈R3是三维基线 b i b_i bi在历元 i i i时的未知增量。
- Λ = d i a g ( λ 1 , ⋯ , λ f ) \Lambda=diag(\lambda_1,\cdots,\lambda_f) Λ=diag(λ1,⋯,λf)， λ j \lambda_j λj为频率 j j j的波长。
- G i = D m T A i ∈ R ( m − 1 ) × 3 G_i=D^T_mA_i\in R^{(m-1)\times3} Gi=DmTAi∈R(m−1)×3,双差观测方程中的设计矩阵，将基线参数 Δ b i \Delta b_i Δbi（接收机位置的增量）与双差相位观测值关联起来。双差观测方程线性化具体细节可以参考GNSS差分定位方程及其线性化
- A i ∈ R m × 3 A_i\in R^{m\times 3} Ai∈Rm×3的每一行为接收机-卫星方向向量，描述卫星几何构型，直接影响基线解算的精度。
- D m ∈ R m × ( m − 1 ) D_m\in R^{m\times(m-1)} Dm∈Rm×(m−1)为卫星间差分矩阵，用于将单差（站间单差）观测值转换为双差观测值。在此将第一颗卫星的观测方程作为参考，其余卫星的观测方程都减去第一个卫星的观测方程。因此其具体形式为：
D m = [ − 1 1 0 ⋯ 0 − 1 0 1 ⋯ 0 ⋮ ⋮ ⋮ ⋱ ⋮ − 1 0 0 ⋯ 1 ] D_m = \begin{bmatrix} -1 & 1 & 0 & \cdots & 0 \\ -1 & 0 & 1 & \cdots & 0 \\ \vdots & \vdots & \vdots & \ddots & \vdots \\ -1 & 0 & 0 & \cdots & 1 \end{bmatrix} Dm= −1−1⋮−110⋮001⋮0⋯⋯⋱⋯00⋮1 - I m − 1 ∈ R ( m − 1 ) × ( m − 1 ) I_{m-1} \in R^{(m-1)\times(m-1)} Im−1∈R(m−1)×(m−1)为单位矩阵， e ∈ R f e \in R^f e∈Rf为全 1 1 1向量，通过克罗内克积将 Λ \Lambda Λ和 G i G_i Gi的维度与 E ( Δ ϕ i ) E(\Delta \phi_i) E(Δϕi)一致。
经过分析可以发现方程的个数为 f ( m − 1 ) f(m-1) f(m−1)，但是未知数一共有 f ( m − 1 ) + 3 f(m-1)+3 f(m−1)+3个，包括 f ( m − 1 ) f(m-1) f(m−1)个模糊度， 3 3 3个位置变量，所以需要两个历元的观测数据构造 2 f ( m − 1 ) 2f(m-1) 2f(m−1)个方程，这时一共有 f ( m − 1 ) + 6 f(m-1)+6 f(m−1)+6个未知量。
为了保证可以准确解算出未知量，所以需要保证 2 f ( m − 1 ) > f ( m − 1 ) + 6 2f(m-1)>f(m-1)+6 2f(m−1)>f(m−1)+6，当只有一个频点时，至少需要7颗卫星。
这时可以列出双历元仅相位观测方程:
E ( Δ ϕ 1 Δ ϕ 2 ) = [ Λ ⊗ I m − 1 e ⊗ G 1 0 Λ ⊗ I m − 1 0 e ⊗ G 2 ] [ a Δ b 1 Δ b 2 ] \begin{equation} \mathbf{E}\begin{pmatrix} \Delta \phi_1 \\ \Delta \phi_2 \end{pmatrix} = \begin{bmatrix} \Lambda \otimes I_{m-1} & e \otimes G_1 & 0 \\ \Lambda \otimes I_{m-1} & 0 & e \otimes G_2 \end{bmatrix} \begin{bmatrix}a\\\Delta b_1\\\Delta b_2\end{bmatrix} \end{equation} E(Δϕ1Δϕ2)=[Λ⊗Im−1Λ⊗Im−1e⊗G100e⊗G2] aΔb1Δb2
为简化后续计算令 A = [ Λ ⊗ I m − 1 e ⊗ G 1 0 Λ ⊗ I m − 1 0 e ⊗ G 2 ] A=\begin{bmatrix} \Lambda \otimes I_{m-1} & e \otimes G_1 & 0 \\ \Lambda \otimes I_{m-1} & 0 & e \otimes G_2 \end{bmatrix} A=[Λ⊗Im−1Λ⊗Im−1e⊗G100e⊗G2] , x = [ a Δ b 1 Δ b 2 ] T x=\begin{bmatrix}a&\Delta b_1&\Delta b_2\end{bmatrix}^{T} x=[aΔb1Δb2]T。
2 观测噪声的协方差矩阵
假设单频未差分相位观测噪声为独立同分布的高斯白噪声，其方差为 σ ϕ 2 \sigma_\phi^2 σϕ2，则站间单差后的噪声方差为：
Var ( ϕ i , r − ϕ i , b ) = Var ( ϕ i , r ) + Var ( ϕ i , b ) = 2 σ ϕ 2 \text{Var}(\phi_{i,r} - \phi_{i,b}) = \text{Var}(\phi_{i,r}) + \text{Var}(\phi_{i,b}) = 2\sigma_\phi^2 Var(ϕi,r−ϕi,b)=Var(ϕi,r)+Var(ϕi,b)=2σϕ2
其中 ϕ i , r \phi_{i,r} ϕi,r 和 ϕ i , b \phi_{i,b} ϕi,b 分别为参考站和流动站的相位观测值。
双差（星间差）操作通过矩阵 D m D_m Dm 实现。设单差观测向量为 Δ ϕ i s ∈ R m \Delta \phi^s_i \in \mathbb{R}^m Δϕis∈Rm，则双差观测向量为：
Δ ϕ i = D m T Δ ϕ i s \Delta \phi_i = D_m^T \Delta \phi^s_i Δϕi=DmTΔϕis
其协方差矩阵为：
Q ϕ i ϕ i = D m T ⋅ Q Δ ϕ i s ⋅ D m = 2 σ ϕ 2 D m T D m Q_{\phi_i \phi_i} = D_m^T \cdot Q_{\Delta \phi^s_i} \cdot D_m = 2\sigma_\phi^2 D_m^T D_m Qϕiϕi=DmT⋅QΔϕis⋅Dm=2σϕ2DmTDm
同时我们知道卫星仰角越低，观测噪声越大。所以定义仰角权重矩阵 C = diag ( w 1 − 1 , … , w m − 1 ) C = \text{diag}(w_1^{-1}, \ldots, w_m^{-1}) C=diag(w1−1,…,wm−1)，其中 w s = sin 2 ( θ s ) w_s = \sin^2(\theta_s) ws=sin2(θs)（ θ s \theta_s θs 为卫星仰角）。
双差协方差矩阵需结合仰角权重：
Q ϕ i ϕ i = 2 σ ϕ 2 D m T C D m Q_{\phi_i \phi_i} = 2\sigma_\phi^2 D_m^T C D_m Qϕiϕi=2σϕ2DmTCDm
令 W i − 1 = D m T C D m W_i^{-1} = D_m^T C D_m Wi−1=DmTCDm，则：
Q ϕ i ϕ i = 2 σ ϕ 2 W i − 1 Q_{\phi_i \phi_i} = 2\sigma_\phi^2 W_i^{-1} Qϕiϕi=2σϕ2Wi−1
对于 f f f 个频率的观测值，假设各频率噪声独立，协方差矩阵扩展为克罗内克积形式：
Q ϕ i ϕ i = 2 σ ϕ 2 [ I f ⊗ W i − 1 ] Q_{\phi_i \phi_i} = 2\sigma_\phi^2 \left[ I_f \otimes W_i^{-1} \right] Qϕiϕi=2σϕ2[If⊗Wi−1]
其中 I f I_f If 是 f × f f \times f f×f 单位矩阵，表示不同频率的独立性。
同时假设不同历元（如 i = 1 i=1 i=1 和 i = 2 i=2 i=2）的观测噪声不相关，因此：
Q ϕ 1 ϕ 2 = 0. Q_{\phi_1 \phi_2} = 0. Qϕ1ϕ2=0.
所以载波相位观测值的协方差矩阵为：
Q ϕ ϕ = 2 σ ϕ 2 [ I f ⊗ W 1 − 1 0 0 I f ⊗ W 2 − 1 ] \begin{equation} Q_{\phi \phi} = 2\sigma_\phi^2\begin{bmatrix} I_f \otimes W_1^{-1}& 0 \\ 0&I_f \otimes W_2^{-1} \end{bmatrix} \end{equation} Qϕϕ=2σϕ2[If⊗W1−100If⊗W2−1]
3 模糊度方差与基线方差
使用加权最小二乘方法求解双历元仅相位观测方程:
x = ( A T Q ϕ ϕ − 1 A ) − 1 A T Q ϕ ϕ − 1 Δ ϕ \begin{equation} x=(A^TQ_{\phi\phi}^{-1}A)^{-1}A^TQ_{\phi\phi}^{-1}\Delta \phi \end{equation} x=(ATQϕϕ−1A)−1ATQϕϕ−1Δϕ
因此解算结果的协方差矩阵为：
E ( X X T ) = E { [ ( A T Q ϕ ϕ − 1 A ) − 1 A T Q ϕ ϕ − 1 Δ ϕ ] [ ( A T Q ϕ ϕ − 1 A ) − 1 A T Q ϕ ϕ − 1 Δ ϕ ] T } = ( A T Q ϕ ϕ − 1 A ) − 1 E(XX^T) =E\left\{[(A^TQ_{\phi\phi}^{-1}A)^{-1}A^TQ_{\phi\phi}^{-1}\Delta \phi][(A^TQ_{\phi\phi}^{-1}A)^{-1}A^TQ_{\phi\phi}^{-1}\Delta \phi]^T\right\} =(A^TQ_{\phi\phi}^{-1}A)^{-1} E(XXT)=E{[(ATQϕϕ−1A)−1ATQϕϕ−1Δϕ][(ATQϕϕ−1A)−1ATQϕϕ−1Δϕ]T}=(ATQϕϕ−1A)−1
在求解该表达式过程中会遇到矩阵求逆，矩阵求逆较为复杂，可以采用先求解最小二乘的正规矩阵the least-squares normal matrix N = A T Q ϕ ϕ − 1 A N=A^TQ_{\phi\phi}^{-1}A N=ATQϕϕ−1A ，然后利用分块矩阵求逆公式进行求解。
因此：
N = A T Q ϕ ϕ − 1 A = 1 2 σ ϕ 2 [ Λ ⊗ I m − 1 Λ ⊗ I m − 1 e T ⊗ G 1 T 0 0 e T ⊗ G 2 T ] [ I f ⊗ W 1 0 0 I f ⊗ W 2 ] [ Λ ⊗ I m − 1 e ⊗ G 1 0 Λ ⊗ I m − 1 0 e ⊗ G 2 ] = 1 2 σ ϕ 2 [ Λ ⊗ W 1 Λ ⊗ W 2 e T ⊗ G 1 T W 1 0 0 e T ⊗ G 2 T W 2 ] [ Λ ⊗ I m − 1 e ⊗ G 1 0 Λ ⊗ I m − 1 0 e ⊗ G 2 ] = 1 2 σ ϕ 2 [ Λ 2 ⊗ ( W 1 + W 2 ) Λ e ⊗ W 1 G 1 Λ e ⊗ W 2 G 2 e T Λ ⊗ G 1 T W 1 f G 1 T W 1 G 1 0 e T Λ ⊗ G 2 T W 2 0 f G 2 T W 2 G 2 ] = [ N a a N a b N a b T N b b ] \begin{aligned} N=&A^TQ_{\phi\phi}^{-1}A\\ =&\frac{1}{2\sigma_\phi^2}\begin{bmatrix}\Lambda \otimes I_{m-1} &\Lambda \otimes I_{m-1} \\ e^T \otimes G_1^T & 0\\0 & e^T \otimes G_2^T\end{bmatrix} \begin{bmatrix} I_f \otimes W_1& 0 \\0&I_f \otimes W_2\end{bmatrix} \begin{bmatrix}\Lambda \otimes I_{m-1} & e \otimes G_1 & 0 \\\Lambda \otimes I_{m-1} & 0 & e \otimes G_2\end{bmatrix}\\ =&\frac{1}{2\sigma_\phi^2}\begin{bmatrix}\Lambda \otimes W_1 &\Lambda \otimes W_2 \\ e^T \otimes G_1^TW_1 & 0\\0 & e^T \otimes G_2^TW_2\end{bmatrix} \begin{bmatrix}\Lambda \otimes I_{m-1} & e \otimes G_1 & 0 \\\Lambda \otimes I_{m-1} & 0 & e \otimes G_2\end{bmatrix}\\ =&\frac{1}{2\sigma_\phi^2}\begin{bmatrix}\Lambda^2 \otimes (W_1+W_2) &\Lambda e \otimes W_1G_1 &\Lambda e \otimes W_2G_2 \\ e^T \Lambda \otimes G_1^TW_1 & fG_1^TW_1G_1&0\\e^T \Lambda \otimes G_2^TW_2 & 0& fG_2^TW_2G_2\end{bmatrix}\\ =&\begin{bmatrix} N_{aa}&N_{ab}\\ N_{ab}^T&N_{bb}\end{bmatrix} \end{aligned} N=====ATQϕϕ−1A2σϕ21 Λ⊗Im−1eT⊗G1T0Λ⊗Im−10eT⊗G2T [If⊗W100If⊗W2][Λ⊗Im−1Λ⊗Im−1e⊗G100e⊗G2]2σϕ21 Λ⊗W1eT⊗G1TW10Λ⊗W20eT⊗G2TW2 [Λ⊗Im−1Λ⊗Im−1e⊗G100e⊗G2]2σϕ21 Λ2⊗(W1+W2)eTΛ⊗G1TW1eTΛ⊗G2TW2Λe⊗W1G1fG1TW1G10Λe⊗W2G20fG2TW2G2 [NaaNabTNabNbb]
其中：
- N a a = 1 2 σ ϕ 2 Λ 2 ⊗ ( W 1 + W 2 ) N_{aa} = \frac{1}{2\sigma_\phi^2} \Lambda^2 \otimes (W_1 + W_2) Naa=2σϕ21Λ2⊗(W1+W2),
- N a b = 1 2 σ ϕ 2 ∑ i = 1 2 u i T ⊗ Λ e ⊗ W i G i N_{ab} = \frac{1}{2\sigma_\phi^2} \sum_{i=1}^2 u_i^T \otimes \Lambda e \otimes W_i G_i Nab=2σϕ21∑i=12uiT⊗Λe⊗WiGi,
- N b b = f 2 σ ϕ 2 ∑ i = 1 2 u i u i T ⊗ G i T W i G i N_{bb} = \frac{f}{2\sigma_\phi^2} \sum_{i=1}^2 u_i u_i^T \otimes G_i^T W_i G_i Nbb=2σϕ2f∑i=12uiuiT⊗GiTWiGi.
且 u 1 = [ 1 0 ] T u_1=\begin{bmatrix}1&0\end{bmatrix}^T u1=[10]T和 u 2 = [ 0 1 ] T u_2=\begin{bmatrix}0&1\end{bmatrix}^T u2=[01]T，基于此，结算结果的协方差矩阵为
N − 1 = [ Q a ^ a ^ Q a ^ b ^ Q b ^ a ^ Q b ^ b ^ ] N^{-1} = \begin{bmatrix} Q_{\hat{a}\hat{a}} & Q_{\hat{a}\hat{b}} \\ Q_{\hat{b}\hat{a}} & Q_{\hat{b}\hat{b}} \end{bmatrix} N−1=[Qa^a^Qb^a^Qa^b^Qb^b^]
其中：
- Q a ^ a ^ = [ N a a − N a b N b b − 1 N a b T ] − 1 Q_{\hat{a}\hat{a}}= [N_{aa} - N_{ab} N_{bb}^{-1} N_{ab}^T]^{-1} Qa^a^=[Naa−NabNbb−1NabT]−1为模糊度方差矩阵；
- Q b ^ b ^ = N b b − 1 + N b b − 1 N a b T Q a ^ a ^ N a b N b b − 1 Q_{\hat{b}\hat{b}} = N_{bb}^{-1} + N_{bb}^{-1}N_{ab}^TQ_{\hat{a}\hat{a}}N_{ab}N_{bb}^{-1} Qb^b^=Nbb−1+Nbb−1NabTQa^a^NabNbb−1为模糊度未固定时的基线向量的方差矩阵；
- Q b ˇ b ˇ = N b b − 1 = 2 σ ϕ 2 f ∑ i = 1 2 u i u i T ⊗ ( G i T W i G i ) − 1 Q_{\check{b}\check{b}}=N_{bb}^{-1}=\frac{2\sigma_{\phi}^2}{f}\sum_{i = 1}^{2}\boldsymbol{u}_i\boldsymbol{u}_i^T\otimes(\boldsymbol{G}_i^T\boldsymbol{W}_i\boldsymbol{G}_i)^{-1} Qbˇbˇ=Nbb−1=f2σϕ2∑i=12uiuiT⊗(GiTWiGi)−1为模糊度固定时的基线向量的方差矩阵。
分块矩阵和矩阵的逆矩阵表达式可以参考矩阵菜谱The Matrix Cookbook，网上也有很多相关资料，可以直接参考分块矩阵怎么求逆？- 知乎 - David Sun的回答中的情况（4）。
3.1 Q a ^ a ^ Q_{\hat{a}\hat{a}} Qa^a^的详细推导
首先求解 N a b N b b − 1 N a b T N_{ab} N_{bb}^{-1} N_{ab}^T NabNbb−1NabT。
N a b N b b − 1 N a b T = 1 2 f σ ϕ 2 [ Λ e ⊗ W 1 G 1 Λ e ⊗ W 2 G 2 ] [ ( G 1 T W 1 G 1 ) − 1 0 0 ( G 2 T W 2 G 2 ) − 1 ] [ e T Λ T ⊗ G 1 T W 1 e T Λ T ⊗ G 2 T W 2 ] = 1 2 f σ ϕ 2 [ Λ e ⊗ W 1 G 1 ( G 1 T W 1 G 1 ) − 1 Λ e ⊗ W 2 G 2 ( G 2 T W 2 G 2 ) − 1 ] [ e T Λ ⊗ G 1 T W 1 e T Λ ⊗ G 2 T W 2 ] = 1 2 f σ ϕ 2 Λ e e T Λ ⊗ ∑ i = 1 2 W i G i ( G i T W i G i ) − 1 G i T W i = 1 2 σ ϕ 2 Λ P e Λ ⊗ ∑ i = 1 2 W i P G i \begin{aligned} N_{ab} N_{bb}^{-1} N_{ab}^T =&\frac{1}{2f\sigma_\phi^2}\begin{bmatrix}\Lambda e \otimes W_1G_1 &\Lambda e \otimes W_2G_2 \end{bmatrix}\begin{bmatrix}(G_1^TW_1G_1)^{-1}&0\\ 0& (G_2^TW_2G_2)^{-1}\end{bmatrix}\begin{bmatrix}e^T\Lambda^T \otimes G_1^TW_1 \\e^T\Lambda^T \otimes G_2^TW_2 \end{bmatrix}\\ =&\frac{1}{2f\sigma_\phi^2}\begin{bmatrix}\Lambda e \otimes W_1G_1(G_1^TW_1G_1)^{-1} &\Lambda e \otimes W_2G_2(G_2^TW_2G_2)^{-1} \end{bmatrix} \begin{bmatrix} e^T\Lambda \otimes G_1^TW_1 \\e^T\Lambda \otimes G_2^TW_2 \end{bmatrix}\\ =&\frac{1}{2f\sigma_\phi^2}\Lambda ee^T\Lambda \otimes \sum_{i=1}^{2}W_iG_i(G_i^TW_iG_i)^{-1}G_i^TW_i \\ =&\frac{1}{2\sigma_\phi^2}\Lambda P_e\Lambda \otimes \sum_{i=1}^{2}W_i\mathcal{P}_{G_i} \end{aligned} NabNbb−1NabT====2fσϕ21[Λe⊗W1G1Λe⊗W2G2][(G1TW1G1)−100(G2TW2G2)−1][eTΛT⊗G1TW1eTΛT⊗G2TW2]2fσϕ21[Λe⊗W1G1(G1TW1G1)−1Λe⊗W2G2(G2TW2G2)−1][eTΛ⊗G1TW1eTΛ⊗G2TW2]2fσϕ21ΛeeTΛ⊗i=1∑2WiGi(GiTWiGi)−1GiTWi2σϕ21ΛPeΛ⊗i=1∑2WiPGi
其中：
- P e = ( 1 / f ) e e T \mathcal{P}_e=(1/f)ee^T Pe=(1/f)eeT为投影矩阵，其正交投影矩阵为 P e ⊥ = I − P e \mathcal{P}_e^{\perp}=I-\mathcal{P}_e Pe⊥=I−Pe
- P G i = G i G i + \mathcal{P}_{G_i}=G_iG_i^+ PGi=GiGi+加权正交投影矩阵，同时 G i + = ( G i T W i G i ) − 1 G i T W i G_i^+=(G_i^TW_iG_i)^{-1}G_i^TW_i Gi+=(GiTWiGi)−1GiTWi为 G i G_i Gi的最小二乘逆，另外 P G i ⊥ = I − P G i \mathcal{P}_{G_i}^{\perp}=I-\mathcal{P}_{G_i} PGi⊥=I−PGi。
投影矩阵的构造与性质可参考正交投影、投影矩阵、最小二乘法LS，而加权正交投影满足幂等性，加权内积下的对称性。
所以：
Q a ^ a ^ − 1 = N a a − N a b N b b − 1 N a b T = 1 2 σ ϕ 2 Λ 2 ⊗ ( W 1 + W 2 ) − 1 2 σ ϕ 2 Λ P e Λ ⊗ ∑ i = 1 2 W i P G i = 1 σ ϕ 2 Λ ( I − P e ) Λ ⊗ W + 1 2 σ ϕ 2 Λ P e Λ ⊗ ∑ i = 1 2 W i − 1 2 σ ϕ 2 Λ P e Λ ⊗ ∑ i = 1 2 W i P G i = 1 σ ϕ 2 { Λ P e ⊥ Λ ⊗ W + 1 2 Λ P e Λ ⊗ ∑ i = 1 2 W i P G i ⊥ } \begin{align*} Q_{\hat{a}\hat{a}}^{-1}=& N_{aa} - N_{ab} N_{bb}^{-1} N_{ab}^T\\ =&\frac{1}{2\sigma_\phi^2} \Lambda^2 \otimes (W_1 + W_2)-\frac{1}{2\sigma_\phi^2}\Lambda P_e\Lambda \otimes \sum_{i=1}^{2}W_i\mathcal{P}_{G_i} \\ =&\frac{1}{\sigma_\phi^2} \Lambda(I-\mathcal{P}_e)\Lambda\otimes W+\frac{1}{2\sigma_\phi^2} \Lambda\mathcal{P}_e\Lambda\otimes \sum_{i=1}^{2}W_i-\frac{1}{2\sigma_\phi^2}\Lambda P_e\Lambda \otimes \sum_{i=1}^{2}W_i\mathcal{P}_{G_i}\\ =&\frac{1}{\sigma_\phi^2}\left\{\Lambda\mathcal{P}_e^{\perp}\Lambda\otimes W +\frac{1}{2} \Lambda\mathcal{P}_e\Lambda\otimes \sum_{i=1}^{2}W_i\mathcal{P}_{G_i}^{\perp}\right \} \end{align*} Qa^a^−1====Naa−NabNbb−1NabT2σϕ21Λ2⊗(W1+W2)−2σϕ21ΛPeΛ⊗i=1∑2WiPGiσϕ21Λ(I−Pe)Λ⊗W+2σϕ21ΛPeΛ⊗i=1∑2Wi−2σϕ21ΛPeΛ⊗i=1∑2WiPGiσϕ21{ΛPe⊥Λ⊗W+21ΛPeΛ⊗i=1∑2WiPGi⊥}
其中 W = ( 1 / 2 ) ∑ i = 1 2 W i W = (1/2)\sum_{i = 1}^{2}W_i W=(1/2)∑i=12Wi，同时不难验证上式为对角阵， Q a ^ a ^ Q_{\hat{a}\hat{a}} Qa^a^直接由上式求逆即可。
Q a ^ a ^ = σ ϕ 2 [ Λ − 1 P e ⊥ Λ − 1 ⊗ W − 1 ] + 2 σ ϕ 2 [ Λ − 1 P e Λ − 1 ⊗ Q ρ ^ ρ ^ ] \begin{equation} Q_{\hat{a}\hat{a}}=\sigma_{\phi}^2 [\Lambda^{-1}\mathcal{P}_e^{\perp}\Lambda^{-1}\otimes W^{-1}] + 2\sigma_{\phi}^2[\Lambda^{-1}\mathcal{P}_e\Lambda^{-1}\otimes Q_{\hat{\rho}\hat{\rho}}] \end{equation} Qa^a^=σϕ2[Λ−1Pe⊥Λ−1⊗W−1]+2σϕ2[Λ−1PeΛ−1⊗Qρ^ρ^]
其中 Q ρ ^ ρ ^ = ( ∑ i = 1 2 W i P G i ⊥ ) − 1 Q_{\hat{\rho}\hat{\rho}}=(\sum_{i=1}^{2}W_i\mathcal{P}_{G_i}^{\perp})^{-1} Qρ^ρ^=(∑i=12WiPGi⊥)−1
3.2 Q b ^ b ^ Q_{\hat{b}\hat{b}} Qb^b^的详细推导
根据前边推导，我们已经得知
Q b ˇ b ˇ = N b b − 1 = 2 σ ϕ 2 f ∑ i = 1 2 u i u i T ⊗ ( G i T W i G i ) − 1 Q a ^ a ^ = σ ϕ 2 [ Λ − 1 P e ⊥ Λ − 1 ⊗ W − 1 ] ⏟ 1 ◯ + 2 σ ϕ 2 [ Λ − 1 P e Λ − 1 ⊗ Q ρ ^ ρ ^ ] ⏟ 2 ◯ Q_{\check{b}\check{b}}=N_{bb}^{-1}=\frac{2\sigma_{\phi}^2}{f}\sum_{i = 1}^{2}\boldsymbol{u}_i\boldsymbol{u}_i^T\otimes(G_i^TW_iG_i)^{-1}\\ Q_{\hat{a}\hat{a}}=\underbrace{\sigma_{\phi}^2 [\Lambda^{-1}\mathcal{P}_e^{\perp}\Lambda^{-1}\otimes W^{-1}]}_{\textcircled{1}}+ \underbrace{2\sigma_{\phi}^2[\Lambda^{-1}\mathcal{P}_e\Lambda^{-1}\otimes Q_{\hat{\rho}\hat{\rho}}]}_{\textcircled{2}} Qbˇbˇ=Nbb−1=f2σϕ2i=1∑2uiuiT⊗(GiTWiGi)−1Qa^a^=1◯ σϕ2[Λ−1Pe⊥Λ−1⊗W−1]+2◯ 2σϕ2[Λ−1PeΛ−1⊗Qρ^ρ^]
模糊度未固定时的基线向量的方差矩阵 Q b ^ b ^ = N b b − 1 + N b b − 1 N a b T Q a ^ a ^ N a b N b b − 1 Q_{\hat{b}\hat{b}} = N_{bb}^{-1} + N_{bb}^{-1}N_{ab}^TQ_{\hat{a}\hat{a}}N_{ab}N_{bb}^{-1} Qb^b^=Nbb−1+Nbb−1NabTQa^a^NabNbb−1，
首先求解 N b b − 1 N a b T Q a ^ a ^ N a b N b b − 1 N_{bb}^{-1}N_{ab}^TQ_{\hat{a}\hat{a}}N_{ab}N_{bb}^{-1} Nbb−1NabTQa^a^NabNbb−1
N b b − 1 N a b T Q a ^ a ^ N a b N b b − 1 = 1 f [ ( G 1 T W 1 G 1 ) − 1 0 0 ( G 2 T W 2 G 2 ) − 1 ] [ e T Λ ⊗ G 1 T W 1 e T Λ ⊗ G 2 T W 2 ] Q a ^ a ^ N a b N b b − 1 = 1 f 2 [ ( G 1 T W 1 G 1 ) − 1 ( e T Λ ⊗ G 1 T W 1 ) ( G 2 T W 2 G 2 ) − 1 ( e T Λ ⊗ G 2 T W 2 ) ] ⏟ 6 × f ( m − 1 ) Q a ^ a ^ [ ( Λ e ⊗ W 1 G 1 ) ( G 1 T W 1 G 1 ) − 1 ( Λ e ⊗ W 2 G 2 ) ( G 2 T W 2 G 2 ) − 1 ] ⏟ f ( m − 1 ) × 6 \begin{aligned} &N_{bb}^{-1}N_{ab}^TQ_{\hat{a}\hat{a}}N_{ab}N_{bb}^{-1}\\ &=\frac{1}{f} \begin{bmatrix}(G_1^TW_1G_1)^{-1}&0\\0&(G_2^TW_2G_2)^{-1}\end{bmatrix} \begin{bmatrix}e^T \Lambda \otimes G_1^TW_1 \\e^T \Lambda \otimes G_2^TW_2\end{bmatrix} Q_{\hat{a}\hat{a}}N_{ab}N_{bb}^{-1}\\ &= \frac{1}{f^2} \underbrace{ \begin{bmatrix}(G_1^TW_1G_1)^{-1}(e^T \Lambda \otimes G_1^TW_1)\\(G_2^TW_2G_2)^{-1}(e^T \Lambda \otimes G_2^TW_2)\end{bmatrix}}_{6\times f(m-1)} Q_{\hat{a}\hat{a}} \underbrace{\begin{bmatrix}(\Lambda e\otimes W_1G_1)(G_1^TW_1G_1)^{-1}&(\Lambda e\otimes W_2G_2)(G_2^TW_2G_2)^{-1}\end{bmatrix}}_{f(m-1)\times 6} \end{aligned} Nbb−1NabTQa^a^NabNbb−1=f1[(G1TW1G1)−100(G2TW2G2)−1][eTΛ⊗G1TW1eTΛ⊗G2TW2]Qa^a^NabNbb−1=f216×f(m−1) [(G1TW1G1)−1(eTΛ⊗G1TW1)(G2TW2G2)−1(eTΛ⊗G2TW2)]Qa^a^f(m−1)×6 [(Λe⊗W1G1)(G1TW1G1)−1(Λe⊗W2G2)(G2TW2G2)−1]
对于 Q a ^ a ^ Q_{\hat{a}\hat{a}} Qa^a^的组成 1 ◯ \textcircled{1} 1◯，我们以 i = 1 i=1 i=1为例证明上式的结果为 0 0 0：
上 式 1 ◯ , i = 1 = σ ϕ 2 f 2 ( G 1 T W 1 G 1 ) − 1 ( e T Λ ⊗ G 1 T W 1 ) ( Λ − 1 P e ⊥ Λ − 1 ⊗ W − 1 ) ( Λ e ⊗ W 1 G 1 ) ( G 1 T W 1 G 1 ) − 1 = σ ϕ 2 f 2 ( G 1 T W 1 G 1 ) − 1 ( e T P e ⊥ Λ − 1 ⊗ G 1 T W 1 W − 1 ) ( Λ e ⊗ W 1 G 1 ) ( G 1 T W 1 G 1 ) − 1 = σ ϕ 2 f 2 ( G 1 T W 1 G 1 ) − 1 ( e T P e ⊥ e ⏟ = 0 ⊗ G 1 T W 1 W − 1 W 1 G 1 ) ( G 1 T W 1 G 1 ) − 1 = 0 \begin{aligned} 上式_{\textcircled{1},i=1}&=\frac{\sigma_{\phi}^2}{f^2}(G_1^TW_1G_1)^{-1}(e^T \Lambda \otimes G_1^TW_1)(\Lambda^{-1}\mathcal{P}_e^{\perp}\Lambda^{-1}\otimes W^{-1})(\Lambda e\otimes W_1G_1)(G_1^TW_1G_1)^{-1}\\ &=\frac{\sigma_{\phi}^2}{f^2}(G_1^TW_1G_1)^{-1}(e^T\mathcal{P}_e^{\perp}\Lambda^{-1} \otimes G_1^TW_1W^{-1})(\Lambda e\otimes W_1G_1)(G_1^TW_1G_1)^{-1}\\ &=\frac{\sigma_{\phi}^2}{f^2}(G_1^TW_1G_1)^{-1}(\underbrace{e^T\mathcal{P}_e^{\perp}e}_{=0}\otimes G_1^TW_1W^{-1}W_1G_1)(G_1^TW_1G_1)^{-1}\\ &=0 \end{aligned} 上式1◯,i=1=f2σϕ2(G1TW1G1)−1(eTΛ⊗G1TW1)(Λ−1Pe⊥Λ−1⊗W−1)(Λe⊗W1G1)(G1TW1G1)−1=f2σϕ2(G1TW1G1)−1(eTPe⊥Λ−1⊗G1TW1W−1)(Λe⊗W1G1)(G1TW1G1)−1=f2σϕ2(G1TW1G1)−1(=0 eTPe⊥e⊗G1TW1W−1W1G1)(G1TW1G1)−1=0
其中 e T P e ⊥ e = 0 e^T\mathcal{P}_e^{\perp}e=0 eTPe⊥e=0，将 P e ⊥ \mathcal{P}_e^{\perp} Pe⊥的表达式代入即可求得。
接下来求解上式关于 Q a ^ a ^ Q_{\hat{a}\hat{a}} Qa^a^的组成 2 ◯ \textcircled{2} 2◯的计算结果：
上 式 2 ◯ = 2 σ ϕ 2 f 2 [ ( G 1 T W 1 G 1 ) − 1 ( e T Λ ⊗ G 1 T W 1 ) ( G 2 T W 2 G 2 ) − 1 ( e T Λ ⊗ G 2 T W 2 ) ] ( Λ − 1 P e Λ − 1 ⊗ Q ρ ^ ρ ^ ) [ ( G 1 T W 1 G 1 ) − 1 ( e T Λ ⊗ G 1 T W 1 ) ( G 2 T W 2 G 2 ) − 1 ( e T Λ ⊗ G 2 T W 2 ) ] T = 2 σ ϕ 2 f 2 [ ( G 1 T W 1 G 1 ) − 1 ( e T P e Λ − 1 ⊗ G 1 T W 1 Q ρ ^ ρ ^ ) ( G 2 T W 2 G 2 ) − 1 ( e T P e Λ − 1 ⊗ G 2 T W 2 Q ρ ^ ρ ^ ) ] [ ( Λ e ⊗ W 1 G 1 ) ( G 1 T W 1 G 1 ) − 1 ( Λ e ⊗ W 2 G 2 ) ( G 2 T W 2 G 2 ) − 1 ] = 2 σ ϕ 2 f [ G 1 + Q ρ ^ ρ ^ G 1 + T G 1 + Q ρ ^ ρ ^ G 2 + T G 2 + Q ρ ^ ρ ^ G 1 + T G 2 + Q ρ ^ ρ ^ G 2 + T ] \begin{aligned} 上式_{\textcircled{2}} &=\frac{2\sigma_{\phi}^2}{f^2} \begin{bmatrix}(G_1^TW_1G_1)^{-1}(e^T \Lambda \otimes G_1^TW_1)\\(G_2^TW_2G_2)^{-1}(e^T \Lambda \otimes G_2^TW_2)\end{bmatrix} (\Lambda^{-1}\mathcal{P}_e\Lambda^{-1}\otimes Q_{\hat{\rho}\hat{\rho}}) \begin{bmatrix}(G_1^TW_1G_1)^{-1}(e^T \Lambda \otimes G_1^TW_1)\\(G_2^TW_2G_2)^{-1}(e^T \Lambda \otimes G_2^TW_2)\end{bmatrix}^T\\ &=\frac{2\sigma_{\phi}^2}{f^2} \begin{bmatrix}(G_1^TW_1G_1)^{-1}(e^T \mathcal{P}_e\Lambda^{-1} \otimes G_1^TW_1Q_{\hat{\rho}\hat{\rho}})\\(G_2^TW_2G_2)^{-1}(e^T \mathcal{P}_e\Lambda^{-1}\otimes G_2^TW_2Q_{\hat{\rho}\hat{\rho}})\end{bmatrix} \begin{bmatrix}(\Lambda e\otimes W_1G_1)(G_1^TW_1G_1)^{-1}&(\Lambda e\otimes W_2G_2)(G_2^TW_2G_2)^{-1}\end{bmatrix}\\ &=\frac{2\sigma_{\phi}^2}{f}\begin{bmatrix} G_1^+Q_{\hat{\rho}\hat{\rho}}G_1^{+T}&G_1^{+}Q_{\hat{\rho}\hat{\rho}}G_2^{+T} \\ G_2^{+}Q_{\hat{\rho}\hat{\rho}}G_1^{+T}&G_2^{+}Q_{\hat{\rho}\hat{\rho}}G_2^{+T} \end{bmatrix} \end{aligned} 上式2◯=f22σϕ2[(G1TW1G1)−1(eTΛ⊗G1TW1)(G2TW2G2)−1(eTΛ⊗G2TW2)](Λ−1PeΛ−1⊗Qρ^ρ^)[(G1TW1G1)−1(eTΛ⊗G1TW1)(G2TW2G2)−1(eTΛ⊗G2TW2)]T=f22σϕ2[(G1TW1G1)−1(eTPeΛ−1⊗G1TW1Qρ^ρ^)(G2TW2G2)−1(eTPeΛ−1⊗G2TW2Qρ^ρ^)][(Λe⊗W1G1)(G1TW1G1)−1(Λe⊗W2G2)(G2TW2G2)−1]=f2σϕ2[G1+Qρ^ρ^G1+TG2+Qρ^ρ^G1+TG1+Qρ^ρ^G2+TG2+Qρ^ρ^G2+T]
所以：
Q b ^ b ^ = Q b ˇ b ˇ + N b b − 1 N a b T Q a ^ a ^ N a b N b b − 1 = 2 σ ϕ 2 f [ ( G 1 T W 1 G 1 ) − 1 0 0 ( G 2 T W 2 G 2 ) − 1 ] + 2 σ ϕ 2 f [ G 1 + Q ρ ^ ρ ^ G 1 + T G 1 + Q ρ ^ ρ ^ G 2 + T G 2 + Q ρ ^ ρ ^ G 1 + T G 2 + Q ρ ^ ρ ^ G 2 + T ] \begin{equation} \begin{aligned} Q_{\hat{b}\hat{b}} &=Q_{\check{b}\check{b}} + N_{bb}^{-1}N_{ab}^TQ_{\hat{a}\hat{a}}N_{ab}N_{bb}^{-1}\\ &=\frac{2\sigma_{\phi}^2}{f}\begin{bmatrix} (G_1^TW_1G_1)^{-1}&0 \\ 0&(G_2^TW_2G_2)^{-1} \end{bmatrix}+\frac{2\sigma_{\phi}^2}{f}\begin{bmatrix} G_1^+Q_{\hat{\rho}\hat{\rho}}G_1^{+T}&G_1^{+}Q_{\hat{\rho}\hat{\rho}}G_2^{+T} \\ G_2^{+}Q_{\hat{\rho}\hat{\rho}}G_1^{+T}&G_2^{+}Q_{\hat{\rho}\hat{\rho}}G_2^{+T} \end{bmatrix} \end{aligned} \end{equation} Qb^b^=Qbˇbˇ+Nbb−1NabTQa^a^NabNbb−1=f2σϕ2[(G1TW1G1)−100(G2TW2G2)−1]+f2σϕ2[G1+Qρ^ρ^G1+TG2+Qρ^ρ^G1+TG1+Qρ^ρ^G2+TG2+Qρ^ρ^G2+T]
由上式可以注意到 Q b ^ b ^ Q_{\hat{b}\hat{b}} Qb^b^由固定解方差 Q b ˇ b ˇ Q_{\check{b}\check{b}} Qbˇbˇ和模糊度相关项组成。同时 Q b ^ b ^ Q_{\hat{b}\hat{b}} Qb^b^存在非对角线部分，体现了两个历元基线参数之间的协方差，原因为即使基线 b 1 b_1 b1和 b 2 b_2 b2被建模为独立的参数，它们仍通过共享的模糊度 a a a间接关联,模糊度的估计误差 ( a ^ − a ) (\hat{a}-a) (a^−a)会被传递到两个历元的基线解算中,导致两者的解之间存在统计相关性。
![问卷调查[mqtt dht]](http://pic.xiahunao.cn/问卷调查[mqtt dht])
