1162 lines
38 KiB
C++
1162 lines
38 KiB
C++
/*
|
|
* ysfv.h
|
|
* Copyright (c) 2016 Chris Johnson
|
|
* Copyright (c) 2025 Christopher Herb
|
|
* Based on work:
|
|
* - YLowpass by Chris Johnson, 2016
|
|
* - YHighpass by Chris Johnson, 2016
|
|
* - YBandpass by Chris Johnson, 2016
|
|
* - YNotch by Chris Johnson, 2016
|
|
* This file is a derivative/major refactor consolidating the above modules.
|
|
*
|
|
* Permission is hereby granted, free of charge, to any person obtaining a copy
|
|
* of this software and associated documentation files (the "Software"), to deal
|
|
* in the Software without restriction, including without limitation the rights
|
|
* to use, copy, modify, merge, publish, distribute, sublicense, and/or sell
|
|
* copies of the Software, and to permit persons to whom the Software is
|
|
* furnished to do so, subject to the following conditions:
|
|
*
|
|
* The above copyright notice and this permission notice shall be included in
|
|
* all copies or substantial portions of the Software.
|
|
*
|
|
* THE SOFTWARE IS PROVIDED "AS IS", WITHOUT WARRANTY OF ANY KIND, EXPRESS OR
|
|
* IMPLIED, INCLUDING BUT NOT LIMITED TO THE WARRANTIES OF MERCHANTABILITY,
|
|
* FITNESS FOR A PARTICULAR PURPOSE AND NONINFRINGEMENT. IN NO EVENT SHALL THE
|
|
* AUTHORS OR COPYRIGHT HOLDERS BE LIABLE FOR ANY CLAIM, DAMAGES OR OTHER
|
|
* LIABILITY, WHETHER IN AN ACTION OF CONTRACT, TORT OR OTHERWISE, ARISING FROM,
|
|
* OUT OF OR IN CONNECTION WITH THE SOFTWARE OR THE USE OR OTHER DEALINGS IN THE
|
|
* SOFTWARE.
|
|
*
|
|
* Changes:
|
|
* - 2025-11-06 Christopher Herb:
|
|
* - Consolidated YLowpass, YHighpass, YBandpass, YNotch into one (state-variable) filter
|
|
* - Templated audio buffer i/o
|
|
* - Converted to procedural programming style
|
|
*/
|
|
|
|
#pragma once
|
|
|
|
#define _USE_MATH_DEFINES
|
|
#include <array>
|
|
#include <math.h>
|
|
#include <stdint.h>
|
|
|
|
using namespace std;
|
|
|
|
namespace trnr {
|
|
|
|
////////////
|
|
// COMMON //
|
|
////////////
|
|
|
|
enum {
|
|
Y_BIQ_FREQ,
|
|
Y_BIQ_RESO,
|
|
Y_BIQ_A0,
|
|
Y_BIQ_A1,
|
|
Y_BIQ_A2,
|
|
Y_BIQ_B1,
|
|
Y_BIQ_B2,
|
|
Y_BIQ_A_A0,
|
|
Y_BIQ_A_A1,
|
|
Y_BIQ_A_A2,
|
|
Y_BIQ_B_A1,
|
|
Y_BIQ_B_A2,
|
|
Y_BIQ_A_B0,
|
|
Y_BIQ_A_B1,
|
|
Y_BIQ_A_B2,
|
|
Y_BIQ_B_B1,
|
|
Y_BIQ_B_B2,
|
|
Y_BIQ_S_L1,
|
|
Y_BIQ_S_L2,
|
|
Y_BIQ_S_R1,
|
|
Y_BIQ_S_R2,
|
|
Y_BIQ_TOTAL
|
|
};
|
|
|
|
enum {
|
|
Y_FIX_FREQ,
|
|
Y_FIX_RESO,
|
|
Y_FIX_A0,
|
|
Y_FIX_A1,
|
|
Y_FIX_A2,
|
|
Y_FIX_B1,
|
|
Y_FIX_B2,
|
|
Y_FIX_S_L1,
|
|
Y_FIX_S_L2,
|
|
Y_FIX_S_R1,
|
|
Y_FIX_S_R2,
|
|
Y_FIX_TOTAL
|
|
};
|
|
|
|
/////////////
|
|
// LOWPASS //
|
|
/////////////
|
|
|
|
struct ylowpass {
|
|
double samplerate;
|
|
array<double, Y_BIQ_TOTAL> biquad;
|
|
|
|
double powFactorA;
|
|
double powFactorB;
|
|
double inTrimA;
|
|
double inTrimB;
|
|
double outTrimA;
|
|
double outTrimB;
|
|
|
|
array<double, Y_FIX_TOTAL> fixA;
|
|
array<double, Y_FIX_TOTAL> fixB;
|
|
|
|
uint32_t fpdL;
|
|
uint32_t fpdR;
|
|
// default stuff
|
|
|
|
float drive;
|
|
float frequency;
|
|
float resonance;
|
|
float edge;
|
|
float output;
|
|
float mix; // parameters. Always 0-1, and we scale/alter them elsewhere.
|
|
};
|
|
|
|
inline void ylowpass_init(ylowpass& y, double samplerate)
|
|
{
|
|
y.samplerate = samplerate;
|
|
y.drive = 0.1f;
|
|
y.frequency = 1.0f;
|
|
y.resonance = 0.0f;
|
|
y.edge = 0.1f;
|
|
y.output = 0.9f;
|
|
y.mix = 1.0f;
|
|
|
|
for (int x = 0; x < Y_BIQ_TOTAL; x++) { y.biquad[x] = 0.0; }
|
|
y.powFactorA = 1.0;
|
|
y.powFactorB = 1.0;
|
|
y.inTrimA = 0.1;
|
|
y.inTrimB = 0.1;
|
|
y.outTrimA = 1.0;
|
|
y.outTrimB = 1.0;
|
|
for (int x = 0; x < Y_FIX_TOTAL; x++) {
|
|
y.fixA[x] = 0.0;
|
|
y.fixB[x] = 0.0;
|
|
}
|
|
|
|
y.fpdL = 1.0;
|
|
while (y.fpdL < 16386) y.fpdL = rand() * UINT32_MAX;
|
|
y.fpdR = 1.0;
|
|
while (y.fpdR < 16386) y.fpdR = rand() * UINT32_MAX;
|
|
}
|
|
|
|
template <typename t_sample>
|
|
inline void ylowpass_process_block(ylowpass& y, t_sample** inputs, t_sample** outputs,
|
|
int blockSize)
|
|
{
|
|
t_sample* in1 = inputs[0];
|
|
t_sample* in2 = inputs[1];
|
|
t_sample* out1 = outputs[0];
|
|
t_sample* out2 = outputs[1];
|
|
|
|
int inFramesToProcess = blockSize;
|
|
double overallscale = 1.0;
|
|
overallscale /= 44100.0;
|
|
overallscale *= y.samplerate;
|
|
|
|
y.inTrimA = y.inTrimB;
|
|
y.inTrimB = y.drive * 10.0;
|
|
|
|
y.biquad[Y_BIQ_FREQ] = pow(y.frequency, 3) * 20000.0;
|
|
if (y.biquad[Y_BIQ_FREQ] < 15.0) y.biquad[Y_BIQ_FREQ] = 15.0;
|
|
y.biquad[Y_BIQ_FREQ] /= y.samplerate;
|
|
y.biquad[Y_BIQ_RESO] = (pow(y.resonance, 2) * 15.0) + 0.5571;
|
|
y.biquad[Y_BIQ_A_A0] = y.biquad[Y_BIQ_A_B0];
|
|
y.biquad[Y_BIQ_A_A1] = y.biquad[Y_BIQ_A_B1];
|
|
y.biquad[Y_BIQ_A_A2] = y.biquad[Y_BIQ_A_B2];
|
|
y.biquad[Y_BIQ_B_A1] = y.biquad[Y_BIQ_B_B1];
|
|
y.biquad[Y_BIQ_B_A2] = y.biquad[Y_BIQ_B_B2];
|
|
// previous run through the buffer is still in the filter, so we move it
|
|
// to the A section and now it's the new starting point.
|
|
double K = tan(M_PI * y.biquad[Y_BIQ_FREQ]);
|
|
double norm = 1.0 / (1.0 + K / y.biquad[Y_BIQ_RESO] + K * K);
|
|
y.biquad[Y_BIQ_A_B0] = K * K * norm;
|
|
y.biquad[Y_BIQ_A_B1] = 2.0 * y.biquad[Y_BIQ_A_B0];
|
|
y.biquad[Y_BIQ_A_B2] = y.biquad[Y_BIQ_A_B0];
|
|
y.biquad[Y_BIQ_B_B1] = 2.0 * (K * K - 1.0) * norm;
|
|
y.biquad[Y_BIQ_B_B2] = (1.0 - K / y.biquad[Y_BIQ_RESO] + K * K) * norm;
|
|
// for the coefficient-interpolated biquad filter
|
|
|
|
y.powFactorA = y.powFactorB;
|
|
y.powFactorB = pow(y.edge + 0.9, 4);
|
|
|
|
// 1.0 == target neutral
|
|
|
|
y.outTrimA = y.outTrimB;
|
|
y.outTrimB = y.output;
|
|
|
|
double wet = y.mix;
|
|
|
|
y.fixA[Y_FIX_FREQ] = y.fixB[Y_FIX_FREQ] = 20000.0 / y.samplerate;
|
|
y.fixA[Y_FIX_RESO] = y.fixB[Y_FIX_RESO] = 0.7071; // butterworth Q
|
|
|
|
K = tan(M_PI * y.fixA[Y_FIX_FREQ]);
|
|
norm = 1.0 / (1.0 + K / y.fixA[Y_FIX_RESO] + K * K);
|
|
y.fixA[Y_FIX_A0] = y.fixB[Y_FIX_A0] = K * K * norm;
|
|
y.fixA[Y_FIX_A1] = y.fixB[Y_FIX_A1] = 2.0 * y.fixA[Y_FIX_A0];
|
|
y.fixA[Y_FIX_A2] = y.fixB[Y_FIX_A2] = y.fixA[Y_FIX_A0];
|
|
y.fixA[Y_FIX_B1] = y.fixB[Y_FIX_B1] = 2.0 * (K * K - 1.0) * norm;
|
|
y.fixA[Y_FIX_B2] = y.fixB[Y_FIX_B2] = (1.0 - K / y.fixA[Y_FIX_RESO] + K * K) * norm;
|
|
// for the fixed-position biquad filter
|
|
|
|
for (int s = 0; s < blockSize; s++) {
|
|
double inputSampleL = *in1;
|
|
double inputSampleR = *in2;
|
|
if (fabs(inputSampleL) < 1.18e-23) inputSampleL = y.fpdL * 1.18e-17;
|
|
if (fabs(inputSampleR) < 1.18e-23) inputSampleR = y.fpdR * 1.18e-17;
|
|
double drySampleL = inputSampleL;
|
|
double drySampleR = inputSampleR;
|
|
|
|
double temp = (double)s / inFramesToProcess;
|
|
y.biquad[Y_BIQ_A0] =
|
|
(y.biquad[Y_BIQ_A_A0] * temp) + (y.biquad[Y_BIQ_A_B0] * (1.0 - temp));
|
|
y.biquad[Y_BIQ_A1] =
|
|
(y.biquad[Y_BIQ_A_A1] * temp) + (y.biquad[Y_BIQ_A_B1] * (1.0 - temp));
|
|
y.biquad[Y_BIQ_A2] =
|
|
(y.biquad[Y_BIQ_A_A2] * temp) + (y.biquad[Y_BIQ_A_B2] * (1.0 - temp));
|
|
y.biquad[Y_BIQ_B1] =
|
|
(y.biquad[Y_BIQ_B_A1] * temp) + (y.biquad[Y_BIQ_B_B1] * (1.0 - temp));
|
|
y.biquad[Y_BIQ_B2] =
|
|
(y.biquad[Y_BIQ_B_A2] * temp) + (y.biquad[Y_BIQ_B_B2] * (1.0 - temp));
|
|
// this is the interpolation code for the biquad
|
|
double powFactor = (y.powFactorA * temp) + (y.powFactorB * (1.0 - temp));
|
|
double inTrim = (y.inTrimA * temp) + (y.inTrimB * (1.0 - temp));
|
|
double outTrim = (y.outTrimA * temp) + (y.outTrimB * (1.0 - temp));
|
|
|
|
inputSampleL *= inTrim;
|
|
inputSampleR *= inTrim;
|
|
|
|
temp = (inputSampleL * y.fixA[Y_FIX_A0]) + y.fixA[Y_FIX_S_L1];
|
|
y.fixA[Y_FIX_S_L1] = (inputSampleL * y.fixA[Y_FIX_A1]) -
|
|
(temp * y.fixA[Y_FIX_B1]) + y.fixA[Y_FIX_S_L2];
|
|
y.fixA[Y_FIX_S_L2] =
|
|
(inputSampleL * y.fixA[Y_FIX_A2]) - (temp * y.fixA[Y_FIX_B2]);
|
|
inputSampleL = temp; // fixed biquad filtering ultrasonics
|
|
temp = (inputSampleR * y.fixA[Y_FIX_A0]) + y.fixA[Y_FIX_S_R1];
|
|
y.fixA[Y_FIX_S_R1] = (inputSampleR * y.fixA[Y_FIX_A1]) -
|
|
(temp * y.fixA[Y_FIX_B1]) + y.fixA[Y_FIX_S_R2];
|
|
y.fixA[Y_FIX_S_R2] =
|
|
(inputSampleR * y.fixA[Y_FIX_A2]) - (temp * y.fixA[Y_FIX_B2]);
|
|
inputSampleR = temp; // fixed biquad filtering ultrasonics
|
|
|
|
// encode/decode courtesy of torridgristle under the MIT license
|
|
if (inputSampleL > 1.0) inputSampleL = 1.0;
|
|
else if (inputSampleL > 0.0)
|
|
inputSampleL = 1.0 - pow(1.0 - inputSampleL, powFactor);
|
|
if (inputSampleL < -1.0) inputSampleL = -1.0;
|
|
else if (inputSampleL < 0.0)
|
|
inputSampleL = -1.0 + pow(1.0 + inputSampleL, powFactor);
|
|
if (inputSampleR > 1.0) inputSampleR = 1.0;
|
|
else if (inputSampleR > 0.0)
|
|
inputSampleR = 1.0 - pow(1.0 - inputSampleR, powFactor);
|
|
if (inputSampleR < -1.0) inputSampleR = -1.0;
|
|
else if (inputSampleR < 0.0)
|
|
inputSampleR = -1.0 + pow(1.0 + inputSampleR, powFactor);
|
|
|
|
temp = (inputSampleL * y.biquad[Y_BIQ_A0]) + y.biquad[Y_BIQ_S_L1];
|
|
y.biquad[Y_BIQ_S_L1] = (inputSampleL * y.biquad[Y_BIQ_A1]) -
|
|
(temp * y.biquad[Y_BIQ_B1]) + y.biquad[Y_BIQ_S_L2];
|
|
y.biquad[Y_BIQ_S_L2] =
|
|
(inputSampleL * y.biquad[Y_BIQ_A2]) - (temp * y.biquad[Y_BIQ_B2]);
|
|
inputSampleL = temp; // coefficient interpolating biquad filter
|
|
temp = (inputSampleR * y.biquad[Y_BIQ_A0]) + y.biquad[Y_BIQ_S_R1];
|
|
y.biquad[Y_BIQ_S_R1] = (inputSampleR * y.biquad[Y_BIQ_A1]) -
|
|
(temp * y.biquad[Y_BIQ_B1]) + y.biquad[Y_BIQ_S_R2];
|
|
y.biquad[Y_BIQ_S_R2] =
|
|
(inputSampleR * y.biquad[Y_BIQ_A2]) - (temp * y.biquad[Y_BIQ_B2]);
|
|
inputSampleR = temp; // coefficient interpolating biquad filter
|
|
|
|
// encode/decode courtesy of torridgristle under the MIT license
|
|
if (inputSampleL > 1.0) inputSampleL = 1.0;
|
|
else if (inputSampleL > 0.0)
|
|
inputSampleL = 1.0 - pow(1.0 - inputSampleL, (1.0 / powFactor));
|
|
if (inputSampleL < -1.0) inputSampleL = -1.0;
|
|
else if (inputSampleL < 0.0)
|
|
inputSampleL = -1.0 + pow(1.0 + inputSampleL, (1.0 / powFactor));
|
|
if (inputSampleR > 1.0) inputSampleR = 1.0;
|
|
else if (inputSampleR > 0.0)
|
|
inputSampleR = 1.0 - pow(1.0 - inputSampleR, (1.0 / powFactor));
|
|
if (inputSampleR < -1.0) inputSampleR = -1.0;
|
|
else if (inputSampleR < 0.0)
|
|
inputSampleR = -1.0 + pow(1.0 + inputSampleR, (1.0 / powFactor));
|
|
|
|
inputSampleL *= outTrim;
|
|
inputSampleR *= outTrim;
|
|
|
|
temp = (inputSampleL * y.fixB[Y_FIX_A0]) + y.fixB[Y_FIX_S_L1];
|
|
y.fixB[Y_FIX_S_L1] = (inputSampleL * y.fixB[Y_FIX_A1]) -
|
|
(temp * y.fixB[Y_FIX_B1]) + y.fixB[Y_FIX_S_L2];
|
|
y.fixB[Y_FIX_S_L2] =
|
|
(inputSampleL * y.fixB[Y_FIX_A2]) - (temp * y.fixB[Y_FIX_B2]);
|
|
inputSampleL = temp; // fixed biquad filtering ultrasonics
|
|
temp = (inputSampleR * y.fixB[Y_FIX_A0]) + y.fixB[Y_FIX_S_R1];
|
|
y.fixB[Y_FIX_S_R1] = (inputSampleR * y.fixB[Y_FIX_A1]) -
|
|
(temp * y.fixB[Y_FIX_B1]) + y.fixB[Y_FIX_S_R2];
|
|
y.fixB[Y_FIX_S_R2] =
|
|
(inputSampleR * y.fixB[Y_FIX_A2]) - (temp * y.fixB[Y_FIX_B2]);
|
|
inputSampleR = temp; // fixed biquad filtering ultrasonics
|
|
|
|
if (wet < 1.0) {
|
|
inputSampleL = (inputSampleL * wet) + (drySampleL * (1.0 - wet));
|
|
inputSampleR = (inputSampleR * wet) + (drySampleR * (1.0 - wet));
|
|
}
|
|
|
|
// begin 32 bit stereo floating point dither
|
|
int expon;
|
|
frexpf((float)inputSampleL, &expon);
|
|
y.fpdL ^= y.fpdL << 13;
|
|
y.fpdL ^= y.fpdL >> 17;
|
|
y.fpdL ^= y.fpdL << 5;
|
|
inputSampleL +=
|
|
((double(y.fpdL) - uint32_t(0x7fffffff)) * 5.5e-36l * pow(2, expon + 62));
|
|
frexpf((float)inputSampleR, &expon);
|
|
y.fpdR ^= y.fpdR << 13;
|
|
y.fpdR ^= y.fpdR >> 17;
|
|
y.fpdR ^= y.fpdR << 5;
|
|
inputSampleR +=
|
|
((double(y.fpdR) - uint32_t(0x7fffffff)) * 5.5e-36l * pow(2, expon + 62));
|
|
// end 32 bit stereo floating point dither
|
|
|
|
*out1 = inputSampleL;
|
|
*out2 = inputSampleR;
|
|
|
|
in1++;
|
|
in2++;
|
|
out1++;
|
|
out2++;
|
|
}
|
|
}
|
|
|
|
//////////////
|
|
// HIGHPASS //
|
|
//////////////
|
|
|
|
struct yhighpass {
|
|
double samplerate;
|
|
array<double, Y_BIQ_TOTAL> biquad;
|
|
|
|
double powFactorA;
|
|
double powFactorB;
|
|
double inTrimA;
|
|
double inTrimB;
|
|
double outTrimA;
|
|
double outTrimB;
|
|
|
|
array<double, Y_FIX_TOTAL> fixA;
|
|
array<double, Y_FIX_TOTAL> fixB;
|
|
|
|
uint32_t fpdL;
|
|
uint32_t fpdR;
|
|
// default stuff
|
|
|
|
float drive;
|
|
float frequency;
|
|
float resonance;
|
|
float edge;
|
|
float output;
|
|
float mix; // parameters. Always 0-1, and we scale/alter them elsewhere.
|
|
};
|
|
|
|
inline void yhighpass_init(yhighpass& y, double samplerate)
|
|
{
|
|
y.samplerate = samplerate;
|
|
y.drive = 0.1f;
|
|
y.frequency = 1.0f;
|
|
y.resonance = 0.0f;
|
|
y.edge = 0.1f;
|
|
y.output = 0.9f;
|
|
y.mix = 1.0f;
|
|
y.fpdL = 0;
|
|
y.fpdR = 0;
|
|
y.biquad.fill(0);
|
|
|
|
y.powFactorA = 1.0;
|
|
y.powFactorB = 1.0;
|
|
y.inTrimA = 0.1;
|
|
y.inTrimB = 0.1;
|
|
y.outTrimA = 1.0;
|
|
y.outTrimB = 1.0;
|
|
for (int x = 0; x < Y_FIX_TOTAL; x++) {
|
|
y.fixA[x] = 0.0;
|
|
y.fixB[x] = 0.0;
|
|
}
|
|
|
|
y.fpdL = 1.0;
|
|
while (y.fpdL < 16386) y.fpdL = rand() * UINT32_MAX;
|
|
y.fpdR = 1.0;
|
|
while (y.fpdR < 16386) y.fpdR = rand() * UINT32_MAX;
|
|
}
|
|
|
|
template <typename t_sample>
|
|
inline void yhighpass_process_block(yhighpass& y, t_sample** inputs, t_sample** outputs,
|
|
int blockSize)
|
|
{
|
|
t_sample* in1 = inputs[0];
|
|
t_sample* in2 = inputs[1];
|
|
t_sample* out1 = outputs[0];
|
|
t_sample* out2 = outputs[1];
|
|
|
|
int inFramesToProcess = blockSize;
|
|
double overallscale = 1.0;
|
|
overallscale /= 44100.0;
|
|
overallscale *= y.samplerate;
|
|
|
|
y.inTrimA = y.inTrimB;
|
|
y.inTrimB = y.drive * 10.0;
|
|
|
|
y.biquad[Y_BIQ_FREQ] = pow(y.frequency, 3) * 20000.0;
|
|
if (y.biquad[Y_BIQ_FREQ] < 15.0) y.biquad[Y_BIQ_FREQ] = 15.0;
|
|
y.biquad[Y_BIQ_FREQ] /= y.samplerate;
|
|
y.biquad[Y_BIQ_RESO] = (pow(y.resonance, 2) * 15.0) + 0.5571;
|
|
y.biquad[Y_BIQ_A_A0] = y.biquad[Y_BIQ_A_B0];
|
|
y.biquad[Y_BIQ_A_A1] = y.biquad[Y_BIQ_A_B1];
|
|
y.biquad[Y_BIQ_A_A2] = y.biquad[Y_BIQ_A_B2];
|
|
y.biquad[Y_BIQ_B_A1] = y.biquad[Y_BIQ_B_B1];
|
|
y.biquad[Y_BIQ_B_A2] = y.biquad[Y_BIQ_B_B2];
|
|
// previous run through the buffer is still in the filter, so we move it
|
|
// to the A section and now it's the new starting point.
|
|
double K = tan(M_PI * y.biquad[Y_BIQ_FREQ]);
|
|
double norm = 1.0 / (1.0 + K / y.biquad[Y_BIQ_RESO] + K * K);
|
|
y.biquad[Y_BIQ_A_B0] = norm;
|
|
y.biquad[Y_BIQ_A_B1] = -2.0 * y.biquad[Y_BIQ_A_B0];
|
|
y.biquad[Y_BIQ_A_B2] = y.biquad[Y_BIQ_A_B0];
|
|
y.biquad[Y_BIQ_B_B1] = 2.0 * (K * K - 1.0) * norm;
|
|
y.biquad[Y_BIQ_B_B2] = (1.0 - K / y.biquad[Y_BIQ_RESO] + K * K) * norm;
|
|
// for the coefficient-interpolated biquad filter
|
|
|
|
y.powFactorA = y.powFactorB;
|
|
y.powFactorB = pow(y.edge + 0.9, 4);
|
|
|
|
// 1.0 == target neutral
|
|
|
|
y.outTrimA = y.outTrimB;
|
|
y.outTrimB = y.output;
|
|
|
|
double wet = y.mix;
|
|
|
|
y.fixA[Y_FIX_FREQ] = y.fixB[Y_FIX_FREQ] = 20000.0 / y.samplerate;
|
|
y.fixA[Y_FIX_RESO] = y.fixB[Y_FIX_RESO] = 0.7071; // butterworth Q
|
|
|
|
K = tan(M_PI * y.fixA[Y_FIX_FREQ]);
|
|
norm = 1.0 / (1.0 + K / y.fixA[Y_FIX_RESO] + K * K);
|
|
y.fixA[Y_FIX_A0] = y.fixB[Y_FIX_A0] = K * K * norm;
|
|
y.fixA[Y_FIX_A1] = y.fixB[Y_FIX_A1] = 2.0 * y.fixA[Y_FIX_A0];
|
|
y.fixA[Y_FIX_A2] = y.fixB[Y_FIX_A2] = y.fixA[Y_FIX_A0];
|
|
y.fixA[Y_FIX_B1] = y.fixB[Y_FIX_B1] = 2.0 * (K * K - 1.0) * norm;
|
|
y.fixA[Y_FIX_B2] = y.fixB[Y_FIX_B2] = (1.0 - K / y.fixA[Y_FIX_RESO] + K * K) * norm;
|
|
// for the fixed-position biquad filter
|
|
|
|
for (int s = 0; s < blockSize; s++) {
|
|
double inputSampleL = *in1;
|
|
double inputSampleR = *in2;
|
|
if (fabs(inputSampleL) < 1.18e-23) inputSampleL = y.fpdL * 1.18e-17;
|
|
if (fabs(inputSampleR) < 1.18e-23) inputSampleR = y.fpdR * 1.18e-17;
|
|
double drySampleL = inputSampleL;
|
|
double drySampleR = inputSampleR;
|
|
|
|
double temp = (double)s / inFramesToProcess;
|
|
y.biquad[Y_BIQ_A0] =
|
|
(y.biquad[Y_BIQ_A_A0] * temp) + (y.biquad[Y_BIQ_A_B0] * (1.0 - temp));
|
|
y.biquad[Y_BIQ_A1] =
|
|
(y.biquad[Y_BIQ_A_A1] * temp) + (y.biquad[Y_BIQ_A_B1] * (1.0 - temp));
|
|
y.biquad[Y_BIQ_A2] =
|
|
(y.biquad[Y_BIQ_A_A2] * temp) + (y.biquad[Y_BIQ_A_B2] * (1.0 - temp));
|
|
y.biquad[Y_BIQ_B1] =
|
|
(y.biquad[Y_BIQ_B_A1] * temp) + (y.biquad[Y_BIQ_B_B1] * (1.0 - temp));
|
|
y.biquad[Y_BIQ_B2] =
|
|
(y.biquad[Y_BIQ_B_A2] * temp) + (y.biquad[Y_BIQ_B_B2] * (1.0 - temp));
|
|
// this is the interpolation code for the biquad
|
|
double powFactor = (y.powFactorA * temp) + (y.powFactorB * (1.0 - temp));
|
|
double inTrim = (y.inTrimA * temp) + (y.inTrimB * (1.0 - temp));
|
|
double outTrim = (y.outTrimA * temp) + (y.outTrimB * (1.0 - temp));
|
|
|
|
inputSampleL *= inTrim;
|
|
inputSampleR *= inTrim;
|
|
|
|
temp = (inputSampleL * y.fixA[Y_FIX_A0]) + y.fixA[Y_FIX_S_L1];
|
|
y.fixA[Y_FIX_S_L1] = (inputSampleL * y.fixA[Y_FIX_A1]) -
|
|
(temp * y.fixA[Y_FIX_B1]) + y.fixA[Y_FIX_S_L2];
|
|
y.fixA[Y_FIX_S_L2] =
|
|
(inputSampleL * y.fixA[Y_FIX_A2]) - (temp * y.fixA[Y_FIX_B2]);
|
|
inputSampleL = temp; // fixed biquad filtering ultrasonics
|
|
temp = (inputSampleR * y.fixA[Y_FIX_A0]) + y.fixA[Y_FIX_S_R1];
|
|
y.fixA[Y_FIX_S_R1] = (inputSampleR * y.fixA[Y_FIX_A1]) -
|
|
(temp * y.fixA[Y_FIX_B1]) + y.fixA[Y_FIX_S_R2];
|
|
y.fixA[Y_FIX_S_R2] =
|
|
(inputSampleR * y.fixA[Y_FIX_A2]) - (temp * y.fixA[Y_FIX_B2]);
|
|
inputSampleR = temp; // fixed biquad filtering ultrasonics
|
|
|
|
// encode/decode courtesy of torridgristle under the MIT license
|
|
if (inputSampleL > 1.0) inputSampleL = 1.0;
|
|
else if (inputSampleL > 0.0)
|
|
inputSampleL = 1.0 - pow(1.0 - inputSampleL, powFactor);
|
|
if (inputSampleL < -1.0) inputSampleL = -1.0;
|
|
else if (inputSampleL < 0.0)
|
|
inputSampleL = -1.0 + pow(1.0 + inputSampleL, powFactor);
|
|
if (inputSampleR > 1.0) inputSampleR = 1.0;
|
|
else if (inputSampleR > 0.0)
|
|
inputSampleR = 1.0 - pow(1.0 - inputSampleR, powFactor);
|
|
if (inputSampleR < -1.0) inputSampleR = -1.0;
|
|
else if (inputSampleR < 0.0)
|
|
inputSampleR = -1.0 + pow(1.0 + inputSampleR, powFactor);
|
|
|
|
temp = (inputSampleL * y.biquad[Y_BIQ_A0]) + y.biquad[Y_BIQ_S_L1];
|
|
y.biquad[Y_BIQ_S_L1] = (inputSampleL * y.biquad[Y_BIQ_A1]) -
|
|
(temp * y.biquad[Y_BIQ_B1]) + y.biquad[Y_BIQ_S_L2];
|
|
y.biquad[Y_BIQ_S_L2] =
|
|
(inputSampleL * y.biquad[Y_BIQ_A2]) - (temp * y.biquad[Y_BIQ_B2]);
|
|
inputSampleL = temp; // coefficient interpolating biquad filter
|
|
temp = (inputSampleR * y.biquad[Y_BIQ_A0]) + y.biquad[Y_BIQ_S_R1];
|
|
y.biquad[Y_BIQ_S_R1] = (inputSampleR * y.biquad[Y_BIQ_A1]) -
|
|
(temp * y.biquad[Y_BIQ_B1]) + y.biquad[Y_BIQ_S_R2];
|
|
y.biquad[Y_BIQ_S_R2] =
|
|
(inputSampleR * y.biquad[Y_BIQ_A2]) - (temp * y.biquad[Y_BIQ_B2]);
|
|
inputSampleR = temp; // coefficient interpolating biquad filter
|
|
|
|
// encode/decode courtesy of torridgristle under the MIT license
|
|
if (inputSampleL > 1.0) inputSampleL = 1.0;
|
|
else if (inputSampleL > 0.0)
|
|
inputSampleL = 1.0 - pow(1.0 - inputSampleL, (1.0 / powFactor));
|
|
if (inputSampleL < -1.0) inputSampleL = -1.0;
|
|
else if (inputSampleL < 0.0)
|
|
inputSampleL = -1.0 + pow(1.0 + inputSampleL, (1.0 / powFactor));
|
|
if (inputSampleR > 1.0) inputSampleR = 1.0;
|
|
else if (inputSampleR > 0.0)
|
|
inputSampleR = 1.0 - pow(1.0 - inputSampleR, (1.0 / powFactor));
|
|
if (inputSampleR < -1.0) inputSampleR = -1.0;
|
|
else if (inputSampleR < 0.0)
|
|
inputSampleR = -1.0 + pow(1.0 + inputSampleR, (1.0 / powFactor));
|
|
|
|
inputSampleL *= outTrim;
|
|
inputSampleR *= outTrim;
|
|
|
|
temp = (inputSampleL * y.fixB[Y_FIX_A0]) + y.fixB[Y_FIX_S_L1];
|
|
y.fixB[Y_FIX_S_L1] = (inputSampleL * y.fixB[Y_FIX_A1]) -
|
|
(temp * y.fixB[Y_FIX_B1]) + y.fixB[Y_FIX_S_L2];
|
|
y.fixB[Y_FIX_S_L2] =
|
|
(inputSampleL * y.fixB[Y_FIX_A2]) - (temp * y.fixB[Y_FIX_B2]);
|
|
inputSampleL = temp; // fixed biquad filtering ultrasonics
|
|
temp = (inputSampleR * y.fixB[Y_FIX_A0]) + y.fixB[Y_FIX_S_R1];
|
|
y.fixB[Y_FIX_S_R1] = (inputSampleR * y.fixB[Y_FIX_A1]) -
|
|
(temp * y.fixB[Y_FIX_B1]) + y.fixB[Y_FIX_S_R2];
|
|
y.fixB[Y_FIX_S_R2] =
|
|
(inputSampleR * y.fixB[Y_FIX_A2]) - (temp * y.fixB[Y_FIX_B2]);
|
|
inputSampleR = temp; // fixed biquad filtering ultrasonics
|
|
|
|
if (wet < 1.0) {
|
|
inputSampleL = (inputSampleL * wet) + (drySampleL * (1.0 - wet));
|
|
inputSampleR = (inputSampleR * wet) + (drySampleR * (1.0 - wet));
|
|
}
|
|
|
|
// begin 32 bit stereo floating point dither
|
|
int expon;
|
|
frexpf((float)inputSampleL, &expon);
|
|
y.fpdL ^= y.fpdL << 13;
|
|
y.fpdL ^= y.fpdL >> 17;
|
|
y.fpdL ^= y.fpdL << 5;
|
|
inputSampleL +=
|
|
((double(y.fpdL) - uint32_t(0x7fffffff)) * 5.5e-36l * pow(2, expon + 62));
|
|
frexpf((float)inputSampleR, &expon);
|
|
y.fpdR ^= y.fpdR << 13;
|
|
y.fpdR ^= y.fpdR >> 17;
|
|
y.fpdR ^= y.fpdR << 5;
|
|
inputSampleR +=
|
|
((double(y.fpdR) - uint32_t(0x7fffffff)) * 5.5e-36l * pow(2, expon + 62));
|
|
// end 32 bit stereo floating point dither
|
|
|
|
*out1 = inputSampleL;
|
|
*out2 = inputSampleR;
|
|
|
|
in1++;
|
|
in2++;
|
|
out1++;
|
|
out2++;
|
|
}
|
|
}
|
|
|
|
//////////////
|
|
// BANDPASS //
|
|
//////////////
|
|
|
|
struct ybandpass {
|
|
double samplerate;
|
|
array<double, Y_BIQ_TOTAL> biquad;
|
|
|
|
double powFactorA;
|
|
double powFactorB;
|
|
double inTrimA;
|
|
double inTrimB;
|
|
double outTrimA;
|
|
double outTrimB;
|
|
|
|
array<double, Y_FIX_TOTAL> fixA;
|
|
array<double, Y_FIX_TOTAL> fixB;
|
|
|
|
uint32_t fpdL;
|
|
uint32_t fpdR;
|
|
// default stuff
|
|
|
|
float drive;
|
|
float frequency;
|
|
float resonance;
|
|
float edge;
|
|
float output;
|
|
float mix; // parameters. Always 0-1, and we scale/alter them elsewhere.
|
|
};
|
|
|
|
inline void ybandpass_init(ybandpass& y, double samplerate)
|
|
{
|
|
y.samplerate = samplerate;
|
|
y.drive = 0.1f;
|
|
y.frequency = 1.0f;
|
|
y.resonance = 0.02f;
|
|
y.edge = 0.1f;
|
|
y.output = 0.9f;
|
|
y.mix = 1.0f;
|
|
y.fpdL = 0;
|
|
y.fpdR = 0;
|
|
y.biquad.fill(0);
|
|
|
|
for (int x = 0; x < Y_BIQ_TOTAL; x++) { y.biquad[x] = 0.0; }
|
|
y.powFactorA = 1.0;
|
|
y.powFactorB = 1.0;
|
|
y.inTrimA = 0.1;
|
|
y.inTrimB = 0.1;
|
|
y.outTrimA = 1.0;
|
|
y.outTrimB = 1.0;
|
|
for (int x = 0; x < Y_FIX_TOTAL; x++) {
|
|
y.fixA[x] = 0.0;
|
|
y.fixB[x] = 0.0;
|
|
}
|
|
|
|
y.fpdL = 1.0;
|
|
while (y.fpdL < 16386) y.fpdL = rand() * UINT32_MAX;
|
|
y.fpdR = 1.0;
|
|
while (y.fpdR < 16386) y.fpdR = rand() * UINT32_MAX;
|
|
}
|
|
|
|
template <typename t_sample>
|
|
inline void ybandpass_process_block(ybandpass& y, t_sample** inputs, t_sample** outputs,
|
|
int blockSize)
|
|
{
|
|
t_sample* in1 = inputs[0];
|
|
t_sample* in2 = inputs[1];
|
|
t_sample* out1 = outputs[0];
|
|
t_sample* out2 = outputs[1];
|
|
|
|
int inFramesToProcess = blockSize;
|
|
double overallscale = 1.0;
|
|
overallscale /= 44100.0;
|
|
overallscale *= y.samplerate;
|
|
|
|
y.inTrimA = y.inTrimB;
|
|
y.inTrimB = y.drive * 10.0;
|
|
|
|
y.biquad[Y_BIQ_FREQ] = pow(y.frequency, 3) * 20000.0;
|
|
if (y.biquad[Y_BIQ_FREQ] < 15.0) y.biquad[Y_BIQ_FREQ] = 15.0;
|
|
y.biquad[Y_BIQ_FREQ] /= y.samplerate;
|
|
y.biquad[Y_BIQ_RESO] = (pow(y.resonance, 2) * 15.0) + 0.5571;
|
|
y.biquad[Y_BIQ_A_A0] = y.biquad[Y_BIQ_A_B0];
|
|
// y.biquad[biq_aA1] = y.biquad[biq_aB1];
|
|
y.biquad[Y_BIQ_A_A2] = y.biquad[Y_BIQ_A_B2];
|
|
y.biquad[Y_BIQ_B_A1] = y.biquad[Y_BIQ_B_B1];
|
|
y.biquad[Y_BIQ_B_A2] = y.biquad[Y_BIQ_B_B2];
|
|
// previous run through the buffer is still in the filter, so we move it
|
|
// to the A section and now it's the new starting point.
|
|
double K = tan(M_PI * y.biquad[Y_BIQ_FREQ]);
|
|
double norm = 1.0 / (1.0 + K / y.biquad[Y_BIQ_RESO] + K * K);
|
|
y.biquad[Y_BIQ_A_B0] = K / y.biquad[Y_BIQ_RESO] * norm;
|
|
// y.biquad[y.biq_aB1] = 0.0; //bandpass can simplify the biquad kernel: leave out
|
|
// this multiply
|
|
y.biquad[Y_BIQ_A_B2] = -y.biquad[Y_BIQ_A_B0];
|
|
y.biquad[Y_BIQ_B_B1] = 2.0 * (K * K - 1.0) * norm;
|
|
y.biquad[Y_BIQ_B_B2] = (1.0 - K / y.biquad[Y_BIQ_RESO] + K * K) * norm;
|
|
// for the coefficient-interpolated biquad filter
|
|
|
|
y.powFactorA = y.powFactorB;
|
|
y.powFactorB = pow(y.drive + 0.9, 4);
|
|
|
|
// 1.0 == target neutral
|
|
|
|
y.outTrimA = y.outTrimB;
|
|
y.outTrimB = y.outTrimB;
|
|
|
|
double wet = y.mix;
|
|
|
|
y.fixA[Y_FIX_FREQ] = y.fixB[Y_FIX_FREQ] = 20000.0 / y.samplerate;
|
|
y.fixA[Y_FIX_RESO] = y.fixB[Y_FIX_RESO] = 0.7071; // butterworth Q
|
|
|
|
K = tan(M_PI * y.fixA[Y_FIX_FREQ]);
|
|
norm = 1.0 / (1.0 + K / y.fixA[Y_FIX_RESO] + K * K);
|
|
y.fixA[Y_FIX_A0] = y.fixB[Y_FIX_A0] = K * K * norm;
|
|
y.fixA[Y_FIX_A1] = y.fixB[Y_FIX_A1] = 2.0 * y.fixA[Y_FIX_A0];
|
|
y.fixA[Y_FIX_A2] = y.fixB[Y_FIX_A2] = y.fixA[Y_FIX_A0];
|
|
y.fixA[Y_FIX_B1] = y.fixB[Y_FIX_B1] = 2.0 * (K * K - 1.0) * norm;
|
|
y.fixA[Y_FIX_B2] = y.fixB[Y_FIX_B2] = (1.0 - K / y.fixA[Y_FIX_RESO] + K * K) * norm;
|
|
// for the fixed-position biquad filter
|
|
|
|
for (int s = 0; s < blockSize; s++) {
|
|
double inputSampleL = *in1;
|
|
double inputSampleR = *in2;
|
|
if (fabs(inputSampleL) < 1.18e-23) inputSampleL = y.fpdL * 1.18e-17;
|
|
if (fabs(inputSampleR) < 1.18e-23) inputSampleR = y.fpdR * 1.18e-17;
|
|
double drySampleL = inputSampleL;
|
|
double drySampleR = inputSampleR;
|
|
|
|
double temp = (double)s / inFramesToProcess;
|
|
y.biquad[Y_BIQ_A0] =
|
|
(y.biquad[Y_BIQ_A_A0] * temp) + (y.biquad[Y_BIQ_A_B0] * (1.0 - temp));
|
|
// biquad[biq_a1] = (biquad[biq_aA1]*temp)+(biquad[biq_aB1]*(1.0-temp));
|
|
y.biquad[Y_BIQ_A2] =
|
|
(y.biquad[Y_BIQ_A_A2] * temp) + (y.biquad[Y_BIQ_A_B2] * (1.0 - temp));
|
|
y.biquad[Y_BIQ_B1] =
|
|
(y.biquad[Y_BIQ_B_A1] * temp) + (y.biquad[Y_BIQ_B_B1] * (1.0 - temp));
|
|
y.biquad[Y_BIQ_B2] =
|
|
(y.biquad[Y_BIQ_B_A2] * temp) + (y.biquad[Y_BIQ_B_B2] * (1.0 - temp));
|
|
// this is the interpolation code for the biquad
|
|
double powFactor = (y.powFactorA * temp) + (y.powFactorB * (1.0 - temp));
|
|
double inTrim = (y.inTrimA * temp) + (y.inTrimB * (1.0 - temp));
|
|
double outTrim = (y.outTrimA * temp) + (y.outTrimB * (1.0 - temp));
|
|
|
|
inputSampleL *= inTrim;
|
|
inputSampleR *= inTrim;
|
|
|
|
temp = (inputSampleL * y.fixA[Y_FIX_A0]) + y.fixA[Y_FIX_S_L1];
|
|
y.fixA[Y_FIX_S_L1] = (inputSampleL * y.fixA[Y_FIX_A1]) -
|
|
(temp * y.fixA[Y_FIX_B1]) + y.fixA[Y_FIX_S_L2];
|
|
y.fixA[Y_FIX_S_L2] =
|
|
(inputSampleL * y.fixA[Y_FIX_A2]) - (temp * y.fixA[Y_FIX_B2]);
|
|
inputSampleL = temp; // fixed biquad filtering ultrasonics
|
|
temp = (inputSampleR * y.fixA[Y_FIX_A0]) + y.fixA[Y_FIX_S_R1];
|
|
y.fixA[Y_FIX_S_R1] = (inputSampleR * y.fixA[Y_FIX_A1]) -
|
|
(temp * y.fixA[Y_FIX_B1]) + y.fixA[Y_FIX_S_R2];
|
|
y.fixA[Y_FIX_S_R2] =
|
|
(inputSampleR * y.fixA[Y_FIX_A2]) - (temp * y.fixA[Y_FIX_B2]);
|
|
inputSampleR = temp; // fixed biquad filtering ultrasonics
|
|
|
|
// encode/decode courtesy of torridgristle under the MIT license
|
|
if (inputSampleL > 1.0) inputSampleL = 1.0;
|
|
else if (inputSampleL > 0.0)
|
|
inputSampleL = 1.0 - pow(1.0 - inputSampleL, powFactor);
|
|
if (inputSampleL < -1.0) inputSampleL = -1.0;
|
|
else if (inputSampleL < 0.0)
|
|
inputSampleL = -1.0 + pow(1.0 + inputSampleL, powFactor);
|
|
if (inputSampleR > 1.0) inputSampleR = 1.0;
|
|
else if (inputSampleR > 0.0)
|
|
inputSampleR = 1.0 - pow(1.0 - inputSampleR, powFactor);
|
|
if (inputSampleR < -1.0) inputSampleR = -1.0;
|
|
else if (inputSampleR < 0.0)
|
|
inputSampleR = -1.0 + pow(1.0 + inputSampleR, powFactor);
|
|
|
|
temp = (inputSampleL * y.biquad[Y_BIQ_A0]) + y.biquad[Y_BIQ_S_L1];
|
|
y.biquad[Y_BIQ_S_L1] = -(temp * y.biquad[Y_BIQ_B1]) + y.biquad[Y_BIQ_S_L2];
|
|
y.biquad[Y_BIQ_S_L2] =
|
|
(inputSampleL * y.biquad[Y_BIQ_A2]) - (temp * y.biquad[Y_BIQ_B2]);
|
|
inputSampleL = temp; // coefficient interpolating biquad filter
|
|
temp = (inputSampleR * y.biquad[Y_BIQ_A0]) + y.biquad[Y_BIQ_S_R1];
|
|
y.biquad[Y_BIQ_S_R1] = -(temp * y.biquad[Y_BIQ_B1]) + y.biquad[Y_BIQ_S_R2];
|
|
y.biquad[Y_BIQ_S_R2] =
|
|
(inputSampleR * y.biquad[Y_BIQ_A2]) - (temp * y.biquad[Y_BIQ_B2]);
|
|
inputSampleR = temp; // coefficient interpolating biquad filter
|
|
|
|
// encode/decode courtesy of torridgristle under the MIT license
|
|
if (inputSampleL > 1.0) inputSampleL = 1.0;
|
|
else if (inputSampleL > 0.0)
|
|
inputSampleL = 1.0 - pow(1.0 - inputSampleL, (1.0 / powFactor));
|
|
if (inputSampleL < -1.0) inputSampleL = -1.0;
|
|
else if (inputSampleL < 0.0)
|
|
inputSampleL = -1.0 + pow(1.0 + inputSampleL, (1.0 / powFactor));
|
|
if (inputSampleR > 1.0) inputSampleR = 1.0;
|
|
else if (inputSampleR > 0.0)
|
|
inputSampleR = 1.0 - pow(1.0 - inputSampleR, (1.0 / powFactor));
|
|
if (inputSampleR < -1.0) inputSampleR = -1.0;
|
|
else if (inputSampleR < 0.0)
|
|
inputSampleR = -1.0 + pow(1.0 + inputSampleR, (1.0 / powFactor));
|
|
|
|
inputSampleL *= outTrim;
|
|
inputSampleR *= outTrim;
|
|
|
|
temp = (inputSampleL * y.fixB[Y_FIX_A0]) + y.fixB[Y_FIX_S_L1];
|
|
y.fixB[Y_FIX_S_L1] = (inputSampleL * y.fixB[Y_FIX_A1]) -
|
|
(temp * y.fixB[Y_FIX_B1]) + y.fixB[Y_FIX_S_L2];
|
|
y.fixB[Y_FIX_S_L2] =
|
|
(inputSampleL * y.fixB[Y_FIX_A2]) - (temp * y.fixB[Y_FIX_B2]);
|
|
inputSampleL = temp; // fixed biquad filtering ultrasonics
|
|
temp = (inputSampleR * y.fixB[Y_FIX_A0]) + y.fixB[Y_FIX_S_R1];
|
|
y.fixB[Y_FIX_S_R1] = (inputSampleR * y.fixB[Y_FIX_A1]) -
|
|
(temp * y.fixB[Y_FIX_B1]) + y.fixB[Y_FIX_S_R2];
|
|
y.fixB[Y_FIX_S_R2] =
|
|
(inputSampleR * y.fixB[Y_FIX_A2]) - (temp * y.fixB[Y_FIX_B2]);
|
|
inputSampleR = temp; // fixed biquad filtering ultrasonics
|
|
|
|
if (wet < 1.0) {
|
|
inputSampleL = (inputSampleL * wet) + (drySampleL * (1.0 - wet));
|
|
inputSampleR = (inputSampleR * wet) + (drySampleR * (1.0 - wet));
|
|
}
|
|
|
|
// begin 32 bit stereo floating point dither
|
|
int expon;
|
|
frexpf((float)inputSampleL, &expon);
|
|
y.fpdL ^= y.fpdL << 13;
|
|
y.fpdL ^= y.fpdL >> 17;
|
|
y.fpdL ^= y.fpdL << 5;
|
|
inputSampleL +=
|
|
((double(y.fpdL) - uint32_t(0x7fffffff)) * 5.5e-36l * pow(2, expon + 62));
|
|
frexpf((float)inputSampleR, &expon);
|
|
y.fpdR ^= y.fpdR << 13;
|
|
y.fpdR ^= y.fpdR >> 17;
|
|
y.fpdR ^= y.fpdR << 5;
|
|
inputSampleR +=
|
|
((double(y.fpdR) - uint32_t(0x7fffffff)) * 5.5e-36l * pow(2, expon + 62));
|
|
// end 32 bit stereo floating point dither
|
|
|
|
*out1 = inputSampleL;
|
|
*out2 = inputSampleR;
|
|
|
|
in1++;
|
|
in2++;
|
|
out1++;
|
|
out2++;
|
|
}
|
|
}
|
|
|
|
///////////
|
|
// NOTCH //
|
|
///////////
|
|
|
|
struct ynotch {
|
|
double samplerate;
|
|
array<double, Y_BIQ_TOTAL> biquad;
|
|
|
|
double powFactorA;
|
|
double powFactorB;
|
|
double inTrimA;
|
|
double inTrimB;
|
|
double outTrimA;
|
|
double outTrimB;
|
|
|
|
array<double, Y_FIX_TOTAL> fixA;
|
|
array<double, Y_FIX_TOTAL> fixB;
|
|
|
|
uint32_t fpdL;
|
|
uint32_t fpdR;
|
|
// default stuff
|
|
|
|
float drive;
|
|
float frequency;
|
|
float resonance;
|
|
float edge;
|
|
float output;
|
|
float mix; // parameters. Always 0-1, and we scale/alter them elsewhere.
|
|
};
|
|
|
|
inline void ynotch_init(ynotch& y, double samplerate)
|
|
{
|
|
y.samplerate = samplerate;
|
|
y.drive = 0.1f;
|
|
y.frequency = 1.0f;
|
|
y.resonance = 0.0f;
|
|
y.edge = 0.1f;
|
|
y.output = 0.9f;
|
|
y.mix = 1.0f;
|
|
y.fpdL = 0;
|
|
y.fpdR = 0;
|
|
y.biquad.fill(0);
|
|
|
|
y.powFactorA = 1.0;
|
|
y.powFactorB = 1.0;
|
|
y.inTrimA = 0.1;
|
|
y.inTrimB = 0.1;
|
|
y.outTrimA = 1.0;
|
|
y.outTrimB = 1.0;
|
|
for (int x = 0; x < Y_FIX_TOTAL; x++) {
|
|
y.fixA[x] = 0.0;
|
|
y.fixB[x] = 0.0;
|
|
}
|
|
|
|
y.fpdL = 1.0;
|
|
while (y.fpdL < 16386) y.fpdL = rand() * UINT32_MAX;
|
|
y.fpdR = 1.0;
|
|
while (y.fpdR < 16386) y.fpdR = rand() * UINT32_MAX;
|
|
}
|
|
|
|
template <typename t_sample>
|
|
inline void ynotch_process_block(ynotch& y, t_sample** inputs, t_sample** outputs,
|
|
int blockSize)
|
|
{
|
|
t_sample* in1 = inputs[0];
|
|
t_sample* in2 = inputs[1];
|
|
t_sample* out1 = outputs[0];
|
|
t_sample* out2 = outputs[1];
|
|
|
|
int inFramesToProcess = blockSize;
|
|
double overallscale = 1.0;
|
|
overallscale /= 44100.0;
|
|
overallscale *= y.samplerate;
|
|
|
|
y.inTrimA = y.inTrimB;
|
|
y.inTrimB = y.drive * 10.0;
|
|
|
|
y.biquad[Y_BIQ_FREQ] = pow(y.frequency, 3) * 20000.0;
|
|
if (y.biquad[Y_BIQ_FREQ] < 15.0) y.biquad[Y_BIQ_FREQ] = 15.0;
|
|
y.biquad[Y_BIQ_FREQ] /= y.samplerate;
|
|
y.biquad[Y_BIQ_RESO] = (pow(y.resonance, 2) * 15.0) + 0.0001;
|
|
y.biquad[Y_BIQ_A_A0] = y.biquad[Y_BIQ_A_B0];
|
|
y.biquad[Y_BIQ_A_A1] = y.biquad[Y_BIQ_A_B1];
|
|
y.biquad[Y_BIQ_A_A2] = y.biquad[Y_BIQ_A_B2];
|
|
y.biquad[Y_BIQ_B_A1] = y.biquad[Y_BIQ_B_B1];
|
|
y.biquad[Y_BIQ_B_A2] = y.biquad[Y_BIQ_B_B2];
|
|
// previous run through the buffer is still in the filter, so we move it
|
|
// to the A section and now it's the new starting point.
|
|
double K = tan(M_PI * y.biquad[Y_BIQ_FREQ]);
|
|
double norm = 1.0 / (1.0 + K / y.biquad[Y_BIQ_RESO] + K * K);
|
|
y.biquad[Y_BIQ_A_B0] = (1.0 + K * K) * norm;
|
|
y.biquad[Y_BIQ_A_B1] = 2.0 * (K * K - 1) * norm;
|
|
y.biquad[Y_BIQ_A_B2] = y.biquad[Y_BIQ_A_B0];
|
|
y.biquad[Y_BIQ_B_B1] = y.biquad[Y_BIQ_A_B1];
|
|
y.biquad[Y_BIQ_B_B2] = (1.0 - K / y.biquad[Y_BIQ_RESO] + K * K) * norm;
|
|
// for the coefficient-interpolated biquad filter
|
|
|
|
y.powFactorA = y.powFactorB;
|
|
y.powFactorB = pow(y.edge + 0.9, 4);
|
|
|
|
// 1.0 == target neutral
|
|
|
|
y.outTrimA = y.outTrimB;
|
|
y.outTrimB = y.output;
|
|
|
|
double wet = y.mix;
|
|
|
|
y.fixA[Y_FIX_FREQ] = y.fixB[Y_FIX_FREQ] = 20000.0 / y.samplerate;
|
|
y.fixA[Y_FIX_RESO] = y.fixB[Y_FIX_RESO] = 0.7071; // butterworth Q
|
|
|
|
K = tan(M_PI * y.fixA[Y_FIX_FREQ]);
|
|
norm = 1.0 / (1.0 + K / y.fixA[Y_FIX_RESO] + K * K);
|
|
y.fixA[Y_FIX_A0] = y.fixB[Y_FIX_A0] = K * K * norm;
|
|
y.fixA[Y_FIX_A1] = y.fixB[Y_FIX_A1] = 2.0 * y.fixA[Y_FIX_A0];
|
|
y.fixA[Y_FIX_A2] = y.fixB[Y_FIX_A2] = y.fixA[Y_FIX_A0];
|
|
y.fixA[Y_FIX_B1] = y.fixB[Y_FIX_B1] = 2.0 * (K * K - 1.0) * norm;
|
|
y.fixA[Y_FIX_B2] = y.fixB[Y_FIX_B2] = (1.0 - K / y.fixA[Y_FIX_RESO] + K * K) * norm;
|
|
// for the fixed-position biquad filter
|
|
|
|
for (int s = 0; s < blockSize; s++) {
|
|
double inputSampleL = *in1;
|
|
double inputSampleR = *in2;
|
|
if (fabs(inputSampleL) < 1.18e-23) inputSampleL = y.fpdL * 1.18e-17;
|
|
if (fabs(inputSampleR) < 1.18e-23) inputSampleR = y.fpdR * 1.18e-17;
|
|
double drySampleL = inputSampleL;
|
|
double drySampleR = inputSampleR;
|
|
|
|
double temp = (double)s / inFramesToProcess;
|
|
y.biquad[Y_BIQ_A0] =
|
|
(y.biquad[Y_BIQ_A_A0] * temp) + (y.biquad[Y_BIQ_A_B0] * (1.0 - temp));
|
|
y.biquad[Y_BIQ_A1] =
|
|
(y.biquad[Y_BIQ_A_A1] * temp) + (y.biquad[Y_BIQ_A_B1] * (1.0 - temp));
|
|
y.biquad[Y_BIQ_A2] =
|
|
(y.biquad[Y_BIQ_A_A2] * temp) + (y.biquad[Y_BIQ_A_B2] * (1.0 - temp));
|
|
y.biquad[Y_BIQ_B1] =
|
|
(y.biquad[Y_BIQ_B_A1] * temp) + (y.biquad[Y_BIQ_B_B1] * (1.0 - temp));
|
|
y.biquad[Y_BIQ_B2] =
|
|
(y.biquad[Y_BIQ_B_A2] * temp) + (y.biquad[Y_BIQ_B_B2] * (1.0 - temp));
|
|
// this is the interpolation code for the biquad
|
|
double powFactor = (y.powFactorA * temp) + (y.powFactorB * (1.0 - temp));
|
|
double inTrim = (y.inTrimA * temp) + (y.inTrimB * (1.0 - temp));
|
|
double outTrim = (y.outTrimA * temp) + (y.outTrimB * (1.0 - temp));
|
|
|
|
inputSampleL *= inTrim;
|
|
inputSampleR *= inTrim;
|
|
|
|
temp = (inputSampleL * y.fixA[Y_FIX_A0]) + y.fixA[Y_FIX_S_L1];
|
|
y.fixA[Y_FIX_S_L1] = (inputSampleL * y.fixA[Y_FIX_A1]) -
|
|
(temp * y.fixA[Y_FIX_B1]) + y.fixA[Y_FIX_S_L2];
|
|
y.fixA[Y_FIX_S_L2] =
|
|
(inputSampleL * y.fixA[Y_FIX_A2]) - (temp * y.fixA[Y_FIX_B2]);
|
|
inputSampleL = temp; // fixed biquad filtering ultrasonics
|
|
temp = (inputSampleR * y.fixA[Y_FIX_A0]) + y.fixA[Y_FIX_S_R1];
|
|
y.fixA[Y_FIX_S_R1] = (inputSampleR * y.fixA[Y_FIX_A1]) -
|
|
(temp * y.fixA[Y_FIX_B1]) + y.fixA[Y_FIX_S_R2];
|
|
y.fixA[Y_FIX_S_R2] =
|
|
(inputSampleR * y.fixA[Y_FIX_A2]) - (temp * y.fixA[Y_FIX_B2]);
|
|
inputSampleR = temp; // fixed biquad filtering ultrasonics
|
|
|
|
// encode/decode courtesy of torridgristle under the MIT license
|
|
if (inputSampleL > 1.0) inputSampleL = 1.0;
|
|
else if (inputSampleL > 0.0)
|
|
inputSampleL = 1.0 - pow(1.0 - inputSampleL, powFactor);
|
|
if (inputSampleL < -1.0) inputSampleL = -1.0;
|
|
else if (inputSampleL < 0.0)
|
|
inputSampleL = -1.0 + pow(1.0 + inputSampleL, powFactor);
|
|
if (inputSampleR > 1.0) inputSampleR = 1.0;
|
|
else if (inputSampleR > 0.0)
|
|
inputSampleR = 1.0 - pow(1.0 - inputSampleR, powFactor);
|
|
if (inputSampleR < -1.0) inputSampleR = -1.0;
|
|
else if (inputSampleR < 0.0)
|
|
inputSampleR = -1.0 + pow(1.0 + inputSampleR, powFactor);
|
|
|
|
temp = (inputSampleL * y.biquad[Y_BIQ_A0]) + y.biquad[Y_BIQ_S_L1];
|
|
y.biquad[Y_BIQ_S_L1] = (inputSampleL * y.biquad[Y_BIQ_A1]) -
|
|
(temp * y.biquad[Y_BIQ_B1]) + y.biquad[Y_BIQ_S_L2];
|
|
y.biquad[Y_BIQ_S_L2] =
|
|
(inputSampleL * y.biquad[Y_BIQ_A2]) - (temp * y.biquad[Y_BIQ_B2]);
|
|
inputSampleL = temp; // coefficient interpolating biquad filter
|
|
temp = (inputSampleR * y.biquad[Y_BIQ_A0]) + y.biquad[Y_BIQ_S_R1];
|
|
y.biquad[Y_BIQ_S_R1] = (inputSampleR * y.biquad[Y_BIQ_A1]) -
|
|
(temp * y.biquad[Y_BIQ_B1]) + y.biquad[Y_BIQ_S_R2];
|
|
y.biquad[Y_BIQ_S_R2] =
|
|
(inputSampleR * y.biquad[Y_BIQ_A2]) - (temp * y.biquad[Y_BIQ_B2]);
|
|
inputSampleR = temp; // coefficient interpolating biquad filter
|
|
|
|
// encode/decode courtesy of torridgristle under the MIT license
|
|
if (inputSampleL > 1.0) inputSampleL = 1.0;
|
|
else if (inputSampleL > 0.0)
|
|
inputSampleL = 1.0 - pow(1.0 - inputSampleL, (1.0 / powFactor));
|
|
if (inputSampleL < -1.0) inputSampleL = -1.0;
|
|
else if (inputSampleL < 0.0)
|
|
inputSampleL = -1.0 + pow(1.0 + inputSampleL, (1.0 / powFactor));
|
|
if (inputSampleR > 1.0) inputSampleR = 1.0;
|
|
else if (inputSampleR > 0.0)
|
|
inputSampleR = 1.0 - pow(1.0 - inputSampleR, (1.0 / powFactor));
|
|
if (inputSampleR < -1.0) inputSampleR = -1.0;
|
|
else if (inputSampleR < 0.0)
|
|
inputSampleR = -1.0 + pow(1.0 + inputSampleR, (1.0 / powFactor));
|
|
|
|
inputSampleL *= outTrim;
|
|
inputSampleR *= outTrim;
|
|
|
|
temp = (inputSampleL * y.fixB[Y_FIX_A0]) + y.fixB[Y_FIX_S_L1];
|
|
y.fixB[Y_FIX_S_L1] = (inputSampleL * y.fixB[Y_FIX_A1]) -
|
|
(temp * y.fixB[Y_FIX_B1]) + y.fixB[Y_FIX_S_L2];
|
|
y.fixB[Y_FIX_S_L2] =
|
|
(inputSampleL * y.fixB[Y_FIX_A2]) - (temp * y.fixB[Y_FIX_B2]);
|
|
inputSampleL = temp; // fixed biquad filtering ultrasonics
|
|
temp = (inputSampleR * y.fixB[Y_FIX_A0]) + y.fixB[Y_FIX_S_R1];
|
|
y.fixB[Y_FIX_S_R1] = (inputSampleR * y.fixB[Y_FIX_A1]) -
|
|
(temp * y.fixB[Y_FIX_B1]) + y.fixB[Y_FIX_S_R2];
|
|
y.fixB[Y_FIX_S_R2] =
|
|
(inputSampleR * y.fixB[Y_FIX_A2]) - (temp * y.fixB[Y_FIX_B2]);
|
|
inputSampleR = temp; // fixed biquad filtering ultrasonics
|
|
|
|
if (wet < 1.0) {
|
|
inputSampleL = (inputSampleL * wet) + (drySampleL * (1.0 - wet));
|
|
inputSampleR = (inputSampleR * wet) + (drySampleR * (1.0 - wet));
|
|
}
|
|
|
|
// begin 32 bit stereo floating point dither
|
|
int expon;
|
|
frexpf((float)inputSampleL, &expon);
|
|
y.fpdL ^= y.fpdL << 13;
|
|
y.fpdL ^= y.fpdL >> 17;
|
|
y.fpdL ^= y.fpdL << 5;
|
|
inputSampleL +=
|
|
((double(y.fpdL) - uint32_t(0x7fffffff)) * 5.5e-36l * pow(2, expon + 62));
|
|
frexpf((float)inputSampleR, &expon);
|
|
y.fpdR ^= y.fpdR << 13;
|
|
y.fpdR ^= y.fpdR >> 17;
|
|
y.fpdR ^= y.fpdR << 5;
|
|
inputSampleR +=
|
|
((double(y.fpdR) - uint32_t(0x7fffffff)) * 5.5e-36l * pow(2, expon + 62));
|
|
// end 32 bit stereo floating point dither
|
|
|
|
*out1 = inputSampleL;
|
|
*out2 = inputSampleR;
|
|
|
|
in1++;
|
|
in2++;
|
|
out1++;
|
|
out2++;
|
|
}
|
|
}
|
|
|
|
/////////
|
|
// SVF //
|
|
/////////
|
|
|
|
enum ysvf_types {
|
|
Y_LOWPASS,
|
|
Y_HIGHPASS,
|
|
Y_BANDPASS,
|
|
Y_NOTCH
|
|
};
|
|
|
|
struct ysvf {
|
|
ysvf_types filter_type;
|
|
ylowpass lowpass;
|
|
yhighpass highpass;
|
|
ybandpass bandpass;
|
|
ynotch notch;
|
|
};
|
|
|
|
enum ysvf_parameters {
|
|
Y_DRIVE,
|
|
Y_FREQUENCY,
|
|
Y_RESONANCE,
|
|
Y_EDGE,
|
|
};
|
|
|
|
inline void ysvf_init(ysvf& y, double samplerate = 44100)
|
|
{
|
|
ylowpass_init(y.lowpass, samplerate);
|
|
yhighpass_init(y.highpass, samplerate);
|
|
ybandpass_init(y.bandpass, samplerate);
|
|
ynotch_init(y.notch, samplerate);
|
|
y.filter_type = ysvf_types::Y_LOWPASS;
|
|
}
|
|
|
|
inline void ysvf_set_param(ysvf& y, ysvf_parameters param, float value)
|
|
{
|
|
switch (param) {
|
|
case Y_DRIVE:
|
|
y.lowpass.drive = value;
|
|
y.highpass.drive = value;
|
|
y.bandpass.drive = value;
|
|
y.notch.drive = value;
|
|
break;
|
|
case Y_FREQUENCY:
|
|
y.lowpass.frequency = value;
|
|
y.highpass.frequency = value;
|
|
y.bandpass.frequency = value;
|
|
y.notch.frequency = value;
|
|
break;
|
|
case Y_RESONANCE:
|
|
y.lowpass.resonance = value;
|
|
y.highpass.resonance = value;
|
|
y.bandpass.resonance = value;
|
|
y.notch.resonance = value;
|
|
break;
|
|
case Y_EDGE:
|
|
y.lowpass.edge = value;
|
|
y.highpass.edge = value;
|
|
y.bandpass.edge = value;
|
|
y.notch.edge = value;
|
|
break;
|
|
}
|
|
}
|
|
|
|
template <typename t_sample>
|
|
inline void ysvf_process_samples(ysvf& y, t_sample** inputs, t_sample** outputs,
|
|
int block_size)
|
|
{
|
|
switch (y.filter_type) {
|
|
case Y_LOWPASS:
|
|
ylowpass_process_block(y.lowpass, inputs, outputs, block_size);
|
|
break;
|
|
case ysvf_types::Y_HIGHPASS:
|
|
yhighpass_process_block(y.highpass, inputs, outputs, block_size);
|
|
break;
|
|
case ysvf_types::Y_BANDPASS:
|
|
ybandpass_process_block(y.bandpass, inputs, outputs, block_size);
|
|
break;
|
|
case ysvf_types::Y_NOTCH:
|
|
ynotch_process_block(y.notch, inputs, outputs, block_size);
|
|
break;
|
|
}
|
|
}
|
|
} // namespace trnr
|