wandering1
40a2816a99
- UN 通道 5 点最小二乘标定完成：r=555.382141, b=-1, 残差 0.038%, R²=0.99999967 - UA/IA/UN 三路 VALID 置位，保护解卡 - ROM 默认值更新为当前标定值（ADC_CAL_R_Q16_VOLT=555.144958） - UB/UC/UN 通过 104 写入 ROM 默认值，config.ini 同步 - cal_fit.py 新增 ROLD 表自动读取，不再需要手传 --r-old - 执行卡精简至 126 行，新增 cal_collect.py 采集脚本与 ROLD 填表 - 标定脚本（cal_fit/cal_write/cal_collect）归集到 test/host/cal/ - 全库文档路径引用同步更新（8 文件） - UN 当前态文档改写：明确零序电压 U0 语义、本台架未接线状态 - 删除已失效的外部知识底座文档链接
559 lines
23 KiB
Python
559 lines
23 KiB
Python
#!/usr/bin/env python3
|
||
# -*- coding: utf-8 -*-
|
||
"""N 点最小二乘标定求解器：**全仓唯一实现**（v2.2 起固件侧拟合已删除）。
|
||
|
||
为什么是唯一实现：标定是台架一次性离线动作——注入已知量 → 读屏上 `meas:` →
|
||
在 PC 上解出 r/b → 经 104 写入设备。设备本体只需要接受校正后的参数，不需要
|
||
自己拟合曲线。固件里曾有过一份 C 实现，生产代码零调用（被 linker 移除），
|
||
留着只会变成第二个真相源，v2.2 已整段删除。所以**本脚本的整数算术就是口径**，
|
||
`--selftest` 是它的回归锚点（原 C 侧金标向量已搬进来）。
|
||
|
||
把台架注入的一串 (屏上读数, 注入已知量) 点拟合成 y = a*x + c（普通最小二乘），
|
||
再折回标定层能落盘的 (r_q16, b)：
|
||
|
||
r_new_q16 = round(a * r_old_q16) # 一次项 = 端子↔引脚变比
|
||
c_pin_mV = c * 65536 / r_new_q16 # 截距从二次值域折回引脚域
|
||
b_new = b_old + round(c_pin_mV * 65536 * 1000 / k_q16)
|
||
|
||
推导（disp 单位 mV）：折算模型 disp = (code + b) × (k_q16/65536) × 1e-3 ×
|
||
(r_q16/65536)，把 disp_new = a×disp_old + c 的 code 一次项与常数项分别对齐，
|
||
一次项给出 r_new = a·r_old，常数项给出上式 b_new（b_old = 写之前设备里的零点，
|
||
必须作为入参传进来，因为读数是在那个零点下读的）。
|
||
|
||
为什么用两点之外的算法：两点没有冗余，单点读数一抖解出的 r 就整条平移；N 点最小二乘
|
||
把残差摊到所有点上，并且能顺手判「非线性/削顶/接线反」——这正是两点的盲区。
|
||
|
||
全程 int64 定点（不用 float 做累加）：结果确定、可跨平台复现，同一组输入永远
|
||
解出同一个 r_q16。整数路径的圆整次序是本脚本语义的一部分，改一处就要同步
|
||
--selftest 的期望值。
|
||
|
||
用法（仓库根目录，解释器用 .venv）：
|
||
|
||
自检（回归锚点，改动算法后必跑）：
|
||
.venv/Scripts/python.exe test/host/cal_fit.py --selftest
|
||
|
||
单路（x = `meas:` 行原样读数 mV/mA；y = 继保仪注入值 V/A，内部 ×1000）：
|
||
.venv/Scripts/python.exe test/host/cal_fit.py --sem ua \\
|
||
--x 5021.5,10041.77,25000,50024.09,99933.02 --y 5,10,25,50,100
|
||
|
||
批量（读执行卡 <!-- CALDATA:BEGIN/END --> 之间的长表，每行一个点）：
|
||
.venv/Scripts/python.exe test/host/cal_fit.py --from-doc \\
|
||
docs/测试说明/专项验证/ADC多点标定执行流程.md
|
||
|
||
公共参数：
|
||
--r-old 本次读数采集时设备里的变比（浮点）：默认 ua/un=1000.0、ia=29.4118；
|
||
也可 ua=1000,un=1000,ia=29.4118 形式逐路给
|
||
--b-old 本次读数采集时设备里的零点（整数码值，默认 0）：裸值或 ua=0,ia=-2 形式
|
||
|
||
⚠ **同一路的 x 必须在同一个 r 下采集**：拿旧 r 时代的读数和新 r 时代的读数混着
|
||
做直线拟合，斜率会被系统性掰弯（表现是解出的 r 差几十个百分点，而合法带仍然放行）。
|
||
|r_new/r_old - 1| > 10% 时打 WARN 就是给这件事提个醒。
|
||
|
||
判据（任一不满足即 REJECT，打印原因，不算可写）：
|
||
n 不在 [3, 8] 冗余不足或超上限
|
||
den <= 0 示值全同 / 点位没拉开
|
||
x 与 y 相邻增量方向相反 接线反、通道串、削顶的典型征兆
|
||
a <= 0 零斜率/负斜率
|
||
r_new_q16 不在 [1.0, 32767.99] 与 run_par_cal_r_set() 的粗守同带
|
||
|b_new| > 8192 零点越带（先查前端，不用系数凑）
|
||
最大残差 > 满量程跨度 2% 非线性/削顶
|
||
|
||
输出纯 ASCII（中文只在本注释里）；exit 0 = 全部解出且合法，2 = 有拒绝/未填行；
|
||
--selftest 的 exit 0 = 全部断言通过，1 = 有失败。
|
||
创建日期：2026-11-13（两点 → N 点最小二乘；同日 v2.2 固件侧拟合删除，本脚本转为唯一实现）
|
||
"""
|
||
import argparse
|
||
import re
|
||
import sys
|
||
|
||
# ---- 判据参数（本脚本是唯一实现，这些值即为口径；改动会使 GOLDEN 期望值失效）----
|
||
# ⚠ 前 5 个是**独立于固件**的算法参数（原先与已删除的 ADC_CAL_FIT_* 同值，
|
||
# 现在以本文件为准）；后 4 个是设备侧合法带，必须与 App/algo/adc_cal.h /
|
||
# run_par.c 保持一致——那两处仍是设备侧安全闸（写参数时的拦截），不随拟合删除。
|
||
K_Q16 = 10000000 # 该路前端斜率（UA/UN/IA 三路同源，同值）
|
||
B_MAX = 8192 # 设备侧 ADC_CAL_B_MAX（零点幅值上限，码值）
|
||
R_MIN = 65536 # 设备侧 ADC_CAL_R_MIN（比值 1.0）
|
||
R_MAX = 2147483647 # 设备侧 ADC_CAL_R_MAX（int32 上限）
|
||
MILLI_MAX = 200000 # 单点值域上界（mV / mA），±200 V / ±200 A
|
||
C_MAX = 100000 # 截距 c 上界（mV），±100 V
|
||
N_MIN = 3 # 点数下限（少于 3 点无冗余）
|
||
N_MAX = 8 # 点数上限
|
||
RES_PCT = 2 # 最大残差 ÷ 满量程跨度上限（百分比）
|
||
# 比值上界的两位小数**向下取整**得 32767.99，不能用 R_MAX/65536.0 直接 %.2f：
|
||
# 2147483647/65536 = 32767.99998474，四舍五入会印成 32768.00，与 adc_cal.h /
|
||
# run_par.c 的"1.0 ~ 32767.99"口径差一个数（提示文案里的带边界必须和真带一致）。
|
||
R_RATIO_MAX = (R_MAX * 100 // 65536) / 100.0 # 32767.99（仅用于提示文案）
|
||
|
||
# 语义槽位（App/algo/adc_input.h 枚举序：UA0/UB1/UC2/UN3/IA4/IB5/IC6/IN7）
|
||
SEM = {
|
||
"ua": dict(slot=0, r_old=1000.0), # 电压族 TPAFE
|
||
"un": dict(slot=3, r_old=1000.0), # 电压族 TPAFE
|
||
"ia": dict(slot=4, r_old=29.4118), # 电流族 TPAFE
|
||
}
|
||
|
||
# 拒绝原因码（负数）：本脚本的对外契约，--selftest 与 REASON 文案都按这张表
|
||
ERR_ARG = -1
|
||
ERR_POINTS = -2
|
||
ERR_DEN = -3
|
||
ERR_MONO = -4
|
||
ERR_SLOPE = -5
|
||
ERR_R_BAND = -6
|
||
ERR_B_BAND = -7
|
||
ERR_RESID = -8
|
||
|
||
REASON = {
|
||
ERR_ARG: "bad argument (|x|,|y| < " + str(MILLI_MAX) +
|
||
", k_q16 in (0,1e8], r_old in band)",
|
||
ERR_POINTS: "n outside [%d,%d]: got %%d point(s)" % (N_MIN, N_MAX),
|
||
ERR_DEN: "x span collapsed (den<=0): readings not spread",
|
||
ERR_MONO: "non-monotonic: adjacent x/y increments oppose",
|
||
ERR_SLOPE: "slope <= 0: readings flat or reversed vs injected value",
|
||
ERR_R_BAND: "r_new out of band [1.0, %.2f]" % R_RATIO_MAX,
|
||
ERR_B_BAND: "|b_new| beyond %d codes" % B_MAX,
|
||
ERR_RESID: "residual > %d%% of y span (non-linearity / clipping)" % RES_PCT,
|
||
}
|
||
|
||
|
||
# ------------------------------------------------------------ 整数算术（唯一实现）
|
||
# v2.2 起这里不再是"镜像"：App/algo/adc_cal.c 的拟合算术已整段删除，本文件就是口径。
|
||
def div_trunc(num, den):
|
||
"""整数除法：向零截断（Python // 向负无穷取整，不能直接用）。"""
|
||
q = abs(num) // abs(den)
|
||
return q if (num < 0) == (den < 0) else -q
|
||
|
||
|
||
def div_round(num, den):
|
||
"""四舍五入（half away from zero）；den 恒 > 0。
|
||
|
||
这是本仓库拟合算术**唯一的圆整口径**，不要换成 `(num + den/2)//den` 那种
|
||
写法：Python 的 // 对负数是向负无穷取整，负数会在 -0.5 处折向 0，
|
||
与期望的 half-away-from-zero 不同。先判符号再取绝对值，不依赖语言行为。
|
||
"""
|
||
if num >= 0:
|
||
return (num + den // 2) // den
|
||
return -((-num + den // 2) // den)
|
||
|
||
|
||
def round_int(v):
|
||
"""浮点 → 整数：二分之一远离零圆整（入参量化口径，与 to_milli 配对）。"""
|
||
return int(v + 0.5) if v >= 0 else -int(-v + 0.5)
|
||
|
||
|
||
def to_milli(v):
|
||
"""读数 → 整数毫单位（mV / mA 域，即 fit() 的直接入参）。"""
|
||
return round_int(v)
|
||
|
||
|
||
def fit(xs, ys, k_q16=K_Q16, r_old_q16=None, b_old=0):
|
||
"""N 点最小二乘标定求解（唯一实现）。返回 (rc, info)；rc == 0 时 info 是解。
|
||
|
||
xs / ys 已是整数 mV / mA；r_old_q16 = 本轮读数采集时的变比 Q16 原值。
|
||
整数路径决定 rc 与 r_q16/b；info 里的 a/res%/R2 是报告用浮点量，不参与判定。
|
||
"""
|
||
n = len(xs)
|
||
if r_old_q16 is None:
|
||
return ERR_ARG, None
|
||
if (k_q16 <= 0) or (k_q16 > 100000000):
|
||
return ERR_ARG, None
|
||
if (r_old_q16 < R_MIN) or (r_old_q16 > R_MAX):
|
||
return ERR_ARG, None
|
||
if (n < N_MIN) or (n > N_MAX):
|
||
return ERR_POINTS, None
|
||
for v in list(xs) + list(ys):
|
||
if abs(v) >= MILLI_MAX:
|
||
return ERR_ARG, None
|
||
|
||
# 一趟：累加 + 单调性（相邻点方向必须一致）
|
||
sx = sy = sxy = sxx = 0
|
||
for i in range(n):
|
||
x = xs[i]
|
||
y = ys[i]
|
||
sx += x
|
||
sy += y
|
||
sxy += x * y
|
||
sxx += x * x
|
||
if i > 0:
|
||
dx = x - xs[i - 1]
|
||
dy = y - ys[i - 1]
|
||
if (dx > 0 and dy < 0) or (dx < 0 and dy > 0):
|
||
return ERR_MONO, None
|
||
|
||
num = n * sxy - sx * sy
|
||
den = n * sxx - sx * sx
|
||
if den <= 0:
|
||
return ERR_DEN, None
|
||
|
||
a_q16 = div_round(num * 65536, den)
|
||
if a_q16 <= 0:
|
||
return ERR_SLOPE, None
|
||
|
||
# 先按除法把 a_q16 夹到"r_new 不会越上界"的范围，再做 a_q16 × r_old——
|
||
# 顺序反了会在极端点位下把 int64 顶穿（a_q16 可到 1e17 量级）
|
||
if a_q16 > (R_MAX * 65536) // r_old_q16:
|
||
return ERR_R_BAND, None
|
||
r_new_q16 = (a_q16 * r_old_q16 + 32768) >> 16
|
||
if (r_new_q16 < R_MIN) or (r_new_q16 > R_MAX):
|
||
return ERR_R_BAND, None
|
||
|
||
# 截距 c：y = a·x + c ⇒ c = (Σy − a·Σx) / n。两遍式（用已算好的 a_q16）。
|
||
# 若改成从 num/den 代数展开的 (Σy·den − num·Σx)/(n·den)，分子峰值到 1.2e19，
|
||
# 超出 int64 上界——Python 的大整数不会溢出，但这条路径要能直接照抄成
|
||
# 定点实现，所以这里**刻意保持 int64 可容纳的两遍式**，别改成展开式
|
||
c_milli = div_trunc(sy * 65536 - a_q16 * sx, n * 65536)
|
||
if c_milli > C_MAX or c_milli < -C_MAX:
|
||
return ERR_B_BAND, None # 截距域离谱 = b 必然越带（且护住下一步的 ×2^32）
|
||
|
||
# 截距折回引脚域，保留 Q16 亚毫伏精度（截到整 mV 会带来约 6.5 码 / mV 的 b
|
||
# 量化，把真正的细分零点修正抹平）；两步走，别一次乘完（会顶穿 int64）
|
||
c_pin_q16 = div_round(c_milli * 65536 * 65536, r_new_q16)
|
||
db = div_round(c_pin_q16 * 1000, k_q16)
|
||
b_new = b_old + db
|
||
if (b_new > B_MAX) or (b_new < -B_MAX):
|
||
return ERR_B_BAND, None
|
||
|
||
# 残差：|y_i − (a·x_i + c)| 与满量程跨度比
|
||
y_min = min(ys)
|
||
y_max = max(ys)
|
||
y_span = y_max - y_min
|
||
res_max = 0
|
||
for i in range(n):
|
||
pred = div_round(a_q16 * xs[i], 65536) + c_milli
|
||
res = ys[i] - pred
|
||
if res < 0:
|
||
res = -res
|
||
if res > res_max:
|
||
res_max = res
|
||
if (res_max * 100) > (y_span * RES_PCT):
|
||
return ERR_RESID, None
|
||
|
||
# 报告用浮点量（不参与判定）
|
||
a_f = a_q16 / 65536.0
|
||
y_mean = float(sy) / n
|
||
ss_res = 0.0
|
||
ss_tot = 0.0
|
||
for i in range(n):
|
||
d = ys[i] - (a_f * xs[i] + c_milli)
|
||
t = ys[i] - y_mean
|
||
ss_res += d * d
|
||
ss_tot += t * t
|
||
r2 = (1.0 - ss_res / ss_tot) if ss_tot > 0.0 else 0.0
|
||
|
||
return 0, dict(n=n, a_q16=a_q16, a=a_f, c_milli=c_milli,
|
||
r_new_q16=r_new_q16, r_new=r_new_q16 / 65536.0,
|
||
c_pin_mV=c_pin_q16 / 65536.0, db=db, b=b_new,
|
||
res_max=res_max,
|
||
res_pct=(100.0 * res_max / y_span) if y_span else 0.0,
|
||
r2=r2, y_span=y_span,
|
||
dev_pct=100.0 * (r_new_q16 / float(r_old_q16) - 1.0))
|
||
|
||
|
||
# ---------------------------------------------------------------- 输入解析
|
||
def parse_floats(text):
|
||
return [float(s) for s in text.split(",") if s.strip() != ""]
|
||
|
||
|
||
def parse_map(text, defaults, conv):
|
||
"""--r-old / --b-old：裸值作用于所有路，或 "ua=…,ia=…" 逐路给。"""
|
||
if text is None:
|
||
return dict(defaults)
|
||
if "=" not in text:
|
||
v = conv(text)
|
||
return {k: v for k in defaults}
|
||
out = dict(defaults)
|
||
for part in text.split(","):
|
||
part = part.strip()
|
||
if not part:
|
||
continue
|
||
k, v = part.split("=", 1)
|
||
k = k.strip().lower()
|
||
if k not in defaults:
|
||
raise ValueError("unknown sem %r (use ua/un/ia)" % k)
|
||
out[k] = conv(v)
|
||
return out
|
||
|
||
|
||
def parse_cal_table(path):
|
||
"""提取 <!-- CALDATA:BEGIN --> ~ <!-- CALDATA:END --> 之间的长表。
|
||
|
||
列序：路 | 序号 | 屏上读数(mV/mA) | 注入已知量(V/A)，每行一个点。
|
||
返回 {sem: [(idx, x_s, y_s), ...]}；标记不存在返回 None。
|
||
两格全空 = 该点本轮没测（调用方跳过）；只填一格 = 半填错误（报 INCOMPLETE）。
|
||
"""
|
||
with open(path, encoding="utf-8") as f:
|
||
text = f.read()
|
||
m = re.search(r"<!--\s*CALDATA:BEGIN\s*-->(.*?)<!--\s*CALDATA:END\s*-->",
|
||
text, re.S)
|
||
if m is None:
|
||
return None
|
||
groups = {}
|
||
for line in m.group(1).splitlines():
|
||
line = line.strip()
|
||
if not line.startswith("|"):
|
||
continue
|
||
cells = [c.strip() for c in line.strip("|").split("|")]
|
||
sem = cells[0].lower()
|
||
if sem not in SEM:
|
||
continue # 表头行 / 分隔行自然跳过
|
||
if len(cells) < 4:
|
||
raise ValueError("row for %r has %d columns; need 4 "
|
||
"[sem|idx|reading|injected]" % (sem, len(cells)))
|
||
if len(cells) > 4:
|
||
# 典型来源：没清空的**旧两点表**（6 列）。静默按前四列解会解出垃圾，
|
||
# 所以这里明说——比让用户去猜"怎么 REJECT 了"省一轮台架
|
||
raise ValueError("row for %r has %d columns (>4): that is the old "
|
||
"two-point table; the N-point table needs exactly "
|
||
"4 columns [sem|idx|reading|injected]"
|
||
% (sem, len(cells)))
|
||
groups.setdefault(sem, []).append((cells[1], cells[2], cells[3]))
|
||
return groups
|
||
|
||
|
||
def parse_rold_table(path):
|
||
"""提取 <!-- ROLD:BEGIN --> ~ <!-- ROLD:END --> 之间的 r/b 旧值表。
|
||
|
||
列序：路 | r_old | b_old [| 备注]，返回 {sem: {"r": float, "b": int}}。
|
||
标记不存在或表为空返回空 dict（调用方回退到 --r-old / SEM 默认）。
|
||
"""
|
||
with open(path, encoding="utf-8") as f:
|
||
text = f.read()
|
||
m = re.search(r"<!--\s*ROLD:BEGIN\s*-->(.*?)<!--\s*ROLD:END\s*-->",
|
||
text, re.S)
|
||
if m is None:
|
||
return {}
|
||
result = {}
|
||
for line in m.group(1).splitlines():
|
||
line = line.strip()
|
||
if not line.startswith("|"):
|
||
continue
|
||
cells = [c.strip() for c in line.strip("|").split("|")]
|
||
sem = cells[0].lower()
|
||
if sem not in SEM:
|
||
continue
|
||
if len(cells) < 3 or not cells[1] or not cells[2]:
|
||
continue
|
||
try:
|
||
result[sem] = {"r": float(cells[1]), "b": int(float(cells[2]))}
|
||
except ValueError:
|
||
continue
|
||
return result
|
||
def parse_cell(s):
|
||
if s in ("", "-", "—"):
|
||
return None
|
||
return float(s)
|
||
|
||
|
||
# ---------------------------------------------------------------- 求解与输出
|
||
def solve_and_print(sem, xs, ys, r_old_q16, b_old):
|
||
"""解一路并打印判读；成功返回 (sem, slot, r_new, b)，拒绝返回 None。"""
|
||
ch = SEM[sem]
|
||
rc, info = fit(xs, ys, K_Q16, r_old_q16, b_old)
|
||
if rc != 0:
|
||
reason = REASON[rc]
|
||
if rc == ERR_POINTS:
|
||
reason = reason % len(xs)
|
||
print("sem=%s slot=%d n=%d REJECT %s" % (sem, ch["slot"], len(xs), reason))
|
||
return None
|
||
print("sem=%s slot=%d n=%d OK r=%.6f r_q16=%d b=%d a=%.9f res%%=%.3f "
|
||
"R2=%.8f" % (sem, ch["slot"], info["n"], info["r_new"],
|
||
info["r_new_q16"], info["b"], info["a"],
|
||
info["res_pct"], info["r2"]))
|
||
if abs(info["dev_pct"]) > 10.0:
|
||
print("sem=%s WARN r_new dev %+.2f%% vs r_old=%.6f - all points must be "
|
||
"read at the same r, and units are mV vs V"
|
||
% (sem, info["dev_pct"], r_old_q16 / 65536.0))
|
||
return (sem, ch["slot"], info["r_new"], info["b"])
|
||
|
||
|
||
def write_all(done):
|
||
"""末尾合并写入命令：先全部 r，再全部 b（cal_write.py 的两段顺序契约）。"""
|
||
parts = ["--r 0x%04X=%.6f" % (0x8050 + s, r) for _, s, r, _ in done]
|
||
parts += ["--b 0x%04X=%d" % (0x8048 + s, b) for _, s, _, b in done]
|
||
print("WRITE-ALL: test/host/cal_write.py " + " ".join(parts))
|
||
|
||
|
||
def run_single(a):
|
||
missing = [n for n, v in (("--sem", a.sem), ("--x", a.x), ("--y", a.y))
|
||
if v is None]
|
||
if missing:
|
||
print("ERROR: single-route mode needs %s (or use --from-doc)"
|
||
% ",".join(missing))
|
||
return 2
|
||
try:
|
||
xv = parse_floats(a.x)
|
||
yv = parse_floats(a.y)
|
||
r_map = parse_map(a.r_old, {s: SEM[s]["r_old"] for s in SEM}, float)
|
||
b_map = parse_map(a.b_old, {s: 0 for s in SEM}, int)
|
||
except ValueError as exc:
|
||
print("ERROR: %s" % exc)
|
||
return 2
|
||
if len(xv) != len(yv):
|
||
print("ERROR: --x has %d points, --y has %d" % (len(xv), len(yv)))
|
||
return 2
|
||
xs = [to_milli(v) for v in xv] # 屏上读数本来就是 mV/mA
|
||
ys = [to_milli(v * 1000.0) for v in yv] # 注入已知量 V/A → mV/mA
|
||
r_old_q16 = round_int(r_map[a.sem] * 65536.0)
|
||
res = solve_and_print(a.sem, xs, ys, r_old_q16, b_map[a.sem])
|
||
if res is None:
|
||
return 2
|
||
write_all([res])
|
||
return 0
|
||
|
||
|
||
def run_from_doc(a):
|
||
try:
|
||
groups = parse_cal_table(a.from_doc)
|
||
except OSError as exc:
|
||
print("BAD-DOC: %s" % exc)
|
||
return 2
|
||
except ValueError as exc:
|
||
print("BAD-TABLE: %s" % exc)
|
||
return 2
|
||
if groups is None:
|
||
print("NO-MARKER: no <!-- CALDATA:BEGIN/END --> table in %s" % a.from_doc)
|
||
return 2
|
||
if not groups:
|
||
print("EMPTY: CALDATA table has no ua/un/ia rows")
|
||
return 2
|
||
# ROLD table in doc overrides --r-old defaults; CLI --r-old still wins if given explicitly
|
||
rold = parse_rold_table(a.from_doc)
|
||
defaults = {s: SEM[s]["r_old"] for s in SEM}
|
||
if rold:
|
||
for s in rold:
|
||
defaults[s] = rold[s]["r"]
|
||
try:
|
||
r_map = parse_map(a.r_old, defaults, float)
|
||
b_defaults = {s: (rold[s]["b"] if s in rold else 0) for s in SEM}
|
||
b_map = parse_map(a.b_old, b_defaults, int)
|
||
except ValueError as exc:
|
||
print("ERROR: %s" % exc)
|
||
return 2
|
||
|
||
done = []
|
||
failed = 0
|
||
for sem in sorted(groups):
|
||
xs, ys = [], []
|
||
for idx, x_s, y_s in groups[sem]:
|
||
xv = parse_cell(x_s)
|
||
yv = parse_cell(y_s)
|
||
if xv is None and yv is None:
|
||
continue # 整行空 = 本轮没测这个点
|
||
if xv is None or yv is None:
|
||
print("sem=%s idx=%s INCOMPLETE: need both reading and injected "
|
||
"value (or leave both empty); skipped" % (sem, idx))
|
||
failed += 1
|
||
continue
|
||
xs.append(to_milli(xv))
|
||
ys.append(to_milli(yv * 1000.0))
|
||
if not xs:
|
||
print("sem=%s SKIP: no filled points" % sem)
|
||
continue
|
||
r_old_q16 = round_int(r_map[sem] * 65536.0)
|
||
res = solve_and_print(sem, xs, ys, r_old_q16, b_map[sem])
|
||
if res is None:
|
||
failed += 1
|
||
else:
|
||
done.append(res)
|
||
if done:
|
||
write_all(done)
|
||
return 0 if (not failed and done) else 2
|
||
|
||
|
||
# ---------------------------------------------------------------- selftest 回归
|
||
# 金标向量：原固件单测（v2.1 的 tc_00b_fit）里的那一组，v2.2 固件侧拟合删除后
|
||
# 搬到这里。输入与期望值**全部写死**，不在这里现算——现算等于自证，锚点就没意义了。
|
||
#
|
||
# 数据怎么来的（推导留档，便于将来重构时重新生成）：
|
||
# code = 500/700/900/1100/1300，真值 r_true=600.0、b_true=+12，
|
||
# k_q16=1e7（TPAFE 152.587890625 µV/码），r_old=1000.0，b_old=0
|
||
# 引脚电压 = code × 0.152587890625 mV → x（旧 r 下的屏上读数）= 该值 × 1000/1
|
||
# y（注入已知量）= (code + 12) × 0.152587890625 × 600.0
|
||
GOLDEN = {
|
||
"x": [76294, 106812, 137329, 167847, 198364],
|
||
"y": [46875, 65186, 83496, 101807, 120117],
|
||
"k_q16": K_Q16,
|
||
"r_old_q16": 1000 * 65536,
|
||
"b_old": 0,
|
||
"expect_rc": 0,
|
||
"expect_r_q16": 39322000,
|
||
"expect_b": 12,
|
||
}
|
||
|
||
# 拒绝路径断言：(名称, xs, ys, r_old_q16, b_old, 期望 rc, 期望原因前缀)
|
||
NG_CASES = [
|
||
("n<3", [1000, 2000], [1000, 2000], 1000 * 65536, 0, ERR_POINTS, "[3,8]"),
|
||
("x all same", [5000, 5000, 5000], [1000, 2000, 3000],
|
||
1000 * 65536, 0, ERR_DEN, "x span collapsed"),
|
||
("mono reversed", [1000, 2000, 3000, 4000, 5000],
|
||
[1000, 2000, 3000, 2000, 1000], 1000 * 65536, 0, ERR_MONO, "non-monotonic"),
|
||
("r out of band", [1000, 101000, 199999], [1000, 1001, 1002],
|
||
1000 * 65536, 0, ERR_R_BAND, "out of band"),
|
||
]
|
||
|
||
|
||
def run_selftest():
|
||
"""跑金标向量 + 拒绝路径断言。exit 0 = 全过，1 = 有失败。"""
|
||
fails = 0
|
||
total = 0
|
||
|
||
def chk(cond, label):
|
||
nonlocal fails, total
|
||
total += 1
|
||
if cond:
|
||
print("[PASS] %s" % label)
|
||
else:
|
||
print("[FAIL] %s" % label)
|
||
fails += 1
|
||
|
||
g = GOLDEN
|
||
rc, info = fit(g["x"], g["y"], g["k_q16"], g["r_old_q16"], g["b_old"])
|
||
chk(rc == g["expect_rc"], "golden: rc == %d (got %d)" % (g["expect_rc"], rc))
|
||
if info is None:
|
||
chk(False, "golden: solver returned no solution")
|
||
else:
|
||
chk(info["r_new_q16"] == g["expect_r_q16"],
|
||
"golden: r_q16 == %d (got %d)" % (g["expect_r_q16"], info["r_new_q16"]))
|
||
chk(info["b"] == g["expect_b"],
|
||
"golden: b == %d (got %d)" % (g["expect_b"], info["b"]))
|
||
|
||
for name, xs, ys, ro, bo, want_rc, want_txt in NG_CASES:
|
||
rc, _ = fit(xs, ys, K_Q16, ro, bo)
|
||
chk(rc == want_rc, "reject %-14s rc == %d (got %d)" % (name, want_rc, rc))
|
||
chk(want_txt in REASON.get(want_rc, ""),
|
||
"reject %-14s message mentions %r" % (name, want_txt))
|
||
|
||
print("=== cal_fit selftest: %d/%d PASSED (failed=%d) ===" %
|
||
(total - fails, total, fails))
|
||
return 0 if fails == 0 else 1
|
||
|
||
|
||
def main():
|
||
ap = argparse.ArgumentParser(
|
||
description="N-point least-squares calibration solver "
|
||
"(sole implementation; firmware no longer fits)")
|
||
ap.add_argument("--selftest", action="store_true",
|
||
help="run the golden-vector regression (exit 0 = all passed)")
|
||
ap.add_argument("--sem", choices=sorted(SEM), default=None)
|
||
ap.add_argument("--x", default=None,
|
||
help="comma-separated on-screen readings (mV or mA)")
|
||
ap.add_argument("--y", default=None,
|
||
help="comma-separated injected values (V or A; x1000 inside)")
|
||
ap.add_argument("--r-old", default=None,
|
||
help="r in the device while readings were taken "
|
||
"(default ua/un=1000.0, ia=29.4118; or ua=..,ia=..)")
|
||
ap.add_argument("--b-old", default=None,
|
||
help="b in the device while readings were taken "
|
||
"(codes, default 0; or ua=..,ia=..)")
|
||
ap.add_argument("--from-doc", default=None,
|
||
help="read the CALDATA-marked long table from a doc")
|
||
a = ap.parse_args()
|
||
|
||
if a.selftest:
|
||
return run_selftest()
|
||
if a.from_doc is not None:
|
||
return run_from_doc(a)
|
||
return run_single(a)
|
||
|
||
|
||
if __name__ == "__main__":
|
||
sys.exit(main())
|