从 Suzaku 指数盘模型到 Figure 3 的方向性条件点 — From the Suzaku exponential-disk model to a directional conditional point on Figure 3

面向物理系本科一年级 · 可执行教程

Author

M31 CGM Team

Published

July 18, 2026

从 Suzaku 指数盘模型到 Figure 3 的方向性条件点

教程目标：理解 Ueda et al. (2022) 用 Suzaku 拟合银河系软 X 射线背景得到的指数电子密度盘模型，如何沿 M31 十四个 XMM 视线投影、吸收转换，并浓缩成 Figure 3 上的一个方向性条件模型点。 Tutorial goal: Understand how the exponential electron-density disk model fitted by Ueda et al. (2022) to the Suzaku soft X-ray background is projected along the fourteen M31 XMM sightlines, absorption-converted, and condensed into one directional conditional model point on Figure 3.

目标读者：物理系本科一年级，已学完普通物理（电磁学/光学），了解基本的原子物理概念（能级、跃迁）和简单的微积分，但不需要天文观测经验。 Target audience: First-year physics undergraduates who have completed general physics (electromagnetism/optics), basic atomic physics, and elementary calculus, without requiring astronomical observing experience.

核心问题：Ueda+22 的盘模型是一个三维密度公式——给定银河系内任意一点 \((R, z)\)，就能算出该处的电子密度 \(n_e\)。但 Figure 3 需要的是沿某条视线方向的 X 射线表面亮度。如何把一个”密度公式”变成”亮度值”？又如何把 14 条视线的预测浓缩成 Figure 3 上的一个点？ Core question: The Ueda+22 disk model is a 3D density formula — given any point \((R, z)\) in the Galaxy, it returns the electron density \(n_e\) there. Figure 3 needs an X-ray surface brightness along a sightline. How do we turn a “density formula” into a “brightness value”? And how do we condense 14 sightline predictions into one point on Figure 3?

0. 准备工作：环境与数据

本教程只需要 Python 标准科学栈。所有数据文件已随教程提供。

Code 
import numpy as np
import pandas as pd
import matplotlib.pyplot as plt
from pathlib import Path

plt.style.use('seaborn-v0_8-whitegrid')
plt.rcParams['font.size'] = 12
plt.rcParams['figure.dpi'] = 150
plt.rcParams['savefig.dpi'] = 300
# NOTE: do NOT add a CJK font fallback here — see SKILL.md
# All matplotlib text labels stay English-only.

# The Ueda disk CSV lives directly under assets/ (symlinked to the paper draft directory).
# The primary measurements CSV (for signed projected radii) lives next to the .qmd file.
DATA = Path("assets")
PRIMARY = Path("m31_cgmsum_v19_primary_measurements_public.csv")
print("Environment ready!")
Environment ready!

1. 背景知识：Suzaku 看到的银河系热气体

1.1 Suzaku 与软 X 射线背景

Suzaku 是日本宇宙航空研究开发机构（JAXA）2005 年发射的 X 射线天文卫星，搭载 X 射线成像光谱仪（XIS）。它的低地球轨道本底极低，特别适合观测弥漫的软 X 射线背景（约 0.2–2 keV）。

Ueda et al. (2022) 用 Suzaku 观测了约 130 个场（parent sample），覆盖银河系 \(75^\circ < l < 285^\circ\)、\(|b| > 15^\circ\) 的大片高银纬天空。他们对每个场的软 X 射线背景谱进行拟合，分离出约 0.22 keV 的银河系晕（Milky-Way halo, MWH）成分——这就是我们要研究的银河系热气体。

1.2 为什么用”盘模型”描述热气体？

银河系的热晕气体并不均匀——离银盘越远，气体越稀薄。Ueda+22 用一个简单的指数盘模型拟合观测到的发射量分布：

\[ n_e(R, z) = n_{e0} \, \exp\!\left(-\frac{R - R_\odot}{R_0}\right) \exp\!\left(-\frac{|z|}{z_0}\right) \]

其中：

符号 含义 Ueda+22 Table 3 (2005–2009, \(|l|>105^\circ\), N=36)
\(n_{e0}\) 太阳位置的电子密度 \(3.4 \times 10^{-3}\) cm\(^{-3}\)
\(R_0\) 径向标长 \(7.0\) kpc
\(z_0\) 垂直标高 \(2.7\) kpc
\(R_\odot\) 太银心距 \(8.2\) kpc（本项目采用值）

