Fisher 信息：相机标定真实数据实操

本文继续相机内参标定的工作, 之前仿真篇在有真值的数字孪生里算了 CRLB、验证了「在估计处算观测信息」这套方法可信（《Fisher 信息：相机标定仿真实操》，理论见《Fisher 信息与精度下界》）。但仿真篇留了一个缺口：真实标定没有真值，第 6 步的 Monte Carlo 验证做不了。

本文把同一套方法跑在真实采的图上，并用 k 折、独立测试、数据一致性把那步验证补回来——回答：这台相机的内参能定到多准、卡在哪里、是算法问题还是硬件问题。

一、概念与理论基础

标定一台相机，本质是：拍一堆棋盘图 → 检测角点 → 解出内参（焦距/主点/畸变）。这套流程里「估得多准」受几个量约束。本节把全文用到的概念一次性讲清——每个都讲「是什么、为什么、怎么用、和下一个什么关系」，后面步骤直接引用，不再回头解释。

1.1 角点检测噪声：σ_u / σ_v / σ_a / σ_det

OpenCV 的 findChessboardCorners + cornerSubPix 把每个角点定到亚像素，但定不完美——真实角点位置与检测位置之间有个误差 $Δp = (Δu, Δv)$（u 水平、v 垂直）。这个误差是随机的，每张图不一样，它的大小就是标定的「噪声地板」：无论算法多好，都不可能比单次检测的位置误差更准。

• σ_u / σ_v：u 方向和 v 方向各自的标准差——检测器在水平方向抖多大（σ_u）、垂直方向抖多大（σ_v）。静态连拍（相机+棋盘都不动）时，同一角点的「真值」不变，跨帧的位置变化全是噪声，跨帧 std 就是它。
• σ_det：径向噪声 $√(E‖Δp‖²) = √(σ_u² + σ_v²)$——把 u/v 两方向误差合成「一个角点的位置误差 RMS」。它和 OpenCV/MATLAB 的重投影 RMS 是同一刻度（都是径向位置误差），可以直接相除比较。
• σ_a：单轴均方根 $σ_det/√2 = √((σ_u²+σ_v²)/2)$——「平均每个方向的噪声」。只在各向同性 （σ_u=σ_v）时它才严格等于 σ_u、σ_v；Fisher 的各向同性公式 $I=(1/σ_a²)JᵀJ$ 用它。
• 各向同性 / 各向异性：σ_u=σ_v 叫各向同性（噪声在水平/垂直一样大，散布是圆）；σ_u≠σ_v 叫各向异性（散布是椭圆）。各向异性时 σ_a 没有唯一定义，Fisher 要用 $Σ=diag(σ_u², σ_v²)$ （u/v 方向分别加权），而不是 σ_a²I（假设两方向同噪声）。

为什么能直接测 σ_det：相机+棋盘静止连拍，同一角点跨帧的位置变化只由检测噪声引起，跨帧 std 即噪声 std。这是全套流程里唯一不依赖标定、直接测噪声的环节，所以 σ_det 是后面所有判据的分母。

1.2 重投影 RMS 与四个误差源

标定解出内参后，用内参+每图外参（位姿）把 3D 棋盘角点「投影」回像素，和检测到的角点比，有个残差。所有角点残差的 RMS 叫重投影 RMS（cv2.calibrateCamera 返回的 rms）：

RMS = √( (1/(M·N)) Σ_i Σ_j ‖proj_ij − detected_ij‖² )

RMS 不是单一来源，是四个误差源叠加：

• e_model（模型误差）：畸变模型（Brown 多项式 k1-3， p1p2）表达不了真实镜头的全部像差，拟合后剩下的部分。
• e_det（检测噪声）：= σ_det，角点检测的随机误差。
• e_target（靶标误差）：标定板不完美——方格边长制造误差、板面翘曲，使 objectPoints（3D 角点）本身带偏差。
• e_num（数值误差）：≈0，浮点/优化残差。

平方叠加（若独立）：$RMS² ≈ e_model² + σ_det² + e_target²$。

RMS 怎么算、参数哪来的——RMS 是「重投影残差」的均方根。算它需要两样：模型预测的角点位置 proj 和 检测到的角点位置 detected。detected 来自 findChessboardCorners；proj 来自下面的正向投影链（cv2.projectPoints）。 ⚠ 关键：这条链里的参数（内参 + 每图外参）不是预先已知的，而是被 calibrateCamera 解出来的——尤其外参（每张图各自的位姿），它确实来自图像角点的拟合（见 ②）。

① 正向投影链（projectPoints：给定参数 → 预测像素 proj；这是 projectPoints 的真实内部过程）：

1. 板坐标：3D 角点 $X_w$ 在板自身坐标系（Z=0 平面，方格间距 25mm，即 $(i·25, j·25, 0)$ mm）。
2. 外参变换（板→相机）：$X_c = R·X_w + t$，R 由旋转向量 rvec 经 Rodrigues 转成 3×3。—— rvec/tvec 是每张图各自的位姿，未知、待解。
3. 透视投影（相机→归一化像面）：$x = X_c.x / X_c.z, y = Y_c.y / X_c.z$（除以深度，得归一化坐标）。
4. 畸变（归一化→畸变后）：径向 $x_d = x(1 + k1·r² + k2·r⁴ + k3·r⁶)$ 加切向 $+ 2p1·x·y + p2·(r²+2x²)$，$r² = x² + y²$。
5. 内参（畸变后→像素）：$u = fx·x_d + cx, v = fy·y_d + cy$ = 模型预测的像素位置 proj。

② 参数（内参 + 外参）怎么解出来（calibrateCamera，这才是「标定」本身）：内参 9 个 + 每图外参 6 个都是未知数，从「检测角点 ↔ 板角点」的对应关系解出来——

• 初值：Zhang 单应矩阵法——每张图的角点对应给一个 2D↔3D 单应矩阵，分解出该图粗略位姿（rvec, tvec）；多张图的位姿约束线性解出粗略内参 K。
• 精化：LM 非线性最小二乘，同时调全部参数（9 内参 + 每图 6 外参），目标就是让下面的 RMS 最小。高斯噪声下最小二乘 = 最大似然（MLE）。

③ 残差 + RMS（用解出的参数走一遍 ①）：
6. 残差：$e = proj − detected$（每角点「模型预测 − 角点检测」，2D 矢量）。
7. RMS：所有角点所有图 $RMS = √( (1/(M·N)) Σ ‖e‖² )$。

一句话：calibrateCamera 迭代地把（内参 + 每图外参）调到让 RMS 最小 → 收敛后，用最终参数走一遍正向投影链 → 得到报告的 RMS。正向链（①）是 projectPoints 的真实计算；外参不是「输入」而是 ② 里被解的未知数——它从图像角点拟合而来（每张图摆位不同，位姿各异）。

关键关系 RMS vs σ_det：若模型完美（e_model=0）+ 靶标完美（e_target=0），RMS→σ_det（只剩检测噪声）。所以 $RMS/σ_det → 1$ 是「到噪声地板」的判据——标定充分利用了数据，只剩不可避免的检测噪声。 RMS>σ_det 说明有系统误差（e_model 或 e_target 没清零）。本例 RMS=0.14、σ_det=0.046，RMS 是 3 倍 → 有系统误差（§六归因）。

1.3 Fisher 信息与 CRLB

Fisher 信息回答：「这批数据能把参数约束到多准？」约束越强（信息越大），参数估得越准。

• 雅可比 J：每个角点观测对每个参数的偏导 $∂f/∂θ$（数值中心差分就能算，不必手推）。J 大→该参数对观测影响大→容易被约束→信息大。J 的形状 $(M·N·2, 参数数)$，每行一个标量观测方程（u 或 v），每列一个参数。
• 信息矩阵：高斯观测下 $I = Jᵀ Σ⁻¹ J$（Σ=噪声协方差）。各向同性 Σ=σ_a²I 退化成 $(1/σ_a²)JᵀJ$。
• CRLB（Cramér-Rao 下界）：任何无偏估计的协方差 $Cov(θ̂) ⪰ I⁻¹$。第 k 内参的 std 下界 $CRLB_k = √[inv(I_β)]_kk$——「不管用什么算法，该参数 1σ 误差都不可能小于它」。

为什么 CRLB 有用：它给「理论天花板」。真实标定的散布（std_emp）应 ≥ CRLB；$std_emp/CRLB ≈ 1$ 说明达到理论极限（统计有效），$≫1$ 说明有 CRLB 没覆盖的额外误差（系统误差）。

1.4 边缘化外参

参数 θ = 内参（9 个）+ 每图外参（6 个，位姿）。外参是 nuisance（干扰参数）：我们不关心它，但它和内参一起从同一批数据估出来，带不确定度，会「吃掉」内参的表观信息——因为标定器分不清「焦距变一点」和「位姿变一点」对投影的影响，外参不确定性会渗进内参。

直接用全 I 求逆会把外参当「已知」，CRLB 假性乐观（约 √487 倍，仿真篇第 4 步实测）。正确做法：分块 $I = [[A, B], [Bᵀ, C]]$（A=内参×内参， C=外参×外参， B=耦合），用 Schur 补 $I_β = A − B·C⁻¹·Bᵀ$ 扣掉外参不确定性，得 9×9 内参信息矩阵。CRLB 从 $inv(I_β)$ 算。（用 solve(C, Bᵀ) 解方程而非显式 inv(C)，数值更稳。）

1.5 模型选择：bias-variance 权衡与交叉验证

畸变模型加参数（k1 → k2 → k3 → p1p2 → …）能降 e_model（拟合更好），但参数多了会过拟合—— 把检测噪声当畸变学进去，验证集反而变差。这就是 bias-variance 权衡：

$$ \mathrm{Error}_{\text{hold-out}}=\underbrace{\mathrm{Bias}^2(\text{model})}_{\text{模型不够,欠拟合}}+\underbrace{\mathrm{Var}(\hat\theta)}_{\text{参数估计抖动}}+\underbrace{\sigma_{\text{det}}^2}_{\text{检测噪声}} $$

• 参数太少 → Bias² 大（欠拟合，真实畸变表达不了）。
• 参数太多 → Variance 大（过拟合，把噪声当信号）。
• 存在最优阶：再加参数 hold-out 反而升。

交叉验证（k 折 hold-out）：把数据分训练集/验证集，训练集拟合内参、验证集（固定内参解外参）看泛化。 hold-out RMS 不再下降（降幅 < 5%）= 模型饱和，该阶是最优——真实畸变已被表达，再加参数只是过拟合。

1.6 地板判据

「地板判据」把三个量串起来：

$$ \mathrm{RMS}_{\text{hold-out}}\approx\mathrm{RMS}_{\text{in-sample}}\approx\sigma_{\text{det}} $$

三个量近似相等，标定才「达到了噪声地板」：

• $RMS_hold-out ≈ RMS_in-sample$：没过拟合（训练/验证同分布，MLE 一致）。
• $RMS_in-sample ≈ σ_det$：到地板（只剩检测噪声，系统误差清零）。

任何一环断开都说明问题：hold-out≫in-sample 是过拟合；两者≫σ_det 是有系统误差。详见 §七。

1.7 仿真 vs 真实

Fisher 信息 $I(θ)$ 是参数的函数，要在某个求值点上算。两种求值方式：

• 仿真（期望信息）：有真值 θ_true，在真值处算 I、CRLB——这是金标准，但只有仿真算得出（只有仿真有真值）。
• 真实（观测信息）：无真值，在标定估计 θ̂ 处算。σ_a 用重投影 RMS/√2（含模型残差）。
• MLE 一致性：θ̂ → θ_true（大样本），故观测信息 ≈ 期望信息——仿真替你验证过「在估计处算观测信息」和「在真值处算金标准」几乎重合（std/CRLB≈1），真实标定里用观测信息法才站得住。

本报告就是真实版：在 θ̂（8 组标定均值）处算观测信息 → CRLB。

1.8 硬件上限：镜头像差与靶标误差

RMS 里的 e_model 和 e_target，根子在硬件：

• 镜头剩余像差（residual aberration）：真实镜头有球差、彗差、像散、场曲、畸变、色差等像差，是 高阶、非对称、非多项式的。Brown 模型（几个多项式系数）只能拟合「多项式能表达」的部分， 剩余像差 = 拟合后还剩下的部分，装不进任何几个系数，变成 RMS 里的 e_model。
• 靶标误差（target error）：标定假设「棋盘是完美等间距方格」，真实靶标方格边长有制造误差、板面翘曲 → objectPoints 带系统误差 → RMS 里的 e_target。材质影响大：纸/亚克力（便宜，面不平、精度差）vs 玻璃镀铬（贵，面平如镜、格准）。

为什么是「硬件上限」：像差和靶标都是物理硬件，软件/算法改不动。要降 RMS 必须换硬件（剩余像差更小的镜头 + 玻璃镀铬靶）。所以「RMS 卡在某值」常常是硬件天花板，不是算法没做好——这正是本报告要区分的核心：算法做到了它能做的最好，剩下是硬件的事。

理论框架到此建立完毕，核心计算式就三条。观测对参数的灵敏度聚成信息矩阵：

$$ \mathcal I=\mathbf J^\top\boldsymbol\Sigma^{-1}\mathbf J $$

边缘化掉外参（nuisance）后得到纯内参信息（Schur 补）：

$$ \mathcal I_\beta=\mathbf A-\mathbf B\mathbf C^{-1}\mathbf B^\top $$

于是每个内参的理论 std 下界：

$$ \mathrm{CRLB}_k=\sqrt{[\mathcal I_\beta^{-1}]_{kk}} $$

（推导——Fisher 信息定义、score 函数、高斯形式、CRLB 不等式证明、MLE 渐近有效性、信息几何——见《Fisher 信息与精度下界》；本文只讲「怎么用、为什么这么算」。）

读完这节，后面的 10 步实操就是把这些概念串成流水线：测 σ_det（1.1）→ 标定得 RMS（1.2）→算 Fisher/CRLB（1.3-1.4）→ 多折定模型阶（1.5）→ 用地板判据（1.6）检验 → 归因到硬件（1.8）。

二、实验目的、内容与解决流程

2.1 实验目的

拿到这台相机的标定结果，表面看 RMS=0.14 px 还不错。但「不错」是相对什么？我们要回答三个递进的科学问题：

1. 精度天花板：在当前拍摄条件（棋盘、位姿、图数）下，内参估计理论上最好能定到多准？有没有一个「物理极限」告诉我们别徒劳？
2. 当前差距：真实标定的散布，达到这个理论极限了吗？差距来自哪里——估计算法不够好（流程问题），还是畸变模型不够（模型问题），还是镜头/靶标本身（硬件问题）？
3. 优化方向：如果没达到极限，该往哪使劲？多拍图？换模型？换镜头？还是固定棋盘？

这三个问题串成一条主线：先用 Fisher 信息算出理论极限（CRLB），再用真实散布和它对比暴露 gap，最后用多折 + 剔除实验归因 gap 的来源。

2.2 实验内容

围绕这三个问题，做 5 个核心实验（对应 §五 的基础 5 步 ①-⑤；⑥-⑩ 是后续的扩展验证）：

实验	做什么	回答哪个问题
① 噪声测量	静态连拍测 σ_det（§1.1）	精度地板——所有判据的分母
② 8 组重复标定	8 组标定 θ̂ + 散布 std_emp（§1.2）	当前精度（经验不确定度）
③ Fisher CRLB	θ̂ 处算观测信息 → CRLB（§1.3-1.4）	理论天花板（问题 1）
④ 多折交叉验证	逐级畸变 + k 折 hold-out（§1.5）	模型阶数——排除「模型不够」（问题 2 的候选 a）
⑤ 综合判据	九条判据 + 地板判据 + 剔微动（§1.6）	归因 + 结论（问题 2、3）

2.3 初步结果与暴露的问题

跑完前两步（①②）就拿到一个反直觉的矛盾，这也是整篇报告要解的谜题：

• 表面看好：RMS=0.14 px，8 组标定极稳（变异 0.15%），fx/fy/k1 估得高度一致——标定质量似乎不错。
• 但暴露问题：σ_det=0.046 px，RMS 是噪声地板的 3 倍。按 §1.2 的 $RMS² ≈ e_model² + σ_det² + e_target²$，这意味着 RMS 里有 $0.14² − 0.046² ≈ 0.017$ 的「额外方差」不是检测噪声——是系统误差。

矛盾点：标定「稳」（8 组一致）却「没到地板」（RMS≫σ_det）。中间这段差距是谁造成的？光看 RMS 这个数分不清，因为有三种可能，且处理方式完全不同：

• （a） 模型不够：畸变阶数太低（只 5 参数 k1-3+p1p2），真实畸变表达不了 → e_model 大。 解法：加畸变参数（k4-6， s1s2），重标定。
• （b） 棋盘微动：采图时棋盘没完全静止，「静止连拍」存在位移，被当成系统误差。 解法：固定棋盘，重新采集。
• （c） 硬件天花板：镜头剩余像差 + 靶标不完美（§1.8）→ 物理决定的 e_model + e_target。 解法：换更好镜头 + 玻璃镀铬靶（换硬件）。

所以必须先归因、再优化，否则徒劳（如果实际是硬件问题，加再多畸变参数也没用）。

2.4 解决流程

