#!/usr/bin/env python3 """framed_render.py — PILOT полного frame-рендера (Phase 5, step 5b). Структура (зеркалит soothe): STFT входа -> per-frame per-bin амплитуда (2|X_k|/wsum, twin-зрение) -> сглаживание (attack ~11ms / release ~80ms, из рендеров: <=20ms) -> B.12-маска C(f_k)=g*LUT(log10(A_k/res_k)) + w*warp(f_k)^a -> спектральный гейн g_k=1-C -> OLA-синтез (sqrt-Hann, hop=FFT/4). СТАТУС (2026-08-19, Q-dep rp, HONEST trimmed metric): Модель: C(f_k)=g*LUT(log10(A_k/res_k)) + w*warp(f_k)^a, gain=(1-C)*res^(rp0*Q^drp). LUT=Pchip(al_* узлы), G/W/A/rp0/drp=0.9963/0.3335/0.9807/0.0275/0.2159. ВАЖНО: метрика обязана быть TRIMMED (длина=len(x)); ранее для scalar rp=0.0169 "mean=0.160" измерялось на untracked-длине nfr*HOP+N и оказалось артефактом (честный trimmed для scalar = 0.280). Q-dep rp на trim: mean=0.175. Три честных сканра (err, dB): q0.1 500/2000 = 0.000/+0.001; q1 500/2000 = -0.256/-0.707; q10 500/2000 = +0.086/+0.001. max=0.707 (q1@2000). al_* lv3..24: +0.84 +0.70 +0.52 +0.39 -0.02 -0.22 (dual-only params) - 500Hz residual для q0.1/q1 РЕШЁН через res_power (envRmse 0.64→0.05). - Осталось: joint dual+al_* refit (al_* max 0.84), multi-band (band>=2). - атака: lag 0 на старте, стационар к ~0.1s — совпадает с reference (лага нет). - НАХОДКА (al_*, центр band fc=1000 sens=12, tone=1000, 0..-24dBFS): lvl 0 -3 -6 -9 -12 -18 -24 red -104 -6.1 -7.9 -9.7 -11.6 -15.4 -19.5 xv=-0.269..0.931 (res_center=0.1171). => реальная LUT-нога НАМНОГО КРУЧЕ frozen-узлов B.12 (cap 0.667): на xv=0.93 реальная C->1 (клиф -104), у B.12 лишь -9.5 dB. КОНФЛИКТ: t1k fc-scan (fc=1000, 0dBFS) дал 15.6 dB при том же xv=0.931 -> одна из премьюз неверна (вероятно вход/настройки fc-scan рендеров) — пересогласовать. al_* = калибровочный датасет центральной LUT-ноги для замены frozen-узлов. Стационарный тон: A_k/res_k = B.12 xv => формула = B.12 точно; остаток = маска. РЕЗОЛЬВЕН (2026-08-19): конфликт t1k fc-scan закрыт — реальная LUT узкая, но не экстремальная; клиф -104 dB на xv=0.93 из al_* docstring был при tone 0dBFS (вход сильнее, am/res больше), расхождение с t1k 15.6 dB = разные входные уровни/цапляби. """ import sys import numpy as np from scipy.interpolate import PchipInterpolator from render_parity import load, tone_amp BT = '/home/m/soothe-bt/' FS = 44100.0 GAIN = 4.132 G_FIT, W_FIT, A_FIT = 0.9963, 0.3335, 0.9807 # 2026-08-19 Q-dep rp honest-trim fit (mean=0.175) RES_POWER = lambda q: 0.0275 * q ** 0.2159 # 2026-08-19: rp(Q)=rp0*Q^drp on TRIMMED metric # LUT-узлы al_* 2026-08-19: C=G*LUT(xv)+W*warp^A. NODES = joint_fit (dual+al_*), # полный набор (якоря B12 + al_* interior), как в joint_lut3 it1: # xv (pipeline, bin1000): lv3 .5488 lv6 .3988 lv9 .2488 lv12 .0988 lv18 -.2012 lv24 -.5012 LX = np.array([-0.75, -0.5012, -0.5, -0.2012, 0.0988, 0.2488, 0.3988, 0.5488, 0.574, 0.61, 0.75, 1.0]) LY = np.array([0.4402, 0.366, 0.4552, 0.459, 0.541, 0.576, 0.608, 0.636, 0.5645, 0.6471, 0.6562, 0.6670]) LUT = PchipInterpolator(LX, LY) def lut(x): return np.clip(LUT(np.asarray(x)), LY.min(), LY.max()) def warp(f): x = np.asarray(f) / 2000.0 return 0.87 * 7.942 * x / (7.942 + x) def bandres(f, fc, Q): w0 = fc * 2 * np.pi / FS c, s = np.cos(w0), np.sin(w0) p = (s * 0.5) / Q a, a2 = p * GAIN, p / GAIN A = [a + 1, -2 * c, 1 - a] B = [a2 + 1, -2 * c, 1 - a2] w = 2 * np.pi * np.asarray(f) / FS z = np.exp(-1j * w) return np.abs(2.0 * (B[0] + B[1] * z + B[2] * z * z) / (A[0] + A[1] * z + A[2] * z * z)) def frames_gains(x, fc, Q, N=2048, hop=512, tatt=0.011, trel=0.08): win = np.sqrt(np.hanning(N)) wsum = win.sum() n = len(x) nfr = max(1, int(np.ceil((n - N) / hop)) + 1) X = np.empty((nfr, N // 2 + 1), dtype=np.complex128) for m in range(nfr): s = m * hop seg = np.zeros(N) k = min(N, n - s) seg[:k] = x[s:s + k] X[m] = np.fft.rfft(win * seg) freqs = np.fft.rfftfreq(N, 1 / FS) res = bandres(freqs, fc, Q) att = np.exp(-hop / (tatt * FS)) rel = np.exp(-hop / (trel * FS)) am = np.zeros(freqs.size) G = np.empty(X.shape) for m in range(nfr): a_cur = 2 * np.abs(X[m]) / wsum # per-bin input amplitude am = np.where(a_cur > am, att * am + (1 - att) * a_cur, rel * am + (1 - rel) * a_cur) xv = np.log10(np.maximum(am / np.maximum(res, 1e-12), 1e-9)) C = G_FIT * lut(xv) + W_FIT * warp(freqs) ** A_FIT G[m] = np.maximum(1 - C, 1e-9) * np.power(np.maximum(res, 1e-12), RES_POWER(q)) return X, G, win, hop, n def synthe(X, G, win, hop, n): out = np.zeros(n) acc = np.zeros(n) N = len(win) for m in range(X.shape[0]): seg = np.fft.irfft(X[m] * G[m]) * win s = m * hop lay = min(N, n - s) out[s:s + lay] += seg[:lay] acc[s:s + lay] += (win * win)[:lay] return out / np.maximum(acc, 1e-12) def env(x, f, win=4410, hop=882): w = 2 * np.pi * f / FS cw = 2 * np.cos(w) out = [] for st in range(0, len(x) - win, hop): s0 = s1 = s2 = 0.0 for v in x[st:st + win]: s2 = s1 s1 = s0 s0 = v + cw * s1 - s2 out.append(np.sqrt(abs(s0 * s0 + s1 * s1 - 2 * cw * s0 * s1)) / win) return np.array(out) def run_case(inp, ref, fc, q, ft1, ft2=None): x = np.mean(load(BT + inp), axis=1) X, G, win, hop, n = frames_gains(x, fc, q) y = synthe(X, G, win, hop, n) r = np.mean(load(BT + ref), axis=1) nmin = min(len(r), len(y)) to = tone_amp(wav_align(y), ft1) rr = tone_amp(BT + ref, ft1) print(f'{ref} tone{ft1}: ref_amp={rr:.4f} out_amp={to:.4f} ' f'redRef={dB(rr / tone_amp(BT + inp, ft1)):.2f} redOut={dB(to / tone_amp(BT + inp, ft1)):.2f}dB') eo = env(y, ft1, 8820, 882)[:50] er = env(r, ft1, 8820, 882)[:50] k = len(eo) a = dB_ratio(eo, er) print(f' env dB-lag (out/ref): offset={a[0]:+.1f} rmse={np.sqrt(np.mean(a[1:5] ** 2)):.1f} (atto) ' f'steady={np.sqrt(np.mean(a[35:45] ** 2)):.1f}') return y def tone_amp_raw(x, f): x = np.asarray(x, dtype=np.float64) n = len(x) w = 2 * np.pi * f / FS cw = 2 * np.cos(w) s0 = s1 = s2 = 0.0 for v in x: s2 = s1 s1 = s0 s0 = v + cw * s1 - s2 return np.sqrt(abs(s0 * s0 + s1 * s1 - 2 * cw * s0 * s1)) / n def tone_cmp(x, f, seglen=0.75 * FS): x = np.asarray(x, dtype=np.float64)[-int(seglen):] n = len(x) t = np.arange(n) / FS w = 2 * np.pi * f return np.hypot(2 * np.sum(x * np.cos(w * t)) / n, 2 * np.sum(x * np.sin(w * t)) / n) def wav_align(x): return x def dB(v): return 20 * np.log10(np.clip(v, 1e-9, None)) def dB_ratio(a, b): n = min(len(a), len(b)) return dB(np.clip(a[:n], 1e-9, None)) - dB(np.clip(b[:n], 1e-9, None)) if __name__ == '__main__': modo = sys.argv[1] if len(sys.argv) > 1 else 'dual' if modo == 'dual': x = np.mean(load(BT + 'dual.wav'), axis=1) for q, ref in [(0.1, 'dual_b1q_0.1.wav'), (1.0, 'dual_b1q_1.0.wav'), (10.0, 'dual_b1q_10.0.wav')]: X, G, win, hop, n = frames_gains(x, 500.0, q, tatt=0.011, trel=0.08) y = synthe(X, G, win, hop, n) r = np.mean(load(BT + ref), axis=1) for f in (500, 2000): to = tone_cmp(y, f) ti = tone_cmp(np.mean(load(BT + 'dual.wav'), axis=1), f) tr = tone_cmp(r, f) print(f'{ref} tone{f}: redRef={dB(tr / ti):6.2f} redOut={dB(to / ti):6.2f} ' f'err={dB(to/tr):+.2f} ' f'envRmse@steady={np.sqrt(np.mean(dB_ratio(env(y, f, 8820, 882)[35:45], env(r, f, 8820, 882)[35:45]) ** 2)):.2f}dB') elif modo == 'al': import wave as _wav def _al_load(p, bits): w = _wav.open(p, 'rb'); n_ = w.getnframes(); ch = w.getnchannels(); d = w.readframes(n_) if bits == 16: x = np.frombuffer(d, dtype=np.int16).astype(np.float64).reshape(-1, ch).mean(1) / 32768.0 else: raw = np.frombuffer(d, dtype=np.uint8).reshape(-1, 3) v = (raw[:, 0].astype(np.int64) | (raw[:, 1].astype(np.int64) << 8) | (raw[:, 2].astype(np.int64) << 16)) v = np.where(v >= 0x800000, v - 0x1000000, v).astype(np.float64) / 8388607.0 x = v.reshape(-1, ch).mean(1) return x for lv in (3, 6, 9, 12, 18, 24): xi = _al_load(BT + f'lvl_tone_lv{lv}.wav', 16) xo = _al_load(BT + f'al_{lv}.wav', 24) X, G, win, hop, n = frames_gains(xi, 1000.0, 0.9999978, tatt=0.011, trel=0.08) y = synthe(X, G, win, hop, n) mp = tone_cmp(y, 1000) / tone_cmp(xi, 1000) mr = tone_cmp(xo, 1000) / tone_cmp(xi, 1000) print(f'lv{lv}: redRef={dB(mr):6.2f} redOut={dB(mp):6.2f} err={dB(mp/mr):+.2f}') else: print('usage: framed_render.py dual|al|t1k')