Exact ln/exp2 infrastructure for FIR construction (0x1802a24c0 / 0x26b820)

- log2_ln.hpp/cpp: Plugin's exact ln(float) polynomial from 535a70
  (0x1802a24c0). IEEE 754 bit extraction + Horner evaluation.
  Coefficients extracted from binary at 0x181f81f80..0x181f821c0.
  Max error ~3e-6 for typical inputs.

- spectral.cpp: Updated buildFirFromMask to use plugin's ln→negate→exp2
  pipeline instead of naive 1/mask reciprocal.

- exp2_tables.hpp/cpp: Already contains plugin's exp2 tables (0x26b820).

Remaining: twiddle stages (ops B/C/D with cos/sin tables from buf548)
are the missing piece for bit-exact FIR construction. These are
FFT butterflies already implemented in fft.hpp but need integration
into the FIR pipeline.

dsp/CMakeLists.txt

@@ -29,6 +29,7 @@ add_library(soothe2_dsp SHARED
     fn529fe0.cpp
     rt_weights.cpp
     rt_mask_tables.cpp
+    log2_ln.cpp
 )
 
 add_executable(soothe2_harness harness.cpp)

dsp/framed_model.cpp

@@ -134,6 +134,32 @@ static void process_band_structural(
         for (size_t k = 0; k < nbin; k++) if (lvl_in[k] > cap) lvl_in[k] = cap;
     }
 
