changeset 107:3e648eec95cb

- Added new helper functions: xtract_windowed() and xtract_features_from_subframes()
- Added windowing functions (window.c)

--- a/src/Makefile.am	Thu Dec 27 20:37:15 2007 +0000
+++ b/src/Makefile.am	Fri Dec 28 19:34:51 2007 +0000
@@ -1,6 +1,7 @@
 MAINTAINERCLEANFILES = Makefile.in
 
-SOURCES = libxtract.c descriptors.c scalar.c vector.c delta.c init.c fini.c
+SOURCES = libxtract.c descriptors.c scalar.c vector.c delta.c init.c\
+	window.c fini.c helper.c
 
 if BUILD_FFT
 FFT_DEFINE = -DXTRACT_FFT

--- a/src/delta.c	Thu Dec 27 20:37:15 2007 +0000
+++ b/src/delta.c	Fri Dec 28 19:34:51 2007 +0000
@@ -32,8 +32,7 @@
 
 int xtract_lnorm(const float *data, const int N, const void *argv , float *result){
 
-    int feature,
-        n,
+    int n,
         type;
 
     float order,
@@ -80,7 +79,7 @@
 
 int xtract_difference_vector(const float *data, const int N, const void *argv, float *result){
 
-    float *frame1,
+    const float *frame1,
           *frame2;
 
     int n;

--- /dev/null	Thu Jan 01 00:00:00 1970 +0000
+++ b/src/helper.c	Fri Dec 28 19:34:51 2007 +0000
@@ -0,0 +1,65 @@
+/* libxtract feature extraction library
+ *  
+ * Copyright (C) 2006 Jamie Bullock
+ *
+ * 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., 51 Franklin Street, Fifth Floor, Boston, MA  02110-1301, 
+ * USA.
+ */
+
+
+/* helper.c: helper functions. */
+
+#include "xtract/libxtract.h"
+
+int xtract_windowed(const float *data, const int N, const void *argv, float *result){
+
+    int n;
+    const float *window;
+
+    n = N;
+    window = (const float *)argv;
+
+    while(n--)
+        result[n] = data[n] * window[n];
+
+    return XTRACT_SUCCESS;
+
+}
+
+int xtract_features_from_subframes(const float *data, const int N, const int feature, const void *argv, float *result){
+
+    const float *frame1,
+          *frame2;
+    float *result1,
+          *result2;
+
+    int n,
+        rv;
+
+    n = N >> 1;
+
+    frame1 = data;
+    frame2 = data + n;
+    result1 = result;
+    result2 = result + n;
+
+    rv = xtract[feature](frame1, n, argv, result1);
+
+    if(rv == XTRACT_SUCCESS)
+        rv = xtract[feature](frame2, n, argv, result2);
+
+    return rv;
+
+}

--- a/src/init.c	Thu Dec 27 20:37:15 2007 +0000
+++ b/src/init.c	Fri Dec 28 19:34:51 2007 +0000
@@ -18,7 +18,7 @@
  * USA.
  */
 
-/* init.c: defines functions that extract a feature as a single value from an input vector */
+/* init.c: defines initialisation and free functions. Also contains library constructor routine. */
 
 #ifdef HAVE_CONFIG_H
 #    include <config.h>
@@ -28,6 +28,7 @@
 #include <stdlib.h>
 
 #include "xtract/libxtract.h"
+#include "xtract_window_private.h"
 #define DEFINE_GLOBALS
 #include "xtract_globals_private.h"
 
@@ -43,7 +44,7 @@
 
     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;
+          freq_bw_mel, *mel_peak, *height_norm, *lin_peak;
 
     mel_peak = height_norm = lin_peak = NULL;
     fft_peak = NULL;
@@ -60,9 +61,9 @@
     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;
-    
+            lin_peak == NULL || fft_peak == NULL)
+        return XTRACT_MALLOC_FAILED;
+
     M = N >> 1;
 
     mel_peak[0] = mel_freq_min;
@@ -91,41 +92,41 @@
     }
 
     i = 0;