要区分这三种原因，需要一个理论基准作参照——这就是 Fisher 信息的价值：它给出「给定数据+模型+噪声，理论最优能到多准」（CRLB），和真实散布一比，gap 立刻暴露。整个流程围绕「理论基准 vs 真实表现」的对比展开，5 步各解决一个子问题，前一步的输出文件喂给后一步，形成完整证据链：

flowchart LR
    subgraph 数据
        S[/"cal_sigama/ 5组×80张静止连拍"/]
        D[/"cal_data/ 80位姿×8帧"/]
    end
    S --> M1["① measure_noise 测 σ_det"] --> N[("noise.json")]
    D --> M2["② repeat_calib 8组标定"] --> R[("repeat.json")]
    D --> M3["③ fisher_eval CRLB"]
    N --> M3
    R --> M3
    M3 --> C[("crlb.json")]
    D --> M4["④ cross_validate 多折定阶"] --> CV[("cv.json")]
    N --> M5["⑤ report 九条判据"]
    R --> M5
    C --> M5
    CV --> M5
    M5 --> OUT[("report.json")]

逐步逻辑（每步解决什么、为什么必须先做它）：

1. 测 σ_det（地板）：独立测出检测噪声，作后面所有判据的分母。没有它，RMS 多少都无从判断好坏—— 0.14 是好是坏？要看它离 σ_det 多远。
2. 8 组重复标定（经验）：拿到真实散布 std_emp（经验不确定度）。这是「真实标定实际有多准」。
3. Fisher CRLB（理论）：在 θ̂ 处算观测信息 → CRLB（理论不确定度）。这是「理论上最好能多准」。 → 把 2 和 3 比：$std_emp/CRLB ≈ 1$ 达到极限；$≫1$ 有 CRLB 没覆盖的额外误差（系统误差）。
4. 多折交叉验证（归因 a）：逐级畸变 + k 折，看 hold-out RMS 平台在哪。排除或确认「模型不够」（候选 a）——若 k1-3 已饱和，加阶降不动 RMS，说明不是模型的锅。
5. 综合判据（归因 b/c + 结论）：九条判据 + 地板判据 + 剔微动实验。排除或确认「微动」（候选 b）—— 若剔微动后 RMS 没降，说明不是微动；剩下的就是硬件（候选 c）。

最终把「RMS 为什么 > σ_det」的三候选逐一排除，定位到唯一原因，给出可靠的优化方向。后面的展开顺序：§三（数据采集）→ §四（先定义一个贯穿全程的独立测试工具）→ §五（10 步实操，把这条流程逐步落地）。

三、数据采集与有效性检查

3.1 采集设备与两套数据的思路

设备：海康工业相机（2560×1440），通过海康 HCNetSDK（libhcnetsdk.so）用 Python ctypes 调用，相机 IP 192.168.2.2:8000，用户 admin。标定板：11×8 内部角点（↔ 12×9 方格，88 角点），方格边长 25 mm（实测准确；只定外参距离尺度，与 fx/fy/cx/cy 等像素量内参无关，见 §9.4）。

针对 §2.1 的科学问题，采三套数据：

数据集	拍法	科学目的	要求
cal_sigama/	相机+棋盘都静止，连拍	测 σ_det（§1.1）	绝对不动，让跨帧位置变化只由检测噪声引起
cal_data/	相机静止、棋盘摆不同位姿	标定 + 8 组重复（§1.7）	位姿多样 + 整板入画 + 每位姿连拍 8 帧
cal_test/	相机静止、棋盘不同位姿（独立采）	独立测试——泛化金标准（§五第 9 步）	不参与训练，仅测试；64 张 SDK 直采（与训练同源）

标定数据为什么这么要求（位姿多样性的几何原因仿真篇第 1 步详述）：

• 位姿多样（位置九宫格、俯仰/偏航/滚转、远近）：打破参数耦合退化（如 fx/cy 耦合），让 Fisher 信息矩阵条件数可控——否则求 CRLB 会出 NaN。
• 整板入画（所有角点在画面内）：角点出画面会把畸变半径撑爆、λ_max 飙到 1e14，求逆爆炸。
• 每位姿连拍 8 帧：8 帧 = 8 个「同位姿、不同噪声实例」的样本，重组为 8 重复组（按第 i 帧抽），做真实版 Monte Carlo（替代仿真换噪声种子）。

3.2 采集脚本

采集分两层：SDK 封装（登录一次连拍多张）+ 批量采集（CLI）。

SDK 抓图封装——把海康 SDK 的「加载库/登录/抓图/登出/清理」封装成 HikGrabber 类， 登录一次、连拍多张（原 demo 每张都登录登出，每次几百 ms；封装后抓一张只 ~150 ms）。用 with 自动登出清理。完整源码：

#!/usr/bin/env python3
"""海康 HCNetSDK 抓图封装(可复用模块)。运行前:export LD_LIBRARY_PATH=<sdk>/lib"""
import ctypes, os, time

SDK_LIB_DIR = ".../EN-HCNetSDK.../lib"
DEFAULT_HOST, DEFAULT_PORT, DEFAULT_USER, DEFAULT_PWD = "192.168.2.2", 8000, "admin", "..."
CHANNEL, PIC_SIZE, PIC_QUALITY = 1, 0xff, 0     # 通道1;0xff=当前码流分辨率;0=最好质量

class _NET_DVR_LOCAL_SDK_PATH(ctypes.Structure):
    _fields_ = [("sPath", ctypes.c_char * 256), ("byRes", ctypes.c_ubyte * 128)]

class _NET_DVR_USER_LOGIN_INFO(ctypes.Structure):
    _fields_ = [("sDeviceAddress", ctypes.c_char * 129), ("byUseTransport", ctypes.c_ubyte),
                ("wPort", ctypes.c_uint16), ("sUserName", ctypes.c_char * 64),
                ("sPassword", ctypes.c_char * 64), ("cbLoginResult", ctypes.c_void_p),
                ("pUser", ctypes.c_void_p), ("bUseAsynLogin", ctypes.c_uint32),
                ("byProxyType", ctypes.c_ubyte), ("byUseUTCTime", ctypes.c_ubyte),
                ("byLoginMode", ctypes.c_ubyte), ("byHttps", ctypes.c_ubyte),
                ("iProxyID", ctypes.c_int32), ("byVerifyMode", ctypes.c_ubyte),
                ("byRes3", ctypes.c_ubyte * 119)]

class _NET_DVR_JPEGPARA(ctypes.Structure):
    _fields_ = [("wPicSize", ctypes.c_uint16), ("wPicQuality", ctypes.c_uint16)]

class HikGrabber:
    """海康相机抓图器:登录一次,可连续抓多张 JPEG。
        with HikGrabber(host="192.168.2.2") as g:
            g.capture("a.jpg"); g.capture("b.jpg")
    """
    def __init__(self, host=DEFAULT_HOST, port=DEFAULT_PORT, user=DEFAULT_USER,
                 pwd=DEFAULT_PWD, channel=CHANNEL, pic_size=PIC_SIZE,
                 pic_quality=PIC_QUALITY, lib_dir=SDK_LIB_DIR):
        self._channel, self._pic_size, self._pic_quality = channel, pic_size, pic_quality
        self._lib_dir, self._uid, self._sdk = lib_dir, -1, None
        self._load(lib_dir); self._proto(); self._set_sdk_paths()
        if not self._sdk.NET_DVR_Init():
            self._raise_err("NET_DVR_Init 失败")
        self._sdk.NET_DVR_SetConnectTime(5000, 3)
        self._login(host, port, user, pwd)

    def _load(self, lib_dir):
        self._sdk = ctypes.CDLL(os.path.join(lib_dir, "libhcnetsdk.so"))

    def _proto(self):
        s = self._sdk
        s.NET_DVR_Init.restype = ctypes.c_bool
        s.NET_DVR_Login_V40.restype = ctypes.c_long
        s.NET_DVR_Login_V40.argtypes = [ctypes.c_void_p, ctypes.c_void_p]
        s.NET_DVR_CaptureJPEGPicture.restype = ctypes.c_bool
        s.NET_DVR_CaptureJPEGPicture.argtypes = [ctypes.c_long, ctypes.c_long,
                                                  ctypes.c_void_p, ctypes.c_char_p]
        s.NET_DVR_GetLastError.restype = ctypes.c_uint32

    def _set_sdk_paths(self):
        """告诉 SDK 组件库(HCNetSDKCom)与加解密库(libcrypto/libssl)所在路径。"""
        d = self._lib_dir
        p = _NET_DVR_LOCAL_SDK_PATH(); p.sPath = d.encode()
        self._sdk.NET_DVR_SetSDKInitCfg(2, ctypes.byref(p))                       # SDK_PATH
        pe = _NET_DVR_LOCAL_SDK_PATH(); pe.sPath = os.path.join(d, "libcrypto.so.1.1").encode()
        self._sdk.NET_DVR_SetSDKInitCfg(3, ctypes.byref(pe))                      # LIBEAY_PATH
        ps = _NET_DVR_LOCAL_SDK_PATH(); ps.sPath = os.path.join(d, "libssl.so.1.1").encode()
        self._sdk.NET_DVR_SetSDKInitCfg(4, ctypes.byref(ps))                      # SSLEAY_PATH

    def _login(self, host, port, user, pwd):
        info = _NET_DVR_USER_LOGIN_INFO()
        info.sDeviceAddress = host.encode(); info.wPort = port
        info.sUserName = user.encode(); info.sPassword = pwd.encode()
        info.byLoginMode = 0; info.byHttps = 0       # 0-Private(8000 私有协议)
        devbuf = ctypes.create_string_buffer(4096)   # NET_DVR_DEVICEINFO_V40,仅接收不解析
        uid = self._sdk.NET_DVR_Login_V40(ctypes.byref(info), devbuf)
        if uid < 0:
            self._sdk.NET_DVR_Cleanup(); self._raise_err("登录失败")
        self._uid = uid

    def _raise_err(self, msg):
        raise RuntimeError(f"{msg}  (NET_DVR_GetLastError={self._sdk.NET_DVR_GetLastError()})")

    def capture(self, out_path) -> float:
        """抓一张 JPEG 到 out_path,返回耗时(ms)。失败抛 RuntimeError。"""
        jpg = _NET_DVR_JPEGPARA(); jpg.wPicSize = self._pic_size; jpg.wPicQuality = self._pic_quality
        t0 = time.monotonic()
        ok = self._sdk.NET_DVR_CaptureJPEGPicture(
            self._uid, self._channel, ctypes.byref(jpg), str(out_path).encode())
        dt = (time.monotonic() - t0) * 1000
        if not ok:
            self._raise_err("抓图失败")
        return dt

    def close(self):
        if self._uid >= 0:
            self._sdk.NET_DVR_Logout(self._uid); self._uid = -1
        self._sdk.NET_DVR_Cleanup()

    def __enter__(self): return self
    def __exit__(self, *exc): self.close(); return False

批量采集——每批一个时间戳文件夹、每张图用时间戳（到微秒）命名，保证同批次唯一、按字典序可排序。支持自动间隔（定时拍，如 -i 1.0）和交互模式（--interactive，手动摆板按回车拍）。完整源码：

#!/usr/bin/env python3
"""批量采集标定图:登录一次连拍 N 张,每批一个时间戳文件夹,每张图时间戳命名。"""
import argparse, os, time
from datetime import datetime
from sdk_grabber import HikGrabber

def ts_folder(): return datetime.now().strftime("%Y%m%d_%H%M%S")          # 批次文件夹(到秒)
def ts_image():  return datetime.now().strftime("%Y%m%d_%H%M%S_%f")       # 图片名(到微秒,唯一)

def main():
    ap = argparse.ArgumentParser(description="批量采集海康相机标定图")
    ap.add_argument("-n", "--count", type=int, default=10, help="采集张数")
    ap.add_argument("-i", "--interval", type=float, default=1.0, help="每张间隔秒数")
    ap.add_argument("--interactive", action="store_true", help="交互模式:每张前按回车(手动摆板)")
    ap.add_argument("-o", "--out-dir", default="captures", help="根输出目录")
    args = ap.parse_args()
    if args.count <= 0: ap.error("--count 必须为正整数")
    if not os.environ.get("LD_LIBRARY_PATH"):
        print("警告:未设置 LD_LIBRARY_PATH,SDK 依赖可能加载失败")

    batch_dir = os.path.join(args.out_dir, ts_folder())
    os.makedirs(batch_dir, exist_ok=True)
    print(f"准备采集 {args.count} 张 -> {batch_dir}")
    print(f"模式:{'交互(按回车拍下一张)' if args.interactive else f'自动间隔 {args.interval}s'}")

    ok = fail = 0
    with HikGrabber() as g:
        print("登录成功,开始抓图 ...")
        for k in range(1, args.count + 1):
            if args.interactive:
                input(f"[{k}/{args.count}] 摆好标定板后按回车抓图 ...")
            path = os.path.join(batch_dir, ts_image() + ".jpg")
            try:
                dt = g.capture(path)
                print(f"  [{k}/{args.count}] OK  {os.path.getsize(path):>7} bytes  {dt:>5.0f} ms"
                      f"  -> {os.path.basename(path)}")
                ok += 1
            except RuntimeError as e:
                print(f"  [{k}/{args.count}] FAIL  {e}"); fail += 1
            if not args.interactive and k < args.count:
                time.sleep(args.interval)
    print(f"\n完成:成功 {ok} / 失败 {fail}")

if __name__ == "__main__":
    main()

采集流程：

采集流程：设好 SDK 库路径 → 在每位姿摆好棋盘、按回车连拍 8 张 → 重复 80 个不同位姿（位置九宫格 + 俯仰/偏航/滚转 + 远近），每次产出一个时间戳文件夹。
噪声数据 cal_sigama 类似，但相机+棋盘都不动，连拍 80 张/组（-n 80 -i 0.5 自动间隔），5 组。

3.3 数据组织与规模

采集后（手动跑了 80 轮标定 + 5 组噪声）：

数据集	结构	总图数	用途
cal_sigama/	5 个时间戳文件夹 × 80 张	400	测 σ_det（每组静止连拍）
cal_data/	80 个时间戳文件夹 × 8 张	640	标定（每位姿连拍 8 帧）
cal_test/	64 张（扁平，不同位姿）	43	独立测试（不参与训练，测泛化精度）

cal_test 测试集说明（独立验证，比交叉验证更严格的金标准）：

• 完全不参与训练：用训练内参在 cal_test 每张图解外参（PnP）+ 重投影 RMS，测泛化精度——完全新的数据，没被模型见过。
• SDK 直采（与训练同源，文件大小~650KB≈cal_data 616KB）：64/65 检测成功（65 张文件，1 张未检出），但测试 RMS=0.27 px（≫ 训练 0.14）→ 绝对值偏高含 cal_test 泛化 gap，不代表内参变差。
• 一致性校验不适用（§3.4）：cal_test 是不同位姿的测试集（非静止连拍），run_detect 的 σ_det=439 是位姿差异误报，不是质量问题。
• 相对比较有效（§9.4）：各组内参（8 组 / 10 折 / FULL）在 cal_test 上 RMS 都 ≈0.27（差在第 4 位）→ 内参可靠，0.27 是 cal_test 泛化 gap（训练0.14→测试0.27）。

目录结构：

capture/
├── cal_sigama/20260717_165330/  (第1组静止连拍 80 张)  ... (5 组)
└── cal_data/20260717_172523/    (位姿1,连拍 8 帧)     ... (80 位姿)

cal_data 重组：80 位姿 × 8 帧 → 按第 i 帧抽成 8 个重复组（第 i 组 = 各文件夹的第 i 帧），每组 80 张覆盖同样 80 位姿、只差噪声实例 → 真实版 Monte Carlo（§1.7）。

3.4 数据有效性检查

拍完不能直接用，要检查每张图棋盘检测成功 + 角点位置合理。用检测脚本（递归检测 + 一致性校验）：

# detect_corners 核心:检测 + 亚像素细化
gray = cv2.cvtColor(img, cv2.COLOR_BGR2GRAY)
found, corners = cv2.findChessboardCorners(gray, board.pattern_size,
                  flags=cv2.CALIB_CB_ADAPTIVE_THRESH | cv2.CALIB_CB_NORMALIZE_IMAGE
                       | cv2.CALIB_CB_FAST_CHECK)     # FAST_CHECK 大图加速(90ms vs 2.5s)
if found:
    corners = cv2.cornerSubPix(gray, corners, (11,11), (-1,-1),
                  (cv2.TERM_CRITERIA_EPS|cv2.TERM_CRITERIA_MAX_ITER, 50, 1e-3))  # 亚像素

# 一致性校验核心(noise_stats):同组角点位置应几乎一致
P = np.stack([c.reshape(-1,2) for c in group_corners])     # (M=8, N=88, 2)
sigma_axis = P.std(axis=0, ddof=1)                          # 每角点每轴 std
sigma_det = np.sqrt(np.mean(sigma_axis[:,0]**2) + np.mean(sigma_axis[:,1]**2))
rms_img = np.sqrt(((P - P.mean(0))**2).sum(2).mean(1))      # 每图 RMS 偏移
outliers = [img for img, rms in zip(imgs, rms_img) if rms > 3*np.median(rms_img)]  # 离群

检查三件事：

1. 检测成功率：每张图能否找到 88 角点（整板入画 + 对焦 + 光照）。
2. 亚像素细化：cornerSubPix 把角点定到亚像素（标定必须，否则掉精度）。
3. 组内一致性：同组（静止连拍）角点位置应几乎一致；离群（RMS>3×median）说明该帧棋盘动/模糊。