关键直觉：\(R_0\) 控制密度在银盘平面内随银心距衰减的快慢，\(z_0\) 控制密度随银盘高度衰减的快慢。两个标长越大，气体分布越”平”。

1.3 本项目的额外假设

Ueda+22 原文只给了 \(n_e\) 公式。要把它转换成 X 射线亮度，本项目补充了三个假设：

  1. \(|z|\) 对称性：盘模型对银盘上下对称，即 \(n_e(R, z) = n_e(R, -z)\)。这与原文的 \(|z|\) 形式一致。
  2. \(n_e / n_H = 1.2\)：电子密度与氢核密度的比值。完全电离的氢+氦等离子体（氦丰度按质量约 25%）给出 \(n_e/n_H \approx 1.2\)。因此 \(n_H = n_e / 1.2\)。
  3. \(R_\odot = 8.2\) kpc：太阳到银心的距离（与项目其他部分保持一致）。

1.4 从密度到亮度：发射量积分

X 射线的热轫致辐射 + 线发射强度正比于发射量（Emission Measure, EM）：

\[ \mathrm{EM} = \int n_e(s)\, n_H(s)\, ds \]

其中积分沿视线方向 \(s\) 进行，单位是 cm\(^{-6}\) pc。直觉上：气体越密、路径越长，发射越强。

1.5 吸收转换与单位

得到 EM 后，还要经过两步：

  1. APEC × phabs 吸收转换：用 \(kT = 0.22\) keV、Lodders 元素丰度的 APEC 等离子体模型生成内禀谱，再用每条视线的 HI4PI 中性氢柱密度 \(N_H\) 做 phabs 光致吸收，最后在 0.5–2.0 keV 波段积分，得到 absorbed flux per EM（单位：flux unit / (cm\(^{-6}\) pc)）。
  2. 乘以 EM 得到 absorbed brightness：每条视线的预测亮度 = EM × flux-per-EM。

这套转换与项目里其他条件模型（如 Locatelli+24、Kaaret+20）使用完全相同的 APEC + phabs 流水线，保证 Figure 3 上各点可比。

Code 
graph TD
    A["Ueda+22 Table 3<br/>n_e0, R0, z0"] --> B["指数盘 n_e(R,z)"]
    B --> C["沿每条 M31 视线积分<br/>EM = ∫ n_e n_H ds"]
    C --> D["APEC(kT=0.22, Lodders)<br/>× phabs(HI4PI N_H)"]
    D --> E["absorbed 0.5-2.0 keV 亮度<br/>(14 条视线)"]
    E --> F["逆方差加权平均<br/>= Figure 3 的一个点"]

graph TD
    A["Ueda+22 Table 3<br/>n_e0, R0, z0"] --> B["指数盘 n_e(R,z)"]
    B --> C["沿每条 M31 视线积分<br/>EM = ∫ n_e n_H ds"]
    C --> D["APEC(kT=0.22, Lodders)<br/>× phabs(HI4PI N_H)"]
    D --> E["absorbed 0.5-2.0 keV 亮度<br/>(14 条视线)"]
    E --> F["逆方差加权平均<br/>= Figure 3 的一个点"]

2. 加载原始数据

