Skip to content

Commit

Permalink
update biquad
Browse files Browse the repository at this point in the history
  • Loading branch information
@jakubfi
jakubfi committed Jun 2, 2024
1 parent 076bf75 commit 95578a1
Show file tree
Hide file tree
Showing 2 changed files with 202 additions and 27 deletions.
219 changes: 194 additions & 25 deletions src/external/biquad/biquad.c
View file
Original file line number Diff line number Diff line change
@@ -1,10 +1,10 @@
// code taken from, stripped down and simplified for EM400:
//
// (c) Copyright 2016, Sean Connelly (@voidqk), http://syntheti.cc
// Simplified for EM400. Taken from:

// (c) Copyright 2016, Sean Connelly (@velipso), https://sean.cm
// MIT License
// Project Home: https://github.com/voidqk/sndfilter
// Project Home: https://github.com/velipso/sndfilter

#include "external/biquad/biquad.h"
#include "biquad.h"
#include <math.h>

// biquad filtering is based on a small sliding window, where the different filters are a result of
Expand All @@ -14,9 +14,8 @@
// b0, b1, b2, a1, a2 transformation coefficients
// xn0, xn1, xn2 the unfiltered sample at position x[n], x[n-1], and x[n-2]
// yn1, yn2 the filtered sample at position y[n-1] and y[n-2]
void sf_biquad_process(sf_biquad_state_st *state, int size, float *input,
float *output){

void sf_biquad_process(sf_biquad_state_st *state, int size, float *input, float *output)
{
// pull out the state into local variables
float b0 = state->b0;
float b1 = state->b1;
Expand All @@ -29,21 +28,18 @@ void sf_biquad_process(sf_biquad_state_st *state, int size, float *input,
float yn2 = state->yn2;

// loop for each sample
for (int n = 0; n < size; n++){
for (int n = 0; n < size; n++) {
// get the current sample
float xn0 = input[n];

// the formula is the same for each channel
float L =
output[n] =
b0 * xn0 +
b1 * xn1 +
b2 * xn2 -
a1 * yn1 -
a2 * yn2;

// save the result
output[n] = (float){ L };

// slide everything down one sample
xn2 = xn1;
xn1 = xn0;
Expand All @@ -70,15 +66,17 @@ void sf_biquad_process(sf_biquad_state_st *state, int size, float *input,
// https://git.io/v10H2

// clear the samples saved across process boundaries
static inline void state_reset(sf_biquad_state_st *state){
static inline void state_reset(sf_biquad_state_st *state)
{
state->xn1 = 0;
state->xn2 = 0;
state->yn1 = 0;
state->yn2 = 0;
}

// set the coefficients so that the output is the input scaled by `amt`
static inline void state_scale(sf_biquad_state_st *state, float amt){
static inline void state_scale(sf_biquad_state_st *state, float amt)
{
state->b0 = amt;
state->b1 = 0.0f;
state->b2 = 0.0f;
Expand All @@ -87,26 +85,29 @@ static inline void state_scale(sf_biquad_state_st *state, float amt){
}

// set the coefficients so that the output is an exact copy of the input
static inline void state_passthrough(sf_biquad_state_st *state){
static inline void state_passthrough(sf_biquad_state_st *state)
{
state_scale(state, 1.0f);
}

// set the coefficients so that the output is zeroed out
static inline void state_zero(sf_biquad_state_st *state){
static inline void state_zero(sf_biquad_state_st *state)
{
state_scale(state, 0.0f);
}

// initialize the biquad state to be a lowpass filter
void sf_lowpass(sf_biquad_state_st *state, int rate, float cutoff, float resonance){
void sf_lowpass(sf_biquad_state_st *state, int rate, float cutoff, float resonance)
{
state_reset(state);
float nyquist = rate * 0.5f;
cutoff /= nyquist;

if (cutoff >= 1.0f)
if (cutoff >= 1.0f) {
state_passthrough(state);
else if (cutoff <= 0.0f)
} else if (cutoff <= 0.0f) {
state_zero(state);
else{
} else {
resonance = powf(10.0f, resonance * 0.05f); // convert resonance from dB to linear
float theta = (float)M_PI * 2.0f * cutoff;
float alpha = sinf(theta) / (2.0f * resonance);
Expand All @@ -121,16 +122,17 @@ void sf_lowpass(sf_biquad_state_st *state, int rate, float cutoff, float resonan
}
}

void sf_highpass(sf_biquad_state_st *state, int rate, float cutoff, float resonance){
void sf_highpass(sf_biquad_state_st *state, int rate, float cutoff, float resonance)
{
state_reset(state);
float nyquist = rate * 0.5f;
cutoff /= nyquist;

if (cutoff >= 1.0f)
if (cutoff >= 1.0f) {
state_zero(state);
else if (cutoff <= 0.0f)
} else if (cutoff <= 0.0f) {
state_passthrough(state);
else{
} else {
resonance = powf(10.0f, resonance * 0.05f); // convert resonance from dB to linear
float theta = (float)M_PI * 2.0f * cutoff;
float alpha = sinf(theta) / (2.0f * resonance);
Expand All @@ -145,3 +147,170 @@ void sf_highpass(sf_biquad_state_st *state, int rate, float cutoff, float resona
}
}

void sf_bandpass(sf_biquad_state_st *state, int rate, float freq, float Q)
{
state_reset(state);
float nyquist = rate * 0.5f;
freq /= nyquist;

if (freq <= 0.0f || freq >= 1.0f) {
state_zero(state);
} else if (Q <= 0.0f) {
state_passthrough(state);
} else {
float w0 = (float)M_PI * 2.0f * freq;
float alpha = sinf(w0) / (2.0f * Q);
float k = cosf(w0);
float a0inv = 1.0f / (1.0f + alpha);
state->b0 = a0inv * alpha;
state->b1 = 0;
state->b2 = a0inv * -alpha;
state->a1 = a0inv * -2.0f * k;
state->a2 = a0inv * (1.0f - alpha);
}
}

void sf_notch(sf_biquad_state_st *state, int rate, float freq, float Q)
{
state_reset(state);
float nyquist = rate * 0.5f;
freq /= nyquist;

if (freq <= 0.0f || freq >= 1.0f) {
state_passthrough(state);
} else if (Q <= 0.0f) {
state_zero(state);
} else {
float w0 = (float)M_PI * 2.0f * freq;
float alpha = sinf(w0) / (2.0f * Q);
float k = cosf(w0);
float a0inv = 1.0f / (1.0f + alpha);
state->b0 = a0inv;
state->b1 = a0inv * -2.0f * k;
state->b2 = a0inv;
state->a1 = a0inv * -2.0f * k;
state->a2 = a0inv * (1.0f - alpha);
}
}

void sf_peaking(sf_biquad_state_st *state, int rate, float freq, float Q, float gain)
{
state_reset(state);
float nyquist = rate * 0.5f;
freq /= nyquist;

if (freq <= 0.0f || freq >= 1.0f) {
state_passthrough(state);
return;
}

float A = powf(10.0f, gain * 0.025f); // square root of gain converted from dB to linear

if (Q <= 0.0f) {
state_scale(state, A * A); // scale by A squared
return;
}

float w0 = (float)M_PI * 2.0f * freq;
float alpha = sinf(w0) / (2.0f * Q);
float k = cosf(w0);
float a0inv = 1.0f / (1.0f + alpha / A);
state->b0 = a0inv * (1.0f + alpha * A);
state->b1 = a0inv * -2.0f * k;
state->b2 = a0inv * (1.0f - alpha * A);
state->a1 = a0inv * -2.0f * k;
state->a2 = a0inv * (1.0f - alpha / A);
}

void sf_allpass(sf_biquad_state_st *state, int rate, float freq, float Q)
{
state_reset(state);
float nyquist = rate * 0.5f;
freq /= nyquist;

if (freq <= 0.0f || freq >= 1.0f) {
state_passthrough(state);
} else if (Q <= 0.0f) {
state_scale(state, -1.0f); // invert the sample
} else {
float w0 = (float)M_PI * 2.0f * freq;
float alpha = sinf(w0) / (2.0f * Q);
float k = cosf(w0);
float a0inv = 1.0f / (1.0f + alpha);
state->b0 = a0inv * (1.0f - alpha);
state->b1 = a0inv * -2.0f * k;
state->b2 = a0inv * (1.0f + alpha);
state->a1 = a0inv * -2.0f * k;
state->a2 = a0inv * (1.0f - alpha);
}
}

// WebAudio hardcodes Q=1
void sf_lowshelf(sf_biquad_state_st *state, int rate, float freq, float Q, float gain)
{
state_reset(state);
float nyquist = rate * 0.5f;
freq /= nyquist;

if (freq <= 0.0f || Q == 0.0f) {
state_passthrough(state);
return;
}

float A = powf(10.0f, gain * 0.025f); // square root of gain converted from dB to linear

if (freq >= 1.0f) {
state_scale(state, A * A); // scale by A squared
return;
}

float w0 = (float)M_PI * 2.0f * freq;
float ainn = (A + 1.0f / A) * (1.0f / Q - 1.0f) + 2.0f;
if (ainn < 0) ainn = 0;
float alpha = 0.5f * sinf(w0) * sqrtf(ainn);
float k = cosf(w0);
float k2 = 2.0f * sqrtf(A) * alpha;
float Ap1 = A + 1.0f;
float Am1 = A - 1.0f;
float a0inv = 1.0f / (Ap1 + Am1 * k + k2);
state->b0 = a0inv * A * (Ap1 - Am1 * k + k2);
state->b1 = a0inv * 2.0f * A * (Am1 - Ap1 * k);
state->b2 = a0inv * A * (Ap1 - Am1 * k - k2);
state->a1 = a0inv * -2.0f * (Am1 + Ap1 * k);
state->a2 = a0inv * (Ap1 + Am1 * k - k2);
}

// WebAudio hardcodes Q=1
void sf_highshelf(sf_biquad_state_st *state, int rate, float freq, float Q, float gain)
{
state_reset(state);
float nyquist = rate * 0.5f;
freq /= nyquist;

if (freq >= 1.0f || Q == 0.0f) {
state_passthrough(state);
return;
}

float A = powf(10.0f, gain * 0.025f); // square root of gain converted from dB to linear

if (freq <= 0.0f) {
state_scale(state, A * A); // scale by A squared
return;
}

float w0 = (float)M_PI * 2.0f * freq;
float ainn = (A + 1.0f / A) * (1.0f / Q - 1.0f) + 2.0f;
if (ainn < 0) ainn = 0;
float alpha = 0.5f * sinf(w0) * sqrtf(ainn);
float k = cosf(w0);
float k2 = 2.0f * sqrtf(A) * alpha;
float Ap1 = A + 1.0f;
float Am1 = A - 1.0f;
float a0inv = 1.0f / (Ap1 - Am1 * k + k2);
state->b0 = a0inv * A * (Ap1 + Am1 * k + k2);
state->b1 = a0inv * -2.0f * A * (Am1 + Ap1 * k);
state->b2 = a0inv * A * (Ap1 + Am1 * k - k2);
state->a1 = a0inv * 2.0f * (Am1 - Ap1 * k);
state->a2 = a0inv * (Ap1 - Am1 * k - k2);
}
10 changes: 8 additions & 2 deletions src/external/biquad/biquad.h
View file
Original file line number Diff line number Diff line change
@@ -1,6 +1,6 @@
// (c) Copyright 2016, Sean Connelly (@voidqk), http://syntheti.cc
// (c) Copyright 2016, Sean Connelly (@velipso), https://sean.cm
// MIT License
// Project Home: https://github.com/voidqk/sndfilter
// Project Home: https://github.com/velipso/sndfilter

// biquad filtering based on WebAudio specification:
// https://webaudio.github.io/web-audio-api/#the-biquadfilternode-interface
Expand Down Expand Up @@ -42,6 +42,12 @@ typedef struct {
// these functions will initialize an sf_biquad_state_st structure based on the desired filter
void sf_lowpass (sf_biquad_state_st *state, int rate, float cutoff, float resonance);
void sf_highpass (sf_biquad_state_st *state, int rate, float cutoff, float resonance);
void sf_bandpass (sf_biquad_state_st *state, int rate, float freq, float Q);
void sf_notch (sf_biquad_state_st *state, int rate, float freq, float Q);
void sf_peaking (sf_biquad_state_st *state, int rate, float freq, float Q, float gain);
void sf_allpass (sf_biquad_state_st *state, int rate, float freq, float Q);
void sf_lowshelf (sf_biquad_state_st *state, int rate, float freq, float Q, float gain);
void sf_highshelf(sf_biquad_state_st *state, int rate, float freq, float Q, float gain);

// this function will process the input sound based on the state passed
// the input and output buffers should be the same size
Expand Down

0 comments on commit 95578a1

Please sign in to comment.