对每个数据目录跑一遍递归角点检测 + 组内一致性校验 + 可视化（cal_data 标定图、cal_sigama 噪声图各一次）。

3.5 检查结果

cal_data（640 张标定图）：

• 检测成功率：640/640 全部成功（88 角点全找到，整板入画）——采集质量好。
• 组内一致性（80 组，每组 8 帧）：σ_det 分布——
  • 56 组（70%）σ_det < 0.08 px（真静止，角点位置几乎不变）
  • 8 组（10%）σ_det 0.08-0.2（轻微动）
  • 16 组（20%）σ_det > 0.2（明显微动——采这些位姿时棋盘没完全稳住）
  • 中位 0.043，均值 0.137（被微动组长尾拉高）
• 离群图：2 张（RMS > 3×median，棋盘整体抖 ~0.8 px；逐角点分析显示 88 角点均匀偏移 → 物理抖动，非检测错位）。

cal_sigama（400 张噪声图）：

• 检测成功率：400/400 全成功。
• 5 组 σ_det：[0.037， 0.036， 0.022， 0.066， 0.048]，pooled=0.046 px（§五第 1 步主值）。
• σ_v/σ_u ≈ 1.9（各向异性，垂直噪声近 2 倍水平）。

cal_test（65 张 → 64 张有效）：

• 检测成功率：64/65（65 张文件，1 张检测失败；64 张全找到 88 角点）。
• 一致性校验不适用：cal_test 是不同位姿测试集（非静止连拍），run_detect 的 σ_det=439 px 是位姿差异误报（64 张位姿各异），不是检测噪声或质量问题。
• 质量评估：64 张有效检测 → 棋盘可见、对焦 OK（检测层面可用）；但测试 RMS=0.27 px（≫训练 0.14）→ 泛化 gap（测试 RMS 0.27 > 训练 0.14）（§五第 9 步详析：各组内参在 cal_test 上 RMS 一致 → 内参可靠，0.27 是 cal_test 自身泛化 gap（训练0.14→测试0.27））。

结论：三套数据（标定 cal_data + 噪声 cal_sigama + 测试 cal_test）检测全部成功（100%），整体有效；cal_data 有 20% 微动组——这本身是 §六归因要验证的（微动是否是 RMS>σ_det 的主因；spoiler：剔微动后 RMS 没降，不是主因，是硬件）。

数据示例（原始棋盘 + 角点检测标注）：

上：原始采集图（cal_data 2560×1440）；下：findChessboardCorners + cornerSubPix 亚像素检测的 88 角点（11×8 棋盘，彩色连线+序号）。cal_data 640/640、cal_test 64/65 检测成功。

四、独立测试方法

本节独立于实验流程，定义一个可复用的评估工具：给定任意内参 θ̂，在 cal_test（64 张，不参与训练）
上算泛化精度 + 多角度指标 + 光流图。实验任何步骤（② 8 组 / ⑥ 10 折 / ⑧ 最终内参）得到 θ̂ 后，都可随时调用它验证泛化——不等流程走完。

工具：评估脚本（evaluate）+ 测试脚本（test_calib） 输入：任意内参 θ̂（json）+ cal_test 64 张 输出：多角度指标（RMS / RMS_u / RMS_v / bias / std / per-image 分布 / 边缘中心 / max）+ 每组光流图

用法：输入任意内参 θ̂ + cal_test，对每组内参（repeat 8 组 / cv5 / full）分别算多角度指标 + 出光流图。

指标解读：

• RMS：泛化精度（各组一致 → 内参可靠；绝对偏高 → 泛化 gap）。$RMS = √(mean‖e‖²)$，无方向的散布度量。
• bias_u/v：残差均值 $bias_u = mean(e_u)$、$bias_v = mean(e_v)$（$e = proj − detected$，所有测试角点）。bias 是有方向的系统偏移：若内参把角点整体往某方向挪（如 fx 偏大 → 投影整体偏右），bias 会是非零常数。本例 bias_u≈+0.0000、bias_v≈−0.0002 px → 正负抵消、均值归零 → 无系统偏置（内参投影居中正确，RMS 0.27 是散布不是偏移）。对比 RMS 0.27 ≫ bias ≈0：随机散布大但系统偏移为零。
• std_u/v：随机散布（v>u → 各向异性）
• RMS/σ_det：≫1 → 泛化 gap（但相对比较仍有效——各组吃同样测试噪声）
• 光流图：系统误差空间分布（边缘 vs 中心、方向性）

本例结果（各组内参在 cal_test 上，64 张）：

内参来源	RMS	bias_u	bias_v	std_u	std_v	RMS/σ	max
repeat 8 组	~0.27	≈0	-0.0001	~0.42	~0.44	~5.9	~0.72
FULL（k=5）	0.27	≈0	-0.0001	0.42	0.44	5.9	3.3

各组 RMS 一致（差第 4 位）→ 内参可靠；0.27 是 cal_test 泛化 gap（训练0.14→测试0.27）（非内参问题）；bias≈0（无系统偏置）。

五、10 步实操

概念（§一）、数据（§三）、独立测试工具（§四）都已就位，下面把整条评估流程拆成 10 步落地： ①-⑤ 回答核心三问（地板在哪 / 差距多少 / 朝哪使劲），⑥-⑩ 从可靠性、数据量、泛化、一致性等更多角度反复验证同一个结论。

概述:流程、目的、输入输出与理论连贯性

怎么做：把「标定能到多准、卡在哪里」拆成 5 个基础步 + 5 个扩展验证步，每步的输出文件喂给下一步，形成完整证据链——不是孤立实验，而是一环扣一环：前一步的结论是后一步的前提。最终把 RMS>σ_det 的三候选（模型/微动/硬件）逐一排除，定位到唯一原因，并以数据一致性收尾确认内参是物理常量。

每步的目的（回答什么）、输入、输出：

步	目的（回答什么问题）	输入	输出
① 测噪声	精度地板在哪？（分母）	cal_sigama 静止连拍	σ_det（noise.json）
② 8 组标定	当前实际有多准？（经验不确定度）	cal_data 80 位姿×8 帧	θ̂ + std_emp（repeat.json）
③ Fisher CRLB	理论上能多准？（精度天花板）	θ̂ + σ_det + 数据几何	CRLB（crlb.json）
④ 多折定阶	模型阶数够不够？（排除模型不够）	cal_data 8 组	hold-out 平台（cv.json）
⑤ 综合判定	归因 + 结论（gap 是谁？）	全部 json	九条判据 + 地板判据（report.json）
⑥ 5 参数 10 折	可靠性（θ̂ 跨折散布）	cal_data 8 组×10 折	θ̂ 散布（cv5.json）
⑦ 精度 vs 数据量	加图能否提升？（硬件 vs 数据）	cal_data k=1-5 帧	rms/CRLB vs k（sub.json）
⑧ 最终内参 + 光流	整合（三路）+ 误差空间分布	全部 json + cal_data	θ̂_final + 光流图（full.json）
⑨ 独立测试	泛化精度（cal_test 金标准）	cal_test 64 张	测试 RMS + 多角度指标 + 光流
⑩ 数据一致性	内参是否数据无关（物理常量？）	train/test/mix 三路	Δθ（consistency.json）

理论怎么连贯（每步用 §一 的哪个概念）：

① σ_det(§1.1 噪声地板)→ 分母
   ↓
② θ̂ + std_emp(§1.2 标定)→ 经验精度
   ↓
③ CRLB(§1.3-1.4 Fisher,θ̂ 处观测信息)→ 理论天花板
   ↓ ② vs ③:std/CRLB≈1 有效;≫1 有额外误差
④ 多折定阶(§1.5 bias-variance)→ 选 5 参数(k1-3+p1p2 范用)
   ↓
⑤ 综合判定(§1.6 地板判据 + 剔微动)→ 归因:模型(排除)→微动(排除)→硬件
   ↓ ============ 扩展验证(⑥-⑩)============
⑥ 5 参数 10 折(可靠性)→ std/CRLB≈5,标定对位姿子集敏感
   ↓
⑦ 精度 vs 数据量(§1.8 硬件)→ rms 不降、CRLB∝1/√k 降→硬件上限
   ↓
⑧ 最终内参 + 光流(三路整合)→ θ̂_final + 系统误差空间分布(边缘/v>u/分块趋势)
   ↓
⑨ 独立测试(cal_test 金标准)→ 泛化可靠(各组 RMS 一致 + bias≈0)
   ↓
⑩ 数据一致性(train/test/mix)→ Δfx<1px → 内参是数据无关的物理常量

每步用前一步的输出文件 + §一的概念，逻辑闭合：没有 σ_det①，RMS② 多少都无从判断好坏；没有 CRLB③，散布② 多大都不知是否到极限；没有多折④，不知 RMS>σ_det 是不是模型问题。

关键设计：所有「重复」实验（②③④⑤）共用固定折法（seed=42，第 4/6 步同法）——8 个重复组用同一套折索引，组间只差噪声实例 → 散布=纯噪声贡献，公平可比（若折法不同，散布混入位姿覆盖差异）。

第 1 步｜测噪声 σ_det

1.1 数据：静态连拍

输入 cal_sigama/，5 组，每组 80 张静止连拍（共 400 张）：

• 拍法：相机 + 棋盘都不动，连拍 80 帧。同一角点的「真值」位置不变 → 跨帧的位置变化全是检测噪声。
• 5 组：不同位姿/区域重复采，看稳定性 + pooled 合并（消除组间位姿差）。
• 为什么静止能测噪声：角点检测器（findChessboardCorners + cornerSubPix 亚像素）每次检测有随机误差 $Δp=(Δu,Δv)$。静止连拍时，同一角点 80 次检测的散布 = 检测器噪声 std。组均值 μ（80 次平均）是真值最佳估计（噪声被 √80 压低），残差 $P−μ$ 即每次的噪声实例。

1.2 原理:σ_det 的连拍计算

对每组（M=80 张 × N=88 角点 × 2 轴）：

1. 堆张量 P，（M， N， 2）：$P[m,j,:] = 第 m 张图第 j 角点的 (u,v)$。
2. 组均值（真值估计）：$μ_j = (1/M) Σ_m P[m,j,:]$，噪声被 √M 压低。
3. 每角点每轴 std（ddof=1）：$sigma_axis[j,k] = std_m(P[m,j,k])$，跨 80 帧 = 该角点该轴噪声。
4. 单轴方差（E[Δu²]）：$var_u = mean_j(sigma_axis[j,0]²)$——mean of variance（平方后平均），不是 mean of std（Jensen 不等式会低估）。
5. 径向噪声：$σ_det = √(var_u + var_v) = √(E[Δu²]+E[Δv²]) = √(E‖Δp‖²)$。

实现细节（§1.1）：

• $σ_det = √(var_u+var_v)$，不是 $(σ_u+σ_v)/2·√2$（算术平均低估 ~5%）。
• $var = mean(std²)$ 不是 $mean(std)²$（Jensen）。
• std 用 ddof=1（除以 M−1）。
• 各向异性（σ_v≠σ_u）时 σ_a=√（（var_u+var_v）/2）=σ_det/√2 只是「单轴均方根」（径向关系）；Fisher 用 $Σ=diag(σ_u²,σ_v²)$。

1.3 完整源码

源码：

def noise_stats(corners_list):
    """一组静止连拍角点 → 噪声统计。"""
    P = np.stack([np.asarray(c).reshape(-1, 2) for c in corners_list])  # (M,N,2)
    M, N, _ = P.shape
    mu = P.mean(axis=0)                    # (N,2) 组均值 = 真值估计(噪声被 √M 压低)
    sigma_axis = P.std(axis=0, ddof=1)     # (N,2) 每角点每轴 std(ddof=1)
    # E[Δu²] = mean of per-corner variance(不是 mean of std —— Jensen 会低估)
    var_u = float(np.mean(sigma_axis[:, 0] ** 2))   # E[Δu²]
    var_v = float(np.mean(sigma_axis[:, 1] ** 2))   # E[Δv²]
    sigma_u = float(np.sqrt(var_u))
    sigma_v = float(np.sqrt(var_v))
    sigma_a_mean = float(np.sqrt((var_u + var_v) / 2))   # 单轴均方根 = σ_det/√2
    sigma_det = float(np.sqrt(var_u + var_v))            # √(E‖Δp‖²)
    diff = P - mu
    rms_img = np.sqrt((diff ** 2).sum(axis=2).mean(axis=1))   # 每图径向 RMS(离群判定用)
    return dict(M=M, N=N, sigma_u=sigma_u, sigma_v=sigma_v,
                sigma_a_mean=sigma_a_mean, sigma_det=sigma_det, ...)

def pooled_noise(groups_corners):
    """多组残差合并(扣各自组均值消除位姿差)。"""
    resid = []
    for cl in groups_corners:
        P = np.stack([np.asarray(c).reshape(-1, 2) for c in cl])   # (M,N,2)
        resid.append(P - P.mean(axis=0))    # 扣该组均值 → 纯噪声残差
    R = np.concatenate(resid, axis=0)        # (ΣM, N, 2)
    var_u = float(np.var(R[..., 0], ddof=1))   # E[Δu²] pooled
    var_v = float(np.var(R[..., 1], ddof=1))
    return dict(sigma_a_pooled=np.sqrt((var_u+var_v)/2),
                sigma_det_pooled=np.sqrt(var_u + var_v),   # √(var_u+var_v)
                sigma_u_pooled=np.sqrt(var_u), sigma_v_pooled=np.sqrt(var_v), ...)

CLI（5 组逐组 + pooled 合并）：

for g in groups:                                     # 5 组静止连拍
    corners = [detect_corners(p, board).corners for p in g]   # 每组 80 张
    st = noise_stats(corners)                        # 逐组 σ_u/σ_v/σ_det
    all_corners.append(corners)
    per_group.append((g.name, st))
pooled = pooled_noise(all_corners)                   # 5 组残差合并 → σ_det_pooled(主值)

1.4 结果：5 组 + pooled

逐组（每组 80 张静止连拍）：

组	σ_u	σ_v	σ_det （px）	u/v 各向异性
165330	0.016	0.034	0.037	⚠ v=2.1×u
165943	0.012	0.036	0.034	⚠ v=3.0×u
171346	0.014	0.015	0.022	✓ 各向同性
171554	0.033	0.054	0.066	⚠ v=1.6×u
172135	0.021	0.043	0.048	⚠ v=2.0×u
pooled（400 张）	0.021	0.040	0.046	σ_v≈1.9×σ_u

组间一致性：σ_det max/min≈3（组 171346 最稳 0.022，组 171554 最乱 0.066），但都静止级（<0.1）。

pooled（5 组残差合并，扣各自组均值，消除位姿差）：

• σ_a=0.032， σ_det=0.046 px（报告主值，喂后续所有判据的分母）
• σ_u=0.021， σ_v=0.040（σ_v≈1.9×σ_u → 各向异性，Fisher 用 Σ=diag）

验证：

• 落工业 C-mount 玻璃靶的常见区间（σ_det 0.02-0.05）→ 检测器正常。
• vs 仿真默认 σ_det=0.1：实测好 2.2×（CRLB ∝ σ_a，精度潜力相应提升）。
• σ_det=0.046 是物理地板（检测器决定的，与标定无关）。

图：6 个角点（中心 #43 + 边缘 #0/#87 等）跨 80 帧的位置散布（单位毫像素 mpx）。每个红点 = 一帧检测的残差 Δp=检测−组均值，散布大小 = 该角点检测噪声。v 方向（纵）散布 > u 方向（横） → 各向异性（σ_v>σ_u）。边缘角点（#0/#87）散布 ≈ 中心（#43），说明检测噪声在视场内较均匀。

→ σ_det 是第 3 步 Fisher 信息（I=（1/σ_a²）JᵀJ）、第 5 步地板判据（RMS/σ_det）的分母。

第 2 步｜8 组重复标定

2.1 原理：8 组重组（真实版 Monte Carlo）

cal_data 是 80 位姿 × 8 帧。重组：第 i 组 = 每个位姿的第 i 帧（每组 80 张）。 8 组覆盖同样 80 位姿，只有「每位姿取了 8 帧里的哪一帧」不同 = 只有噪声实例不同 → 8 组同分布。各自独立标定（每组 80 张）→ $θ̂_1…θ̂_8$，散布 $std_emp = std(θ̂, ddof=1)$ = 经验不确定度（真实版 Monte Carlo，替代仿真换噪声种子）。 为什么合法：8 组 rms 变异 0.15%（极稳）→ 同位姿假设成立 → 散布只反映噪声，不含位姿差异。

2.2 源码（重组 + 标定）

# 1. 重组:第 i 组 = 各文件夹的第 i 帧
per_dir = [sorted(p for p in g.iterdir() if ...) for g in groups]   # 每文件夹 8 帧
repeats = [[per_dir[d][i] for d in range(n_groups)] for i in range(n_per)]  # 8 组 × 80 张

# 2. 逐组标定
thetas = []
for rep in repeats:                                   # 8 组
    corners = [detect_corners(p, board).corners for p in rep]
    cal = calibrate_camera(board, corners, shape)     # cv2.calibrateCamera → K, dist, rms
    thetas.append(cal.theta)                          # [fx,fy,cx,cy,k1,k2,p1,p2,k3]

# 3. 散布(经验不确定度)
theta_mean = np.array(thetas).mean(axis=0)
theta_std  = np.array(thetas).std(axis=0, ddof=1)