Code 
ueda = pd.read_csv(DATA / "m31_cgmsum_ueda2022_disk_m31_footprint_absorbed_0p5_2p0.csv")
print(f"Loaded {len(ueda)} M31 footprint fields")
print(f"Columns: {list(ueda.columns)}")
ueda.head(8)
Loaded 14 M31 footprint fields
Columns: ['obsid', 'galactic_l_deg', 'galactic_b_deg', 'nh_hi4pi_1e22_cm-2', 'nominal_em_cm-6_pc', 'absorbed_0p5_2p0_flux_per_em_fluxunit', 'nominal_absorbed_0p5_2p0_fluxunit', 'n0_minus_1sigma_absorbed_0p5_2p0_fluxunit', 'n0_plus_1sigma_absorbed_0p5_2p0_fluxunit', 'measurement_staterr_absorbed_0p5_2p0_fluxunit']
obsid galactic_l_deg galactic_b_deg nh_hi4pi_1e22_cm-2 nominal_em_cm-6_pc absorbed_0p5_2p0_flux_per_em_fluxunit nominal_absorbed_0p5_2p0_fluxunit n0_minus_1sigma_absorbed_0p5_2p0_fluxunit n0_plus_1sigma_absorbed_0p5_2p0_fluxunit measurement_staterr_absorbed_0p5_2p0_fluxunit
0 800730201 119.831752 -19.888005 0.057791 0.002190 146.437298 0.320667 0.302081 0.339807 0.136534
1 800730301 119.508203 -20.173004 0.057602 0.002180 146.551086 0.319424 0.300910 0.338490 0.122502
2 800730501 119.959338 -20.297171 0.055286 0.002166 147.951714 0.320496 0.301920 0.339625 0.141052
3 800730601 120.221879 -19.978270 0.057930 0.002179 146.354000 0.318906 0.300422 0.337941 0.187279
4 800730701 120.355845 -20.368759 0.057273 0.002157 146.749162 0.316472 0.298130 0.335362 0.175246
5 800730801 120.090295 -20.632593 0.055956 0.002147 147.545133 0.316771 0.298411 0.335678 0.435051
6 800730901 119.868483 -20.906513 0.055921 0.002136 147.566082 0.315247 0.296976 0.334064 0.121189
7 800731101 120.520910 -20.803429 0.056700 0.002132 147.094590 0.313600 0.295424 0.332318 0.631616

数据包含 14 行，每行对应一个 M31 XMM 观测场。列含义如下：

列名 含义
obsid XMM-Newton 观测号
galactic_l_deg, galactic_b_deg 场中心的银道坐标
nh_hi4pi_1e22_cm-2 HI4PI 中性氢柱密度（\(10^{22}\) cm\(^{-2}\)）
nominal_em_cm-6_pc 名义发射量 EM（cm\(^{-6}\) pc）
absorbed_0p5_2p0_flux_per_em_fluxunit absorbed 0.5–2.0 keV flux per EM（flux unit / (cm\(^{-6}\) pc)）
nominal_absorbed_0p5_2p0_fluxunit 名义 absorbed 亮度（Figure 3 单位：\(10^{-15}\) erg cm\(^{-2}\) s\(^{-1}\) arcmin\(^{-2}\)）
n0_minus_1sigma_absorbed_0p5_2p0_fluxunit \(n_{e0}\) 下沿（\(-1\sigma\)）对应预测
n0_plus_1sigma_absorbed_0p5_2p0_fluxunit \(n_{e0}\) 上沿（\(+1\sigma\)）对应预测
measurement_staterr_absorbed_0p5_2p0_fluxunit 该场的实测 CGMsum 统计误差（用于加权）

3. 第一步：理解密度 → EM → 亮度的转换链

3.1 验证 nominal = EM × flux-per-EM

每条视线的名义预测亮度由两步乘法得到：

\[ S_{\mathrm{abs},i} = \mathrm{EM}_i \times F_{\mathrm{per\,EM}}(N_{H,i}, kT=0.22) \]

其中 \(F_{\mathrm{per\,EM}}\) 依赖该视线的吸收柱密度 \(N_H\)（HI4PI），但与 EM 无关——它是一个纯光谱转换因子。

Code 
em = ueda["nominal_em_cm-6_pc"].to_numpy()
flux_per_em = ueda["absorbed_0p5_2p0_flux_per_em_fluxunit"].to_numpy()
nominal = ueda["nominal_absorbed_0p5_2p0_fluxunit"].to_numpy()

manual = em * flux_per_em
print("EM      × flux-per-EM  =  nominal (CSV)")
for i in range(min(5, len(ueda))):
    print(f"  {em[i]:.6e} × {flux_per_em[i]:.6f} = {manual[i]:.6f}   (CSV: {nominal[i]:.6f})")

np.testing.assert_allclose(manual, nominal, rtol=1e-6)
print(f"\n✓ All {len(ueda)} fields: manual EM × flux-per-EM matches CSV nominal to < 1e-6")
EM      × flux-per-EM  =  nominal (CSV)
  2.189788e-03 × 146.437298 = 0.320667   (CSV: 0.320667)
  2.179606e-03 × 146.551086 = 0.319424   (CSV: 0.319424)
  2.166217e-03 × 147.951714 = 0.320496   (CSV: 0.320496)
  2.179003e-03 × 146.354000 = 0.318906   (CSV: 0.318906)
  2.156553e-03 × 146.749162 = 0.316472   (CSV: 0.316472)

✓ All 14 fields: manual EM × flux-per-EM matches CSV nominal to < 1e-6

