Ignore:
File:
1 edited

Legend:

Unmodified
Added
Removed
  • src/filterbank.c

    r45134c5 r7c6c806d  
    11/*
    2    Copyright (C) 2007 Amaury Hazan <ahazan@iua.upf.edu>
    3                   and Paul Brossier <piem@piem.org>
     2   Copyright (C) 2007 Amaury Hazan
     3   Ported to aubio from LibXtract
     4   http://libxtract.sourceforge.net/
     5   
    46
    57   This program is free software; you can redistribute it and/or modify
     
    1921*/
    2022
    21 /* part of this mfcc implementation were inspired from LibXtract
    22    http://libxtract.sourceforge.net/
    23 */
    24 
    2523#include "aubio_priv.h"
    26 #include "sample.h"
    2724#include "filterbank.h"
    2825
    29 #define USE_EQUAL_GAIN 1
    30 #define VERY_SMALL_NUMBER 2e-42
    3126
    32 /** \brief A structure to store a set of n_filters filters of lenghts win_s */
    33 struct aubio_filterbank_t_ {
    34     uint_t win_s;
    35     uint_t n_filters;
    36     fvec_t **filters;
     27// Struct Declaration
     28
     29/** \brief A structure to store a set of n_filters Mel filters */
     30typedef struct aubio_mel_filter_ {
     31    int n_filters;
     32    smpl_t **filters;
    3733};
    3834
    39 aubio_filterbank_t * new_aubio_filterbank(uint_t n_filters, uint_t win_s){
    40   /** allocating space for filterbank object */
    41   aubio_filterbank_t * fb = AUBIO_NEW(aubio_filterbank_t);
    42   uint_t filter_cnt;
    43   fb->win_s=win_s;
    44   fb->n_filters=n_filters;
     35// Initialization
    4536
    46   /** allocating filter tables */
    47   fb->filters=AUBIO_ARRAY(fvec_t*,n_filters);
    48   for (filter_cnt=0; filter_cnt<n_filters; filter_cnt++)
    49     /* considering one-channel filters */
    50     fb->filters[filter_cnt]=new_fvec(win_s, 1);
     37int aubio_mfcc_init(int N, smpl_t nyquist, int style, smpl_t freq_min, smpl_t freq_max, int freq_bands, smpl_t **fft_tables){
    5138
    52   return fb;
    53 }
     39    int n, i, k, *fft_peak, M, next_peak;
     40    smpl_t norm, mel_freq_max, mel_freq_min, norm_fact, height, inc, val,
     41        freq_bw_mel, *mel_peak, *height_norm, *lin_peak;
    5442
    55 aubio_filterbank_t * new_aubio_filterbank_mfcc(uint_t n_filters, uint_t win_s, smpl_t samplerate, smpl_t freq_min, smpl_t freq_max){
    56   smpl_t nyquist = samplerate/2.;
    57   uint_t style = 1;
    58   aubio_filterbank_t * fb = new_aubio_filterbank(n_filters, win_s);
     43    mel_peak = height_norm = lin_peak = NULL;
     44    fft_peak = NULL;
     45    norm = 1;
    5946
    60   uint_t n, i, k, *fft_peak, M, next_peak;
    61   smpl_t norm, mel_freq_max, mel_freq_min, norm_fact, height, inc, val,
    62          freq_bw_mel, *mel_peak, *height_norm, *lin_peak;
     47    mel_freq_max = 1127 * log(1 + freq_max / 700);
     48    mel_freq_min = 1127 * log(1 + freq_min / 700);
     49    freq_bw_mel = (mel_freq_max - mel_freq_min) / freq_bands;
    6350
    64   mel_peak = height_norm = lin_peak = NULL;
    65   fft_peak = NULL;
    66   norm = 1;
     51    mel_peak = (smpl_t *)malloc((freq_bands + 2) * sizeof(smpl_t));
     52    /* +2 for zeros at start and end */
     53    lin_peak = (smpl_t *)malloc((freq_bands + 2) * sizeof(smpl_t));
     54    fft_peak = (int *)malloc((freq_bands + 2) * sizeof(int));
     55    height_norm = (smpl_t *)malloc(freq_bands * sizeof(smpl_t));
    6756
    68   mel_freq_max = 1127 * log(1 + freq_max / 700);
    69   mel_freq_min = 1127 * log(1 + freq_min / 700);
    70   freq_bw_mel = (mel_freq_max - mel_freq_min) / fb->n_filters;
     57    if(mel_peak == NULL || height_norm == NULL ||
     58                    lin_peak == NULL || fft_peak == NULL)
     59                    return XTRACT_MALLOC_FAILED;
     60   
     61    M = N >> 1;
    7162
    72   mel_peak = (smpl_t *)malloc((fb->n_filters + 2) * sizeof(smpl_t));
    73   /* +2 for zeros at start and end */
    74   lin_peak = (smpl_t *)malloc((fb->n_filters + 2) * sizeof(smpl_t));
    75   fft_peak = (uint_t *)malloc((fb->n_filters + 2) * sizeof(uint_t));
    76   height_norm = (smpl_t *)malloc(fb->n_filters * sizeof(smpl_t));
     63    mel_peak[0] = mel_freq_min;
     64    lin_peak[0] = 700 * (exp(mel_peak[0] / 1127) - 1);
     65    fft_peak[0] = lin_peak[0] / nyquist * M;
    7766
    78   if(mel_peak == NULL || height_norm == NULL ||
    79       lin_peak == NULL || fft_peak == NULL)
    80     return NULL;
    8167
    82   M = fb->win_s >> 1;
    83 
    84   mel_peak[0] = mel_freq_min;
    85   lin_peak[0] = 700 * (exp(mel_peak[0] / 1127) - 1);
    86   fft_peak[0] = lin_peak[0] / nyquist * M;
    87 
    88   for (n = 1; n <= fb->n_filters; n++){ 
     68    for (n = 1; n <= freq_bands; n++){ 
    8969    /*roll out peak locations - mel, linear and linear on fft window scale */
    90     mel_peak[n] = mel_peak[n - 1] + freq_bw_mel;
    91     lin_peak[n] = 700 * (exp(mel_peak[n] / 1127) -1);
    92     fft_peak[n] = lin_peak[n] / nyquist * M;
    93   }
    94 
    95   for (n = 0; n < fb->n_filters; n++){
    96     /*roll out normalised gain of each peak*/
    97     if (style == USE_EQUAL_GAIN){
    98       height = 1;
    99       norm_fact = norm;
    100     }
    101     else{
    102       height = 2 / (lin_peak[n + 2] - lin_peak[n]);
    103       norm_fact = norm / (2 / (lin_peak[2] - lin_peak[0]));
    104     }
    105     height_norm[n] = height * norm_fact;
    106   }
    107 
    108   i = 0;
    109 
    110   for(n = 0; n < fb->n_filters; n++){
    111 
    112     /*calculate the rise increment*/
    113     if(n > 0)
    114       inc = height_norm[n] / (fft_peak[n] - fft_peak[n - 1]);
    115     else
    116       inc = height_norm[n] / fft_peak[n];
    117     val = 0; 
    118 
    119     /*zero the start of the array*/
    120     for(k = 0; k < i; k++)
    121       fb->filters[n]->data[0][k]=0.f;
    122 
    123     /*fill in the rise */
    124     for(; i <= fft_peak[n]; i++){
    125       fb->filters[n]->data[0][k]=val;
    126       val += inc;
     70        mel_peak[n] = mel_peak[n - 1] + freq_bw_mel;
     71        lin_peak[n] = 700 * (exp(mel_peak[n] / 1127) -1);
     72        fft_peak[n] = lin_peak[n] / nyquist * M;
    12773    }
    12874
    129     /*calculate the fall increment */
    130     inc = height_norm[n] / (fft_peak[n + 1] - fft_peak[n]);
    131 
    132     val = 0;
    133     next_peak = fft_peak[n + 1];
    134 
    135     /*reverse fill the 'fall' */
    136     for(i = next_peak; i > fft_peak[n]; i--){
    137       fb->filters[n]->data[0][k]=val;
    138       val += inc;
     75    for (n = 0; n < freq_bands; n++){
     76        /*roll out normalised gain of each peak*/
     77        if (style == XTRACT_EQUAL_GAIN){
     78            height = 1;
     79            norm_fact = norm;
     80        }
     81        else{
     82            height = 2 / (lin_peak[n + 2] - lin_peak[n]);
     83            norm_fact = norm / (2 / (lin_peak[2] - lin_peak[0]));
     84        }
     85        height_norm[n] = height * norm_fact;
    13986    }
    14087
    141     /*zero the rest of the array*/
    142     for(k = next_peak + 1; k < fb->win_s; k++)
    143       fb->filters[n]->data[0][k]=0.f;
    144   }
     88    i = 0;
     89   
     90    for(n = 0; n < freq_bands; n++){
     91 
     92  /*calculate the rise increment*/
     93        if(n > 0)
     94            inc = height_norm[n] / (fft_peak[n] - fft_peak[n - 1]);
     95        else
     96            inc = height_norm[n] / fft_peak[n];
     97        val = 0; 
     98 
     99  /*zero the start of the array*/
     100  for(k = 0; k < i; k++)
     101     fft_tables[n][k] = 0.f;
     102 
     103  /*fill in the rise */
     104        for(; i <= fft_peak[n]; i++){
     105            fft_tables[n][i] = val;
     106            val += inc;
     107        }
     108 
     109        /*calculate the fall increment */
     110        inc = height_norm[n] / (fft_peak[n + 1] - fft_peak[n]);
     111 
     112        val = 0;
     113  next_peak = fft_peak[n + 1];
     114 
     115  /*reverse fill the 'fall' */
     116        for(i = next_peak; i > fft_peak[n]; i--){
     117            fft_tables[n][i] = val;
     118            val += inc;
     119        }
    145120
    146   free(mel_peak);
    147   free(lin_peak);
    148   free(height_norm);
    149   free(fft_peak);
     121  /*zero the rest of the array*/
     122  for(k = next_peak + 1; k < N; k++)
     123      fft_tables[n][k] = 0.f;
     124    }
    150125
    151   return fb;
     126    free(mel_peak);
     127    free(lin_peak);
     128    free(height_norm);
     129    free(fft_peak);
     130
     131    return XTRACT_SUCCESS;
    152132
    153133}
    154 
    155 
    156 void del_aubio_filterbank(aubio_filterbank_t * fb){
    157   uint_t filter_cnt;
    158   /** deleting filter tables first */
    159   for (filter_cnt=0; filter_cnt<fb->n_filters; filter_cnt++)
    160     del_fvec(fb->filters[filter_cnt]);
    161   AUBIO_FREE(fb->filters);
    162   AUBIO_FREE(fb);
    163 }
    164 
    165 void aubio_filterbank_do(aubio_filterbank_t * f, cvec_t * in, fvec_t *out) {
    166   uint_t n, filter_cnt;
    167   for(filter_cnt = 0; (filter_cnt < f->n_filters)
    168     && (filter_cnt < out->length); filter_cnt++){
    169       out->data[0][filter_cnt] = 0.f;
    170       for(n = 0; n < in->length; n++){
    171           out->data[0][filter_cnt] += in->norm[0][n]
    172             * f->filters[filter_cnt]->data[0][n];
    173       }
    174       out->data[0][filter_cnt] =
    175         LOG(out->data[0][filter_cnt] < VERY_SMALL_NUMBER ?
    176             VERY_SMALL_NUMBER : out->data[0][filter_cnt]);
    177   }
    178 
    179   return;
    180 }
Note: See TracChangeset for help on using the changeset viewer.