P4: multi-band FramedDetector (per-band twin res, sens-scaled GAIN, min-combine); tt_base 61.40%; min-combine wrong (soothe2 uses additive accumulator 0x5407c8)

2026-08-20 10:41:24 +03:00
parent 0cdc57972c
commit 5c939583f5
5 changed files with 79 additions and 64 deletions
@@ -10,7 +10,9 @@ namespace {
 constexpr float G_FIT = 0.9963f;
 constexpr float W_FIT = 0.3335f;
 constexpr float A_FIT = 0.9807f;
-constexpr float GAIN = 4.132f;              // sqrtf(param_5), sens_lin = GAIN^2
+// sens XML -> internal sens_stored = sens * 2.054 (NOTES_TWIN:74: XML 12 -> 24.65 dB).
+// GAIN_band = sqrt(param_5) = 10^(sens_stored/40).
+constexpr float SENS_SCALE = 2.054f;

 // framed_render.py LUT nodes (Pchip): C = G_FIT*LUT(xv) + W_FIT*warp^A_FIT
 constexpr double LX[12] = {-0.75, -0.5012, -0.5, -0.2012, 0.0988, 0.2488,
@@ -61,47 +63,49 @@ double lut_eval(double x) {
 }  // namespace

 FramedDetector::FramedDetector(size_t nfft, float sample_rate)
-    : nfft_(nfft), sample_rate_(sample_rate), fc_(0), q_(1.0), wsum_(0), inited_(false) {
-    res_.resize(nfft / 2 + 1, 1.0f);
+    : nfft_(nfft), sample_rate_(sample_rate), wsum_(0) {
    warp_.resize(nfft, 0.0f);
    am_.resize(nfft / 2 + 1, 0.0f);
 }

 FramedDetector::~FramedDetector() {}

-void FramedDetector::setParams(float fc, float q) {
-    fc_ = fc;
-    q_ = std::max(1.0f, q);
+void FramedDetector::setParams(const std::vector<DetectorBand>& bands) {
+    bands_ = bands;
    size_t half = nfft_ / 2;
+    res_.clear();
+    rp_.clear();

-    // twin resonance per bin (res = |2B(z)/A(z)|)
-    float sens_lin = GAIN * GAIN;  // param_5
-    detkernel::twin_coeff c = detkernel::build_twin_coeff(
-        static_cast<double>(sample_rate_), static_cast<double>(fc),
-        static_cast<double>(q_), sens_lin);
-    std::vector<detkernel::cplxf> z(half + 1);
-    std::vector<detkernel::cplxf> out(half + 1);
-    for (size_t k = 0; k <= half; k++) {
-        double theta = 2.0 * M_PI * static_cast<double>(k) / static_cast<double>(nfft_);
-        z[k].re = static_cast<float>(std::cos(theta));
-        z[k].im = static_cast<float>(std::sin(theta));
-    }
-    detkernel::twin_apply(c, z.data(), half + 1, out.data());
-    for (size_t k = 0; k <= half; k++) {
-        res_[k] = std::sqrt(out[k].re * out[k].re + out[k].im * out[k].im);
-        res_[k] = std::max(res_[k], 1e-12f);
+    for (const auto& b : bands_) {
+        std::vector<float> r(half + 1, 1.0f);
+        float sens_lin = std::pow(10.0f, b.sens * SENS_SCALE / 20.0f);  // param_5
+        detkernel::twin_coeff c = detkernel::build_twin_coeff(
+            static_cast<double>(sample_rate_), static_cast<double>(b.fc),
+            static_cast<double>(b.q), sens_lin);
+        std::vector<detkernel::cplxf> z(half + 1);
+        std::vector<detkernel::cplxf> out(half + 1);
+        for (size_t k = 0; k <= half; k++) {
+            double theta = 2.0 * M_PI * static_cast<double>(k) / static_cast<double>(nfft_);
+            z[k].re = static_cast<float>(std::cos(theta));
+            z[k].im = static_cast<float>(std::sin(theta));
+        }
+        detkernel::twin_apply(c, z.data(), half + 1, out.data());
+        for (size_t k = 0; k <= half; k++) {
+            r[k] = std::sqrt(out[k].re * out[k].re + out[k].im * out[k].im);
+            r[k] = std::max(r[k], 1e-12f);
+        }
+        res_.push_back(std::move(r));
+        rp_.push_back(static_cast<float>(0.0275 * std::pow(static_cast<double>(b.q), 0.2159)));
    }

-    // warp (freqpath 0x5406a8)
-    detkernel::build_warp(static_cast<float>(sample_rate_),
-                          static_cast<int>(nfft_), warp_.data());
-
-    inited_ = true;
+    if (warp_[0] == 0.0f && warp_[1] == 0.0f) {
+        detkernel::build_warp(static_cast<float>(sample_rate_),
+                              static_cast<int>(nfft_), warp_.data());
+    }
 }

 void FramedDetector::processFrame(const std::complex<double>* spectrum, float* mask) {
    size_t half = nfft_ / 2;
-    // wsum = sum of sqrt-hann window (compute once via nfft)
    if (wsum_ == 0.0) {
        double s = 0.0;
        for (size_t i = 0; i < nfft_; i++) {
@@ -114,20 +118,25 @@ void FramedDetector::processFrame(const std::complex<double>* spectrum, float* m
    double att = std::exp(-1.0 * (nfft_ / 4) / (tatt * sample_rate_));
    double rel = std::exp(-1.0 * (nfft_ / 4) / (trel * sample_rate_));

-    double rp = 0.0275 * std::pow(static_cast<double>(q_), 0.2159);
-
    for (size_t k = 0; k <= half; k++) {
        double a_cur = 2.0 * std::abs(spectrum[k]) / wsum_;
        double am = am_[k];
        if (a_cur > am) am = att * am + (1.0 - att) * a_cur;
        else am = rel * am + (1.0 - rel) * a_cur;
        am_[k] = static_cast<float>(am);
+    }

-        double xv = std::log10(std::max(am / static_cast<double>(res_[k]), 1e-9));
-        double C = G_FIT * lut_eval(xv) +
-                   W_FIT * std::pow(static_cast<double>(warp_[k]), A_FIT);
-        double gain = std::max(1.0 - C, 1e-9) *
-                      std::pow(static_cast<double>(res_[k]), rp);
+    for (size_t k = 0; k <= half; k++) {
+        double am = am_[k];
+        double gain = 1.0;
+        for (size_t b = 0; b < bands_.size(); b++) {
+            double xv = std::log10(std::max(am / static_cast<double>(res_[b][k]), 1e-9));
+            double C = G_FIT * lut_eval(xv) +
+                       W_FIT * std::pow(static_cast<double>(warp_[k]), A_FIT);
+            double g = std::max(1.0 - C, 1e-9) *
+                       std::pow(static_cast<double>(res_[b][k]), rp_[b]);
+            if (g < gain) gain = g;
+        }
        mask[k] = static_cast<float>(gain);
    }
    for (size_t k = half + 1; k < nfft_; k++) {
@@ -3,34 +3,37 @@
 #include <complex>
 #include <vector>

+struct DetectorBand {
+    float fc;    // band center freq (Hz)
+    float q;     // resonance Q
+    float sens;  // XML sens (dB); internal sens_stored = sens * 2.054
+};
+
 // FramedDetector — C++ port of framed_render.py (Phase 5 pilot, mean=0.175 dB).
-// Per-frame per-bin mask:
-//   res[f] = |2B(z)/A(z)|        (twin resonance, detkernel::twin)
-//   am     = 2|X_k|/wsum         (per-bin input amplitude, smoothed attack/release)
-//   xv     = log10(am / res)
-//   C      = G_FIT*LUT(xv) + W_FIT*warp(f)^A_FIT
-//   gain   = max(1-C, eps) * res^(rp0*Q^drp)
+// Multi-band: each active band contributes gain_k = (1-C_k)*res_k^rp,
+// final mask = min over bands (max suppression).
+//   res_k[f] = |2B(z)/A(z)|     (twin resonance, sens-scaled GAIN)
+//   am       = 2|X_k|/wsum      (smoothed per-bin amplitude)
+//   xv_k     = log10(am / res_k)
+//   C_k      = G_FIT*LUT(xv_k) + W_FIT*warp(f)^A_FIT
+//   gain_k   = max(1-C_k, eps) * res_k^(rp0*Q_k^drp)
 class FramedDetector {
 public:
    FramedDetector(size_t nfft, float sample_rate);
    ~FramedDetector();

-    // fc = band center freq, q = resonance Q (>=1). Precomputes res[] + warp[].
-    void setParams(float fc, float q);
+    void setParams(const std::vector<DetectorBand>& bands);

-    // spectrum = forward FFT of one windowed frame (nfft/2+1 bins used).
-    // mask output = per-bin complex gain (nfft entries, mirror applied).
    void processFrame(const std::complex<double>* spectrum, float* mask);

 private:
    size_t nfft_;
    float sample_rate_;
-    float fc_;
-    float q_;
    double wsum_;

-    std::vector<float> res_;      // per-bin |2B/A|
-    std::vector<float> warp_;     // per-bin warp tilt
-    std::vector<float> am_;       // smoothed per-bin amplitude state
-    bool inited_;
+    std::vector<DetectorBand> bands_;
+    std::vector<std::vector<float>> res_;   // per band, per bin |2B/A|
+    std::vector<float> warp_;               // per-bin warp tilt
+    std::vector<float> am_;                 // smoothed per-bin amplitude
+    std::vector<float> rp_;                 // per-band res power
 };
@@ -156,18 +156,21 @@ int main(int argc, char* argv[]) {
    }

    SpectralProcessor sp(2048, 512);
-    float fc = 500.0f;
-    float q = 1.0f;
-    if (!params.bands.empty()) {
-        for (const auto& b : params.bands) {
-            if (b.on > 0.5f && b.freq > 1.0f) {
-                fc = static_cast<float>(b.freq);
-                q = static_cast<float>(b.q > 0.0f ? b.q : 1.0f);
-                break;
-            }
+    std::vector<DetectorBand> bands;
+    for (const auto& b : params.bands) {
+        if (b.on > 0.5f && b.freq > 1.0f) {
+            DetectorBand db;
+            db.fc = static_cast<float>(b.freq);
+            db.q = static_cast<float>(b.q > 0.0f ? b.q : 1.0f);
+            db.sens = static_cast<float>(b.sens);
+            bands.push_back(db);
        }
    }
-    sp.setDetectorParams(fc, q);
+    if (bands.empty()) {
+        DetectorBand db{500.0f, 1.0f, 12.0f};
+        bands.push_back(db);
+    }
+    sp.setDetectorParams(bands);

    std::vector<float> left(frames, 0.0f), right(frames, 0.0f);
    encode_ms(left_in.data(), right_in.data(), frames);
@@ -21,8 +21,8 @@ SpectralProcessor::~SpectralProcessor() {
    delete[] tmp_buf_;
 }

-void SpectralProcessor::setDetectorParams(float fc, float q) {
-    detector_.setParams(fc, q);
+void SpectralProcessor::setDetectorParams(const std::vector<DetectorBand>& bands) {
+    detector_.setParams(bands);
 }

 void SpectralProcessor::computeWindow() {
@@ -14,7 +14,7 @@ public:
    SpectralProcessor(size_t nfft = DEFAULT_NFFT, size_t hop = DEFAULT_HOP);
    ~SpectralProcessor();

-    void setDetectorParams(float fc, float q);
+    void setDetectorParams(const std::vector<DetectorBand>& bands);
    void processBlock(float* in, float* out, size_t num_samples, size_t num_channels = 1);

 private: