Ignore:
File:
1 edited

Legend:

Unmodified
Added
Removed
  • src/filterbank.c

    rbc4ba75 r45134c5  
    1919*/
    2020
     21/* part of this mfcc implementation were inspired from LibXtract
     22   http://libxtract.sourceforge.net/
     23*/
    2124
    2225#include "aubio_priv.h"
    2326#include "sample.h"
    2427#include "filterbank.h"
    25 
    26 #include "stdio.h"
    2728
    2829#define USE_EQUAL_GAIN 1
     
    5253}
    5354
    54 aubio_filterbank_t * new_aubio_filterbank_mfcc(uint_t n_filters, uint_t win_s, uint_t samplerate, smpl_t freq_min, smpl_t freq_max){
    55  
     55aubio_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){
    5656  smpl_t nyquist = samplerate/2.;
    5757  uint_t style = 1;
     
    142142    for(k = next_peak + 1; k < fb->win_s; k++)
    143143      fb->filters[n]->data[0][k]=0.f;
    144 
    145 
    146144  }
    147145
     
    151149  free(fft_peak);
    152150
    153 
    154151  return fb;
    155152
    156153}
    157154
    158 /*
    159 FB initialization based on Slaney's auditory toolbox
    160 TODO:
    161   *solve memory leak problems while
    162   *solve quantization issues when constructing signal:
    163     *bug for win_s=512
    164     *corrections for win_s=1024 -> why even filters with smaller amplitude
    165 
    166 */
    167 
    168 aubio_filterbank_t * new_aubio_filterbank_mfcc2(uint_t n_filters, uint_t win_s, uint_t samplerate, smpl_t freq_min, smpl_t freq_max){
    169  
    170   aubio_filterbank_t * fb = new_aubio_filterbank(n_filters, win_s);
    171  
    172  
    173   //slaney params
    174   smpl_t lowestFrequency = 133.3333;
    175   smpl_t linearSpacing = 66.66666666;
    176   smpl_t logSpacing = 1.0711703;
    177 
    178   uint_t linearFilters = 13;
    179   uint_t logFilters = 27;
    180   uint_t allFilters = linearFilters + logFilters;
    181  
    182   //buffers for computing filter frequencies
    183   fvec_t * freqs=new_fvec(allFilters+2 , 1);
    184  
    185   fvec_t * lower_freqs=new_fvec( allFilters, 1);
    186   fvec_t * upper_freqs=new_fvec( allFilters, 1);
    187   fvec_t * center_freqs=new_fvec( allFilters, 1);
    188 
    189  
    190   fvec_t * triangle_heights=new_fvec( allFilters, 1);
    191   //lookup table of each bin frequency in hz
    192   fvec_t * fft_freqs=new_fvec(win_s, 1);
    193 
    194   uint_t filter_cnt, bin_cnt;
    195  
    196   //first step: filling all the linear filter frequencies
    197   for(filter_cnt=0; filter_cnt<linearFilters; filter_cnt++){
    198     freqs->data[0][filter_cnt]=lowestFrequency+ filter_cnt*linearSpacing;
    199   }
    200   smpl_t lastlinearCF=freqs->data[0][filter_cnt-1];
    201  
    202   //second step: filling all the log filter frequencies
    203   for(filter_cnt=0; filter_cnt<logFilters+2; filter_cnt++){
    204     freqs->data[0][filter_cnt+linearFilters]=lastlinearCF*(pow(logSpacing,filter_cnt+1));
    205   }
    206 
    207   //Option 1. copying interesting values to lower_freqs, center_freqs and upper freqs arrays
    208   //TODO: would be nicer to have a reference to freqs->data, anyway we do not care in this init step
    209    
    210   for(filter_cnt=0; filter_cnt<allFilters; filter_cnt++){
    211     lower_freqs->data[0][filter_cnt]=freqs->data[0][filter_cnt];
    212     center_freqs->data[0][filter_cnt]=freqs->data[0][filter_cnt+1];
    213     upper_freqs->data[0][filter_cnt]=freqs->data[0][filter_cnt+2];
    214   }
    215 
    216 
    217   //computing triangle heights so that each triangle has unit area
    218   for(filter_cnt=0; filter_cnt<allFilters; filter_cnt++){
    219     triangle_heights->data[0][filter_cnt]=2./(upper_freqs->data[0][filter_cnt]-lower_freqs->data[0][filter_cnt]);
    220   }
    221  
    222  
    223   //AUBIO_DBG
    224   AUBIO_DBG("filter tables frequencies\n");
    225   for(filter_cnt=0; filter_cnt<allFilters; filter_cnt++)
    226     AUBIO_DBG("filter n. %d %f %f %f %f\n",filter_cnt, lower_freqs->data[0][filter_cnt], center_freqs->data[0][filter_cnt], upper_freqs->data[0][filter_cnt], triangle_heights->data[0][filter_cnt]);
    227 
    228   //filling the fft_freqs lookup table, which assigns the frequency in hz to each bin
    229   for(bin_cnt=0; bin_cnt<win_s; bin_cnt++){
    230     //TODO: check the formula!
    231     fft_freqs->data[0][bin_cnt]= (smpl_t)samplerate* (smpl_t)bin_cnt/ (smpl_t)win_s;
    232   }
    233 
    234   //building each filter table
    235   for(filter_cnt=0; filter_cnt<allFilters; filter_cnt++){
    236 
    237     //TODO:check special case : lower freq =0
    238     //calculating rise increment in mag/Hz
    239     smpl_t riseInc= triangle_heights->data[0][filter_cnt]/(center_freqs->data[0][filter_cnt]-lower_freqs->data[0][filter_cnt]);
    240    
    241 
    242     AUBIO_DBG("\nfilter %d",filter_cnt);
    243     //zeroing begining of filter
    244     AUBIO_DBG("\nzero begin\n");
    245     for(bin_cnt=0; bin_cnt<win_s-1; bin_cnt++){
    246       //zeroing beigining of array
    247       fb->filters[filter_cnt]->data[0][bin_cnt]=0.f;
    248       AUBIO_DBG(".");
    249       //AUBIO_DBG("%f %f %f\n", fft_freqs->data[0][bin_cnt], fft_freqs->data[0][bin_cnt+1], lower_freqs->data[0][filter_cnt]);
    250       if(fft_freqs->data[0][bin_cnt]<= lower_freqs->data[0][filter_cnt] && fft_freqs->data[0][bin_cnt+1]> lower_freqs->data[0][filter_cnt]){
    251         break;
    252       }
    253     }
    254     bin_cnt++;
    255    
    256     AUBIO_DBG("\npos slope\n");
    257     //positive slope
    258     for(; bin_cnt<win_s-1; bin_cnt++){
    259       AUBIO_DBG(".");
    260       fb->filters[filter_cnt]->data[0][bin_cnt]=(fft_freqs->data[0][bin_cnt]-lower_freqs->data[0][filter_cnt])*riseInc;
    261       //if(fft_freqs->data[0][bin_cnt]<= center_freqs->data[0][filter_cnt] && fft_freqs->data[0][bin_cnt+1]> center_freqs->data[0][filter_cnt])
    262       if(fft_freqs->data[0][bin_cnt+1]> center_freqs->data[0][filter_cnt])
    263         break;
    264     }
    265     //bin_cnt++;
    266    
    267    
    268     //negative slope
    269     AUBIO_DBG("\nneg slope\n");   
    270     for(; bin_cnt<win_s-1; bin_cnt++){
    271       //AUBIO_DBG(".");
    272      
    273       //checking whether last value is less than 0...
    274       smpl_t val=triangle_heights->data[0][filter_cnt]-(fft_freqs->data[0][bin_cnt]-center_freqs->data[0][filter_cnt])*riseInc;
    275       if(val>=0)
    276         fb->filters[filter_cnt]->data[0][bin_cnt]=val;
    277       else fb->filters[filter_cnt]->data[0][bin_cnt]=0.f;
    278      
    279       //if(fft_freqs->data[0][bin_cnt]<= upper_freqs->data[0][bin_cnt] && fft_freqs->data[0][bin_cnt+1]> upper_freqs->data[0][filter_cnt])
    280       //TODO: CHECK whether bugfix correct
    281       if(fft_freqs->data[0][bin_cnt+1]> upper_freqs->data[0][filter_cnt])
    282         break;
    283     }
    284     //bin_cnt++;
    285    
    286     AUBIO_DBG("\nzero end\n");
    287     //zeroing tail
    288     for(; bin_cnt<win_s; bin_cnt++)
    289       //AUBIO_DBG(".");
    290       fb->filters[filter_cnt]->data[0][bin_cnt]=0.f;
    291 
    292   }
    293  
    294  
    295 
    296   del_fvec(freqs);
    297   del_fvec(lower_freqs);
    298   del_fvec(upper_freqs);
    299   del_fvec(center_freqs);
    300 
    301   del_fvec(triangle_heights);
    302   del_fvec(fft_freqs);
    303 
    304   return fb;
    305 
    306 }
    307 
    308 void aubio_dump_filterbank(aubio_filterbank_t * fb){
    309 
    310   FILE * mlog;
    311   mlog=fopen("filterbank.txt","w");
    312  
    313   int k,n;
    314   //dumping filter values
    315   //smpl_t area_tmp=0.f;
    316   for(n = 0; n < fb->n_filters; n++){
    317     for(k = 0; k < fb->win_s; k++){
    318       fprintf(mlog,"%f ",fb->filters[n]->data[0][k]);
    319     }
    320     fprintf(mlog,"\n");
    321   }
    322  
    323   if(mlog) fclose(mlog);
    324 }
    325155
    326156void del_aubio_filterbank(aubio_filterbank_t * fb){
Note: See TracChangeset for help on using the changeset viewer.