ref: b01bd4a8331c5e24f88bb0659a1f723fe7bf12ca
parent: 59c046dc1d15355548e24d6e29bae383782e656a
author: Paul Brossier <piem@piem.org>
date: Mon Oct 19 10:45:13 EDT 2009
src/temporal: derive biquad from filter, use in peakpicker
--- a/src/aubio.h
+++ b/src/aubio.h
@@ -71,8 +71,8 @@
#if HAVE_SAMPLERATE
#include "temporal/resampler.h"
#endif /* HAVE_SAMPLERATE */
-#include "temporal/biquad.h"
#include "temporal/filter.h"
+#include "temporal/biquad.h"
#include "temporal/a_weighting.h"
#include "temporal/c_weighting.h"
#include "spectral/filterbank.h"
--- a/src/onset/peakpick.c
+++ b/src/onset/peakpick.c
@@ -20,6 +20,8 @@
#include "aubio_priv.h"
#include "fvec.h"
#include "mathutils.h"
+#include "lvec.h"
+#include "temporal/filter.h"
#include "temporal/biquad.h"
#include "onset/peakpick.h"
@@ -42,7 +44,7 @@
aubio_pickerfn_t pickerfn;
/** biquad lowpass filter */
- aubio_biquad_t * biquad;
+ aubio_filter_t * biquad;
/** original onsets */
fvec_t * onset_keep;
/** modified onsets */
@@ -94,7 +96,7 @@
/* filter onset_proc */
/** \bug filtfilt calculated post+pre times, should be only once !? */
- //aubio_biquad_do_filtfilt(p->biquad,onset_proc,scratch);
+ aubio_filter_do_filtfilt (p->biquad, onset_proc, scratch);
for (i = 0; i < p->channels; i++) {
/* calculate mean and median for onset_proc */
@@ -159,14 +161,17 @@
t->onset_proc = new_fvec(t->win_post+t->win_pre+1, channels);
t->onset_peek = new_fvec(3, channels);
- /* cutoff: low-pass filter cutoff [0.34, 1] */
- /* t->cutoff=0.34; */
- t->biquad = new_aubio_biquad(0.1600,0.3200,0.1600,-0.5949,0.2348);
+ /* cutoff: low-pass filter with cutoff reduced frequency at 0.34
+ generated with octave butter function: [b,a] = butter(2, 0.34);
+ */
+ t->biquad = new_aubio_filter_biquad (0.15998789, 0.31997577, 0.15998789,
+ -0.59488894, 0.23484048, channels);
+
return t;
}
void del_aubio_peakpicker(aubio_peakpicker_t * p) {
- del_aubio_biquad(p->biquad);
+ del_aubio_filter(p->biquad);
del_fvec(p->onset_keep);
del_fvec(p->onset_proc);
del_fvec(p->onset_peek);
--- a/src/temporal/biquad.c
+++ b/src/temporal/biquad.c
@@ -19,96 +19,34 @@
#include "aubio_priv.h"
#include "fvec.h"
-#include "mathutils.h"
-#include "temporal/biquad.h"
+#include "lvec.h"
+#include "temporal/filter.h"
-/** \note this file needs to be in double or more less precision would lead to large
- * errors in the output
- */
+uint_t
+aubio_filter_set_biquad (aubio_filter_t * f, lsmp_t b0, lsmp_t b1, lsmp_t b2,
+ lsmp_t a1, lsmp_t a2) {
+ uint_t order = aubio_filter_get_order (f);
+ lvec_t *bs = aubio_filter_get_feedforward (f);
+ lvec_t *as = aubio_filter_get_feedback (f);
-struct _aubio_biquad_t {
- lsmp_t a2;
- lsmp_t a3;
- lsmp_t b1;
- lsmp_t b2;
- lsmp_t b3;
- lsmp_t o1;
- lsmp_t o2;
- lsmp_t i1;
- lsmp_t i2;
-};
-
-void aubio_biquad_do(aubio_biquad_t * b, fvec_t * in) {
- uint_t i,j;
- lsmp_t i1 = b->i1;
- lsmp_t i2 = b->i2;
- lsmp_t o1 = b->o1;
- lsmp_t o2 = b->o2;
- lsmp_t a2 = b->a2;
- lsmp_t a3 = b->a3;
- lsmp_t b1 = b->b1;
- lsmp_t b2 = b->b2;
- lsmp_t b3 = b->b3;
-
- i=0; // works in mono only !!!
- //for (i=0;i<in->channels;i++) {
- for (j = 0; j < in->length; j++) {
- lsmp_t i0 = in->data[i][j];
- lsmp_t o0 = b1 * i0 + b2 * i1 + b3 * i2
- - a2 * o1 - a3 * o2;// + 1e-37;
- in->data[i][j] = o0;
- i2 = i1;
- i1 = i0;
- o2 = o1;
- o1 = o0;
+ if (order != 3) {
+ AUBIO_ERROR ("order of biquad filter must be 3, not %d\n", order);
+ return AUBIO_FAIL;
}
- b->i2 = i2;
- b->i1 = i1;
- b->o2 = o2;
- b->o1 = o1;
- //}
+ bs->data[0][0] = b0;
+ bs->data[0][1] = b1;
+ bs->data[0][2] = b2;
+ as->data[0][0] = 1.;
+ as->data[0][1] = a1;
+ as->data[0][1] = a2;
+ return AUBIO_OK;
}
-void aubio_biquad_do_filtfilt(aubio_biquad_t * b, fvec_t * in, fvec_t * tmp) {
- uint_t j,i=0;
- uint_t length = in->length;
- lsmp_t mir;
- /* mirroring */
- mir = 2*in->data[i][0];
- b->i1 = mir - in->data[i][2];
- b->i2 = mir - in->data[i][1];
- /* apply filtering */
- aubio_biquad_do(b,in);
- /* invert */
- for (j = 0; j < length; j++)
- tmp->data[i][length-j-1] = in->data[i][j];
- /* mirror again */
- mir = 2*tmp->data[i][0];
- b->i1 = mir - tmp->data[i][2];
- b->i2 = mir - tmp->data[i][1];
- /* apply filtering */
- aubio_biquad_do(b,tmp);
- /* invert back */
- for (j = 0; j < length; j++)
- in->data[i][j] = tmp->data[i][length-j-1];
-}
-
-aubio_biquad_t * new_aubio_biquad(
- lsmp_t b1, lsmp_t b2, lsmp_t b3,
- lsmp_t a2, lsmp_t a3) {
- aubio_biquad_t * b = AUBIO_NEW(aubio_biquad_t);
- b->a2 = a2;
- b->a3 = a3;
- b->b1 = b1;
- b->b2 = b2;
- b->b3 = b3;
- b->i1 = 0.;
- b->i2 = 0.;
- b->o1 = 0.;
- b->o2 = 0.;
- return b;
-}
-
-void del_aubio_biquad(aubio_biquad_t * b) {
- AUBIO_FREE(b);
+aubio_filter_t *
+new_aubio_filter_biquad (lsmp_t b0, lsmp_t b1, lsmp_t b2, lsmp_t a1, lsmp_t a2,
+ uint_t channels)
+{
+ aubio_filter_t *f = new_aubio_filter (3, channels);
+ aubio_filter_set_biquad (f, b0, b1, b2, a1, a2);
+ return f;
}
--- a/src/temporal/biquad.h
+++ b/src/temporal/biquad.h
@@ -26,11 +26,14 @@
This file implements a normalised biquad filter (second order IIR):
- \f$ y[n] = b_1 x[n] + b_2 x[n-1] + b_3 x[n-2] - a_2 y[n-1] - a_3 y[n-2] \f$
+ \f$ y[n] = b_0 x[n] + b_1 x[n-1] + b_2 x[n-2] - a_1 y[n-1] - a_2 y[n-2] \f$
The filtfilt version runs the filter twice, forward and backward, to
compensate the phase shifting of the forward operation.
+ See also <a href="http://en.wikipedia.org/wiki/Digital_biquad_filter">Digital
+ biquad filter</a> on wikipedia.
+
*/
#ifdef __cplusplus
@@ -37,41 +40,31 @@
extern "C" {
#endif
-/** biquad filter object */
-typedef struct _aubio_biquad_t aubio_biquad_t;
+/** set coefficients of a biquad filter
-/** filter input vector
+ \param b0 forward filter coefficient
+ \param b1 forward filter coefficient
+ \param b2 forward filter coefficient
+ \param a1 feedback filter coefficient
+ \param a2 feedback filter coefficient
- \param b biquad object as returned by new_aubio_biquad
- \param in input vector to filter
-
*/
-void aubio_biquad_do(aubio_biquad_t * b, fvec_t * in);
-/** filter input vector forward and backward
+uint_t
+aubio_filter_set_biquad (aubio_filter_t * f, lsmp_t b0, lsmp_t b1, lsmp_t b2,
+ lsmp_t a1, lsmp_t a2);
- \param b biquad object as returned by new_aubio_biquad
- \param in input vector to filter
- \param tmp memory space to use for computation
-
-*/
-void aubio_biquad_do_filtfilt(aubio_biquad_t * b, fvec_t * in, fvec_t * tmp);
/** create new biquad filter
+ \param b0 forward filter coefficient
\param b1 forward filter coefficient
\param b2 forward filter coefficient
- \param b3 forward filter coefficient
+ \param a1 feedback filter coefficient
\param a2 feedback filter coefficient
- \param a3 feedback filter coefficient
*/
-aubio_biquad_t * new_aubio_biquad(lsmp_t b1, lsmp_t b2, lsmp_t b3, lsmp_t a2, lsmp_t a3);
-
-/** delete biquad filter
-
- \param b biquad object to delete
-
-*/
-void del_aubio_biquad(aubio_biquad_t * b);
+aubio_filter_t *
+new_aubio_filter_biquad (lsmp_t b0, lsmp_t b1, lsmp_t b2, lsmp_t a1, lsmp_t a2,
+ uint_t channels);
#ifdef __cplusplus
}
--- a/swig/aubio.i
+++ b/swig/aubio.i
@@ -84,10 +84,8 @@
uint_t aubio_filter_set_c_weighting (aubio_filter_t * b, uint_t samplerate);
/* biquad */
-aubio_biquad_t * new_aubio_biquad(lsmp_t b1, lsmp_t b2, lsmp_t b3, lsmp_t a2, lsmp_t a3);
-void aubio_biquad_do(aubio_biquad_t * b, fvec_t * in);
-void aubio_biquad_do_filtfilt(aubio_biquad_t * b, fvec_t * in, fvec_t * tmp);
-void del_aubio_biquad(aubio_biquad_t * b);
+aubio_filter_t * new_aubio_filter_biquad(lsmp_t b1, lsmp_t b2, lsmp_t b3, lsmp_t a2, lsmp_t a3, uint_t channels);
+uint_t aubio_filter_set_biquad (aubio_filter_t * b, lsmp_t b1, lsmp_t b2, lsmp_t b3, lsmp_t a2, lsmp_t a3);
/* hist */
aubio_hist_t * new_aubio_hist(smpl_t flow, smpl_t fhig, uint_t nelems, uint_t channels);
--- a/tests/src/test-biquad.c
+++ b/tests/src/test-biquad.c
@@ -5,12 +5,12 @@
uint_t win_s = 1024; /* window size */
uint_t channels = 1; /* number of channel */
fvec_t * in = new_fvec (win_s, channels); /* input buffer */
- aubio_biquad_t * o = new_aubio_biquad(0.3,0.2,0.1,0.2,0.3);
+ aubio_filter_t * o = new_aubio_filter_biquad(0.3,0.2,0.1,0.2,0.3, channels);
- aubio_biquad_do_filtfilt(o,in,in);
- aubio_biquad_do(o,in);
+ aubio_filter_do_filtfilt(o,in,in);
+ aubio_filter_do(o,in);
- del_aubio_biquad(o);
+ del_aubio_filter(o);
del_fvec(in);
return 0;
}