diff --git a/src/ChoppingKinky.cpp b/src/ChoppingKinky.cpp index 52c11d3..9d00858 100644 --- a/src/ChoppingKinky.cpp +++ b/src/ChoppingKinky.cpp @@ -67,19 +67,14 @@ struct ChoppingKinky : Module { dsp::SchmittTrigger trigger; bool output_a_to_chopp; - float previous_a = 0.0; - VariableOversampling<> oversampler[NUM_CHANNELS]; int oversamplingIndex = 2; - - dsp::RCFilter choppDC; + + dsp::BiquadFilter blockDCFilter; bool blockDC = false; - unsigned int currentFrame = 0; - - ChoppingKinky() { config(NUM_PARAMS, NUM_INPUTS, NUM_OUTPUTS, NUM_LIGHTS); configParam(FOLD_A_PARAM, 0.f, 2.f, 0.f, ""); @@ -89,14 +84,13 @@ struct ChoppingKinky : Module { // calculate up/downsampling rates onSampleRateChange(); - } void onSampleRateChange() override { // SampleRateConverter needs integer value int sampleRate = APP->engine->getSampleRate(); - - choppDC.setCutoff(10.3f / sampleRate); + + blockDCFilter.setParameters(dsp::BiquadFilter::HIGHPASS, 10.3f / sampleRate, M_SQRT1_2, 1.0f); for (int channel_idx = 0; channel_idx < NUM_CHANNELS; channel_idx++) { oversampler[channel_idx].setOversamplingIndex(oversamplingIndex); @@ -150,33 +144,46 @@ struct ChoppingKinky : Module { } } + bool chopp_is_required = outputs[OUT_CHOPP_OUTPUT].isConnected(); + bool a_is_required = outputs[OUT_A_OUTPUT].isConnected() || chopp_is_required; + bool b_is_required = outputs[OUT_B_OUTPUT].isConnected() || chopp_is_required; + in_a = in_a * gain_a; in_b = in_b * gain_b; float in_chopp = output_a_to_chopp ? 1.f : 0.f; - oversampler[CHANNEL_A].upsample(in_a); - oversampler[CHANNEL_B].upsample(in_b); - oversampler[CHANNEL_CHOPP].upsample(in_chopp); + if (a_is_required) { + oversampler[CHANNEL_A].upsample(in_a); + } + if (b_is_required) { + oversampler[CHANNEL_B].upsample(in_b); + } + if (chopp_is_required) { + oversampler[CHANNEL_CHOPP].upsample(in_chopp); + } float* osBufferA = oversampler[CHANNEL_A].getOSBuffer(); float* osBufferB = oversampler[CHANNEL_B].getOSBuffer(); float* osBufferChopp = oversampler[CHANNEL_CHOPP].getOSBuffer(); for (int i = 0; i < oversampler[0].getOversamplingRatio(); i++) { - // TOO SLOW - //break; - osBufferA[i] = foldResponse(4.0f * osBufferA[i], -0.2); - osBufferB[i] = foldResponseSin(4.0f * osBufferB[i]); - osBufferChopp[i] = osBufferChopp[i] * osBufferA[i] + (1.f - osBufferChopp[i]) * osBufferB[i]; + if (a_is_required) { + osBufferA[i] = foldResponse(4.0f * osBufferA[i], -0.2); + } + if (b_is_required) { + osBufferB[i] = foldResponseSin(4.0f * osBufferB[i]); + } + if (chopp_is_required) { + osBufferChopp[i] = osBufferChopp[i] * osBufferA[i] + (1.f - osBufferChopp[i]) * osBufferB[i]; + } } - float out_a = oversampler[CHANNEL_A].downsample(); - float out_b = oversampler[CHANNEL_B].downsample(); - float out_chopp = oversampler[CHANNEL_CHOPP].downsample(); + float out_a = a_is_required ? oversampler[CHANNEL_A].downsample() : 0.f; + float out_b = b_is_required ? oversampler[CHANNEL_B].downsample() : 0.f; + float out_chopp = chopp_is_required ? oversampler[CHANNEL_CHOPP].downsample() : 0.f; if (blockDC) { - choppDC.process(out_chopp); - out_chopp -= choppDC.lowpass(); + out_chopp = blockDCFilter.process(out_chopp); } previous_a = in_a; diff --git a/src/VariableOversampling.hpp b/src/VariableOversampling.hpp index 3220a6e..d796624 100644 --- a/src/VariableOversampling.hpp +++ b/src/VariableOversampling.hpp @@ -2,336 +2,86 @@ #include #include -#include -/** Digital IIR filter processor. Using TDF-II structure: - * https://ccrma.stanford.edu/~jos/filters/Transposed_Direct_Forms.html -*/ -template -struct IIRFilter { - /** transfer function numerator coefficients: b_0, b_1, etc.*/ - T b[ORDER] = {}; - - /** transfer function denominator coefficients: a_0, a_1, etc.*/ - T a[ORDER] = {}; - /** filter state */ - T z[ORDER]; +/** + High-order filter to be used for anti-aliasing or anti-imaging. + The template parameter N should be 1/2 the desired filter order. - IIRFilter() { - reset(); - } + Currently uses an 2*N-th order Butterworth filter. + @TODO: implement Chebyshev, Elliptic filter options. +*/ +template +class AAFilter { +public: + AAFilter() = default; - void reset() { - std::fill(z, &z[ORDER], 0.0f); - } + /** Calculate Q values for a Butterworth filter of a given order */ + static std::vector calculateButterQs(int order) { + const int lim = int (order / 2); + std::vector Qs; - void setCoefficients(const T* b, const T* a) { - for (int i = 0; i < ORDER; i++) { - this->b[i] = b[i]; + for (int k = 1; k <= lim; ++k) { + auto b = -2.0f * std::cos((2.0f * k + order - 1) * 3.14159 / (2.0f * order)); + Qs.push_back(1.0f / b); } - for (int i = 1; i < ORDER; i++) { - this->a[i] = a[i]; - } - } - template - inline typename std::enable_if ::type - process(T x) noexcept { - T y = z[1] + x * b[0]; - z[1] = x * b[1] - y * a[1]; - return y; - } - - template - inline typename std::enable_if ::type - process(T x) noexcept { - T y = z[1] + x * b[0]; - z[1] = z[2] + x * b[1] - y * a[1]; - z[2] = x * b[2] - y * a[2]; - return y; - } - - template - inline typename std::enable_if <(N > 3), T>::type - process(T x) noexcept { - T y = z[1] + x * b[0]; - - for (int i = 1; i < ORDER-1; ++i) - z[i] = z[i+1] + x * b[i] - y * a[i]; - - z[ORDER-1] = x * b[ORDER-1] - y * a[ORDER-1]; - - return y; - } - - /** Computes the complex transfer function $H(s)$ at a particular frequency - s: normalized angular frequency equal to $2 \pi f / f_{sr}$ ($\pi$ is the Nyquist frequency) - */ - std::complex getTransferFunction(T s) { - // Compute sum(a_k z^-k) / sum(b_k z^-k) where z = e^(i s) - std::complex bSum(b[0], 0); - std::complex aSum(1, 0); - for (int i = 1; i < ORDER; i++) { - T p = -i * s; - std::complex z(simd::cos(p), simd::sin(p)); - bSum += b[i] * z; - aSum += a[i - 1] * z; - } - return bSum / aSum; + std::reverse(Qs.begin(), Qs.end()); + return Qs; } - T getFrequencyResponse(T f) { - return simd::abs(getTransferFunction(2 * M_PI * f)); + /** + * Resets the filter to process at a new sample rate. + * + * @param sampleRate: The base (i.e. pre-oversampling) sample rate of the audio being processed + * @param osRatio: The oversampling ratio at which the filter is being used + */ + void reset(float sampleRate, int osRatio) { + float fc = 0.98f * (sampleRate / 2.0f); + auto Qs = calculateButterQs(2 * N); + + for (int i = 0; i < N; ++i) + filters[i].setParameters(dsp::BiquadFilter::Type::LOWPASS, fc / (osRatio * sampleRate), Qs[i], 1.0f); } - T getFrequencyPhase(T f) { - return simd::arg(getTransferFunction(2 * M_PI * f)); - } -}; - -template -struct TOnePole : IIRFilter<2, T> { - enum Type { - LOWPASS, - HIGHPASS, - NUM_TYPES - }; - - TOnePole() { - setParameters(LOWPASS, 0.0f); - } - - /** Calculates and sets the one-pole transfer function coefficients. - f: normalized frequency (cutoff frequency / sample rate), must be less than 0.5 - */ - void setParameters(Type type, float f) { - switch (type) { - case LOWPASS: { - this->a[1] = -std::exp(-2.f * M_PI * f); - this->b[0] = 1.f + this->a[0]; - } break; - - case HIGHPASS: { - this->a[1] = std::exp(-2.f * M_PI * (0.5f - f)); - this->b[0] = 1.f - this->a[0]; - } break; - - default: break; - } - } -}; + inline float process(float x) noexcept { + for (int i = 0; i < N; ++i) + x = filters[i].process(x); -typedef TOnePole<> OnePole; - -template -struct TBiquadFilter : IIRFilter<3, T> { - enum Type { - LOWPASS, - HIGHPASS, - LOWSHELF, - HIGHSHELF, - BANDPASS, - PEAK, - NOTCH, - NUM_TYPES - }; - - TBiquadFilter() { - setParameters(LOWPASS, 0.f, 0.f, 1.f); + return x; } - /** Calculates and sets the biquad transfer function coefficients. - f: normalized frequency (cutoff frequency / sample rate), must be less than 0.5 - Q: quality factor - V: gain - */ - void setParameters(Type type, float f, float Q, float V) { - float K = std::tan(M_PI * f); - switch (type) { - case LOWPASS: { - float norm = 1.f / (1.f + K / Q + K * K); - this->b[0] = K * K * norm; - this->b[1] = 2.f * this->b[0]; - this->b[2] = this->b[0]; - this->a[1] = 2.f * (K * K - 1.f) * norm; - this->a[2] = (1.f - K / Q + K * K) * norm; - } break; - - case HIGHPASS: { - float norm = 1.f / (1.f + K / Q + K * K); - this->b[0] = norm; - this->b[1] = -2.f * this->b[0]; - this->b[2] = this->b[0]; - this->a[1] = 2.f * (K * K - 1.f) * norm; - this->a[2] = (1.f - K / Q + K * K) * norm; - - } break; - - case LOWSHELF: { - float sqrtV = std::sqrt(V); - if (V >= 1.f) { - float norm = 1.f / (1.f + M_SQRT2 * K + K * K); - this->b[0] = (1.f + M_SQRT2 * sqrtV * K + V * K * K) * norm; - this->b[1] = 2.f * (V * K * K - 1.f) * norm; - this->b[2] = (1.f - M_SQRT2 * sqrtV * K + V * K * K) * norm; - this->a[1] = 2.f * (K * K - 1.f) * norm; - this->a[2] = (1.f - M_SQRT2 * K + K * K) * norm; - } - else { - float norm = 1.f / (1.f + M_SQRT2 / sqrtV * K + K * K / V); - this->b[0] = (1.f + M_SQRT2 * K + K * K) * norm; - this->b[1] = 2.f * (K * K - 1) * norm; - this->b[2] = (1.f - M_SQRT2 * K + K * K) * norm; - this->a[1] = 2.f * (K * K / V - 1.f) * norm; - this->a[2] = (1.f - M_SQRT2 / sqrtV * K + K * K / V) * norm; - } - } break; - - case HIGHSHELF: { - float sqrtV = std::sqrt(V); - if (V >= 1.f) { - float norm = 1.f / (1.f + M_SQRT2 * K + K * K); - this->b[0] = (V + M_SQRT2 * sqrtV * K + K * K) * norm; - this->b[1] = 2.f * (K * K - V) * norm; - this->b[2] = (V - M_SQRT2 * sqrtV * K + K * K) * norm; - this->a[1] = 2.f * (K * K - 1.f) * norm; - this->a[2] = (1.f - M_SQRT2 * K + K * K) * norm; - } - else { - float norm = 1.f / (1.f / V + M_SQRT2 / sqrtV * K + K * K); - this->b[0] = (1.f + M_SQRT2 * K + K * K) * norm; - this->b[1] = 2.f * (K * K - 1.f) * norm; - this->b[2] = (1.f - M_SQRT2 * K + K * K) * norm; - this->a[1] = 2.f * (K * K - 1.f / V) * norm; - this->a[2] = (1.f / V - M_SQRT2 / sqrtV * K + K * K) * norm; - } - } break; - - case BANDPASS: { - float norm = 1.f / (1.f + K / Q + K * K); - this->b[0] = K / Q * norm; - this->b[1] = 0.f; - this->b[2] = -this->b[0]; - this->a[1] = 2.f * (K * K - 1.f) * norm; - this->a[2] = (1.f - K / Q + K * K) * norm; - } break; - - case PEAK: { - float c = 1.0f / K; - float phi = c * c; - float Knum = c / Q; - float Kdenom = Knum; - - if(V > 1.0f) - Knum *= V; - else - Kdenom /= V; - - float norm = phi + Kdenom + 1.0; - this->b[0] = (phi + Knum + 1.0f) / norm; - this->b[1] = 2.0f * (1.0f - phi) / norm; - this->b[2] = (phi - Knum + 1.0f) / norm; - this->a[1] = 2.0f * (1.0f - phi) / norm; - this->a[2] = (phi - Kdenom + 1.0f) / norm; - } break; - - case NOTCH: { - float norm = 1.f / (1.f + K / Q + K * K); - this->b[0] = (1.f + K * K) * norm; - this->b[1] = 2.f * (K * K - 1.f) * norm; - this->b[2] = this->b[0]; - this->a[1] = this->b[1]; - this->a[2] = (1.f - K / Q + K * K) * norm; - } break; - - default: break; - } - } -}; - -typedef TBiquadFilter<> BiquadFilter; - - -/** - High-order filter to be used for anti-aliasing or anti-imaging. - The template parameter N should be 1/2 the desired filter order. - - Currently uses an 2*N-th order Butterworth filter. - @TODO: implement Chebyshev, Elliptic filter options. -*/ -template -class AAFilter -{ -public: - AAFilter() = default; - - /** Calculate Q values for a Butterworth filter of a given order */ - static std::vector calculateButterQs(int order) { - const int lim = int (order / 2); - std::vector Qs; - - for(int k = 1; k <= lim; ++k) { - auto b = -2.0f * std::cos((2.0f * k + order - 1) * 3.14159 / (2.0f * order)); - Qs.push_back(1.0f / b); - } - - std::reverse(Qs.begin(), Qs.end()); - return Qs; - } - - /** - * Resets the filter to process at a new sample rate. - * - * @param sampleRate: The base (i.e. pre-oversampling) sample rate of the audio being processed - * @param osRatio: The oversampling ratio at which the filter is being used - */ - void reset(float sampleRate, int osRatio) { - float fc = 0.98f * (sampleRate / 2.0f); - auto Qs = calculateButterQs(2*N); - - for(int i = 0; i < N; ++i) - filters[i].setParameters(BiquadFilter::Type::LOWPASS, fc / (osRatio * sampleRate), Qs[i], 1.0f); - } - - inline float process(float x) noexcept { - for(int i = 0; i < N; ++i) - x = filters[i].process(x); - - return x; - } - private: - BiquadFilter filters[N]; + dsp::BiquadFilter filters[N]; }; /** * Base class for oversampling of any order - */ + */ class BaseOversampling { public: - BaseOversampling() = default; - virtual ~BaseOversampling() {} + BaseOversampling() = default; + virtual ~BaseOversampling() {} - /** Resets the oversampler for processing at some base sample rate */ - virtual void reset(float /*baseSampleRate*/) = 0; + /** Resets the oversampler for processing at some base sample rate */ + virtual void reset(float /*baseSampleRate*/) = 0; - /** Upsample a single input sample and update the oversampled buffer */ - virtual void upsample(float) noexcept = 0; + /** Upsample a single input sample and update the oversampled buffer */ + virtual void upsample(float) noexcept = 0; - /** Output a downsampled output sample from the current oversampled buffer */ - virtual float downsample() noexcept = 0; + /** Output a downsampled output sample from the current oversampled buffer */ + virtual float downsample() noexcept = 0; - /** Returns a pointer to the oversampled buffer */ - virtual float* getOSBuffer() noexcept = 0; + /** Returns a pointer to the oversampled buffer */ + virtual float* getOSBuffer() noexcept = 0; }; -/** +/** Class to implement an oversampled process. To use, create an object and prepare using `reset()`. - + Then use the following code to process samples: @code oversample.upsample(x); @@ -341,51 +91,50 @@ public: @endcode */ template -class Oversampling : public BaseOversampling -{ +class Oversampling : public BaseOversampling { public: - Oversampling() = default; - virtual ~Oversampling() {} + Oversampling() = default; + virtual ~Oversampling() {} - void reset(float baseSampleRate) override { - aaFilter.reset(baseSampleRate, ratio); - aiFilter.reset(baseSampleRate, ratio); - std::fill(osBuffer, &osBuffer[ratio], 0.0f); - } + void reset(float baseSampleRate) override { + aaFilter.reset(baseSampleRate, ratio); + aiFilter.reset(baseSampleRate, ratio); + std::fill(osBuffer, &osBuffer[ratio], 0.0f); + } - inline void upsample(float x) noexcept override { - osBuffer[0] = ratio * x; - std::fill(&osBuffer[1], &osBuffer[ratio], 0.0f); + inline void upsample(float x) noexcept override { + osBuffer[0] = ratio * x; + std::fill(&osBuffer[1], &osBuffer[ratio], 0.0f); - for(int k = 0; k < ratio; k++) - osBuffer[k] = aiFilter.process(osBuffer[k]); - } + for (int k = 0; k < ratio; k++) + osBuffer[k] = aiFilter.process(osBuffer[k]); + } - inline float downsample() noexcept override { - float y = 0.0f; - for(int k = 0; k < ratio; k++) - y = aaFilter.process(osBuffer[k]); + inline float downsample() noexcept override { + float y = 0.0f; + for (int k = 0; k < ratio; k++) + y = aaFilter.process(osBuffer[k]); - return y; - } + return y; + } - inline float* getOSBuffer() noexcept override { - return osBuffer; - } + inline float* getOSBuffer() noexcept override { + return osBuffer; + } - float osBuffer[ratio]; + float osBuffer[ratio]; private: - AAFilter aaFilter; // anti-aliasing filter - AAFilter aiFilter; // anti-imaging filter + AAFilter aaFilter; // anti-aliasing filter + AAFilter aiFilter; // anti-imaging filter }; -/** +/** Class to implement an oversampled process, with variable oversampling factor. To use, create an object, set the oversampling factor using `setOversamplingindex()` and prepare using `reset()`. - + Then use the following code to process samples: @code oversample.upsample(x); @@ -398,45 +147,57 @@ private: template class VariableOversampling { public: - VariableOversampling() = default; + VariableOversampling() = default; - /** Prepare the oversampler to process audio at a given sample rate */ - void reset(float sampleRate) { - for(auto* os : oss) - os->reset(sampleRate); - } + /** Prepare the oversampler to process audio at a given sample rate */ + void reset(float sampleRate) { + for (auto* os : oss) + os->reset(sampleRate); + } - /** Sets the oversampling factor as 2^idx */ - void setOversamplingIndex (int newIdx) { osIdx = newIdx; } + /** Sets the oversampling factor as 2^idx */ + void setOversamplingIndex(int newIdx) { + osIdx = newIdx; + } + + /** Returns the oversampling index */ + int getOversamplingIndex() const noexcept { + return osIdx; + } - /** Returns the oversampling index */ - int getOversamplingIndex() const noexcept { return osIdx; } + /** Upsample a single input sample and update the oversampled buffer */ + inline void upsample(float x) noexcept { + oss[osIdx]->upsample(x); + } - /** Upsample a single input sample and update the oversampled buffer */ - inline void upsample(float x) noexcept { oss[osIdx]->upsample(x); } + /** Output a downsampled output sample from the current oversampled buffer */ + inline float downsample() noexcept { + return oss[osIdx]->downsample(); + } - /** Output a downsampled output sample from the current oversampled buffer */ - inline float downsample() noexcept { return oss[osIdx]->downsample(); } + /** Returns a pointer to the oversampled buffer */ + inline float* getOSBuffer() noexcept { + return oss[osIdx]->getOSBuffer(); + } - /** Returns a pointer to the oversampled buffer */ - inline float* getOSBuffer() noexcept { return oss[osIdx]->getOSBuffer(); } + /** Returns the current oversampling factor */ + int getOversamplingRatio() const noexcept { + return 1 << osIdx; + } - /** Returns the current oversampling factor */ - int getOversamplingRatio() const noexcept { return 1 << osIdx; } - private: - enum { - NumOS = 5, // number of oversampling options - }; - - int osIdx = 0; - - Oversampling<1 << 0, filtN> os0; // 1x - Oversampling<1 << 1, filtN> os1; // 2x - Oversampling<1 << 2, filtN> os2; // 4x - Oversampling<1 << 3, filtN> os3; // 8x - Oversampling<1 << 4, filtN> os4; // 16x - BaseOversampling* oss[NumOS] = { &os0, &os1, &os2, &os3, &os4 }; + enum { + NumOS = 5, // number of oversampling options + }; + + int osIdx = 0; + + Oversampling < 1 << 0, filtN > os0; // 1x + Oversampling < 1 << 1, filtN > os1; // 2x + Oversampling < 1 << 2, filtN > os2; // 4x + Oversampling < 1 << 3, filtN > os3; // 8x + Oversampling < 1 << 4, filtN > os4; // 16x + BaseOversampling* oss[NumOS] = { &os0, &os1, &os2, &os3, &os4 }; };