Matching: split the envelope match into a static equalizer and a remainder
2 changed files+69−0
matching/README.mdmodified+2−0View file
@@ -7,6 +7,8 @@ This directory holds scripts for imposing the frequency-dependent decay of a rec
77 - `analyze.py WAV ...`: per-band levels, fitted decay times and the strongest partials of one or more clips.
88 - `match.py REFERENCE.wav SIMULATED.{wav,f32} OUT`: aligns both sounds at their onsets, measures band envelopes on a third-octave grid with a 4096-sample STFT, and writes `OUT_decay.wav` (the simulation with the recording's per-band decay rates, its own initial balance kept), `OUT_envelope.wav` (the recording's whole band envelope, balance and decay, imposed on the simulation), and the aligned inputs `OUT_reference_aligned.wav` and `OUT_sim_aligned.wav`. All outputs are level-matched on the body of the sound (RMS over 0.1 to 1 s) and the attack is soft-clipped.
99 - `plot.py OUT`: spectrograms of the four clips and band decay curves of the recording against the simulation.
10+- `eq.py REFERENCE.wav SIMULATED.f32 OUT`: splits the envelope-mode gain table into a static per-band equalizer (the balance at 0.1 s) and a time-varying remainder, reports both, and writes the simulation through the static equalizer (`OUT_eq.wav`, 28 band gains) and through a cubic-in-log-frequency approximation of it (`OUT_eq4.wav`, 4 coefficients).
11+- `partials.py CLIP ...`: the dominant partials of several clips side by side.
1012 - `modes.ts OUT.bin [key=value ...]` and `modes.py OUT.bin OUT.json [fmax]`: linearized mode frequencies of the model sheet (all three displacement fields, including the stiffening by the initial bow), for comparing the model's partials with a recording without running the simulator. Any parameter can be overridden, for example `bowAmplitude=0` or `E=178e9`.
1113
1214 `match.py` takes `--bands N` to set the number of bands per octave (default 3). Twelfth-octave bands are narrow enough to weight individual partials, which makes the result closer to the recording but says less about the loss model.
matching/eq.pyadded+67−0View file
@@ -0,0 +1,67 @@
1+"""How much of the envelope match is a static equalizer?
2+
3+Usage: python3 eq.py REFERENCE.wav SIMULATED.f32 OUT_PREFIX
4+
5+Recomputes the envelope-mode gain table of match.py, splits it into a static part
6+(the per-band gain at 0.1 s) and a time-varying remainder, reports the size of the
7+remainder, and writes
8+ OUT_eq.wav the simulation through the static per-band equalizer (one gain per band)
9+ OUT_eq4.wav the same with the equalizer smoothed to a cubic in log frequency (4 numbers)
10+"""
11+import sys, numpy as np
12+from scipy.signal import stft, istft
13+import match as M
14+
15+ref_path, sim_path, out = sys.argv[1:4]
16+FS, NPERSEG, HOP = M.FS, M.NPERSEG, M.HOP
17+ref = M.load(ref_path)
18+sim = M.load(sim_path)
19+o_ref, o_sim = M.onset(ref), M.onset(sim)
20+noise = ref[: max(int(0.2 * FS), o_ref - int(0.05 * FS))]
21+ref, sim = ref[o_ref:], sim[o_sim:]
22+ref_full = ref
23+n = min(len(ref), len(sim))
24+ref, sim = ref[:n], sim[:n]
25+f, t, Zr = stft(ref, FS, nperseg=NPERSEG, noverlap=NPERSEG - HOP)
26+_, _, Zs = stft(sim, FS, nperseg=NPERSEG, noverlap=NPERSEG - HOP)
27+Er, Es = M.band_envelopes(np.abs(Zr) ** 2, f), M.band_envelopes(np.abs(Zs) ** 2, f)
28+_, _, Zn = stft(noise, FS, nperseg=NPERSEG, noverlap=NPERSEG - HOP)
29+En = M.band_envelopes(np.abs(Zn) ** 2, f)
30+_, _, Zt = stft(ref_full[-int(0.5 * FS):], FS, nperseg=NPERSEG, noverlap=NPERSEG - HOP)
31+Et = M.band_envelopes(np.abs(Zt) ** 2, f)
32+floor = np.maximum(np.median(En, axis=1) + 3, np.median(Et, axis=1) - 3)
33+i_ref = np.searchsorted(t, 0.1)
34+C = M.CENTERS
35+
36+G0 = (Er - Es)[:, i_ref] # static part: balance at 0.1 s
37+print("static equalizer (dB gain applied to the simulation per third-octave band):")
38+for c, g in zip(C, G0):
39+ print(f" {c:7.0f} Hz {g:+6.1f}")
40+# Time-varying remainder over the frames where the recording is above its floor.
41+print("\nremainder after the static part (recording above floor + 10 dB):")
42+print(" band trusted until remainder RMS remainder at end of trusted range")
43+for b in range(len(C)):
44+ below = np.where(Er[b, i_ref:] < floor[b] + 10)[0]
45+ last = i_ref + (below[0] if len(below) else len(t) - i_ref) - 1
46+ if last <= i_ref + 5:
47+ continue
48+ r = (Er - Es)[b, i_ref:last] - G0[b]
49+ print(f" {C[b]:6.0f} Hz {t[last]:5.2f} s {np.sqrt(np.mean(r ** 2)):5.1f} dB {r[-1]:+6.1f} dB")
50+
51+# Static equalizer applied to the simulation, full and smoothed.
52+lf = np.log(np.maximum(f, 1.0))
53+lc = np.log(C)
54+usable = (C > 100) & (C < 6000) # where the recording informs the balance
55+coef = np.polyfit(lc[usable], G0[usable], 3)
56+G4 = np.polyval(coef, np.clip(lc, lc[usable].min(), lc[usable].max()))
57+print("\ncubic in log frequency (4 coefficients) versus the 28 band gains, RMS difference over 100 Hz to 6 kHz:",
58+ f"{np.sqrt(np.mean((G4 - G0)[usable] ** 2)):.1f} dB")
59+i_a, i_b = int(0.1 * FS), int(1.0 * FS)
60+body = lambda x: np.sqrt(np.mean(x[i_a:i_b] ** 2))
61+target = body(ref)
62+for name, G in (("eq", G0), ("eq4", G4)):
63+ Gf = np.interp(lf, lc, G)
64+ _, y = istft(Zs * 10 ** (Gf / 20)[:, None], FS, nperseg=NPERSEG, noverlap=NPERSEG - HOP)
65+ y = y[:n]
66+ M.save(f"{out}_{name}.wav", np.tanh(y * (target / body(y))))
67+print(f"wrote {out}_eq.wav and {out}_eq4.wav")