一台钟静止在 \(K'\) 系的原点 \(x'=0\)。它在 \(K\) 系中的位置是 \(x=vt\)。代入式 (3.3) 的第二式： \[t'=\gamma\Big(t-\frac{v\cdot vt}{c^2}\Big)=\gamma t\Big(1-\frac{v^2}{c^2}\Big)=\frac{t}{\gamma}.\] 钟显示的是 \(t'\)。于是在 \(K\) 系看来，运动的钟走过的时间只有 \(t/\gamma\)。把这个结论写成微分形式，并推广到速度随时间变化的情形： \[\tag{3.4} \dd\tau=\sqrt{1-\frac{v^2(t)}{c^2}}\,\dd t .\] \(\tau\) 叫原时（proper time），即钟自己走过的时间；\(t\) 是惯性系 \(K\) 的时间。

波源在 \(K\) 系中以速度 \(\vect{v}\) 运动，在它自己的静止系中以频率 \(\nu_0\) 发出波；观察者静止在 \(K\) 系中远处。设波源到观察者的方向与 \(\vect{v}\) 的夹角为 \(\theta\)（在 \(K\) 系中量度）。用“波峰计数”的办法求观察者测到的频率：

1. 相邻两个波峰，在波源的原时中间隔 \(1/\nu_0\)；由式 (3.4)，在 \(K\) 系中间隔 \(\Delta t_{\mathrm e}=\gamma/\nu_0\)。

2. 在这段时间内波源朝观察者方向前进了 \(v\cos\theta\,\Delta t_{\mathrm e}\)，第二个波峰少走这段路，早到 \(v\cos\theta\,\Delta t_{\mathrm e}/c\)。观察者处两波峰的到达间隔为 \(\Delta t_{\mathrm r}=\Delta t_{\mathrm e}(1-\beta\cos\theta)\)，\(\beta=v/c\)。（观察者在远处，\(\theta\) 在一个周期内的变化可忽略。）

3. 观察者测得的频率 \(\nu=1/\Delta t_{\mathrm r}\)：

\[\tag{3.6} \nu=\frac{\nu_0}{\gamma(1-\beta\cos\theta)} .\] 对 \(\beta\) 展开到二阶：\(1/(1-\beta\cos\theta)=1+\beta\cos\theta+\beta^2\cos^2\theta+\order{\beta^3}\)，\(1/\gamma=1-\beta^2/2+\order{\beta^4}\)，相乘得 \begin{align*} \frac{\nu-\nu_0}{\nu_0}=\;&\underbrace{\beta\cos\theta}_{\text{一阶多普勒}}\\ &+\underbrace{\beta^2\cos^2\theta-\frac{\beta^2}{2}}_{\text{二阶项}}+\order{\beta^3}.\tag{3.7} \end{align*}

第一步：加速火箭。 在惯性系中，一枚长 \(h\) 的火箭从静止开始以加速度 \(g\) 沿 \(+z\) 方向加速（图 1(a)）。\(t=0\) 时底部光源（此刻速度为零）发出频率 \(\nu_{\mathrm b}\) 的光。光到达顶部用时 \(\Delta t\approx h/c\)（火箭在这段时间内移动 \(\frac12g\Delta t^2\)，相对修正为 \(gh/c^2\)，属于更高阶）。此时顶部接收器的速度 \(u=g\Delta t=gh/c\)，方向与光的传播方向相同，即背离光源。

第二步：一阶多普勒。 接收器以 \(u\ll c\) 背离光源运动，相邻两个波峰追上它需要多走 \(u\Delta t_{\mathrm r}\) 的路程，测到的频率降低为 \(\nu_{\mathrm t}=\nu_{\mathrm b}(1-u/c)\)（对运动接收器的一阶多普勒公式；二阶项 \(u^2/c^2=(gh/c^2)^2\) 可忽略）。所以 \[\frac{\nu_{\mathrm t}-\nu_{\mathrm b}}{\nu_{\mathrm b}}=-\frac{u}{c}=-\frac{gh}{c^2}.\] 第三步：等效原理。 把火箭换成静止在均匀引力场 \(g\) 中、高为 \(h\) 的塔（图 1(b)），结果相同：从塔底向上发出的光，在塔顶测得的频率低了 \[\tag{3.9} \frac{\Delta\nu}{\nu}=-\frac{gh}{c^2}.\]

把从 \(\vect{r}_1\) 到 \(\vect{r}_2\) 的光路分成许多小段，每段内场近似均匀。沿一小段上升 \(\dd z\)，局部重力加速度为 \(g\)，引力势的增量 \(\dd\Phi=g\,\dd z\)（\(\Phi\) 是单位质量的引力势能，向上增大）。由式 (3.9)，这一段造成 \(\dd\nu/\nu=-g\,\dd z/c^2=-\dd\Phi/c^2\)。各段相乘（取对数后相加）： \begin{gather*} \ln\frac{\nu_2}{\nu_1}=-\frac{\Phi_2-\Phi_1}{c^2}\quad\Longrightarrow\\ \frac{\nu_2}{\nu_1}=\exp\!\Big(-\frac{\Phi_2-\Phi_1}{c^2}\Big)\approx1-\frac{\Phi_2-\Phi_1}{c^2}. \end{gather*} 光路的水平部分与 \(\Phi\) 无关，不贡献。注意：每一段用的式 (3.9) 本身只准到 \(gh/c^2\) 的一阶，所以这一推导只在 \(\Delta\Phi/c^2\) 的一阶成立，\(\exp\) 形式不应当作二阶的结果（二阶项与“局部的 \(g\) 用谁的尺和钟量”有关，要用广义相对论处理）。地面附近 \(|\Phi|/c^2\approx7\times10^{-10}\)，二阶项约 \(10^{-19}\)，远小于本讲讨论的任何效应。

