mirror of
https://gitlab.freedesktop.org/pipewire/pipewire.git
synced 2025-10-31 22:25:38 -04:00
7785f61677
The current biquad calculations are based on RBJ's cookbook [1], except for low-/highpass. Since the filter configuration is also based on using the definition of Q, it makes sense to also align the remaining calculations to use the same filter cookbook instead of using resonance which doesn't result in the same coefficients as when using Q. [1] = https://www.w3.org/TR/audio-eq-cookbook/
374 lines
9.5 KiB
C
374 lines
9.5 KiB
C
/* Copyright (c) 2013 The Chromium OS Authors. All rights reserved.
|
|
* Use of this source code is governed by a BSD-style license that can be
|
|
* found in the LICENSE file.
|
|
*/
|
|
|
|
/* Copyright (C) 2010 Google Inc. All rights reserved.
|
|
* Use of this source code is governed by a BSD-style license that can be
|
|
* found in the LICENSE.WEBKIT file.
|
|
*/
|
|
|
|
#include <math.h>
|
|
#include "biquad.h"
|
|
|
|
#ifndef M_PI
|
|
#define M_PI 3.14159265358979323846
|
|
#endif
|
|
|
|
/* Q = 1 / sqrt(2), also resulting Q value when S = 1 */
|
|
#define BIQUAD_DEFAULT_Q 0.707106781186548
|
|
|
|
static void set_coefficient(struct biquad *bq, double b0, double b1, double b2,
|
|
double a0, double a1, double a2)
|
|
{
|
|
double a0_inv = 1 / a0;
|
|
bq->b0 = (float)(b0 * a0_inv);
|
|
bq->b1 = (float)(b1 * a0_inv);
|
|
bq->b2 = (float)(b2 * a0_inv);
|
|
bq->a1 = (float)(a1 * a0_inv);
|
|
bq->a2 = (float)(a2 * a0_inv);
|
|
}
|
|
|
|
static void biquad_lowpass(struct biquad *bq, double cutoff, double Q)
|
|
{
|
|
/* Limit cutoff to 0 to 1. */
|
|
cutoff = fmax(0.0, fmin(cutoff, 1.0));
|
|
|
|
if (cutoff == 1 || cutoff == 0) {
|
|
/* When cutoff is 1, the z-transform is 1.
|
|
* When cutoff is zero, nothing gets through the filter, so set
|
|
* coefficients up correctly.
|
|
*/
|
|
set_coefficient(bq, cutoff, 0, 0, 1, 0, 0);
|
|
return;
|
|
}
|
|
|
|
/* Set Q to a sane default value if not set */
|
|
if (Q <= 0)
|
|
Q = BIQUAD_DEFAULT_Q;
|
|
|
|
/* Compute biquad coefficients for lowpass filter */
|
|
/* H(s) = 1 / (s^2 + s/Q + 1) */
|
|
double w0 = M_PI * cutoff;
|
|
double alpha = sin(w0) / (2 * Q);
|
|
double k = cos(w0);
|
|
|
|
double b0 = (1 - k) / 2;
|
|
double b1 = 1 - k;
|
|
double b2 = (1 - k) / 2;
|
|
double a0 = 1 + alpha;
|
|
double a1 = -2 * k;
|
|
double a2 = 1 - alpha;
|
|
|
|
set_coefficient(bq, b0, b1, b2, a0, a1, a2);
|
|
}
|
|
|
|
static void biquad_highpass(struct biquad *bq, double cutoff, double Q)
|
|
{
|
|
/* Limit cutoff to 0 to 1. */
|
|
cutoff = fmax(0.0, fmin(cutoff, 1.0));
|
|
|
|
if (cutoff == 1 || cutoff == 0) {
|
|
/* When cutoff is one, the z-transform is 0. */
|
|
/* When cutoff is zero, we need to be careful because the above
|
|
* gives a quadratic divided by the same quadratic, with poles
|
|
* and zeros on the unit circle in the same place. When cutoff
|
|
* is zero, the z-transform is 1.
|
|
*/
|
|
set_coefficient(bq, 1 - cutoff, 0, 0, 1, 0, 0);
|
|
return;
|
|
}
|
|
|
|
/* Set Q to a sane default value if not set */
|
|
if (Q <= 0)
|
|
Q = BIQUAD_DEFAULT_Q;
|
|
|
|
/* Compute biquad coefficients for highpass filter */
|
|
/* H(s) = s^2 / (s^2 + s/Q + 1) */
|
|
double w0 = M_PI * cutoff;
|
|
double alpha = sin(w0) / (2 * Q);
|
|
double k = cos(w0);
|
|
|
|
double b0 = (1 + k) / 2;
|
|
double b1 = -(1 + k);
|
|
double b2 = (1 + k) / 2;
|
|
double a0 = 1 + alpha;
|
|
double a1 = -2 * k;
|
|
double a2 = 1 - alpha;
|
|
|
|
set_coefficient(bq, b0, b1, b2, a0, a1, a2);
|
|
}
|
|
|
|
static void biquad_bandpass(struct biquad *bq, double frequency, double Q)
|
|
{
|
|
/* No negative frequencies allowed. */
|
|
frequency = fmax(0.0, frequency);
|
|
|
|
/* Don't let Q go negative, which causes an unstable filter. */
|
|
Q = fmax(0.0, Q);
|
|
|
|
if (frequency <= 0 || frequency >= 1) {
|
|
/* When the cutoff is zero, the z-transform approaches 0, if Q
|
|
* > 0. When both Q and cutoff are zero, the z-transform is
|
|
* pretty much undefined. What should we do in this case?
|
|
* For now, just make the filter 0. When the cutoff is 1, the
|
|
* z-transform also approaches 0.
|
|
*/
|
|
set_coefficient(bq, 0, 0, 0, 1, 0, 0);
|
|
return;
|
|
}
|
|
if (Q <= 0) {
|
|
/* When Q = 0, the above formulas have problems. If we
|
|
* look at the z-transform, we can see that the limit
|
|
* as Q->0 is 1, so set the filter that way.
|
|
*/
|
|
set_coefficient(bq, 1, 0, 0, 1, 0, 0);
|
|
return;
|
|
}
|
|
|
|
double w0 = M_PI * frequency;
|
|
double alpha = sin(w0) / (2 * Q);
|
|
double k = cos(w0);
|
|
|
|
double b0 = alpha;
|
|
double b1 = 0;
|
|
double b2 = -alpha;
|
|
double a0 = 1 + alpha;
|
|
double a1 = -2 * k;
|
|
double a2 = 1 - alpha;
|
|
|
|
set_coefficient(bq, b0, b1, b2, a0, a1, a2);
|
|
}
|
|
|
|
static void biquad_lowshelf(struct biquad *bq, double frequency, double Q,
|
|
double db_gain)
|
|
{
|
|
/* Clip frequencies to between 0 and 1, inclusive. */
|
|
frequency = fmax(0.0, fmin(frequency, 1.0));
|
|
|
|
double A = pow(10.0, db_gain / 40);
|
|
|
|
if (frequency == 1) {
|
|
/* The z-transform is a constant gain. */
|
|
set_coefficient(bq, A * A, 0, 0, 1, 0, 0);
|
|
return;
|
|
}
|
|
if (frequency <= 0) {
|
|
/* When frequency is 0, the z-transform is 1. */
|
|
set_coefficient(bq, 1, 0, 0, 1, 0, 0);
|
|
return;
|
|
}
|
|
|
|
/* Set Q to an equivalent value to S = 1 if not specified */
|
|
if (Q <= 0)
|
|
Q = BIQUAD_DEFAULT_Q;
|
|
|
|
double w0 = M_PI * frequency;
|
|
double alpha = sin(w0) / (2 * Q);
|
|
double k = cos(w0);
|
|
double k2 = 2 * sqrt(A) * alpha;
|
|
double a_plus_one = A + 1;
|
|
double a_minus_one = A - 1;
|
|
|
|
double b0 = A * (a_plus_one - a_minus_one * k + k2);
|
|
double b1 = 2 * A * (a_minus_one - a_plus_one * k);
|
|
double b2 = A * (a_plus_one - a_minus_one * k - k2);
|
|
double a0 = a_plus_one + a_minus_one * k + k2;
|
|
double a1 = -2 * (a_minus_one + a_plus_one * k);
|
|
double a2 = a_plus_one + a_minus_one * k - k2;
|
|
|
|
set_coefficient(bq, b0, b1, b2, a0, a1, a2);
|
|
}
|
|
|
|
static void biquad_highshelf(struct biquad *bq, double frequency, double Q,
|
|
double db_gain)
|
|
{
|
|
/* Clip frequencies to between 0 and 1, inclusive. */
|
|
frequency = fmax(0.0, fmin(frequency, 1.0));
|
|
|
|
double A = pow(10.0, db_gain / 40);
|
|
|
|
if (frequency == 1) {
|
|
/* The z-transform is 1. */
|
|
set_coefficient(bq, 1, 0, 0, 1, 0, 0);
|
|
return;
|
|
}
|
|
if (frequency <= 0) {
|
|
/* When frequency = 0, the filter is just a gain, A^2. */
|
|
set_coefficient(bq, A * A, 0, 0, 1, 0, 0);
|
|
return;
|
|
}
|
|
|
|
/* Set Q to an equivalent value to S = 1 if not specified */
|
|
if (Q <= 0)
|
|
Q = BIQUAD_DEFAULT_Q;
|
|
|
|
double w0 = M_PI * frequency;
|
|
double alpha = sin(w0) / (2 * Q);
|
|
double k = cos(w0);
|
|
double k2 = 2 * sqrt(A) * alpha;
|
|
double a_plus_one = A + 1;
|
|
double a_minus_one = A - 1;
|
|
|
|
double b0 = A * (a_plus_one + a_minus_one * k + k2);
|
|
double b1 = -2 * A * (a_minus_one + a_plus_one * k);
|
|
double b2 = A * (a_plus_one + a_minus_one * k - k2);
|
|
double a0 = a_plus_one - a_minus_one * k + k2;
|
|
double a1 = 2 * (a_minus_one - a_plus_one * k);
|
|
double a2 = a_plus_one - a_minus_one * k - k2;
|
|
|
|
set_coefficient(bq, b0, b1, b2, a0, a1, a2);
|
|
}
|
|
|
|
static void biquad_peaking(struct biquad *bq, double frequency, double Q,
|
|
double db_gain)
|
|
{
|
|
/* Clip frequencies to between 0 and 1, inclusive. */
|
|
frequency = fmax(0.0, fmin(frequency, 1.0));
|
|
|
|
/* Don't let Q go negative, which causes an unstable filter. */
|
|
Q = fmax(0.0, Q);
|
|
|
|
double A = pow(10.0, db_gain / 40);
|
|
|
|
if (frequency <= 0 || frequency >= 1) {
|
|
/* When frequency is 0 or 1, the z-transform is 1. */
|
|
set_coefficient(bq, 1, 0, 0, 1, 0, 0);
|
|
return;
|
|
}
|
|
if (Q <= 0) {
|
|
/* When Q = 0, the above formulas have problems. If we
|
|
* look at the z-transform, we can see that the limit
|
|
* as Q->0 is A^2, so set the filter that way.
|
|
*/
|
|
set_coefficient(bq, A * A, 0, 0, 1, 0, 0);
|
|
return;
|
|
}
|
|
|
|
double w0 = M_PI * frequency;
|
|
double alpha = sin(w0) / (2 * Q);
|
|
double k = cos(w0);
|
|
|
|
double b0 = 1 + alpha * A;
|
|
double b1 = -2 * k;
|
|
double b2 = 1 - alpha * A;
|
|
double a0 = 1 + alpha / A;
|
|
double a1 = -2 * k;
|
|
double a2 = 1 - alpha / A;
|
|
|
|
set_coefficient(bq, b0, b1, b2, a0, a1, a2);
|
|
}
|
|
|
|
static void biquad_notch(struct biquad *bq, double frequency, double Q)
|
|
{
|
|
/* Clip frequencies to between 0 and 1, inclusive. */
|
|
frequency = fmax(0.0, fmin(frequency, 1.0));
|
|
|
|
/* Don't let Q go negative, which causes an unstable filter. */
|
|
Q = fmax(0.0, Q);
|
|
|
|
if (frequency <= 0 || frequency >= 1) {
|
|
/* When frequency is 0 or 1, the z-transform is 1. */
|
|
set_coefficient(bq, 1, 0, 0, 1, 0, 0);
|
|
return;
|
|
}
|
|
if (Q <= 0) {
|
|
/* When Q = 0, the above formulas have problems. If we
|
|
* look at the z-transform, we can see that the limit
|
|
* as Q->0 is 0, so set the filter that way.
|
|
*/
|
|
set_coefficient(bq, 0, 0, 0, 1, 0, 0);
|
|
return;
|
|
}
|
|
|
|
double w0 = M_PI * frequency;
|
|
double alpha = sin(w0) / (2 * Q);
|
|
double k = cos(w0);
|
|
|
|
double b0 = 1;
|
|
double b1 = -2 * k;
|
|
double b2 = 1;
|
|
double a0 = 1 + alpha;
|
|
double a1 = -2 * k;
|
|
double a2 = 1 - alpha;
|
|
|
|
set_coefficient(bq, b0, b1, b2, a0, a1, a2);
|
|
}
|
|
|
|
static void biquad_allpass(struct biquad *bq, double frequency, double Q)
|
|
{
|
|
/* Clip frequencies to between 0 and 1, inclusive. */
|
|
frequency = fmax(0.0, fmin(frequency, 1.0));
|
|
|
|
/* Don't let Q go negative, which causes an unstable filter. */
|
|
Q = fmax(0.0, Q);
|
|
|
|
if (frequency <= 0 || frequency >= 1) {
|
|
/* When frequency is 0 or 1, the z-transform is 1. */
|
|
set_coefficient(bq, 1, 0, 0, 1, 0, 0);
|
|
return;
|
|
}
|
|
|
|
if (Q <= 0) {
|
|
/* When Q = 0, the above formulas have problems. If we
|
|
* look at the z-transform, we can see that the limit
|
|
* as Q->0 is -1, so set the filter that way.
|
|
*/
|
|
set_coefficient(bq, -1, 0, 0, 1, 0, 0);
|
|
return;
|
|
}
|
|
|
|
double w0 = M_PI * frequency;
|
|
double alpha = sin(w0) / (2 * Q);
|
|
double k = cos(w0);
|
|
|
|
double b0 = 1 - alpha;
|
|
double b1 = -2 * k;
|
|
double b2 = 1 + alpha;
|
|
double a0 = 1 + alpha;
|
|
double a1 = -2 * k;
|
|
double a2 = 1 - alpha;
|
|
|
|
set_coefficient(bq, b0, b1, b2, a0, a1, a2);
|
|
}
|
|
|
|
void biquad_set(struct biquad *bq, enum biquad_type type, double freq, double Q,
|
|
double gain)
|
|
{
|
|
/* Clear history values. */
|
|
bq->type = type;
|
|
bq->x1 = 0;
|
|
bq->x2 = 0;
|
|
|
|
switch (type) {
|
|
case BQ_LOWPASS:
|
|
biquad_lowpass(bq, freq, Q);
|
|
break;
|
|
case BQ_HIGHPASS:
|
|
biquad_highpass(bq, freq, Q);
|
|
break;
|
|
case BQ_BANDPASS:
|
|
biquad_bandpass(bq, freq, Q);
|
|
break;
|
|
case BQ_LOWSHELF:
|
|
biquad_lowshelf(bq, freq, Q, gain);
|
|
break;
|
|
case BQ_HIGHSHELF:
|
|
biquad_highshelf(bq, freq, Q, gain);
|
|
break;
|
|
case BQ_PEAKING:
|
|
biquad_peaking(bq, freq, Q, gain);
|
|
break;
|
|
case BQ_NOTCH:
|
|
biquad_notch(bq, freq, Q);
|
|
break;
|
|
case BQ_ALLPASS:
|
|
biquad_allpass(bq, freq, Q);
|
|
break;
|
|
case BQ_NONE:
|
|
case BQ_RAW:
|
|
/* Default is an identity filter. */
|
|
set_coefficient(bq, 1, 0, 0, 1, 0, 0);
|
|
break;
|
|
}
|
|
}
|