3.2 可视化：EM、N_H 与 flux-per-EM 的关系

Code 
fig, (ax1, ax2) = plt.subplots(1, 2, figsize=(10, 4.5), constrained_layout=True)

ax1.scatter(ueda["galactic_b_deg"], em * 1e3, s=50, color="#2c8a4a", edgecolor="white", alpha=0.85)
ax1.set_xlabel("Galactic latitude b (deg)")
ax1.set_ylabel("Emission measure EM (10^-3 cm^-6 pc)")
ax1.set_title("EM vs Galactic latitude (14 M31 fields)")

sc = ax2.scatter(ueda["nh_hi4pi_1e22_cm-2"], flux_per_em, s=50, color="#2676b8", edgecolor="white", alpha=0.85)
ax2.set_xlabel("N_H (HI4PI, 10^22 cm^-2)")
ax2.set_ylabel("Absorbed flux per EM (flux unit / (cm^-6 pc))")
ax2.set_title("Flux-per-EM vs N_H")

plt.suptitle("Ueda+22 disk: EM and spectral conversion across the M31 footprint", fontsize=13, fontweight="bold")
plt.show()

14 个 M31 场的 EM、N_H 与 flux-per-EM。左：EM 随银纬变化；右：flux-per-EM 随 N_H 变化（N_H 越大，吸收越强，flux-per-EM 越低）。

观察： - 14 个场的 EM 都在约 \(2.0\)–\(2.2 \times 10^{-3}\) cm\(^{-6}\) pc 之间——盘模型在 M31 方向的预测相当均匀，因为 M31 的 14 个场在银道上很集中（\(l \approx 120^\circ\)–\(122^\circ\)，\(b \approx -20^\circ\)–\(-23^\circ\)）。 - \(N_H\) 越大的视线，flux-per-EM 越低——吸收更强，软 X 射线被吃掉更多。这是 phabs 吸收的直接体现。

4. 第二步：14 条视线的预测分布

4.1 为什么 14 条视线的预测不一样？

虽然盘模型在 M31 方向相当均匀，但每条视线的 \(N_H\) 略有不同（\(0.048\)–\(0.069 \times 10^{22}\) cm\(^{-2}\)），EM 也随银纬微弱变化。这导致 14 条视线的名义预测亮度有一个确定性的散布——这不是统计误差，而是盘模型在不同视线上的几何投影差异。

4.2 加入投影半径，看空间结构

Code 
# Merge with primary measurements to get signed projected M31 radius
primary = pd.read_csv(PRIMARY)[["obsid", "signed_rproj_kpc", "side"]]
merged = ueda.merge(primary, on="obsid", how="inner", validate="one_to_one")
assert len(merged) == 14, "directional diagnostic requires all fourteen M31 fields"
merged = merged.sort_values("signed_rproj_kpc").reset_index(drop=True)
merged[["obsid", "signed_rproj_kpc", "side", "nominal_absorbed_0p5_2p0_fluxunit"]]
obsid signed_rproj_kpc side nominal_absorbed_0p5_2p0_fluxunit
0 800732801 -26.096167 South/SE 0.304488
1 800732301 -16.734316 South/SE 0.302719
2 800732001 -13.206633 South/SE 0.289468
3 800731901 -11.173314 South/SE 0.292670
4 800731501 10.522662 North/NW 0.307894
5 800731601 13.073245 North/NW 0.300974
6 800731101 13.451164 North/NW 0.313600
7 800730801 18.923548 North/NW 0.316771
8 800730901 19.014166 North/NW 0.315247
9 800730701 19.544786 North/NW 0.316472
10 800730501 23.396669 North/NW 0.320496
11 800730601 25.028293 North/NW 0.318906
12 800730301 28.687554 North/NW 0.319424
13 800730201 28.793728 North/NW 0.320667

signed_rproj_kpc 是该场相对于 M31 中心的带符号投影半径：北/西北为正，南/东南为负。这样我们能看到盘模型预测在 M31 两侧是否有差异。

4.3 可视化：14 条视线的预测与加权

Code 
radius = merged["signed_rproj_kpc"].to_numpy()
nominal = merged["nominal_absorbed_0p5_2p0_fluxunit"].to_numpy()
low = merged["n0_minus_1sigma_absorbed_0p5_2p0_fluxunit"].to_numpy()
high = merged["n0_plus_1sigma_absorbed_0p5_2p0_fluxunit"].to_numpy()
errors = merged["measurement_staterr_absorbed_0p5_2p0_fluxunit"].to_numpy()