光在路上真的“丢了频率”吗？考虑一个静态的引力场：场不随时间变化，光源和接收器都静止（图 1(c)）。设想用某个公共的时间标度 \(t\) 来描述这一过程（具体选法见第 3.4 节），由于一切都不随 \(t\) 变化，第 \(1\) 个波峰从底部传到顶部所用的 \(t\) 时间与第 \(n\) 个波峰相同。因此在任意一段 \(t\) 时间 \(\Delta t\) 内，底部发出的波峰数等于顶部收到的波峰数： \[\nu_{\mathrm b}\,\Delta\tau_{\mathrm b}=\nu_{\mathrm t}\,\Delta\tau_{\mathrm t},\] 其中 \(\Delta\tau_{\mathrm b}\)、\(\Delta\tau_{\mathrm t}\) 是在同一段 \(\Delta t\) 中底部钟、顶部钟各自走过的时间，\(\nu_{\mathrm b}\)、\(\nu_{\mathrm t}\) 是用当地的钟测得的频率。于是 \[\tag{3.10} \frac{\Delta\tau_{\mathrm t}}{\Delta\tau_{\mathrm b}}=\frac{\nu_{\mathrm b}}{\nu_{\mathrm t}}=1+\frac{\Phi_{\mathrm t}-\Phi_{\mathrm b}}{c^2} .\] 顶部的钟走得快（\(\Phi_{\mathrm t}>\Phi_{\mathrm b}\)）。“光的红移”和“低处的钟走得慢”是同一件事的两种说法；能直接比较的只有两个钟的读数。

设一台钟在 ECI 系中位于 \(\vect{r}(t)\)，速度 \(\vect{v}=\dd\vect{r}/\dd t\)，当地引力势 \(\Phi(\vect r)\)（\(\Phi<0\)）。

第一步：静止在当地的“辅助钟”。 在 \(\vect r\) 处放一台静止的辅助钟。由式 (3.10)（取“底部”为 \(\vect r\) 处、“顶部”为 \(\Phi=0\) 的无穷远处），它与坐标时的速率之比为 \[\frac{\dd\tau_{\mathrm s}}{\dd t}=1+\frac{\Phi(\vect r)}{c^2}+\order{c^{-4}} .\] 第二步：运动的钟与辅助钟比较。 在辅助钟附近一个足够小的区域里，引力的影响由等效原理可以消去，狭义相对论成立。运动的钟以速度 \(v_{\mathrm{loc}}\)（辅助钟所在处的观测者用当地的尺和钟测得）掠过辅助钟，由式 (3.4) 与钟假设， \[\frac{\dd\tau}{\dd\tau_{\mathrm s}}=1-\frac{v_{\mathrm{loc}}^2}{2c^2}+\order{c^{-4}} .\] 第三步：\(v_{\mathrm{loc}}\) 与坐标速度 \(v\) 的差别。 当地的钟比坐标时慢一个因子 \(1+\Phi/c^2\)；当地的尺与坐标长度之间也只差 \(\order{\Phi/c^2}\) 的因子（弱场度规的空间部分为 \((1-2\Phi/c^2)\dd\vect r^2\)，见本节检验框），所以 \(v_{\mathrm{loc}}^2=v^2[1+\order{\Phi/c^2}]\)，代入后差别是 \(\order{v^2\Phi/c^4}\)，与第二步中略去的项同阶。

第四步：相乘。 \(\dd\tau/\dd t=(\dd\tau/\dd\tau_{\mathrm s})(\dd\tau_{\mathrm s}/\dd t)\)，展开后乘积项 \(\Phi v^2/(2c^4)\) 同样是 \(\order{c^{-4}}\)： \[\tag{3.11} \boxed{\frac{\dd\tau}{\dd t}=1+\frac{\Phi(\vect r)}{c^2}-\frac{v^2}{2c^2}+\order{c^{-4}}}\]

