From 175389b149226271dae86b4b477b9d4bed741452 Mon Sep 17 00:00:00 2001 From: Josh Gao Date: Mon, 5 Oct 2026 22:03:36 +0800 Subject: [PATCH] =?UTF-8?q?content(app-G):=20=E6=98=8E=E7=A1=AE=E4=BD=8D?= =?UTF-8?q?=E7=BD=AE=E9=80=9F=E5=BA=A6=E6=B7=B7=E5=90=88=E5=8D=8F=E8=AE=AE?= =?UTF-8?q?=E4=B8=8E=E6=96=B9=E5=90=91=E7=AC=A6=E5=8F=B7?= MIME-Version: 1.0 Content-Type: text/plain; charset=UTF-8 Content-Transfer-Encoding: 8bit 区分候选角导数和位置轨迹的时间导数，给出频移近似误差。 明确保留混合协议及统一可微轨迹时各接口的修改条件。 补充独立差分、固定节点 OLS 与方向符号检查，全部通过。 XeLaTeX 多遍与引用检查通过；PDF 转换受外部进程文件占用阻塞。 --- check_finite_range_bearing.py | 56 +++++++++++++++++++++++ parts/G-userPositionCrb.tex | 85 ++++++++++++++++++++++++++++++++++- 2 files changed, 139 insertions(+), 2 deletions(-) diff --git a/check_finite_range_bearing.py b/check_finite_range_bearing.py index 8a805c9..6b39f0e 100644 --- a/check_finite_range_bearing.py +++ b/check_finite_range_bearing.py @@ -74,6 +74,35 @@ def geometry(x, t, ell): return w, W, wd, WD, b, bx, bd, bdx +def mixed_motion(x, t, ell): + """Synthetic non-spin rotation: deliberately not an astronomical model.""" + p0, P0 = intersection(x) + a, k = 7.292115e-5*t, 3e-5 + z = np.array([[np.cos(a), -np.sin(a), 0.], + [np.sin(a), np.cos(a), 0.], [0., 0., 1.]]) + a = k*t + q = np.array([[1., 0., 0.], [0., np.cos(a), -np.sin(a)], + [0., np.sin(a), np.cos(a)]]) + generator = np.array([[0., 0., 0.], [0., 0., -k], [0., k, 0.]]) + B = z@q + Bd = SPIN@B + B@generator + C, Cd = SPIN@B, SPIN@Bd + p, pd, vu = B@p0, Bd@p0, C@p0 + s = SATELLITES[ell] + t*VELOCITIES[ell] + .5*t*t*ACCELERATIONS[ell] + vs = VELOCITIES[ell] + t*ACCELERATIONS[ell] + r = np.linalg.norm(s-p) + w = (s-p)/r + S = np.eye(3)-np.outer(w, w) + wd = S@(vs-pd)/r + proxy = S@(vs-vu)/r + W = -S@B@P0/r + b = -(vs-vu)@w/WAVELENGTH + bg = -(vs-pd)@w/WAVELENGTH + bd = -((ACCELERATIONS[ell]-Cd@p0)@w + (vs-vu)@wd)/WAVELENGTH + bx = -(-w@C@P0+(vs-vu)@W)/WAVELENGTH + return p, pd, vu, w, wd, proxy, b, bg, bd, bx + + def coefficients(x, ell): data = [geometry(x, t, ell) for t in TIMES] ref = [geometry(x, t, 0) for t in TIMES] @@ -168,6 +197,33 @@ def main(): Tref[:2, :2] = finite_difference(lambda x: bearing(intersection(x)[0], SATELLITES[1]), x) Tref[2:, :2] = Cx check("reference-change full Jacobian", Tref, finite_difference(switch, z), 2e-6) + for ell in range(2): + p, pd, vu, w, wd, proxy, b, bg, bd, bx = mixed_motion(x, .41, ell) + step = 1e-3 + check("mixed position velocity", pd, + (mixed_motion(x, .41+step, ell)[0]-mixed_motion(x, .41-step, ell)[0])/(2*step), + 2e-6) + check("mixed true LOS time derivative", wd, + (mixed_motion(x, .41+step, ell)[3]-mixed_motion(x, .41-step, ell)[3])/(2*step)) + check("mixed node Doppler Jacobian", bx, + finite_difference(lambda x: mixed_motion(x, .41, ell)[6], x)) + check("mixed Doppler time derivative", bd, + (mixed_motion(x, .41+step, ell)[6]-mixed_motion(x, .41-step, ell)[6])/(2*step)) + check("mixed Doppler discrepancy", b-bg, (vu-pd)@w/WAVELENGTH) + check("mixed proxy LOS discrepancy", proxy-wd, + -(np.eye(3)-np.outer(w, w))@(vu-pd)/np.linalg.norm( + SATELLITES[ell]+.41*VELOCITIES[ell]+.5*.41**2*ACCELERATIONS[ell]-p)) + assert np.linalg.norm(vu-pd) > 1, "Fixture must expose the mixed protocol" + e = -w + check("orientation steering identity", + np.exp(-2j*np.pi*ARRAY@w/WAVELENGTH), + np.exp(2j*np.pi*ARRAY@e/WAVELENGTH)) + check("orientation Doppler identity", b, (VELOCITIES[ell]+.41*ACCELERATIONS[ell]-vu)@e/WAVELENGTH) + mixed_coeff = lambda x: np.linalg.lstsq(DESIGN, [ + mixed_motion(x, t, 1)[6]-mixed_motion(x, t, 0)[6] for t in TIMES], rcond=None)[0] + mixed_jac = np.linalg.lstsq(DESIGN, [ + mixed_motion(x, t, 1)[9]-mixed_motion(x, t, 0)[9] for t in TIMES], rcond=None)[0] + check("mixed fixed-design OLS Jacobian", mixed_jac, finite_difference(mixed_coeff, x)) print("PASS: synthetic local formula checks; code-side geometry and CP/IP acceptance still required.") diff --git a/parts/G-userPositionCrb.tex b/parts/G-userPositionCrb.tex index d31c15a..16e0362 100644 --- a/parts/G-userPositionCrb.tex +++ b/parts/G-userPositionCrb.tex @@ -290,7 +290,67 @@ Doppler 部分，由 $\nu_\ell=\nu_{\mathrm{ref}}+\Delta\nu_\ell(\mathbf p)$ 与 % 这里 $\mathbf R(t)$ 为地固系到惯性系的已知旋转；导数按物理秒计算。$\mathbf V_\xi$、$\mathbf A_\xi$ 分别是卫星--源相对速度与相对加速度对两角的 Jacobian，卫星星历不随候选角变化。匀速自转近似下 $\dot{\mathbf B}=\boldsymbol\Omega_{\mathrm E}\mathbf B$、$\ddot{\mathbf B}=\boldsymbol\Omega_{\mathrm E}^2\mathbf B$。若选用惯性静止源特例，则令 $\mathbf B=\mathbf I_3$ 且两阶导数为零；两种源运动约定须在生成与估计中一致。 -令 $\mathbf v_\ell^{\mathrm{rel}}=\mathbf v_\ell-\mathbf v_u$、$\mathbf a_\ell^{\mathrm{rel}}=\mathbf a_\ell^{\mathrm{orb}}-\mathbf a_u$。各星传播方向、时间导数与候选角导数分别为 +\paragraph{位置映射与速度近似分别给定时} +部分实现采用逐历元坐标旋转给出位置，同时另用简化自转速度给出 Doppler。可将这一确定性协议写为 +% +\begin{equation} + \label{eq:G_bearing_mixed_motion} + \begin{aligned} + \mathbf p(t)&=\mathbf B(t)\mathbf p_0,\qquad + \widetilde{\mathbf v}_u(t)=\mathbf C(t)\mathbf p_0,\\ + \frac{\partial\mathbf p}{\partial\boldsymbol\xi^\Tran} + &=\mathbf B(t)\mathbf P_\xi,\qquad + \frac{\partial\widetilde{\mathbf v}_u}{\partial\boldsymbol\xi^\Tran} + =\mathbf C(t)\mathbf P_\xi,\\ + \boldsymbol\epsilon_v(t) + &=\widetilde{\mathbf v}_u(t)-\dot{\mathbf p}(t) + =(\mathbf C(t)-\dot{\mathbf B}(t))\mathbf p_0 . + \end{aligned} +\end{equation} +% +这里 $\mathbf B(t)$ 与 $\mathbf C(t)$ 为已知且不随候选角变化的矩阵。例如，$\widetilde{\mathbf v}_u=\boldsymbol\Omega_{\mathrm E}\mathbf p$ 对应 $\mathbf C=\boldsymbol\Omega_{\mathrm E}\mathbf B$。只有当 $\dot{\mathbf B}=\boldsymbol\Omega_{\mathrm E}\mathbf B$ 时，它才是该位置轨迹的精确速度；包含其他时间变化旋转的坐标转换须使用其实际导数检查这一等式。候选角导数与时间导数是不同操作：即使 $\mathbf C\ne\dot{\mathbf B}$，\eqref{eq:G_bearing_mixed_motion} 的两个角 Jacobian 仍分别成立。 + +令 $\mathbf v_\ell^{\,g}=\mathbf v_\ell-\dot{\mathbf B}\mathbf p_0$ 为位置轨迹的实际相对速度，$\widetilde{\mathbf v}_\ell^{\mathrm{rel}}=\mathbf v_\ell-\mathbf C\mathbf p_0$ 为 Doppler 近似中采用的相对速度。沿用逐历元位置计算的 LOS，有 +% +\begin{equation} + \label{eq:G_bearing_mixed_time_derivatives} + \begin{aligned} + \dot{\mathbf w}_\ell + &=\frac{\mathbf S_\ell\mathbf v_\ell^{\,g}}{r_\ell},\\ + \widetilde b_\ell + &=-\frac{(\widetilde{\mathbf v}_\ell^{\mathrm{rel}})^\Tran\mathbf w_\ell} + {\lambda_{\mathrm c}},\\ + \frac{\mathrm d\widetilde b_\ell}{\mathrm dt} + &=-\frac{ + (\mathbf a_\ell^{\mathrm{orb}}-\dot{\mathbf C}\mathbf p_0)^\Tran\mathbf w_\ell + +(\widetilde{\mathbf v}_\ell^{\mathrm{rel}})^\Tran\dot{\mathbf w}_\ell} + {\lambda_{\mathrm c}},\\ + \frac{\partial\widetilde b_\ell}{\partial\boldsymbol\xi^\Tran} + &=-\frac{-\mathbf w_\ell^\Tran\mathbf C\mathbf P_\xi+ + (\widetilde{\mathbf v}_\ell^{\mathrm{rel}})^\Tran\mathbf W_\ell} + {\lambda_{\mathrm c}}. + \end{aligned} +\end{equation} +% +其中 $\mathbf S_\ell$、$\mathbf W_\ell$ 的定义见 \eqref{eq:G_bearing_los_derivatives}。若 $\mathbf C=\boldsymbol\Omega_{\mathrm E}\mathbf B$，则 $\dot{\mathbf C}=\dot{\boldsymbol\Omega}_{\mathrm E}\mathbf B+\boldsymbol\Omega_{\mathrm E}\dot{\mathbf B}$，一般不能换成 $\boldsymbol\Omega_{\mathrm E}^2\mathbf B$。时间导数第一行使用实际位置轨迹的速度；将 $\widetilde{\mathbf v}_\ell^{\mathrm{rel}}$ 代入所得量只是代理 LOS rate，与真实导数的差为 $-\mathbf S_\ell\boldsymbol\epsilon_v/r_\ell$。角度相关的时间导数可直接对 \eqref{eq:G_bearing_mixed_time_derivatives} 求导，使用 $\mathbf B\mathbf P_\xi$、$\dot{\mathbf B}\mathbf P_\xi$、$\mathbf C\mathbf P_\xi$ 和 $\dot{\mathbf C}\mathbf P_\xi$，不混用两种速度。 + +保留该近似时，\eqref{eq:G_bearing_window_ols} 的节点向量改用 $\widetilde b_\ell-\widetilde b_r$，节点 Jacobian 改用 \eqref{eq:G_bearing_mixed_time_derivatives} 的最后一行，再乘固定的 $\mathbf X^\dagger$。这不要求 $\mathbf C=\dot{\mathbf B}$，也不要求计算瞬时 rate；OLS rate 是所选节点频率直线的斜率。由同一节点协议生成、拟合和求导的均值模型可在数学上自洽，其参考 bearing 与经纬度的坐标变换、实信息矩阵和干扰消元仍成立。物理上，该频移近似偏离位置轨迹的距离变化率，须另行量化： +% +\begin{equation} + \label{eq:G_bearing_motion_discrepancy} + \widetilde b_\ell-b_\ell^{\,g} + =\frac{\boldsymbol\epsilon_v^\Tran\mathbf w_\ell}{\lambda_{\mathrm c}}, + \qquad + (\widetilde b_\ell-\widetilde b_r)-(b_\ell^{\,g}-b_r^{\,g}) + =\frac{\boldsymbol\epsilon_v^\Tran(\mathbf w_\ell-\mathbf w_r)} + {\lambda_{\mathrm c}}. +\end{equation} +% +有限距离的不同视线使该误差通常不能全部吸收到公共参考频移。数值检查应分别差分位置轨迹、近似速度和节点 Doppler：前者用于检验 LOS 时间导数，后两者用于检验近似 Doppler 的角导数及瞬时导数。微分步长采用物理秒或弧度，须检查步长收敛，并区分解析旋转、历元插值与时间转换在工作点的可微性。 + +若要求同一个可微物理源轨迹，可选择原有 $\mathbf B(t)$ 并令 $\mathbf C=\dot{\mathbf B}$，或明确采用匀速自转近似 $\mathbf B(t)=\exp((t-t_0)\boldsymbol\Omega_{\mathrm E})$。前者更新源速度和加速度，后者同时更新逐历元源位置；两者均不增加未知源运动参数。生成、候选、OLS 节点与 Jacobian、瞬时 rate、信息矩阵及高保真积分相位须使用同一选择，重建受影响的参考系数与目标真值。仅改 CRB 侧速度导数会破坏同模关系；这种物理协议更新也不再是旧样本的纯坐标转换。是否更新生产协议由场景与近似误差的证据决定，局部自洽性检查本身不要求静默更换协议。 + +以下 LOS 与 Doppler rate 公式采用 \eqref{eq:G_bearing_source_motion} 的一致轨迹；保留分别给定的位置和速度时，使用 \eqref{eq:G_bearing_mixed_time_derivatives}。令 $\mathbf v_\ell^{\mathrm{rel}}=\mathbf v_\ell-\mathbf v_u$、$\mathbf a_\ell^{\mathrm{rel}}=\mathbf a_\ell^{\mathrm{orb}}-\mathbf a_u$。各星传播方向、时间导数与候选角导数分别为 % \begin{equation} \label{eq:G_bearing_los_derivatives} @@ -327,6 +387,27 @@ Doppler 部分，由 $\nu_\ell=\nu_{\mathrm{ref}}+\Delta\nu_\ell(\mathbf p)$ 与 % 几何 Doppler 的时间计算仍保留实际 LOS 变化。空间导向变化很小并不足以证明 Doppler rate 中的 LOS 导数项可忽略。 +\paragraph{来向、传播方向与空间相位符号} +在第 $\ell$ 颗卫星定义星到源的来向 $\mathbf e_\ell=-\mathbf w_\ell$，本地方向由同一姿态旋转给出 $\mathbf e_\ell^{\mathrm{loc}}=\boldsymbol\Omega_\ell^\Tran\mathbf e_\ell$、$\mathbf w_\ell^{\mathrm{loc}}=-\mathbf e_\ell^{\mathrm{loc}}$。沿用时间相位 $e^{j2\pi f t}$ 的约定，本小节模板满足 +% +\begin{equation} + \label{eq:G_bearing_orientation_signs} + \begin{aligned} + {}[\mathbf a_\ell]_m + &=\exp\!\left(-j\frac{2\pi}{\lambda_{\mathrm c}} + (\mathbf r_{\ell,m}^{\mathrm{loc}})^\Tran\mathbf w_\ell^{\mathrm{loc}}\right) + =\exp\!\left(+j\frac{2\pi}{\lambda_{\mathrm c}} + (\mathbf r_{\ell,m}^{\mathrm{loc}})^\Tran\mathbf e_\ell^{\mathrm{loc}}\right),\\ + b_\ell + &=-\frac{(\mathbf v_\ell^{\mathrm{rel}})^\Tran\mathbf w_\ell}{\lambda_{\mathrm c}} + =+\frac{(\mathbf v_\ell^{\mathrm{rel}})^\Tran\mathbf e_\ell}{\lambda_{\mathrm c}} . + \end{aligned} +\end{equation} +% +同样的方向替换适用于 \eqref{eq:G_bearing_mixed_time_derivatives} 的近似相对速度。对应角坐标可写为 $\alpha_w=\operatorname{wrap}(\alpha_e+\pi)$、$\beta_w=-\beta_e$；空间方向 Jacobian 同时满足 $\partial\mathbf w_\ell/\partial\boldsymbol\xi^\Tran=-\partial\mathbf e_\ell/\partial\boldsymbol\xi^\Tran$。因此，改用相反方向作为模板输入时，空间指数符号和 Doppler 投影符号须配套变换。两种输入约定表示同一个模板，而不是改变接收信号。 + +实际阵列实现可能采用不同的阵元位移定义、基带共轭约定或角输入约定，应从其模板定义确定对应关系。在真值及有效域内候选点，分别比较由 $\mathbf w_\ell$ 和 $-\mathbf w_\ell$ 两种写法得到的 steering 与均值，并在固定增益／相位后检查一致性；还应核对角导数和随时间变化的 Doppler 相位。单独翻转方向或单独翻转 Doppler 不能据此判定等价，实际空间指数也不能仅由函数名或 bearing 朝向猜测。 + \paragraph{有限距离 Doppler、rate 与窗口拟合} 沿用 \eqref{eq:A6_latlon_doppler_jac} 的传播方向与 Doppler 符号，几何频移 $b_\ell$（Hz）及其瞬时 rate（Hz/s）为 % @@ -355,7 +436,7 @@ Doppler 部分，由 $\nu_\ell=\nu_{\mathrm{ref}}+\Delta\nu_\ell(\mathbf p)$ 与 \begin{equation} \label{eq:G_bearing_window_ols} \begin{aligned} - [\mathbf d_\ell]_p&=b_\ell(t_0+s_p)-b_r(t_0+s_p),\\ + {}[\mathbf d_\ell]_p&=b_\ell(t_0+s_p)-b_r(t_0+s_p),\\ \mathbf c_\ell&=\begin{bmatrix}\Delta f_\ell\\\Delta\dot f_\ell\end{bmatrix} =\mathbf X^\dagger\mathbf d_\ell,\\ \frac{\partial\mathbf c_\ell}{\partial\boldsymbol\xi^\Tran}