From e85dbd43b845d913be8adf006b8083bf388ab502 Mon Sep 17 00:00:00 2001 From: Holger Wirtz Date: Fri, 21 Jun 2019 08:39:39 +0200 Subject: [PATCH] Removed Splie interpolation. Tried another calculation for filter coefficients. --- config.h | 2 +- effect_modulated_delay.cpp | 116 +++++++++++++---------------------- effect_modulated_delay.h | 7 ++- spline.cpp | 120 ------------------------------------- spline.h | 42 ------------- 5 files changed, 47 insertions(+), 240 deletions(-) delete mode 100644 spline.cpp delete mode 100644 spline.h diff --git a/config.h b/config.h index dd3383f..5cb6ea7 100644 --- a/config.h +++ b/config.h @@ -57,7 +57,7 @@ #define USE_XFADE_DATA 1 // CHORUS parameters #define INTERPOLATION_WINDOW_SIZE 5 // use only odd numbers!!! -//#define INTERPOLATE_MODE 11 +#define INTERPOLATE_MODE 11 #define CHORUS_WAVEFORM WAVEFORM_TRIANGLE // WAVEFORM_SINE WAVEFORM_TRIANGLE WAVEFORM_SAWTOOTH WAVEFORM_SAWTOOTH_REVERSE #define CHORUS_DELAY_LENGTH_SAMPLES (15*AUDIO_BLOCK_SAMPLES) // one AUDIO_BLOCK_SAMPLES = 2.902ms; you need doubled length, e.g. delay point is 20ms, so you need up to 40ms delay! diff --git a/effect_modulated_delay.cpp b/effect_modulated_delay.cpp index 4e4df89..f1cad50 100644 --- a/effect_modulated_delay.cpp +++ b/effect_modulated_delay.cpp @@ -25,7 +25,6 @@ #include #include "limits.h" #include "effect_modulated_delay.h" -#include "spline.h" #include "config.h" /******************************************************************/ @@ -62,7 +61,7 @@ boolean AudioEffectModulatedDelay::begin(short *delayline, int d_length) filter.numStages = 1; filter.pState = filter_state; filter.pCoeffs = filter_coeffs; - calcModFilterCoeff(5000.0, 0.0, 5.0); + calcModFilterCoeff(500.0); return (true); } @@ -86,13 +85,6 @@ void AudioEffectModulatedDelay::update(void) float mod_number; float mod_fraction; -#ifdef INTERPOLATE_MODE - int8_t j; - float x[INTERPOLATION_WINDOW_SIZE]; - float y[INTERPOLATION_WINDOW_SIZE]; - Spline s(x, y, INTERPOLATION_WINDOW_SIZE, INTERPOLATE_MODE); -#endif - // (Filter implementation: https://web.fhnw.ch/technik/projekte/eit/Fruehling2016/MuelZum/html/parametric__equalizer__example_8c_source.html) arm_q15_to_float(modulation->data, modulation_f32, AUDIO_BLOCK_SAMPLES); arm_biquad_cascade_df1_f32(&filter, modulation_f32, modulation_f32, AUDIO_BLOCK_SAMPLES); @@ -111,23 +103,8 @@ void AudioEffectModulatedDelay::update(void) mod_idx = *mp * float(_delay_length >> 1); mod_fraction = modff(mod_idx, &mod_number); -#ifdef INTERPOLATE_MODE - // Generate a an array with the size of INTERPOLATION_WINDOW_SIZE of x/y values around mod_idx for interpolation - uint8_t c; - int16_t c_mod_idx = _circ_idx - int(round(mod_idx)); // This is the pointer to the value in the circular buffer at the current modulation index - for (j = ~(INTERPOLATION_WINDOW_SIZE >> 1) | 0x01, c = 0; j <= INTERPOLATION_WINDOW_SIZE >> 1; j++, c++) // only another way to say: from -INTERPOLATION_WINDOW_SIZE/2 to INTERPOLATION_WINDOW_SIZE/2 - { - int16_t jc_mod_idx = (c_mod_idx + j) % _delay_length; // The modulation index pointer plus the value of the current window pointer - if (jc_mod_idx < 0) - y[c] = float(_delayline[_delay_length + jc_mod_idx]); - else - y[c] = float(_delayline[jc_mod_idx]); - x[c] = float(j); - } - *bp = int(round(s.value(mod_fraction))); -#else // Simple interpolation - int16_t c_mod_idx = (_circ_idx - int(round(mod_idx))) % _delay_length; + int16_t c_mod_idx = (_circ_idx + int(round(mod_idx))) % _delay_length; float value1, value2; if (c_mod_idx < 0) @@ -140,7 +117,15 @@ void AudioEffectModulatedDelay::update(void) value1 = _delayline[c_mod_idx - 1]; value2 = _delayline[c_mod_idx]; } - *bp = mod_fraction * value1 + (1.0 - mod_fraction) * value2; + *bp = int(round(mod_fraction * value1 + (1.0 - mod_fraction) * value2)); + +#ifdef DEBUG + float m = (value2 - value1) / (SHRT_MAX >> 1); + if (m > 1.0 || m < -1.0) + { + Serial.print(F("WARNING m=")); + Serial.println(m, 4); + } #endif bp++; // next audio data @@ -164,56 +149,39 @@ void AudioEffectModulatedDelay::setDelay(float milliseconds) _delay_length = min(AUDIO_SAMPLE_RATE * milliseconds / 500, _max_delay_length); } -void AudioEffectModulatedDelay::calcModFilterCoeff(float32_t cFrq, float32_t gain, float32_t width) +void AudioEffectModulatedDelay::calcModFilterCoeff(float32_t cFrq) { - /* Calculate intermediate values */ - // float32_t A = sqrt(pow(10, gain / 20.0f)); - // float32_t w0 = 2.0f * PI * cFrq / ((float32_t)AUDIO_SAMPLE_RATE_EXACT); - // float32_t cosw0 = cos(w0); - // float32_t sinw0 = sin(w0); - // float32_t alpha = sinw0 / (2.0f * width); - - /* Calculate coefficients */ - // float32_t b0 = 1.0f + alpha * A; - // float32_t b1 = -2.0f * cosw0; - // float32_t b2 = 1.0f - alpha * A; - // float32_t a0 = 1.0f + alpha / A; - // float32_t a1 = -2.0f * cosw0; - // float32_t a2 = 1.0f - alpha / A; - -/* https://stackoverflow.com/questions/20924868/calculate-coefficients-of-2nd-order-butterworth-low-pass-filter/20932062 - ff=cutoff_frq/sample_rate=AUDIO_SAMPLE_RATE_EXACT/1000 - const double ita =1.0/ tan(M_PI*ff); - const double q=sqrt(2.0); - b0 = 1.0 / (1.0 + q*ita + ita*ita); - b1= 2*b0; - b2= b0; - a1 = 2.0 * (ita*ita - 1.0) * b0; - a2 = -(1.0 - q*ita + ita*ita) * b0; - - ff=1000/44117.64706=0.02266666666 - ita=804.60898525 - q=1.414213 - - */ - // 1kHz 2nd order Butterworth lowpass filter coefficients - // calculated with Iowa IIR FIlter Designer 6.5 - float32_t b0 = 0.124589380980617656; - float32_t b1 = 0.124589380980617656; - float32_t b2 = 0.0; - float32_t a0 = 1.000000000000000000; - float32_t a1 = -0.750821238038764660; - float32_t a2 = 0.0; - - /* Normalize so a0 = 1 */ - filter_coeffs[0] = b0 / a0; - filter_coeffs[1] = b1 / a0; - filter_coeffs[2] = b2 / a0; - filter_coeffs[3] = -a1 / a0; - filter_coeffs[4] = -a2 / a0; + const float sqrt2 = 1.4142135623730950488; + + float QcRaw = (2 * PI * cFrq) / AUDIO_SAMPLE_RATE_EXACT; // Find cutoff frequency in [0..PI] + float QcWarp = tan(QcRaw); // Warp cutoff frequency + + float gain = 1 / (1 + sqrt2 / QcWarp + 2 / (QcWarp * QcWarp)); + filter_coeffs[2] = (1 - sqrt2 / QcWarp + 2 / (QcWarp * QcWarp)) * gain; + filter_coeffs[1] = (2 - 2 * 2 / (QcWarp * QcWarp)) * gain; + filter_coeffs[0] = 1; + filter_coeffs[3] = 1 * gain; + filter_coeffs[4] = 2 * gain; + + /* + // 1.1kHz 2nd order Butterworth lowpass filter coefficients + // calculated with Iowa IIR FIlter Designer 6.5 + float32_t b0 = 0.072959657268266670; + float32_t b1 = 0.072959657268266670; + float32_t b2 = 0.0; + float32_t a0 = 1.000000000000000000; + float32_t a1 = -0.854080685463466605; + float32_t a2 = 0.0; + + // Normalize so a0 = 1 + filter_coeffs[0] = b0 / a0; + filter_coeffs[1] = b1 / a0; + filter_coeffs[2] = b2 / a0; + filter_coeffs[3] = -a1 / a0; + filter_coeffs[4] = -a2 / a0; */ } -void AudioEffectModulatedDelay::setModFilter(float cFrq, float gain, float width) +void AudioEffectModulatedDelay::setModFilter(float cFrq) { - calcModFilterCoeff(cFrq, gain, width); + calcModFilterCoeff(cFrq); } diff --git a/effect_modulated_delay.h b/effect_modulated_delay.h index 450a14a..c4a2373 100644 --- a/effect_modulated_delay.h +++ b/effect_modulated_delay.h @@ -26,9 +26,9 @@ #include "Arduino.h" #include "AudioStream.h" +#include "config.h" /*************************************************************************/ - // A u d i o E f f e c t M o d u l a t e d D e l a y // Written by Pete (El Supremo) Jan 2014 // 140219 - correct storage class (not static) @@ -45,16 +45,17 @@ class AudioEffectModulatedDelay : boolean begin(short *delayline, int delay_length); virtual void update(void); virtual void setDelay(float milliseconds); - virtual void setModFilter(float cFrq, float gain, float width); + virtual void setModFilter(float cFrq); private: - virtual void calcModFilterCoeff(float32_t cFrq, float32_t gain, float32_t width); + virtual void calcModFilterCoeff(float32_t cFrq); audio_block_t *inputQueueArray[2]; int16_t *_delayline; uint16_t _circ_idx; uint16_t _max_delay_length; uint16_t _delay_length; + // filter data arm_biquad_casd_df1_inst_f32 filter; float32_t modulation_f32[AUDIO_BLOCK_SAMPLES]; diff --git a/spline.cpp b/spline.cpp deleted file mode 100644 index 47e1599..0000000 --- a/spline.cpp +++ /dev/null @@ -1,120 +0,0 @@ -#include "Arduino.h" -#include "spline.h" -#include - -Spline::Spline(void) { - _prev_point = 0; -} - -Spline::Spline( float x[], float y[], int numPoints, int degree ) -{ - setPoints(x, y, numPoints); - setDegree(degree); - _prev_point = 0; -} - -Spline::Spline( float x[], float y[], float m[], int numPoints ) -{ - setPoints(x, y, m, numPoints); - setDegree(Hermite); - _prev_point = 0; -} - -void Spline::setPoints( float x[], float y[], int numPoints ) { - _x = x; - _y = y; - _length = numPoints; -} - -void Spline::setPoints( float x[], float y[], float m[], int numPoints ) { - _x = x; - _y = y; - _m = m; - _length = numPoints; -} - -void Spline::setDegree( int degree ) { - _degree = degree; -} - -float Spline::value( float x ) -{ - if ( _x[0] > x ) { - return _y[0]; - } - else if ( _x[_length - 1] < x ) { - return _y[_length - 1]; - } - else { - for (int i = 0; i < _length; i++ ) - { - int index = ( i + _prev_point ) % _length; - - if ( _x[index] == x ) { - _prev_point = index; - return _y[index]; - } else if ( (_x[index] < x) && (x < _x[index + 1]) ) { - _prev_point = index; - return calc( x, index ); - } - } - } - return (0.0); -} - -float Spline::calc( float x, int i ) -{ - switch ( _degree ) { - case 0: - return _y[i]; - case 1: - if ( _x[i] == _x[i + 1] ) { - // Avoids division by 0 - return _y[i]; - } else { - return _y[i] + (_y[i + 1] - _y[i]) * ( x - _x[i]) / ( _x[i + 1] - _x[i] ); - } - case Hermite: - return hermite( ((x - _x[i]) / (_x[i + 1] - _x[i])), _y[i], _y[i + 1], _m[i], _m[i + 1], _x[i], _x[i + 1] ); - case Catmull: - if ( i == 0 ) { - // x prior to spline start - first point used to determine tangent - return _y[1]; - } else if ( i == _length - 2 ) { - // x after spline end - last point used to determine tangent - return _y[_length - 2]; - } else { - float t = (x - _x[i]) / (_x[i + 1] - _x[i]); - float m0 = (i == 0 ? 0 : catmull_tangent(i) ); - float m1 = (i == _length - 1 ? 0 : catmull_tangent(i + 1) ); - return hermite( t, _y[i], _y[i + 1], m0, m1, _x[i], _x[i + 1]); - } - } - return (0.0); -} - -float Spline::hermite( float t, float p0, float p1, float m0, float m1, float x0, float x1 ) { - return (hermite_00(t) * p0) + (hermite_10(t) * (x1 - x0) * m0) + (hermite_01(t) * p1) + (hermite_11(t) * (x1 - x0) * m1); -} -float Spline::hermite_00( float t ) { - return (2 * pow(t, 3)) - (3 * pow(t, 2)) + 1; -} -float Spline::hermite_10( float t ) { - return pow(t, 3) - (2 * pow(t, 2)) + t; -} -float Spline::hermite_01( float t ) { - return (3 * pow(t, 2)) - (2 * pow(t, 3)); -} -float Spline::hermite_11( float t ) { - return pow(t, 3) - pow(t, 2); -} - -float Spline::catmull_tangent( int i ) -{ - if ( _x[i + 1] == _x[i - 1] ) { - // Avoids division by 0 - return 0; - } else { - return (_y[i + 1] - _y[i - 1]) / (_x[i + 1] - _x[i - 1]); - } -} diff --git a/spline.h b/spline.h deleted file mode 100644 index 563c742..0000000 --- a/spline.h +++ /dev/null @@ -1,42 +0,0 @@ -/* - Library for 1-d splines - Copyright Ryan Michael - Licensed under the GPLv3 -*/ -#ifndef spline_h -#define spline_h - -#include "Arduino.h" - -#define Hermite 10 -#define Catmull 11 - -class Spline -{ - public: - Spline( void ); - Spline( float x[], float y[], int numPoints, int degree = 1 ); - Spline( float x[], float y[], float m[], int numPoints ); - float value( float x ); - void setPoints( float x[], float y[], int numPoints ); - void setPoints( float x[], float y[], float m[], int numPoints ); - void setDegree( int degree ); - - private: - float calc( float, int); - float* _x; - float* _y; - float* _m; - int _degree; - int _length; - int _prev_point; - - float hermite( float t, float p0, float p1, float m0, float m1, float x0, float x1 ); - float hermite_00( float t ); - float hermite_10( float t ); - float hermite_01( float t ); - float hermite_11( float t ); - float catmull_tangent( int i ); -}; - -#endif