取以 \(\omega\uvec z\) 转动的参考系，令转动系的 \(x'\) 轴始终指向驱动场方向 \((\cos\omega t,\sin\omega t,0)\)。由式 (3.15)，转动系中的进动矢量为 \(\vect\Omega-\omega\uvec z\)： \begin{gather*} \Big(\dv{\vect S}{t}\Big)_{\text{转}}=\vect\Omega_{\mathrm{eff}}\times\vect S, \\ \vect\Omega_{\mathrm{eff}}=\Omega_{\mathrm R}\uvec x'-\Delta\,\uvec z,\qquad \Delta\equiv\omega-\omega_0 . \tag{3.19} \end{gather*} 在转动系中驱动场不再转动，\(\vect\Omega_{\mathrm{eff}}\) 是常矢量，大小 \(W=\sqrt{\Omega_{\mathrm R}^2+\Delta^2}\)，与 \(z\) 轴的夹角 \(\vartheta\) 满足 \(\tan\vartheta=\Omega_{\mathrm R}/(-\Delta)\)。运动就是绕固定轴的匀速转动。\(P_{\mathrm e}\) 只与 \(S_z\) 有关，而 \(z\) 轴在两个参考系中相同，所以转动系中算出的 \(P_{\mathrm e}\) 就是实验室中的 \(P_{\mathrm e}\)。

原子开始时在基态，\(\vect S(0)=-\uvec z\)。在转动系中绕单位矢量 \(\uvec n=\vect\Omega_{\mathrm{eff}}/W=(\Omega_{\mathrm R},0,-\Delta)/W\) 以角速度 \(W\) 转动。把 \(\vect S(0)\) 分解为平行和垂直于 \(\uvec n\) 的部分： \[\vect S_\parallel=(\vect S(0)\cdot\uvec n)\uvec n=\frac{\Delta}{W}\uvec n,\qquad \vect S_\perp=\vect S(0)-\vect S_\parallel .\] 平行部分不变；垂直部分在垂直于 \(\uvec n\) 的平面内转动：\(\vect S_\perp(t)=\vect S_\perp\cos Wt+(\uvec n\times\vect S_\perp)\sin Wt\)。取 \(z\) 分量：\(S_{\parallel z}=\frac{\Delta}{W}\cdot\frac{-\Delta}{W}=-\frac{\Delta^2}{W^2}\)，\(S_{\perp z}=-1+\frac{\Delta^2}{W^2}=-\frac{\Omega_{\mathrm R}^2}{W^2}\)；\((\uvec n\times\vect S_\perp)_z=(\uvec n\times\vect S(0))_z=-(n_x\uvec x+n_z\uvec z)\times\uvec z\cdot\uvec z=0\)。于是 \[S_z(t)=-\frac{\Delta^2}{W^2}-\frac{\Omega_{\mathrm R}^2}{W^2}\cos Wt ,\] 代入式 (3.17)，用 \(1-\cos Wt=2\sin^2(Wt/2)\)： \[\tag{3.20} P_{\mathrm e}(t)=\frac{\Omega_{\mathrm R}^2}{\Omega_{\mathrm R}^2+\Delta^2}\sin^2\frac{\sqrt{\Omega_{\mathrm R}^2+\Delta^2}\;t}{2} .\]

对一个原子测量一次：以概率 \(p\) 发现它在上能级（记 \(X=1\)），以概率 \(1-p\) 在下能级（\(X=0\)）。\(\avg X=p\)，\(\avg{X^2}=p\)，所以 \(\mathrm{Var}(X)=p-p^2=p(1-p)\)。\(N\) 个独立原子中处于上能级的个数 \(N_{\mathrm e}=\sum X_k\) 服从二项分布，由上一框 \begin{gather*} \avg{N_{\mathrm e}}=Np, \\ \mathrm{Var}(N_{\mathrm e})=Np(1-p), \\ \sigma_{\hat p}=\sqrt{\frac{p(1-p)}{N}}, \tag{3.21} \end{gather*} 其中 \(\hat p=N_{\mathrm e}/N\) 是测得的激发比例。即使每个原子的状态完全相同、仪器毫无噪声，\(\hat p\) 仍有这个起伏，它来自测量结果的随机性，叫量子投影噪声[29]。\(p=0\) 或 \(1\) 时噪声为零，\(p=1/2\) 时最大。

设钟每隔 \(T_{\mathrm c}\) 给出一个频率测量值 \(y_j\)，各次相互独立、方差为 \(\sigma_1^2\)（称为白频率噪声，例如量子投影噪声）。一段 \(\tau_{\mathrm a}\) 内有 \(n=\tau_{\mathrm a}/T_{\mathrm c}\) 个值，平均值的方差为 \(\sigma_1^2/n\)。相邻两段独立，差的方差为 \(2\sigma_1^2/n\)，代入式 (3.23)： \[\tag{3.24} \sigma_y(\tau_{\mathrm a})=\frac{\sigma_1}{\sqrt n}=\sigma_1\sqrt{\frac{T_{\mathrm c}}{\tau_{\mathrm a}}} .\] 在双对数图上是斜率 \(-1/2\) 的直线。这时 Allan 偏差就等于“平均 \(\tau_{\mathrm a}\) 后的频率测量的标准差”。

\(P_{\mathrm e}(\Delta)\) 在共振处取极大、对 \(\Delta\) 偶对称，在峰顶 \(\dd P_{\mathrm e}/\dd\Delta=0\)，测不出 \(\Delta\) 的符号。办法是交替地把频率放在 \(\omega+\Delta_{\mathrm m}\) 与 \(\omega-\Delta_{\mathrm m}\)（调制深度 \(\Delta_{\mathrm m}\) 取在半高处），两次测得的概率相减： \begin{align*} \xi(\Delta)&\equiv P_{\mathrm e}(\Delta+\Delta_{\mathrm m})-P_{\mathrm e}(\Delta-\Delta_{\mathrm m})\\ &\approx 2P_{\mathrm e}'(\Delta_{\mathrm m})\,\Delta\qquad(|\Delta|\ll\Delta_{\mathrm m}), \end{align*} 这里对 \(\Delta\) 作了一阶泰勒展开，并用了 \(P_{\mathrm e}'(-\Delta_{\mathrm m})=-P_{\mathrm e}'(\Delta_{\mathrm m})\)（偶函数的导数是奇函数）。\(\xi\) 正比于 \(\Delta\)，符号指示本振偏高还是偏低，叫误差信号。斜率 \(|P_{\mathrm e}'|\) 在半高处最大，所以调制深度取半宽。