RMS（§1.2）：标定后 $RMS=√(mean‖proj−detected‖²)$，含四误差源（e_model+e_det+e_target+e_num）。

2.3 结果：8 组标定 θ̂ 逐组

组	fx	fy	cx	cy	k1	k2	p1	p2	k3	rms
0	1855.058	1845.492	1269.49	725.02	-0.40515	0.23024	7.4e-05	-2.1e-05	-0.08118	0.1379
1	1855.030	1845.464	1269.56	725.08	-0.40501	0.22989	6.8e-05	-2.6e-05	-0.08091	0.1376
2	1855.090	1845.506	1269.63	725.13	-0.40512	0.23013	6.9e-05	-2.8e-05	-0.08106	0.1383
3	1855.156	1845.575	1269.58	725.10	-0.40514	0.23021	7.0e-05	-2.8e-05	-0.08115	0.1381
4	1855.157	1845.577	1269.61	725.05	-0.40508	0.23005	7.4e-05	-2.8e-05	-0.08102	0.1379
5	1854.971	1845.404	1269.62	725.11	-0.40508	0.23013	6.9e-05	-2.8e-05	-0.08113	0.1380
6	1854.979	1845.406	1269.58	725.04	-0.40509	0.23009	7.3e-05	-2.5e-05	-0.08105	0.1377
7	1854.978	1845.402	1269.66	725.26	-0.40509	0.23014	6.3e-05	-2.7e-05	-0.08111	0.1380
mean	1855.053	1845.478	1269.59	725.10	-0.40509	0.23011	7.0e-05	-2.6e-05	-0.08108	0.1379
std（ddof=1）	0.077	0.072	0.05	0.07	4.4e-5	1.1e-4	3.6e-6	2.3e-6	8.8e-5	2e-4

解读：fx/fy/k1 八组极一致（std/fx=0.004%）→ 标定高度可重复；rms 变异 0.15% → 同位姿假设成立 → std_emp 合法。

2.4 图表：8 组 θ̂ 散布

图：8 组标定的 fx/fy/cx/cy/k1（柱状，相对均值偏差 %）。柱高一致 → 可重复；微小差异 = 噪声实例造成的散布（std_emp）。

→ θ_mean（8 组均值）是第 3 步 Fisher 求值点；theta_std 是经验不确定度，等下和 CRLB 比看有效性。

第 3 步｜Fisher CRLB

3.1 原理:理论精度下界的计算

CRLB 回答「仅凭检测噪声，这组数据最多能把内参定多准」。在标定估计 θ̂ 处：

1. 数值雅可比 J：投影函数 f（θ，外参） 对全部参数（9 内参 + 每图 6 外参）的中心差分偏导。 J: (M·N·2, 9+6·M) = (80×88×2, 489) = (14080, 489)。步长 h=1e-6·max(1,|p|)。
2. 信息矩阵 I = Jᵀ Σ⁻¹ J：各向同性 Σ=σ_a²·I；各向异性时用 Σ=diag(σ_u²,σ_v²)，即 Jᵀ·diag(1/σ_u²,1/σ_v² 交错)·J（实测 σ_v≈2σ_u，见第 1 步）。
3. Schur 边缘化外参（ nuisance 参数，§1.4）：外参每图都不同、不关心，但要扣掉它对内参不确定度的"稀释"。分块 I=[[A,B],[Bᵀ,C]]（A=内参块，C=外参块）： I_β = A − B·C⁻¹·Bᵀ（9×9），这才是纯内参的有效信息。
4. CRLB = √|diag（I_β⁻¹）|：每个内参的理论 std 下界。

两个易错点（踩过坑）：（a） 求值点用 theta_mean（8 组均值）并在该点重解外参 （solve_extrinsics_for_theta），不用第 0 组标定的外参——否则求值点偏；（b） 各向异性 CRLB 的 σ_u/σ_v 用标定残差的实际 u/v std（从 projectPoints−detected 实算），不是把噪声 σ 按 rms/noise 缩放。 J 可信性交叉验证：板角点 ∂u/∂k3 的解析式 fx·x·r⁶ 与数值 J 对应元素比对，相对差<5% 才信 CRLB。

3.2 源码（雅可比 + 信息矩阵 + J 验证，与 pipeline/step2-5 同构）

def numerical_jacobian(theta_intr, rvecs, tvecs, board):
    """数值雅可比(中心差分)。J: (n_img*n_pts*2, 9+6*n_img)。"""
    n_img, n_pts = len(rvecs), board.shape[0]
    n_param = 9 + 6 * n_img
    params = np.empty(n_param)
    params[:9] = theta_intr
    for i in range(n_img):                       # 外参拼进参数向量
        params[9+6*i:9+6*i+3] = np.asarray(rvecs[i]).flatten()
        params[9+6*i+3:9+6*i+6] = np.asarray(tvecs[i]).flatten()
    def f_all(p):                                # 投影函数:全图全角点 → (n_img,n_pts,2)
        fx, fy, cx, cy = p[:4]; dist = p[4:9]
        K = np.array([[fx,0,cx],[0,fy,cy],[0,0,1.]])
        out = np.empty((n_img, n_pts, 2))
        for i in range(n_img):
            rv, tv = p[9+6*i:9+6*i+3], p[9+6*i+3:9+6*i+6]
            pts, _ = cv2.projectPoints(board.reshape(-1,1,3), rv.reshape(3,1),
                                       tv.reshape(3,1), K, dist)
            out[i] = pts.reshape(-1, 2)
        return out
    J = np.empty((n_img*n_pts*2, n_param))
    for k in range(n_param):
        h = 1e-6 * max(1.0, abs(params[k]))
        pp, pm = params.copy(), params.copy(); pp[k] += h; pm[k] -= h
        J[:, k] = (f_all(pp) - f_all(pm)).reshape(-1) / (2*h)   # 中心差分
    return J

def crlb_from_J(J, sigma_a=None, sigma_u=None, sigma_v=None):
    """I = JᵀΣ⁻¹J,Schur 边缘化外参 → CRLB=√|diag(inv(I_β))|。"""
    M = J.shape[0]
    if sigma_u is not None:                      # 各向异性:u/v 行交替加权
        w = np.empty(M); w[0::2] = 1/sigma_u**2; w[1::2] = 1/sigma_v**2
    else:                                        # 各向同性
        w = np.full(M, 1/sigma_a**2)
    I = J.T @ (w[:, None] * J)
    A, B, C = I[:9,:9], I[:9,9:], I[9:,9:]
    Ibeta = A - B @ np.linalg.solve(C, B.T)      # Schur 补:9×9 纯内参信息
    return np.sqrt(np.abs(np.diag(np.linalg.inv(Ibeta)))), float(np.linalg.cond(Ibeta))

def jacobian_check(J, theta, rvecs, tvecs, board, img_idx=0, corner_idx=None):
    """交叉验证:角点 ∂u/∂k3,解析 fx·x·r⁶ vs 数值 J[row, k3列]。"""
    if corner_idx is None: corner_idx = board.shape[0]-1      # 板角(r 大→灵敏)
    R0, _ = cv2.Rodrigues(np.asarray(rvecs[img_idx]).reshape(3,1))
    Xc = R0 @ board[corner_idx] + np.asarray(tvecs[img_idx]).flatten()
    x, y = Xc[0]/Xc[2], Xc[1]/Xc[2]; r2 = x*x + y*y
    ana = theta[0] * x * r2**3                   # fx·x·r⁶
    row = (img_idx*board.shape[0] + corner_idx)*2 + 0
    return ana, float(J[row, 8])                 # 解析 vs 数值

3.3 结果（cal_data/ + noise.json + repeat.json → crlb.json）

J 交叉验证 ∂u/∂k3：解析 fx·x·r⁶=-464.6 vs 数值=-465.4，相对差 0.17% → J 可信。

9 个内参的 CRLB vs 实测散布（iso_rms 主视角，σ_a=rms/√2；CRLB=仅检测噪声决定的理论 std 下界）：

内参	CRLB（理论噪声底）	std_emp（8 组散布）	std/CRLB	解读
fx	0.0543	0.0767	1.41	✓ 有效
fy	0.0544	0.0725	1.33	✓
cx	0.0091	0.0524	5.76	⚠ 超 CRLB
cy	0.0072	0.0746	10.3	⚠
k1	3.0e-5	4.4e-5	1.43	✓
k2	4.0e-5	1.08e-4	2.69	接近有效（略受位姿影响）
p1	6.8e-6	3.6e-6	0.53	≈0，小数相除不稳
p2	5.5e-6	2.3e-6	0.41	≈0，小数相除不稳
k3	4.5e-5	8.8e-5	1.97	✓

有效性主项均值（fx,fy,cx,cy,k1,k3）：iso_rms=3.70，aniso_rms=3.85（p1/p2≈0 跳过）。

9×9 内参信息矩阵 I_β（Schur 边缘化外参后；对角=各内参信息量，越大越准；非对角=参数间耦合，负耦合会互相"吃掉"信息）：

	fx	fy	cx	cy	k1	k2	p1	p2	k3
fx	1.7e+03	-1.4e+03	4e+02	-4.3e+02	1.2e+06	7.7e+05	-1.4e+05	9.4e+05	4.2e+05
fy	-1.4e+03	1.5e+03	-2.7e+02	5.1e+02	-7.1e+05	-4.4e+05	1e+05	-3.1e+05	-2.5e+05
cx	4e+02	-2.7e+02	1.2e+04	-2.2e+02	3.2e+05	2.1e+05	3.4e+05	1.9e+06	1.7e+05
cy	-4.3e+02	5.1e+02	-2.2e+02	2.2e+04	-2.9e+05	-2.2e+05	7.2e+06	-1.4e+06	-1.3e+05
k1	1.2e+06	-7.1e+05	3.2e+05	-2.9e+05	3.7e+09	2.3e+09	-3e+08	2.9e+09	1.5e+09
k2	7.7e+05	-4.4e+05	2.1e+05	-2.2e+05	2.3e+09	2.2e+09	-2.2e+08	2.4e+09	1.2e+09
p1	-1.4e+05	1e+05	3.4e+05	7.2e+06	-3e+08	-2.2e+08	2.4e+10	1.1e+09	-1.1e+08
p2	9.4e+05	-3.1e+05	1.9e+06	-1.4e+06	2.9e+09	2.4e+09	1.1e+09	3.6e+10	1.7e+09
k3	4.2e+05	-2.5e+05	1.7e+05	-1.3e+05	1.5e+09	1.2e+09	-1.1e+08	1.7e+09	1.2e+09

读法：对角 fx-fx=1.7e3、cy-cy=2.2e4、p2-p2=3.6e10（主点/切向信息量很大→CRLB 很小）；fx-fy 强负耦合（-1.4e3）解释了 fx、fy 不能无限独立变准；cond(I_β)=2.0e8、λ_min=181.6 ≫ 浮点地板 8e-6 → 数值稳定无 NaN。

3.4 解读

• fx/fy/k1/k3 达理论有效（std/CRLB≈1.3–2）：这四个参数的数据已被用尽，实测散布就是噪声决定的理论下界附近。
• cx/cy 超 CRLB 5–10 倍——但不是因为主点散布大。关键看绝对值：cx 的 std_emp=0.052 和 fx 的 0.077 同量级；ratio 大（5.8、10.3）纯粹是因为 CRLB_cx/cy 极小（0.009/0.007）。 88 角点×80 图=7040 观测，主点由棋盘多姿态的对称中心强约束，信息量巨大 → 理论预测能定到 0.009px。
• 实测达不到 0.009px 的根因是 CRLB 只管随机误差，实测散布含系统误差（棋盘微动 → 整板平移 ≡ 主点反向平移，cx/cy 对此最敏感）。这在 §6.1（剔微动）与 §6.2 闭环验证。

→ 衔接第 4 步：cx/cy 超界可能来自模型阶数不够（高阶畸变没建好→残差→主点偏），先做多折定阶排除； → 衔接第 5 步：把 std/CRLB 作为九条判据之一（第 9 条：fx std/CRLB≈1 ✓）。

第 4 步｜多折模型阶数

4.1 原理：逐级加畸变，找 hold-out RMS 平台

bias-variance 权衡（§1.5）：加畸变参数降 Bias²（欠拟合残差）但升 Variance（过拟合噪声），有一个最优阶。做法：逐级加畸变（k1 → k1-2 → k1-3 → +p1p2 → +s1s2 → +k4-6），每级做 k 折 hold-out：

• 训练集（75 张）拟合内参（该级 flags）→ 固定内参对验证集（5 张）解外参（solvePnP）→ 验证集重投影 RMS。
• hold-out RMS 衡量泛化，in-sample RMS 衡量拟合；两者差<20% = 无过拟合。
• 平台判据（九条第 7）：相邻阶数 hold-out 降幅 <5% = 模型饱和，再加参数是过拟合（降反升）。

8 组固定折法（seed=42，8 组共用同一套折索引）→ 组间只有噪声实例不同，跨组 hold-out 可叠加比较。这是真实数据版的 Monte Carlo：用 8 组重复连拍代替换噪声种子。

4.2 源码（阶梯 flags + hold-out + 固定折法）

# 阶数序列(文档 §3.3):逐级加畸变参数
DIST_LADDERS = [
    ("k1",        cv2.CALIB_FIX_K2|cv2.CALIB_FIX_K3|cv2.CALIB_ZERO_TANGENT_DIST),
    ("k1-2",      cv2.CALIB_FIX_K3|cv2.CALIB_ZERO_TANGENT_DIST),
    ("k1-3",      cv2.CALIB_ZERO_TANGENT_DIST),
    ("k1-3+p1p2", 0),                       # 默认 Brown 5 系数
    ("+s1s2",     cv2.CALIB_THIN_PRISM_MODEL),
    ("+k4-6",     cv2.CALIB_RATIONAL_MODEL),
]

def holdout_fold(board, train_corners, val_corners, shape, flags):
    """训练拟合内参(flags 控制阶数)→ 固定内参对验证集解外参 → 验证集重投影 RMS。"""
    cal = calibrate_camera(board, train_corners, shape, flags=flags)
    objp = board.object_points().astype(np.float32).reshape(-1,1,3)
    sq = []
    for c in val_corners:
        pts = np.asarray(c).astype(np.float32).reshape(-1,1,2)
        ok, r, t = cv2.solvePnP(objp, pts, cal.K, cal.dist, flags=cv2.SOLVEPNP_ITERATIVE)
        proj, _ = cv2.projectPoints(objp, r, t, cal.K, cal.dist)
        sq.append(((proj.reshape(-1,2) - pts.reshape(-1,2))**2).sum(1))
    return float(np.sqrt(np.concatenate(sq).mean())), float(cal.rms)   # hold-out, in-sample

# 固定折索引(8 组共用 → 组间折法相同,只差噪声)
rng = np.random.default_rng(42)
folds = [set(rng.permutation(n_poses)[:5].tolist()) for _ in range(40)]   # 40 折,每折留 5

4.3 结果（cal_data/ → cv.json，8 组×40 折=320 点/级）

[         k1]  320 点  hold-out RMS mean=1.175  (欠拟合:只 1 阶径向不够)
[       k1-2]  320 点  hold-out RMS mean=0.258
[       k1-3]  320 点  hold-out RMS mean=0.142  ← 饱和平台
[ k1-3+p1p2]  320 点  hold-out RMS mean=0.142  (+p1p2 不降:切向≈0)
[      +s1s2]  320 点  hold-out RMS mean=0.318  (过拟合:薄棱镜多余)
[      +k4-6]   ~点    hold-out RMS mean=13.85  (崩:8 参数不收敛/数值病态)

平台判据:k1-2→k1-3 降 45%(<饱和前);k1-3→+p1p2 降 0.3% ← 饱和;之后降反升 = 过拟合
in-sample vs hold-out(k1-3):gap≈+1% ✓ 无过拟合

4.4 解读

• 最优畸变模型 = k1-3（3 阶径向，hold-out 在此饱和=0.142）。p1p2≈0（第 2 步）切向可去； s1s2/k4-6 过拟合或不收敛。模型阶数不是 RMS>σ_det 的原因——加阶不降反升。
• 这排除归因候选 （a） 模型不够：RMS=0.142 ≫ σ_det=0.046 不是因为畸变建得不够。
• → 衔接第 5 步：平台值 0.142 进地板判据 RMS_hold-out≈RMS_in-sample≫σ_det； → 衔接第 6 步：定阶后固定 5 参数（k1-3+p1p2，范用）做 10 折可靠性验证。

第 5 步｜综合报告

5.1 原理：九条判据 + 地板判据，一次性定性

把 1-4 步的数字汇总成可勾选的判据（九条），把"这组数据标定得好不好"变成客观打钩表。关键是地板判据（§1.6、§七）：RMS_hold-out ≈ RMS_in-sample ≫ σ_det 三量各检验一件事——

• RMS_hold-out ≈ RMS_in-sample：训练泛化一致 → 无过拟合（模型没记住噪声）。
• 两者都 ≫ σ_det：RMS 没压到噪声地板 → 卡硬件上限（系统误差，非随机噪声）。

判据里红绿分明：绿（✓2 无过拟合、✓6 约束充分、✓7 模型饱和、✓9 有效）说明算法/数据无问题； 红（✗1 RMS 未到 1.5σ_det、✗3 RMS/σ_det∉[0.8，1.2]）都指向同一根因：RMS 卡硬件，不是算法问题。 求值点：CRLB 的求值点用 theta_mean 并重解外参（同第 3 步），不沿用第 0 组标定外参。

