diff --git a/model_fir.py b/model_fir.py
new file mode 100644
index 0000000..469f0ec
--- /dev/null
+++ b/model_fir.py
@@ -0,0 +1,128 @@
+#!/usr/bin/env python3
+"""model_fir.py — РЕЗУЛЬТАТ диагностического моста (2026-08-18).
+
+Вопрос: даёт ли РЕАЛЬНАЯ цепочка FUN_180529fe0 (masкa -> warp -> FIR/DC) глубину
+dual_b1q с реDD-константами из декомпа? Ответ: ДА, через АДДИТИВНЫЙ аккумулятор.
+
+МОДЕЛЬ (rmse=0.236 dB на 36 точках):
+  C(f) = g·LUT(log10(L0/res_band(f)))  +  w·warp(f)^a
+  red  = -20·log10(1-C)
+  g=1.221  w=0.358  a=3.143   (фит, bounds [0.9..1.5],[0.2..1.2],[2..6])
+  LUT   = PCHIP B.11 (заморожен), warp = 0.87·Kx/(K+x), K=exp(2.0723), x=f/2000
+
+КЛЮЧЕВЫЕ ВЫВОДЫ:
+1. dual500-константа объяснена БЕЗ tilt: res_band(500; fc=500)=0.117 не зависит от Q
+   (твин в центре полосы Q-независим); C500 = g·LUT(0.574)·warp^a(500≈0) = 0.692 -> 10.23 dB.
+2. Глубина dual2000 НЕ мультипликативный warp·LUT (ratio warp 2000/500=3.66 -> rmse>10 dB);
+   её даёт добавка w·warp(f)^3.14: на 2000 Гц вклад 0.358·0.773^3.14=0.149, на 1000<=0.02
+   (иначе t1k/t1kq рушатся) => точка фита фиксирована единственной. Экзотика a≈π: кандидат —
+   кратный каскад warp (freq-axis 0x540698 / ∏0x540688 / двойной FFT-проход).
+3. Остаток dual2000 при Q=0.1 (14.4 vs 15.2) = тонкая форма LUT-колена 0.574 (0.7 dB) —
+   не закрывается per-bin формами; требует точного ур-я импульса (FFT-уровень 0x535a70).
+4. Порядок измерений в фите: DUAL.ravel() (500,2000 попарно) | T1KQ | T1K.
+5. res(2000) растёт 0.216(Q=0.1)..1.988(Q=10) и входит ЧЕРЕЗ LUT(xv), где xv=log10(L0/res):
+   xv(2000)=0.308(Q=0.1)..-0.656(Q=10) -> C2000-Q-зависимость совпадает.
+"""
+import numpy as np
+from scipy.optimize import least_squares
+from scipy.interpolate import PchipInterpolator
+
+FS = 44100.0
+GAIN = 4.132
+
+QS = [0.1, 0.2, 0.3, 0.5, 0.7, 1.0, 1.5, 2.0, 3.0, 5.0, 10.0]
+DUAL = np.array([(10.220, 15.224), (10.219, 13.963), (10.219, 12.947),
+                 (10.219, 11.725), (10.218, 11.119), (10.216, 10.689),
+                 (10.210, 10.412), (10.203, 10.305), (10.182, 10.225),
+                 (10.113, 10.183), (9.822, 10.165)])
+FCS = [800.0, 900.0, 950.0, 1000.0, 1050.0, 1100.0, 1200.0]
+T1KQ = np.array([7.868, 8.536, 8.726, 8.788, 8.729, 8.575, 8.115])
+T1K = np.array([14.548, 15.332, 15.553, 15.626, 15.557, 15.378, 14.840])
+L0_DUAL = 10 ** (-7.142 / 20)
+L0_T1KQ = 10 ** (-18.063 / 20)
+L0_T1K = 1.0
+
+LUT_X = np.array([-0.750, -0.500, -0.250, 0.000, 0.250, 0.500, 0.574, 0.610, 0.750, 1.000])
+LUT_Y = np.array([0.4402, 0.4552, 0.4813, 0.5072, 0.5329, 0.5332, 0.5645, 0.6471, 0.6562, 0.6670])
+_LUT = PchipInterpolator(LUT_X, LUT_Y)
+
+
+def lut(x):
+    return np.float64(np.clip(_LUT(float(x)), LUT_Y[0], LUT_Y[-1]))
+
+
+def res_band(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 * 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 warp(f):
+    x = f / 2000.0
+    return 0.87 * 7.942 * x / (7.942 + x)
+
+
+def make_combos():
+    combos = []
+    for q in QS:
+        combos += [(500.0, 500.0, q, L0_DUAL), (2000.0, 500.0, q, L0_DUAL)]
+    for fc in FCS:
+        combos.append((1000.0, fc, 0.9999978, L0_T1KQ))
+    for fc in FCS:
+        combos.append((1000.0, fc, 0.9999978, L0_T1K))
+    return combos
+
+
+COMBOS = make_combos()
+
+
+def pred(P):
+    g, w, a = P
+    out = []
+    for ft, fc, Q, L0 in COMBOS:
+        xv = np.log10(L0 / res_band(ft, fc, Q))
+        m = g * lut(xv) + w * warp(ft) ** a
+        out.append(-20 * np.log10(1 - np.clip(m, 0, 0.999)))
+    return np.array(out)
+
+
+def run():
+    meas = np.concatenate([DUAL.ravel(), T1KQ, T1K])
+    rng = np.random.default_rng(2)
+    best = None
+    for trial in range(40):
+        try:
+            p0 = rng.uniform([0.9, 0.2, 2.0], [1.5, 1.2, 6.0])
+            r = least_squares(lambda p: pred(p) - meas, p0,
+                              bounds=([0.9, 0.2, 2.0], [1.5, 1.2, 6.0]),
+                              max_nfev=400, ftol=1e-11, x_scale='jac')
+            pr = pred(r.x)
+            rm = np.sqrt(np.mean((pr - meas) ** 2))
+            if best is None or rm < best[0]:
+                best = (rm, pr, r.x)
+        except Exception:
+            pass
+    rm, pr, P = best
+    seg = lambda i0, i1: np.sqrt(np.mean((pr[i0:i1] - meas[i0:i1]) ** 2))
+    print(f'rmse={rm:.4f} dual={seg(0,22):.3f} t1kq={seg(22,29):.3f} t1k={seg(29,36):.3f} '
+          f'P(g,w,a)={np.round(P,4)}  [warp^a: 2000={warp(2000)**P[2]:.3f}, 1000={warp(1000)**P[2]:.3f}]')
+    print('d500=' + ' '.join(f'{pr[2*i]:5.2f}/{meas[2*i]:5.2f}' for i in range(11)))
+    print('d2k =' + ' '.join(f'{pr[2*i+1]:5.2f}/{meas[2*i+1]:5.2f}' for i in range(11)))
+    print('t1kq=' + ' '.join(f'{pr[22+j]:5.2f}/{meas[22+j]:5.2f}' for j in range(7)))
+    print('t1k =' + ' '.join(f'{pr[29+j]:5.2f}/{meas[29+j]:5.2f}' for j in range(7)))
+    # числовая справка res/warp
+    print('\nres(500;fc500)=%.3f (Q-независим, объясняет C500-константу)' % res_band(500, 500, 1.0))
+    print('res(2000): ' + ' '.join('%.3f' % res_band(2000, 500, q) for q in QS))
+    print('xv(2000)=log10(L0D/res): ' +
+          ' '.join('%.3f' % np.log10(L0_DUAL / res_band(2000, 500, q)) for q in QS))
+
+
+if __name__ == '__main__':
+    run()
\ No newline at end of file