src/temporal: derive biquad from filter, use in peakpicker
authorPaul Brossier <piem@piem.org>
Mon, 19 Oct 2009 12:45:13 +0000 (14:45 +0200)
committerPaul Brossier <piem@piem.org>
Mon, 19 Oct 2009 12:45:13 +0000 (14:45 +0200)
src/aubio.h
src/onset/peakpick.c
src/temporal/biquad.c
src/temporal/biquad.h
swig/aubio.i
tests/src/test-biquad.c

index 379d6ef1a3babfe5a2b2254cd3ece7e75d147cbd..6e30008257a62938940c82d3b0a4ff3d069cb67e 100644 (file)
@@ -71,8 +71,8 @@ extern "C"
 #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"
index 17b9a6a6e211f160ee866d90ec7fa5297281fc89..fefe23597d32d81e3599e3136262f269bb8c00c3 100644 (file)
@@ -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 @@ struct _aubio_peakpicker_t {
        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 @@ aubio_peakpicker_do (aubio_peakpicker_t * p, fvec_t * onset, fvec_t * out)
 
   /* 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 @@ aubio_peakpicker_t * new_aubio_peakpicker(uint_t channels) {
        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);
index 7674ccaae0604ad324158649884eddb89b013abc..d715b41211217d8bec00366d09c0fff75124f9f2 100644 (file)
 
 #include "aubio_priv.h"
 #include "fvec.h"
-#include "mathutils.h"
-#include "temporal/biquad.h"
-
-/** \note this file needs to be in double or more less precision would lead to large
- * errors in the output 
- */
-
-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;
+#include "lvec.h"
+#include "temporal/filter.h"
+
+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);
+
+  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;
-  //}
-}
-
-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;
+  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 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;
 }
index 2ee3a810e7abdd8f38530438ce90b52cef86df8d..94f6fd9a333a77d716396ede53ec07b1c5d61ab5 100644 (file)
 
   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
 extern "C" {
 #endif
 
-/** biquad filter object */
-typedef struct _aubio_biquad_t aubio_biquad_t;
-
-/** filter input vector
+/** set coefficients of a biquad filter
 
-  \param b biquad object as returned by new_aubio_biquad
-  \param in input vector to filter
+  \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
 
 */
-void aubio_biquad_do(aubio_biquad_t * b, fvec_t * in);
-/** filter input vector forward and backward
-
-  \param b biquad object as returned by new_aubio_biquad
-  \param in input vector to filter
-  \param tmp memory space to use for computation
+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);
 
-*/
-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
 }
index a314ae5cd62d873b9549263c4ad79a8888a102ab..3f87c8c2f27b6623bae36e689ce3223561fe41c4 100644 (file)
@@ -84,10 +84,8 @@ aubio_filter_t * new_aubio_filter_c_weighting (uint_t channels, uint_t samplerat
 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);
index 319e2307bb7527e3b3ccf672bad806e2d3fd8863..78d73188f9a16067cc59ca024d43a3b71d54c627 100644 (file)
@@ -5,12 +5,12 @@ int main(){
         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;
 }