added mffc + bark
authorahazan <ahazan@mtg55>
Wed, 5 Sep 2007 16:46:29 +0000 (18:46 +0200)
committerahazan <ahazan@mtg55>
Wed, 5 Sep 2007 16:46:29 +0000 (18:46 +0200)
examples/aubiomfcc.c [new file with mode: 0644]
src/aubiofilterbank.c [new file with mode: 0644]
src/aubiofilterbank.h [new file with mode: 0644]
src/bark.c [new file with mode: 0644]
src/bark.h [new file with mode: 0644]
src/mfcc.c [new file with mode: 0644]
src/mfcc.h [new file with mode: 0644]

diff --git a/examples/aubiomfcc.c b/examples/aubiomfcc.c
new file mode 100644 (file)
index 0000000..cca6605
--- /dev/null
@@ -0,0 +1,82 @@
+/*
+   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);
+      
+      //TODO: extract Magnitude buffer f_fvec from fftgrain cvec
+
+      //compute mfccs
+      aubio_mffc_do (magbuf, nframes, filterbank, outbuf);
+      
+      
+
+      /* 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);
+  examples_common_process(aubio_process,process_print);
+  examples_common_del();
+  debug("End of program.\n");
+  fflush(stderr);
+  return 0;
+}
+
diff --git a/src/aubiofilterbank.c b/src/aubiofilterbank.c
new file mode 100644 (file)
index 0000000..a4fa199
--- /dev/null
@@ -0,0 +1,123 @@
+/*
+   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 "aubiofilterbank.h"
+
+// Initialization
+
+int aubio_mfcc_init(int N, float nyquist, int style, float freq_min, float freq_max, int freq_bands, float **fft_tables){
+
+    int n, i, k, *fft_peak, M, next_peak; 
+    float 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 = (float *)malloc((freq_bands + 2) * sizeof(float)); 
+    /* +2 for zeros at start and end */
+    lin_peak = (float *)malloc((freq_bands + 2) * sizeof(float));
+    fft_peak = (int *)malloc((freq_bands + 2) * sizeof(int));
+    height_norm = (float *)malloc(freq_bands * sizeof(float));
+
+    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;
+
+}
\ No newline at end of file
diff --git a/src/aubiofilterbank.h b/src/aubiofilterbank.h
new file mode 100644 (file)
index 0000000..acd32e5
--- /dev/null
@@ -0,0 +1,44 @@
+/*
+   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.
+
+*/
+
+#ifndef AUBIOFILTERBANK_H
+#include AUBIOFILTERBANK_H
+
+#define NYQUIST 22050.f
+
+// Struct Declaration
+
+/** \brief A structure to store a set of n_filters Mel filters */
+typedef struct aubio_mel_filter_ {
+    int n_filters;
+    float **filters;
+} aubio_mel_filter;
+
+// 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, float nyquist, int style, float freq_min, float freq_max, int freq_bands, float **fft_tables);
+
+#endif
\ No newline at end of file
diff --git a/src/bark.c b/src/bark.c
new file mode 100644 (file)
index 0000000..652f125
--- /dev/null
@@ -0,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
diff --git a/src/bark.h b/src/bark.h
new file mode 100644 (file)
index 0000000..7c03394
--- /dev/null
@@ -0,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
diff --git a/src/mfcc.c b/src/mfcc.c
new file mode 100644 (file)
index 0000000..4205b2c
--- /dev/null
@@ -0,0 +1,62 @@
+/*
+   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 "mffc.h"
+#include "aubiofilterbank.h"
+
+// Computation
+
+int aubio_mfcc_do(const float *data, const int N, const void *argv, float *result){
+
+    aubio_mel_filter *f;
+    int n, filter;
+
+    f = (aubio_mel_filter *)argv;
+    
+    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]);
+    }
+
+    for(n = filter + 1; n < N; n++) result[n] = 0; 
+    
+    aubio_dct_do(result, f->n_filters, NULL, result);
+    
+    return XTRACT_SUCCESS;
+}
+
+int aubio_dct_do(const float *data, const int N, const void *argv, float *result){
+    
+    fftwf_plan plan;
+    
+    plan = 
+        fftwf_plan_r2r_1d(N, (float *) data, result, FFTW_REDFT00, FFTW_ESTIMATE);
+    
+    fftwf_execute(plan);
+    fftwf_destroy_plan(plan);
+
+    return XTRACT_SUCCESS;
+}
\ No newline at end of file
diff --git a/src/mfcc.h b/src/mfcc.h
new file mode 100644 (file)
index 0000000..03a9c5d
--- /dev/null
@@ -0,0 +1,145 @@
+/*
+   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 
+
+#define NYQUIST 22050.f
+
+//libXtract enums
+// TODO: remove them 
+
+/** \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;
+
+
+
+
+// Computation
+
+/** \brief Extract Mel Frequency Cepstral Coefficients based on a method described by Rabiner
+ * 
+ * \param *data: a pointer to the first element in an array of spectral magnitudes, e.g. the first half of the array pointed to by *resul from xtract_spectrum()
+ * \param N: the number of array elements to be considered
+ * \param *argv: a pointer to a data structure of type xtract_mel_filter, containing n_filters coefficient tables to make up a mel-spaced filterbank
+ * \param *result: a pointer to an array containing the resultant MFCC
+ * 
+ * The data structure pointed to by *argv must be obtained by first calling xtract_init_mfcc
+ */
+int aubio_mfcc_do(const float *data, const int N, const void *argv, float *result);
+
+/** \brief Extract the Discrete Cosine transform of a time domain signal
+ * \param *data: a pointer to the first element in an array of floats representing an audio vector
+ * \param N: the number of array elements to be considered
+ * \param *argv: a pointer to NULL 
+ * \param *result: a pointer to an array containing resultant dct coefficients
+ */
+int aubio_dct_do(const float *data, const int N, const void *argv, float *result);
+
+
+#endif
\ No newline at end of file