[pulseaudio-discuss] [PATCH 1/8] lfe-filter: Import code from the Chrome OS audio server

David Henningsson david.henningsson at canonical.com
Tue Mar 24 02:29:12 PDT 2015


The chrome OS audio server has some already existing code, which
has been made available under a BSD-style license, which should be
safe to import by us.

Signed-off-by: David Henningsson <david.henningsson at canonical.com>
---
 LICENSE                             |   3 +
 src/pulsecore/filter/LICENSE.WEBKIT |  27 +++
 src/pulsecore/filter/biquad.c       | 368 ++++++++++++++++++++++++++++++++++++
 src/pulsecore/filter/biquad.h       |  57 ++++++
 src/pulsecore/filter/crossover.c    | 188 ++++++++++++++++++
 src/pulsecore/filter/crossover.h    |  70 +++++++
 6 files changed, 713 insertions(+)
 create mode 100644 src/pulsecore/filter/LICENSE.WEBKIT
 create mode 100644 src/pulsecore/filter/biquad.c
 create mode 100644 src/pulsecore/filter/biquad.h
 create mode 100644 src/pulsecore/filter/crossover.c
 create mode 100644 src/pulsecore/filter/crossover.h

diff --git a/LICENSE b/LICENSE
index 226c4ce..6932317 100644
--- a/LICENSE
+++ b/LICENSE
@@ -29,6 +29,9 @@ considered too small and stable to be considered as an external library) use the
 more permissive MIT license. This include the device reservation DBus protocol
 and realtime kit implementations.
 
+A more permissive BSD-style license is used for LFE filters, see
+src/pulsecore/filter/LICENSE.WEBKIT for details.
+
 Additionally, a more permissive Sun license is used for code that performs
 u-law, A-law and linear PCM conversions.
 