-   
+
     for(n = 0; n < freq_bands; n++){
-	
-		// calculate the rise increment
+
+        // calculate the rise increment
         if(n==0)
             inc = height_norm[n] / fft_peak[n];
         else
             inc = height_norm[n] / (fft_peak[n] - fft_peak[n - 1]);
         val = 0;	
-	
-		// zero the start of the array
-		for(k = 0; k < i; k++)
-		   fft_tables[n][k] = 0.f;
-	
-		// fill in the rise
+
+        // 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' 
+        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;
+        // zero the rest of the array
+        for(k = next_peak + 1; k < N; k++)
+            fft_tables[n][k] = 0.f;
     }
 
 
@@ -145,9 +146,9 @@
 
     float *input, *output;
     int optimisation;
-    
+
     input = output = NULL;
-    
+
     fprintf(stderr, "Optimisation level: %d\n", XTRACT_FFT_OPTIMISATION_LEVEL);
 
     if(XTRACT_FFT_OPTIMISATION_LEVEL == 0)
@@ -209,14 +210,63 @@
     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 */
+    /*FIX shohuld use rounding, but couldn't get it to work */
 
     return XTRACT_SUCCESS;
 }
 
+float *xtract_init_window(const int N, const int type){
+
+    float *window;
+
+    window = malloc(N * sizeof(float));
+
+    switch (type) {
+        case XTRACT_GAUSS:
+            gauss(window, N, 0.4);
+            break;
+        case XTRACT_HAMMING:
+            hamming(window, N);
+            break;
+        case XTRACT_HANN:
+            hann(window, N);
+            break;
+        case XTRACT_BARTLETT:
+            bartlett(window, N);
+            break;
+        case XTRACT_TRIANGULAR:
+            triangular(window, N);
+            break;
+        case XTRACT_BARTLETT_HANN:
+            bartlett_hann(window, N);
+            break;
+        case XTRACT_BLACKMAN:
+            blackman(window, N);
+            break;
+        case XTRACT_KAISER:
+            kaiser(window, N, 3 * PI);
+            break;
+        case XTRACT_BLACKMAN_HARRIS:
+            blackman_harris(window, N);
+            break;
+        default:
+            hann(window, N);
+            break;
+    }
+
+    return window;
+
+}
+
+void xtract_free_window(float *window){
+
+    free(window);
+
+}
+
 #ifdef __GNUC__
 __attribute__((constructor)) void init()
 #else

--- a/src/libxtract.c	Thu Dec 27 20:37:15 2007 +0000
+++ b/src/libxtract.c	Fri Dec 28 19:34:51 2007 +0000
@@ -64,13 +64,13 @@
     xtract_hps,
     xtract_f0,
     xtract_failsafe_f0,
-/* xtract_delta.h */
+    /* xtract_delta.h */
     xtract_lnorm,
     xtract_flux,
     xtract_attack_time,
     xtract_decay_time,
     xtract_difference_vector,
-/* xtract_vector.h */
+    /* xtract_vector.h */
     xtract_autocorrelation,
     xtract_amdf,
     xtract_asdf,
@@ -82,6 +82,8 @@
     xtract_dct,
     xtract_harmonic_spectrum,
     xtract_lpc,
-    xtract_lpcc
+    xtract_lpcc,
+    /* xtract_helper.h */
+    xtract_windowed
 };
 

--- a/src/scalar.c	Thu Dec 27 20:37:15 2007 +0000
+++ b/src/scalar.c	Fri Dec 28 19:34:51 2007 +0000
@@ -19,7 +19,7 @@
  */
 
 
-/* xtract_scalar.c: defines functions that extract a feature as a single value from an input vector */
+/* scalar.c: defines functions that extract a feature as a single value from an input vector */
 
 #include "xtract/libxtract.h"
 #include "xtract_macros_private.h"