第 \(k\) 个循环本振的失谐为 \(\Delta_k\)，测得误差信号 \(\xi_k=-\kappa\Delta_k+(\text{噪声})\)，\(\kappa=2|P_{\mathrm e}'(\Delta_{\mathrm m})|>0\)（式中取 \(P_{\mathrm e}'(\Delta_{\mathrm m})<0\)，峰的右侧斜率为负）。修正规则 \(\Delta_{k+1}=\Delta_k+G\xi_k/\kappa\)（\(G\) 为无量纲增益）给出 \(\Delta_{k+1}=(1-G)\Delta_k\)：\(0<G<2\) 时失谐按几何级数衰减到零，\(G=1\) 一步到位，\(G\) 小则响应慢但平均掉更多噪声。锁定后本振的平均频率就是原子频率，本振的频率起伏在比若干个循环更长的时间尺度上被原子“纠正”。

在转动系中逐段计算（图 3）。

(1) 初态 \(\vect S=-\uvec z\)。第一个 \(\pi/2\) 脉冲绕 \(x'\) 轴转 \(90^\circ\)：由 \(\dd\vect S/\dd t=\Omega_{\mathrm R}\uvec x'\times\vect S\)，起始时 \(\uvec x'\times(-\uvec z)=+\uvec y'\)，\(\vect S\) 朝 \(+y'\) 转去，转过 \(90^\circ\) 后 \(\vect S=\uvec y'\)。

(2) 自由演化 \(T\)：\(\dd\vect S/\dd t=-\Delta\uvec z\times\vect S\)，\(\vect S\) 在赤道面内绕 \(z\) 轴转过角度 \(-\Delta T\)。因为 \(-\Delta\uvec z\times\uvec y'=+\Delta\uvec x'\)，\(\vect S=\sin(\Delta T)\,\uvec x'+\cos(\Delta T)\,\uvec y'\)。物理意义：原子的叠加态以 \(\omega_0\) 自由演化，而参考系（即振荡器的相位）以 \(\omega\) 转动；\(\Delta T\) 是 \(T\) 内两者积累的相位差。

(3) 第二个 \(\pi/2\) 脉冲绕 \(x'\) 轴再转 \(90^\circ\)：\(\uvec x'\) 分量不变，\(\uvec y'\to\uvec z\)（由 \(\uvec x'\times\uvec y'=\uvec z\)）。末态 \(\vect S=\sin(\Delta T)\uvec x'+\cos(\Delta T)\uvec z\)。

代入式 (3.17)： \[\tag{4.3} P_{\mathrm e}=\frac{1+\cos\Delta T}{2}=\cos^2\frac{\Delta T}{2} .\]

设最高点比微波腔高 \(H\)。原子在腔处的速度 \(v_0=\sqrt{2gH}\)（能量守恒），从向上穿过腔到向下穿过腔的时间是竖直上抛的往返时间 \[\tag{4.5} T=\frac{2v_0}{g}=2\sqrt{\frac{2H}{g}} .\] \(H={1\,}\mathrm{m}\) 时 \(T={0.90\,}\mathrm{s}\)，由式 (4.4) 条纹宽 \(\Delta\nu={0.55\,}\mathrm{Hz}\)，\(Q=2\nu_0T=1.7\times10^{10}\)，比原子束钟窄约 200 倍。

设锁定点在 Ramsey 条纹的半高处（第 4.2 节），即 \(\Delta T=\pm\pi/2\)。由式 (4.3)， \begin{gather*} p=\frac{1+\cos\Delta T}{2},\qquad \dv{p}{\Delta}=-\frac{T}{2}\sin\Delta T ,\\ \Big|\dv{p}{\Delta}\Big|_{\Delta T=\pi/2}=\frac{T}{2},\qquad p=\frac12 . \end{gather*} 一次循环中测量 \(N\) 个原子，测得的激发比例 \(\hat p\) 有量子投影噪声 \(\sigma_{\hat p}=\sqrt{p(1-p)/N}=1/(2\sqrt N)\)（式 (3.21)）。把它换算成对 \(\Delta\) 的估计误差：概率的起伏除以斜率（局部一阶泰勒展开，要求 \(\sigma_{\hat p}\ll1\)）， \[\sigma_\Delta=\frac{\sigma_{\hat p}}{|\dd p/\dd\Delta|}=\frac{1/(2\sqrt N)}{T/2}=\frac{1}{T\sqrt N}\quad(\mathrm{rad/s}).\] 换成分数频率（除以 \(\omega_0\)），一次循环的频率误差为 \(\sigma_1=\dfrac{1}{\omega_0T\sqrt N}\)。

实际伺服在两侧交替探测（第 4.2 节）：误差信号 \(\xi=\hat p_+-\hat p_-\) 的标准差为 \(\sqrt2\,\sigma_{\hat p}\)（两次独立测量的方差相加），斜率为 \(|\dd\xi/\dd\Delta|=T\)，所以每两个循环给出一个方差为 \(2\sigma_{\hat p}^2/T^2=1/(2NT^2)\) 的 \(\Delta\) 估计。时间 \(\tau_{\mathrm a}\) 内有 \(\tau_{\mathrm a}/(2T_{\mathrm c})\) 个这样的估计，平均后方差为 \(\frac{1}{2NT^2}\cdot\frac{2T_{\mathrm c}}{\tau_{\mathrm a}}=\frac{T_{\mathrm c}}{NT^2\tau_{\mathrm a}}\)，与“每个循环一个 \(\sigma_\Delta=1/(T\sqrt N)\) 的独立测量”的算法结果相同。

