Changeset b01bd4a


Ignore:
Timestamp:
Oct 19, 2009, 2:45:13 PM (15 years ago)
Author:
Paul Brossier <piem@piem.org>
Branches:
feature/autosink, feature/cnn, feature/cnn_org, feature/constantq, feature/crepe, feature/crepe_org, feature/pitchshift, feature/pydocstrings, feature/timestretch, fix/ffmpeg5, master, pitchshift, sampler, timestretch, yinfft+
Children:
b849106
Parents:
59c046d
Message:

src/temporal: derive biquad from filter, use in peakpicker

Files:
6 edited

Legend:

Unmodified
Added
Removed
  • src/aubio.h

    r59c046d rb01bd4a  
    7272#include "temporal/resampler.h"
    7373#endif /* HAVE_SAMPLERATE */
     74#include "temporal/filter.h"
    7475#include "temporal/biquad.h"
    75 #include "temporal/filter.h"
    7676#include "temporal/a_weighting.h"
    7777#include "temporal/c_weighting.h"
  • src/onset/peakpick.c

    r59c046d rb01bd4a  
    2121#include "fvec.h"
    2222#include "mathutils.h"
     23#include "lvec.h"
     24#include "temporal/filter.h"
    2325#include "temporal/biquad.h"
    2426#include "onset/peakpick.h"
     
    4345
    4446        /** biquad lowpass filter */
    45         aubio_biquad_t * biquad;
     47        aubio_filter_t * biquad;
    4648        /** original onsets */
    4749        fvec_t * onset_keep;
     
    9597  /* filter onset_proc */
    9698  /** \bug filtfilt calculated post+pre times, should be only once !? */
    97   //aubio_biquad_do_filtfilt(p->biquad,onset_proc,scratch);
     99  aubio_filter_do_filtfilt (p->biquad, onset_proc, scratch);
    98100
    99101  for (i = 0; i < p->channels; i++) {
     
    160162        t->onset_peek = new_fvec(3, channels);
    161163
    162         /* cutoff: low-pass filter cutoff [0.34, 1] */
    163         /* t->cutoff=0.34; */
    164         t->biquad = new_aubio_biquad(0.1600,0.3200,0.1600,-0.5949,0.2348);
     164        /* cutoff: low-pass filter with cutoff reduced frequency at 0.34
     165   generated with octave butter function: [b,a] = butter(2, 0.34);
     166   */
     167  t->biquad = new_aubio_filter_biquad (0.15998789, 0.31997577, 0.15998789,
     168    -0.59488894, 0.23484048, channels);
     169
    165170        return t;
    166171}
    167172
    168173void del_aubio_peakpicker(aubio_peakpicker_t * p) {
    169         del_aubio_biquad(p->biquad);
     174        del_aubio_filter(p->biquad);
    170175        del_fvec(p->onset_keep);
    171176        del_fvec(p->onset_proc);
  • src/temporal/biquad.c

    r59c046d rb01bd4a  
    2020#include "aubio_priv.h"
    2121#include "fvec.h"
    22 #include "mathutils.h"
    23 #include "temporal/biquad.h"
     22#include "lvec.h"
     23#include "temporal/filter.h"
    2424
    25 /** \note this file needs to be in double or more less precision would lead to large
    26  * errors in the output
    27  */
     25uint_t
     26aubio_filter_set_biquad (aubio_filter_t * f, lsmp_t b0, lsmp_t b1, lsmp_t b2,
     27    lsmp_t a1, lsmp_t a2) {
     28  uint_t order = aubio_filter_get_order (f);
     29  lvec_t *bs = aubio_filter_get_feedforward (f);
     30  lvec_t *as = aubio_filter_get_feedback (f);
    2831
    29 struct _aubio_biquad_t {
    30   lsmp_t a2;
    31   lsmp_t a3;
    32   lsmp_t b1;
    33   lsmp_t b2;
    34   lsmp_t b3;
    35   lsmp_t o1;
    36   lsmp_t o2;
    37   lsmp_t i1;
    38   lsmp_t i2;
    39 };
    40 
    41 void aubio_biquad_do(aubio_biquad_t * b, fvec_t * in) {
    42   uint_t i,j;
    43   lsmp_t i1 = b->i1;
    44   lsmp_t i2 = b->i2;
    45   lsmp_t o1 = b->o1;
    46   lsmp_t o2 = b->o2;
    47   lsmp_t a2 = b->a2;
    48   lsmp_t a3 = b->a3;
    49   lsmp_t b1 = b->b1;
    50   lsmp_t b2 = b->b2;
    51   lsmp_t b3 = b->b3;
    52 
    53   i=0; // works in mono only !!!
    54   //for (i=0;i<in->channels;i++) {
    55   for (j = 0; j < in->length; j++) {
    56     lsmp_t i0 = in->data[i][j];
    57     lsmp_t o0 = b1 * i0 + b2 * i1 + b3 * i2
    58       - a2 * o1 - a3 * o2;// + 1e-37;
    59     in->data[i][j] = o0;
    60     i2 = i1;
    61     i1 = i0;
    62     o2 = o1;
    63     o1 = o0;
     32  if (order != 3) {
     33    AUBIO_ERROR ("order of biquad filter must be 3, not %d\n", order);
     34    return AUBIO_FAIL;
    6435  }
    65   b->i2 = i2;
    66   b->i1 = i1;
    67   b->o2 = o2;
    68   b->o1 = o1;
    69   //}
     36  bs->data[0][0] = b0;
     37  bs->data[0][1] = b1;
     38  bs->data[0][2] = b2;
     39  as->data[0][0] = 1.;
     40  as->data[0][1] = a1;
     41  as->data[0][1] = a2;
     42  return AUBIO_OK;
    7043}
    7144
    72 void aubio_biquad_do_filtfilt(aubio_biquad_t * b, fvec_t * in, fvec_t * tmp) {
    73   uint_t j,i=0;
    74   uint_t length = in->length;
    75   lsmp_t mir;
    76   /* mirroring */
    77   mir = 2*in->data[i][0];
    78   b->i1 = mir - in->data[i][2];
    79   b->i2 = mir - in->data[i][1];
    80   /* apply filtering */
    81   aubio_biquad_do(b,in);
    82   /* invert  */
    83   for (j = 0; j < length; j++)
    84     tmp->data[i][length-j-1] = in->data[i][j];
    85   /* mirror again */
    86   mir = 2*tmp->data[i][0];
    87   b->i1 = mir - tmp->data[i][2];
    88   b->i2 = mir - tmp->data[i][1];
    89   /* apply filtering */
    90   aubio_biquad_do(b,tmp);
    91   /* invert back */
    92   for (j = 0; j < length; j++)
    93     in->data[i][j] = tmp->data[i][length-j-1];
     45aubio_filter_t *
     46new_aubio_filter_biquad (lsmp_t b0, lsmp_t b1, lsmp_t b2, lsmp_t a1, lsmp_t a2,
     47    uint_t channels)
     48{
     49  aubio_filter_t *f = new_aubio_filter (3, channels);
     50  aubio_filter_set_biquad (f, b0, b1, b2, a1, a2);
     51  return f;
    9452}
    95 
    96 aubio_biquad_t * new_aubio_biquad(
    97     lsmp_t b1, lsmp_t b2, lsmp_t b3,
    98     lsmp_t a2, lsmp_t a3) {
    99   aubio_biquad_t * b = AUBIO_NEW(aubio_biquad_t);
    100   b->a2 = a2;
    101   b->a3 = a3;
    102   b->b1 = b1;
    103   b->b2 = b2;
    104   b->b3 = b3;
    105   b->i1 = 0.;
    106   b->i2 = 0.;
    107   b->o1 = 0.;
    108   b->o2 = 0.;
    109   return b;
    110 }
    111 
    112 void del_aubio_biquad(aubio_biquad_t * b) {
    113   AUBIO_FREE(b);
    114 }
  • src/temporal/biquad.h

    r59c046d rb01bd4a  
    2727  This file implements a normalised biquad filter (second order IIR):
    2828 
    29   \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$
     29  \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$
    3030
    3131  The filtfilt version runs the filter twice, forward and backward, to
    3232  compensate the phase shifting of the forward operation.
     33
     34  See also <a href="http://en.wikipedia.org/wiki/Digital_biquad_filter">Digital
     35  biquad filter</a> on wikipedia.
    3336
    3437*/
     
    3841#endif
    3942
    40 /** biquad filter object */
    41 typedef struct _aubio_biquad_t aubio_biquad_t;
     43/** set coefficients of a biquad filter
    4244
    43 /** filter input vector
    44 
    45   \param b biquad object as returned by new_aubio_biquad
    46   \param in input vector to filter
     45  \param b0 forward filter coefficient
     46  \param b1 forward filter coefficient
     47  \param b2 forward filter coefficient
     48  \param a1 feedback filter coefficient
     49  \param a2 feedback filter coefficient
    4750
    4851*/
    49 void aubio_biquad_do(aubio_biquad_t * b, fvec_t * in);
    50 /** filter input vector forward and backward
     52uint_t
     53aubio_filter_set_biquad (aubio_filter_t * f, lsmp_t b0, lsmp_t b1, lsmp_t b2,
     54    lsmp_t a1, lsmp_t a2);
    5155
    52   \param b biquad object as returned by new_aubio_biquad
    53   \param in input vector to filter
    54   \param tmp memory space to use for computation
     56/** create new biquad filter
     57
     58  \param b0 forward filter coefficient
     59  \param b1 forward filter coefficient
     60  \param b2 forward filter coefficient
     61  \param a1 feedback filter coefficient
     62  \param a2 feedback filter coefficient
    5563
    5664*/
    57 void aubio_biquad_do_filtfilt(aubio_biquad_t * b, fvec_t * in, fvec_t * tmp);
    58 /** create new biquad filter
    59 
    60   \param b1 forward filter coefficient
    61   \param b2 forward filter coefficient
    62   \param b3 forward filter coefficient
    63   \param a2 feedback filter coefficient
    64   \param a3 feedback filter coefficient
    65 
    66 */
    67 aubio_biquad_t * new_aubio_biquad(lsmp_t b1, lsmp_t b2, lsmp_t b3, lsmp_t a2, lsmp_t a3);
    68 
    69 /** delete biquad filter
    70  
    71   \param b biquad object to delete
    72 
    73 */
    74 void del_aubio_biquad(aubio_biquad_t * b);
     65aubio_filter_t *
     66new_aubio_filter_biquad (lsmp_t b0, lsmp_t b1, lsmp_t b2, lsmp_t a1, lsmp_t a2,
     67    uint_t channels);
    7568
    7669#ifdef __cplusplus
  • swig/aubio.i

    r59c046d rb01bd4a  
    8585
    8686/* biquad */
    87 aubio_biquad_t * new_aubio_biquad(lsmp_t b1, lsmp_t b2, lsmp_t b3, lsmp_t a2, lsmp_t a3);
    88 void aubio_biquad_do(aubio_biquad_t * b, fvec_t * in);
    89 void aubio_biquad_do_filtfilt(aubio_biquad_t * b, fvec_t * in, fvec_t * tmp);
    90 void del_aubio_biquad(aubio_biquad_t * b);
     87aubio_filter_t * new_aubio_filter_biquad(lsmp_t b1, lsmp_t b2, lsmp_t b3, lsmp_t a2, lsmp_t a3, uint_t channels);
     88uint_t aubio_filter_set_biquad (aubio_filter_t * b, lsmp_t b1, lsmp_t b2, lsmp_t b3, lsmp_t a2, lsmp_t a3);
    9189
    9290/* hist */
  • tests/src/test-biquad.c

    r59c046d rb01bd4a  
    66        uint_t channels   = 1;                          /* number of channel */
    77        fvec_t * in       = new_fvec (win_s, channels); /* input buffer */
    8         aubio_biquad_t * o = new_aubio_biquad(0.3,0.2,0.1,0.2,0.3);
     8        aubio_filter_t * o = new_aubio_filter_biquad(0.3,0.2,0.1,0.2,0.3, channels);
    99
    10         aubio_biquad_do_filtfilt(o,in,in);
    11         aubio_biquad_do(o,in);
     10        aubio_filter_do_filtfilt(o,in,in);
     11        aubio_filter_do(o,in);
    1212
    13         del_aubio_biquad(o);
     13        del_aubio_filter(o);
    1414        del_fvec(in);
    1515        return 0;
Note: See TracChangeset for help on using the changeset viewer.