--- /dev/null	Thu Jan 01 00:00:00 1970 +0000
+++ b/src/window.c	Fri Dec 28 19:34:51 2007 +0000
@@ -0,0 +1,182 @@
+/* libxtract feature extraction library
+ *  
+ * Copyright (C) 2006 Jamie Bullock
+ *
+ * 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., 51 Franklin Street, Fifth Floor, Boston, MA  02110-1301, 
+ * USA.
+ */
+
+/* window.c: defines window generation functions (formulae courtesy of Wikipedia (http://en.wikipedia.org/wiki/Window_function) */
+
+#include <math.h>
+
+#include "xtract_window_private.h"
+
+void gauss(float *window, const int N, const float sd){
+
+    int n;
+    const float M = N - 1;
+    float num,
+          den,
+          exponent;
+
+    for (n = 0; n < N; n++) {
+
+        num = n - M / 2.f;
+        den = sd * M / 2.f;
+        
+        exponent = -0.5 * powf(num / den, 2);
+
+        window[n] = exp(exponent);
+
+    }
+}
+
+void hamming(float *window, const int N){
+
+    int n;
+    const float M = N - 1;
+
+    for (n = 0; n < N; n++)
+        window[n] = 0.53836 - (0.46164 * cosf(2.0 * PI * (float)n / M));
+
+}
+
+void hann(float *window, const int N){
+
+    int n;
+    const float M = N - 1;
+
+    for (n = 0; n < N; n++)
+        window[n] = 0.5 * (1.0 - cosf(2.0 * PI * (float)n / M));
+
+}
+
+void bartlett(float *window, const int N){
+
+    int n;
+    const float M = N - 1;
+
+    for (n = 0; n < N; n++)
+        window[n] = 2.f / M * (M / 2.f - fabsf(n - M / 2.f));
+
+}
+
+void triangular(float *window, const int N){
+
+    int n;
+    const float M = N - 1;
+
+    for (n = 0; n < N; n++)
+        window[n] = 2.f / N * (N / 2.f - fabsf(n - M / 2.f));
+}
+
+void bartlett_hann(float *window, const int N){
+
+    int n;
+    const float M = N - 1,
+          a0 = 0.62,
+          a1 = 0.5,
+          a2 = 0.38;
+    float term1 = 0.f,
+          term2 = 0.f;
+
+    for (n = 0; n < N; n++){
+
+        term1 = a1 * fabsf(n / M - 0.5);
+        term2 = a2 * cosf(2.0 * PI * (float)n / M);
+
+        window[n] = a0 - term1 - term2;
+    }
+}
+
+void blackman(float *window, const int N){
+
+    int n;
+    const float M = N - 1,
+          a0 = 0.42,
+          a1 = 0.5,
+          a2 = 0.08;
+    float term1 = 0.f,
+          term2 = 0.f;
+
+    for (n = 0; n < N; n++) {
+    
+        term1 = a1 * cosf(2.0 * PI * (float)n / M);
+        term2 = a2 * cosf(4.0 * PI * (float)n / M);
+
+        window[n] = a0 - term1 + term2;
+    }
+}
+
+#define BIZ_EPSILON 1E-21 // Max error acceptable 
+
+/* Based on code from mplayer window.c, and somewhat beyond me */
+float besselI0(float x){
+
+  float temp;
+  float sum   = 1.0;
+  float u     = 1.0;
+  float halfx = x/2.0;
+  int      n     = 1;
+
+  do {
+
+    temp = halfx/(float)n;
+    u *=temp * temp;
+    sum += u;
+    n++;
+
+  } while (u >= BIZ_EPSILON * sum);
+
+  return(sum);
+
+}
+
+void kaiser(float *window, const int N, const float alpha){
+
+    int n;
+    const float M = N - 1;
+    float num;
+
+    for (n = 0; n < N; n++) {
+
+        num = besselI0(alpha * sqrtf(1.0 - powf((2.0 * n / M - 1), 2)));
+        window[n] = num / besselI0(alpha);
+        
+    }
+}
+
+void blackman_harris(float *window, const int N){
+
+    int n;
+    const float M = N - 1,
+          a0 = 0.35875,
+          a1 = 0.48829,
+          a2 = 0.14128,
+          a3 = 0.01168;
+    float term1 = 0.f,
+          term2 = 0.f,
+          term3 = 0.f;
+
+    for (n = 0; n < N; n++) {
+
+        term1 = a1 * cosf(2.0 * PI * n / M);
+        term2 = a2 * cosf(4.0 * PI * n / M);
+        term3 = a3 * cosf(6.0 * PI * n / M);
+
+        window[n] = a0 - term1 + term2 - term3;
+    }
+}