5.2 源码（判据计算 + 求值点处理）

# 求值点:theta_mean 处重解外参(非第 0 组 cal 外参)
rvecs, tvecs = solve_extrinsics_for_theta(board, corners, theta_mean)
J = numerical_jacobian(theta_mean[:9], rvecs, tvecs, board.object_points())
crlb, cond, _ = crlb_from_J(J, sigma_a=cal.rms / np.sqrt(2))

rms, sd = cal.rms, noise["sigma_det"]
checks = [
    ("1. in-sample RMS ≤ 1.5σ_det",  rms <= 1.5*sd, f"{rms:.3f}≤{1.5*sd:.3f}?"),
    ("3. RMS/σ_det ∈ [0.8,1.2]",      0.8<=rms/sd<=1.2, f"{rms/sd:.2f} 到地板?"),
    ("6. fx 跨子集 cv < 0.3%",        fx_cv<0.3, f"{fx_cv:.3f}% 约束充分?"),
    ("7. 加畸变 hold-out 降 <5%",     drop<5, f"k1-3→+p1p2 降{drop:.1f}% 饱和?"),
    ("2. hold-out vs in-sample <20%", gap<20, f"gap {gap:+.1f}% 无过拟合?"),
    ("9. fx std/CRLB ≈ 1",            0.5<=ratio<=2, f"{ratio:.2f} 有效?"),
]
# 地板判据:RMS_hold-out=0.142 ≈ RMS_in-sample=0.138 ≫ σ_det=0.046 (3.1×)

5.3 结果（全部 json → report.json）

九条判据:
  ✓ 2. hold-out vs in-sample <20%:k1-3 gap +1.4%  ✓ 无过拟合
  ✓ 6. fx 跨子集 cv <0.3%:0.0039%  ✓ 约束充分
  ✓ 7. 加畸变 hold-out 降<5%:k1-3→+p1p2 降 0.3%  ✓ 饱和
  ✓ 9. fx std/CRLB ≈1:1.41  ✓ 有效
  ✗ 1. in-sample RMS ≤1.5σ_det:0.138≤0.069? ✗ 卡硬件
  ✗ 3. RMS/σ_det ∈[0.8,1.2]:3.0 ✗ 未到地板(卡硬件上限)
地板判据:RMS_hold-out 0.142 ≈ RMS_in-sample 0.138 ≫ σ_det 0.046 (3.1×)
  → hold-out≈in-sample ✓ 无过拟合;两者≫σ_det → 卡硬件上限(镜头像差+靶标)

5.4 解读

绿判据全过（无过拟合/约束充分/饱和/有效）→ 算法与模型没问题；红判据（✗1✗3）只反映"RMS 没压到噪声地板"，根因是硬件系统误差。三条 RMS（0.142/0.138/0.046）的层级关系是整个评估的核心结论：标定达到泛化一致（泛化一致），但精度被硬件封顶。→ 衔接第 6-9 步：从更多角度（可靠性/数据量/ 独立测试/数据一致性）反复验证这个结论。

第 6 步｜5 参数 10 折验证（可靠性）

6.1 原理：固定模型，看标定数值稳不稳

第 4 步定了阶（k1-3 饱和），但实际范用选 5 参数（k1-3+p1p2）——切向 p1p2 虽≈0，但不崩溃、对一般镜头更通用。本步固定 5 参数，8 组各 10 折（固定折法 seed=42，组间只差噪声），每折训练 72 张标定 → θ̂_fold，跨 80 折 的 std = 标定的经验可靠性。

这是比 CRLB 更"粗"但更"真"的不确定度：CRLB 只含随机噪声（预测 0.054），而跨折 std 还包含 “用哪些位姿/哪批图"带来的散布（子集变化 + 组间微动叠加）。两者之比 std/CRLB 直接反映"实测比理论差几倍”。

6.2 源码（固定 5 参数 + 80 折）

FLAGS_5PARAM = 0   # k1-3 + p1p2(默认 Brown 5 参数)

def holdout_fold(board, train_corners, val_corners, shape):
    cal = calibrate_camera(board, train_corners, shape, flags=FLAGS_5PARAM)  # 训练→θ̂
    objp = board.object_points().astype(np.float32).reshape(-1,1,3)
    sq = []
    for c in val_corners:
        pts = np.asarray(c).astype(np.float32).reshape(-1,1,2)
        ok, r, t = cv2.solvePnP(objp, pts, cal.K, cal.dist)          # 固定内参解外参
        proj, _ = cv2.projectPoints(objp, r, t, cal.K, cal.dist)
        sq.append(((proj.reshape(-1,2) - pts.reshape(-1,2))**2).sum(1))
    return cal.theta, float(np.sqrt(np.concatenate(sq).mean())), float(cal.rms)

rng = np.random.default_rng(42)
folds = [set(rng.permutation(n_poses)[:8].tolist()) for _ in range(10)]   # 10 折,每折留 8
thetas = []
for rep_i in range(8):                      # 8 组 × 10 折 = 80 折
    for val_set in folds:
        theta, ho, ins = holdout_fold(board, train_72, val_8, shape)
        thetas.append(theta)
theta_std_80fold = np.array(thetas).std(0, ddof=1)   # 跨 80 折散布 = 经验可靠性

6.3 结果（cal_data/ + crlb.json → cv5.json）

θ̂ 跨 80 折（8 组 × 10 折，固定 5 参数）的 mean / std / cv，与理论 CRLB 对照：

内参	mean	std（跨 80 折）	cv（%）	CRLB（理论噪声底）	std/CRLB
fx	1854.9	±0.279	0.0151	0.0543	5.15
fy	1845.3	±0.273	0.0148	0.0544	5.02
cx	1269.6	±0.333	0.0262	0.0091	36.5
cy	725.12	±0.174	0.0240	0.0072	24.1
k1	−0.40502	±2.1e-4	0.0510	3.0e-5	6.78
k2	0.2299	±5.4e-4	0.2345	4.0e-5	13.5
p1	6.9e-5	±1.5e-5	22.0	6.8e-6	2.24
p2	−2.3e-5	±9.4e-6	40.3	5.5e-6	1.69
k3	−0.08090	±4.7e-4	0.583	4.5e-5	10.6

hold-out RMS（80 折验证集）：mean=0.1445，median=0.1411，p95=0.168，max=0.169，min=0.127（相对 in-sample 0.138 gap +5% → 无过拟合）。

结果分析：

• θ̂ 跨 80 折极稳：fx=1854.9±0.28（cv 0.015%），换不同位姿子集 / 噪声实例，内参几乎不变 → 标定可靠。
• 但跨折 std ≫ CRLB：fx std/CRLB≈5、cx≈37、cy≈24。实测散布是"纯噪声理论"的几倍到几十倍——差额来自"用哪些位姿"的子集敏感性 + 组间微动，不是纯随机噪声。
• p1/p2 的 std/CRLB<1（2.2、1.7）是因为真值≈0、信息少、小数相除不稳（非异常）；k2 的 cv 偏大（0.23%）是高阶径向弱约束。
• 这解释了第 3 步 cx/cy 超 CRLB：理论只算噪声，实测还吃系统 / 位姿贡献。

6.4 解读

• → 衔接第 7 步：跨折散布 > 理论 CRLB，加更多数据能否压下去？第 7 步看 k 帧 vs 精度。
• → 衔接第 8 步：80 折 std 是最终不确定度的三路之一（取保守上界）。

第 7 步｜精度 vs 数据量（硬件 vs 数据）

7.1 原理：加图降的是理论精度，降不动实际 RMS

每位姿取前 k 帧（k=1，2，3，4，5 → 80/160/240/320/400 张），看加图能否提升精度。两个量要分开看：

• CRLB ∝ 1/√k：随机误差随数据量下降（理论精度变好）。
• RMS（实测残差）卡硬件：系统误差（镜头像差+靶标）与数据量无关，基本不动。

所以预期是"加图 → CRLB 降、RMS 不降"。这把"硬件上限"和"数据不足"区分：若 RMS 随 k 降，说明之前是数据不足；若不降，说明已是硬件封顶。工程约束：640 全量无初值会卡（外参 6×640=3844， LM 解 3853×3853 系统超慢），所以给 theta_mean 初值（USE_INTRINSIC_GUESS），并用 k=5（400 张）作为"最大可用数据"。

7.2 源码（给初值的 k 帧标定 + Fisher CRLB）

def calibrate_with_guess(board, corners, shape, theta0):
    """给初值的 5 参数标定(USE_INTRINSIC_GUESS,加速大图数收敛)。"""
    K0 = np.array([[theta0[0],0,theta0[2]],[0,theta0[1],theta0[3]],[0,0,1.]])
    dist0 = np.array(theta0[4:9])
    objp = board.object_points().astype(np.float32).reshape(-1,1,3)
    rms, K, dist, rvecs, tvecs = cv2.calibrateCamera(
        [objp]*len(corners), [c.astype(np.float32).reshape(-1,1,2) for c in corners],
        (shape[1], shape[0]), K0, dist0, flags=cv2.CALIB_USE_INTRINSIC_GUESS)
    return CalibrationResult(K, dist[:5], rvecs, tvecs, rms, len(corners), shape)

for k in [1, 2, 3, 4, 5]:
    corners_k = [前 k 帧的所有 80 位姿]                 # 80×k 张
    cal = calibrate_with_guess(board, corners_k, shape, theta_mean)
    J = numerical_jacobian(cal.theta[:9], cal.rvecs, cal.tvecs, board.object_points())
    crlb, _, _ = crlb_from_J(J, sigma_a=cal.rms/np.sqrt(2))   # 该 k 下的理论精度

7.3 结果（cal_data/ + repeat.json → sub.json）

k	图数	rms	CRLB_fx
1	80	0.1379	0.0545
2	160	0.1379	0.0385
3	240	0.1380	0.0314
4	320	0.1380	0.0272
5	400	0.1380（+0.1%）	0.0241（−55.8%）

解读:
  - rms:k=1→5  0.1379→0.1380 (+0.1%)。基本不动 = 系统误差(硬件)与数据量无关,加图降不动。
  - CRLB_fx:k=1→5  0.0545→0.0241 (−55.8%)。理论精度随数据 ∝1/√k 改善(随机误差降)。
  - θ̂ 随 k 趋稳(数据越多估计越稳)。

补充：文件夹级交叉验证（训练每位姿全 8 帧 = 432 张）：

• 3 折 hold-out：[0.1415， 0.1399， 0.1372]，mean=0.1396（std 0.0018，极稳）
• 对比图片级（每位姿 1 帧，75 张）0.1424 → 数据 ×6，hold-out 只降 2%
• 再次印证：rms 卡硬件，加数据（×6）降不动。

7.4 解读

结论（硬件 vs 数据的关键判决）：加 5 倍数据，RMS 几乎不动（0.1379→0.1380），只有 CRLB 降。说明 RMS=0.138 不是"图不够多"，而是镜头剩余像差 + 靶标不平的硬件天花板——排除归因候选（b）运动（微动）、印证候选（c）硬件。要更准须换硬件（玻璃镀铬靶 + 更平镜头），不是加图。

实践含义：没必要在同一位置重复拍摄标定图。同一位姿的连拍帧高度相似（只差噪声实例），对内参估计是重复信息——k=1→5 帧（80→400 张）RMS 不动、CRLB 虽降但远大于实际散布（std/CRLB≈5），多拍只会无谓增加标定求解的参数规模（外参 6×N 增大、LM 变慢），不改善实际精度。每位姿拍 1 张、把精力放在位姿多样性（位置九宫格 + 俯仰/偏航/滚转 + 远近）上才是有效信息——多样位姿才真正增加 Fisher 信息、改善参数解耦（否则 fx/cy 等耦合退化、CRLB 爆 NaN）。

重复结构（8 帧/位姿）的唯一价值是重复性分析（8 组重组 = 真实版 Monte Carlo），不是"更多标定数据"。 → 衔接第 8 步：每位姿取 1 张（80 张扁平）给出 θ̂_final。

第 8 步｜最终内参 + 误差光流

8.1 原理：三路整合给最终内参，光流看系统误差在哪

最终内参 = calibrate_full 对扁平 80 张（每位姿 1 张，无重复帧）的标定结果 θ̂_full，不确定度取 三路保守上界 max(CRLB, std_repeat, std_cv5)——任一条路给出的不确定度都不会比这更乐观，稳妥。

光流图回答"RMS=0.138 这个系统误差在图像什么位置、朝什么方向"：固定最终内参，逐图解外参→投影→ 每角点残差（proj−detected）按像素位置画箭头。看三件事：边缘是否更大（镜头边缘像差）、某方向是否系统偏（各向异性/装配偏心）、残差是随机还是有结构（剩噪声 vs 模型/硬件系统误差）。

8.2 源码（三路整合 + 逐角点残差光流）

# 整合三路
theta_final = np.array(full["theta_full"])                  # 扁平 80 张标定(calibrate_full 产出)
crlb_full   = np.array(full["crlb"])                        # CRLB(理论噪声底)
std_repeat  = np.array(repeat["theta_std"])                 # 8 组散布(经验)
std_cv5     = np.array(cv5["thetas"]).std(0, ddof=1)        # 80 折散布(经验)
for k, name in enumerate(PARAM_NAMES):
    unc = max(crlb_full[k], std_repeat[k], std_cv5[k])      # 三路取保守上界

# 逐角点残差 → 按像素位置画箭头
for i, c in enumerate(all_corners):
    proj, _ = cv2.projectPoints(objp, cal.rvecs[i], cal.tvecs[i], K, dist)
    det = np.asarray(c).reshape(-1,2); prj = proj.reshape(-1,2)
    dx = prj[:,0]-det[:,0]; dy = prj[:,1]-det[:,1]          # 残差(预测-检测)
ax.quiver(det_x, det_y, dx*scale, dy*scale, ...)            # 按像素位置画箭头

8.3 结果（最终内参表 + 训练光流）

最终内参 θ̂_final（扁平 80 张标定，每位姿挑 1 张），不确定度取三路保守上界：

内参	θ̂_final	± 不确定度	（CRLB	8组std	10折std）
fx	1855.058	±0.279	0.054	0.077	0.28）
fy	1845.492	±0.275	0.054	0.073	0.27）
cx	1269.49	±0.300	0.009	0.052	0.30）
cy	725.02	±0.300	0.007	0.070	0.30）
k1	-0.40515	±4e-5	3e-5	4e-5	2e-5）
k2	0.23024	—		p1≈0/p2≈0（切向可忽略）
k3	-0.08118	±9e-5	5e-5	9e-5	—）

θ̂_final 现由 calibrate_full 对扁平 80 张（每位姿 1 张，无重复帧）标定得出——同一位姿的连拍帧是重复信息，加图只降 CRLB 不降 RMS（第 7 步），80 张足够。若用 400 张（每位姿 5 张， run_pipeline --final-frames 5）得 fx=1855.098，二者在 ±0.4px 重复性内一致。

• fx 相对精度 0.0041%（不确定度 0.279/fx），远优于工程判据 0.1-0.3%。
• 训练光流：RMS=0.169，边缘/中心≈1.06（边缘略大，镜头边缘像差），v>u（各向异性，印证 σ_v>σ_u）。

逐角点光流箭头多而杂（系统趋势 + 检测噪声混在一起）。标准测试流程里另算一份分块平均残差场 （按像素位置分 8×10 块，块内取平均，滤掉随机噪声显系统趋势）+ 它的向量相关性，作为第 9 步测试评估的固定输出（见 9.4），不在此单独展开。

→ 衔接第 9 步：用 cal_test 独立测试泛化，并给出分块残差场的向量相关性（系统误差稳定性）。

第 9 步｜独立测试（cal_test 泛化）

9.1 原理：金标准——完全不参与训练的新数据

交叉验证（第 4、6 步）的验证集终究来自同一批采集；cal_test（64 张 SDK 直采，完全不参与训练，不同时段）是更可靠的金标准。把各组训练内参（repeat 8 组 / cv5 80 折 / full）拿到 cal_test 上： 固定内参解每图外参（solvePnP）→ 重投影 RMS = 泛化精度。多角度指标不只看 RMS： bias（系统偏置，应≈0）、std_u/v（各向异性）、E/C（边缘/中心）、per-image 分布、max 残差。各组吃同样的测试噪声，故相对排名可信（质量差只抬高绝对值，不改变结论）。

9.2 源码（固定内参解外参 + 多角度指标）

def test_intrinsics(board, test_corners, theta):
    """用内参 theta 在 test 角点上解 PnP(固定内参保外参)→ 每图 + 总重投影 RMS。"""
    K = np.array([[theta[0],0,theta[2]],[0,theta[1],theta[3]],[0,0,1.]])
    dist = theta[4:9]
    objp = board.object_points().astype(np.float32).reshape(-1,1,3)
    per_img, all_sq = [], []
    for c in test_corners:
        pts = np.asarray(c).astype(np.float32).reshape(-1,1,2)
        ok, r, t = cv2.solvePnP(objp, pts, K, dist, flags=cv2.SOLVEPNP_ITERATIVE)
        proj, _ = cv2.projectPoints(objp, r, t, K, dist)
        d2 = ((proj.reshape(-1,2) - pts.reshape(-1,2))**2).sum(1)
        per_img.append(float(np.sqrt(d2.mean()))); all_sq.append(d2)
    return float(np.sqrt(np.concatenate(all_sq).mean())), per_img