每次循环用时 \(T_{\mathrm c}\)（含原子制备、飞行与探测，\(T_{\mathrm c}>T\)），各次的投影噪声相互独立，是白频率噪声。由式 (3.24)： \begin{align*} \sigma_y(\tau_{\mathrm a})&=\frac{1}{\omega_0T\sqrt N}\sqrt{\frac{T_{\mathrm c}}{\tau_{\mathrm a}}}\\ &=\frac{1}{\pi Q}\sqrt{\frac{T_{\mathrm c}}{N\tau_{\mathrm a}}},\qquad Q=2\nu_0T .\tag{4.6} \end{align*} 第二个等号用了 \(\omega_0T=2\pi\nu_0T=\pi Q\)。这正是文献中量子投影噪声限制下的稳定度公式[32]。

令 \(\vect r=\vect r_0+\delta\vect r\)，\(b=b_0+\delta b\)。几何距离对 \(\vect r\) 的梯度是 \(\nabla_{\vect r}|\vect r-\vect r_i|=(\vect r-\vect r_i)/|\vect r-\vect r_i|\)，即从卫星指向接收机的单位矢量。记由 \(\vect r_0\) 指向卫星 \(i\) 的单位矢量为 \(\uvec e_i=(\vect r_i-\vect r_0)/d_i^{(0)}\)，\(d_i^{(0)}=|\vect r_0-\vect r_i|\)。一阶泰勒展开： \[|\vect r-\vect r_i|=d_i^{(0)}-\uvec e_i\cdot\delta\vect r+\order{\frac{|\delta\vect r|^2}{d_i^{(0)}}} .\] （物理意义：接收机朝卫星移动 \(\delta\vect r\)，距离缩短 \(\uvec e_i\cdot\delta\vect r\)。）定义残差 \(l_i\equiv\rho_i-d_i^{(0)}-cb_0\)，代入式 (5.2)： \[\tag{5.3} l_i=-\uvec e_i\cdot\delta\vect r+c\,\delta b+\varepsilon_i .\] 略去的二阶项约 \(|\delta\vect r|^2/(2d_i)\)：初值偏 10 km 时约 2.5 m，偏 1 km 时约 2.5 cm。注意 \(b\) 在式 (5.2) 中本来就是线性的，所以 \(b_0\) 可以取任意值（例如 0）。

\(n>4\) 时方程多于未知量，又有噪声，一般没有精确解。取使残差平方和最小的 \(\vect u\)： \[J(\vect u)=\sum_{i=1}^n\Big(l_i-\sum_jG_{ij}u_j\Big)^2 .\] 对每个 \(u_k\) 求偏导并令其为零： \begin{gather*} \pdv{J}{u_k}=-2\sum_iG_{ik}\Big(l_i-\sum_jG_{ij}u_j\Big)=0\\ \Longrightarrow\quad \sum_j\Big(\sum_iG_{ik}G_{ij}\Big)u_j=\sum_iG_{ik}l_i . \end{gather*} 左边括号里的 \(\sum_iG_{ik}G_{ij}\) 就是矩阵 \(\mathsf G^{\mathsf T}\mathsf G\) 的 \((k,j)\) 元，右边是 \(\mathsf G^{\mathsf T}\vect l\) 的第 \(k\) 元。写成矩阵形式，即正规方程： \begin{gather*} (\mathsf G^{\mathsf T}\mathsf G)\,\vect u=\mathsf G^{\mathsf T}\vect l\qquad\Longrightarrow \\ \hat{\vect u}=\mathsf C\,\mathsf G^{\mathsf T}\vect l,\qquad \mathsf C\equiv(\mathsf G^{\mathsf T}\mathsf G)^{-1}. \tag{5.4} \end{gather*} \(\mathsf G^{\mathsf T}\mathsf G\) 是 \(4\times4\) 的对称方阵。

设各伪距的误差 \(\varepsilon_i\) 相互独立、均值为零、标准差都是 \(\sigma_\rho\)。\(\hat{\vect u}\) 是 \(\vect l\) 的线性函数，其误差为 \(\delta\hat{\vect u}=\mathsf C\mathsf G^{\mathsf T}\vect\varepsilon\)。由“独立量之和的方差相加”，可证（附录 D） \[\tag{5.5} \avg{\delta\hat u_j\,\delta\hat u_k}=\sigma_\rho^2\,C_{jk} .\] 对角元 \(\sigma_\rho^2C_{jj}\) 是第 \(j\) 个未知量的方差，非对角元描述两个未知量的误差是否相关。