--- /dev/null	Thu Jan 01 00:00:00 1970 +0000
+++ b/src/xtract_window_private.h	Fri Dec 28 19:34:51 2007 +0000
@@ -0,0 +1,98 @@
+/* libxtract feature extraction library
+ *  
+ * Copyright (C) 2006 Jamie Bullock
+ *
+ * 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., 51 Franklin Street, Fifth Floor, Boston, MA  02110-1301, 
+ * USA.
+ */
+
+/* xtract_window_private.h: declares window generation functions */
+
+#define PI 3.1415926535897931
+
+/** \brief generate a Gaussian window
+ *
+ * \param *window a pointer to an array to contain the window data
+ * \param N the number of elements in the array pointed to by *window
+ * \param sd the standard deviation of the "distribution" represented by the Gaussian curve. The higher the value of sd, the wider the curve. Generally sd <= 0.5
+ *
+ */
+void gauss(float *window, const int N, const float sd);
+
+/** \brief generate a Hamming window
+ *
+ * \param *window a pointer to an array to contain the window data
+ * \param N the number of elements in the array pointed to by *window
+ *
+ */
+void hamming(float *window, const int N);
+
+/** \brief generate a Hann window
+ *
+ * \param *window a pointer to an array to contain the window data
+ * \param N the number of elements in the array pointed to by *window
+ *
+ */
+void hann(float *window, const int N);
+
+/** \brief generate a Bartlett window
+ *
+ * \param *window a pointer to an array to contain the window data
+ * \param N the number of elements in the array pointed to by *window
+ *
+ */
+void bartlett(float *window, const int N);
+
+/** \brief generate a Triangular window
+ *
+ * \param *window a pointer to an array to contain the window data
+ * \param N the number of elements in the array pointed to by *window
+ *
+ */
+void triangular(float *window, const int N);
+
+/** \brief generate a Bartlett-Hann  window
+ *
+ * \param *window a pointer to an array to contain the window data
+ * \param N the number of elements in the array pointed to by *window
+ *
+ */
+void bartlett_hann(float *window, const int N);
+
+/** \brief generate a Blackman window
+ *
+ * \param *window a pointer to an array to contain the window data
+ * \param N the number of elements in the array pointed to by *window
+ *
+ */
+void blackman(float *window, const int N);
+
+/** \brief generate a Kaiser window
+ *
+ * \param *window a pointer to an array to contain the window data
+ * \param N the number of elements in the array pointed to by *window
+ * \param alpha The larger the value of |alpha|, the narrower the window becomes
+ *
+ */
+void kaiser(float *window, const int N, const float alpha);
+
+/** \brief generate a Blackman-Harris window
+ *
+ * \param *window a pointer to an array to contain the window data
+ * \param N the number of elements in the array pointed to by *window
+ *
+ */
+void blackman_harris(float *window, const int N);
+

--- a/swig/xtract.i	Thu Dec 27 20:37:15 2007 +0000
+++ b/swig/xtract.i	Fri Dec 28 19:34:51 2007 +0000
@@ -5,6 +5,7 @@
 %{
 #include "xtract/xtract_scalar.h"
 #include "xtract/xtract_vector.h"
+#include "xtract/xtract_helper.h"
 #include "xtract/libxtract.h"
 %}
 
@@ -19,4 +20,5 @@
 %clear float *result;
 
 %include "xtract/xtract_vector.h"
+%include "xtract/xtract_helper.h"
 %include "xtract/libxtract.h"

--- a/xtract/Makefile.am	Thu Dec 27 20:37:15 2007 +0000
+++ b/xtract/Makefile.am	Fri Dec 28 19:34:51 2007 +0000
@@ -3,5 +3,5 @@
 libxtractdir = $(includedir)/xtract
 
 libxtract_HEADERS = libxtract.h xtract_macros.h xtract_types.h xtract_delta.h \
-		    xtract_scalar.h  xtract_vector.h
+		    xtract_scalar.h  xtract_vector.h xtract_helper.h
 

--- a/xtract/libxtract.h	Thu Dec 27 20:37:15 2007 +0000
+++ b/xtract/libxtract.h	Fri Dec 28 19:34:51 2007 +0000
@@ -60,6 +60,7 @@
 #include "xtract_delta.h"
 #include "xtract_types.h"
 #include "xtract_macros.h"
+#include "xtract_helper.h"
 
 /** \defgroup libxtract API
   *
@@ -67,7 +68,7 @@
   * @{
   */
 