# evaluate_one 另算:bias/std_u/std_v/E/C/per-img p95/max(多角度)

9.3 结果（cal_test 64 张，多角度指标）

指标	值	解读
测试 RMS	0.27（各组一致到第 4 位）	内参可靠
bias_u/v	≈0 / −0.0001	无系统偏置 ✓
std_u/v	0.42 / 0.44	v>u（各向异性，印证 σ_v>σ_u）
RMS/σ_det	5.9	≫1（泛化 gap）
max 残差	0.72 px	个别差图
各组测试 RMS	0.27±<0.01	repeat/cv5/full 几乎相同

== 相对比较(各组吃同样测试噪声,排名可信)==
   泛化最好:FULL(400 张训练)= 0.272
   各组(repeat#0-7 / cv5 / full)测试 RMS 全部落在 0.27±0.01 → 内参一致、可靠

9.4 分块残差场 + 向量相关性（标准测试输出）

逐角点光流太杂（系统趋势 + 检测噪声混在一起）。标准测试流程固定再算一份分块平均残差场：残差按像素位置分 8×10 块、块内取平均 → 零均值随机噪声被 √N 压低，剩下系统趋势。再把每个块的平均 （du，dv） 摊成一个向量，对 8 个训练重复组 + 测试集（都用 θ̂_final）各算一个向量，看它们的相关。

结果一：分块残差场（去噪声后的系统趋势）

数据	块平均径向 mean	max	全局 RMS
训练 8 组（各自）	≈0.030 px	0.15–0.20	0.138
测试 cal_test	0.120 px	0.941（左下边缘）	0.272

块平均（0.03/0.12）远小于全局 RMS（0.14/0.27）：大部分残差是随机噪声，被平均压掉。测试集左下角 max≈0.94 是 k3 表达不了的边缘像差——硬件上限的空间位置。

结果二：残差向量相关性（核心）

组间相关(8 组,同位姿不同帧):  mean=0.945   min=0.910   max=0.984   ← 高度可复现
组-vs-测试(8 组各对测试):      mean=0.012   (逐组 -0.02~0.04)        ← 不相关
训练平均场 vs 测试场:          0.012

怎么读（两层结论，都重要）：

• 组间相关 0.94 ≫ 0：同位姿跨噪声实例（8 组只差取了哪一帧）残差场几乎一样 → 残差是 系统性/确定性的（由位姿决定），不是随机噪声。这直接证明 RMS>σ_det 是系统误差，而非噪声涨落（若纯噪声，组间相关应≈0，RMS 也就≈σ_det）。
• 组-vs-测试 ≈ 0：换成测试集的位姿后残差场与训练场不相关 → 系统残差场是位姿相关的（镜头畸变残余按不同位姿投影出不同图案），不是一张固定的像素像差图。说明"硬件上限"是按位姿表现的镜头/模型残余，加数据（第 7 步）降不动它的根源在此——它不是某张静态像差，而是模型（k1-3）对任意位姿都留一点吃不掉的残余。

train-vs-test ≈ 0 的含义（详解）

字面：Pearson 相关衡量两个向量"去均值、归一化后图案是否匹配"。≈0 = 训练残差场的空间图案与测试残差场在统计上无关（图案对不上）。两层原因：

• （a） 训练残差是"拟合后残差"（post-fit）：θ̂_final 本就是在训练数据上拟合来的，优化器已把训练残差压到极小（均值 0.03px）；测试残差是"泛化残差"（0.12px），模型没见过这些位姿。一个是被抑制的样本内残差，一个是样本外残差，结构本就不同。
• （b） 残差按"图像像素 （u，v）"分块，但物理来源是"场角/角点"相关：镜头径向畸变残余 + 靶标方格不均匀都随场角变，不随像素位置变。同一像素块在不同位姿下被不同场角的角点命中 → 块内均值依赖位姿。训练位姿 ≠ 测试位姿 → 分块图案不同 → 相关 ≈ 0。对照：8 个训练组位姿相同 → 命中混合相同 → 相关 0.94。

它说明什么（可下结论）：

• 系统残差场不是一张"固定的像素像差图"。若镜头存在简单的"每像素固定偏移"，则 train/test 按像素分块的残差应是同一张图 → 相关高；≈0 否定了这一点。
• 系统残差是位姿相关的：随标定板怎么摆（俯仰/偏航/远近/位置）而变。

它不说明什么（避免误读）：

• ≠ “镜头像差是随机的”。像差是固定的物理量，只是不表现为"固定像素偏移图"（它场角相关，投影到像素后随位姿变）。
• ≠ “测试用了不同相机/内参”。同一台相机、同一个 θ̂_final。
• 不否定"硬件上限"。根因仍是镜头剩余像差 + 靶标不均匀，只是按位姿调制表现。

方法论注意：

• 训练侧是 post-fit 残差（被压小），故 train-vs-test 是不对称比较；但组间 0.94 证明本方法能检出真实的图案相关（图案真匹配时会高），所以 ≈0 是"图案真不匹配"，不是低信号伪影。
• 若要更干净地检验"有没有固定像差"，应比较两个位姿覆盖相似的独立数据集，或直接看畸变模型自身预测的残差场——但定性结论（位姿相关、非固定像素图）不变。

那镜头的固定像差在哪——模型已经把它捕获了

train-vs-test≈0 说明残差不是固定像差图。但镜头确实有固定像差——它被畸变模型（k1-3+p1p2）捕获了。下面把"模型捕获的畸变"这张图算出来。

δ（u，v）是什么、怎么算：对图像每个像素位置 （u，v），先用内参反算归一化坐标 $x = (u−cx)/fx, y = (v−cy)/fy$（假设没畸变时该像素对应的方向），再代入畸变模型算出畸变后的 $（x_d, y_d）$（径向 $x_d ← x·(1+k1·r²+k2·r⁴+k3·r⁶)$ + 切向，$r²=x²+y²$），最后换回像素位移 $δ = ( fx·(x_d−x), fy·(y_d−y) )$，模长 $|δ| = √(δu²+δv²)$。

275 px 是什么意思：δ 就是"畸变把角点挪了多远"——理想针孔（无畸变）下角点该落在像素 P₀，实际有畸变落在 P₁，$δ = P₁ − P₀$。它只依赖像素位置（给定 θ̂ 即固定，与棋盘怎么摆无关），所以是一张固定的像差图。中心 $r≈0$ 处 $|δ|≈0$；越往边缘 $r$ 越大，桶形畸变（k1<0）把角点往中心拉，到图像边角 $|δ|$ 最大 ≈ 275 px——即边缘的角点被畸变挪动了约 275 个像素。把它和实测残差场并排：

场	边缘 max	性质
模型像差场 δ（u，v）	~275 px	固定像素图（纯 （u，v） 函数），模型已捕获
实测残差场（test）	0.94 px	位姿相关，模型吃不掉的零头

残差 max / 像差 max = 0.94/275 = 0.34% → 模型捕获了 99.9%+ 的镜头畸变。

这把 train-vs-test≈0 的疑问说清：固定像差是存在的（275px，模型已建成 δ（u，v） 这张固定图）；残差（0.94px）之所以位姿相关、不构成固定图，是因为它是"δ 之外的高阶残余 + 靶标方格不均匀 + 外参传播"的混合，按像素分块时被位姿打乱。所以"硬件上限"的真相是：镜头主畸变（275px 量级）已被模型吃掉，剩下的 <1px 残差是靶标/高阶/位姿残余——这部分加图降不动（第 7 步），换镜头+玻璃靶才能降。

9.5 解读

• 各组内参在独立测试集上泛化一致（0.27±0.01）→ 内参不是过拟合训练集，是真实物理量。
• 0.27 > 训练 0.14 是 cal_test 的泛化 gap：测试集位姿覆盖/难度不同（且含测试集自身噪声）， 不是内参问题。bias≈0 说明无系统偏置。
• → 衔接第 10 步：泛化一致还不够，最后用"数据一致性"正面回答"内参是不是数据无关的物理常量"。

9.6 测试误差 vs 标定板大小：近拍脱焦（非畸变）

测试 RMS 0.27 比训练 0.14 高，除位姿覆盖不同外，还有一个结构性来源：标定板在画面里越大（离相机越近）的图，误差越大。把 64 张测试图按各自的板视大小（角点到质心均距）与每图 RMS 算相关：

相关	系数	含义
RMS vs 板视大小	+0.68	强正相关：板越大/越近，误差越大
RMS vs 均归一半径（边缘程度）	−0.14	几乎无关 → 不是边缘畸变
板视大小 vs 边缘半径	−0.13	大板并不更靠边缘

按板大小分组：大板（>293 px）平均 RMS 0.317，小板 0.173，约 2 倍。

判定为近距离脱焦，而非畸变或模型问题：

• 若是边缘畸变，误差应随角点离图像中心的归一半径增大（畸变 $r²$ 项）；但这里与边缘半径几乎无关（−0.14），与板视大小强相关（+0.68）。
• 板大意味着每个特征占的像素更多，亚像素检测本应更准（误差更小），实测反而更大 → 排除像素尺度因素。
• 剩下最合理的解释：板离相机越近（画面里越大），越偏离镜头焦平面 → 角点越糊 → $cornerSubPix$ 检测噪声越大 → 重投影 RMS 越大。即测试集里的近拍大板脱焦，抬高了测试 RMS。

这不影响内参可靠性（各组内参在这些图上 RMS 仍一致），但解释了泛化 gap 的一部分来源；要降低需在采集端覆盖相近的距离段、或收光圈增大景深。

9.7 每步内参 → cal_test 泛化 RMS 汇总

把每个产出内参的步骤（②⑥⑦⑧）的内参，都拿到同一份 cal_test 上评估泛化 RMS，横向对比：

步骤	内参来源	训练图数	cal_test 泛化 RMS
② repeat	8 组重复标定（各 80）	80	0.2717 ± 0.0000（8 组几乎相同）
⑥ cv5	5 参数 10 折（各 72）	72	0.2718 ± 0.0002 [0.272， 0.272]
⑦ subsample	k=1…5 帧	80–400	≈0.2717（与数据量无关，见第 7 步）
⑧ final	θ̂_final（扁平 80 张）	80	0.2717

结论：不管用多少数据（80→400）、哪种方法（单组/10 折/全量），内参在独立测试集上的泛化 RMS 都≈0.2717，差异在小数点后第 4 位。这把散落在 ②⑥⑦⑧ 的内参用"泛化 RMS"这把统一尺子量了一遍，结果高度一致——印证第 7 步"RMS 卡硬件、与数据量/方法无关"，也支撑 §九"内参是数据无关的物理常量"。

第 10 步｜数据一致性验证

10.1 原理：换一批数据，内参变不变

第 9 步证明内参泛化可靠，但"可靠"是"各组估的差不多"。本步正面回答更强的命题：内参是不是相机的物理常量、与具体训练数据无关？ 用三套不同来源的数据各自独立标定，比 θ̂：

• train：cal_data 第 0 组（80 个位姿各 1 帧）
• test：cal_test 全部（64 张，SDK 直采、不同时段）
• mix：train + test 拼一起（144 张）

若三套 Δθ 极小（如 Δfx<1px，即 fx 的 0.1‰ 量级），说明换数据内参几乎不变 → 物理常量。 初值处理：三套都用同一 theta_mean 作初值（USE_INTRINSIC_GUESS）——既保证小样本（64 张）不陷局部极小，又让三路差异只来自数据本身、不含优化器噪声（无初值时 64 张 test 会发散：fx=1853.67 rms=0.64；有初值 fx=1854.27 rms=0.26）。

10.2 源码（三路同初值标定 + Δθ）

def calib_theta(board, corners, shape, theta0):
    """给初值标定(USE_INTRINSIC_GUESS):三路同初值 → 差异只来自数据。"""
    K0 = np.array([[theta0[0],0,theta0[2]],[0,theta0[1],theta0[3]],[0,0,1.]])
    dist0 = np.array(theta0[4:9])
    rms, K, dist, _, _ = cv2.calibrateCamera(
        [objp]*len(corners), [c.astype(np.float32).reshape(-1,1,2) for c in corners],
        (W,H), K0, dist0, flags=cv2.CALIB_USE_INTRINSIC_GUESS)
    return theta, rms, len(corners)

theta0 = repeat["theta_mean"]                  # 三路共同初值
th_tr, _, _ = calib_theta(board, c_train80,  shape, theta0)
th_te, _, _ = calib_theta(board, c_test64,   shape, theta0)
th_mx, _, _ = calib_theta(board, c_train80+c_test64, shape, theta0)
T = np.array([th_tr, th_te, th_mx])            # (3,9)
delta = T.max(0) - T.min(0)                    # 三路 Δ(max-min)

10.3 结果（--train-dir cal_data/ --test-dir cal_test/ → consistency.json）

全部 9 个内参的三路标定结果（train/test/mix 各自独立标定，三路同 theta_mean 初值）：

内参	train （80张）	test （64张）	mix （144张）	Δ（max−min）	Δ/ref
fx	1855.0583	1854.2717	1854.3535	0.787 px	0.0424%
fy	1845.4919	1844.4019	1844.6214	1.090 px	0.0591%
cx	1269.4919	1267.6240	1268.5471	1.868 px	0.1472%
cy	725.0157	724.1036	723.6907	1.325 px	0.1829%
k1	−0.405145	−0.412421	−0.405636	0.0073	1.78%
k2	0.230242	0.26993	0.233813	0.0397	16.22%
p1	7.36e−05	5.78e−04	3.20e−04	5.05e−04	155.8%
p2	−2.14e−05	−1.15e−04	5.08e−07	1.16e−04	255.1%
k3	−0.08118	−0.12970	−0.08507	0.0485	49.18%

	train	test	mix
图数	80	64	144
RMS	0.1379	0.2589	0.2061

• 主内参 fx/fy/cx/cy 极稳：绝对 Δ 只有 0.8–1.9 px（=fx/fy 的万分之几），这才是"内参稳不稳" 的硬指标。Δ/ref <0.2%。
• k2/k3/p1/p2 的相对 Δ 看起来很大（16%–255%），但这不是内参不稳，两个原因叠加：
  • 分母小：p1/p2 本身 ≈1e−4 量级（切向畸变几乎为 0），绝对 Δ 才 1e−4 px 量级——这是浮点噪声/ 小数相除，第 4 步多折已证 p1p2≈0 可去掉。
  • 高阶项弱约束：k2/k3 是高阶径向，64 张 test 单独标定时约束弱、估得偏（test 的 k3=−0.13 偏离 train/mix 的 −0.085）；但 mix（144 张）把它们拉回 train 附近（k3 mix=−0.085 ≈ train=−0.081）—— 数据一多就收敛到同一个值，恰恰证明是同一个物理量。
• 结论不变：核心几何参数（fx/fy/cx/cy）跨三套独立数据 Δ<2px（<0.2%），内参是数据无关的物理常量；高阶畸变项的表观大散布是小样本弱约束 + 小分母的假象，mix 已收敛。

10.4 解读

Δfx=0.787px（=fx 的 0.042%）：三套独立来源的数据（80/64/144 张，含完全不同的 cal_test）估出的 fx 仅差万分之零点几 → 内参是相机的物理常量、与具体训练数据无关。这是对"最终结果可靠、可泛化"最强有力的正面证据——不是"估得稳"，而是"换数据也估出同一个值"。

• RMS 差异（train 0.14 < test 0.26）来自位姿覆盖/难度不同，不是内参不稳。
• 残余 Δθ 与各小样本自身的标定噪声（CRLB，64 张 fx~0.3-0.5px）同量级，属正常统计涨落。
• capstone 结论：整条评估链（噪声→标定→CRLB→定阶→可靠性→数据量→独立测试→数据一致性）收敛到同一个 θ̂_final，它就是这台相机的内参物理常量，fx 相对精度 0.0041%，精度被硬件封顶在 RMS=0.14。

六、结果分析与归因

前 5 节（概念 → 目的 → 数据 → 测试工具 → 10 步实操）已把流程跑完、把数字摆齐，各步的「解读」里也零散点过归因。本节把这些零散判断收拢成一条闭环：集中回答最关键的一个问题——RMS（0.14）为什么是 σ_det（0.046）的 3 倍？ 把三个候选（模型 / 微动 / 硬件）逐一证伪排除，定位到唯一原因，并与仿真对照，框定系统误差的下限。

6.1 RMS>σ_det（3.1×）归因 —— 三候选逐一排除

• （a） 模型不够？ 第 4 步证伪：k1-3 饱和，加阶不降反升。
• （b） 棋盘微动？ fisher_eval --max-sigma-det 0.08 剔微动位姿：

	全 80 位姿	静止位姿
rms	0.138	0.141（未降）
fx std/CRLB	1.41	→ 有效
cx std/CRLB	5.76	仍超

rms 没降 → 微动对 rms 贡献小；但 fx/fy 剔微动后达有效，cx/cy 仍超。

• （c） 镜头像差 + 靶标？ 剩下即真相。
  • 镜头剩余像差：真实镜头像差（球差/彗差/像散/场曲/畸变/色差）是高阶非对称非多项式的，Brown 模型（k1-3，p1p2）只拟合多项式可表达部分，剩余像差=拟合后剩下的，装不进任何几个系数。
  • 靶标误差：标定假设「完美等间距方格」，真实靶标边长有制造误差、板面翘曲 → objectPoints 带系统误差。纸/亚克力靶（面不平）vs 玻璃镀铬靶（面平如镜、格准）。
  • 硬件上限：这俩是物理硬件，算法改不动。rms=0.14 是这台镜头+靶的天花板，落工业锚点（0.1-0.25）→ 正常。

6.2 cx/cy 超 CRLB 的原因

CRLB 只管随机误差；cx/cy 额外散布来自系统误差（主点对靶标对称性/检测偏置最敏感——整板平移≡主点反向平移）。剔微动后 fx/fy 达有效，cx/cy 仍超 → 主点对其他系统误差敏感。std/CRLB>1 是预期，非算法失败。

6.3 对照仿真（gap=系统误差，这是仿真的价值）

	仿真（step6）	真实（本文）
std/CRLB 主项	≈1.0（模型/靶标/数值都清零，只剩噪声）	3.70

gap=被仿真清零的误差源（硬件像差+靶标）。仿真定义算法金标准，真实偏离它的部分=真实系统误差地板。

6.4 各向同性 vs 各向异性（σ_v≈2σ_u）

σ_v/σ_u≈1.9 → 各向异性。Fisher 应 $Σ=diag(σ_u²,σ_v²)$（非 σ_a²I）。本例 aniso CRLB 略大（v 方向噪声大→信息少），更，但不改结论。 为什么 σ_v>σ_u：可能棋盘水平放置时垂直方向重力微振、传感器读出方向性、镜头某方向像质差。

七、地板判据（详解）

§六 的归因反复用到「地板判据」这条总判据（RMS_hold-out ≈ RMS_in-sample ≫ σ_det），这里把它单独展开说透：是什么、每个等式各检验什么、本例落在哪里。

$$ \mathrm{RMS}_{\text{hold-out}}\approx\mathrm{RMS}_{\text{in-sample}}\approx\sigma_{\text{det}} $$

三个量近似相等，才算标定达到了噪声地板。

用法

等式	检验	不成立说明
$RMS_hold-out ≈ RMS_in-sample$	过拟合	hold-out≫in-sample → 参数太多，把噪声当信号学
$RMS_in-sample ≈ σ_det$	到地板	≫σ_det → 有系统误差（模型/靶标/像差）
$RMS_hold-out ≈ σ_det$	泛化到地板	同上，但在未见数据上

理论原理

RMS² = e_model² + σ_det² + e_target²（四误差源，§一）。若模型完美（e_model=0）+ 靶标完美（e_target=0），只剩检测噪声 → RMS=σ_det。hold-out≈in-sample 因训练/验证同分布（MLE 一致）。

本例结论

RMS_hold-out=0.142  ≈  RMS_in-sample=0.138   ✓ 无过拟合
                 ≫  σ_det=0.046              (3.1×) ✗ 未到地板

hold-out≈in-sample（green）：k1-3 模型没过拟合，参数阶选对了。 两者≫σ_det（red）：有系统误差。结合第 4 步（模型饱和）+ §6.1（剔微动 rms 没降），定位系统误差=镜头剩余像差+靶标（硬件上限），非算法/模型。 结论：标定（没过拟合、参数阶对、算法充分利用了数据），但没到地板（卡在硬件）。这是「算法做到了它能做的最好，剩下是硬件的事」。

八、九条判据（可勾选的达标表）

地板判据是一条「总判据」；工业上还有一组可逐条勾选的达标判据（九条）。下面把本例结果填进表里，红绿分明——绿的代表算法/数据没问题，红的都指向同一个硬件根因。

#	判据	结果	达标
1	in-sample RMS ≤ 1.5σ_det	0.138 ≤ 0.069？	✗（系统误差）
2	hold-out vs in-sample <20%	k1-3 gap +3.3%	✓
3	RMS/σ_det ∈ [0.8，1.2]	3.0	✗（卡硬件）
4	残差 Moran’s I 不显著	（需残差图）	—
5	χ²/DOF ∈ [0.8，1.5]	≫1（用静态 σ_a）	✗（系统误差）
6	fx 跨子集 cv <0.3%	0.004%	✓
7	加畸变 hold-out 降 <5%	k1-3→+p1p2： -0.1%	✓（饱和）
8	独立度量 Δ ∈ ±δ	（无 CMM）	—
9	std/CRLB ≈1	fx 1.4 / cx 5.8	部分（fx✓ cx✗）

2/6/7/9（fx）绿 = 无过拟合/约束充分/模型饱和/fx 有效。1/3/5 红，根因不是算法或模型（多折已证），是硬件像差+靶标——这是硬件问题，不是流程问题。

判据表给的是定性结论（算法没问题、卡硬件）；下一节给定量落点——最终内参的数值与可保证精度。

九、最终内参与可保证精度（综合全部证据，核心结论）

这一节是全文落点：综合前 10 步 + 畸变场分析，明确给出最终内参数值与能保证到什么精度，每个数字都有出处。关键是把"精度"拆成两件不同的事——重复性（可保证）与准确性（有边界）。

9.1 最终内参 θ̂_final

内参	θ̂_final	含义
fx	1855.058 px	水平焦距
fy	1845.492 px	垂直焦距（fx≠fy → 像素微非正方/镜头微非对称）
cx	1269.49 px	主点 u（图宽 2560，近中心 1280）
cy	725.02 px	主点 v（图高 1440，近中心 720）
k1	−0.40515	径向 1 阶（主畸变，桶形）
k2	0.23024	径向 2 阶
k3	−0.08118	径向 3 阶
p1	7.4e−5	切向（≈0，镜头对中良好）
p2	−2.1e−5	切向（≈0）

为什么是这个值：θ̂_final = full.json 的 theta_full，由 calibrate_full 对扁平 80 张（每位姿挑 1 张，无重复帧）标定得出。三条独立证据都指向它：① 数据一致性 train/test/mix 收敛 Δfx<1px（第 10 步）；② 三路不确定度（8 组/10 折/CRLB）同量级；③ 独立测试各组内参泛化一致（第 9 步）。它就是这台相机的内参物理常量。

θ̂_final 具体怎么算出来（逐步）：

1. 输入：cal_data_final/ 里 80 张扁平图（每位姿挑第 1 帧，无重复帧；由 run_pipeline 自动从 cal_data 构建）。每张图 detect_corners → 88 个亚像素角点。
2. 初值：从 repeat.json 读 theta_mean（第 2 步 8 组重复标定的均值，fx≈1855.053）。给初值是为了开 CALIB_USE_INTRINSIC_GUESS——大图数无初值时 LM 容易发散或超慢。
3. 联合估计（核心）：cv2.calibrateCamera(objpoints, imgpoints, (W,H), K0, dist0, flags=CALIB_USE_INTRINSIC_GUESS) 同时估 9 个内参 [fx,fy,cx,cy,k1,k2,p1,p2,k3]（5 参数 Brown 模型：径向 k1-3 + 切向 p1p2）+ 每图 6 个外参 [rvec,tvec]（80 图 × 6 = 480 个），用 LM 迭代最小化全部角点重投影残差平方和 Σ‖proj(θ, 外参) − detected‖²（高斯噪声下 = 最大似然估计）。共 9+480=489 个参数。
4. 结果：LM 收敛 → θ̂_final = [fx=1855.058, fy=1845.492, cx=1269.49, cy=725.02, k1=−0.40515, k2=0.23024, p1=7.4e−5, p2=−2.1e−5, k3=−0.08118]，重投影 rms=0.1379 px。
5. CRLB（理论噪声底，θ̂_final 处算）：数值雅可比 J（投影函数对 9 内参 + 480 外参的中心差分偏导，步长 1e-6·max(1,|p|)）→ 信息矩阵 I = Jᵀ Σ⁻¹ J（各向同性 Σ=σ_a²·I，σ_a=rms/√2）→ Schur 补 I_β = A − B·C⁻¹·Bᵀ 把 480 个外参 nuisance 边缘化掉，得 9×9 纯内参信息 → CRLB_k = √[inv(I_β)]_kk。 fx 的 CRLB=0.054 px（“仅检测噪声决定的理论下界”，实测重复性 0.28 px 是它的 ~5×，差额是位姿/系统贡献）。

一句话：80 张图 + theta_mean 初值 → calibrateCamera 联合解 9 内参 + 480 外参（最小二乘=最大似然）→ θ̂_final；再在 θ̂_final 处算雅可比 → Schur 边缘化外参 → CRLB。下游 final_result 把它和 8 组 std、10 折 std 三路合成"保证 ±0.4 px"的精度结论。

9.2 先把"精度"拆成两件不同的事

精度有两个含义，混为一谈会过度承诺：

• 重复性 （repeatability / precision）：重新标定一遍，数值变多少？——能保证，9.3 给四条证据。
• 准确性 （accuracy）：估值离物理"真值"多远？——部分保证，受靶标品质等系统因素约束（9.4）。

9.3 重复性——四条独立证据，从乐观到保守

内参	CRLB（噪声理论底）	8 组 std（同位姿）	10 折 std（换位姿子集）	跨数据集 Δ（换整批）	保证 ±（Δ/2）
fx	0.054	0.077	0.279	0.79	±0.4 px （0.02%）
fy	0.054	0.072	0.273	1.09	±0.55 px （0.03%）
cx	0.009	0.052	0.333	1.87	±0.9 px （0.07%）
cy	0.007	0.075	0.174	1.33	±0.7 px （0.09%）
k1	3.0e−5	4.4e−5	2.1e−4	0.0073	±0.004 （0.9%）
k2	3.9e−5	1.1e−4	5.4e−4	0.0397	±0.020 （8.6%）
k3	4.6e−5	8.8e−5	4.7e−4	0.0485	±0.024 （30%）

四条证据从窄到宽，每条多算一类误差源：

• CRLB（crlb.json/full.json）：仅检测噪声决定的理论 std 下界（第 3 步）。乐观——只算随机噪声，实际达不到。
• 8 组 std（repeat.json）：同 80 位姿、换不同帧（只换噪声实例）。重复性最窄。
• 10 折 std（cv5.json）：换位姿子集（72 张/折 × 80 折）。含"用哪些位姿"的敏感性 → 比 8 组大 ~4×。
• 跨数据集 Δ（consistency.json）：train/test/mix 三套独立来源（含不同时段 cal_test）的 max−min。最宽、最。

保证声明（取最保守的跨数据集证据）：

fx = 1855.058 ± 0.4 px（重复性，0.02%） —— 任意一次独立重标定（不同帧/位姿/时段）的 fx 都落在 1855.058 ± 0.4 内；1σ ≈ 0.28 px（10 折 std）。fx 相对精度 0.02%，远优于工程判据 0.1–0.3%。

怎么理解 CRLB（0.054）与实际（0.28）差 ~5×：CRLB 只算噪声，实际重复性还含"位姿子集变化 + 组间微动"贡献（第 6 步的 std/CRLB≈5 同源）。这 5× gap 是"系统/位姿贡献"，不是算法不行，是数据决定（第 7 步加图降不动）。

9.4 准确性——能保证什么，不能保证什么

• 投影/使用精度（可保证，有证据）：用 θ̂_final 投影，独立测试集 cal_test 上重投影 RMS = 0.27 px （第 9 步）→ 任意新场景的点投影落在真实位置 ~0.27 px RMS 内（p95 ~0.5 px）。这是应用（3D 测量/视觉/AR）直接关心的"用起来多准"，可保证。
• 绝对准确性（部分保证）：
  • bias_u/v ≈ 0（第 9 步）→ 无系统投影偏置；
  • 数据无关性（train/test/mix 收敛，第 10 步）→ 估值稳定，是物理真值的最佳估计；
  • 但绝对准确性受 靶标平整度（板非理想共面 → 内参偏）+ 镜头模型残余（k1-3 之外的高阶，第 9 步 9.4 节残差 0.27 px）约束，残留 ~0.27 px 量级系统误差。
• 不能从数据内完全保证的：没有计量级认证参考靶（已知精确角度/尺寸）→ 无法分离"稳定估值" 与"物理真值"。注：fx/fy/cx/cy 是像素量，与方格标称边长无关（边长只定外参距离尺度），但受 方格不均匀 + 板不平整 bias。要更高绝对精度须换认证靶 + 剩余像差更小的镜头。

9.5 一句话最终结论

最终内参 θ̂_final = {fx=1855.058， fy=1845.492， cx=1269.49， cy=725.02， k1=−0.4051， k2=0.2302， k3=−0.0812， p1≈p2≈0}；重复性 fx ±0.4 px（0.02%，跨独立数据集保证，1σ≈0.28 px）；投影精度 0.27 px RMS （独立测试集）。精度被硬件封顶（镜头 275 px 主畸变已被模型吃掉，剩 <1 px 是靶标/高阶/位姿残余），算法已随机误差达到下界——要更准须换硬件（认证靶 + 更平镜头），不是加图或改代码。

9.6 进一步提升精度的途径

结论既然是「卡硬件」，优化就只能朝硬件使劲。下列方向按预期收益从大到小排：

1. 换玻璃镀铬靶（收益最大）：纸/亚克力靶面不平、格不准，是 e_target 的主因；玻璃镀铬面平如镜、方格精度高，可直接砍掉靶标系统误差这一项。
2. 采集时固定棋盘：消除「静止连拍存在微动」（§三.5 的 20% 微动组），让 σ_det 测得更真、 cx/cy 散布更接近 CRLB。
3. 各向异性 Σ：Fisher 用 $Σ=diag(σ_u²,σ_v²)$（σ_v≈2σ_u），CRLB 略升但更（第 3 步已支持）。
4. 量准方格边长：边长只定外参尺度、不影响 fx/fy/cx/cy（像素量），但量准能改善 3D 测量的绝对尺度。
5. 去掉切向 p1p2：多折已证 p1p2≈0（第 4 步），CALIB_ZERO_TANGENT_DIST 减 2 个弱约束参数。

注意：加畸变参数（k4-6/s1s2）不在列——第 4 步多折已证会过拟合或不收敛，RMS 不降反升。

十、完整标定流程

把 §五 的各环节收敛成一个可运行的完整流程：在训练图上标定内参并估计不确定度，在独立测试集上评估泛化精度，并给出训练与测试的残差图。测试集不参与标定，仅检验内参对未见数据的预测能力。

10.1 算法

整个流程分七步，每步有明确的输入/输出与原理：

1. 角点检测与亚像素细化：findChessboardCorners 在二值化图上找黑白方格交汇的角点拓扑 → cornerSubPix 在灰度图上拟合二次曲面把角点精化到亚像素（精度 ≈ σ_det，§1.1）。
2. 3D–2D 对应：标定板内角点按已知拓扑生成 3D 物点（板坐标 Z=0 网格，方格边长 scale），与检测到的 2D 角点按相同行列配对。
3. 内参 + 外参联合估计：cv2.calibrateCamera（LM）联合解 9 内参 + 每图 6 外参，最小化重投影残差平方和（高斯噪声下 = 最大似然，§1.2）。外参并非预先已知，而是和内参一起被解出（初值由 Zhang 单应矩阵法给出，再 LM 精化）；图多（≥60）时先用一部分估初值、再用 CALIB_USE_INTRINSIC_GUESS 全量优化。
4. 不确定度（k 折交叉标定散布）：把训练图随机分 k 折，每折用其余 k−1 份标定 → θ̂_fold，跨折 std = 经验不确定度。它不依赖噪声模型，对照 Fisher CRLB（§1.3 理论值）。
5. 训练重投影 RMS + 残差图：用标定内参 + 各图外参投影回像素、与检测角点比，算每图 RMS；逐角点残差按像素位置分块平均 → 误差级别图（看中心 vs 边缘）。
6. 测试集泛化评估：固定训练得到的内参，对每张测试图解 PnP（只解外参）→ 重投影 RMS = 泛化精度；统计 bias（系统偏置）、std_u/v（各向异性）。
7. 测试残差光流图：测试集逐角点残差按像素位置画箭头，看系统误差的空间结构。

10.2 输入与输出

输入：训练图目录（标定用）+ 标定板配置（棋盘规格 + 方格边长）+（可选）独立测试图目录（不参与标定，仅评估泛化）。

输出：内参 θ̂ + 训练 RMS；不确定度（k 折 ±std）+ 每图 RMS 分布；训练残差图（误差级别）、测试残差光流图；测试泛化 RMS、bias、std；结果 JSON。

10.3 完整源码

#!/usr/bin/env python3
# -*- coding: utf-8 -*-
"""完整标定流程:训练图标定内参 + 不确定度估计 + 独立测试集泛化评估 + 残差图。

输入:
  dir            标定图目录(训练,递归读图)
  --board        标定板配置 YAML
  --test-dir     独立测试图目录(可选;不参与标定,仅评估泛化)

输出:
  1. 内参 θ̂(fx,fy,cx,cy,k1,k2,p1,p2,k3)+ 重投影 RMS
  2. 不确定度:k 折交叉标定散布(±std)、每图 RMS 分布(p50/p95/max)
  3. 训练残差图:全图误差级别示意图(残差按像素分块平均 → 分级热力图)
  4. 测试评估(有 --test-dir 时):固定内参对测试集解 PnP → 泛化 RMS、bias、测试残差光流图
  5. 结果 JSON(calibration_result.json)

图多(≥60)时先用一部分估初值,再用 USE_INTRINSIC_GUESS 全量优化(避免大图数 LM 发散/超慢)。
"""
from __future__ import annotations

if __name__ == "__main__" and __package__ is None:
    import sys
    from pathlib import Path
    sys.path.insert(0, str(Path(__file__).resolve().parent.parent.parent))

import argparse
import json
import time
from pathlib import Path

import cv2
import numpy as np
import matplotlib
matplotlib.use("Agg")
import matplotlib.pyplot as plt
from matplotlib.colors import ListedColormap, BoundaryNorm

from core.base.plot_style import use_cjk
use_cjk()

from core.base.calibrate import PARAM_NAMES, CalibrationResult, calibrate_camera
from core.base.config import load_board
from core.base.detect import detect_corners

_IMG_EXT = {".jpg", ".jpeg", ".png", ".bmp", ".tif", ".tiff"}


def calibrate_with_guess(board, corners, shape, theta0):
    """给初值的 5 参数标定(USE_INTRINSIC_GUESS):大图数收敛快、不发散。"""
    K0 = np.array([[theta0[0], 0, theta0[2]], [0, theta0[1], theta0[3]], [0, 0, 1.]],
                  dtype=np.float64)
    dist0 = np.array(theta0[4:9], dtype=np.float64)
    objp = board.object_points().astype(np.float32).reshape(-1, 1, 3)
    objpoints = [objp for _ in corners]
    imgpoints = [np.asarray(c).astype(np.float32).reshape(-1, 1, 2) for c in corners]
    H, W = shape
    rms, K, dist, rvecs, tvecs = cv2.calibrateCamera(
        objpoints, imgpoints, (W, H), K0, dist0, flags=cv2.CALIB_USE_INTRINSIC_GUESS)
    return CalibrationResult(K=K, dist=np.asarray(dist).flatten()[:5], rvecs=list(rvecs),
                             tvecs=list(tvecs), rms=float(rms), n_images=len(corners),
                             image_shape=tuple(shape), message=f"rms={rms:.4f}px")


def kfold_uncertainty(board, corners, shape, theta0, k=5, seed=42):
    """k 折交叉标定:每折用 (k-1)/k 的图标定 → θ̂;跨折 std(ddof=1)= 经验不确定度。"""
    n = len(corners)
    if n < k * 6:
        return None, None
    rng = np.random.default_rng(seed)
    idx = rng.permutation(n)
    folds = np.array_split(idx, k)
    thetas = []
    for fi in range(k):
        val = set(folds[fi].tolist())
        train = [corners[i] for i in range(n) if i not in val]
        if len(train) < 6:
            continue
        try:
            cal = calibrate_with_guess(board, train, shape, theta0)
            thetas.append(cal.theta)
        except Exception:
            continue
    if len(thetas) < 2:
        return None, None
    T = np.array(thetas)
    return T.std(0, ddof=1), T.mean(0)


def per_image_rms(board, corners, cal):
    """每图重投影 RMS(用标定内参 + 该图外参)。"""
    K, dist = cal.K, cal.dist
    objp = board.object_points().astype(np.float32).reshape(-1, 1, 3)
    out = []
    for i, c in enumerate(corners):
        proj, _ = cv2.projectPoints(objp, cal.rvecs[i], cal.tvecs[i], K, dist)
        d2 = ((proj.reshape(-1, 2) - np.asarray(c).reshape(-1, 2)) ** 2).sum(1)
        out.append(float(np.sqrt(d2.mean())))
    return np.array(out)


def error_level_map(board, corners, cal, shape, save, grid=(20, 20)):
    """训练残差图:残差按像素位置分块平均 → 分级热力图。

    级别(相对全局 RMS r):优秀 <0.5r | 良好 0.5-1.5r | 一般 1.5-3r | 差 >3r。
    """
    K, dist = cal.K, cal.dist
    objp = board.object_points().astype(np.float32).reshape(-1, 1, 3)
    px, py, mag = [], [], []
    for i, c in enumerate(corners):
        proj, _ = cv2.projectPoints(objp, cal.rvecs[i], cal.tvecs[i], K, dist)
        e = proj.reshape(-1, 2) - np.asarray(c).reshape(-1, 2)
        det = np.asarray(c).reshape(-1, 2)
        px.append(det[:, 0]); py.append(det[:, 1]); mag.append(np.sqrt((e ** 2).sum(1)))
    px = np.concatenate(px); py = np.concatenate(py); mag = np.concatenate(mag)
    H, W = shape
    nr, nc = grid
    cell = np.full((nr, nc), np.nan)
    iu = np.clip(np.digitize(px, np.linspace(0, W, nc + 1)) - 1, 0, nc - 1)
    iv = np.clip(np.digitize(py, np.linspace(0, H, nr + 1)) - 1, 0, nr - 1)
    r = cal.rms
    for i in range(nr):
        for j in range(nc):
            m = (iv == i) & (iu == j)
            if m.sum():
                cell[i, j] = mag[m].mean()
    lvl = np.full((nr, nc), -1)
    lvl[(cell >= 0) & (cell < 0.5 * r)] = 0
    lvl[(cell >= 0.5 * r) & (cell < 1.5 * r)] = 1
    lvl[(cell >= 1.5 * r) & (cell < 3 * r)] = 2
    lvl[cell >= 3 * r] = 3

    fig, axes = plt.subplots(1, 2, figsize=(15, 6))
    ax = axes[0]
    im = ax.imshow(cell, extent=[0, W, H, 0], aspect='equal', cmap='RdYlGn_r',
                   vmin=0, vmax=max(np.nanmax(cell), 3 * r))
    ax.set_title(f"训练残差热力图(块平均,px)\n全局 RMS={r:.4f}px  块 max={np.nanmax(cell):.3f}px")
    ax.set_xlabel("u (px)"); ax.set_ylabel("v (px)")
    plt.colorbar(im, ax=ax, fraction=0.046, label="块平均残差 (px)")
    ax = axes[1]
    cmap = ListedColormap(["#2ecc71", "#f1c40f", "#e67e22", "#e74c3c"])
    norm = BoundaryNorm([-0.5, 0.5, 1.5, 2.5, 3.5], cmap.N)
    ax.imshow(np.where(lvl >= 0, lvl, np.nan), extent=[0, W, H, 0], aspect='equal',
              cmap=cmap, norm=norm)
    tot = (lvl >= 0).sum()
    pcts = [(lvl == lv).sum() / tot * 100 for lv in range(4)]
    names = ["优秀(<0.5r)", "良好(0.5-1.5r)", "一般(1.5-3r)", "差(>3r)"]
    legend = "\n".join(f"{cmap.colors[i]} ■ {names[i]}: {pcts[i]:.1f}%" for i in range(4))
    ax.set_title("误差级别分布\n" + legend, fontsize=10, loc='left')
    ax.set_xlabel("u (px)"); ax.set_ylabel("v (px)")
    ax.set_xticks([]); ax.set_yticks([])
    plt.tight_layout()
    plt.savefig(save, dpi=110)
    plt.close()
    return dict(cell_max=float(np.nanmax(cell)), pcts={names[i]: float(pcts[i]) for i in range(4)})


def evaluate_test(board, test_corners, cal, shape, save_flow=None, scale=50.0):
    """测试集泛化评估:固定内参 cal,对每张测试图解 PnP(solvePnP)→ 重投影残差。

    返回泛化 RMS、每图 RMS 分布、bias/std(u/v);可选画测试残差光流图(quiver)。
    """
    K, dist = cal.K, cal.dist
    objp = board.object_points().astype(np.float32).reshape(-1, 1, 3)
    all_e, per, px, py, dx, dy = [], [], [], [], [], []
    for c in test_corners:
        pts = np.asarray(c).astype(np.float32).reshape(-1, 1, 2)
        ok, r, t = cv2.solvePnP(objp, pts, K, dist, flags=cv2.SOLVEPNP_ITERATIVE)
        proj, _ = cv2.projectPoints(objp, r, t, K, dist)
        e = proj.reshape(-1, 2) - pts.reshape(-1, 2)
        all_e.append(e); per.append(float(np.sqrt((e ** 2).sum(1).mean())))
        pts2 = pts.reshape(-1, 2)
        px.append(pts2[:, 0]); py.append(pts2[:, 1]); dx.append(e[:, 0]); dy.append(e[:, 1])
    E = np.concatenate(all_e, 0)
    rms = float(np.sqrt((E ** 2).sum(1).mean()))
    if save_flow:
        px = np.concatenate(px); py = np.concatenate(py)
        dx = np.concatenate(dx); dy = np.concatenate(dy)
        mag = np.sqrt(dx ** 2 + dy ** 2)
        H, W = shape
        fig, ax = plt.subplots(figsize=(10, 5.6))
        ax.quiver(px, py, dx * scale, dy * scale, mag, angles='xy', scale_units='xy',
                  scale=1, cmap='viridis', width=0.003, alpha=0.85)
        ax.set_xlim(0, W); ax.set_ylim(H, 0); ax.set_aspect('equal')
        ax.set_title(f"测试集残差光流(箭头 ×{scale:.0f},泛化 RMS={rms:.3f}px)")
        ax.set_xlabel("u (px)"); ax.set_ylabel("v (px)")
        plt.tight_layout(); plt.savefig(save_flow, dpi=110); plt.close()
    return dict(rms=rms, n=len(test_corners),
                per_mean=float(np.mean(per)), per_p50=float(np.percentile(per, 50)),
                per_p95=float(np.percentile(per, 95)), per_max=float(np.max(per)),
                bias_u=float(E[:, 0].mean()), bias_v=float(E[:, 1].mean()),
                std_u=float(E[:, 0].std()), std_v=float(E[:, 1].std()))


def main():
    ap = argparse.ArgumentParser(description="完整标定流程:标定 + 不确定度 + 测试评估 + 残差图")
    ap.add_argument("dir", help="标定图目录(训练,递归读图)")
    ap.add_argument("--board", default=None, help="标定板配置 YAML")
    ap.add_argument("--test-dir", default=None, help="独立测试图目录(不参与标定,仅评估泛化)")
    ap.add_argument("--save-dir", default="outputs", help="输出目录")
    ap.add_argument("--n-folds", type=int, default=5, help="不确定度的折数")
    ap.add_argument("--grid-rows", type=int, default=20)
    ap.add_argument("--grid-cols", type=int, default=20)
    args = ap.parse_args()

    def detect_all(d):
        imgs = sorted(p for p in Path(d).rglob("*") if p.suffix.lower() in _IMG_EXT)
        cs, sh = [], None
        for p in imgs:
            r = detect_corners(str(p), board)
            if r.success:
                cs.append(r.corners); sh = r.image_shape
        return cs, sh, len(imgs)

    board = load_board(args.board)
    print(f"板子:{board.board_type} {board.cols}×{board.rows}  方格 {board.square_size_mm}{board.unit}")
    corners, shape, n_img = detect_all(args.dir)
    print(f"训练:读 {n_img} 张 → 有效 {len(corners)} 张  shape={shape}")
    if len(corners) < 8 or shape is None:
        ap.error("有效训练图太少(<8)或无 shape")

    test_corners = []
    if args.test_dir:
        test_corners, tsh, n_test = detect_all(args.test_dir)
        print(f"测试:读 {n_test} 张 → 有效 {len(test_corners)} 张(不参与标定)")
    print()

    # 1. 内参标定
    print("== 1. 内参标定(5 参数 Brown:fx,fy,cx,cy,k1,k2,p1,p2,k3)==")
    t0 = time.monotonic()
    if len(corners) >= 60:
        cal0 = calibrate_camera(board, corners[:60], shape)
        cal = calibrate_with_guess(board, corners, shape, cal0.theta)
        print(f"  图较多({len(corners)}),用前 60 张估初值后全量标定(USE_INTRINSIC_GUESS)")
    else:
        cal = calibrate_camera(board, corners, shape)
    theta = cal.theta
    print(f"  完成 {time.monotonic()-t0:.0f}s  训练 RMS={cal.rms:.4f}px")
    print(f"  θ̂ = {dict(zip(PARAM_NAMES, [f'{x:.6g}' for x in theta]))}")

    # 2. 不确定度(k 折交叉标定散布)
    print(f"\n== 2. 不确定度({args.n_folds} 折交叉标定散布)==")
    std, _ = kfold_uncertainty(board, corners, shape, theta, k=args.n_folds)
    if std is None:
        print("  图太少,无法算 k 折不确定度")
    else:
        print(f"   {'内参':>4}  {'θ̂':>14}  {'±std(k折)':>12}  {'cv(%)':>8}")
        for k_, name in enumerate(PARAM_NAMES):
            cv = std[k_] / abs(theta[k_]) * 100 if abs(theta[k_]) > 1e-12 else float('nan')
            print(f"   {name:>4}  {theta[k_]:>14.6g}  {std[k_]:>12.4g}  {cv:>8.4f}")

    # 3. 训练每图 RMS 分布 + 残差图
    per = per_image_rms(board, corners, cal)
    print(f"\n== 3. 训练重投影 RMS(每图)==")
    print(f"  mean={per.mean():.4f}  p50={np.percentile(per,50):.4f}  "
          f"p95={np.percentile(per,95):.4f}  max={per.max():.4f}px")

    save_dir = Path(args.save_dir); save_dir.mkdir(parents=True, exist_ok=True)
    train_png = save_dir / "error_level_map.png"
    lvl_info = error_level_map(board, corners, cal, shape, train_png,
                               grid=(args.grid_rows, args.grid_cols))
    print(f"  训练残差图:{train_png.name}(块 max={lvl_info['cell_max']:.3f}px)")

    # 4. 测试集泛化评估
    out = {
        "param_names": PARAM_NAMES, "theta": theta.tolist(),
        "theta_std_kfold": std.tolist() if std is not None else None,
        "rms_train": cal.rms, "n_train": len(corners), "shape": list(shape),
        "per_image_rms_train": {"mean": float(per.mean()), "p50": float(np.percentile(per, 50)),
                                "p95": float(np.percentile(per, 95)), "max": float(per.max())},
        "error_level": lvl_info,
    }
    if test_corners:
        print(f"\n== 4. 测试集泛化评估(固定内参 → solvePnP → 重投影)==")
        test_png = save_dir / "test_residual_flow.png"
        te = evaluate_test(board, test_corners, cal, shape, test_png)
        print(f"  泛化 RMS={te['rms']:.4f}px  (每图 mean={te['per_mean']:.4f}, "
              f"p50={te['per_p50']:.4f}, p95={te['per_p95']:.4f}, max={te['per_max']:.4f})")
        print(f"  bias: u={te['bias_u']:+.4f}  v={te['bias_v']:+.4f}px  "
              f"std: u={te['std_u']:.4f}  v={te['std_v']:.4f}px")
        print(f"  测试残差光流图:{test_png.name}")
        print(f"  → 训练 RMS {cal.rms:.4f} < 测试 {te['rms']:.4f}:差值是泛化 gap(测试集位姿/难度不同),")
        print(f"    bias≈0 说明无系统投影偏置;std_v>std_u 印证各向异性。")
        out["test_eval"] = te
    else:
        print("\n(未给 --test-dir,跳过测试集评估)")

    (save_dir / "calibration_result.json").write_text(
        json.dumps(out, indent=2, ensure_ascii=False), encoding="utf-8")
    print(f"\n已保存 → {save_dir}/calibration_result.json")


if __name__ == "__main__":
    main()

10.4 结果

训练（80 张，每位姿 1 张）：

内参	θ̂	±std（5 折）	cv（%）
fx	1855.058	±0.483	0.0260
fy	1845.492	±0.462	0.0250
cx	1269.49	±0.420	0.0331
cy	725.02	±0.329	0.0453
k1	−0.40515	±4.5e−4	0.110

训练 RMS = 0.1379 px（每图 p50=0.133，p95=0.176，max=0.218）。

测试（64 张独立位姿，固定内参 → solvePnP）：泛化 RMS = 0.2717 px（每图 p50=0.207，p95=0.438，max=0.720）；bias u=+0.0000、v=−0.0002 px（无系统偏置）；std u=0.179、v=0.205 px（v>u，各向异性，印证 σ_v>σ_u）。

训练 RMS 0.138 < 测试 0.272：差值是泛化 gap（测试集位姿/难度不同，并非内参变差）；bias≈0 说明无系统投影偏置。结果与 §九 一致——内参是可靠、可泛化的物理量，精度被硬件封顶。

小结

整条评估链收敛到一句话：这台相机+靶，fx 相对精度 0.02%（跨独立数据集可保证），重投影 RMS 0.14 px 是硬件天花板——随机误差已达下界，内参经 10 步验证是可靠的物理常量，要更准须换硬件，不是加图或改代码。

回到开头的叙事：仿真篇留的缺口是「真实标定没有真值，第 6 步的 Monte Carlo 验证做不了」。本文用三件事把它补上——k 折交叉标定（替代换噪声种子，给经验不确定度）、cal_test 独立测试（完全不参与训练的泛化金标准）、train/test/mix 数据一致性（证明换整批数据内参也不变）。没有真值，但换几个角度交叉验证，依然能把「这台相机估得多准、卡在哪」说清楚：随机精度可保证，绝对精度被硬件封顶。

参考资料

• 完整代码仓库：github.com/zywvvd/camera_calibration（采图 + 标定 + Fisher 信息精度分析，一条龙）
• OpenCV calibrateCamera / solvePnP / projectPoints 文档
• Zhang Z. A Flexible New Technique for Camera Calibration， IEEE TPAMI 2000
• 理论推导：《Fisher 信息与精度下界》
• 仿真验证：《Fisher 信息：相机标定仿真实操》

文章链接：
https://www.zywvvd.com/notes/study/camera-imaging/fisher-information-real/fisher-information-real/