由式 (3.11)，卫星钟相对坐标时： \[\frac{\dd\tau_{\mathrm{sat}}}{\dd t}=1-\frac{\mu}{ac^2}-\frac{\mu}{2ac^2}=1-\frac{3\mu}{2ac^2}.\] 先用最简单的地面钟模型：静止在不转动的地球表面 \(r=R_{\mathrm E}\) 处，\(\dd\tau_{\mathrm{gnd}}/\dd t=1-\mu/(R_{\mathrm E}c^2)\)。两者相除并保留到 \(1/c^2\) 阶（\((1+\epsilon_1)/(1+\epsilon_2)\approx1+\epsilon_1-\epsilon_2\)）： \begin{align*} y_{\mathrm{sat}}\equiv\frac{\dd\tau_{\mathrm{sat}}}{\dd\tau_{\mathrm{gnd}}}-1=\;&\underbrace{\frac{\mu}{c^2}\Big(\frac{1}{R_{\mathrm E}}-\frac1a\Big)}_{\text{引力（广义相对论）项}}\\ &\underbrace{-\;\frac{\mu}{2ac^2}}_{\text{速度（狭义相对论）项}} .\tag{6.1} \end{align*}

令式 (6.1) 为零：\(\frac{1}{R_{\mathrm E}}-\frac1a=\frac{1}{2a}\)，即 \[\tag{6.2} a=\frac32R_{\mathrm E}\approx{9567\,}\mathrm{km},\] 即轨道高度约 3190 km。用完整模型 \(L_{\mathrm G}=3\mu/(2ac^2)\) 得 \(a={9546\,}\mathrm{km}=1.497R_{\mathrm E}\)[24]。低于此高度（如 400 km 高的低轨卫星），速度项占优，卫星钟变慢；高于此高度，引力项占优，卫星钟变快。\(a\to\infty\) 时 \(y_{\mathrm{sat}}\to\mu/(R_{\mathrm E}c^2)\)：远离地球的静止钟比地面钟快 \(7\times10^{-10}\)。

把式 (6.5) 代入式 (3.11)： \begin{align*} \frac{\dd\tau}{\dd t}&=1-\frac{\mu}{rc^2}-\frac{\mu}{2c^2}\Big(\frac2r-\frac1a\Big)\\ &=1-\frac{3\mu}{2ac^2}-\frac{2\mu}{c^2}\Big(\frac1r-\frac1a\Big). \end{align*} 第二项是常值（与圆轨道相同），第三项随轨道位置变化。对时间积分。由开普勒方程，\(n\,\dd t=(1-e\cos E)\,\dd E\)，即 \(\dd t=\frac{r}{an}\dd E\)（用了 \(r=a(1-e\cos E)\)）。于是 \begin{align*} \Big(\frac1r-\frac1a\Big)\dd t&=\Big(\frac1r-\frac1a\Big)\frac{r}{an}\dd E\\ &=\frac{1}{an}\Big(1-\frac ra\Big)\dd E=\frac{e\cos E}{an}\,\dd E . \end{align*} 从近地点积分：\(\int_0^t(\frac1r-\frac1a)\dd t=\frac{e\sin E}{an}\)。又 \(an=\sqrt{\mu/a}\)，所以 \begin{gather*} \tau(t)=\Big(1-\frac{3\mu}{2ac^2}\Big)t+\Delta t_{\mathrm r},\\ \Delta t_{\mathrm r}=-\frac{2\sqrt{\mu a}}{c^2}\,e\sin E .\tag{6.6} \end{gather*}

接收机计算位置时用的是随地球转动的“地固坐标系”（ECEF）：卫星星历给出的是发射时刻卫星的地固坐标，接收机自己在地固系中静止。但光速不变只在惯性系中成立（Ashby 强调，若在转动系中直接假定光速为 \(c\) 将导致严重错误[24]），所以传播时间必须在 ECI 系中计算。

取 ECI 系与地固系在发射时刻 \(t_{\mathrm e}\) 重合。发射时卫星在 \(\vect r_{\mathrm S}\)，接收机在 \(\vect r_{\mathrm R}\)，记 \(\vect D=\vect r_{\mathrm R}-\vect r_{\mathrm S}\)，\(D=|\vect D|\)。在 ECI 系中接收机以 \(\vect v_{\mathrm R}=\vect\omega_{\mathrm E}\times\vect r_{\mathrm R}\) 运动，接收时刻 \(t_{\mathrm r}=t_{\mathrm e}+\Delta t\) 位于 \(\vect r_{\mathrm R}+\vect v_{\mathrm R}\Delta t\)（略去 \(\order{\omega_{\mathrm E}^2\Delta t^2}\) 的向心项）。光速不变给出 \[c\Delta t=|\vect D+\vect v_{\mathrm R}\Delta t|=D+\uvec D\cdot\vect v_{\mathrm R}\Delta t+\order{\frac{v_{\mathrm R}^2\Delta t^2}{D}} .\] 右边第二项已是小量，其中的 \(\Delta t\) 可用零阶值 \(D/c\) 代替： \[c\Delta t=D+\frac{\vect D\cdot\vect v_{\mathrm R}}{c} .\] 用矢量混合积的轮换性 \(\vect D\cdot(\vect\omega_{\mathrm E}\times\vect r_{\mathrm R})=\vect\omega_{\mathrm E}\cdot(\vect r_{\mathrm R}\times\vect D)\)，以及 \(\vect r_{\mathrm R}\times\vect D=\vect r_{\mathrm R}\times(\vect r_{\mathrm R}-\vect r_{\mathrm S})=\vect r_{\mathrm S}\times\vect r_{\mathrm R}\)： \begin{gather*} \Delta t=\frac Dc+\Delta t_{\mathrm{S}}, \\ \Delta t_{\mathrm S}=\frac{\vect\omega_{\mathrm E}\cdot(\vect r_{\mathrm S}\times\vect r_{\mathrm R})}{c^2}=\frac{\omega_{\mathrm E}(x_{\mathrm S}y_{\mathrm R}-y_{\mathrm S}x_{\mathrm R})}{c^2} \\ =\frac{2\omega_{\mathrm E}A_z}{c^2}, \tag{6.7} \end{gather*} \(A_z=\frac12(\vect r_{\mathrm S}\times\vect r_{\mathrm R})_z\) 是三角形 O–S–R 在赤道面上投影的带符号面积（从 S 到 R 为逆时针，即信号向东传播时为正）。

