diff --git a/check_finite_range_bearing.py b/check_finite_range_bearing.py index 3fe6d3b..5290254 100644 --- a/check_finite_range_bearing.py +++ b/check_finite_range_bearing.py @@ -1,4 +1,4 @@ -"""Deterministic formula checks for Appendix G; run: python check_finite_range_bearing.py. +"""Deterministic formula checks for Appendices F/G; run: python check_finite_range_bearing.py. Synthetic sphere and known trajectories only: no constellation or MC validation. """ @@ -182,6 +182,53 @@ def main(): assert np.linalg.matrix_rank(projected) == 4 cov = np.linalg.inv(info) check("full covariance reparameterization", T@np.linalg.inv(T.T@info@T)@T.T, cov) + # Appendix F: structured CP has real frame amplitudes and one phase per star. + for selected in [[0], [0, 1]]: + Jcp = np.vstack([jac[ell] for ell in selected]) + N = np.zeros((len(Jcp), len(selected)*(len(TIMES)+1)), dtype=complex) + block_size = len(TIMES)*len(ARRAY) + for star, ell in enumerate(selected): + rows = slice(star*block_size, (star+1)*block_size) + N[rows, star*(len(TIMES)+1)] = 1j*means[ell] + for frame in range(len(TIMES)): + block = slice(star*block_size+frame*len(ARRAY), + star*block_size+(frame+1)*len(ARRAY)) + N[block, star*(len(TIMES)+1)+frame+1] = means[ell][frame*len(ARRAY):(frame+1)*len(ARRAY)] + real_N = np.vstack((N.real, N.imag)) + scores = np.vstack((Jcp.real, Jcp.imag))*np.sqrt(2) + residual_scores = scores-real_N@np.linalg.lstsq(real_N, scores, rcond=None)[0] + I = residual_scores.T@residual_scores + K, coupling, rate_info = I[:3, :3], I[:3, 3], I[3, 3] + CK = np.linalg.inv(K) + CU = np.linalg.inv(K-np.outer(coupling, coupling)/rate_info) + v = CK@coupling + h = rate_info-coupling@v + assert h > 0, "K/U covariance identities require identifiable rate" + check("CP K/U full covariance update", CU-CK, np.outer(v, v)/h, 1e-9) + q = coupling[:2]-K[:2, 2]*coupling[2]/K[2, 2] + G = K[:2, :2]-np.outer(K[:2, 2], K[:2, 2])/K[2, 2] + check("partial angle-rate coupling", v[:2], np.linalg.solve(G, q), 1e-9) + metric = np.diag([np.cos(x[1])**2, 1.]) + check("CP K/U angular metric update", np.trace(metric@(CU[:2, :2]-CK[:2, :2])), + v[:2]@metric@v[:2]/h, 1e-9) + if len(selected) == 1: + check("SS spatial-temporal orthogonality", I[:2, 2:], np.zeros((2, 2)), 1e-10) + check("SS exact angular subblock equality", CU[:2, :2], CK[:2, :2], 1e-10) + assert CU[2, 2] > CK[2, 2], "First-sample time origin couples frequency and rate" + # Nonzero EFIM loss can affect only Doppler after the covariance update. + K = np.array([[3., .3, .5], [.3, 2., -.4], [.5, -.4, 4.]]) + CK = np.linalg.inv(K) + for v in [np.array([0., 0., .7]), np.array([.2, -.3, .7])]: + coupling = K@v + rate_info = coupling@v+2. + CU = np.linalg.inv(K-np.outer(coupling, coupling)/rate_info) + check("synthetic rate-only covariance update", CU-CK, np.outer(v, v)/2., 1e-12) + assert np.linalg.norm(coupling[:2]) > 0 + if np.linalg.norm(v[:2]) == 0: + check("nonzero EFIM loss with equal angle bounds", CU[:2, :2], CK[:2, :2], 1e-12) + assert CU[2, 2] > CK[2, 2] + else: + assert np.trace(metric@(CU[:2, :2]-CK[:2, :2])) > 0 delta = .19 shifted = TIMES-delta for ell in range(2): diff --git a/parts/F-multiFrameCrb.tex b/parts/F-multiFrameCrb.tex index ca59eca..58e4a10 100644 --- a/parts/F-multiFrameCrb.tex +++ b/parts/F-multiFrameCrb.tex @@ -547,8 +547,14 @@ 因此 $\dfrac{\partial \nu_{\ell_{\mathrm{ref}}}}{\partial x}=0$ 且 $\dfrac{\partial \dot\nu_{\ell_{\mathrm{ref}}}}{\partial x}=0$ 对 $x\in\{\alpha,\beta\}$ 成立。 由式~\eqref{eq:A7_qx_cp} 可知，reference-link 的 \gls{doa} 导数不含 Doppler-mapping-induced temporal phase derivative，而只含空间阵列导向导数。 在幅度和公共相位干扰方向被投影消去后，该 \gls{doa} 得分的残差部分是与 $\mathbf a_{\ell_{\mathrm{ref}}}(\boldsymbol\xi)$ 正交的空间切向量；相反，$\nu_{\mathrm{ref}}$ 和 $\dot\nu_{\mathrm{ref}}$ 的残差得分仍具有沿参考导向矢量的时间相位结构。 -在标准窄带阵列流形下，二者的 Fisher 内积为零，因此式~\eqref{eq:A7_SS_xi_rate_orthogonality_condition} 成立。 +这里还要求同一星的阵列响应在各帧冻结，已知导频与角度无关，噪声采用本节的时空白噪声， +参考 Doppler 和 rate 为独立自由标量，并按结构化幅度及公共相位消元。 +在这些时空可分条件下，二者的 Fisher 内积为零，因此式~\eqref{eq:A7_SS_xi_rate_orthogonality_condition} 成立。 这解释了单星多帧下 Doppler rate 已知与未知 \gls{doa} \gls{crb} 的严格重合。 +参考相对映射为零排除了角度的时间相位导数，但这一点本身不替代上述空间--时间正交条件。 +这一结论也适用于\appref{subsec:G_reference_bearing} 的有限距离参考来向参数化， +只需该角坐标为正则局部坐标；已知高度不把两个参考标量改为角度的已知函数。 +观测时间不必居中，首帧首样本为时间原点同样成立。 在多星情形中，非参考星通常满足 % @@ -566,8 +572,84 @@ % 此时 \gls{doa} 导数包含由 $\dfrac{\partial\Delta\nu_\ell}{\partial x}$ 与 $\dfrac{\partial\Delta\dot\nu_\ell}{\partial x}$ 诱导的时间相位分量。 因此，$\mathbf i_{\xi\dot\nu}$ 和 $\mathbf a_{\xi\nu}^{\mathrm K}$ 一般不再严格为零，Doppler rate 已知与未知 \gls{doa} \gls{crb} 的完全相等性也不应被作为一般结论。 -若数值上多星 \gls{doa} \gls{crb} 仍与 Doppler rate 已知基准基本重合，更稳妥的解释是：在所选代表性几何下，\gls{doa} 残差得分与残差化的 Doppler rate 干扰方向的耦合较弱。 -这一定量现象适合作为 \gls{efim} 度量的自洽性检查，而不应升级为主要理论贡献。 +角相关的时间映射不自动保证角度界有可见变化。需进一步消去参考 Doppler， +再比较角度的有效信息及协方差，而不是仅检查 \eqref{eq:A7_CP_U_EFIM_from_K} 是否有非零扣减。 + +\paragraph{消去参考 Doppler 后的角度条件} +沿用 \eqref{eq:A7_CP_K_partition_xi_nu_rate} 的分块，定义 +% +\begin{equation} + \label{eq:A7_angle_rate_partial_coupling} + \begin{aligned} + \mathbf G&=\mathbf A_{\xi\xi}^{\mathrm K} + -\frac{\mathbf a_{\xi\nu}^{\mathrm K}\mathbf a_{\nu\xi}^{\mathrm K}} + {a_{\nu\nu}^{\mathrm K}},\\ + \mathbf q&=\mathbf i_{\xi\dot\nu} + -\frac{\mathbf a_{\xi\nu}^{\mathrm K}i_{\nu\dot\nu}} + {a_{\nu\nu}^{\mathrm K}},\qquad + d=i_{\dot\nu\dot\nu}-\frac{i_{\nu\dot\nu}^2}{a_{\nu\nu}^{\mathrm K}}. + \end{aligned} +\end{equation} +% +假设扩展信息矩阵正定，因而 $a_{\nu\nu}^{\mathrm K}>0$、$d>0$。 +先消去 $\nu_{\mathrm{ref}}$ 后，关于角度和 rate 的矩阵为 +$\left[\begin{smallmatrix}\mathbf G&\mathbf q\\\mathbf q^\Tran&d\end{smallmatrix}\right]$。 +因此，两个角度协方差子块为 +% +\begin{equation} + \label{eq:A7_angle_covariance_KU} + \mathbf C_{\xi}^{\mathrm K}=\mathbf G^{-1},\qquad + \mathbf C_{\xi}^{\mathrm U} + =\left(\mathbf G-\frac{\mathbf q\mathbf q^\Tran}{d}\right)^{-1},\qquad + \mathbf C_{\xi}^{\mathrm U}=\mathbf C_{\xi}^{\mathrm K} + \ \Longleftrightarrow\ \mathbf q=\mathbf0. +\end{equation} +% +这由 Schur 消元的顺序等价性直接得到。$\mathbf q$ 是消去幅相和参考 Doppler 后的角度--rate 耦合； +它可以因投影抵消而为零，即使 $\mathbf i_{\xi\dot\nu}$ 和 +$\mathbf a_{\xi\nu}^{\mathrm K}$ 各自不为零。 +式~\eqref{eq:A7_SS_xi_rate_orthogonality_condition} 是足够条件， +而 \eqref{eq:A7_angle_covariance_KU} 给出这一正定工作点上更一般的判据。 + +\paragraph{完整协方差更新与合成角度指标} +记 $\mathbf C_{\eta}^{\mathrm K}=(\mathcal I_{\eta\eta}^{\mathrm K})^{-1}$，令 +% +\begin{equation} + \label{eq:A7_rate_covariance_update} + \begin{aligned} + \mathbf v&=\mathbf C_{\eta}^{\mathrm K}\mathbf i_{\eta\dot\nu},\qquad + h=i_{\dot\nu\dot\nu} + -\mathbf i_{\dot\nu\eta}\mathbf C_{\eta}^{\mathrm K}\mathbf i_{\eta\dot\nu}>0,\\ + \mathbf C_{\eta}^{\mathrm U}-\mathbf C_{\eta}^{\mathrm K} + &=\frac{\mathbf v\mathbf v^\Tran}{h},\qquad + \mathbf C_{\xi}^{\mathrm U}-\mathbf C_{\xi}^{\mathrm K} + =\frac{\mathbf v_\xi\mathbf v_\xi^\Tran}{h},\qquad + \mathbf v_\xi=\mathbf G^{-1}\mathbf q. + \end{aligned} +\end{equation} +% +将 \eqref{eq:A7_rate_covariance_update} 代回两个矩阵的乘积，或使用单秩逆更新，可验证该恒等式。 +完整 \gls{efim} 相同要求 $\mathbf i_{\eta\dot\nu}=\mathbf0$； +角度子块相同只要求 $\mathbf v_\xi=\mathbf0$，允许参考 Doppler 方差增大。 +对 \eqref{eq:A7_direction_mse_bound} 的同一局部角度指标，有 +% +\begin{equation} + \label{eq:A7_angle_metric_KU_difference} + B_\theta^{\mathrm U}-B_\theta^{\mathrm K} + =\frac{\mathbf v_\xi^\Tran\mathbf W_\xi\mathbf v_\xi}{h},\qquad + \sqrt{B_\theta^{\mathrm U}}-\sqrt{B_\theta^{\mathrm K}} + =\frac{B_\theta^{\mathrm U}-B_\theta^{\mathrm K}} + {\sqrt{B_\theta^{\mathrm U}}+\sqrt{B_\theta^{\mathrm K}}}. +\end{equation} +% +在 $\cos\beta\ne0$ 的正则角坐标中，$\mathbf W_\xi\succ0$， +所以该合成指标严格相等与角度协方差子块严格相等等价； +这来自本模型的半正定单秩更新，不是任意两个协方差具有相同加权迹的一般性质。 +对非满秩的其他投影指标，指标相等只要求 $\mathbf W_\xi^{1/2}\mathbf v_\xi=\mathbf0$。 +经纬度实现转为参考来向时，须先按 \eqref{eq:G_bearing_covariance_transform} 变换完整协方差，保留交叉项。 +多星指标近似一致的直接条件是 \eqref{eq:A7_angle_metric_KU_difference} 相对于报告尺度足够小， +不要求完整 \gls{efim} 扣减矩阵很小。时间、几何、空间信息和参考 Doppler 的消元共同决定这个量。 +若信息奇异，上述普通逆不适用，应先检查可估角度子空间；不能由完整矩阵的伪逆子块判定 K/U 相等。 \subsubsection{放松连续相位模型下的 CRB}