diff --git a/src/pulsecore/filter/LICENSE.WEBKIT b/src/pulsecore/filter/LICENSE.WEBKIT
new file mode 100644
index 0000000..2f69d9f
--- /dev/null
+++ b/src/pulsecore/filter/LICENSE.WEBKIT
@@ -0,0 +1,27 @@
+/*
+ * Copyright (C) 2010 Google Inc. All rights reserved.
+ *
+ * Redistribution and use in source and binary forms, with or without
+ * modification, are permitted provided that the following conditions
+ * are met:
+ *
+ * 1.  Redistributions of source code must retain the above copyright
+ *     notice, this list of conditions and the following disclaimer.
+ * 2.  Redistributions in binary form must reproduce the above copyright
+ *     notice, this list of conditions and the following disclaimer in the
+ *     documentation and/or other materials provided with the distribution.
+ * 3.  Neither the name of Apple Computer, Inc. ("Apple") nor the names of
+ *     its contributors may be used to endorse or promote products derived
+ *     from this software without specific prior written permission.
+ *
+ * THIS SOFTWARE IS PROVIDED BY APPLE AND ITS CONTRIBUTORS "AS IS" AND ANY
+ * EXPRESS OR IMPLIED WARRANTIES, INCLUDING, BUT NOT LIMITED TO, THE IMPLIED
+ * WARRANTIES OF MERCHANTABILITY AND FITNESS FOR A PARTICULAR PURPOSE ARE
+ * DISCLAIMED. IN NO EVENT SHALL APPLE OR ITS CONTRIBUTORS BE LIABLE FOR ANY
+ * DIRECT, INDIRECT, INCIDENTAL, SPECIAL, EXEMPLARY, OR CONSEQUENTIAL DAMAGES
+ * (INCLUDING, BUT NOT LIMITED TO, PROCUREMENT OF SUBSTITUTE GOODS OR SERVICES;
+ * LOSS OF USE, DATA, OR PROFITS; OR BUSINESS INTERRUPTION) HOWEVER CAUSED AND
+ * ON ANY THEORY OF LIABILITY, WHETHER IN CONTRACT, STRICT LIABILITY, OR TORT
+ * (INCLUDING NEGLIGENCE OR OTHERWISE) ARISING IN ANY WAY OUT OF THE USE OF
+ * THIS SOFTWARE, EVEN IF ADVISED OF THE POSSIBILITY OF SUCH DAMAGE.
+ */
diff --git a/src/pulsecore/filter/biquad.c b/src/pulsecore/filter/biquad.c
new file mode 100644
index 0000000..b28256d
--- /dev/null
+++ b/src/pulsecore/filter/biquad.c
@@ -0,0 +1,368 @@
+/* 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 max
+#define max(a, b) ({ __typeof__(a) _a = (a);	\
+			__typeof__(b) _b = (b);	\
+			_a > _b ? _a : _b; })
+#endif
+
+#ifndef min
+#define min(a, b) ({ __typeof__(a) _a = (a);	\
+			__typeof__(b) _b = (b);	\
+			_a < _b ? _a : _b; })
+#endif
+
+#ifndef M_PI
+#define M_PI 3.14159265358979323846
+#endif
+
+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 = b0 * a0_inv;
+	bq->b1 = b1 * a0_inv;
+	bq->b2 = b2 * a0_inv;
+	bq->a1 = a1 * a0_inv;
+	bq->a2 = a2 * a0_inv;
+}
+
+static void biquad_lowpass(struct biquad *bq, double cutoff, double resonance)
+{
+	/* Limit cutoff to 0 to 1. */
+	cutoff = max(0.0, min(cutoff, 1.0));
+
+	if (cutoff == 1) {
+		/* When cutoff is 1, the z-transform is 1. */
+		set_coefficient(bq, 1, 0, 0, 1, 0, 0);
+	} else if (cutoff > 0) {
+		/* Compute biquad coefficients for lowpass filter */
+		resonance = max(0.0, resonance); /* can't go negative */
+		double g = pow(10.0, 0.05 * resonance);
+		double d = sqrt((4 - sqrt(16 - 16 / (g * g))) / 2);
+
+		double theta = M_PI * cutoff;
+		double sn = 0.5 * d * sin(theta);
+		double beta = 0.5 * (1 - sn) / (1 + sn);
+		double gamma = (0.5 + beta) * cos(theta);
+		double alpha = 0.25 * (0.5 + beta - gamma);
+
+		double b0 = 2 * alpha;
+		double b1 = 2 * 2 * alpha;
+		double b2 = 2 * alpha;
+		double a1 = 2 * -gamma;
+		double a2 = 2 * beta;
+
+		set_coefficient(bq, b0, b1, b2, 1, a1, a2);
+	} else {
+		/* When cutoff is zero, nothing gets through the filter, so set
+		 * coefficients up correctly.
+		 */
+		set_coefficient(bq, 0, 0, 0, 1, 0, 0);
+	}
+}
+
+static void biquad_highpass(struct biquad *bq, double cutoff, double resonance)
+{
+	/* Limit cutoff to 0 to 1. */
+	cutoff = max(0.0, min(cutoff, 1.0));
+
+	if (cutoff == 1) {
+		/* The z-transform is 0. */
+		set_coefficient(bq, 0, 0, 0, 1, 0, 0);
+	} else if (cutoff > 0) {
+		/* Compute biquad coefficients for highpass filter */
+		resonance = max(0.0, resonance); /* can't go negative */
+		double g = pow(10.0, 0.05 * resonance);
+		double d = sqrt((4 - sqrt(16 - 16 / (g * g))) / 2);
+
+		double theta = M_PI * cutoff;
+		double sn = 0.5 * d * sin(theta);
+		double beta = 0.5 * (1 - sn) / (1 + sn);
+		double gamma = (0.5 + beta) * cos(theta);
+		double alpha = 0.25 * (0.5 + beta + gamma);
+
+		double b0 = 2 * alpha;
+		double b1 = 2 * -2 * alpha;
+		double b2 = 2 * alpha;
+		double a1 = 2 * -gamma;
+		double a2 = 2 * beta;
+
+		set_coefficient(bq, b0, b1, b2, 1, a1, a2);
+	} else {
+		/* 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, 0, 0, 1, 0, 0);
+	}
+}
+
+static void biquad_bandpass(struct biquad *bq, double frequency, double Q)
+{
+	/* No negative frequencies allowed. */
+	frequency = max(0.0, frequency);
+
+	/* Don't let Q go negative, which causes an unstable filter. */
+	Q = max(0.0, Q);
+
+	if (frequency > 0 && frequency < 1) {
+		double w0 = M_PI * frequency;
+		if (Q > 0) {
+			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);
+		} else {
+			/* 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);
+		}
+	} else {
+		/* 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);
+	}
+}
+
+static void biquad_lowshelf(struct biquad *bq, double frequency, double db_gain)
+{
+	/* Clip frequencies to between 0 and 1, inclusive. */
+	frequency = max(0.0, min(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);
+	} else if (frequency > 0) {
+		double w0 = M_PI * frequency;
+		double S = 1; /* filter slope (1 is max value) */
+		double alpha = 0.5 * sin(w0) *
+			sqrt((A + 1 / A) * (1 / S - 1) + 2);
+		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);
+	} else {
+		/* When frequency is 0, the z-transform is 1. */
+		set_coefficient(bq, 1, 0, 0, 1, 0, 0);
+	}
+}
+
+static void biquad_highshelf(struct biquad *bq, double frequency,
+			     double db_gain)
+{
+	/* Clip frequencies to between 0 and 1, inclusive. */
+	frequency = max(0.0, min(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);
+	} else if (frequency > 0) {
+		double w0 = M_PI * frequency;
+		double S = 1; /* filter slope (1 is max value) */
+		double alpha = 0.5 * sin(w0) *
+			sqrt((A + 1 / A) * (1 / S - 1) + 2);
+		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);
+	} else {
+		/* When frequency = 0, the filter is just a gain, A^2. */
+		set_coefficient(bq, A * A, 0, 0, 1, 0, 0);
+	}
+}
+
+static void biquad_peaking(struct biquad *bq, double frequency, double Q,
+			   double db_gain)
+{
+	/* Clip frequencies to between 0 and 1, inclusive. */
+	frequency = max(0.0, min(frequency, 1.0));
+
+	/* Don't let Q go negative, which causes an unstable filter. */
+	Q = max(0.0, Q);
+
+	double A = pow(10.0, db_gain / 40);
+
+	if (frequency > 0 && frequency < 1) {
+		if (Q > 0) {
+			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);
+		} else {
+			/* 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);
+		}
+	} else {
+		/* When frequency is 0 or 1, the z-transform is 1. */
+		set_coefficient(bq, 1, 0, 0, 1, 0, 0);
+	}
+}
+
+static void biquad_notch(struct biquad *bq, double frequency, double Q)
+{
+	/* Clip frequencies to between 0 and 1, inclusive. */
+	frequency = max(0.0, min(frequency, 1.0));
+
+	/* Don't let Q go negative, which causes an unstable filter. */
+	Q = max(0.0, Q);
+
+	if (frequency > 0 && frequency < 1) {
+		if (Q > 0) {
+			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);
+		} else {
+			/* 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);
+		}
+	} else {
+		/* When frequency is 0 or 1, the z-transform is 1. */
+		set_coefficient(bq, 1, 0, 0, 1, 0, 0);
+	}
+}
+
+static void biquad_allpass(struct biquad *bq, double frequency, double Q)
+{
+	/* Clip frequencies to between 0 and 1, inclusive. */
+	frequency = max(0.0, min(frequency, 1.0));
+
+	/* Don't let Q go negative, which causes an unstable filter. */
+	Q = max(0.0, Q);
+
+	if (frequency > 0 && frequency < 1) {
+		if (Q > 0) {
+			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);
+		} else {
+			/* 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);
+		}
+	} else {
+		/* When frequency is 0 or 1, the z-transform is 1. */
+		set_coefficient(bq, 1, 0, 0, 1, 0, 0);
+	}
+}
+
+void biquad_set(struct biquad *bq, enum biquad_type type, double freq, double Q,
+		double gain)
+{
+	/* Default is an identity filter. Also clear history values. */
+	set_coefficient(bq, 1, 0, 0, 1, 0, 0);
+	bq->x1 = 0;
+	bq->x2 = 0;
+	bq->y1 = 0;
+	bq->y2 = 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, gain);
+		break;
+	case BQ_HIGHSHELF:
+		biquad_highshelf(bq, freq, 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:
+		break;
+	}
+}
diff --git a/src/pulsecore/filter/biquad.h b/src/pulsecore/filter/biquad.h
new file mode 100644
index 0000000..c584aa9
--- /dev/null
+++ b/src/pulsecore/filter/biquad.h
@@ -0,0 +1,57 @@
+/* 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.
+ */
+
+#ifndef BIQUAD_H_
+#define BIQUAD_H_
+
+#ifdef __cplusplus
+extern "C" {
+#endif
+
+/* The biquad filter parameters. The transfer function H(z) is (b0 + b1 * z^(-1)
+ * + b2 * z^(-2)) / (1 + a1 * z^(-1) + a2 * z^(-2)).  The previous two inputs
+ * are stored in x1 and x2, and the previous two outputs are stored in y1 and
+ * y2.
+ *
+ * We use double during the coefficients calculation for better accurary, but
+ * float is used during the actual filtering for faster computation.
+ */
+struct biquad {
+	float b0, b1, b2;
+	float a1, a2;
+	float x1, x2;
+	float y1, y2;
+};
+
+/* The type of the biquad filters */
+enum biquad_type {
+	BQ_NONE,
+	BQ_LOWPASS,
+	BQ_HIGHPASS,
+	BQ_BANDPASS,
+	BQ_LOWSHELF,
+	BQ_HIGHSHELF,
+	BQ_PEAKING,
+	BQ_NOTCH,
+	BQ_ALLPASS
+};
+
+/* Initialize a biquad filter parameters from its type and parameters.
+ * Args:
+ *    bq - The biquad filter we want to set.
+ *    type - The type of the biquad filter.
+ *    frequency - The value should be in the range [0, 1]. It is relative to
+ *        half of the sampling rate.
+ *    Q - Quality factor. See Web Audio API for details.
+ *    gain - The value is in dB. See Web Audio API for details.
+ */
+void biquad_set(struct biquad *bq, enum biquad_type type, double freq, double Q,
+		double gain);
+
+#ifdef __cplusplus
+} /* extern "C" */
+#endif
+
+#endif /* BIQUAD_H_ */
diff --git a/src/pulsecore/filter/crossover.c b/src/pulsecore/filter/crossover.c
new file mode 100644
index 0000000..11a8c6e
--- /dev/null
+++ b/src/pulsecore/filter/crossover.c
@@ -0,0 +1,188 @@
+/* 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.
+ */
+
+#include "crossover.h"
+#include "biquad.h"
+
+static void lr4_set(struct lr4 *lr4, enum biquad_type type, float freq)
+{
+	struct biquad q;
+	biquad_set(&q, type, freq, 0, 0);
+	lr4->b0 = q.b0;
+	lr4->b1 = q.b1;
+	lr4->b2 = q.b2;
+	lr4->a1 = q.a1;
+	lr4->a2 = q.a2;
+	lr4->x1 = 0;
+	lr4->x2 = 0;
+	lr4->y1 = 0;
+	lr4->y2 = 0;
+	lr4->z1 = 0;
+	lr4->z2 = 0;
+}
+
+/* Split input data using two LR4 filters, put the result into the input array
+ * and another array.
+ *
+ * data0 --+-- lp --> data0
+ *         |
+ *         \-- hp --> data1
+ */
+static void lr4_split(struct lr4 *lp, struct lr4 *hp, int count, float *data0,
+		      float *data1)
+{
+	float lx1 = lp->x1;
+	float lx2 = lp->x2;
+	float ly1 = lp->y1;
+	float ly2 = lp->y2;
+	float lz1 = lp->z1;
+	float lz2 = lp->z2;
+	float lb0 = lp->b0;
+	float lb1 = lp->b1;
+	float lb2 = lp->b2;
+	float la1 = lp->a1;
+	float la2 = lp->a2;
+
+	float hx1 = hp->x1;
+	float hx2 = hp->x2;
+	float hy1 = hp->y1;
+	float hy2 = hp->y2;
+	float hz1 = hp->z1;
+	float hz2 = hp->z2;
+	float hb0 = hp->b0;
+	float hb1 = hp->b1;
+	float hb2 = hp->b2;
+	float ha1 = hp->a1;
+	float ha2 = hp->a2;
+
+	int i;
+	for (i = 0; i < count; i++) {
+		float x, y, z;
+		x = data0[i];
+		y = lb0*x + lb1*lx1 + lb2*lx2 - la1*ly1 - la2*ly2;
+		z = lb0*y + lb1*ly1 + lb2*ly2 - la1*lz1 - la2*lz2;
+		lx2 = lx1;
+		lx1 = x;
+		ly2 = ly1;
+		ly1 = y;
+		lz2 = lz1;
+		lz1 = z;
+		data0[i] = z;
+
+		y = hb0*x + hb1*hx1 + hb2*hx2 - ha1*hy1 - ha2*hy2;
+		z = hb0*y + hb1*hy1 + hb2*hy2 - ha1*hz1 - ha2*hz2;
+		hx2 = hx1;
+		hx1 = x;
+		hy2 = hy1;
+		hy1 = y;
+		hz2 = hz1;
+		hz1 = z;
+		data1[i] = z;
+	}
+
+	lp->x1 = lx1;
+	lp->x2 = lx2;
+	lp->y1 = ly1;
+	lp->y2 = ly2;
+	lp->z1 = lz1;
+	lp->z2 = lz2;
+
+	hp->x1 = hx1;
+	hp->x2 = hx2;
+	hp->y1 = hy1;
+	hp->y2 = hy2;
+	hp->z1 = hz1;
+	hp->z2 = hz2;
+}
+
+/* Split input data using two LR4 filters and sum them back to the original
+ * data array.
+ *
+ * data --+-- lp --+--> data
+ *        |        |
+ *        \-- hp --/
+ */
+static void lr4_merge(struct lr4 *lp, struct lr4 *hp, int count, float *data)
+{
+	float lx1 = lp->x1;
+	float lx2 = lp->x2;
+	float ly1 = lp->y1;
+	float ly2 = lp->y2;
+	float lz1 = lp->z1;
+	float lz2 = lp->z2;
+	float lb0 = lp->b0;
+	float lb1 = lp->b1;
+	float lb2 = lp->b2;
+	float la1 = lp->a1;
+	float la2 = lp->a2;
+
+	float hx1 = hp->x1;
+	float hx2 = hp->x2;
+	float hy1 = hp->y1;
+	float hy2 = hp->y2;
+	float hz1 = hp->z1;
+	float hz2 = hp->z2;
+	float hb0 = hp->b0;
+	float hb1 = hp->b1;
+	float hb2 = hp->b2;
+	float ha1 = hp->a1;
+	float ha2 = hp->a2;
+
+	int i;
+	for (i = 0; i < count; i++) {
+		float x, y, z;
+		x = data[i];
+		y = lb0*x + lb1*lx1 + lb2*lx2 - la1*ly1 - la2*ly2;
+		z = lb0*y + lb1*ly1 + lb2*ly2 - la1*lz1 - la2*lz2;
+		lx2 = lx1;
+		lx1 = x;
+		ly2 = ly1;
+		ly1 = y;
+		lz2 = lz1;
+		lz1 = z;
+
+		y = hb0*x + hb1*hx1 + hb2*hx2 - ha1*hy1 - ha2*hy2;
+		z = hb0*y + hb1*hy1 + hb2*hy2 - ha1*hz1 - ha2*hz2;
+		hx2 = hx1;
+		hx1 = x;
+		hy2 = hy1;
+		hy1 = y;
+		hz2 = hz1;
+		hz1 = z;
+		data[i] = z + lz1;
+	}
+
+	lp->x1 = lx1;
+	lp->x2 = lx2;
+	lp->y1 = ly1;
+	lp->y2 = ly2;
+	lp->z1 = lz1;
+	lp->z2 = lz2;
+
+	hp->x1 = hx1;
+	hp->x2 = hx2;
+	hp->y1 = hy1;
+	hp->y2 = hy2;
+	hp->z1 = hz1;
+	hp->z2 = hz2;
+}
+
+void crossover_init(struct crossover *xo, float freq1, float freq2)
+{
+	int i;
+	for (i = 0; i < 3; i++) {
+		float f = (i == 0) ? freq1 : freq2;
+		lr4_set(&xo->lp[i], BQ_LOWPASS, f);
+		lr4_set(&xo->hp[i], BQ_HIGHPASS, f);
+	}
+}
+
+void crossover_process(struct crossover *xo, int count, float *data0,
+		       float *data1, float *data2)
+{
+	lr4_split(&xo->lp[0], &xo->hp[0], count, data0, data1);
+	lr4_merge(&xo->lp[1], &xo->hp[1], count, data0);
+	lr4_split(&xo->lp[2], &xo->hp[2], count, data1, data2);
+}
diff --git a/src/pulsecore/filter/crossover.h b/src/pulsecore/filter/crossover.h
new file mode 100644
index 0000000..99a601c
--- /dev/null
+++ b/src/pulsecore/filter/crossover.h
@@ -0,0 +1,70 @@
+/* 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.
+ */
+
+#ifndef CROSSOVER_H_
+#define CROSSOVER_H_
+
+#ifdef __cplusplus
+extern "C" {
+#endif
+
+/* An LR4 filter is two biquads with the same parameters connected in series:
+ *
+ * x -- [BIQUAD] -- y -- [BIQUAD] -- z
+ *
+ * Both biquad filter has the same parameter b[012] and a[12],
+ * The variable [xyz][12] keep the history values.
+ */
+struct lr4 {
+	float b0, b1, b2;
+	float a1, a2;
+	float x1, x2;
+	float y1, y2;
+	float z1, z2;
+};
+
+/* Three bands crossover filter:
+ *
+ * INPUT --+-- lp0 --+-- lp1 --+---> LOW (0)
+ *         |         |         |
+ *         |         \-- hp1 --/
+ *         |
+ *         \-- hp0 --+-- lp2 ------> MID (1)
+ *                   |
+ *                   \-- hp2 ------> HIGH (2)
+ *
+ *            [f0]       [f1]
+ *
+ * Each lp or hp is an LR4 filter, which consists of two second-order
+ * lowpass or highpass butterworth filters.
+ */
+struct crossover {
+	struct lr4 lp[3], hp[3];
+};
+
+/* Initializes a crossover filter
+ * Args:
+ *    xo - The crossover filter we want to initialize.
+ *    freq1 - The normalized frequency splits low and mid band.
+ *    freq2 - The normalized frequency splits mid and high band.
+ */
+void crossover_init(struct crossover *xo, float freq1, float freq2);
+
+/* Splits input samples to three bands.
+ * Args:
+ *    xo - The crossover filter to use.
+ *    count - The number of input samples.
+ *    data0 - The input samples, also the place to store low band output.
+ *    data1 - The place to store mid band output.
+ *    data2 - The place to store high band output.
+ */
+void crossover_process(struct crossover *xo, int count, float *data0,
+		       float *data1, float *data2);
+
+#ifdef __cplusplus
+} /* extern "C" */
+#endif
+
+#endif /* CROSSOVER_H_ */
-- 
1.9.1



More information about the pulseaudio-discuss mailing list