\(\frac{\dd}{\dd t}|c_{\mathrm e}|^2=2\Re(c_{\mathrm e}^*\dot c_{\mathrm e})\)。由式 (A.1)，\(c_{\mathrm e}^*\dot c_{\mathrm e}=-\frac{\ii}{2}(-\Delta|c_{\mathrm e}|^2+\Omega_{\mathrm R}c_{\mathrm e}^*c_{\mathrm g})\)。第一项是纯虚数，实部为零；第二项的实部为 \(\frac{\Omega_{\mathrm R}}{2}\Im(c_{\mathrm e}^*c_{\mathrm g})\)（用了 \(\Re(-\ii z)=\Im z\)）。所以 \(\frac{\dd}{\dd t}|c_{\mathrm e}|^2=\Omega_{\mathrm R}\Im(c_{\mathrm e}^*c_{\mathrm g})=\frac{\Omega_{\mathrm R}}{2}S_y\)。由归一化，\(\frac{\dd}{\dd t}|c_{\mathrm g}|^2=-\frac{\Omega_{\mathrm R}}{2}S_y\)。相减： \[\dot S_z=\Omega_{\mathrm R}S_y .\]

\(\frac{\dd}{\dd t}(c_{\mathrm e}^*c_{\mathrm g})=\dot c_{\mathrm e}^*c_{\mathrm g}+c_{\mathrm e}^*\dot c_{\mathrm g}\)。由式 (A.1) 及其复共轭 \(\dot c_{\mathrm e}^*=\frac{\ii}{2}(-\Delta c_{\mathrm e}^*+\Omega_{\mathrm R}c_{\mathrm g}^*)\)： \begin{align*} \frac{\dd}{\dd t}(c_{\mathrm e}^*c_{\mathrm g})&=\frac{\ii}{2}\big(-\Delta c_{\mathrm e}^*c_{\mathrm g}+\Omega_{\mathrm R}|c_{\mathrm g}|^2\big)-\frac{\ii}{2}\big(\Omega_{\mathrm R}|c_{\mathrm e}|^2+\Delta c_{\mathrm e}^*c_{\mathrm g}\big)\\ &=-\ii\Delta\,c_{\mathrm e}^*c_{\mathrm g}-\frac{\ii\Omega_{\mathrm R}}{2}S_z . \end{align*} 记 \(c_{\mathrm e}^*c_{\mathrm g}=\frac12(S_x+\ii S_y)\)，取实部与虚部： \[\dot S_x=\Delta S_y,\qquad \dot S_y=-\Delta S_x-\Omega_{\mathrm R}S_z .\]

椭圆的半长轴 \(a\)、半短轴 \(b_{\mathrm e}=a\sqrt{1-e^2}\)，焦点（地心）在距中心 \(ae\) 处。以椭圆中心为原点、长轴为 \(x\) 轴，椭圆上的点可写成 \(x=a\cos E\)，\(y=b_{\mathrm e}\sin E\)——这就是偏近点角 \(E\) 的定义（\((a\cos E,a\sin E)\) 是半径 \(a\) 的辅助圆上横坐标相同的点）。到焦点 \((ae,0)\) 的距离： \begin{align*} r^2&=(a\cos E-ae)^2+b_{\mathrm e}^2\sin^2E\\ &=a^2\big[\cos^2E-2e\cos E+e^2+(1-e^2)(1-\cos^2E)\big]\\ &=a^2(1-e\cos E)^2 . \end{align*}

由开普勒第二定律，焦点到卫星的连线单位时间扫过的面积为常数 \(\pi ab_{\mathrm e}/P\)（\(P\) 为周期）。从近地点到 \(E\) 扫过的面积：先在辅助圆上算“焦点—近地点—投影点”围成的面积，等于圆心扇形 \(\frac12a^2E\) 减去三角形（圆心、焦点、投影点）\(\frac12\cdot ae\cdot a\sin E\)；椭圆是辅助圆在 \(y\) 方向压缩 \(b_{\mathrm e}/a\) 倍，面积也按此比例缩小： \[\mathcal{A}(E)=\frac{b_{\mathrm e}}{a}\cdot\frac{a^2}{2}(E-e\sin E)=\frac{ab_{\mathrm e}}{2}(E-e\sin E).\] 令 \(\mathcal A=\pi ab_{\mathrm e}\,t/P\)，得 \(E-e\sin E=\frac{2\pi}{P}t=nt\)；由开普勒第三定律 \(P^2=4\pi^2a^3/\mu\)，\(n=\sqrt{\mu/a^3}\)。