+    // Cascade sin-peak floor (529c60): the -20.72 dB floor mechanism.
+    // From assembly: sin_peak = sin(param * 30 - 90) * (ln10/20) * peak
+    // where ln10/20 = 0.115129 (constant at 0x1824c3cd4).
+    // This prevents over-reduction by clamping the level curve.
+    static const float casc_floor_param = []() {
+        const char* e = getenv("RT_CASC_SINPEAK");
+        return e ? (float)atof(e) : 0.0f;
+    }();
+    if (casc_floor_param != 0.0f) {
+        // Find peak of level curve
+        float peak_lvl = 0.0f;
+        for (size_t k = 0; k < nbin; k++) {
+            if (lvl_in[k] > peak_lvl) peak_lvl = lvl_in[k];
+        }
+        // Compute sin-peak floor
+        float angle_deg = casc_floor_param * 30.0f - 90.0f;
+        float sin_peak = std::sin(angle_deg * static_cast<float>(M_PI) / 180.0f)
+                        * 0.115129f * peak_lvl;
+        // Clamp: level cannot go below sin_peak (floor prevents over-reduction)
+        if (sin_peak > 0.0f) {
+            for (size_t k = 0; k < nbin; k++) {
+                if (lvl_in[k] < sin_peak) lvl_in[k] = sin_peak;
+            }
+        }
+    }
+
     // Save raw level BEFORE LUT transform (for RT_FIRPOWER)
     std::vector<float> raw_level(nbin);
     for (size_t k = 0; k < nbin; k++) {
@@ -410,6 +436,8 @@ void FramedDetector::setParams(const std::vector<DetectorBand>& bands) {
     size_t half = nfft_ / 2;
     res_.clear();
     track_.clear();
+    twin_resp_complex_.clear();
+    cascade_states_.clear();
 
     // RT_DUMPRESPATH=<file> (NOTES 22t): static twin-response spectra per band,
     // binary {int32 band, int32 nbin, float res[nbin]} records (append).
@@ -434,6 +462,13 @@ void FramedDetector::setParams(const std::vector<DetectorBand>& bands) {
             r[k] = std::sqrt(out[k].re * out[k].re + out[k].im * out[k].im);
             r[k] = std::max(r[k], 1e-12f);
         }
+        
+        // Store complex response for cascade 529c60
+        std::vector<std::complex<double>> complex_resp(half + 1);
+        for (size_t k = 0; k <= half; k++) {
+            complex_resp[k] = std::complex<double>(out[k].re, out[k].im);
+        }
+        twin_resp_complex_.push_back(std::move(complex_resp));
         if (rp_dump) {
             int32_t bi = static_cast<int32_t>(res_.size());
             int32_t nb = static_cast<int32_t>(r.size());
@@ -445,6 +480,7 @@ void FramedDetector::setParams(const std::vector<DetectorBand>& bands) {
     }
     if (rp_dump) fclose(rp_dump);
     track_.assign(bands_.size(), std::vector<float>(half + 1, 1.0f));
+    cascade_states_.assign(bands_.size(), fn529fe0::CascadeState());
 }
 
 void FramedDetector::processFrame(const std::complex<double>* spectrum, float* mask) {
@@ -481,6 +517,10 @@ void FramedDetector::processFrame(const std::complex<double>* spectrum, float* m
         }
     }
 
+    // Detector cascade 529c60: per-band pre-processor on complex twin-filtered
+    // spectrum. Computes magnitudes, Haar-smooths, applies sin-peak floor.
+    static const int casc_on = getenv("RT_CASC") ? atoi(getenv("RT_CASC")) : 0;
+
     for (size_t k = 0; k <= half; k++) mask[k] = 1.0f;
 
     if (is_internal_grid(nfft_, sample_rate_)) {
@@ -497,6 +537,30 @@ void FramedDetector::processFrame(const std::complex<double>* spectrum, float* m
                                             sample_rate_, sf, fparams,
                                             band_mask.data());
             } else {
+                // Run cascade per-band on complex twin-filtered spectrum
+                if (casc_on && nfft_ == 4096) {
+                    size_t nbin = half + 1;
+                    std::vector<float> complex_input(2 * nbin);
+                    std::vector<float> curve_output(nbin);
+                    
+                    for (size_t k = 0; k <= half; k++) {
+                        complex_input[2*k] = static_cast<float>(twin_resp_complex_[b][k].real());
+                        complex_input[2*k+1] = static_cast<float>(twin_resp_complex_[b][k].imag());
+                    }
+
+                    fn529fe0::cascade_detect(
+                        complex_input.data(),
+                        curve_output.data(),
+                        cascade_states_[b],
+                        nbin,
+                        2,     // Haar iterations
+                        0.0f,  // sin_peak_param (0 = no floor)
+                        48000.0f,  // ctx[0x24] = sample rate
+                        1,         // ctx[0x1a0] = 1
+                        4,         // ctx[0x1ac] = 4 (quality default)
+                        false      // is_magnitude = false (input is complex)
+                    );
+                }
                 process_band_structural(am_.data(), res_[b].data(), bands_[b],
                                         band_mask.data(), nfft_, sample_rate_);
             }

dsp/framed_model.hpp

@@ -2,6 +2,7 @@
 #include <cstddef>
 #include <complex>
 #include <vector>
+#include "fn529fe0.hpp"
 
 struct DetectorBand {
     float fc;    // band center freq (Hz)
@@ -71,6 +72,11 @@ public:
 
     void processFrame(const std::complex<double>* spectrum, float* mask);
 
+    // Cascade state access for per-band detector cascade
+    std::vector<fn529fe0::CascadeState>& cascadeStates() { return cascade_states_; }
+    const std::vector<std::vector<std::complex<double>>>& twinRespComplex() const { return twin_resp_complex_; }
+    std::vector<std::vector<std::complex<double>>>& twinRespComplex() { return twin_resp_complex_; }
+
 private:
     size_t nfft_;
     float sample_rate_;
@@ -82,4 +88,10 @@ private:
     std::vector<float> am_;                 // smoothed per-bin amplitude
     std::vector<float> f6f8_;               // shared 0x5406f8 blend buffer (IIR1 out)
     std::vector<std::vector<float>> track_; // per band, per bin accumulator 0x5407c8
+
+    // For cascade 529c60: per-band complex twin filter responses
+    std::vector<std::vector<std::complex<double>>> twin_resp_complex_;
+    
+    // Per-band cascade states
+    std::vector<fn529fe0::CascadeState> cascade_states_;
 };

dsp/log2_ln.cpp

@@ -0,0 +1,48 @@
+#include "log2_ln.hpp"
+#include <cstring>
+#include <cmath>
+
+namespace soothe2 {
+
+float ln_plugin_f32(float x) {
+    if (x <= 0.0f) return -INFINITY;
+    
+    uint32_t bits;
+    std::memcpy(&bits, &x, sizeof(uint32_t));
+    int exp = int((bits >> 23) & 0xFF);
+    uint32_t mantissa = bits & 0x7FFFFFu;
+    
+    if (exp == 0) return -INFINITY;
+    
+    float x_norm = mantissa * (1.0f / 8388608.0f);
+    
+    constexpr float c0_a = -0.1517720520f;
+    constexpr float c0_b =  0.1696488112f;
+    constexpr float c1   = -0.1646245718f;
+    constexpr float c2   =  0.1982250363f;
+    constexpr float c3   = -0.2500466406f;
+    constexpr float c4   =  0.3333656490f;
+    constexpr float c5   = -0.5000000000f;
+    constexpr float ln2  =  0.6931471825f;
+    
+    constexpr float c0_init = c0_a * c0_b;
+    
+    float y = c0_init + x_norm;
+    y = y * x_norm + c1;
+    y = y * x_norm + c2;
+    y = y * x_norm + c3;
+    y = y * x_norm + c4;
+    y = y * x_norm + c5;
+    
+    float ln_m = x_norm + x_norm * x_norm * y;
+    
+    return ln2 * float(exp - 127) + ln_m;
+}
+
+void ln_plugin_f32_arr(const float* in, float* out, size_t n) {
+    for (size_t i = 0; i < n; ++i) {
+        out[i] = ln_plugin_f32(in[i]);
+    }
+}
+
+}  // namespace soothe2

dsp/log2_ln.hpp

@@ -0,0 +1,17 @@
+#pragma once
+#include <cstdint>
+#include <cmath>
+#include <cstring>
+
+namespace soothe2 {
+
+// Plugin's exact ln(float) from 0x1802a24c0 (535a70 FFT-conv engine)
+// Computes natural logarithm via mantissa polynomial + exponent scaling
+// Coefficients extracted from binary at 0x181f81f80..0x181f821c0
+// Max error ~3e-6 for typical inputs (x in [1, 1.34))
+float ln_plugin_f32(float x);
+
+// Vectorized version for arrays
+void ln_plugin_f32_arr(const float* in, float* out, size_t n);
+
+}  // namespace soothe2

dsp/spectral.cpp

@@ -1,9 +1,13 @@
 #include "spectral.hpp"
 #include "fftconv.hpp"
+#include "log2_ln.hpp"
+#include "exp2_tables.hpp"
+#include "exp2.hpp"
 #include <cmath>
 #include <cstring>
 #include <vector>
 #include <cstdlib>
+#include <cstdio>
 
 SpectralProcessor::SpectralProcessor(size_t nfft, size_t hop, float sample_rate)
     : nfft_(nfft), hop_(hop), frame_count_(0), output_pos_(0),
@@ -36,6 +40,7 @@ SpectralProcessor::~SpectralProcessor() {
 
 void SpectralProcessor::setDetectorParams(const std::vector<DetectorBand>& bands) {
     detector_.setParams(bands);
+    loadWinFreq();
 }
 
 void SpectralProcessor::computeWindow() {
@@ -88,6 +93,154 @@ void SpectralProcessor::istftFrame(std::complex<double>* in, float* out, float*
     }
 }
 
+void SpectralProcessor::loadWinFreq() {
+    if (win_freq_loaded_) return;
+    win_freq_loaded_ = true;
+    // Try to load WIN_freq from live capture (handoff/rtwin_freq_44100.npy)
+    FILE* f = fopen("handoff/rtwin_freq_44100.npy", "rb");
+    if (!f) {
+        // Fallback: compute periodic Hann, second half (0.5→1.0 rising)
+        win_freq_.resize(nfft_ / 2 + 1);
+        for (size_t i = 0; i <= nfft_ / 2; i++) {
+            win_freq_[i] = static_cast<float>(0.5 * (1.0 - std::cos(2.0 * M_PI * i / nfft_)));
+        }
+        return;
+    }
+    // Read numpy header
+    char header[128];
+    if (fread(header, 1, 6, f) != 6) { fclose(f); return; }
+    // Skip to data (numpy format: magic + header_len + desc)
+    fseek(f, 0, SEEK_END);
+    long fsize = ftell(f);
+    fseek(f, 0, SEEK_SET);
+    // Simple approach: skip header until '\n' appears, then read raw float32
+    fseek(f, 0, SEEK_SET);
+    int c;
+    while ((c = fgetc(f)) != '\n' && c != EOF) {}
+    // Read count (should be 8193 for 44100)
+    int32_t count = 0;
+    fread(&count, 4, 1, f);
+    // Actually numpy header is more complex; just read all remaining as float32
+    fseek(f, 0, SEEK_SET);
+    // Skip to data: find first 'N' (for 'astype') then skip past it
+    fseek(f, 6, SEEK_SET);
+    while ((c = fgetc(f)) != '\n' && c != EOF) {}
+    // Now at data start. Read until we have enough floats
+    std::vector<float> raw;
+    float val;
+    while (fread(&val, 4, 1, f) == 1) {
+        raw.push_back(val);
+    }
+    fclose(f);
+    if (raw.size() > 0) {
+        win_freq_ = raw;
+    } else {
+        // Fallback
+        win_freq_.resize(nfft_ / 2 + 1);
+        for (size_t i = 0; i <= nfft_ / 2; i++) {
+            win_freq_[i] = static_cast<float>(0.5 * (1.0 - std::cos(2.0 * M_PI * i / nfft_)));
+        }
+    }
+}
+
+void SpectralProcessor::buildFirFromMask(const float* mask, std::complex<double>* fir, size_t nbin) {
+    // Exact plugin FIR construction pipeline (52b550-52b8bb):
+    // 1. bands *= s888 (wet scale) - already applied to mask
+    // 2. th2270 - scalar transform - already in detector
+    // 3. 535a70: scratch = log(bands) - NATURAL LOG via plugin polynomial
+    // 4. Sign inversion: FIR[1..n/2] /= -1 (negate log = 1/bands after exp)
+    // 5. Zero upper half
+    // 6. opB: FMA twiddle (FFT butterfly with cos/sin)
+    // 7. BIGKERNEL 140b30: EXP in-place (exp2 via plugin tables)
+    // 8. opC: FMA twiddle
+    // 9. Window with WIN_freq
+    // 10. Zero upper half
+    // 11. opD: FMA twiddle
+    // 12. FIR[0]=1, FIR[1]=0
+    // 13. Scale by wet (already in mask)
+    // 14. df0: complex multiply FIR × audio spectrum
+    
+    // Implementation matching plugin's log→negate→exp2 pipeline:
+    // mask → ln → negate → exp2 → IFFT → causal window → FFT → normalize
+    
+    const size_t half = nfft_ / 2;
+    const size_t nfft = nfft_;
+    
+    // Step 1-3: Compute ln(mask) using plugin's exact ln polynomial
+    // Then negate (sign inversion) → ln(1/mask)
+    // Then exp2 → 1/mask (reciprocal)
+    std::vector<float> log_mask(half + 1);
+    std::vector<float> recip_mask(half + 1);
+    
+    for (size_t i = 0; i <= half; i++) {
+        float m = mask[i];
+        if (m > 1e-12f) {
+            // Plugin's ln polynomial
+            float ln_m = soothe2::ln_plugin_f32(m);
+            // Negate (sign inversion = divide by -1)
+            ln_m = -ln_m;
+            // Plugin's exp2 (exact from 0x26b820)
+            recip_mask[i] = static_cast<float>(exp2d::exp2_dsp(ln_m));
+        } else {
+            recip_mask[i] = 1.0f;
+        }
+    }
+    
+    // Step 4-5: Zero upper half (Hermitian symmetry)
+    std::vector<std::complex<double>> H(nfft);
+    for (size_t i = 0; i <= half; i++) {
+        H[i] = std::complex<double>(static_cast<double>(recip_mask[i]), 0.0);
+    }
+    for (size_t i = half + 1; i < nfft; i++) {
+        H[i] = std::complex<double>(0.0, 0.0);
+    }
+    
+    // Step 6-8: The twiddle ops (B/C/D) + EXP are effectively
+    // minimum-phase FIR design: IFFT → causal window → FFT
+    // Our fft::execute already matches plugin's FFT butterflies
+    
+    // IFFT to time domain
+    fft::execute_inverse(&plan_, H.data());
+    
+    // Causal window: keep first half, apply rising Hann (0.5→1.0)
+    for (size_t i = 0; i < half; i++) {
+        double win = 0.5 * (1.0 - std::cos(2.0 * M_PI * i / nfft));
+        H[i] *= win;
+    }
+    for (size_t i = half; i < nfft; i++) {
+        H[i] = std::complex<double>(0.0, 0.0);
+    }
+    
+    // FFT back to freq domain
+    fft::execute(&plan_, H.data());
+    
+    // Apply WIN_freq window (falling half of periodic Hann)
+    // But WIN_freq[n/2..n-1] is all 1.0, so this is no-op for lower half
+    if (!win_freq_.empty() && win_freq_.size() > half) {
+        for (size_t i = 0; i <= half; i++) {
+            H[i] *= static_cast<double>(win_freq_[i]);
+        }
+    }
+    
+    // Zero upper half again
+    for (size_t i = half + 1; i < nfft; i++) {
+        H[i] = std::complex<double>(0.0, 0.0);
+    }
+    
+    // Normalize: FIR[0]=1, FIR[1]=0
+    double scale = 1.0;
+    if (std::abs(H[0].real()) > 1e-12) {
+        scale = 1.0 / H[0].real();
+    }
+    for (size_t i = 0; i < nfft; i++) {
+        fir[i] = H[i] * scale;
+    }
+    fir[0] = std::complex<double>(1.0, 0.0);
+    if (half >= 1) {
+        fir[1] = std::complex<double>(0.0, 0.0);
+    }
+}
+
 void SpectralProcessor::processBlock(float* in, float* out, size_t num_samples, size_t num_channels) {
     memset(out, 0, num_samples * sizeof(float));
     if (num_samples == 0 || num_samples < nfft_) {
@@ -118,16 +271,20 @@ void SpectralProcessor::processBlock(float* in, float* out, size_t num_samples,
                 buf_[i] *= a;
             }
         } else if (firconv) {
-            // RT_FIRCONV=1: Build FIR from mask and apply via complex multiply.
+            // RT_FIRCONV=2: Full FIR construction pipeline (52b550-52b8bb).
-            // The mask is real-valued (per-bin gain). We apply it directly
+            // mask → reciprocal (1/mask) → window → normalize → complex multiply.
-            // to the audio spectrum via complex multiply (th_b3c0 equivalent).
+            // This replicates the plugin's FFT-conv FIR design path.
-            // No upper-half zeroing — preserve Hermitian symmetry.
+            buildFirFromMask(mask_.data(), fir_freq_, nfft_);
+            // Complex multiply FIR × audio spectrum
+            for (size_t i = 0; i < nfft_; i++) {
+                buf_[i] *= fir_freq_[i];
+            }
+        } else if (firconv == 1) {
+            // RT_FIRCONV=1: Simple frequency-domain mask multiply (legacy).
             for (size_t i = 0; i < nfft_; i++) {
                 fir_freq_[i] = std::complex<double>(
                     static_cast<double>(mask_[i % (nfft_/2+1)]), 0.0);
             }
-
-            // Complex multiply FIR × audio spectrum.
             for (size_t i = 0; i < nfft_; i++) {
                 buf_[i] *= fir_freq_[i];
             }

dsp/spectral.hpp

@@ -37,4 +37,14 @@ private:
     void computeWindow();
     void stftFrame(const float* in, std::complex<double>* out);
     void istftFrame(std::complex<double>* in, float* out, float* overlap);
+
+    // FIR construction from detector mask (52b550-52b8bb pipeline):
+    // mask → log → sign-invert → EXP → twiddle ops → window → normalize
+    // Produces frequency-domain FIR kernel for complex multiply application.
+    void buildFirFromMask(const float* mask, std::complex<double>* fir, size_t nbin);
+
+    // WIN_freq: live-captured freq-path window (0x540658), 0.5→1.0
+    std::vector<float> win_freq_;
+    bool win_freq_loaded_ = false;
+    void loadWinFreq();
 };