ref: 21bd43cd56c206a650d48e79d966bd4c40dff3dd
parent: cb0415db6f21cc0df0293b63e15f6a17b71457ef
parent: 870855670c7d903af86e10ca5e5051a30e000f90
author: Paul Brossier <piem@piem.org>
date: Sat Sep 8 09:52:48 EDT 2007
merge from Amaury's changes
--- /dev/null
+++ b/examples/aubiomfcc.c
@@ -1,0 +1,113 @@
+/*
+ Copyright (C) 2007 Amaury Hazan
+
+ This program is free software; you can redistribute it and/or modify
+ it under the terms of the GNU General Public License as published by
+ the Free Software Foundation; either version 2 of the License, or
+ (at your option) any later version.
+
+ This program is distributed in the hope that it will be useful,
+ but WITHOUT ANY WARRANTY; without even the implied warranty of
+ MERCHANTABILITY or FITNESS FOR A PARTICULAR PURPOSE. See the
+ GNU General Public License for more details.
+
+ You should have received a copy of the GNU General Public License
+ along with this program; if not, write to the Free Software
+ Foundation, Inc., 675 Mass Ave, Cambridge, MA 02139, USA.
+*/
+
+#include "utils.h"
+
+unsigned int pos = 0; /*frames%dspblocksize*/
+uint_t usepitch = 0;
+
+int aubio_process(float **input, float **output, int nframes);
+int aubio_process(float **input, float **output, int nframes) {+ unsigned int i; /*channels*/
+ unsigned int j; /*frames*/
+
+ for (j=0;j<(unsigned)nframes;j++) {+ if(usejack) {+ for (i=0;i<channels;i++) {+ /* write input to datanew */
+ fvec_write_sample(ibuf, input[i][j], i, pos);
+ /* put synthnew in output */
+ output[i][j] = fvec_read_sample(obuf, i, pos);
+ }
+ }
+ /*time for fft*/
+ if (pos == overlap_size-1) { + /* block loop */
+
+ //compute mag spectrum
+ aubio_pvoc_do (pv,ibuf, fftgrain);
+
+ uint_t n_coefs= n_filters/2 +1;
+ uint_t coef_cnt;
+
+
+ for (coef_cnt=0; coef_cnt<n_coefs ; coef_cnt++)
+ mfcc_outbuf[coef_cnt]=0.f;
+
+ //compute mfccs
+ aubio_mffc_do(fftgrain->norm, nframes, mf, mfcc_outbuf, fft_dct, fftgrain_dct);
+
+ for (coef_cnt=0; coef_cnt<n_coefs ; coef_cnt++)
+ outmsg("%f ",mfcc_outbuf[coef_cnt]);+ outmsg("\n");+
+
+
+ /* end of block loop */
+ pos = -1; /* so it will be zero next j loop */
+ }
+ pos++;
+ }
+ return 1;
+}
+
+void process_print (void);
+void process_print (void) {+ /* output times in seconds
+ write extracted mfccs
+ */
+
+ if (output_filename == NULL) {+ if(frames >= 4) {+ outmsg("%f\n",(frames-4)*overlap_size/(float)samplerate);+ } else if (frames < 4) {+ outmsg("%f\n",0.);+ }
+ }
+}
+
+int main(int argc, char **argv) {+ examples_common_init(argc,argv);
+
+ //allocate and initialize mel filter bank
+
+
+ //allocating global mf (in utils.c)
+ uint_t banksize = (uint) ( sizeof(aubio_mel_filter));
+ mf = (aubio_mel_filter *)getbytes(banksize);
+
+ mf->n_filters = 20;
+ mf->filters = (smpl_t **)getbytes(mf->n_filters * sizeof(smpl_t *));
+ for(n = 0; n < mf->n_filters; n++)
+ mf->filters[n] = (smpl_t *)getbytes((buffer_size/2+1) * sizeof(smpl_t));
+
+ //populating the filter
+ aubio_mfcc_init(buffer_size, nyquist, XTRACT_EQUAL_GAIN, lowfreq, highfreq, mf->n_filters, mf->filters);
+
+ //process
+ examples_common_process(aubio_process,process_print);
+ examples_common_del();
+ debug("End of program.\n");+ fflush(stderr);
+
+ //destroying filterbank
+ free(mf);
+
+ return 0;
+}
+
--- a/examples/utils.c
+++ b/examples/utils.c
@@ -60,7 +60,21 @@
int isonset = 0;
aubio_pickpeak_t * parms;
+/* mfcc objects */
+//parameters
+uint_t n_filters=20;
+uint_t nyquist= samplerate / 2.;
+smpl_t lowfreq=80.f;
+smpl_t highfreq=18000.f;
+// filterbank object
+aubio_mel_filter * mf;
+// DCT mfft and result storage
+aubio_mfft * fft_dct;
+cvec_t * fftgrain_dct;
+smpl_t mfcc_outbuf[11];
+
+
/* pitch objects */
smpl_t pitch = 0.;
aubio_pitchdetection_t * pitchdet;
@@ -300,6 +314,9 @@
obuf = new_fvec(overlap_size, channels);
fftgrain = new_cvec(buffer_size, channels);
+ //init for mfcc process
+ fftgrain_dct= new_cvec(n_filters, channels);
+
if (usepitch) {pitchdet = new_aubio_pitchdetection(buffer_size*4,
overlap_size, channels, samplerate, type_pitch, mode_pitch);
@@ -312,6 +329,11 @@
}
/* phase vocoder */
pv = new_aubio_pvoc(buffer_size, overlap_size, channels);
+
+ // dct phase vocoder
+ //TODO: check size
+ fft_dct = new_aubio_mfft(n_filters, channels);
+
/* onsets */
parms = new_aubio_peakpicker(threshold);
o = new_aubio_onsetdetection(type_onset,buffer_size,channels);
@@ -345,6 +367,11 @@
del_cvec(fftgrain);
del_fvec(onset);
del_fvec(woodblock);
+
+ //mffc related
+ del_aubio_mfft(fft_dct);
+ del_cvec(fftgrain_dct);
+
aubio_cleanup();
}
--- a/examples/utils.h
+++ b/examples/utils.h
@@ -97,6 +97,18 @@
extern int isonset;
extern aubio_pickpeak_t * parms;
+/* mfcc objects */
+// params
+extern uint_t n_filters;
+extern uint_t nyquist;
+extern smpl_t lowfreq;
+extern smpl_t highfreq;
+// filterbank object
+extern aubio_mel_filter * mf;
+// DCT pvoc and result storage
+extern aubio_mfft_t * fft_dct;
+extern cvec_t * fftgrain_dct;
+extern smpl_t mfcc_outbuf[20];
/* pitch objects */
extern smpl_t pitch;
--- a/src/Makefile.am
+++ b/src/Makefile.am
@@ -21,7 +21,9 @@
beattracking.h \
onset.h \
tempo.h \
- filter.h
+ filter.h \
+ mfcc.h
+
nodist_pkginclude_HEADERS = config.h
lib_LTLIBRARIES = libaubio.la
@@ -68,7 +70,9 @@
tempo.c \
tempo.h \
filter.c \
- filter.h
+ filter.h \
+ mfcc.h \
+ mfcc.c
AM_CFLAGS = @AUBIO_CFLAGS@ @FFTWLIB_CFLAGS@ @SAMPLERATE_CFLAGS@
libaubio_la_LIBADD = @FFTWLIB_LIBS@ @SAMPLERATE_LIBS@ @LTLIBOBJS@
--- a/src/aubio.h
+++ b/src/aubio.h
@@ -79,6 +79,7 @@
#include "beattracking.h"
#include "onset.h"
#include "tempo.h"
+#include "mfcc.h"
#ifdef __cplusplus
} /* extern "C" */
--- /dev/null
+++ b/src/bark.c
@@ -1,0 +1,54 @@
+/*
+ Copyright (C) 2006 Amaury Hazan
+ Ported to aubio from LibXtract
+ http://libxtract.sourceforge.net/
+
+
+ This program is free software; you can redistribute it and/or modify
+ it under the terms of the GNU General Public License as published by
+ the Free Software Foundation; either version 2 of the License, or
+ (at your option) any later version.
+
+ This program is distributed in the hope that it will be useful,
+ but WITHOUT ANY WARRANTY; without even the implied warranty of
+ MERCHANTABILITY or FITNESS FOR A PARTICULAR PURPOSE. See the
+ GNU General Public License for more details.
+
+ You should have received a copy of the GNU General Public License
+ along with this program; if not, write to the Free Software
+ Foundation, Inc., 675 Mass Ave, Cambridge, MA 02139, USA.
+
+*/
+
+#define XTRACT_BARK_BANDS 26
+#include "bark.c"
+
+int xtract_init_bark(int N, float sr, int *band_limits){+
+ float edges[] = {0, 100, 200, 300, 400, 510, 630, 770, 920, 1080, 1270, 1480, 1720, 2000, 2320, 2700, 3150, 3700, 4400, 5300, 6400, 7700, 9500, 12000, 15500, 20500, 27000}; /* Takes us up to sr = 54kHz (CCRMA: JOS)*/+
+ int bands = XTRACT_BARK_BANDS;
+
+ while(bands--)
+ band_limits[bands] = edges[bands] / sr * N;
+ /*FIX shohuld use rounding, but couldn't get it to work */
+
+ return XTRACT_SUCCESS;
+}
+
+
+
+
+int xtract_bark_coefficients(const float *data, const int N, const void *argv, float *result){+
+ int *limits, band, n;
+
+ limits = (int *)argv;
+
+ for(band = 0; band < XTRACT_BARK_BANDS - 1; band++){+ for(n = limits[band]; n < limits[band + 1]; n++)
+ result[band] += data[n];
+ }
+
+ return XTRACT_SUCCESS;
+}
\ No newline at end of file
--- /dev/null
+++ b/src/bark.h
@@ -1,0 +1,53 @@
+/*
+ Copyright (C) 2006 Amaury Hazan
+ Ported to aubio from LibXtract
+ http://libxtract.sourceforge.net/
+
+
+ This program is free software; you can redistribute it and/or modify
+ it under the terms of the GNU General Public License as published by
+ the Free Software Foundation; either version 2 of the License, or
+ (at your option) any later version.
+
+ This program is distributed in the hope that it will be useful,
+ but WITHOUT ANY WARRANTY; without even the implied warranty of
+ MERCHANTABILITY or FITNESS FOR A PARTICULAR PURPOSE. See the
+ GNU General Public License for more details.
+
+ You should have received a copy of the GNU General Public License
+ along with this program; if not, write to the Free Software
+ Foundation, Inc., 675 Mass Ave, Cambridge, MA 02139, USA.
+
+*/
+
+#ifndef BARK_H
+#define BARK_H
+
+
+// Initalization
+
+/** \brief A function to initialise bark filter bounds
+ *
+ * A pointer to an array of BARK_BANDS ints most be passed in, and is populated with BARK_BANDS fft bin numbers representing the limits of each band
+ *
+ * \param N: the audio block size
+ * \param sr: The sample audio sample rate
+ * \param *band_limits: a pointer to an array of BARK_BANDS ints
+ */
+int xtract_init_bark(int N, float sr, int *band_limits);
+
+// Computation
+
+/** \brief Extract Bark band coefficients based on a method
+ * \param *data: a pointer to the first element in an array of floats representing the magnitude coefficients from the magnitude spectrum of an audio vector, (e.g. the first half of the array pointed to by *result from xtract_spectrum().
+ * \param N: the number of array elements to be considered
+ * \param *argv: a pointer to an array of ints representing the limits of each bark band. This can be obtained by calling xtract_init_bark.
+ * \param *result: a pointer to an array containing resultant bark coefficients
+ *
+ * The limits array pointed to by *argv must be obtained by first calling xtract_init_bark
+ *
+ */
+int xtract_bark_coefficients(const float *data, const int N, const void *argv, float *result);
+
+
+#endif
\ No newline at end of file
--- /dev/null
+++ b/src/filterbank.c
@@ -1,0 +1,266 @@
+/*
+ Copyright (C) 2007 Amaury Hazan
+ Ported to aubio from LibXtract
+ http://libxtract.sourceforge.net/
+
+
+ This program is free software; you can redistribute it and/or modify
+ it under the terms of the GNU General Public License as published by
+ the Free Software Foundation; either version 2 of the License, or
+ (at your option) any later version.
+
+ This program is distributed in the hope that it will be useful,
+ but WITHOUT ANY WARRANTY; without even the implied warranty of
+ MERCHANTABILITY or FITNESS FOR A PARTICULAR PURPOSE. See the
+ GNU General Public License for more details.
+
+ You should have received a copy of the GNU General Public License
+ along with this program; if not, write to the Free Software
+ Foundation, Inc., 675 Mass Ave, Cambridge, MA 02139, USA.
+
+*/
+
+#include "aubio_priv.h"
+#include "filterbank.h"
+
+
+// Struct Declaration
+
+/** \brief A structure to store a set of n_filters filters of lenghts win_s */
+struct aubio_filterbank_t_ {+ uint_t win_s;
+ uint_t n_filters;
+ fvec_t *filters;
+};
+
+aubio_filterbank_t * new_aubio_filterbank(uint_t n_filters, uint_t win_s){+
+ int filter_cnt;
+ /** allocating space for filterbank object */
+ aubio_filterbank_t * fb = AUBIO_NEW(aubio_filterbank_t);
+ fb->win_s=win_s;
+ fb->n_filters=n_filters;
+
+ /** allocating filter tables */
+ fb->filters=AUBIO_ARRAY(n_filters,f_vec_t);
+ for (filter_cnt=0; filter_cnt<n_filters; filter_cnt++)
+ /* considering one-channel filters */
+ filters[filter_cnt]=new_fvec(win_s, 1);
+
+}
+
+void del_aubio_filterbank(aubio_filterbank_t * fb){+
+ int filter_cnt;
+ /** deleting filter tables first */
+ for (filter_cnt=0; filter_cnt<fb->n_filters; filter_cnt++)
+ del_fvec(fb->filters[filter_cnt]);
+ AUBIO_FREE(fb->filters);
+ AUBIO_FREE(fb);
+
+}
+
+// Initialization
+
+void aubio_filterbank_mfcc_init(aubio_filterbank_t * fb, smpl_t nyquist, int style, smpl_t freq_min, smpl_t freq_max){+
+ int n, i, k, *fft_peak, M, next_peak;
+ smpl_t norm, mel_freq_max, mel_freq_min, norm_fact, height, inc, val,
+ freq_bw_mel, *mel_peak, *height_norm, *lin_peak;
+
+ mel_peak = height_norm = lin_peak = NULL;
+ fft_peak = NULL;
+ norm = 1;
+
+ mel_freq_max = 1127 * log(1 + freq_max / 700);
+ mel_freq_min = 1127 * log(1 + freq_min / 700);
+ freq_bw_mel = (mel_freq_max - mel_freq_min) / fb->n_filters;
+
+ mel_peak = (smpl_t *)malloc((fb->n_filters + 2) * sizeof(smpl_t));
+ /* +2 for zeros at start and end */
+ lin_peak = (smpl_t *)malloc((fb->n_filters + 2) * sizeof(smpl_t));
+ fft_peak = (int *)malloc((fb->n_filters + 2) * sizeof(int));
+ height_norm = (smpl_t *)malloc(fb->n_filters * sizeof(smpl_t));
+
+ if(mel_peak == NULL || height_norm == NULL ||
+ lin_peak == NULL || fft_peak == NULL)
+ return XTRACT_MALLOC_FAILED;
+
+ M = fb->win_s >> 1;
+
+ mel_peak[0] = mel_freq_min;
+ lin_peak[0] = 700 * (exp(mel_peak[0] / 1127) - 1);
+ fft_peak[0] = lin_peak[0] / nyquist * M;
+
+
+ for (n = 1; n <= fb->n_filters; n++){ + /*roll out peak locations - mel, linear and linear on fft window scale */
+ mel_peak[n] = mel_peak[n - 1] + freq_bw_mel;
+ lin_peak[n] = 700 * (exp(mel_peak[n] / 1127) -1);
+ fft_peak[n] = lin_peak[n] / nyquist * M;
+ }
+
+ for (n = 0; n < fb->n_filters; n++){+ /*roll out normalised gain of each peak*/
+ if (style == XTRACT_EQUAL_GAIN){+ height = 1;
+ norm_fact = norm;
+ }
+ else{+ height = 2 / (lin_peak[n + 2] - lin_peak[n]);
+ norm_fact = norm / (2 / (lin_peak[2] - lin_peak[0]));
+ }
+ height_norm[n] = height * norm_fact;
+ }
+
+ i = 0;
+
+ for(n = 0; n < fb->n_filters; n++){+
+ /*calculate the rise increment*/
+ if(n > 0)
+ inc = height_norm[n] / (fft_peak[n] - fft_peak[n - 1]);
+ else
+ inc = height_norm[n] / fft_peak[n];
+ val = 0;
+
+ /*zero the start of the array*/
+ for(k = 0; k < i; k++)
+ //fft_tables[n][k] = 0.f;
+ fb->filters[n]->data[0][k]=0.f;
+
+ /*fill in the rise */
+ for(; i <= fft_peak[n]; i++){ + // fft_tables[n][i] = val;
+ fb->filters[n]->data[0][k]=val;
+ val += inc;
+ }
+
+ /*calculate the fall increment */
+ inc = height_norm[n] / (fft_peak[n + 1] - fft_peak[n]);
+
+ val = 0;
+ next_peak = fft_peak[n + 1];
+
+ /*reverse fill the 'fall' */
+ for(i = next_peak; i > fft_peak[n]; i--){ + //fft_tables[n][i] = val;
+ fb->filters[n]->data[0][k]=val;
+ val += inc;
+ }
+
+ /*zero the rest of the array*/
+ for(k = next_peak + 1; k < fb->win_s; k++)
+ //fft_tables[n][k] = 0.f;
+ fb->filters[n]->data[0][k]=0.f;
+ }
+
+ free(mel_peak);
+ free(lin_peak);
+ free(height_norm);
+ free(fft_peak);
+
+ //return XTRACT_SUCCESS;
+
+}
+
+//to be deleted code
+
+
+// int 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){+//
+// int n, i, k, *fft_peak, M, next_peak;
+// smpl_t norm, mel_freq_max, mel_freq_min, norm_fact, height, inc, val,
+// freq_bw_mel, *mel_peak, *height_norm, *lin_peak;
+//
+// mel_peak = height_norm = lin_peak = NULL;
+// fft_peak = NULL;
+// norm = 1;
+//
+// mel_freq_max = 1127 * log(1 + freq_max / 700);
+// mel_freq_min = 1127 * log(1 + freq_min / 700);
+// freq_bw_mel = (mel_freq_max - mel_freq_min) / freq_bands;
+//
+// mel_peak = (smpl_t *)malloc((freq_bands + 2) * sizeof(smpl_t));
+// /* +2 for zeros at start and end */
+// lin_peak = (smpl_t *)malloc((freq_bands + 2) * sizeof(smpl_t));
+// fft_peak = (int *)malloc((freq_bands + 2) * sizeof(int));
+// height_norm = (smpl_t *)malloc(freq_bands * sizeof(smpl_t));
+//
+// if(mel_peak == NULL || height_norm == NULL ||
+// lin_peak == NULL || fft_peak == NULL)
+// return XTRACT_MALLOC_FAILED;
+//
+// M = N >> 1;
+//
+// mel_peak[0] = mel_freq_min;
+// lin_peak[0] = 700 * (exp(mel_peak[0] / 1127) - 1);
+// fft_peak[0] = lin_peak[0] / nyquist * M;
+//
+//
+// for (n = 1; n <= freq_bands; n++){ +// /*roll out peak locations - mel, linear and linear on fft window scale */
+// mel_peak[n] = mel_peak[n - 1] + freq_bw_mel;
+// lin_peak[n] = 700 * (exp(mel_peak[n] / 1127) -1);
+// fft_peak[n] = lin_peak[n] / nyquist * M;
+// }
+//
+// for (n = 0; n < freq_bands; n++){+// /*roll out normalised gain of each peak*/
+// if (style == XTRACT_EQUAL_GAIN){+// height = 1;
+// norm_fact = norm;
+// }
+// else{+// height = 2 / (lin_peak[n + 2] - lin_peak[n]);
+// norm_fact = norm / (2 / (lin_peak[2] - lin_peak[0]));
+// }
+// height_norm[n] = height * norm_fact;
+// }
+//
+// i = 0;
+//
+// for(n = 0; n < freq_bands; n++){+//
+// /*calculate the rise increment*/
+// if(n > 0)
+// inc = height_norm[n] / (fft_peak[n] - fft_peak[n - 1]);
+// else
+// inc = height_norm[n] / fft_peak[n];
+// val = 0;
+//
+// /*zero the start of the array*/
+// for(k = 0; k < i; k++)
+// fft_tables[n][k] = 0.f;
+//
+// /*fill in the rise */
+// for(; i <= fft_peak[n]; i++){ +// fft_tables[n][i] = val;
+// val += inc;
+// }
+//
+// /*calculate the fall increment */
+// inc = height_norm[n] / (fft_peak[n + 1] - fft_peak[n]);
+//
+// val = 0;
+// next_peak = fft_peak[n + 1];
+//
+// /*reverse fill the 'fall' */
+// for(i = next_peak; i > fft_peak[n]; i--){ +// fft_tables[n][i] = val;
+// val += inc;
+// }
+//
+// /*zero the rest of the array*/
+// for(k = next_peak + 1; k < N; k++)
+// fft_tables[n][k] = 0.f;
+// }
+//
+// free(mel_peak);
+// free(lin_peak);
+// free(height_norm);
+// free(fft_peak);
+//
+// return XTRACT_SUCCESS;
+//
+// }
--- /dev/null
+++ b/src/filterbank.h
@@ -1,0 +1,80 @@
+/*
+ Copyright (C) 2007 Amaury Hazan
+ adapted to aubio from LibXtract
+ http://libxtract.sourceforge.net/
+
+
+ This program is free software; you can redistribute it and/or modify
+ it under the terms of the GNU General Public License as published by
+ the Free Software Foundation; either version 2 of the License, or
+ (at your option) any later version.
+
+ This program is distributed in the hope that it will be useful,
+ but WITHOUT ANY WARRANTY; without even the implied warranty of
+ MERCHANTABILITY or FITNESS FOR A PARTICULAR PURPOSE. See the
+ GNU General Public License for more details.
+
+ You should have received a copy of the GNU General Public License
+ along with this program; if not, write to the Free Software
+ Foundation, Inc., 675 Mass Ave, Cambridge, MA 02139, USA.
+*/
+
+/** \file
+
+ Filterbank object
+
+ General-purpose spectral filterbank object. Comes with mel-filter initialization function.
+
+*/
+
+#ifndef AUBIOFILTERBANK_H
+#define AUBIOFILTERBANK_H
+
+#ifdef __cplusplus
+extern "C" {+#endif
+
+
+
+typedef struct aubio_filterbank_t_ aubio_filterbank_t;
+
+/** create filterbank object
+
+ \param win_s size of analysis buffer (and length the FFT transform)
+ \param n_filters number of filters to create
+
+*/
+
+aubio_filterbank_t * new_aubio_filterbank(uint_t n_filters, uint_t win_s);
+
+/** destroy filterbank object
+
+ \param fb filterbank, as returned by new_aubio_filterbank method
+
+*/
+void del_aubio_filterbank(aubio_filterbank_t * fb);
+
+/** filterbank initialization for mel filters
+
+ \param fb filterbank, as returned by new_aubio_filterbank method
+ \param nyquist nyquist frequency, i.e. half of the sampling rate
+ \param style libxtract style
+ \param freqmin lowest filter frequency
+ \param freqmax highest filter frequency
+
+*/
+void aubio_filterbank_mfcc_init(aubio_filterbank_t * fb, smpl_t nyquist, int style, smpl_t freq_min, smpl_t freq_max);
+
+// Initialization
+
+/** \brief A function to initialise a mel filter bank
+ *
+ * It is up to the caller to pass in a pointer to memory allocated for freq_bands arrays of length N. This function populates these arrays with magnitude coefficients representing the mel filterbank on a linear scale
+ */
+int 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);
+
+#ifdef __cplusplus
+}
+#endif
+
+#endif
--- /dev/null
+++ b/src/mfcc.c
@@ -1,0 +1,215 @@
+/*
+ Copyright (C) 2006 Amaury Hazan
+ Ported to aubio from LibXtract
+ http://libxtract.sourceforge.net/
+
+
+ This program is free software; you can redistribute it and/or modify
+ it under the terms of the GNU General Public License as published by
+ the Free Software Foundation; either version 2 of the License, or
+ (at your option) any later version.
+
+ This program is distributed in the hope that it will be useful,
+ but WITHOUT ANY WARRANTY; without even the implied warranty of
+ MERCHANTABILITY or FITNESS FOR A PARTICULAR PURPOSE. See the
+ GNU General Public License for more details.
+
+ You should have received a copy of the GNU General Public License
+ along with this program; if not, write to the Free Software
+ Foundation, Inc., 675 Mass Ave, Cambridge, MA 02139, USA.
+
+*/
+
+#include "aubio_priv.h"
+#include "sample.h"
+#include "fft.h"
+#include "filterbank.h"
+#include "mfcc.h"
+#include "math.h"
+
+
+
+/** Internal structure for mfcc object **/
+
+struct aubio_mfcc_t_{+
+ /** grain length */
+ uint_t win_s;
+
+ /** sample rate (needed?) */
+ uint_t samplerate;
+
+ /** number of channels */
+ uint_t channels;
+
+ /** filter bank */
+ aubio_filterbank_t * fb;
+
+ /** number of coefficients (= fb->n_filters/2 +1) */
+ uint_t n_coefs;
+
+ /** lowest frequency for filters */
+ smpl_t lowfreq;
+
+ /** highest frequency for filters */
+ smpl_t highfreq;
+
+ /** input buffer for dct * [fb->n_filters] */
+ fvec_t * in_dct;
+
+ /** fft object for dct */
+ aubio_mfft_t * fft_dct;
+
+ /** output buffer for dct */
+ cvec_t * fftgrain_dct;
+
+};
+
+
+aubio_mfcc_t * new_aubio_mfcc (uint_t win_s, uint_t samplerate ,uint_t n_coefs, smpl_t lowfreq, smpl_t highfreq, uint_t channels){+
+
+ /** allocating space for mfcc object */
+
+ aubio_mfcc_t * mfcc = AUBIO_NEW(aubio_mfcc_t);
+
+ mfcc->win_s=win_s;
+ mfcc->samplerate=samplerate;
+ mfcc->channels=channels;
+ mfcc->n_coefs=n_coefs;
+ mfcc->lowfreq=lowfreq;
+ mfcc->highfreq=highfreq;
+
+ /** filterbank allocation */
+ //we need (n_coefs-1)*2 filters to obtain n_coefs coefficients after dct
+ mfcc->fb=new_aubio_filterbank((n_coefs-1)*2, mfcc->win_s);
+
+ /** allocating space for fft object (used for dct) */
+ mfcc->fft_dct=new_aubio_mfft(mfcc->win_s, 1);
+
+ /** allocating buffers */
+
+ mfcc->in_dct=new_fvec(mfcc->win_s, 1);
+
+ mfcc->fftgrain_dct=new_cvec(mfcc->fb->n_filters, 1);
+
+ /** populating the filterbank */
+
+ aubio_filterbank_mfcc_init(mfcc->fb, (mfcc->samplerate)/2, XTRACT_EQUAL_GAIN, mfcc->lowfreq, mfcc->highfreq);
+
+ return mfcc;
+
+};
+
+
+void del_aubio_mfcc(aubio_mfcc_t *mf){+
+ /** deleting filterbank */
+ del_aubio_filterbank(mf->fb);
+ /** deleting mfft object */
+ del_aubio_mfft(mf->fft_dct);
+ /** deleting buffers */
+ del_fvec(mf->in_dct);
+ del_cvec(mf->fftgrain_dct);
+
+ /** deleting mfcc object */
+ AUBIO_FREE(mf);
+
+}
+
+
+// Computation
+
+void aubio_mfcc_do(aubio_mfcc_t * mf, cvec_t *in, fvec_t *out){+
+ aubio_filterbank_t *f = mf->fb;
+ uint_t n, filter_cnt;
+
+ for(filter_cnt = 0; filter_cnt < f->n_filters; filter_cnt++){+ mf->in_dct->data[0][filter_cnt] = 0.f;
+ for(n = 0; n < mf->win_s; n++){+ mf->in_dct->data[0][filter_cnt] += in->norm[0][n] * f->filters[filter_cnt]->data[0][n];
+ }
+ mf->in_dct->data[0][filter_cnt] = LOG(mf->in_dct->data[0][filter_cnt] < XTRACT_LOG_LIMIT ? XTRACT_LOG_LIMIT : mf->in_dct->data[0][filter_cnt]);
+ }
+
+ //TODO: check that zero padding
+ // the following line seems useless since the in_dct buffer has the correct size
+ //for(n = filter + 1; n < N; n++) result[n] = 0;
+
+ aubio_dct_do(mf, mf->in_dct, out);
+
+ //return XTRACT_SUCCESS;
+}
+
+void aubio_dct_do(aubio_mfcc_t * mf, fvec_t *in, fvec_t *out){+
+
+
+ //fvec_t * momo = new_fvec(20, 1);
+ //momo->data = data;
+
+ //compute mag spectrum
+ aubio_mfft_do (mf->fft_dct, in, mf->fftgrain_dct);
+
+ int i;
+ //extract real part of fft grain
+ for(i=0; i<mf->n_coefs ;i++){+ out->data[0][i]= mf->fftgrain_dct->norm[0][i]*COS(mf->fftgrain_dct->phas[0][i]);
+ }
+
+
+ //return XTRACT_SUCCESS;
+}
+
+
+///////// OLD CODE
+
+// int aubio_mfcc_do(const float *data, const int N, const void *argv, float *result, aubio_mfft_t * fft_dct, cvec_t * fftgrain_dct){+//
+// aubio_mel_filter *f;
+// uint_t n, filter;
+//
+// f = (aubio_mel_filter *)argv;
+// printf("%d",f->n_filters);+//
+// for(filter = 0; filter < f->n_filters; filter++){+// result[filter] = 0.f;
+// for(n = 0; n < N; n++){+// result[filter] += data[n] * f->filters[filter][n];
+// }
+// result[filter] = LOG(result[filter] < XTRACT_LOG_LIMIT ? XTRACT_LOG_LIMIT : result[filter]);
+// }
+//
+// //TODO: check that zero padding
+// for(n = filter + 1; n < N; n++) result[n] = 0;
+//
+// aubio_dct_do(result, f->n_filters, NULL, result, fft_dct, fftgrain_dct);
+//
+// return XTRACT_SUCCESS;
+// }
+
+// Added last two arguments to be able to pass from example
+
+// int aubio_dct_do(const float *data, const int N, const void *argv, float *result, aubio_mfft_t * fft_dct, cvec_t * fftgrain_dct){+//
+//
+// //call aubio p_voc in dct setting
+//
+// //TODO: fvec as input? Remove data length, N?
+//
+// fvec_t * momo = new_fvec(20, 1);
+// momo->data = data;
+//
+// //compute mag spectrum
+// aubio_mfft_do (fft_dct, data, fftgrain_dct);
+//
+// int i;
+// //extract real part of fft grain
+// for(i=0; i<N ;i++){+// result[i]= fftgrain_dct->norm[0][i]*COS(fftgrain_dct->phas[0][i]);
+// }
+//
+//
+// return XTRACT_SUCCESS;
+// }
--- /dev/null
+++ b/src/mfcc.h
@@ -1,0 +1,215 @@
+/*
+ Copyright (C) 2006 Amaury Hazan
+ Ported to aubio from LibXtract
+ http://libxtract.sourceforge.net/
+
+
+ This program is free software; you can redistribute it and/or modify
+ it under the terms of the GNU General Public License as published by
+ the Free Software Foundation; either version 2 of the License, or
+ (at your option) any later version.
+
+ This program is distributed in the hope that it will be useful,
+ but WITHOUT ANY WARRANTY; without even the implied warranty of
+ MERCHANTABILITY or FITNESS FOR A PARTICULAR PURPOSE. See the
+ GNU General Public License for more details.
+
+ You should have received a copy of the GNU General Public License
+ along with this program; if not, write to the Free Software
+ Foundation, Inc., 675 Mass Ave, Cambridge, MA 02139, USA.
+
+*/
+
+#ifndef MFCC_H
+#define MFCC_H
+
+#ifdef __cplusplus
+extern "C" {+#endif
+
+#include "filterbank.h"
+
+//libXtract constants and enums
+// TODO: remove them
+
+#define XTRACT_SQ(a) ((a) * (a))
+#define XTRACT_MIN(a, b) ((a) < (b) ? (a) : (b))
+#define XTRACT_MAX(a, b) ((a) > (b) ? (a) : (b))
+#define XTRACT_NEEDS_FFTW printf("LibXtract must be compiled with fftw support to use this function.\n")+#define XTRACT_VERY_SMALL_NUMBER 2e-42
+#define XTRACT_LOG_LIMIT XTRACT_VERY_SMALL_NUMBER
+#define XTRACT_LOG_LIMIT_DB -96.0
+#define XTRACT_DB_SCALE_OFFSET 96.0
+#define XTRACT_VERY_BIG_NUMBER 2e42
+#define XTRACT_SR_UPPER_LIMIT 192000.0
+#define XTRACT_SR_LOWER_LIMIT 22050.0
+#define XTRACT_SR_DEFAULT 44100.0
+#define XTRACT_FUNDAMENTAL_DEFAULT 440.0
+#define XTRACT_CHECK_nyquist if(!nyquist) nyquist = XTRACT_SR_DEFAULT / 2
+#define XTRACT_CHECK_q if(!q) q = XTRACT_SR_DEFAULT / N
+#define XTRACT_IS_ODD(x) (x % 2 != 0 ? 1 : 0)
+#define XTRACT_SR_LIMIT SR_UPPER_LIMIT
+#define XTRACT_FFT_BANDS_MIN 16
+#define XTRACT_FFT_BANDS_MAX 65536
+#define XTRACT_FFT_BANDS_DEF 1024
+#define XTRACT_SPEC_BW_MIN 0.168 /* Minimum spectral bandwidth \
+ (= SR_LOWER_LIMIT / FFT_BANDS_MAX*/
+#define XTRACT_SPEC_BW_MAX 12000.0 /* SR_UPPER_LIMIT / FFT_BANDS_MIN */
+#define XTRACT_SPEC_BW_DEF 43.066 /* SR_DEFAULT / FFT_BANDS_DEF */
+
+/** \brief Enumeration of feature initialisation functions */
+enum xtract_feature_init_ {+ XTRACT_INIT_MFCC = 100,
+ XTRACT_INIT_BARK
+};
+
+/** \brief Enumeration of feature types */
+enum xtract_feature_types_ {+ XTRACT_SCALAR,
+ XTRACT_VECTOR,
+ XTRACT_DELTA
+};
+
+/** \brief Enumeration of mfcc types */
+enum xtract_mfcc_types_ {+ XTRACT_EQUAL_GAIN,
+ XTRACT_EQUAL_AREA
+};
+
+/** \brief Enumeration of return codes */
+enum xtract_return_codes_ {+ XTRACT_SUCCESS,
+ XTRACT_MALLOC_FAILED,
+ XTRACT_BAD_ARGV,
+ XTRACT_BAD_VECTOR_SIZE,
+ XTRACT_NO_RESULT,
+ XTRACT_FEATURE_NOT_IMPLEMENTED
+};
+
+/** \brief Enumeration of spectrum types */
+enum xtract_spectrum_ {+ XTRACT_MAGNITUDE_SPECTRUM,
+ XTRACT_LOG_MAGNITUDE_SPECTRUM,
+ XTRACT_POWER_SPECTRUM,
+ XTRACT_LOG_POWER_SPECTRUM
+};
+
+/** \brief Enumeration of data types*/
+typedef enum type_ {+ XTRACT_FLOAT,
+ XTRACT_FLOATARRAY,
+ XTRACT_INT,
+ XTRACT_MEL_FILTER
+} xtract_type_t;
+
+/** \brief Enumeration of units*/
+typedef enum unit_ {+ /* NONE, ANY */
+ XTRACT_HERTZ = 2,
+ XTRACT_ANY_AMPLITUDE_HERTZ,
+ XTRACT_DBFS,
+ XTRACT_DBFS_HERTZ,
+ XTRACT_PERCENT,
+ XTRACT_SONE
+} xtract_unit_t;
+
+/** \brief Boolean */
+typedef enum {+ XTRACT_FALSE,
+ XTRACT_TRUE
+} xtract_bool_t;
+
+/** \brief Enumeration of vector format types*/
+typedef enum xtract_vector_ {+ /* N/2 magnitude/log-magnitude/power/log-power coeffs and N/2 frequencies */
+ XTRACT_SPECTRAL,
+ /* N spectral amplitudes */
+ XTRACT_SPECTRAL_MAGNITUDES,
+ /* N/2 magnitude/log-magnitude/power/log-power peak coeffs and N/2
+ * frequencies */
+ XTRACT_SPECTRAL_PEAKS,
+ /* N spectral peak amplitudes */
+ XTRACT_SPECTRAL_PEAKS_MAGNITUDES,
+ /* N spectral peak frequencies */
+ XTRACT_SPECTRAL_PEAKS_FREQUENCIES,
+ /* N/2 magnitude/log-magnitude/power/log-power harmonic peak coeffs and N/2
+ * frequencies */
+ XTRACT_SPECTRAL_HARMONICS,
+ /* N spectral harmonic amplitudes */
+ XTRACT_SPECTRAL_HARMONICS_MAGNITUDES,
+ /* N spectral harmonic frequencies */
+ XTRACT_SPECTRAL_HARMONICS_FREQUENCIES,
+ XTRACT_ARBITRARY_SERIES,
+ XTRACT_AUDIO_SAMPLES,
+ XTRACT_MEL_COEFFS,
+ XTRACT_BARK_COEFFS,
+ XTRACT_NO_DATA
+} xtract_vector_t;
+
+
+typedef struct aubio_mfcc_t_ aubio_mfcc_t;
+
+// Creation
+
+/** create mfcc object
+
+ \param win_s size of analysis buffer (and length the FFT transform)
+ \param samplerate
+ \param n_coefs: number of desired coefs
+ \param lowfreq: lowest frequency to use in filterbank
+ \param highfreq highest frequency to use in filterbank
+ \param channels number of channels
+
+*/
+aubio_mfcc_t * new_aubio_mfcc (uint_t win_s, uint_t samplerate ,uint_t n_coefs, smpl_t lowfreq, smpl_t highfreq, uint_t channels);
+
+// Deletion
+
+/** delete mfcc object
+
+ \param mf mfcc object as returned by new_aubio_mfcc
+
+*/
+void del_aubio_mfcc(aubio_mfcc_t *mf);
+
+// Process
+
+/** mfcc object processing
+
+ \param mf mfcc object as returned by new_aubio_mfcc
+ \param in input spectrum (win_s long)
+ \param out output mel coefficients buffer (n_filters/2 +1 long)
+
+*/
+
+void aubio_mfcc_do(aubio_mfcc_t * mf, cvec_t *in, fvec_t *out);
+
+/** intermediate dct involved in aubio_mfcc_do
+
+ \param mf mfcc object as returned by new_aubio_mfcc
+ \param in input spectrum (n_filters long)
+ \param out output mel coefficients buffer (n_filters/2 +1 long)
+
+*/
+
+void aubio_dct_do(aubio_mfcc_t * mf, fvec_t *in, fvec_t *out);
+
+
+
+
+//old code
+
+
+/*
+int aubio_mfcc_do(const float *data, const int N, const void *argv, float *result, aubio_mfft_t *fft_dct, cvec_t *fftgrain_dct);
+
+int aubio_dct_do(const float *data, const int N, const void *argv, float *result, aubio_mfft_t *fft_dct, cvec_t *fftgrain_dct);*/
+
+
+
+
+#ifdef __cplusplus
+}
+#endif
+
+#endif
--
⑨