在近地点 \(r_{\mathrm p}=a(1-e)\) 与远地点 \(r_{\mathrm a}=a(1+e)\) 速度都垂直于矢径，角动量守恒给出 \(v_{\mathrm p}r_{\mathrm p}=v_{\mathrm a}r_{\mathrm a}\)；能量守恒 \(\frac12v_{\mathrm p}^2-\frac{\mu}{r_{\mathrm p}}=\frac12v_{\mathrm a}^2-\frac{\mu}{r_{\mathrm a}}\)。消去 \(v_{\mathrm a}=v_{\mathrm p}r_{\mathrm p}/r_{\mathrm a}\)：\(\frac12v_{\mathrm p}^2(1-\frac{r_{\mathrm p}^2}{r_{\mathrm a}^2})=\mu\frac{r_{\mathrm a}-r_{\mathrm p}}{r_{\mathrm p}r_{\mathrm a}}\)，得 \(v_{\mathrm p}^2=\frac{2\mu r_{\mathrm a}}{r_{\mathrm p}(r_{\mathrm p}+r_{\mathrm a})}\)。比能量 \(\mathcal E=\frac12v_{\mathrm p}^2-\frac{\mu}{r_{\mathrm p}}=\frac{\mu}{r_{\mathrm p}}\Big(\frac{r_{\mathrm a}}{r_{\mathrm p}+r_{\mathrm a}}-1\Big)=-\frac{\mu}{r_{\mathrm p}+r_{\mathrm a}}=-\frac{\mu}{2a}\)。于是任一点 \(v^2=2(\mathcal E+\mu/r)=\mu(2/r-1/a)\)，即式 (6.5)。

由式 (5.4)，\(\hat u_j=\sum_kM_{jk}l_k\)，\(\mathsf M=\mathsf C\mathsf G^{\mathsf T}\)。把 \(\vect l=\mathsf G\vect u_{\text{真}}+\vect\varepsilon\) 代入，\(\mathsf M\mathsf G=\mathsf C\mathsf G^{\mathsf T}\mathsf G=\mathsf I\)，所以 \(\delta\hat u_j=\sum_kM_{jk}\varepsilon_k\)。于是 \begin{align*} \avg{\delta\hat u_j\delta\hat u_l}&=\sum_{k,m}M_{jk}M_{lm}\avg{\varepsilon_k\varepsilon_m}\\ &=\sigma_\rho^2\sum_kM_{jk}M_{lk}=\sigma_\rho^2(\mathsf M\mathsf M^{\mathsf T})_{jl}, \end{align*} 用了 \(\avg{\varepsilon_k\varepsilon_m}=\sigma_\rho^2\)（\(k=m\)）或 0（\(k\ne m\)）。而 \(\mathsf M\mathsf M^{\mathsf T}=\mathsf C\mathsf G^{\mathsf T}\mathsf G\mathsf C^{\mathsf T}=\mathsf C\)（\(\mathsf C\) 对称）。故 \(\avg{\delta\hat u_j\delta\hat u_l}=\sigma_\rho^2C_{jl}\)。\(\sigma_\rho\) 未知时用残差估计 \(s_0^2=\sum v_i^2/(n-m)\)，\(m\) 为未知量个数；除以 \(n-m\) 而不是 \(n\)，是因为拟合用掉了 \(m\) 个自由度，残差平方和的期望为 \((n-m)\sigma_\rho^2\)。

数据 \((h_i,y_i)\)，\(y_i\) 的标准不确定度 \(\sigma_i\)，模型 \(y=y_0+kh\)。最小化 \(\chi^2=\sum w_i(y_i-y_0-kh_i)^2\)，\(w_i=1/\sigma_i^2\)。对 \(y_0\)、\(k\) 求偏导并令其为零： \[y_0\Sigma_0+k\Sigma_h=\Sigma_y,\qquad y_0\Sigma_h+k\Sigma_{hh}=\Sigma_{hy},\] 其中 \(\Sigma_0=\sum w_i\)，\(\Sigma_h=\sum w_ih_i\)，\(\Sigma_{hh}=\sum w_ih_i^2\)，\(\Sigma_y=\sum w_iy_i\)，\(\Sigma_{hy}=\sum w_ih_iy_i\)。解得（克莱姆法则，\(D_w=\Sigma_0\Sigma_{hh}-\Sigma_h^2\)） \[k=\frac{\Sigma_0\Sigma_{hy}-\Sigma_h\Sigma_y}{D_w},\qquad y_0=\frac{\Sigma_{hh}\Sigma_y-\Sigma_h\Sigma_{hy}}{D_w}.\] 这是本附录第一框在 \(\mathsf G\) 第 \(i\) 行为 \((1,h_i)/\sigma_i\)、\(l_i=y_i/\sigma_i\) 时的特例，\(\mathsf G^{\mathsf T}\mathsf G=\begin{pmatrix}\Sigma_0&\Sigma_h\\\Sigma_h&\Sigma_{hh}\end{pmatrix}\)，\(\sigma_\rho=1\)。其逆给出 \(\sigma_{y_0}^2=\Sigma_{hh}/D_w\)，\(\sigma_k^2=\Sigma_0/D_w\)，\(\mathrm{cov}(y_0,k)=-\Sigma_h/D_w\)。若模型正确且 \(\sigma_i\) 估计正确，\(\chi^2\) 的期望为 \(n-2\)。