From c222e054ca39599bd837d3a524fa4e4ae320279a Mon Sep 17 00:00:00 2001 From: Matiq Date: Sat, 22 Aug 2026 12:55:02 +0300 Subject: [PATCH] =?UTF-8?q?docs:=20Phase=20B=20recovery=20=E2=80=94=204=20?= =?UTF-8?q?offline=20detector=20hypotheses=20refuted;=20scripts=20into=20s?= =?UTF-8?q?cripts/?= MIME-Version: 1.0 Content-Type: text/plain; charset=UTF-8 Content-Transfer-Encoding: 8bit --- handoff/NOTES_LEVEL.md | 29 +++++++++ scripts/phaseA_grid_fast.py | 115 ++++++++++++++++++++++++++++++++++++ scripts/phaseB_pool2.py | 56 ++++++++++++++++++ scripts/phaseB_pool3.py | 55 +++++++++++++++++ scripts/phaseB_pooling.py | 110 ++++++++++++++++++++++++++++++++++ scripts/phaseB_rho.py | 71 ++++++++++++++++++++++ scripts/phaseB_temporal.py | 68 +++++++++++++++++++++ 7 files changed, 504 insertions(+) create mode 100644 scripts/phaseA_grid_fast.py create mode 100644 scripts/phaseB_pool2.py create mode 100644 scripts/phaseB_pool3.py create mode 100644 scripts/phaseB_pooling.py create mode 100644 scripts/phaseB_rho.py create mode 100644 scripts/phaseB_temporal.py diff --git a/handoff/NOTES_LEVEL.md b/handoff/NOTES_LEVEL.md index c441851..8957035 100644 --- a/handoff/NOTES_LEVEL.md +++ b/handoff/NOTES_LEVEL.md @@ -1362,3 +1362,32 @@ member -> committed canon stays HEAD (LUT g=0.344/m=4.2). 2. Step 7 live capture ctx+0x188 (constants + possible content-dependent branch), NOTES_CAPTURE.md method. 3. combine consumer hunt (multiband cascade) for dual/comb. + +## ============ UPDATE 2026-08-22a: PHASE B RECOVERY — offline detector hypotheses REFUTED ============ + +Контекст: Phase B сессия 2026-08-21 (вечер) + утро 08-22 осталась незакоммиченной — +скрипты в /tmp/opencode, результаты не сохранялись. Все 4 эксперимента перепрогнаны, +stdout зафиксирован (`phaseB_*.out`), скрипты перенесены в `scripts/phaseB_*.py`. +Движок: phaseA_grid_fast.py (валидированный офлайн-тректор, погрешность 0.09..0.19 dB +на анкорах). Анкоры res_500/al_12/al_24/t1k_1000, критерий pred_err ~ 0 на ВСЕХ. + +| # | Гипотеза | Скрипт | Результат | Статус | +|---|----------|--------|-----------|--------| +| 1b | slide_max pooling lvl (w=3..33) | phaseB_pool2.py | w=3: tones −2.7/−2.3/−4.9, rms 3.02 (base 1.04) | **REFUTED** | +| 1c | neighborhood mean/RMS pooling | phaseB_pool3.py | best mean w=3 rms 1.35; res +1.9..+5.1 vs t1k −1.9..−4.4 tradeoff, ни одна точка не закрывает все 4 | **REFUTED** | +| 2 | temporal dynamics (hold b/dbdecay r/ema a на полной траектории) | phaseB_temporal.py | ВСЕ варианты = baseline (rms 1.03–1.05): метрика в steady-state, состояние успевает устояться до окна | **REFUTED (no-op)** | +| 3 | ρ(IIR1)+Δ joint scan (ρ∈[0.30,0.95], Δ∈[−1.5,+1.5]) | phaseB_rho.py | best ρ=0.830 D=+1.40 → rms 0.654; НО al_24 стабильно −1.19..−1.25, ρ без источника в декомпе (канон 0.692) | **REJECTED** (нарушает golden rule #1, тот же scalar-family тупик) | + +Выводы: +1. Пространственный pooling ЛЮБОГО вида (max/mean/rms) не объясняет content-gap. +2. Temporal-класс гипотез НЕПРОВЕРЯЕМ на steady-state анкорах — нужен переходный + контент (burst-рефы уже есть в soothe-bt) или другой анкорный набор. +3. Лучший (ρ,D) = репараметризация affine-семейства Phase A → упирается в тот же + KEY NEGATIVE RESULT (scalar-семейство не закрывает тон+шум одновременно). + Канон НЕ сменён: HEAD (LUT γ=0.344/MULT=4.2, ρ=0.692). + +NEXT (приоритеты без изменений): +1. Step 7 live capture ctx+0x188 (NOTES_CAPTURE.md метод) — константы A/B/γ + + возможная контент-зависимая ветка. Требует REAPER+плагин. +2. combine consumer hunt (межполосный каскад) для dual/comb. +3. Опционально: transient-анкоры для проверяемости temporal-класса. diff --git a/scripts/phaseA_grid_fast.py b/scripts/phaseA_grid_fast.py new file mode 100644 index 0000000..37d4efc --- /dev/null +++ b/scripts/phaseA_grid_fast.py @@ -0,0 +1,115 @@ +#!/usr/bin/env python3 +"""Phase A step 3 (fast): combo-vectorized law grid-search.""" +import re, sys +import numpy as np + +SRC = '/home/m/re-tools/dsp/rt_mask_tables.cpp' +src = open(SRC).read() +def tab(name): + m = re.search(r'const double %s\[\] = \{(.*?)\};' % name, src, re.S) + return np.array([float(x) for x in re.findall(r'[-+0-9.eE]+', m.group(1))]) +A1, B1 = tab('kIIR_A1'), tab('kIIR_B1') +A2, B2 = tab('kIIR_A2'), tab('kIIR_B2') +A3, B3 = tab('kIIR_A3'), tab('kIIR_B3') +DUMPDIR = {'dump_res_new.bin': '/tmp/', 'dump_t1k.bin': '/tmp/'} + +def load_dump(p): + d = np.loadtxt(p, skiprows=1); return d[:, 2], d[:, 6] +def load_traj(p): + b = open(p, 'rb').read(); off = 0; fr = [] + while off < len(b): + _, nb = np.frombuffer(b, dtype=np.int32, count=2, offset=off); off += 8 + fr.append(np.frombuffer(b, dtype='= lv.max()-6]; keep_n = (keep_dB.shape[0],) + else: + keep_dB = dB; keep_n = None + DATA[name] = (np.ascontiguousarray(keep_dB, dtype=np.float64), W.astype(np.float64), bm, keep_n) + +def gains_batch(dB, W, bm, X0s, SLs, CMs, C_pre=None, FLs=None): + if FLs is None: FLs = np.zeros_like(X0s) + """dB [T,nbin]; returns G [C,T] gain at bm for each combo.""" + import sys + print('gains_batch shapes:', dB.shape, W.shape, bm, X0s.shape, file=sys.stderr) + C, T, N = len(X0s), dB.shape[0], dB.shape[1] + if C_pre is not None: + c = np.broadcast_to(C_pre, (C, T, N)) + else: + X0 = X0s[:, None, None]; SL = SLs[:, None, None]; CM = CMs[:, None, None] + FL = FLs[:, None, None] + c = np.clip(X0 + SL*dB, FL, CM) # [C,T,N] + acc = np.zeros((C, T)) + y = np.empty_like(c) + for i in range(N): + acc = A1[i]*acc + B1[i]*c[:, :, i] + y[:, :, i] = acc + acc = np.zeros((C, T)) + for i in range(N): + acc = A2[i]*acc + B2[i]*y[:, :, i] + y[:, :, i] = 0.8*np.exp2(-acc)*W[i] + # IIR3 bidi x2 on y + for _ in range(2): + st = np.zeros((C, T)) + for i in range(N): + st = y[:, :, i]*B3[i] + st*A3[i] + y[:, :, i] = st + st = y[:, :, -1].copy() + for i in range(N-2, 0, -1): + st = y[:, :, i]*B3[i] + st*A3[i] + y[:, :, i] = st + return y[:, :, bm] + +# old-law reference gains +GO = {} +for name,(dB,W,bm,_) in DATA.items(): + c = np.clip((dB+13.78)/82.07, 0, 1)**0.344*4.2 + c_old = np.clip((dB+13.78)/82.07, 0, 1)**0.344*4.2 + g = gains_batch(dB, W, bm, np.array([0.]), np.array([0.]), np.array([99.]), C_pre=c_old) + GO[name] = float(g[0].mean()) if name=='res_500' else float(np.median(g[0])) + +def evaluate(X0s, SLs, CMs, FLs=None): + out = {} + for name,(dB,W,bm,_) in DATA.items(): + g = gains_batch(dB, W, bm, X0s, SLs, CMs, FLs=FLs) # [C,T] + agg = np.sqrt(np.mean(g**2, axis=1)) if name=='res_500' else np.median(g, axis=1) + out[name] = MEAS_OLD[name] + 20*np.log10(agg/GO[name]) + return out + +X0g = np.arange(1.85, 2.35, 0.05); SLg = np.arange(0.065, 0.102, 0.0025); CMg = np.array([99.]) +FLg = np.array([0., 0.15, 0.3, 0.45, 0.6]) +X0f, SLf, CMf, FLf = [j.ravel() for j in np.meshgrid(X0g, SLg, CMg, FLg, indexing='ij')] +names = list(DATA) + +recs = [] +CH = 120 +for s in range(0, len(X0f), CH): + sl = slice(s, s+CH) + ev = evaluate(X0f[sl], SLf[sl], CMf[sl], FLf[sl]) + for j in range(len(X0f[sl])): + e = {n: ev[n][j] for n in names} + recs.append((sum(v*v for v in e.values())/len(names), + X0f[sl][j], SLf[sl][j], CMf[sl][j], e, FLf[sl][j])) +recs.sort(key=lambda r: r[0]) +print('refined top-12:') +for tot, X0, SL, CM, e, FL in recs[:12]: + print(f' X0={X0:.2f} S={SL:.4f} FL={FL:.2f} rms={np.sqrt(tot):.3f} ' + + ' '.join(f'{n[:5]}:{v:+.2f}' for n,v in e.items())) +MEAS_NEW={'res_500':2.019,'al_12':0.244,'al_24':-0.069,'t1k_1000':-0.321} +ev18=evaluate(np.array([1.8]),np.array([0.11]),np.array([99.])) +print('new(1.8,.11) model-pred vs measured:') +for n in names: + print(f' {n:10} pred{ev18[n][0]:+.3f} meas{MEAS_NEW[n]:+.3f}') diff --git a/scripts/phaseB_pool2.py b/scripts/phaseB_pool2.py new file mode 100644 index 0000000..3e4cc28 --- /dev/null +++ b/scripts/phaseB_pool2.py @@ -0,0 +1,56 @@ +#!/usr/bin/env python3 +"""Phase B step 1b: spatial max-pooling scan on TOP of validated grid engine. + +Only deviation from phaseA_grid_fast.py: trajectory transform before law. +pool_lvl w: sliding max over bins (width w). pool_db == pool_lvl (monotone), +pool_am ~ pool_lvl near flat res -> skip both. +""" +import numpy as np + +exec(open('/tmp/opencode/phaseA_grid_fast.py').read().split("# old-law reference gains")[0]) + +def slide_max(x, w): + if w <= 1: return x + h = w // 2 + xp = np.pad(x, ((0, 0), (h, h)), mode='edge') + win = np.lib.stride_tricks.sliding_window_view(xp, w, axis=1) + return np.ascontiguousarray(win.max(axis=-1)) + +RAW = {} +for name, traj, dump, bm in [ + ('res_500','traj_res500.bin','dump_res_new.bin',85), + ('al_12','traj_al12.bin','dump_t1k.bin',85), + ('al_24','traj_al24.bin','dump_t1k.bin',85), + ('t1k_1000','traj_t1k.bin','dump_t1k.bin',85)]: + RAW[name] = load_traj('/tmp/opencode/'+traj)[WIN[name][0]:WIN[name][1]] + +GO = {} +for name,(dB,W,bm,_) in DATA.items(): + c_old = np.clip((dB+13.78)/82.07, 0, 1)**0.344*4.2 + g = gains_batch(dB, W, bm, np.array([0.]), np.array([0.]), np.array([99.]), C_pre=c_old) + GO[name] = float(g[0].mean()) if name=='res_500' else float(np.median(g[0])) + +def eval_variant(w): + global DATA + saved = {n: DATA[n] for n in DATA} + for name,(dB,W,bm,_) in DATA.items(): + T = slide_max(RAW[name], w) + dBp = np.log10(np.maximum(T, 1e-12))*20.0 + if name != 'res_500': + lv = dBp[:, bm]; keep = dBp[lv >= lv.max()-6] + else: + keep = dBp + DATA[name] = (np.ascontiguousarray(keep), W, bm, None) + out = {} + for name,(dB,W,bm,_) in DATA.items(): + g = gains_batch(dB, W, bm, np.array([1.8]), np.array([0.11]), np.array([99.])) + agg = np.sqrt(np.mean(g**2, axis=1)) if name=='res_500' else np.median(g, axis=1) + out[name] = MEAS_OLD[name] + 20*np.log10(float(agg[0])/GO[name]) + DATA.update(saved) + return out + +print(f'{"w":>3} ' + ' '.join(f'{n:>9}' for n in DATA) + ' rms') +for w in [1, 3, 5, 9, 17, 33]: + e = eval_variant(w) + tot = np.sqrt(sum(v*v for v in e.values())/len(e)) + print(f'{w:>3} ' + ' '.join(f'{v:+9.2f}' for v in e.values()) + f' {tot:.2f}') diff --git a/scripts/phaseB_pool3.py b/scripts/phaseB_pool3.py new file mode 100644 index 0000000..ceb7033 --- /dev/null +++ b/scripts/phaseB_pool3.py @@ -0,0 +1,55 @@ +#!/usr/bin/env python3 +"""Phase B step 1c: neighborhood MEAN/RMS pooling scan (max already refuted). + +Rationale: slide_max raises lvl at tone bins via sidelobes -> over-reduction +(tones broke -2.7..-4.9). Mean/RMS pooling does the opposite for an isolated +narrow peak among quiet neighbours -> less reduction on tones, ~neutral on +wide noise. Pooling in LINEAR lvl domain (am ~ lvl near flat res). +""" +import numpy as np + +exec(open('/tmp/opencode/phaseA_grid_fast.py').read().split("# old-law reference gains")[0]) + +def slide(x, w, op): + if w <= 1: return x + h = w // 2 + xp = np.pad(x, ((0, 0), (h, h)), mode='edge') + win = np.lib.stride_tricks.sliding_window_view(xp, w, axis=1) + if op == 'mean': return win.mean(axis=-1) + return np.sqrt((win ** 2).mean(axis=-1)) + +RAW = {} +for name, traj in [('res_500','traj_res500.bin'), ('al_12','traj_al12.bin'), + ('al_24','traj_al24.bin'), ('t1k_1000','traj_t1k.bin')]: + RAW[name] = load_traj('/tmp/opencode/'+traj)[WIN[name][0]:WIN[name][1]] + +GO = {} +for name,(dB,W,bm,_) in DATA.items(): + c_old = np.clip((dB+13.78)/82.07, 0, 1)**0.344*4.2 + g = gains_batch(dB, W, bm, np.array([0.]), np.array([0.]), np.array([99.]), C_pre=c_old) + GO[name] = float(g[0].mean()) if name=='res_500' else float(np.median(g[0])) + +def eval_variant(op, w): + saved = {n: DATA[n] for n in DATA} + for name,(dB,W,bm,_) in DATA.items(): + T = slide(RAW[name], w, op) + dBp = np.log10(np.maximum(T, 1e-12))*20.0 + if name != 'res_500': + lv = dBp[:, bm]; keep = dBp[lv >= lv.max()-6] + else: + keep = dBp + DATA[name] = (np.ascontiguousarray(keep), W, bm, None) + out = {} + for name,(dB,W,bm,_) in DATA.items(): + g = gains_batch(dB, W, bm, np.array([1.8]), np.array([0.11]), np.array([99.])) + agg = np.sqrt(np.mean(g**2, axis=1)) if name=='res_500' else np.median(g, axis=1) + out[name] = MEAS_OLD[name] + 20*np.log10(float(agg[0])/GO[name]) + DATA.update(saved) + return out + +print(f'{"op":>4} {"w":>3} ' + ' '.join(f'{n:>9}' for n in DATA) + ' rms') +for op in ['mean', 'rms']: + for w in [3, 5, 9, 17, 33]: + e = eval_variant(op, w) + tot = np.sqrt(sum(v*v for v in e.values())/len(e)) + print(f'{op:>4} {w:>3} ' + ' '.join(f'{v:+9.2f}' for v in e.values()) + f' {tot:.2f}') diff --git a/scripts/phaseB_pooling.py b/scripts/phaseB_pooling.py new file mode 100644 index 0000000..947d7ff --- /dev/null +++ b/scripts/phaseB_pooling.py @@ -0,0 +1,110 @@ +#!/usr/bin/env python3 +"""Phase B step 1: spatial max-pooling hypothesis scan (offline, no rebuilds). + +lvl'(t,k) = pool(lvl)(t,k) with width w, then law -> chain -> median/rms gain +at metric bin. Variants: pool_am (pool raw am then /res), pool_lvl, pool_db. +""" +import re +import numpy as np + +SRC = '/home/m/re-tools/dsp/rt_mask_tables.cpp' +src = open(SRC).read() +def tab(name): + m = re.search(r'const double %s\[\] = \{(.*?)\};' % name, src, re.S) + return np.array([float(x) for x in re.findall(r'[-+0-9.eE]+', m.group(1))]) +A1, B1 = tab('kIIR_A1'), tab('kIIR_B1') +A2, B2 = tab('kIIR_A2'), tab('kIIR_B2') +A3, B3 = tab('kIIR_A3'), tab('kIIR_B3') +DUMPDIR = {'dump_res_new.bin': '/tmp/', 'dump_t1k.bin': '/tmp/'} + +def load_dump(p): + d = np.loadtxt(p, skiprows=1); return d[:, 1], d[:, 2], d[:, 6] # am,res,W +def load_traj(p): + b = open(p, 'rb').read(); off = 0; fr = [] + while off < len(b): + _, nb = np.frombuffer(b, dtype=np.int32, count=2, offset=off); off += 8 + fr.append(np.frombuffer(b, dtype='18} {"w":>3} ' + ' '.join(f'{n:>9}' for n in CASES) + ' (pred err, dB)') +# baselines on CORRECT lvl (traj stores lvl_raw already) +GBASE = {} +for name,(T,res_k,W,bm) in CASES.items(): + dB = np.log10(np.maximum(T, 1e-12)) + GBASE[name] = agg(gains_batch(OLD(dB)[None], W, bm)[0], name) +print('sanity new@w=1 (vs validated):') +row=[] +for name,(T,res_k,W,bm) in CASES.items(): + dB = np.log10(np.maximum(T, 1e-12)) + gn = agg(gains_batch(LAW['new'](dB)[None], W, bm)[0], name) + row.append(MEAS_OLD[name] + 20*np.log10(gn/GBASE[name])) +print(f'{"new":>18} {1:>3} ' + ' '.join(f'{v:+9.2f}' for v in row)) + +for variant in ['pool_lvl', 'pool_am', 'pool_db']: + for w in [3, 5, 9, 17]: + row = [] + for name,(T,res_k,W,bm) in CASES.items(): + am = T # stored lvl_raw = am/res*scale -> recover am = lvl/res*scale... careful + # stored lvl_raw = am/res_k * SCALE => am = lvl_raw * res_k / SCALE + am_abs = T * res_k[None,:] / SCALE + if variant == 'pool_am': + lv = slide_max(am_abs[None], w)[0] / res_k[None,:] * SCALE + elif variant == 'pool_lvl': + lv = slide_max(T[None], w)[0] + else: + db_ = np.log10(np.maximum(T, 1e-12))*20 + lv = 10**(slide_max(db_[None], w)[0]/20) + dB = np.log10(np.maximum(lv, 1e-12))*20 + g = gains_batch(LAW['new'](dB)[None], W, bm) + row.append(MEAS_OLD[name] + 20*np.log10(agg(g[0],name)/GBASE[name])) + print(f'{variant:>18} {w:>3} ' + ' '.join(f'{v:+9.2f}' for v in row)) diff --git a/scripts/phaseB_rho.py b/scripts/phaseB_rho.py new file mode 100644 index 0000000..ab1a3ef --- /dev/null +++ b/scripts/phaseB_rho.py @@ -0,0 +1,71 @@ +#!/usr/bin/env python3 +"""Phase B step 3: (rho, Delta) joint scan. + +Hypothesis: content gap lives in IIR1 spike attenuation vs law level. +rho = IIR1 pole (DC-normalized: y = rho*acc + (1-rho)*x), canon rho=0.692. +Delta = additive shift of affine law c = max(1.8+D+0.11*dB, 0). +Anchored at canon old-law gains GO. Criterion: pred_err ~ 0 on ALL 4 anchors. +""" +import numpy as np + +exec(open('/tmp/opencode/phaseA_grid_fast.py').read().split("# old-law reference gains")[0]) + +GO = {} +for name, (dB, W, bm, _) in DATA.items(): + c_old = np.clip((dB + 13.78) / 82.07, 0, 1) ** 0.344 * 4.2 + g = gains_batch(dB, W, bm, np.array([0.]), np.array([0.]), np.array([99.]), C_pre=c_old) + GO[name] = float(g[0].mean()) if name == 'res_500' else float(np.median(g[0])) +print('GO:', {k: round(v, 3) for k, v in GO.items()}) + +def gains_rho(dB, W, bm, rho, DLs): + """dB [T,N]; law c=max(1.8+Dl+0.11*dB,0); IIR1 pole=rho (DC-norm). + returns [C,T] gain at bm for each Delta in DLs.""" + C, T, N = len(DLs), dB.shape[0], dB.shape[1] + c = np.maximum(1.8 + DLs[:, None, None] + 0.11 * dB[None], 0.0) + acc = np.zeros((C, T)); y = np.empty_like(c) + for i in range(N): + acc = rho * acc + (1 - rho) * c[:, :, i] + y[:, :, i] = acc + acc = np.zeros((C, T)) + for i in range(N): + acc = A2[i] * acc + B2[i] * y[:, :, i] + y[:, :, i] = 0.8 * np.exp2(-acc) * W[i] + for _ in range(2): + st = np.zeros((C, T)) + for i in range(N): + st = y[:, :, i] * B3[i] + st * A3[i] + y[:, :, i] = st + st = y[:, :, -1].copy() + for i in range(N - 2, 0, -1): + st = y[:, :, i] * B3[i] + st * A3[i] + y[:, :, i] = st + return y[:, :, bm] + +RHOS = np.linspace(0.30, 0.95, 131) +DLS = np.linspace(-1.5, 1.5, 121) +names = list(DATA) +best = [] +for rho in RHOS: + ev = {} + for name, (dB, W, bm, _) in DATA.items(): + g = gains_rho(dB, W, bm, rho, DLS) + agg = np.sqrt(np.mean(g ** 2, axis=1)) if name == 'res_500' else np.median(g, axis=1) + ev[name] = MEAS_OLD[name] + 20 * np.log10(agg / GO[name]) + E = np.stack([ev[n] for n in names]) # [4, C] + rms = np.sqrt((E ** 2).mean(axis=0)) # per Delta + j = int(rms.argmin()) + best.append((rms[j], rho, DLS[j], E[:, j])) +best.sort() +print('\ntop-10 (rms over 4 anchors):') +for rms, rho, dl, e in best[:10]: + print(f' rho={rho:.3f} D={dl:+.3f} rms={rms:.3f} ' + + ' '.join(f'{n[:5]}:{v:+.2f}' for n, v in zip(names, e))) +print(f'\ncanon rho=0.692 D=0 reference:') +j0 = int(np.argmin(np.abs(DLS))) +for rho in [0.692]: + ev = {} + for name, (dB, W, bm, _) in DATA.items(): + g = gains_rho(dB, W, bm, rho, DLS[j0:j0+1]) + agg = np.sqrt(np.mean(g[0] ** 2)) if name == 'res_500' else np.median(g[0]) + ev[name] = MEAS_OLD[name] + 20 * np.log10(agg / GO[name]) + print(' ' + ' '.join(f'{n[:5]}:{v:+.2f}' for n, v in ev.items())) diff --git a/scripts/phaseB_temporal.py b/scripts/phaseB_temporal.py new file mode 100644 index 0000000..39150dc --- /dev/null +++ b/scripts/phaseB_temporal.py @@ -0,0 +1,68 @@ +#!/usr/bin/env python3 +"""Phase B step 2: temporal detector-dynamics scan (offline). + +Variants with time memory applied to FULL trajectory (state settles before +metric window), then validated pipeline (window, -6dB core, median/rms). +""" +import numpy as np + +exec(open('/tmp/opencode/phaseA_grid_fast.py').read().split("# old-law reference gains")[0]) + +FULL = {} +for name, traj, dump, bm in [ + ('res_500','traj_res500.bin','dump_res_new.bin',85), + ('al_12','traj_al12.bin','dump_t1k.bin',85), + ('al_24','traj_al24.bin','dump_t1k.bin',85), + ('t1k_1000','traj_t1k.bin','dump_t1k.bin',85)]: + FULL[name] = load_traj('/tmp/opencode/'+traj) + +GO = {} +for name,(dB,W,bm,_) in DATA.items(): + c_old = np.clip((dB+13.78)/82.07, 0, 1)**0.344*4.2 + g = gains_batch(dB, W, bm, np.array([0.]), np.array([0.]), np.array([99.]), C_pre=c_old) + GO[name] = float(g[0].mean()) if name=='res_500' else float(np.median(g[0])) + +def t_hold(T, b): # linear peak-hold decay + out = T.copy() + for t in range(1, len(T)): + out[t] = np.maximum(T[t], b*out[t-1]) + return out + +def t_dbdecay(T, r): # dB-domain peak decay r dB/frame + db = np.log10(np.maximum(T, 1e-12))*20.0 + for t in range(1, len(db)): + db[t] = np.maximum(db[t], db[t-1]-r) + return 10**(db/20) + +def t_ema(T, a): # EMA in dB domain + db = np.log10(np.maximum(T, 1e-12))*20.0 + out = db.copy() + for t in range(1, len(db)): + out[t] = a*db[t] + (1-a)*out[t-1] + return 10**(out/20) + +def eval_tf(fn): + out = {} + for name,(dB,W,bm,_) in DATA.items(): + Tt = fn(FULL[name])[WIN[name][0]:WIN[name][1]] + dBp = np.log10(np.maximum(Tt, 1e-12))*20.0 + if name != 'res_500': + lv = dBp[:, bm]; keep = dBp[lv >= lv.max()-6] + else: + keep = dBp + g = gains_batch(np.ascontiguousarray(keep), W, bm, + np.array([1.8]), np.array([0.11]), np.array([99.])) + agg = np.sqrt(np.mean(g**2, axis=1)) if name=='res_500' else np.median(g, axis=1) + out[name] = MEAS_OLD[name] + 20*np.log10(float(agg[0])/GO[name]) + return out + +VARS = [('none', lambda T: T)] +for b in [0.8, 0.9, 0.95, 0.99]: VARS.append((f'hold b={b}', lambda T, b=b: t_hold(T, b))) +for r in [0.25, 0.5, 1.0, 2.0]: VARS.append((f'dbdec r={r}', lambda T, r=r: t_dbdecay(T, r))) +for a in [0.3, 0.5, 0.7]: VARS.append((f'ema a={a}', lambda T, a=a: t_ema(T, a))) + +print(f'{"variant":>12} ' + ' '.join(f'{n:>9}' for n in DATA) + ' rms') +for lbl, fn in VARS: + e = eval_tf(fn) + tot = np.sqrt(sum(v*v for v in e.values())/len(e)) + print(f'{lbl:>12} ' + ' '.join(f'{v:+9.2f}' for v in e.values()) + f' {tot:.2f}')