ECCV 2026 3D 视觉自动驾驶单光子 LiDAR SPAD 死时间马尔可夫更新过程光子计数统计直方图立方体仿真

Markov-Renewal Single-Photon LiDAR Simulator¶

会议: ECCV 2026
论文: ECCV 2026 / 项目页
领域: 3D 视觉 / 自动驾驶
关键词: 单光子 LiDAR、SPAD 死时间、马尔可夫更新过程、光子计数统计、直方图立方体仿真

一句话总结¶

把 SPAD 死时间下的光子注册过程写成马尔可夫更新过程（MRP），用闭式公式解析预测每个像素的光子计数分布（高斯均值与方差）和时间 PDF，再靠「延迟不变性 + 预计算查表」一次并行采样出整幅 H×W×n_b 直方图立方体：统计上与逐光子拒绝采样的 gold standard 无法区分，但生成 128×128×1024 立方体只要 1.01 s 对 40.8 s（>40×）。

研究背景与动机¶

单光子 LiDAR（SP-LiDAR）用 SPAD 阵列配合 TCSPC 计时，把每个激光周期内探测到的光子时间戳累积成直方图，从而以单光子灵敏度做高分辨率 3D 成像。但单光子灵敏度也带来了麻烦的测量统计：光子到达本身是随机点过程，而探测器的状态又反过来改变这个过程——SPAD 每注册一个光子就会进入长达 \(t_d\) 的死时间，期间到达的光子被直接忽略（异步工作模式下要等 \(t_d\) 结束才重新"睁眼"）。这个效应是非线性的：它破坏检测事件之间的统计独立性，产生时间相关性，同时压低光子注册率并且让直方图形状发生畸变（经典 pile-up 非线性）。因此任何忽略死时间的 Poisson / 多项分布仿真器在物理上都是错的。

现有仿真器被卡在一个两难里。一侧是 Poisson 类模型：采样快，但物理上错，生成的高通量直方图既没有 pile-up 畸变，背景光子也是均匀铺满整周期的。另一侧是顺序式「gold standard」仿真器：先按多项泊松采样出所有光子到达，再逐光子判断是否落在死时间窗口内并拒绝，物理正确，但需要 per-photon、per-pixel 的循环，通量、空间分辨率和实现数（realization）任何一个上去就算不动——128×128×1024 的立方体要 40.8 s，而深度学习训练需要成千上万个这样的立方体。这正是本文要打掉的瓶颈：因为缺少可规模化的高保真仿真，SP-LiDAR 的重建网络只能拿简化的低通量仿真器或物理错误的模型来训练。近期的 Zhang 等人试图用神经网络直接学被死时间畸变的统计量绕开逐光子循环，但他们的时间 PDF 估计器每换一个死时间设置就要重训，计数预测在方差上明显不准，整个流程仍是逐像素而非完全并行。

关键矛盾在于：仿真一个直方图需要两样东西——时间 PDF 和光子计数——而前者已经被解决、后者一直没人做对。Rapp 等人的 Markov 链模型（后经 Zhang 等人用"死时间即行移位"的改写加速）能准确给出平稳时间 PDF，但它只描述"下一个注册落在哪个 bin"，把绝对时刻 \(T_k\) 丢掉了，因此回答不了"这段曝光时间里到底注册了多少个光子"，也就定不出直方图的尺度。而标准更新理论（renewal theory）要求间隔独立同分布，在这里不成立——死时间让相邻注册的间隔强相关；自激点过程则计算上不可行。本文的核心 idea 是把注册过程重新表述为马尔可夫更新过程（MRP）：状态转移给出相对时刻，保持时间（holding time）补回绝对时间，于是"计数 \(\times\) 形状"可以分开解析求解——计数是累积保持时间的逆过程，混合性保证它满足中心极限定理，从而得到高斯的极限分布与闭式均值/方差；再由延迟不变性把三维查找表压成二维，整幅立方体用批量内存读取加廉价平移一次采样出来。

方法详解¶

整体框架¶