# Inverse-variance weights (using the measured-total statistical errors)
weights = np.square(errors) ** -1
weights_norm = weights / weights.sum()

fig, (ax1, ax2) = plt.subplots(2, 1, sharex=True, figsize=(9, 6),
                                constrained_layout=True,
                                gridspec_kw={"height_ratios": [1.6, 1]})

ax1.vlines(radius, low, high, color="#2c8a4a", alpha=0.45, linewidth=2,
           label="Marginal n0 sensitivity (not joint posterior)")
ax1.scatter(radius, nominal, color="#2c8a4a", s=45, zorder=3, edgecolor="white",
            label="Nominal field prediction")
ax1.axhspan(nominal.min(), nominal.max(), color="#2c8a4a", alpha=0.08,
            label="Deterministic footprint min-max")
mean = float(np.sum(nominal * weights) / np.sum(weights))
ax1.axhline(mean, color="#d47a2c", linewidth=1.5, label=f"Weighted mean = {mean:.6f}")
ax1.set_ylabel("Predicted brightness (flux unit)")
ax1.set_title("14/14 fields inside adopted angular domain", loc="left", fontsize=10, fontweight="bold")
ax1.legend(loc="upper left", fontsize=8)

ax2.bar(radius, weights_norm, width=1.5, color="#2c8a4a", alpha=0.75)
ax2.set_xlabel("Signed projected M31 radius (kpc)")
ax2.set_ylabel("Normalized\nstatistical weight")

plt.suptitle("Ueda+22 directional footprint and estimator", fontsize=13, fontweight="bold")
plt.show()

14 个 M31 场的 Ueda+22 盘模型预测。绿点：名义预测；绿色竖线：n0 ±1σ 边缘灵敏度（不是联合后验）；橙色横线：逆方差加权平均 = Figure 3 的点。

关键观察： - 14 条视线的名义预测在约 \(0.289\)–\(0.321\) flux unit 之间，散布很窄——盘模型在 M31 方向几乎是”一个值”。 - 绿色竖线（\(n_0 \pm 1\sigma\) 边缘灵敏度）比散布宽——这说明参数不确定性比视线几何散布更大，但两者都不是联合后验区间。 - 北/西北侧（正半径）的预测略高于南/东南侧（负半径），因为盘模型在 M31 两侧的路径长度和 \(N_H\) 略有不同。 - 下图的权重分布不均——某些场（如 OBSID 0800730801，统计误差 0.435）权重很低，因为它们的实测 CGMsum 统计误差很大；误差小的场（如 0800731601，0.094）权重大。

5. 第三步：从 14 个值到一个加权平均点

5.1 为什么用逆方差加权平均？

Figure 3 上的点要代表”14 条视线的整体预测”。我们用逆方差加权平均（inverse-variance weighted mean）：

\[ \bar{S} = \frac{\sum_i w_i S_i}{\sum_i w_i}, \qquad w_i = \frac{1}{\sigma_i^2} \]

其中 \(\sigma_i\) 是第 \(i\) 个场的实测 CGMsum 统计误差（不是模型预测的散布）。直觉上：实测误差小的场”更可信”，应该在平均中占更大比重。

这与 Henley & Shelton 2013 教程里用 median + IQR 不同——H&S13 的 86 条视线是独立的实测样本，取 median 抗离群；而 Ueda+22 的 14 条视线是同一个模型的不同投影，用加权平均更自然，权重来自实测质量。

5.2 计算与验证

Code 
# Inverse-variance weighted mean (matches the frozen Figure 3 ledger)
mean = float(np.sum(nominal * weights) / np.sum(weights))
print(f"14-field inverse-variance weighted mean = {mean:.15f}")

# Verify against the frozen Figure 3 ledger
LEDGER_CENTRAL = 0.311904286700407
LEDGER_LOW  = 0.2894675687845744   # nominal footprint min
LEDGER_HIGH = 0.3206666288460409   # nominal footprint max
LEDGER_CONTEXT_LOW  = 0.27269046920968987   # n0 -1sigma min
LEDGER_CONTEXT_HIGH = 0.339806764996886     # n0 +1sigma max

