脏波束与不唯一的图像

把测得的波朝每个像素加回去

令 $N$ 为独立复样本数量， $y_k$ 为标定后的可见度， $(l,m)$ 为试探的天空位置。给每条信号乘回相反相位，再对实部取平均。试探位置碰到一个亮点时，各条贡献容易对齐；其他位置则部分抵消。这样得到脏图 $D$ 。

把所有可见度都换成一，再做同样的和，就得到脏波束。它是原点处单位点源的响应。把前向模型代入脏图公式，可以看到：每个真实亮点都变成一份移位后的波束。采样之间的空隙，便化为波束中的旁瓣，其中也有负旁瓣。

D(l,m)=\frac1N\sum_k\operatorname{Re}\!\left[y_k e^{2\pi i(u_kl+v_km)}\right],\qquad B(l,m)=\frac1N\sum_k\cos[2\pi(u_kl+v_km)]

采用均匀权重，除以样本数后，原点的波束值等于一。若把共轭样本也明确加入，权重随样本翻倍，得到同一个实结果。

D(l,m)=\sum_p x_p B(l-l_p,m-m_p)

每个源像素的通量，乘一份平移后的脏波束，再全部相加。这说明脏图为什么带着阵列自己的模糊图案。

对应阅读：NRAO Essential Radio Astronomy, Chapter 3, aperture synthesis · EHT Collaboration, First M87 Results IV, arXiv:1906.11241v1

每波束与每像素，数值含义不同

中心 $0.60\,\mathrm{Jy}$ 点源的脏图峰值是 $0.60\,\mathrm{Jy}$ 每单位峰值波束。把脏图所有像素直接相加，通常得不回 $0.60\,\mathrm{Jy}$ ：旁瓣和显示窗口都会影响这个和。重建图则为每个像素保存一份通量，像素通量之和被限制为独立测得的总通量。

默认网格边长为偶数，原点位于四个像素之间。连续波束满足 $B(0,0)=1$ ，采样在像素中心的峰值却可能略低。每幅显示图各用自己的峰值拉伸颜色，并打印该峰值。比较时要同时看形状和单位，相同颜色未必代表相同强度。

亲手构造两幅信号完全相同的天空

假设所有基线都满足 $v=0$ 。把一个点向北移动，东西坐标和通量不变，前向求和里的每个指数都不变。两个北坐标不同的正亮度图像，便给出完全相同的无噪声观测。继续加长东西基线，也分不开这一次移动。

默认 $24\times24$ 图像含 $576$ 个实像素通量。五十四条复观测最多给出 $108$ 条实线性方程；独立总通量再加一条。非负条件能缩小候选范围，但稀疏采样仍留下很多方向，沿着它们改变像素，信号变化很小，甚至完全不变。

缺少短基线时，延展辐射尤其容易漏掉。一个宽而平滑的成分可以带着很多总通量，却在长基线频率上贡献很少。补上总通量有帮助，重建仍要把它分配到那些尚未分清的低频结构里。

对应阅读：NRAO Essential Radio Astronomy, Chapter 3, aperture synthesis · EHT Collaboration, First M87 Results IV, arXiv:1906.11241v1

明确写出两种重建各自喜欢什么

把离散傅里叶算子记为 $A$ ，像素通量向量记为 $x$ 。最小化平均复信号残差平方，再加正则项 $R$ ，同时要求像素非负、通量和等于独立测量。正则项就是观测留有选择时，程序用来选图像的偏好。

平滑选项惩罚上下左右相邻像素之差，每对邻居只数一次。集中选项加入半径加权通量，再加一小项像素平方。 $\rho_p$ 为像素到中心的距离，以整个视场宽度为单位。下面的系数与浏览器完全相同，通量统一用央斯基。

J(x)=\frac1{2N}\sum_k|(Ax)_k-y_k|^2+R(x),\qquad x_p\geq0,\quad\sum_px_p=F_{\mathrm{total}}

数据项单位为央斯基平方。保存的合成信号没有加入随机噪声；每条复样本的五毫央斯基只作参考不确定度，本式没有再用它加权。

R_{\mathrm{smooth}}=\frac{0.05}{2}\sum_{\langle p,q\rangle}(x_p-x_q)^2

邻居差让快速的像素起伏代价更高。即使更细的结构也能预测观测，程序仍可能偏好较宽的结构。

R_{\mathrm{compact}}=(0.004\,\mathrm{Jy})\sum_p\rho_p^2x_p+\frac{0.03}{2}\sum_px_p^2

第一项偏好靠近中心的辐射，第二项不鼓励把通量全塞进极少数像素。它们与邻居平滑是不同选择。

对应阅读：EHT Collaboration, First M87 Results IV, arXiv:1906.11241v1

推导同时保留非负和总通量的一步

对实像素微扰求导，复残差平方给出 $\operatorname{Re}[A^*(Ax-y)]/N$ 。加上正则项梯度，先沿下坡方向走一步，再投影回非负且总通量固定的单纯形。投影的做法是减去同一个阈值，把负值截到零，并调整阈值，让剩余通量恰好为 $F_{\mathrm{total}}$ 。

为什么每个像素减同一个阈值？把试探向量与允许图像的平方距离最小化，再为通量和加入一个拉格朗日乘子。对未碰到零边界的像素求导，得到同一个平移量；本来会变负的像素停在零。将试探值排序，就能找出哪些像素仍为正，并算出阈值。

程序做九十步加速投影梯度。回溯检查局部二次模型是否盖住试探目标值；若完整目标上升，便重启加速动量。这是有限次计算。除了看图片，也要看目标历史和预测信号残差。增加迭代次数与改变正则项，回答的是不同问题。

x^{\mathrm{next}}=\Pi\!\left(x-\eta\nabla J(x)\right),\qquad [\Pi(z)]_p=\max(0,z_p-\tau),\quad \sum_p\max(0,z_p-\tau)=F_{\mathrm{total}}

这个投影把两项条件同时满足。先截负值，再把正值按比例归一化，通常会得到另一个结果。

先读固定的默认结果，再改条件

固定默认值含五十四条样本、三小时观测、 $0.600\,\mathrm{Jy}$ 总通量，以及有限像素双分量源。九十步以后，平滑重建的复信号均方根残差为 $7.041\,\mathrm{mJy}$ ，集中偏好为 $5.934\,\mathrm{mJy}$ 。这里对复数模长平方取平均再开方，衡量的是信号误差，不是图像每个像素的误差。

无噪声数据仍有非零残差，来自有限迭代和明确加入的偏好，不是估计出的天体噪声。两种反演都只读取已标定信号、网格和独立总通量；源参考只送到单独标记的参考图。导出观测后，即使不再运行光源生成器，也能重算这两幅图。

\mathrm{RMS}_{V}=\sqrt{\frac1N\sum_k|(Ax)_k-y_k|^2}

信号残差较小，不能单独保证每个图像细节正确。还要换偏好比较，并在相同采样下检验已知模拟源能恢复到什么程度。

像素网格与视场	$24\times24,\quad120\,\mu\mathrm{as}$
中央格林尼治时角	$H=105^{\circ}$
总通量	$F=0.600\,\mathrm{Jy}$
独立复样本数	54
第一组有向闭合相位	$0.287631^{circ}$
平滑偏好信号均方根残差	$7.0408,mathrm{mJy}$
集中偏好信号均方根残差	$5.9340,mathrm{mJy}$