MaRS 要解的输入是：系统配置 \(\theta=(t_r, N_r, \sigma_t, t_d, n_b, \tau, S, B, N_{\text{iter}})\)，其中 \(t_r\) 是激光重复周期、\(N_r\) 是脉冲数（总曝光时间 \(t=N_r t_r\)）、\(\sigma_t\) 是脉冲宽度、\(t_d\) 是死时间、\(n_b\) 是直方图 bin 数；再加上场景参数——深度图 \(Z\)、反射率图 \(R\) 与背景水平 \(B\)。输出是 \(N_{\text{iter}}\) 个 \(H\times W\times n_b\) 的直方图立方体。与顺序仿真器"先生成每个光子再逐个拒绝"不同，MaRS 的核心策略是把直方图拆成"计数 \(\times\) 形状"：形状由 Markov 链的平稳时间 PDF 给出，计数由本文新提出的 MRP 计数模型给出；两者都只依赖像素级的 \((S,B,t_d)\) 与延迟 \(\tau\)，因此可以预计算、查表、并行采样。整条管线是：统计引擎先算出每像素的计数分布与时间 PDF，谱截断让其中的方差计算变得可算，再用延迟不变性把查找表从三维降到二维，最后批量查表 + 平移 + 两级采样一次性得到整个立方体。

%%{init: {'flowchart': {'rankSpacing': 24, 'nodeSpacing': 28, 'padding': 6, 'wrappingWidth': 400}}}%%
flowchart TD
    A["输入：深度 Z / 反射率 R / 背景 B<br/>+ 系统参数 θ"] --> B["统计引擎：形状 + 计数<br/>Markov 链平稳 PDF ＋ MRP 计数"]
    B --> C["谱截断协方差<br/>长程相关性闭式补齐"]
    C --> D["查表 + 延迟不变性<br/>统计量与延迟解耦"]
    D --> E["并行立方体采样<br/>计数 × 形状一次生成"]
    E --> F["直方图立方体<br/>→ 深度 / 反射率重建"]

关键设计¶

1. MRP 光子计数模型：把 Markov 链丢掉的绝对时间找回来

前面说的 Markov 链模型能把时间 PDF 解得很准，但它的状态是相对时刻 \(X_k = T_k \bmod t_r\)，绝对时刻 \(T_k\) 被完全丢弃，"总共来了多少个光子"就无从谈起——直方图形状对了，辐射尺度还是错的，而尺度直接决定反射率估计的偏置。本文的做法是给这条链补上一个绝对时间维度：把注册过程写成 MRP \(\{(X_k, W_k)\}\)，\(X_k\) 仍是 Markov 状态（相对时刻），\(W_k = T_k - T_{k-1}\) 是绝对尺度上的保持时间，两者联合建模。转移核只比原来的转移矩阵多一个条件保持时间分布：

\[Q_{ij}(t) = \Pr\!\left(X_{k+1}=j,\; W_{k+1}\le t \mid X_k = i\right) = \tilde{P}_{ij}\,F_{ij}(t)\]

其中 \(\tilde{P}_{ij}\) 就是已有的 Markov 转移概率，\(F_{ij}(t)\) 是转移 \(i\to j\) 的保持时间条件 CDF。这一步的巧妙之处在于：计数 \(N(t)\) 是累积保持时间 \(T_n=\sum_{k\le n}W_k\) 的逆过程——注册越密，每次的增量 \(W_k\) 越小。虽然 \(W_k\) 之间因死时间而相关、不满足经典更新理论的 i.i.d. 要求，但底层的 Markov 状态混合得足够快，\(T_n\) 仍满足中心极限定理，把它反解过来就得到光子计数的渐近高斯分布：

\[N(t) \xrightarrow{\;d\;} \mathcal{N}\!\left(\frac{t}{\mu},\; \frac{t\,\sigma^2}{\mu^3}\right)\quad (t\to\infty),\qquad \mu \triangleq \lim_{n\to\infty}\frac{\mathbb{E}[T_n]}{n},\quad \sigma^2 \triangleq \lim_{n\to\infty}\frac{\mathrm{Var}[T_n]}{n}\]

（方差项中 \(\mu^3\) 的分母来自"逆过程"展开，⚠️ 缓存文本此处公式有损坏，以原文为准。）接下来要把 \(\mu,\sigma^2\) 写成可算的闭式：条件于一次转移 \(i\to j\)，保持时间由两部分组成——从探测器重新激活的时刻 \(x_i'\) 到下一个注册时刻 \(x_j\) 的最小清醒间隔 \(t_{i'j}\)，加上若干个完整激光周期的等待；而"一个周期内至少来一个光子"的概率是 \(1-e^{-\Lambda}\)（\(\Lambda=S+B\) 是每周期总光子通量），所以等待的完整周期数服从几何分布，于是

