Comparison of Map-Making Algorithms for CMB Experiments — 故事版¶
arXiv: astro-ph/0501504 ｜ 作者: Poutanen et al. ｜ 年份: 2006
从一条数据流到一张天图¶
想象你坐在 Planck 卫星上。卫星每分钟自转一圈，望远镜沿着天空画出近乎大圆的扫描环。探测器每秒采样约 108 次（\(f_s = 108.3\) Hz），一年下来产出约 34 亿个时间序列样本。这条数据流——时间有序数据（TOD）——就是一切的起点。
你的任务是从这条 TOD 中提取出一张 CMB 天图。
问题在于 TOD 不只有 CMB 信号。仪器噪声——尤其是低频的 \(1/f\) 噪声——像慢慢漂移的水位一样叠加在信号上。如果你只是简单地把落在同一个像素里的样本取平均（naive binning），这些低频漂移就会留在图中，变成条纹状的伪影。 [补充]
制图算法的核心任务就是：在 bin 图之前，尽可能干净地去掉相关噪声。
两条路线：ML 与 Destriping¶
最大似然（ML）路线¶
ML 的思路是数学上最直接的：写下数据的似然函数
然后对天图 \(\mathbf{m}\) 求极值。这里 \(\mathbf{y}\) 是 TOD，\(\mathbf{P}\) 是指向矩阵（每一行告诉你这个样本对应天空的哪个像素），\(\mathbf{N}\) 是噪声协方差矩阵。极值条件给出正规方程
这是一个 \(N_{\text{pix}}\) 维的线性方程组。对 Planck 来说 \(N_{\text{pix}} \sim 10^6 \text{–} 10^8\)，直接求逆不现实。ROMA 和 MapCUMBA 都用预条件共轭梯度法迭代求解，关键技巧是利用噪声协方差矩阵近似为循环矩阵，在频域中变成对角矩阵，从而高效计算 \(\mathbf{N}^{-1}\mathbf{y}\)。 [原文]
代价：需要知道噪声功率谱（即 \(\mathbf{N}\) 的结构），且计算量巨大——本文中每张图需要 192–256 个 CPU 并行计算约 10 分钟。 [原文]
Destriping 路线¶
Destriping 利用了 Planck 扫描策略的特殊性：卫星在每次重指向（repointing）之间的 60 圈扫描画的是同一条环，可以 coadd 成一条"ring"。\(1/f\) 噪声在一个 ring 上的效果可以近似为一个均匀偏移（baseline）。 [原文]
于是问题变成：找到每个 ring 的 baseline 振幅 \(\mathbf{a}\)，从 TOD 中减掉，再 bin 成图。
怎么找 \(\mathbf{a}\)？利用不同 ring 之间的交叉点——同一个像素被不同 ring 扫过时应该给出相同的信号。如果不同 ring 在同一像素处的读数不同，差别就来自它们各自的 baseline。这就是 destriping 方程
的物理含义。矩阵 \(\mathbf{Z}\) 的作用是"对每个样本减去同像素内的平均值"——这恰好消除了信号，只留下需要拟合的 baseline。 [原文]
优势：不需要噪声功率谱，方程维度只有 ring 的数量（约几千），单 CPU 7 分钟就能出一张图。 [原文]
比赛开始：谁的图更好？¶
作者用 Level S 软件仿真了 Planck LFI 100 GHz 一个探测器的一年数据，分 4 种情况（不同膝频、波束形状、有无前景），三个代码分别出图。
回合一：图域比较¶
先看直觉最简单的指标——输出图和"理想图"（binned noiseless map）之间的偏差，即重建误差（reconstruction error）\(\boldsymbol{\varepsilon}\)。
| 指标 | ML（ROMA ≈ MapCUMBA） | Destriping |
|---|---|---|
| 误差图 std | 138.25 μK | 138.46 μK |
| 白噪声下限 | 137.22 μK | 137.22 μK |
ML 的噪声残留更低——意料之中，因为它使用了完整的噪声协方差信息。但差距很小，只有约 0.2 μK。 [原文]
反转来了。作者把误差分解为噪声分量 \(\boldsymbol{\varepsilon}_n\) 和信号分量 \(\boldsymbol{\varepsilon}_p\)：
| 指标 | ML | Destriping |
|---|---|---|
| 信号误差 \(\boldsymbol{\varepsilon}_p\) std（Case 1） | 0.88 μK | 0.28 μK |
ML 的信号误差反而是 destriping 的 3 倍！ 而且 ML 的信号误差在前景辐射梯度最强的地方（银道面附近）尤其大。 [原文]
信号误差从何而来？¶
根源是像素化噪声（pixelisation noise）。信号在一个像素内并不均匀——同一个像素中不同样本的信号值略有不同（亚像素结构）。把它们平均到同一个像素就会引入误差。
关键在于两种算法对像素化噪声的"响应"不同：
- ML：噪声滤波器 \(\mathbf{N}^{-1}\) 会作用于整条 TOD，包括信号部分。像素化噪声的频谱几乎是平坦的（白噪声型），在膝频 \(f_k\) 以下与真正的 \(1/f\) 噪声长得一样——ML 无法区分两者，把像素化噪声也当成 \(1/f\) 噪声来"校正"，反而引入了信号误差。 [原文]
- Destriping：只拟合每个 ring 的均匀 baseline（零频模式），像素化噪声中只有零频分量会影响结果。高次模式全部被忽略，反而避免了 ML 的放大效应。 [原文]
形象地说：ML 拿着一把精细的滤波器去清理噪声，但这把滤波器分不清"真正的 \(1/f\) 噪声"和"像素化噪声"——两者在膝频以下的频段看起来一样平坦，结果 ML 把像素化噪声也当成 \(1/f\) 噪声来"纠正"，反而弄巧成拙。destriping 用的是粗粒度的 baseline 模型，粗到只拟合零频模式，像素化噪声中只有零频分量会影响结果，其他频率分量全部被忽略——反而更安全。 [补充]
比赛继续：谱域比较¶
图域的差异会怎样传递到角功率谱 \(C_\ell\) 估计？
作者用 MASTER 方法估计功率谱，关键公式是
其中 \(F_\ell\) 是滤波函数（制图对信号的畸变），\(\langle \widetilde{N}_\ell \rangle\) 是噪声偏差（用 100 次 MC 估计）。
- Destriping 的 \(F_\ell \approx 1\)——制图几乎不畸变信号。 [原文]
- ML 的 \(F_\ell\) 在 \(\ell = 800 \sim 1000\) 处偏离最大，约 0.6%。 [原文]
如果忽略 \(F_\ell\)（设为 1），ML 的功率谱估计会有一个约 0.6% 的偏差。对 Planck 级精度的实验，这不可忽略。不过，\(F_\ell\) 可以通过纯信号 MC 模拟来估计并校正——校正后两种方法的功率谱估计趋于一致。 [原文]
误差条方面，destriping 比 ML 大约 5%（在 \(\ell \lesssim 1000\)），主要来自更高的噪声残留。在低 \(\ell\) 端被宇宙方差（cosmic variance）主导后，两者差异消失。 [原文]
计算成本：不可忽略的现实考量¶
| ML（ROMA / MapCUMBA） | Destriping | |
|---|---|---|
| 每张图 | ~10 min, 192–256 CPU 并行 | ~7 min, 单 CPU |
| MC 100 次 | ~1000 min × 200 CPU | ~400 min × 1 CPU |
在需要大量 MC 模拟来估计噪声偏差和误差条的现实中，destriping 的计算优势是压倒性的。 [原文]
最终图景¶
这篇文章的核心发现可以用一句话概括：
ML 制图赢在噪声更低，但输在对信号的畸变更大。两者在功率谱层面的差异可以校正，校正后差距很小（误差条约 5%）。考虑到 destriping 的计算优势，它是 Planck 级 CMB 实验中一个极具竞争力的选择。
这个结论打破了"ML 全面最优"的直觉，揭示了制图问题中噪声与信号误差的微妙权衡——在大数据时代的 CMB 实验中，"最优"不只是统计意义上的，还必须考虑计算可行性和系统性偏差。 [补充]