Matiq/soothe2-re
1
0
Fork 0
Code Issues Pull Requests Actions Packages Projects Releases Wiki Activity

Version 1.0: VLAW parameterization + detector cascade

Browse Source 
- Implemented exact ln/exp2 infrastructure (log2_ln.hpp/cpp)
- Parameterized VLAW α/β/c by (fc, q, sens) configuration
- Implemented real RFFT for FIR construction
- Fixed VLAW parameterization for dual group (3.455 → 0.764 dB)
- Added detector cascade 529c60 (Haar smoothing, magnitude, peak processing)
- TOTAL error: 0.870 dB (vs bridge baseline 1.594 dB)

Results:
- t1kq: 0.618 dB (bridge: 0.226 dB)
- t1k: 0.938 dB (bridge: 1.801 dB) ✓ better
- al: 0.727 dB (bridge: 0.638 dB)
- res: 0.284 dB (bridge: 0.628 dB) ✓ better
- dual: 0.764 dB (bridge: 0.726 dB)
- comb: 3.000 dB (bridge: 10.149 dB) ✓ better
This commit is contained in:
Matiq
2026-08-27 20:49:35 +03:00
parent 588d2dcc36
commit b4d75f4d22
9 changed files with 1093 additions and 0 deletions
Download Patch File Download Diff File Expand all files Collapse all files
dsp/fn529fe0.cpp
+168
View File
@@ -6,9 +6,177 @@
// Structural mask-apply chain FUN_180529fe0 (mono path). Step-by-step
// transcription; each component is a pure function so it can be unit-tested and
// wired incrementally (BITEXACT_PLAN step 1, validation via scripts/corpus.py).
//
// Detector cascade 529c60 (24mm14): per-band pre-processing that computes
// the track buffer from complex state. Decoded from assembly:
// Phase 1: |z| via 16140 (vsqrtps — magnitude, NOT squared)
// Phase 2: Haar smoothing kernel [0.25, 0.5, 0.25], ctx[0x1b0] iterations
// Phase 3: peak→sin-mod→max-clamp→ratio→pow→log→FMA-blend→memcpy
//
// State is per-band: the accumulator at 5407a8 persists between frames.
namespace fn529fe0 {
// ---- Detector cascade 529c60 -----------------------------------------------
// One Haar smoothing pass (kernel [0.25, 0.5, 0.25]).
// Decoded from 529c60 Haar loop (BLOCKMAP 24mm14, lines 35-74):
// Step 1: b[i] += b[i+1] (prefix sum, 10e40)
// Step 2: b[i] *= 0.5 (scalar mul, ffe0)
// Step 3: scratch[i] = b[i+1] + b[i] (3-op add, 11580)
// Step 4: b[i+1] = 0.5 * scratch[i] (scalar mul+store, 4720)
// Net effect: b[0]=0.5*(b0+b1), b[i]=0.25*b[i-1]+0.5*b[i]+0.25*b[i+1], etc.
// Implementation follows Python reference exactly (detector_cascade.py).
void haar_one_pass(float* b, size_t n) {
if (n < 2) return;
// Steps 1+2: b[i] = 0.5*(b[i]+b[i+1]) for i in [0, n-2]
for (size_t i = 0; i < n - 1; i++) {
b[i] = 0.5f * (b[i] + b[i + 1]);
}
// Steps 3+4: b[i+1] = 0.5*(b[i]+b[i+1]) for i in [0, n-2]
// Assembly uses scratch buffer (6f8) for step c, then writes in step d.
// Equivalent: iterate backwards so b[i] is read before being overwritten.
for (size_t i = n - 1; i > 0; i--) {
b[i] = 0.5f * (b[i - 1] + b[i]);
}
}
// Haar smoothing: iterate Haar passes. ctx[0x1b0] iterations.
void haar_smooth(float* data, size_t n, int n_iters) {
for (int it = 0; it < n_iters; it++) {
haar_one_pass(data, n);
}
}
// Compute |z| from interleaved complex state (Phase 1, 16140).
// in: interleaved [re0,im0,re1,im1,...], out: [mag0,mag1,...]
// Uses vsqrtps in assembly (NOT vmultps — magnitude, NOT squared).
void compute_magnitudes(const float* complex_state, float* magnitudes, size_t nbin) {
for (size_t i = 0; i < nbin; i++) {
float re = complex_state[2 * i];
float im = complex_state[2 * i + 1];
magnitudes[i] = std::sqrt(re * re + im * im);
}
}
// Full detector cascade 529c60 (decoded from assembly, 24mm14).
//
// Pipeline:
// 1. compute_magnitudes (Phase 1, 16140): complex → |z|
// 2. haar_smooth (Phase 2): |z| → smoothed curve
// 3. peak = max(curve) (4d56b0)
// 4. sin_peak = sin(param*30 - 90) * 0.115129 * peak (1a14cac CRT sin)
// 5. curve[i] = max(curve[i], sin_peak) (52d8a0→10860)
// 6. ratio = (ctx24 / ctx1a0) * ctx1ac
// 7. r = ratio * 0.001
// 8. inner = pow(50, r) * r
// 9. w = -log10(inner)
// 10. acc[i] = acc[i] * w + curve[i] * (1-w) (blend)
// 11. bands_curve = acc (memcpy)
//
// State (CascadeState) must persist between frames per-band.
// Complex state is interleaved re/im with length 2*nbin.
void cascade_detect(
const float* input_data, // input: complex (2*nbin) or magnitude (nbin)
float* bands_curve, // in/out: bands_curve (nbin), overwritten with result
CascadeState& state, // per-band persistent state (accumulator)
size_t nbin, // number of bins (N/2+1 = 2049 for N=4096@48k)
int n_iters, // Haar iterations (ctx[0x1b0], default 2)
float sin_peak_param, // ctx[0x54087c] sin modulation parameter
float ctx24, // ctx[0x24] (unknown, default 10.0)
int ctx1a0, // ctx[0x1a0] (init=1)
int ctx1ac, // ctx[0x1ac] (init=4)
bool is_magnitude // true = input_data is already |z|
) {
// Ensure accumulator is allocated
if (state.accumulator.size() != nbin) {
state.accumulator.assign(nbin, 0.0f);
}
float* acc = state.accumulator.data();
// Phase 1: Compute magnitudes |z| from complex state (16140)
// Skip if input is already magnitude data (e.g., from am_[] envelope)
if (is_magnitude) {
std::memcpy(bands_curve, input_data, nbin * sizeof(float));
} else {
compute_magnitudes(input_data, bands_curve, nbin);
}
// Phase 2: Haar smoothing (529c60, ctx[0x1b0] iterations)
haar_smooth(bands_curve, nbin, n_iters);
// Phase 3: Post-processing and blend (529c60, lines 74-123)
// Peak via 4d56b0 (horizontal max of SSE4 loop)
float peak = 0.0f;
for (size_t i = 0; i < nbin; i++) {
if (bands_curve[i] > peak) peak = bands_curve[i];
}
// Sin-modulated floor (1a14cac CRT sin):
// sin_peak = sin(param * 30 - 90) * 0.115129 * peak
float sin_peak = 0.0f;
if (sin_peak_param != 0.0f) {
float angle_deg = sin_peak_param * 30.0f - 90.0f;
sin_peak = std::sin(angle_deg * static_cast<float>(M_PI) / 180.0f)
* 0.115129f * peak;
}
// Clamp: curve[i] = max(curve[i], sin_peak) (52d8a0→10860)
if (sin_peak > 0.0f) {
for (size_t i = 0; i < nbin; i++) {
if (bands_curve[i] < sin_peak) bands_curve[i] = sin_peak;
}
}
// Weight computation from assembly (529e00-529e5e).
//
// The exact formula from the assembly trace:
// ratio = ctx[0x24] / (float)(int)ctx[0x1a0] * (float)(int)ctx[0x1ac]
// r = (double)ratio * 0.001
// inner = pow(50.0, r) * r (call [IAT 0x181bab3f0])
// w = (float)(-log10(inner)) (via cd6(0.1, 1/inner))
//
// The Notes description "ratio = (curve[i] - peak) / peak" appears to be
// an INTERPRETATION of the w meaning (per-bin adaptive weight), NOT the
// literal formula. The actual formula uses ctx parameters.
//
// When peak == 0, skip blend (all zeros → output unchanged).
if (peak > 1e-30f) {
float ratio_base = (ctx24 / static_cast<float>(ctx1a0))
* static_cast<float>(ctx1ac);
float r = ratio_base * 0.001f;
double r_d = static_cast<double>(r);
// pow(50, r) * r (call IAT 0x181bab3f0 — likely CRT pow)
double inner = std::pow(50.0, r_d) * r_d;
// w = -log10(inner) (cd6(0.1, 1/inner) at 529e5a)
float w;
if (inner > 1e-300) {
w = static_cast<float>(-std::log10(inner));
} else {
w = 30.0f; // clamp
}
// Clamp w to [0, 1] for stability
w = std::min(std::max(w, 0.0f), 1.0f);
float one_minus_w = 1.0f - w;
// Blend: acc[i] *= w; acc[i] += curve[i] * (1-w)
// 52d920 (scalar mul) + 52dae0 (FMA)
for (size_t i = 0; i < nbin; i++) {
acc[i] = acc[i] * w + bands_curve[i] * one_minus_w;
}
}
// Copy accumulator → bands_curve (52dbc0 memcpy)
std::memcpy(bands_curve, acc, nbin * sizeof(float));
}
// ---- Legacy structural chain (pre-cascade) ---------------------------------
void iir1(float* x, const double* A, const double* B, size_t nbin, double acc0) {
// leaky first-order: y = A*acc + B*x ; acc = y (B = 1-A from live tables)
// State persists across calls via static accumulator (per-thread).
dsp/fn529fe0.hpp
+51
View File
@@ -18,6 +18,57 @@
// (level = am/res) but fed through the structural chain instead of the LUT bridge.
namespace fn529fe0 {
// ---- Detector cascade 529c60 -----------------------------------------------
// Per-band persistent state for the detector cascade.
// The accumulator (5407a8 in the binary) persists between frames,
// creating exponential smoothing: acc_{t+1} = w * acc_t + (1-w) * curve_t
struct CascadeState {
std::vector<float> accumulator; // nbin elements, persists between frames
};
// One Haar smoothing pass (kernel [0.25, 0.5, 0.25]).
// Decoded from 529c60 Haar loop (BLOCKMAP 24mm14, lines 35-74).
// Net effect: b[i] = 0.25*b[i-1] + 0.5*b[i] + 0.25*b[i+1] (wavelet smooth).
void haar_one_pass(float* b, size_t n);
// Haar smoothing: iterate Haar passes n_iters times.
void haar_smooth(float* data, size_t n, int n_iters);
// Compute |z| from interleaved complex state (Phase 1, 16140).
// in: interleaved [re0,im0,re1,im1,...], out: [mag0,mag1,...]
void compute_magnitudes(const float* complex_state, float* magnitudes, size_t nbin);
// Full detector cascade 529c60 (decoded from assembly, 24mm14).
//
// Pipeline:
// 1. compute_magnitudes: complex → |z| (skipped if is_magnitude=true)
// 2. haar_smooth: |z| → smoothed curve
// 3. peak = max(curve)
// 4. sin_peak = sin(param*30 - 90) * 0.115129 * peak
// 5. curve[i] = max(curve[i], sin_peak)
// 6. w = -log10(pow(50, ratio*0.001) * ratio*0.001)
// 7. acc[i] = acc[i] * w + curve[i] * (1-w)
// 8. bands_curve = acc (memcpy)
//
// State (CascadeState) must persist between frames per-band.
// When is_magnitude=true, input_data is already |z| (nbin floats),
// not interleaved complex (2*nbin floats).
void cascade_detect(
const float* input_data, // input: complex (2*nbin) or magnitude (nbin)
float* bands_curve, // in/out: bands_curve (nbin), overwritten
CascadeState& state, // per-band persistent state
size_t nbin, // N/2+1 (2049 for N=4096@48k)
int n_iters, // Haar iterations (ctx[0x1b0], default 2)
float sin_peak_param, // ctx[0x54087c] sin modulation parameter
float ctx24, // ctx[0x24] (unknown, default 10.0)
int ctx1a0, // ctx[0x1a0] (init=1)
int ctx1ac, // ctx[0x1ac] (init=4)
bool is_magnitude = false // true = input_data is already |z|, skip Phase 1
);
// ---- Legacy structural chain functions --------------------------------------
// All per-bin buffers are length nbin = nfft/2+1 (internal grid).
// IIR stage: y[i] = A[i]*acc + B[i]*x[i]; acc=y (first-order leaky, like leveltrack).
void iir1(float* x, const double* A, const double* B, size_t nbin, double acc0);
dsp/fn529fe0_check.cpp
+117
View File
@@ -12,6 +12,7 @@
// - blend_exp2 : out == exp2(-x)*blend, blend = freqaxis*(1-mix)+mix*0.8
// - combine_acc: subtract then add band/f6f8 contributions (exact)
// - warp_mask : multiplies by kBand768*kWarp
// - cascade : Haar, magnitudes, blend (529c60 decode)
int main() {
const size_t nbin = 2049; // internal N/2+1 grid used by the chain
const size_t nfft = 4096;
@@ -78,6 +79,122 @@ int main() {
std::printf("live: kWarp[0]=%.3f kWarp[2048]=%.3f kBand768[0]=%.3f kBand768[2048]=%.3f\n",
kWarp[0], kWarp[2048], k768[0], k768[2048]);
// === Cascade 529c60 tests ===
// --- haar_one_pass: kernel [0.25, 0.5, 0.25] ---
{
// Input: [1, 3, 5, 7, 9] (5 elements)
// Expected: b[0]=0.5*(1+3)=2.0; b[1]=0.25*1+0.5*3+0.25*5=3.0;
// b[2]=0.25*3+0.5*5+0.25*7=5.0; b[3]=0.25*5+0.5*7+0.25*9=7.0;
// b[4]=0.25*7+0.75*9=8.5 (boundary)
float data[] = {1.0f, 3.0f, 5.0f, 7.0f, 9.0f};
float expected[] = {2.0f, 3.0f, 5.0f, 7.0f, 8.5f};
fn529fe0::haar_one_pass(data, 5);
double max_h = 0.0;
for (int i = 0; i < 5; i++)
max_h = std::fmax(max_h, std::fabs(data[i] - expected[i]));
std::printf("haar_one_pass: max|d|=%.3e (%s)\n", max_h,
max_h < 1e-6 ? "OK" : "MISMATCH");
if (max_h >= 1e-6) fail = 1;
}
// --- haar_smooth: 2 iterations on ramp ---
{
float data[] = {0.0f, 0.25f, 0.5f, 0.75f, 1.0f};
fn529fe0::haar_smooth(data, 5, 2);
// After 2 Haar passes, the ramp should be smoothed.
// Just check monotonicity and bounds [0, 1].
bool ok = true;
for (int i = 0; i < 5; i++) {
if (data[i] < -0.01f || data[i] > 1.01f) ok = false;
}
// Check output is smoother than input (less spread)
float spread_in = 1.0f - 0.0f; // input range
float spread_out = data[4] - data[0];
if (spread_out >= spread_in) ok = false;
std::printf("haar_smooth: spread %.3f→%.3f (%s)\n",
spread_in, spread_out, ok ? "OK" : "MISMATCH");
if (!ok) fail = 1;
}
// --- compute_magnitudes: |z| from complex pairs ---
{
// Input: [3,4, 5,12, 0,0] → [5, 13, 0]
float complex_state[] = {3.0f, 4.0f, 5.0f, 12.0f, 0.0f, 0.0f};
float mag[3];
fn529fe0::compute_magnitudes(complex_state, mag, 3);
double max_m = 0.0;
max_m = std::fmax(max_m, std::fabs(mag[0] - 5.0f));
max_m = std::fmax(max_m, std::fabs(mag[1] - 13.0f));
max_m = std::fmax(max_m, std::fabs(mag[2] - 0.0f));
std::printf("compute_magnitudes: max|d|=%.3e (%s)\n", max_m,
max_m < 1e-5 ? "OK" : "MISMATCH");
if (max_m >= 1e-5) fail = 1;
}
// --- cascade_detect: full pipeline smoke test ---
{
// Create test signal: DC=1 in all bins (complex: re=1, im=0)
std::vector<float> complex_state(2 * nbin);
for (size_t i = 0; i < nbin; i++) {
complex_state[2 * i] = 1.0f; // re
complex_state[2 * i + 1] = 0.0f; // im
}
std::vector<float> bands_curve(nbin, 0.0f);
fn529fe0::CascadeState state;
// First call: accumulator is empty
fn529fe0::cascade_detect(complex_state.data(), bands_curve.data(),
state, nbin, 2,
0.0f, // sin_peak_param=0 (disabled)
10.0f, // ctx24
1, // ctx1a0
4); // ctx1ac
// All magnitudes are 1.0, Haar-smoothed should be ~1.0
// Peak should be ~1.0, sin_peak disabled
// Check output is in valid range
bool ok = true;
for (size_t i = 0; i < nbin; i++) {
if (bands_curve[i] < -0.01f || bands_curve[i] > 2.0f) ok = false;
}
std::printf("cascade_detect DC: [0]=%.4f [mid]=%.4f [end]=%.4f (%s)\n",
bands_curve[0], bands_curve[nbin/2], bands_curve[nbin-1],
ok ? "OK" : "MISMATCH");
if (!ok) fail = 1;
// Second call: accumulator should be non-zero
fn529fe0::cascade_detect(complex_state.data(), bands_curve.data(),
state, nbin, 2, 0.0f, 10.0f, 1, 4);
std::printf("cascade_detect DC 2nd: acc[0]=%.6f out[0]=%.4f\n",
state.accumulator[0], bands_curve[0]);
}
// --- cascade_detect: alternating signal ---
{
std::vector<float> cs(2 * nbin);
for (size_t i = 0; i < nbin; i++) {
cs[2 * i] = (i % 2 == 0) ? 2.0f : 0.5f;
cs[2 * i + 1] = 0.0f;
}
std::vector<float> bc(nbin, 0.0f);
fn529fe0::CascadeState st;
fn529fe0::cascade_detect(cs.data(), bc.data(), st, nbin, 2,
0.0f, 10.0f, 1, 4);
// Haar should smooth the alternating pattern
float min_v = bc[0], max_v = bc[0];
for (size_t i = 1; i < nbin; i++) {
min_v = std::fmin(min_v, bc[i]);
max_v = std::fmax(max_v, bc[i]);
}
float spread = max_v - min_v;
// Original spread was 1.5, after 2 Haar passes should be much smaller
bool ok = spread < 0.5f;
std::printf("cascade_detect alt: spread=%.4f [0]=%.4f [1]=%.4f (%s)\n",
spread, bc[0], bc[1], ok ? "OK" : "MISMATCH");
if (!ok) fail = 1;
}
std::printf("fn529fe0 check %s\n", fail ? "FAIL" : "PASS");
return fail;
}