\[\mu_{ij} = t_{i'j} + t_r\frac{e^{-\Lambda}}{1-e^{-\Lambda}},\qquad \sigma^2_{ij} = t_r^2\frac{e^{-\Lambda}}{\left(1-e^{-\Lambda}\right)^2}\]

（⚠️ 缓存文本中该引理的两式 OCR 有损坏，此处按几何等待的物理含义复原，以原文为准。）有了这些，构造一阶矩核 \(Q_\mu=\tilde P\odot M\)、二阶矩核 \(Q_{\mu^2}=\tilde P\odot M\odot M\) 与方差核 \(Q_{\sigma^2}=\tilde P\odot\Sigma\)（\(\odot\) 为逐元素乘，\(M,\Sigma\) 收集所有 \(\mu_{ij},\sigma^2_{ij}\)），有效均值就是 \(\mu=\pi Q_\mu \mathbf{1}\)（\(\pi\) 是 Markov 链的平稳分布，由 \(\pi \tilde P=\pi\) 解出），有效方差拆成稳态项加协方差项 \(\sigma^2=\sigma^2_{ss}+\gamma_{ss}\)，其中 \(\sigma^2_{ss}=\pi(Q_{\mu^2}+Q_{\sigma^2})\mathbf{1}\)。至此"计数"第一次有了解析答案，而且它与"形状"共用同一个 \(\tilde P\)、同一套平稳分布，两条线天然自洽。

2. 谱截断协方差：给无穷长的相关性一个收敛的闭式尾巴

\(\gamma_{ss}=\sum_{l\ge1}\gamma_l\) 是个无穷级数，必须截断才好算，但不同 \((S,B,t_d)\) 下 \(\gamma_l\) 的衰减形态差异极大（Fig. 4(a)）：固定截断 c 要么取大导致算得慢，要么取小引入偏差，没有普适的 c。本文借用了 Zhang 等人关于"\(\tilde P\) 的长期动态由少数主导特征值控制"的观察：短程的协方差项依赖全部特征模、但项数少、直接精确求和最划算；长程的相关性则由衰减最慢的几个特征模决定，可以用低秩近似 \(\tilde P^{\,l-1}\approx V_p\Lambda_p^{\,l-1}U_p^\top\) 解析地把尾巴补上（\(\{(\lambda_i,v_i,u_i)\}_{i=1}^p\) 是前 \(p\) 个特征对，满足双正交 \(U_p^\top V_p=I_p\)）：

\[\hat\gamma_{ss}(L,p) = \sum_{l=1}^{L}\gamma_l \;+\; \sum_{i=2}^{p}\frac{\alpha_i\beta_i\,\lambda_i^{\,L}}{1-\lambda_i},\qquad \alpha_i = (\pi Q_\mu V_p)_i,\quad \beta_i=(U_p^\top Q_{\mu^2})_i\]

第一项是精确计算的短程相关性，第二项是用主导特征模写的闭式尾项。这样做的收益不是"近似得更准"，而是把不可算的无穷和变成了常数级开销：实验发现第 5 个之后的特征值模长普遍降到 \(\lesssim 0.7\)，超过 \(L\approx 6\) 步后贡献可忽略，所以统一取 \(L=6\)、\(p=5\) 就能覆盖所有测试过的通量/死时间区间。作者给出的对照很有说服力：在若干困难设置下，即使老老实实累加 20 项协方差也仍然抓不到长程尾巴（方差不收敛），而这条谱规则用几乎为零的额外算力就把它补了回来。

3. 查表 + 延迟不变性：用一条物理对称性把三维表压成二维

即便前两个设计让统计量都可算了，直接算仍然不可 scale：在原生时间分辨率（例如 \(n_b=2^{14}\)）上实时求计数统计很贵，而每个像素的 \((S,B,\tau)\) 都不同，逐像素重算等于把顺序仿真器的病换个地方再犯一次。最自然的想法是预计算查找表（LUT），但把三个变量都表化会直接爆炸：取 \(n_S=256\)、\(n_B=100\)、\(n_\tau=1024\)、\(n_b=1024\)，每个表项要存 \((\mu,\sigma^2)\) 加一个长度 \(n_b\) 的 PDF \(\pi\)（\(n_b+2=1026\) 个 float），总存储 \(256\cdot100\cdot1024\cdot1026\cdot4\approx 108\) GB，完全不实用。

转折点是一条物理对称性——延迟只平移时间，不改变更新结构。把每像素的到达率写成 \(\lambda_\tau(t)=S\,s(t-\tau)+B/t_r\)，改变 \(\tau\) 只是把通量的时间原点挪一挪：光子到达的随机机制、死时间的作用方式都没变，所以注册计数的均值与方差与 \(\tau\) 无关，时间 PDF 也只是循环平移：

\[\mu_\tau=\mu_{\tau'},\qquad \sigma^2_\tau=\sigma^2_{\tau'},\qquad \boldsymbol{\pi}_{\tau'}=\mathrm{shift}_{(\tau'-\tau)}\left(\boldsymbol{\pi}_\tau\right)\]

（距离引起的回波衰减被吸收进 \(S\)（或 \(B\)），作为应用该定理前的一步预处理。）于是 LUT 只需在规范延迟（脉冲居中于 \(t_r/2\)）上对 \((S,B)\) 建表，\(\tau\) 推迟到采样时用一次平移解决：存储从 108 GB 降到 \(256\cdot100\cdot1026\cdot4\text{ B}\approx100\) MB，可以整表驻留、批量访问。整个立方体生成就成了三步——批量查表取 \((\mu,\sigma^2,\pi_0)\)、按 \(Z-t_r/2\) 平移 PDF、先抽计数图 \(n_i\sim\mathcal{N}(\mu,\sigma^2)\) 再抽直方图 \(H_i\sim\text{Multinomial}(n_i,\Pi)\)——全部是向量化内存读取、平移和采样，没有任何逐光子或逐像素循环。注意这张表不是通用的：它绑定在给定的系统配置 \((t_r,n_b,\sigma_t,t_d)\) 与脉冲形状 \(s(\cdot)\) 上，换到别的工况（或别的 \((S,B)\) 范围与离散化）需要重建。

一个完整示例¶

设一块 128×128 的深度图与反射率图，系统参数取 \(t_r\) 对应 7.5 MHz 脉冲、\(n_b=1024\)，查找表按 \(S\in[0.01,10]\)（256 点）、\(B\in[0.01,10]\)（100 点）建好，整表 100.4 MB。对某个像素：它的反射率 \(R\) 与背景 \(B\) 决定 \((S,B)\)，直接落进表网格取表项，或落在相邻格点之间取值，拿到这一像素的 \((\mu,\sigma^2,\pi_0)\)——此时 \(\pi_0\) 还是"脉冲居中"的形状；这个像素的深度 \(Z\) 换算成延迟 \(\tau=2Z/c\) 后，\(\pi_0\) 被循环平移 \(\tau-t_r/2\)，变成该像素真正的回波时间分布。要生成一个 realization，就先从这个像素的高斯计数分布里抽一个整数光子数 \(n\)，再按平移后的 \(\pi\) 做一次多项分布抽样，得到这条 1024 维的直方图。整张图 16384 个像素同时做这件事，重复 \(N_{\text{iter}}\) 次。实测三段耗时分别为：批量查表 \(0.020\pm0.005\) s、时间对齐 \(0.144\pm0.008\) s、直方图采样 \(0.849\pm0.070\) s，合计 \(1.013\pm0.071\) s；同一立方体交给顺序 gold standard 要 \(40.831\pm8.045\) s。

损失函数 / 训练策略¶

MaRS 本身是仿真器，没有可训练参数，也没有损失函数；它的"训练策略"意义在于为下游学习式重建供数。实验中使用一个 DDPM U-Net（Ho 等人），改造后直接吃完整的时空直方图立方体，同一架构分别训练两个模型做深度恢复与反射率恢复，训练到收敛为止，四个版本（Poisson / Renewal / Zhang 等人 / MaRS）共享同一架构与同一张 LUT 规格，区别只在训练数据的仿真器来源。评估分两轮：先在各自仿真器的数据上测试以确认训练本身稳定，再把所有模型统一放到 gold standard 仿真器的数据上测，用后者衡量跨域泛化。

实验关键数据¶

主实验¶

光子计数分布本身的精度（在 \((S,B,t_d)\) 网格上取平均，ground truth 是顺序仿真器给出的经验分布）。这是全篇最核心的一张表，因为它直接检验"计数"这一被前人忽略的统计量。

方法	Wasserstein 距离 ↓	KL 散度 ↓	均值误差 ↓	方差误差 ↓
Poisson	7966.130	24.208	7966.130	9257.464
Renewal（等能量常速率近似）	207.565	18.821	207.564	45.927
Zhang 等人 [39]	2.337	0.409	0.282	64.558
MaRS	0.309	0.013	0.269	1.588

下游深度与反射率重建（同一 U-Net 架构，分别用四个仿真器的数据训练；「自身」= 在自己仿真器的数据上测，「GS」= 在 gold standard 仿真器数据上测，后者最接近真实传感器行为）。

训练用仿真器	测试域	反射率 PSNR ↑	反射率 RMSE ↓	深度 PSNR ↑	深度 RMSE ↓
Poisson	自身	36.4946	0.0722	33.6827	0.0606
Poisson	GS	7.0237	1.6537	7.2490	0.7853
Renewal	自身	29.7230	0.0543	30.9704	0.0514
Renewal	GS	12.1210	0.0959	11.6988	0.0967
Zhang 等人 [39]	自身	22.7444	0.0526	21.9506	0.0506
Zhang 等人 [39]	GS	22.5421	0.0820	22.2137	0.1193
MaRS	自身	24.0321	0.0532	24.2821	0.0553
MaRS	GS	24.2488	0.0739	23.5776	0.0949

（表中两列组对应原文的 Reflectivity 与 Depth 两组指标；原文表格同时给出 SBR=0.6 与 SBR=1.0 两档设置，缓存文本只保留了每行 10 个数值，⚠️ 完整分组与 SBR=1.0 的数值以原文 Table 5 为准。）

消融实验¶

四个仿真器各自"做对了什么统计量"的对照（\(\mu\) 计数均值、\(\sigma^2\) 计数方差、\(\pi\) 时间 PDF、LUT 使用）。这张表解释了后面所有性能差异的来源。

仿真器	死时间校正的 \(\mu\)	死时间校正的 \(\sigma^2\)	死时间校正的 \(\pi\)	使用 LUT
Poisson	✗	✗	✗	✗
Renewal	✗	✗	✓	✓
Zhang 等人 [39]	✓	✗	✓	✓
MaRS	✓	✓	✓	✓

生成 128×128×1024 立方体的耗时分解，以及 LUT 与顺序采样的对照。

阶段 / 方法	时间 (s)
批量查表	0.020 ± 0.005
时间对齐（按 \(\tau\) 平移 PDF）	0.144 ± 0.008
直方图采样	0.849 ± 0.070
MaRS 总计	1.013 ± 0.071
顺序 gold standard	40.831 ± 8.045

关键发现¶

计数方差的缺失是 Zhang 等人 [39] 的主要短板，而不在均值：他们的均值误差 0.282 与 MaRS 的 0.269 基本持平，但方差误差 64.558 对 1.588，差了约 40 倍——说明"把均值做对"远不足以确定直方图的辐射尺度。
在物理错误的仿真器上训练出的网络完全不迁移：Poisson 在自身数据上反射率 PSNR 高达 36.49，换到 GS 域直接掉到 7.02、SSIM 变成负值（−0.075）；Renewal 也从 29.72 掉到 12.12。而 MaRS 的跨域差距几乎为零（反射率 24.03 → 24.25，深度 24.28 → 23.58），在 GS 域上取得最好的深度与反射率恢复。这直接把"死时间校正的光子计数"从分布层面的改进坐实为下游任务层面的收益。
单看域内分数会得出错误结论：Poisson 的域内 PSNR 是全场最高，但那是它自己的仿真数据太"干净"所致，只有跨仿真器评测才暴露真实差距。
LUT 是效率的关键：带 LUT 的 MaRS 在通量、realization 数、像素数三个维度上都近似常数时间；不带 LUT 的版本因反复做大矩阵运算而明显退化；顺序仿真器在三个维度上都急剧变慢。单像素场景下三者的量级也很直观：Poisson 0.2 s、MaRS 0.3 s、gold standard 16.1 s。
谱截断的参数选择有依据：跨 \((S,B,t_d)\) 设置统计发现第 5 个之后的特征值模长 \(\lesssim0.7\)，故统一取 \(L=6\)、\(p=5\)；对照实验显示困难设置下即使累加 20 项协方差也补不上长程尾巴，而谱规则以可忽略的开销精确恢复。此外，场景的 \((S,B,\tau)\) 未必落在 LUT 网格点上，但由此产生的量化误差很小，能被重建网络吸收。
真实硬件验证超出了"与仿真器一致"：作者用 7.5 MHz 脉冲激光 + 死时间 85 ns 的 SPAD 实测，以 9975 个激光周期内的探测数记作一个 realization、重复 10000 次，按 Rapp 等人的标定流程估计信号与背景水平后预测计数分布；MaRS 在固定距离与环境光下的亮目标与暗目标两个像素上都与实测分布吻合。这说明 MRP 计数模型抓的是真实探测器的死时间统计，而不只是复现了另一个仿真器。

亮点与洞察¶

"直方图 = 计数 × 形状"的解耦是全篇的杠杆。前人（Rapp、Zhang 等）把全部精力花在形状（时间 PDF）上并且做得很好，本文指出没被解决的其实是尺度（计数），而尺度错了会以辐射偏置的形式污染反射率估计。把问题重新切一刀，比在旧切法上继续优化更有效。
MRP 的用法很巧：Markov 链负责相对状态、更新过程负责绝对时间，两者用一个条件保持时间分布缝合；再借"计数是累积保持时间的逆过程"把原本难处理的 \(N(t)\) 转成"CLT + 反解"。经典更新理论在这里失效（间隔不 i.i.d.），但混合性足够强，结果依然是高斯的——这个"退一步仍然成立"的判断是理论部分最漂亮的地方。
用一条物理对称性换三个数量级的存储：延迟不变性把一个 108 GB 的三维表变成 100 MB 的二维表加一次平移。任何被"全局平移参数"参数化的仿真（不同距离、不同触发时刻）都可以照搬这个套路：把平移不变的部分预计算，把平移留到采样时做。
可迁移的模板：凡"状态相关的点过程 + 需要绝对计数"的场景——死时间受限的探测器阵列、通信中的泊松体制、单光子成像以外的计数型传感器——都可以套用"状态链 + 保持时间 + CLT 反解 + 低秩谱尾"这条链路。

局限与展望¶

作者明确承认 LUT 不通用：它绑定在给定的系统配置 \((t_r,n_b,\sigma_t,t_d)\) 与脉冲形状上，换工况（或换 \((S,B)\) 的取值范围与离散化）需重新生成。
建表的通量范围被限定在 \(S,B\in[0.01,10]\)（偏室内量级），超出该范围的强反射或远距低通量场景需要重建表或外推，而论文未报告外推行为。
高斯近似来自 CLT 的渐近性，其有效性依赖 \(t=N_r t_r\) 足够大（\(\Lambda\) 足够高）。实验网格集中在高通量死时间区间，低通量/极少光子（\(N(t)\) 很小时）下"高斯 \(\times\) 多项"的组合是否仍成立、失效边界在哪，论文没有给出。
效率实验全部在 NumPy CPU 上完成（AMD Ryzen 7 9700X，32 GB DDR5），没有 GPU 实现，对照基准也是同机的顺序实现；与高度优化的 GPU 光子级仿真相比优势有多大仍是开放问题。
真机验证只覆盖固定距离、固定环境光下的亮/暗两个像素分布，未覆盖距离扫描、多目标、极端动态范围等情形。
可以想到的改进：把 LUT 做成可微或带自适应插值的版本以减小量化误差；显式建模高斯近似的失效区间并在低通量时回退到精确分布；把谱截断的 \((L,p)\) 从经验常数改成按 \(\tilde P\) 谱间隙自适应的选择。

评分¶

新颖性: ⭐⭐⭐⭐ 首次给出死时间下光子计数分布的解析刻画（MRP + CLT + 谱截断尾项），真正拆掉了仿真器的速度-精度两难。
实验充分度: ⭐⭐⭐⭐ 仿真精度对比、真实 SPAD 硬件验证、下游深度/反射率、三维效率 scaling 都覆盖了；但真机验证规模偏小，效率对比未含 GPU 基线。
写作质量: ⭐⭐⭐⭐ "计数 × 形状"的问题重构讲得很清楚，物理直觉与公式衔接顺畅；代价是数学密度偏高，主要推导都推给了补充材料。
价值: ⭐⭐⭐⭐ 为学习式 SP-LiDAR 提供了可规模化、物理忠实的数据生成器，且对已有 Markov 链工具是即插即用式的补全。