Matiq
847725f5fc
- Fixed execute_inverse: removed conj bug, now uses positive twiddle only - Added WOLA normalization factor (wola_sum/hop_ for Hann+hop=N/4) - New detector model: floor(level) + bell_curve * boost, calibrated from measured data (summary.md level sweep, 7 data points) - Transcribed twiddle loader (FUN_18014ec20): Cody-Waite 4-level reduction with minimax sin/cos polynomial, constants from Frida memory dump - Added soothe_constants.hpp with extracted polynomial coefficients - HARNESS parameters updated to match burst500_b1.rpp (depth=0.864) Results: burst500.wav reduction now -5.8 dB vs Ref -6.8 dB (was -14.4 dB)
90 lines
2.7 KiB
C++
90 lines
2.7 KiB
C++
#include "spectral.hpp"
|
|
#include <cmath>
|
|
#include <cstring>
|
|
#include <vector>
|
|
|
|
SpectralProcessor::SpectralProcessor(size_t nfft, size_t hop)
|
|
: nfft_(nfft), hop_(hop), frame_count_(0), output_pos_(0),
|
|
detector_(nfft, 44100.0f) {
|
|
window_ = new double[nfft_];
|
|
computeWindow();
|
|
fft::init_plan(&plan_, static_cast<uint32_t>(std::log2(nfft_)));
|
|
buf_ = new std::complex<double>[nfft_];
|
|
tmp_buf_ = new std::complex<double>[nfft_];
|
|
overlap_.resize(nfft_, 0.0f);
|
|
mask_.resize(nfft_, 1.0f);
|
|
}
|
|
|
|
SpectralProcessor::~SpectralProcessor() {
|
|
delete[] window_;
|
|
delete[] buf_;
|
|
delete[] tmp_buf_;
|
|
}
|
|
|
|
void SpectralProcessor::setDetectorParams(float sharpness, float selectivity, float depth) {
|
|
detector_.setParams(sharpness, selectivity, depth);
|
|
}
|
|
|
|
void SpectralProcessor::computeWindow() {
|
|
for (size_t i = 0; i < nfft_; i++) {
|
|
window_[i] = 0.5 * (1.0 - std::cos(2.0 * M_PI * i / (nfft_ - 1)));
|
|
}
|
|
}
|
|
|
|
void SpectralProcessor::stftFrame(const float* in, std::complex<double>* out) {
|
|
for (size_t i = 0; i < nfft_; i++) {
|
|
out[i] = std::complex<double>(static_cast<double>(in[i]) * window_[i], 0.0);
|
|
}
|
|
fft::execute(&plan_, out);
|
|
}
|
|
|
|
void SpectralProcessor::istftFrame(std::complex<double>* in, float* out, float* overlap) {
|
|
memcpy(tmp_buf_, in, nfft_ * sizeof(std::complex<double>));
|
|
fft::execute_inverse(&plan_, tmp_buf_);
|
|
static bool wola_computed = false;
|
|
static float wola_norm = 1.0f;
|
|
if (!wola_computed) {
|
|
double wola_sum = 0.0;
|
|
for (size_t i = 0; i < nfft_; i++) {
|
|
wola_sum += window_[i] * window_[i];
|
|
}
|
|
wola_norm = static_cast<float>(wola_sum / hop_);
|
|
wola_computed = true;
|
|
}
|
|
for (size_t i = 0; i < nfft_; i++) {
|
|
overlap[i] += static_cast<float>(tmp_buf_[i].real() * window_[i]);
|
|
}
|
|
for (size_t i = 0; i < hop_; i++) {
|
|
out[i] = overlap[i] / wola_norm;
|
|
}
|
|
for (size_t i = 0; i < nfft_ - hop_; i++) {
|
|
overlap[i] = overlap[i + hop_];
|
|
}
|
|
for (size_t i = nfft_ - hop_; i < nfft_; i++) {
|
|
overlap[i] = 0.0f;
|
|
}
|
|
}
|
|
|
|
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_) {
|
|
return;
|
|
}
|
|
|
|
size_t nframes = (num_samples - nfft_) / hop_ + 1;
|
|
|
|
for (size_t f = 0; f < nframes; f++) {
|
|
size_t offset = f * hop_;
|
|
if (offset + nfft_ > num_samples) break;
|
|
stftFrame(in + offset, buf_);
|
|
|
|
detector_.processFrame(buf_, mask_.data());
|
|
|
|
for (size_t i = 0; i < nfft_; i++) {
|
|
buf_[i] *= mask_[i];
|
|
}
|
|
|
|
istftFrame(buf_, out + offset, overlap_.data());
|
|
}
|
|
}
|