np.testing.assert_allclose(mean, LEDGER_CENTRAL, rtol=0, atol=5e-13)
np.testing.assert_allclose(nominal.min(), LEDGER_LOW,  rtol=0, atol=1e-10)
np.testing.assert_allclose(nominal.max(), LEDGER_HIGH, rtol=0, atol=1e-10)
np.testing.assert_allclose(low.min(),  LEDGER_CONTEXT_LOW,  rtol=0, atol=1e-10)
np.testing.assert_allclose(high.max(), LEDGER_CONTEXT_HIGH, rtol=0, atol=1e-10)
print("✓ Matches frozen Figure 3 ledger to < 5e-13")
print(f"  Central:           {mean:.6f}")
print(f"  Footprint min-max: [{nominal.min():.6f}, {nominal.max():.6f}]")
print(f"  n0 context bounds: [{low.min():.6f}, {high.max():.6f}]")
14-field inverse-variance weighted mean = 0.311904286700407
✓ Matches frozen Figure 3 ledger to < 5e-13
  Central:           0.311904
  Footprint min-max: [0.289468, 0.320667]
  n0 context bounds: [0.272690, 0.339807]

5.3 可视化：14 条视线预测与 Figure 3 点

Code 
fig, ax = plt.subplots(figsize=(8, 4.5))
ax.hist(nominal, bins=8, color="#2c8a4a", alpha=0.7, edgecolor="white", linewidth=0.5)
ax.axvline(mean, color="#d47a2c", linewidth=2, label=f"Weighted mean = {mean:.3f}")
ax.axvspan(nominal.min(), nominal.max(), color="#6f7782", alpha=0.12,
           label=f"Footprint min-max [{nominal.min():.3f}, {nominal.max():.3f}]")
ax.axvspan(low.min(), high.max(), color="#2c8a4a", alpha=0.12,
           label=f"n0 context [{low.min():.3f}, {high.max():.3f}]")
ax.axvline(0.964727, color="#b74842", linestyle="--", linewidth=1.5, label="Observed CGMsum = 0.965 (M31)")
ax.set_xlabel("Absorbed 0.5-2.0 keV brightness (10^-15 erg s^-1 cm^-2 arcmin^-2)")
ax.set_ylabel("Number of fields")
ax.set_title("Ueda+22 disk: 14-field prediction distribution")
ax.legend(fontsize=8, loc="upper right")
plt.show()

14 条 M31 视线的 Ueda+22 盘模型预测分布。橙色竖线为逆方差加权平均（Figure 3 的点），灰色带为 footprint min-max，绿色带为 n0 ±1σ context bounds。

关键观察： - 14 条视线的预测集中在 \(0.29\)–\(0.32\) 附近，加权平均 \(0.312\)——这是 Ueda+22 盘模型对 M31 方向银河系前景的预测。 - M31 方向实测 CGMsum（\(0.965\)）远高于这个预测——这说明 Ueda+22 的盘模型只能解释 M31 方向前景的一部分，其余来自 M31 自身的 CGM 或其他未建模成分。 - footprint min-max（灰色）比 \(n_0\) context bounds（绿色）窄——前者是视线几何散布，后者是参数边缘灵敏度，两者不是同一个不确定性的两种估计。

6. 关键假设清单（必须记住！）

步骤 假设 来源
密度模型 指数电子密度盘 \(n_e(R,z)\) Ueda+22 Table 3 (2005–2009, \(|l|>105^\circ\), N=36)
盘参数 \(n_{e0}=3.4\times10^{-3}\) cm\(^{-3}\), \(R_0=7.0\) kpc, \(z_0=2.7\) kpc Ueda+22 Table 3
\(|z|\) 对称 \(n_e(R,z) = n_e(R,-z)\) 本项目（与原文 \(|z|\) 形式一致）
电子/氢比 \(n_e / n_H = 1.2\) 本项目（完全电离 H+He 等离子体）
太银心距 \(R_\odot = 8.2\) kpc 本项目（与项目其他部分一致）
发射量 \(\mathrm{EM} = \int n_e n_H ds\) 本项目积分
光谱模型 APEC(\(kT=0.22\) keV, Lodders 丰度) × phabs(HI4PI) 本项目转换流水线
吸收柱密度 HI4PI 巡天 \(N_H\) 每场独立
统计量 逆方差加权平均（权重来自实测 CGMsum 统计误差） 本项目选择
区间含义 footprint min-max 是几何散布；\(n_0 \pm 1\sigma\) 是参数边缘灵敏度；两者不是联合后验区间 本项目决定
M31 关系 14 个场均在采用 Table 3 行的角域内（\(|l|>105^\circ\)） 本项目验证

7. 动手练习

练习 1：验证一条视线的乘法链

选择第 1 条视线（OBSID 0800730201），用它的 EM 和 flux-per-EM 手动计算 nominal absorbed brightness，并与 CSV 比较。

OBSID 800730201.0:
  EM           = 2.189788e-03 cm^-6 pc
  flux-per-EM  = 146.437298
  Manual nominal = 0.320667
  CSV nominal     = 0.320667
  ✓ Match

练习 2：不加权 vs 加权平均

计算 14 条视线名义预测的简单算术平均和逆方差加权平均，比较两者差异。为什么加权平均更接近”高质量场”的预测？

Simple mean:   0.309985
Weighted mean: 0.311904
Difference:    0.001919
The weighted mean is pulled toward fields with smaller statistical errors,
which represent higher-quality CGMsum measurements.

练习 3：N_H 与 flux-per-EM 的相关性

计算 14 个场的 \(N_H\) 与 flux-per-EM 之间的 Pearson 相关系数。这个相关性是正还是负？为什么？

Pearson r = -0.9998
Negative — larger N_H means stronger photoelectric absorption,
so fewer soft X-ray photons survive, lowering flux-per-EM.

8. 总结：从原始数据到 Figure 3 的完整链路

Code 
graph TD
    A["Ueda+22 Table 3<br/>2005-2009, |l|>105, N=36"] --> B["指数盘 n_e0=3.4e-3<br/>R0=7.0, z0=2.7 kpc"]
    B --> C["本项目补充：<br/>|z| 对称, n_e/n_H=1.2, R_sun=8.2"]
    C --> D["沿 14 条 M31 视线积分<br/>EM = ∫ n_e n_H ds"]
    D --> E["APEC(0.22 keV, Lodders)<br/>× phabs(HI4PI N_H)"]
    E --> F["14 个 absorbed 亮度<br/>(flux unit)"]
    F --> G["逆方差加权平均<br/>= 0.3119"]
    G --> H["Figure 3: Ueda+22 disk<br/>directional conditional point"]

graph TD
    A["Ueda+22 Table 3<br/>2005-2009, |l|>105, N=36"] --> B["指数盘 n_e0=3.4e-3<br/>R0=7.0, z0=2.7 kpc"]
    B --> C["本项目补充：<br/>|z| 对称, n_e/n_H=1.2, R_sun=8.2"]
    C --> D["沿 14 条 M31 视线积分<br/>EM = ∫ n_e n_H ds"]
    D --> E["APEC(0.22 keV, Lodders)<br/>× phabs(HI4PI N_H)"]
    E --> F["14 个 absorbed 亮度<br/>(flux unit)"]
    F --> G["逆方差加权平均<br/>= 0.3119"]
    G --> H["Figure 3: Ueda+22 disk<br/>directional conditional point"]

一句话总结：Ueda et al. (2022) 用 Suzaku 拟合银河系软 X 射线背景，得到约 0.22 keV MWH 成分的指数电子密度盘。本项目沿 M31 十四条 XMM 视线积分 \(\mathrm{EM} = \int n_e n_H ds\)，用 APEC + HI4PI phabs 转 absorbed 0.5–2.0 keV，再用实测 CGMsum 统计权重做逆方差加权平均，得到 \(0.312\)——这是 Ueda+22 盘模型对 M31 方向银河系前景的方向性条件预测，不是 M31 CGM 本身的测量。Footprint min-max 是视线几何散布，\(n_0 \pm 1\sigma\) 是参数边缘灵敏度，两者都不是联合后验区间。

参考资料

  1. Ueda, Y. et al. (2022). “The soft X-ray background with Suzaku. I. Milky Way halo.” PASJ, 74, 1094–1111. doi:10.1093/pasj/psac077
  2. HI4PI 全天中性氢巡天 — HI4PI Collaboration et al. (2016), A&A, 594, A116
  3. Lodders, K. (2003). “Solar System Abundances and Condensation Temperatures of the Elements.” ApJ, 591, 1220
  4. APEC 等离子体发射模型 — Smith et al. (2001), ApJ, 556, L91

教程结束 🎓 下一步：继续阅读 Kaaret+2020 HaloSat tutorial，了解另一个银河系晕方向性模型如何投影到 M31。