-#define XTRACT_FEATURES 57
+#define XTRACT_FEATURES 58
     
 /** \brief Enumeration of features, elements are used as indixes to an array of pointers to feature extracton functions */
 enum xtract_features_ {
@@ -127,7 +128,9 @@
     XTRACT_DCT,
     XTRACT_HARMONIC_SPECTRUM,
     XTRACT_LPC,
-    XTRACT_LPCC
+    XTRACT_LPCC,
+    /* Helper functions */
+    XTRACT_WINDOWED
 };
 
 /** \brief Enumeration of feature initialisation functions */
@@ -197,6 +200,19 @@
     XTRACT_TRUE
 } xtract_bool_t;
 
+/** \brief Window types */
+enum xtract_window_types_ {
+    XTRACT_GAUSS,
+    XTRACT_HAMMING,
+    XTRACT_HANN,
+    XTRACT_BARTLETT,
+    XTRACT_TRIANGULAR,
+    XTRACT_BARTLETT_HANN,
+    XTRACT_BLACKMAN,
+    XTRACT_KAISER,
+    XTRACT_BLACKMAN_HARRIS
+};
+
 /** \brief Enumeration of vector format types*/
 typedef enum xtract_vector_ {
     /* N/2 magnitude/log-magnitude/power/log-power coeffs and N/2 frequencies */
@@ -359,6 +375,21 @@
  */
 int xtract_init_fft(int N, int feature_name);
 
+/** \brief Make a window of a given type and return a pointer to it
+ *
+ * \param N: the size of the window
+ * \param type: the type of the window as given in the enumeration window_types_
+ *
+ */
+float *xtract_init_window(const int N, const int type);
+
+/** \brief Free a window as allocated by xtract_make_window() 
+ * 
+ * \param *window: a pointer to an array of floats as allocated by xtract_make_window()
+ *
+ */
+void xtract_free_window(float *window);
+
 /* \brief A function to build an array of function descriptors */
 void *xtract_make_descriptors();
 

--- /dev/null	Thu Jan 01 00:00:00 1970 +0000
+++ b/xtract/xtract_helper.h	Fri Dec 28 19:34:51 2007 +0000
@@ -0,0 +1,72 @@
+/* libxtract feature extraction library
+ *  
+ * Copyright (C) 2006 Jamie Bullock
+ *
+ * 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., 51 Franklin Street, Fifth Floor, Boston, MA  02110-1301, 
+ * USA.
+ */
+
+/** \file xtract_helper.h: helper functions for making life with libxtract a bit more bearable */
+
+#ifndef XTRACT_HELPER_H
+#define XTRACT_HELPER_H
+
+#ifdef __cplusplus
+extern "C" {
+#endif
+
+/**
+  * \defgroup helper helper functions
+  *
+  * Declares helper functions, and their parameters. 
+  *
+  * \note These functions don't necessarily conform to the prototype used in xtract_scalar.h and xtract_vector.h etc, and as such are intended to be called 'directly' rather than via the xtract[] function pointer array (libxtract.h)
+  *
+  * @{
+  */
+
+/** \brief Apply a window function to an array of length N
+ *
+ * \param *data a pointer to an array of floats
+ * \param N the number of elements in the array pointed to by *data
+ * \param *argv a pointer to a window function as returned by xtract_make_window()
+ * \param *result a pointer to the first element an array containing the windowed data
+ *
+ * It is up to the caller to generate and free the array containing the window, and to allocate and free memory of size N to hold the data pointed to by *result
+ *
+ */
+int xtract_windowed(const float *data, const int N, const void *argv, float *result);
+
+/** \brief Divides the array pointed to by *data into two subframes, and applies a given feature to each subframe, returning them in a single array pointed to by result
+ *
+ * \param *data an array of floats
+ * \param N the number of elements in the array pointed by *data
+ * \param feature an integer representing the feature to be applied to each subframe in data. This will be a value as given in the enumeration xtract_features_ (libxtract.h)
+ * \param *argv a pointer to the argument vector to be passed to the feature extraction function as determined by feature
+ * \param *result a pointer to the 'packed' results of the feature calculation. This may be passed in as *data to xtract_features_from_subframes() to calculate further features on the subframes, or xtract_difference_vector(), to get the difference between the subframes.
+ *
+ */
+int xtract_features_from_subframes(const float *data, const int N, const int feature, const void *argv, float *result);
+
+/** @} */
+
+#ifdef __cplusplus
+}
+#endif
+
+#endif
+
+
+