BogaudioModules

BogaudioModules for VCV Rack
Log | Files | Refs | README | LICENSE

commit 98cb5a4e4558754b748ee5cad0d351f960fd8514
parent 95b491117c96576a24150ef5463486f2368dd784
Author: Matt Demanett <matt@demanett.net>
Date:   Sat, 21 Apr 2018 20:23:18 -0400

Implement fast/clean CIC decimator.

Diffstat:
Mbenchmarks/filter_benchmark.cpp | 22++++++++++++++++++----
Msrc/FMOp.cpp | 2+-
Msrc/FMOp.hpp | 2+-
Msrc/Test.cpp | 22++++++++++++++++++++--
Msrc/Test.hpp | 29++++++++++++++++++++++++-----
Msrc/VCO.cpp | 6+++---
Msrc/VCO.hpp | 6+++---
Msrc/XCO.cpp | 8++++----
Msrc/XCO.hpp | 8++++----
Msrc/dsp/filter.cpp | 48++++++++++++++++++++++++++++++++++++++++++++----
Msrc/dsp/filter.hpp | 33+++++++++++++++++++++++++++++----
11 files changed, 151 insertions(+), 35 deletions(-)

diff --git a/benchmarks/filter_benchmark.cpp b/benchmarks/filter_benchmark.cpp @@ -7,19 +7,33 @@ using namespace bogaudio::dsp; -static void BM_Decimator(benchmark::State& state) { +static void BM_LPFDecimator(benchmark::State& state) { WhiteNoiseGenerator r; const int n = 8; float buf[n]; for (int i = 0; i < n; ++i) { buf[i] = r.next(); } - Decimator d(44100.0, n); + LPFDecimator d(44100.0, n); for (auto _ : state) { - benchmark::DoNotOptimize(d.next(n, buf)); + benchmark::DoNotOptimize(d.next(buf)); } } -BENCHMARK(BM_Decimator); +BENCHMARK(BM_LPFDecimator); + +static void BM_CICDecimator(benchmark::State& state) { + WhiteNoiseGenerator r; + const int n = 8; + float buf[n]; + for (int i = 0; i < n; ++i) { + buf[i] = r.next(); + } + CICDecimator d(4, 8); + for (auto _ : state) { + benchmark::DoNotOptimize(d.next(buf)); + } +} +BENCHMARK(BM_CICDecimator); static void BM_RackDecimator(benchmark::State& state) { WhiteNoiseGenerator r; diff --git a/src/FMOp.cpp b/src/FMOp.cpp @@ -116,7 +116,7 @@ void FMOp::step() { _phasor.advancePhase(); _buffer[i] = _sineTable.nextFromPhasor(_phasor, Phasor::radiansToPhase(offset)); } - float out = _decimator.next(oversample, _buffer); + float out = _decimator.next(_buffer); out *= _level; if (_levelEnvelopeOn) { out *= envelope; diff --git a/src/FMOp.hpp b/src/FMOp.hpp @@ -65,7 +65,7 @@ struct FMOp : Module { ADSR _envelope; Phasor _phasor; SineTableOscillator _sineTable; - Decimator _decimator; + LPFDecimator _decimator; SchmittTrigger _gateTrigger; FMOp() : Module(NUM_PARAMS, NUM_INPUTS, NUM_OUTPUTS, NUM_LIGHTS) { diff --git a/src/Test.cpp b/src/Test.cpp @@ -235,10 +235,10 @@ void Test::step() { squareBuf[i] = _square.nextFromPhasor(_oversampledPhasor); } - out += _sawDecimator.next(OVERSAMPLEN, sawBuf) * oMix; + out += _sawDecimator.next(sawBuf) * oMix; // out += _sawRackDecimator.process(sawBuf) * oMix; - out2 += _squareDecimator.next(OVERSAMPLEN, squareBuf) * oMix; + out2 += _squareDecimator.next(squareBuf) * oMix; // out2 += _squareRackDecimator.process(squareBuf) * oMix; } else { @@ -250,6 +250,24 @@ void Test::step() { outputs[OUT_OUTPUT].value = out * 5.0f; outputs[OUT2_OUTPUT].value = out2 * 5.0f; +#elif DECIMATORS + const int quality = 12; + float sampleRate = engineGetSampleRate(); + float frequency = oscillatorPitch(15000.0); + _saw.setSampleRate(sampleRate); + _saw.setFrequency(frequency / (float)OVERSAMPLEN); + _saw.setQuality(quality); + _cicDecimator.setParams(sampleRate, OVERSAMPLEN); + _lpfDecimator.setParams(sampleRate, OVERSAMPLEN); + + float buf[OVERSAMPLEN] {}; + for (int i = 0; i < OVERSAMPLEN; ++i) { + buf[i] = _saw.next(); + } + outputs[OUT_OUTPUT].value = _cicDecimator.next(buf) * 5.0f; + // outputs[OUT2_OUTPUT].value = _lpfDecimator.next(buf) * 5.0f; + outputs[OUT2_OUTPUT].value = _rackDecimator.process(buf) * 5.0f; + #elif FM const float amplitude = 5.0f; float baseHz = oscillatorPitch(); diff --git a/src/Test.hpp b/src/Test.hpp @@ -16,9 +16,10 @@ extern Model* modelTest; // #define OVERSAMPLING 1 // #define OVERSAMPLED_BL 1 // #define ANTIALIASING 1 +#define DECIMATORS 1 // #define FM 1 // #define PM 1 -#define FEEDBACK_PM 1 +// #define FEEDBACK_PM 1 // #define EG 1 // #define TABLES 1 @@ -54,6 +55,10 @@ extern Model* modelTest; #include "dsp/oscillator.hpp" #include "dsp/decimator.hpp" // rack #include "dsp/filter.hpp" +#elif DECIMATORS +#include "dsp/oscillator.hpp" +#include "dsp/filter.hpp" +#include "dsp/decimator.hpp" // rack #elif FM #include "dsp/oscillator.hpp" #elif PM @@ -105,7 +110,8 @@ struct Test : Module { LowPassFilter _lpf; #elif SINE SineOscillator _sine; - SineTableOscillator _sine2; + SineTable _table; + TablePhasor _sine2; #elif SQUARE SquareOscillator _square; BandLimitedSquareOscillator _square2; @@ -141,10 +147,17 @@ struct Test : Module { Phasor _oversampledPhasor; BandLimitedSawOscillator _saw; BandLimitedSquareOscillator _square; - bogaudio::dsp::Decimator _sawDecimator; - bogaudio::dsp::Decimator _squareDecimator; + bogaudio::dsp::LPFDecimator _sawDecimator; + bogaudio::dsp::LPFDecimator _squareDecimator; rack::Decimator<OVERSAMPLEN, OVERSAMPLEN> _sawRackDecimator; rack::Decimator<OVERSAMPLEN, OVERSAMPLEN> _squareRackDecimator; +#elif DECIMATORS + #define OVERSAMPLEN 8 + #define STAGES 4 + BandLimitedSawOscillator _saw; + bogaudio::dsp::CICDecimator _cicDecimator; + bogaudio::dsp::LPFDecimator _lpfDecimator; + rack::Decimator<OVERSAMPLEN, OVERSAMPLEN> _rackDecimator; #elif FM float _baseHz = 0.0f; float _ratio = 0.0f; @@ -169,13 +182,19 @@ struct Test : Module { Test() : Module(NUM_PARAMS, NUM_INPUTS, NUM_OUTPUTS, NUM_LIGHTS) -#if TABLES +#if SINE + , _table(12) + , _sine2(_table) +#elif DECIMATORS + , _cicDecimator(STAGES) +#elif TABLES , _table(StaticBlepTable::table(), 44100.0, 1000.0) #endif { onReset(); #ifdef SINE + _table.generate(); _sine2.setPhase(M_PI); #elif SAW diff --git a/src/VCO.cpp b/src/VCO.cpp @@ -102,13 +102,13 @@ void VCO::step() { } } if (outputs[SQUARE_OUTPUT].active) { - squareOut += oMix * amplitude * _squareDecimator.next(oversample, _squareBuffer); + squareOut += oMix * amplitude * _squareDecimator.next(_squareBuffer); } if (outputs[SAW_OUTPUT].active) { - sawOut += oMix * amplitude * _sawDecimator.next(oversample, _sawBuffer); + sawOut += oMix * amplitude * _sawDecimator.next(_sawBuffer); } if (outputs[TRIANGLE_OUTPUT].active) { - triangleOut += oMix * amplitude * _triangleDecimator.next(oversample, _triangleBuffer); + triangleOut += oMix * amplitude * _triangleDecimator.next(_triangleBuffer); } } else { diff --git a/src/VCO.hpp b/src/VCO.hpp @@ -61,9 +61,9 @@ struct VCO : Module { BandLimitedSawOscillator _saw; TriangleOscillator _triangle; SineTableOscillator _sine; - Decimator _squareDecimator; - Decimator _sawDecimator; - Decimator _triangleDecimator; + CICDecimator _squareDecimator; + CICDecimator _sawDecimator; + CICDecimator _triangleDecimator; float _squareBuffer[oversample]; float _sawBuffer[oversample]; float _triangleBuffer[oversample]; diff --git a/src/XCO.cpp b/src/XCO.cpp @@ -153,19 +153,19 @@ void XCO::step() { } } if (squareOversample) { - squareOut += oMix * amplitude * _squareDecimator.next(oversample, _squareBuffer); + squareOut += oMix * amplitude * _squareDecimator.next(_squareBuffer); } if (sawOversample) { - sawOut += oMix * amplitude * _sawDecimator.next(oversample, _sawBuffer); + sawOut += oMix * amplitude * _sawDecimator.next(_sawBuffer); } if (triangleOversample) { - triangleOut += amplitude * _triangleDecimator.next(oversample, _triangleBuffer); + triangleOut += amplitude * _triangleDecimator.next(_triangleBuffer); if (!triangleSample) { triangleOut *= oMix; } } if (sineOversample) { - sineOut += amplitude * _sineDecimator.next(oversample, _sineBuffer); + sineOut += amplitude * _sineDecimator.next(_sineBuffer); } } else { diff --git a/src/XCO.hpp b/src/XCO.hpp @@ -96,10 +96,10 @@ struct XCO : Module { BandLimitedSawOscillator _saw; TriangleOscillator _triangle; SineTableOscillator _sine; - Decimator _squareDecimator; - Decimator _sawDecimator; - Decimator _triangleDecimator; - Decimator _sineDecimator; + CICDecimator _squareDecimator; + CICDecimator _sawDecimator; + CICDecimator _triangleDecimator; + CICDecimator _sineDecimator; float _squareBuffer[oversample]; float _sawBuffer[oversample]; float _triangleBuffer[oversample]; diff --git a/src/dsp/filter.cpp b/src/dsp/filter.cpp @@ -177,20 +177,60 @@ float MultipoleFilter::next(float sample) { } -void Decimator::setParams(float sampleRate, int oversample) { +void LPFDecimator::setParams(float sampleRate, int factor) { + _factor = factor; _filter.setParams( MultipoleFilter::LP_TYPE, 4, - oversample * sampleRate, + factor * sampleRate, 0.45f * sampleRate, 0 ); } -float Decimator::next(int n, const float* buf) { +float LPFDecimator::next(const float* buf) { float s = 0.0f; - for (int i = 0; i < n; ++i) { + for (int i = 0; i < _factor; ++i) { s = _filter.next(buf[i]); } return s; } + + +// cascaded integrator–comb decimator: https://en.wikipedia.org/wiki/Cascaded_integrator%E2%80%93comb_filter +CICDecimator::CICDecimator(int stages, int factor) { + assert(stages > 0); + _stages = stages; + _integrators = new T[_stages + 1] {}; + _combs = new T[_stages] {}; + setParams(0.0f, factor); +} + +CICDecimator::~CICDecimator() { + delete[] _integrators; + delete[] _combs; +} + +void CICDecimator::setParams(float _sampleRate, int factor) { + assert(factor > 0); + if (_factor != factor) { + _factor = factor; + _gainCorrection = 1.0f / (float)(pow(_factor, _stages)); + } +} + +float CICDecimator::next(const float* buf) { + for (int i = 0; i < _factor; ++i) { + _integrators[0] = buf[i] * scale; + for (int j = 1; j <= _stages; ++j) { + _integrators[j] += _integrators[j - 1]; + } + } + T s = _integrators[_stages]; + for (int i = 0; i < _stages; ++i) { + T t = s; + s -= _combs[i]; + _combs[i] = t; + } + return _gainCorrection * (s / (float)scale); +} diff --git a/src/dsp/filter.hpp b/src/dsp/filter.hpp @@ -1,5 +1,6 @@ #pragma once +#include <stdint.h> #include <math.h> #include "buffer.hpp" @@ -127,14 +128,38 @@ struct MultipoleFilter : Filter { }; struct Decimator { + Decimator() {} + virtual ~Decimator() {} + + virtual void setParams(float sampleRate, int factor) = 0; + virtual float next(const float* buf) = 0; +}; + +struct LPFDecimator : Decimator { + int _factor; MultipoleFilter _filter; - Decimator(float sampleRate = 1000.0f, int oversample = 4) { - setParams(sampleRate, oversample); + LPFDecimator(float sampleRate = 1000.0f, int factor = 8) { + setParams(sampleRate, factor); } + virtual void setParams(float sampleRate, int factor) override; + virtual float next(const float* buf) override; +}; + +struct CICDecimator : Decimator { + typedef int64_t T; + static constexpr T scale = 1l << 32; + int _stages; + T* _integrators; + T* _combs; + int _factor = 0; + float _gainCorrection; + + CICDecimator(int stages = 4, int factor = 8); + virtual ~CICDecimator(); - void setParams(float sampleRate, int oversample); - float next(int n, const float* buf); + virtual void setParams(float sampleRate, int factor) override; + virtual float next(const float* buf) override; }; } // namespace dsp