From c277634b13117e721e43f34a09cafb93c725fa3f Mon Sep 17 00:00:00 2001 From: Jamie Bullock Date: Wed, 9 Jan 2013 12:45:29 +0000 Subject: switched from single to double precision througout. closes #9 --- src/delta.c | 28 ++-- src/descriptors.c | 44 ++--- src/helper.c | 12 +- src/init.c | 26 +-- src/libxtract.c | 2 +- src/scalar.c | 387 +++++++++++++++++++++----------------------- src/vector.c | 162 ++++++++----------- src/window.c | 96 +++++------ src/xtract_window_private.h | 18 +-- 9 files changed, 361 insertions(+), 414 deletions(-) (limited to 'src') diff --git a/src/delta.c b/src/delta.c index eea91d2..8e6f611 100644 --- a/src/delta.c +++ b/src/delta.c @@ -27,7 +27,7 @@ #include "xtract/libxtract.h" -int xtract_flux(const float *data, const int N, const void *argv , float *result) +int xtract_flux(const double *data, const int N, const void *argv , double *result) { /* FIX: don't be lazy -- take the lnorm of the difference vector! */ @@ -35,20 +35,20 @@ int xtract_flux(const float *data, const int N, const void *argv , float *result } -int xtract_lnorm(const float *data, const int N, const void *argv , float *result) +int xtract_lnorm(const double *data, const int N, const void *argv , double *result) { int n, type; - float order; + double order; - order = *(float *)argv; - type = *((float *)argv+1); + order = *(double *)argv; + type = *((double *)argv+1); - order = order > 0 ? order : 2.f; + order = order > 0 ? order : 2.0; - *result = 0.f; + *result = 0.0; switch(type) { @@ -57,40 +57,40 @@ int xtract_lnorm(const float *data, const int N, const void *argv , float *resul for(n = 0; n < N; n++) { if(data[n] > 0) - *result += powf(data[n], order); + *result += pow(data[n], order); } break; default: for(n = 0; n < N; n++) - *result += powf(data[n], order); + *result += pow(data[n], order); break; } - *result = powf(*result, 1.f / order); + *result = pow(*result, 1.0 / order); return XTRACT_SUCCESS; } -int xtract_attack_time(const float *data, const int N, const void *argv , float *result) +int xtract_attack_time(const double *data, const int N, const void *argv , double *result) { return XTRACT_FEATURE_NOT_IMPLEMENTED; } -int xtract_decay_time(const float *data, const int N, const void *argv, float *result) +int xtract_decay_time(const double *data, const int N, const void *argv, double *result) { return XTRACT_FEATURE_NOT_IMPLEMENTED; } -int xtract_difference_vector(const float *data, const int N, const void *argv, float *result) +int xtract_difference_vector(const double *data, const int N, const void *argv, double *result) { - const float *frame1, + const double *frame1, *frame2; int n; diff --git a/src/descriptors.c b/src/descriptors.c index 6fe7942..15b790b 100644 --- a/src/descriptors.c +++ b/src/descriptors.c @@ -32,7 +32,7 @@ xtract_function_descriptor_t *xtract_make_descriptors(void) int f , F; char *name, *p_name, *desc, *p_desc, *author; - float *argv_min, *argv_max, *argv_def, *result_min, *result_max; + double *argv_min, *argv_max, *argv_def, *result_min, *result_max; int *argc, *year, *argv_donor; xtract_vector_t *data_format, *result_format; xtract_unit_t *data_unit, *argv_unit, *result_unit; @@ -80,8 +80,8 @@ xtract_function_descriptor_t *xtract_make_descriptors(void) *argv_unit = XTRACT_DBFS; break; case XTRACT_SPECTRAL_INHARMONICITY: - *argv_min = 0.f; - *argv_max = XTRACT_SR_UPPER_LIMIT / 2; + *argv_min = 0.0; + *argv_max = XTRACT_SR_UPPER_LIMIT / 2.0; *argv_def = XTRACT_FUNDAMENTAL_DEFAULT; *argv_unit = XTRACT_HERTZ; break; @@ -104,29 +104,29 @@ xtract_function_descriptor_t *xtract_make_descriptors(void) *argv_max = XTRACT_FFT_BANDS_MAX; *argv_def = XTRACT_SPEC_BW_DEF ; *argv_unit = XTRACT_HERTZ; - *(argv_min + 1) = 0.f; - *(argv_max + 1) = 100.f; - *(argv_def + 1) = 95.f; + *(argv_min + 1) = 0.0; + *(argv_max + 1) = 100.0; + *(argv_def + 1) = 95.0; *(argv_unit + 1) = XTRACT_PERCENT; break; case XTRACT_PEAK_SPECTRUM: - *argv_min = XTRACT_SR_LOWER_LIMIT / 2; - *argv_max = XTRACT_SR_UPPER_LIMIT / 2; - *argv_def = XTRACT_SR_DEFAULT / 2; + *argv_min = XTRACT_SR_LOWER_LIMIT / 2.0; + *argv_max = XTRACT_SR_UPPER_LIMIT / 2.0; + *argv_def = XTRACT_SR_DEFAULT / 2.0; *argv_unit = XTRACT_HERTZ; - *(argv_min + 1) = 0.f; - *(argv_max + 1) = 100.f ; - *(argv_def + 1) = 10.f ; + *(argv_min + 1) = 0.0; + *(argv_max + 1) = 100.0 ; + *(argv_def + 1) = 10.0 ; *(argv_unit + 1) = XTRACT_PERCENT; break; case XTRACT_HARMONIC_SPECTRUM: - *argv_min = 0.f; - *argv_max = XTRACT_SR_UPPER_LIMIT / 2; + *argv_min = 0.0; + *argv_max = XTRACT_SR_UPPER_LIMIT / 2.0; *argv_def = XTRACT_FUNDAMENTAL_DEFAULT; *argv_unit = XTRACT_HERTZ; - *(argv_min + 1) = 0.f; - *(argv_max + 1) = 1.f ; - *(argv_def + 1) = .1f ; + *(argv_min + 1) = 0.0; + *(argv_max + 1) = 1.0 ; + *(argv_def + 1) = .1 ; *(argv_unit + 1) = XTRACT_NONE; break; case XTRACT_NOISINESS: @@ -1251,18 +1251,18 @@ xtract_function_descriptor_t *xtract_make_descriptors(void) case XTRACT_HPS: case XTRACT_ROLLOFF: *result_unit = XTRACT_HERTZ; - *result_min = 0.f; - *result_max = XTRACT_SR_UPPER_LIMIT / 2; + *result_min = 0.0; + *result_max = XTRACT_SR_UPPER_LIMIT / 2.0; break; case XTRACT_ZCR: *result_unit = XTRACT_HERTZ; - *result_min = 0.f; + *result_min = 0.0; *result_max = XTRACT_ANY; break; case XTRACT_ODD_EVEN_RATIO: *result_unit = XTRACT_NONE; - *result_min = 0.f; - *result_max = 1.f; + *result_min = 0.0; + *result_max = 1.0; break; case XTRACT_FLATNESS_DB: *result_unit = XTRACT_DBFS; diff --git a/src/helper.c b/src/helper.c index 2a050e9..de33ffa 100644 --- a/src/helper.c +++ b/src/helper.c @@ -35,14 +35,14 @@ #define INDEX 1 #endif -int xtract_windowed(const float *data, const int N, const void *argv, float *result) +int xtract_windowed(const double *data, const int N, const void *argv, double *result) { int n; - const float *window; + const double *window; n = N; - window = (const float *)argv; + window = (const double *)argv; while(n--) result[n] = data[n] * window[n]; @@ -51,12 +51,12 @@ int xtract_windowed(const float *data, const int N, const void *argv, float *res } -int xtract_features_from_subframes(const float *data, const int N, const int feature, const void *argv, float *result) +int xtract_features_from_subframes(const double *data, const int N, const int feature, const void *argv, double *result) { - const float *frame1, + const double *frame1, *frame2; - float *result1, + double *result1, *result2; int n, rv; diff --git a/src/init.c b/src/init.c index 7b932da..87f8caf 100644 --- a/src/init.c +++ b/src/init.c @@ -39,11 +39,11 @@ #include "xtract_globals_private.h" -int xtract_init_mfcc(int N, float nyquist, int style, float freq_min, float freq_max, int freq_bands, float **fft_tables) +int xtract_init_mfcc(int N, double nyquist, int style, double freq_min, double freq_max, int freq_bands, double **fft_tables) { int n, i, k, *fft_peak, M, next_peak; - float norm, mel_freq_max, mel_freq_min, norm_fact, height, inc, val, + double 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; @@ -54,11 +54,11 @@ int xtract_init_mfcc(int N, float nyquist, int style, float freq_min, float freq 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)); + mel_peak = (double *)malloc((freq_bands + 2) * sizeof(double)); /* +2 for zeros at start and end */ - lin_peak = (float *)malloc((freq_bands + 2) * sizeof(float)); + lin_peak = (double *)malloc((freq_bands + 2) * sizeof(double)); fft_peak = (int *)malloc((freq_bands + 2) * sizeof(int)); - height_norm = (float *)malloc(freq_bands * sizeof(float)); + height_norm = (double *)malloc(freq_bands * sizeof(double)); if(mel_peak == NULL || height_norm == NULL || lin_peak == NULL || fft_peak == NULL) @@ -109,7 +109,7 @@ int xtract_init_mfcc(int N, float nyquist, int style, float freq_min, float freq // zero the start of the array for(k = 0; k < i; k++) - fft_tables[n][k] = 0.f; + fft_tables[n][k] = 0.0; // fill in the rise for(; i <= fft_peak[n]; i++) @@ -133,7 +133,7 @@ int xtract_init_mfcc(int N, float nyquist, int style, float freq_min, float freq // zero the rest of the array for(k = next_peak + 1; k < N; k++) - fft_tables[n][k] = 0.f; + fft_tables[n][k] = 0.0; } @@ -239,10 +239,10 @@ void xtract_free_fft(void) } } -int xtract_init_bark(int N, float sr, int *band_limits) +int xtract_init_bark(int N, double 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)*/ + double 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; @@ -253,11 +253,11 @@ int xtract_init_bark(int N, float sr, int *band_limits) return XTRACT_SUCCESS; } -float *xtract_init_window(const int N, const int type) +double *xtract_init_window(const int N, const int type) { - float *window; + double *window; - window = malloc(N * sizeof(float)); + window = malloc(N * sizeof(double)); switch (type) { @@ -296,7 +296,7 @@ float *xtract_init_window(const int N, const int type) return window; } -void xtract_free_window(float *window) +void xtract_free_window(double *window) { free(window); } diff --git a/src/libxtract.c b/src/libxtract.c index 43a7d16..3bc744e 100644 --- a/src/libxtract.c +++ b/src/libxtract.c @@ -24,7 +24,7 @@ #include "xtract/libxtract.h" #define XTRACT_H -int(*xtract[])(const float *, const int, const void *, float *) = +int(*xtract[])(const double *, const int, const void *, double *) = { /* xtract_scalar.h */ xtract_mean, diff --git a/src/scalar.c b/src/scalar.c index 1ef3839..80f0058 100644 --- a/src/scalar.c +++ b/src/scalar.c @@ -34,37 +34,12 @@ #include "xtract/xtract_helper.h" #include "xtract_macros_private.h" -#ifndef powf -#define powf pow -#endif - -#ifndef expf -#define expf exp -#endif - -#ifndef sqrtf -#define sqrtf sqrt -#endif - -#ifndef fabsf -#define fabsf fabs -#endif - - -void test(void) -{ - printf("Hello world\n"); -#ifdef WORDS_BIGENDIAN - printf("Big endian!\n"); -#endif -} - -int xtract_mean(const float *data, const int N, const void *argv, float *result) +int xtract_mean(const double *data, const int N, const void *argv, double *result) { int n = N; - *result = 0.f; + *result = 0.0; while(n--) *result += data[n]; @@ -74,57 +49,57 @@ int xtract_mean(const float *data, const int N, const void *argv, float *result) return XTRACT_SUCCESS; } -int xtract_variance(const float *data, const int N, const void *argv, float *result) +int xtract_variance(const double *data, const int N, const void *argv, double *result) { int n = N; - *result = 0.f; + *result = 0.0; while(n--) - *result += powf(data[n] - *(float *)argv, 2); + *result += pow(data[n] - *(double *)argv, 2); *result = *result / (N - 1); return XTRACT_SUCCESS; } -int xtract_standard_deviation(const float *data, const int N, const void *argv, float *result) +int xtract_standard_deviation(const double *data, const int N, const void *argv, double *result) { - *result = sqrtf(*(float *)argv); + *result = sqrt(*(double *)argv); return XTRACT_SUCCESS; } -int xtract_average_deviation(const float *data, const int N, const void *argv, float *result) +int xtract_average_deviation(const double *data, const int N, const void *argv, double *result) { int n = N; - *result = 0.f; + *result = 0.0; while(n--) - *result += fabsf(data[n] - *(float *)argv); + *result += fabs(data[n] - *(double *)argv); *result /= N; return XTRACT_SUCCESS; } -int xtract_skewness(const float *data, const int N, const void *argv, float *result) +int xtract_skewness(const double *data, const int N, const void *argv, double *result) { int n = N; - float temp = 0.f; + double temp = 0.0; - *result = 0.f; + *result = 0.0; while(n--) { - temp = (data[n] - ((float *)argv)[0]) / ((float *)argv)[1]; - *result += powf(temp, 3); + temp = (data[n] - ((double *)argv)[0]) / ((double *)argv)[1]; + *result += pow(temp, 3); } *result /= N; @@ -133,34 +108,34 @@ int xtract_skewness(const float *data, const int N, const void *argv, float *re return XTRACT_SUCCESS; } -int xtract_kurtosis(const float *data, const int N, const void *argv, float *result) +int xtract_kurtosis(const double *data, const int N, const void *argv, double *result) { int n = N; - float temp = 0.f; + double temp = 0.0; - *result = 0.f; + *result = 0.0; while(n--) { - temp = (data[n] - ((float *)argv)[0]) / ((float *)argv)[1]; - *result += powf(temp, 4); + temp = (data[n] - ((double *)argv)[0]) / ((double *)argv)[1]; + *result += pow(temp, 4); } *result /= N; - *result -= 3.0f; + *result -= 3.0; return XTRACT_SUCCESS; } -int xtract_spectral_centroid(const float *data, const int N, const void *argv, float *result) +int xtract_spectral_centroid(const double *data, const int N, const void *argv, double *result) { int n = (N >> 1); - const float *freqs, *amps; - float FA = 0.f, A = 0.f; + const double *freqs, *amps; + double FA = 0.0, A = 0.0; amps = data; freqs = data + n; @@ -171,39 +146,39 @@ int xtract_spectral_centroid(const float *data, const int N, const void *argv, A += amps[n]; } - if(A == 0.f) - *result = 0.f; + if(A == 0.0) + *result = 0.0; else *result = FA / A; return XTRACT_SUCCESS; } -int xtract_spectral_mean(const float *data, const int N, const void *argv, float *result) +int xtract_spectral_mean(const double *data, const int N, const void *argv, double *result) { return xtract_spectral_centroid(data, N, argv, result); } -int xtract_spectral_variance(const float *data, const int N, const void *argv, float *result) +int xtract_spectral_variance(const double *data, const int N, const void *argv, double *result) { int m; - float A = 0.f; - const float *freqs, *amps; + double A = 0.0; + const double *freqs, *amps; m = N >> 1; amps = data; freqs = data + m; - *result = 0.f; + *result = 0.0; while(m--) { A += amps[m]; - *result += powf(freqs[m] - ((float *)argv)[0], 2) * amps[m]; + *result += pow(freqs[m] - ((double *)argv)[0], 2) * amps[m]; } *result = *result / A; @@ -211,30 +186,30 @@ int xtract_spectral_variance(const float *data, const int N, const void *argv, f return XTRACT_SUCCESS; } -int xtract_spectral_standard_deviation(const float *data, const int N, const void *argv, float *result) +int xtract_spectral_standard_deviation(const double *data, const int N, const void *argv, double *result) { - *result = sqrtf(*(float *)argv); + *result = sqrt(*(double *)argv); return XTRACT_SUCCESS; } -/*int xtract_spectral_average_deviation(const float *data, const int N, const void *argv, float *result){ +/*int xtract_spectral_average_deviation(const double *data, const int N, const void *argv, double *result){ int m; - float A = 0.f; - const float *freqs, *amps; + double A = 0.0; + const double *freqs, *amps; m = N >> 1; amps = data; freqs = data + m; - *result = 0.f; + *result = 0.0; while(m--){ A += amps[m]; - *result += fabsf((amps[m] * freqs[m]) - *(float *)argv); + *result += fabs((amps[m] * freqs[m]) - *(double *)argv); } *result /= A; @@ -242,89 +217,89 @@ int xtract_spectral_standard_deviation(const float *data, const int N, const voi return XTRACT_SUCCESS; }*/ -int xtract_spectral_skewness(const float *data, const int N, const void *argv, float *result) +int xtract_spectral_skewness(const double *data, const int N, const void *argv, double *result) { int m; - const float *freqs, *amps; + const double *freqs, *amps; m = N >> 1; amps = data; freqs = data + m; - *result = 0.f; + *result = 0.0; while(m--) - *result += powf(freqs[m] - ((float *)argv)[0], 3) * amps[m]; + *result += pow(freqs[m] - ((double *)argv)[0], 3) * amps[m]; - *result /= powf(((float *)argv)[1], 3); + *result /= pow(((double *)argv)[1], 3); return XTRACT_SUCCESS; } -int xtract_spectral_kurtosis(const float *data, const int N, const void *argv, float *result) +int xtract_spectral_kurtosis(const double *data, const int N, const void *argv, double *result) { int m; - const float *freqs, *amps; + const double *freqs, *amps; m = N >> 1; amps = data; freqs = data + m; - *result = 0.f; + *result = 0.0; while(m--) - *result += powf(freqs[m] - ((float *)argv)[0], 4) * amps[m]; + *result += pow(freqs[m] - ((double *)argv)[0], 4) * amps[m]; - *result /= powf(((float *)argv)[1], 4); - *result -= 3.0f; + *result /= pow(((double *)argv)[1], 4); + *result -= 3.0; return XTRACT_SUCCESS; } -int xtract_irregularity_k(const float *data, const int N, const void *argv, float *result) +int xtract_irregularity_k(const double *data, const int N, const void *argv, double *result) { int n, M = N - 1; - *result = 0.f; + *result = 0.0; for(n = 1; n < M; n++) - *result += fabsf(data[n] - (data[n-1] + data[n] + data[n+1]) / 3.f); + *result += fabs(data[n] - (data[n-1] + data[n] + data[n+1]) / 3.0); return XTRACT_SUCCESS; } -int xtract_irregularity_j(const float *data, const int N, const void *argv, float *result) +int xtract_irregularity_j(const double *data, const int N, const void *argv, double *result) { int n = N; - double num = 0.f, den = 0.f; + double num = 0.0, den = 0.0; while(n--) { - num += powf(data[n] - data[n+1], 2); - den += powf(data[n], 2); + num += pow(data[n] - data[n+1], 2); + den += pow(data[n], 2); } - *result = (float)(num / den); + *result = (double)(num / den); return XTRACT_SUCCESS; } -int xtract_tristimulus_1(const float *data, const int N, const void *argv, float *result) +int xtract_tristimulus_1(const double *data, const int N, const void *argv, double *result) { int n = N; - float den, p1, temp; + double den, p1, temp; - den = p1 = temp = 0.f; + den = p1 = temp = 0.0; for(n = 0; n < N; n++) { @@ -336,9 +311,9 @@ int xtract_tristimulus_1(const float *data, const int N, const void *argv, float } } - if(den == 0.f || p1 == 0.f) + if(den == 0.0 || p1 == 0.0) { - *result = 0.f; + *result = 0.0; return XTRACT_NO_RESULT; } else @@ -348,14 +323,14 @@ int xtract_tristimulus_1(const float *data, const int N, const void *argv, float } } -int xtract_tristimulus_2(const float *data, const int N, const void *argv, float *result) +int xtract_tristimulus_2(const double *data, const int N, const void *argv, double *result) { int n = N; - float den, p2, p3, p4, ps, temp; + double den, p2, p3, p4, ps, temp; - den = p2 = p3 = p4 = ps = temp = 0.f; + den = p2 = p3 = p4 = ps = temp = 0.0; for(n = 0; n < N; n++) { @@ -373,9 +348,9 @@ int xtract_tristimulus_2(const float *data, const int N, const void *argv, float ps = p2 + p3 + p4; - if(den == 0.f || ps == 0.f) + if(den == 0.0 || ps == 0.0) { - *result = 0.f; + *result = 0.0; return XTRACT_NO_RESULT; } else @@ -386,14 +361,14 @@ int xtract_tristimulus_2(const float *data, const int N, const void *argv, float } -int xtract_tristimulus_3(const float *data, const int N, const void *argv, float *result) +int xtract_tristimulus_3(const double *data, const int N, const void *argv, double *result) { int n = N, count = 0; - float den, num, temp; + double den, num, temp; - den = num = temp = 0.f; + den = num = temp = 0.0; for(n = 0; n < N; n++) { @@ -406,9 +381,9 @@ int xtract_tristimulus_3(const float *data, const int N, const void *argv, float } } - if(den == 0.f || num == 0.f) + if(den == 0.0 || num == 0.0) { - *result = 0.f; + *result = 0.0; return XTRACT_NO_RESULT; } else @@ -418,15 +393,15 @@ int xtract_tristimulus_3(const float *data, const int N, const void *argv, float } } -int xtract_smoothness(const float *data, const int N, const void *argv, float *result) +int xtract_smoothness(const double *data, const int N, const void *argv, double *result) { int n, M; - float *input; + double *input; - input = (float *)malloc(N * sizeof(float)); - memcpy(input, data, N * sizeof(float)); + input = (double *)malloc(N * sizeof(double)); + memcpy(input, data, N * sizeof(double)); if (input[0] <= 0) input[0] = XTRACT_LOG_LIMIT; @@ -439,8 +414,8 @@ int xtract_smoothness(const float *data, const int N, const void *argv, float *r { if(input[n+1] <= 0) input[n+1] = XTRACT_LOG_LIMIT; - *result += fabsf(20.f * logf(input[n]) - (20.f * logf(input[n-1]) + - 20.f * logf(input[n]) + 20.f * logf(input[n+1])) / 3.f); + *result += fabs(20.0 * log(input[n]) - (20.0 * log(input[n-1]) + + 20.0 * log(input[n]) + 20.0 * log(input[n+1])) / 3.0); } free(input); @@ -448,13 +423,13 @@ int xtract_smoothness(const float *data, const int N, const void *argv, float *r return XTRACT_SUCCESS; } -int xtract_spread(const float *data, const int N, const void *argv, float *result) +int xtract_spread(const double *data, const int N, const void *argv, double *result) { return xtract_spectral_variance(data, N, argv, result); } -int xtract_zcr(const float *data, const int N, const void *argv, float *result) +int xtract_zcr(const double *data, const int N, const void *argv, double *result) { int n = N; @@ -462,39 +437,39 @@ int xtract_zcr(const float *data, const int N, const void *argv, float *result) for(n = 1; n < N; n++) if(data[n] * data[n-1] < 0) (*result)++; - *result /= (float)N; + *result /= (double)N; return XTRACT_SUCCESS; } -int xtract_rolloff(const float *data, const int N, const void *argv, float *result) +int xtract_rolloff(const double *data, const int N, const void *argv, double *result) { int n = N; - float pivot, temp, percentile; + double pivot, temp, percentile; - pivot = temp = 0.f; - percentile = ((float *)argv)[1]; + pivot = temp = 0.0; + percentile = ((double *)argv)[1]; while(n--) pivot += data[n]; - pivot *= percentile / 100.f; + pivot *= percentile / 100.0; for(n = 0; temp < pivot; n++) temp += data[n]; - *result = n * ((float *)argv)[0]; - /* *result = (n / (float)N) * (((float *)argv)[1] * .5); */ + *result = n * ((double *)argv)[0]; + /* *result = (n / (double)N) * (((double *)argv)[1] * .5); */ return XTRACT_SUCCESS; } -int xtract_loudness(const float *data, const int N, const void *argv, float *result) +int xtract_loudness(const double *data, const int N, const void *argv, double *result) { int n = N, rv; - *result = 0.f; + *result = 0.0; if(n > XTRACT_BARK_BANDS) { @@ -505,27 +480,27 @@ int xtract_loudness(const float *data, const int N, const void *argv, float *res rv = XTRACT_SUCCESS; while(n--) - *result += powf(data[n], 0.23); + *result += pow(data[n], 0.23); return rv; } -int xtract_flatness(const float *data, const int N, const void *argv, float *result) +int xtract_flatness(const double *data, const int N, const void *argv, double *result) { int n, count, denormal_found; double num, den, temp; - num = 1.f; - den = temp = 0.f; + num = 1.0; + den = temp = 0.0; denormal_found = 0; count = 0; for(n = 0; n < N; n++) { - if((temp = data[n]) != 0.f) + if((temp = data[n]) != 0.0) { if (xtract_is_denormal(num)) { @@ -540,15 +515,15 @@ int xtract_flatness(const float *data, const int N, const void *argv, float *res if(!count) { - *result = 0.f; + *result = 0.0; return XTRACT_NO_RESULT; } - num = powf(num, 1.f / (float)N); - den /= (float)N; + num = pow(num, 1.0 / (double)N); + den /= (double)N; - *result = (float) (num / den); + *result = (double) (num / den); if(denormal_found) return XTRACT_DENORMAL_FOUND; @@ -557,12 +532,12 @@ int xtract_flatness(const float *data, const int N, const void *argv, float *res } -int xtract_flatness_db(const float *data, const int N, const void *argv, float *result) +int xtract_flatness_db(const double *data, const int N, const void *argv, double *result) { - float flatness; + double flatness; - flatness = *(float *)argv; + flatness = *(double *)argv; if (flatness <= 0) flatness = XTRACT_LOG_LIMIT; @@ -573,27 +548,27 @@ int xtract_flatness_db(const float *data, const int N, const void *argv, float * } -int xtract_tonality(const float *data, const int N, const void *argv, float *result) +int xtract_tonality(const double *data, const int N, const void *argv, double *result) { - float sfmdb; + double sfmdb; - sfmdb = *(float *)argv; + sfmdb = *(double *)argv; - *result = XTRACT_MIN(sfmdb / -60.f, 1); + *result = XTRACT_MIN(sfmdb / -60.0, 1); return XTRACT_SUCCESS; } -int xtract_crest(const float *data, const int N, const void *argv, float *result) +int xtract_crest(const double *data, const int N, const void *argv, double *result) { - float max, mean; + double max, mean; - max = mean = 0.f; + max = mean = 0.0; - max = *(float *)argv; - mean = *((float *)argv+1); + max = *(double *)argv; + mean = *((double *)argv+1); *result = max / mean; @@ -601,15 +576,15 @@ int xtract_crest(const float *data, const int N, const void *argv, float *result } -int xtract_noisiness(const float *data, const int N, const void *argv, float *result) +int xtract_noisiness(const double *data, const int N, const void *argv, double *result) { - float h, i, p; /*harmonics, inharmonics, partials */ + double h, i, p; /*harmonics, inharmonics, partials */ - i = p = h = 0.f; + i = p = h = 0.0; - h = *(float *)argv; - p = *((float *)argv+1); + h = *(double *)argv; + p = *((double *)argv+1); i = p - h; @@ -619,28 +594,28 @@ int xtract_noisiness(const float *data, const int N, const void *argv, float *re } -int xtract_rms_amplitude(const float *data, const int N, const void *argv, float *result) +int xtract_rms_amplitude(const double *data, const int N, const void *argv, double *result) { int n = N; - *result = 0.f; + *result = 0.0; while(n--) *result += XTRACT_SQ(data[n]); - *result = sqrtf(*result / (float)N); + *result = sqrt(*result / (double)N); return XTRACT_SUCCESS; } -int xtract_spectral_inharmonicity(const float *data, const int N, const void *argv, float *result) +int xtract_spectral_inharmonicity(const double *data, const int N, const void *argv, double *result) { int n = N >> 1; - float num = 0.f, den = 0.f, fund; - const float *freqs, *amps; + double num = 0.0, den = 0.0, fund; + const double *freqs, *amps; - fund = *(float *)argv; + fund = *(double *)argv; amps = data; freqs = data + n; @@ -648,7 +623,7 @@ int xtract_spectral_inharmonicity(const float *data, const int N, const void *ar { if(amps[n]) { - num += fabsf(freqs[n] - n * fund) * XTRACT_SQ(amps[n]); + num += fabs(freqs[n] - n * fund) * XTRACT_SQ(amps[n]); den += XTRACT_SQ(amps[n]); } } @@ -659,19 +634,19 @@ int xtract_spectral_inharmonicity(const float *data, const int N, const void *ar } -int xtract_power(const float *data, const int N, const void *argv, float *result) +int xtract_power(const double *data, const int N, const void *argv, double *result) { return XTRACT_FEATURE_NOT_IMPLEMENTED; } -int xtract_odd_even_ratio(const float *data, const int N, const void *argv, float *result) +int xtract_odd_even_ratio(const double *data, const int N, const void *argv, double *result) { int M = (N >> 1), n; - float odd = 0.f, even = 0.f, temp; + double odd = 0.0, even = 0.0, temp; for(n = 0; n < M; n++) { @@ -688,9 +663,9 @@ int xtract_odd_even_ratio(const float *data, const int N, const void *argv, floa } } - if(odd == 0.f || even == 0.f) + if(odd == 0.0 || even == 0.0) { - *result = 0.f; + *result = 0.0; return XTRACT_NO_RESULT; } else @@ -700,15 +675,15 @@ int xtract_odd_even_ratio(const float *data, const int N, const void *argv, floa } } -int xtract_sharpness(const float *data, const int N, const void *argv, float *result) +int xtract_sharpness(const double *data, const int N, const void *argv, double *result) { int n = N, rv; - float sl, g; /* sl = specific loudness */ + double sl, g; /* sl = specific loudness */ double temp; - sl = g = 0.f; - temp = 0.f; + sl = g = 0.0; + temp = 0.0; if(n > XTRACT_BARK_BANDS) rv = XTRACT_BAD_VECTOR_SIZE; @@ -718,27 +693,27 @@ int xtract_sharpness(const float *data, const int N, const void *argv, float *re while(n--) { - sl = powf(data[n], 0.23); - g = (n < 15 ? 1.f : 0.066 * expf(0.171 * n)); + sl = pow(data[n], 0.23); + g = (n < 15 ? 1.0 : 0.066 * exp(0.171 * n)); temp += n * g * sl; } - temp = 0.11 * temp / (float)N; - *result = (float)temp; + temp = 0.11 * temp / (double)N; + *result = (double)temp; return rv; } -int xtract_spectral_slope(const float *data, const int N, const void *argv, float *result) +int xtract_spectral_slope(const double *data, const int N, const void *argv, double *result) { - const float *freqs, *amps; - float f, a, + const double *freqs, *amps; + double f, a, F, A, FA, FXTRACT_SQ; /* sums of freqs, amps, freq * amps, freq squared */ int n, M; - F = A = FA = FXTRACT_SQ = 0.f; + F = A = FA = FXTRACT_SQ = 0.0; n = M = N >> 1; amps = data; @@ -754,30 +729,30 @@ int xtract_spectral_slope(const float *data, const int N, const void *argv, floa FXTRACT_SQ += f * f; } - *result = (1.f / A) * (M * FA - F * A) / (M * FXTRACT_SQ - F * F); + *result = (1.0 / A) * (M * FA - F * A) / (M * FXTRACT_SQ - F * F); return XTRACT_SUCCESS; } -int xtract_lowest_value(const float *data, const int N, const void *argv, float *result) +int xtract_lowest_value(const double *data, const int N, const void *argv, double *result) { int n = N; - float temp; + double temp; *result = data[--n]; while(n--) { - if((temp = data[n]) > *(float *)argv) + if((temp = data[n]) > *(double *)argv) *result = XTRACT_MIN(*result, data[n]); } return XTRACT_SUCCESS; } -int xtract_highest_value(const float *data, const int N, const void *argv, float *result) +int xtract_highest_value(const double *data, const int N, const void *argv, double *result) { int n = N; @@ -791,12 +766,12 @@ int xtract_highest_value(const float *data, const int N, const void *argv, float } -int xtract_sum(const float *data, const int N, const void *argv, float *result) +int xtract_sum(const double *data, const int N, const void *argv, double *result) { int n = N; - *result = 0.f; + *result = 0.0; while(n--) *result += *data++; @@ -805,12 +780,12 @@ int xtract_sum(const float *data, const int N, const void *argv, float *result) } -int xtract_nonzero_count(const float *data, const int N, const void *argv, float *result) +int xtract_nonzero_count(const double *data, const int N, const void *argv, double *result) { int n = N; - *result = 0.f; + *result = 0.0; while(n--) *result += (*data++ ? 1 : 0); @@ -819,35 +794,35 @@ int xtract_nonzero_count(const float *data, const int N, const void *argv, float } -int xtract_hps(const float *data, const int N, const void *argv, float *result) +int xtract_hps(const double *data, const int N, const void *argv, double *result) { int n = N, M, m, l, peak_index, position1_lwr; - float *coeffs2, *coeffs3, *product, L, + double *coeffs2, *coeffs3, *product, L, largest1_lwr, peak, ratio1, sr; - sr = *(float*)argv; + sr = *(double*)argv; if(sr == 0) - sr = 44100.f; + sr = 44100.0; - coeffs2 = (float *)malloc(N * sizeof(float)); - coeffs3 = (float *)malloc(N * sizeof(float)); - product = (float *)malloc(N * sizeof(float)); + coeffs2 = (double *)malloc(N * sizeof(double)); + coeffs3 = (double *)malloc(N * sizeof(double)); + product = (double *)malloc(N * sizeof(double)); while(n--) coeffs2[n] = coeffs3[n] = 1; M = N >> 1; - L = N / 3.f; + L = N / 3.0; while(M--) { m = M << 1; - coeffs2[M] = (data[m] + data[m+1]) * 0.5f; + coeffs2[M] = (data[m] + data[m+1]) * 0.5; if(M < L) { l = M * 3; - coeffs3[M] = (data[l] + data[l+1] + data[l+2]) / 3.f; + coeffs3[M] = (data[l] + data[l+1] + data[l+2]) / 3.0; } } @@ -880,7 +855,7 @@ int xtract_hps(const float *data, const int N, const void *argv, float *result) peak_index * 0.6 && ratio1 > 0.1) peak_index = position1_lwr; - *result = sr / (float)peak_index; + *result = sr / (double)peak_index; free(coeffs2); free(coeffs3); @@ -890,25 +865,25 @@ int xtract_hps(const float *data, const int N, const void *argv, float *result) } -int xtract_f0(const float *data, const int N, const void *argv, float *result) +int xtract_f0(const double *data, const int N, const void *argv, double *result) { int M, tau, n; - float sr; + double sr; size_t bytes; - float f0, err_tau_1, err_tau_x, array_max, + double f0, err_tau_1, err_tau_x, array_max, threshold_peak, threshold_centre, *input; - sr = *(float *)argv; + sr = *(double *)argv; if(sr == 0) - sr = 44100.f; + sr = 44100.0; - input = (float *)malloc(bytes = N * sizeof(float)); + input = (double *)malloc(bytes = N * sizeof(double)); input = memcpy(input, data, bytes); - /* threshold_peak = *((float *)argv+1); - threshold_centre = *((float *)argv+2); - printf("peak: %.2f\tcentre: %.2f\n", threshold_peak, threshold_centre);*/ + /* threshold_peak = *((double *)argv+1); + threshold_centre = *((double *)argv+2); + printf("peak: %.2\tcentre: %.2\n", threshold_peak, threshold_centre);*/ /* add temporary dynamic control over thresholds to test clipping effects */ /* FIX: tweak and make into macros */ @@ -949,14 +924,14 @@ int xtract_f0(const float *data, const int N, const void *argv, float *result) /* Estimate fundamental freq */ for (n = 1; n < M; n++) - err_tau_1 = err_tau_1 + fabsf(input[n] - input[n+1]); + err_tau_1 = err_tau_1 + fabs(input[n] - input[n+1]); /* FIX: this doesn't pose too much load if it returns 'early', but if it can't find f0, load can be significant for larger block sizes M^2 iterations! */ for (tau = 2; tau < M; tau++) { err_tau_x = 0; for (n = 1; n < M; n++) { - err_tau_x = err_tau_x + fabsf(input[n] - input[n+tau]); + err_tau_x = err_tau_x + fabs(input[n] - input[n+tau]); } if (err_tau_x < err_tau_1) { @@ -971,27 +946,27 @@ int xtract_f0(const float *data, const int N, const void *argv, float *result) return XTRACT_NO_RESULT; } -int xtract_failsafe_f0(const float *data, const int N, const void *argv, float *result) +int xtract_failsafe_f0(const double *data, const int N, const void *argv, double *result) { - float *spectrum = NULL, argf[2], *peaks = NULL, return_code, sr; + double *spectrum = NULL, argf[2], *peaks = NULL, return_code, sr; return_code = xtract_f0(data, N, argv, result); if(return_code == XTRACT_NO_RESULT) { - sr = *(float *)argv; + sr = *(double *)argv; if(sr == 0) - sr = 44100.f; - spectrum = (float *)malloc(N * sizeof(float)); - peaks = (float *)malloc(N * sizeof(float)); + sr = 44100.0; + spectrum = (double *)malloc(N * sizeof(double)); + peaks = (double *)malloc(N * sizeof(double)); argf[0] = sr; argf[1] = XTRACT_MAGNITUDE_SPECTRUM; xtract_spectrum(data, N, argf, spectrum); - argf[1] = 10.f; + argf[1] = 10.0; xtract_peak_spectrum(spectrum, N >> 1, argf, peaks); - argf[0] = 0.f; + argf[0] = 0.0; xtract_lowest_value(peaks+(N >> 1), N >> 1, argf, result); free(spectrum); diff --git a/src/vector.c b/src/vector.c index f396745..5639374 100644 --- a/src/vector.c +++ b/src/vector.c @@ -31,47 +31,19 @@ #include "xtract/libxtract.h" #include "xtract_macros_private.h" - -#ifndef roundf -float roundf(float f) -{ - if (f - (int)f >= 0.5) - return (float)((int)f + 1); - else - return (float)((int)f); -} -#endif - -#ifndef powf -#define powf pow -#endif - -#ifndef expf -#define expf exp -#endif - -#ifndef sqrtf -#define sqrtf sqrt -#endif - -#ifndef fabsf -#define fabsf fabs -#endif - #include "xtract_globals_private.h" -#include "xtract_macros_private.h" -int xtract_spectrum(const float *data, const int N, const void *argv, float *result) +int xtract_spectrum(const double *data, const int N, const void *argv, double *result) { int vector = 0; int withDC = 0; int normalise = 0; - float q = 0.f; - float temp = 0.f; - float max = 0.f; - float NxN = XTRACT_SQ(N); - float *marker = NULL; + double q = 0.0; + double temp = 0.0; + double max = 0.0; + double NxN = XTRACT_SQ(N); + double *marker = NULL; size_t bytes = N * sizeof(double); double *rfft = NULL; unsigned int n = 0; @@ -86,10 +58,10 @@ int xtract_spectrum(const float *data, const int N, const void *argv, float *res rfft[n] = (double)data[n]; } - q = *(float *)argv; - vector = (int)*((float *)argv+1); - withDC = (int)*((float *)argv+2); - normalise = (int)*((float *)argv+3); + q = *(double *)argv; + vector = (int)*((double *)argv+1); + withDC = (int)*((double *)argv+2); + normalise = (int)*((double *)argv+3); XTRACT_CHECK_q; @@ -122,7 +94,7 @@ int xtract_spectrum(const float *data, const int N, const void *argv, float *res temp = XTRACT_SQ(rfft[nx2]) + XTRACT_SQ(rfft[nx2+1]); if (temp > XTRACT_LOG_LIMIT) { - temp = logf(sqrtf(temp) / (float)N); + temp = log(sqrt(temp) / (double)N); } else { @@ -160,7 +132,7 @@ int xtract_spectrum(const float *data, const int N, const void *argv, float *res } if ((temp = XTRACT_SQ(rfft[n]) + XTRACT_SQ(rfft[N - n])) > XTRACT_LOG_LIMIT) - temp = logf(temp / NxN); + temp = log(temp / NxN); else temp = XTRACT_LOG_LIMIT_DB; @@ -187,7 +159,7 @@ int xtract_spectrum(const float *data, const int N, const void *argv, float *res ++n; } - *marker = (float)(sqrt(XTRACT_SQ(rfft[n*2]) + + *marker = (double)(sqrt(XTRACT_SQ(rfft[n*2]) + XTRACT_SQ(rfft[n*2+1])) / (double)N); XTRACT_SET_FREQUENCY; @@ -208,7 +180,7 @@ int xtract_spectrum(const float *data, const int N, const void *argv, float *res return XTRACT_SUCCESS; } -int xtract_autocorrelation_fft(const float *data, const int N, const void *argv, float *result) +int xtract_autocorrelation_fft(const double *data, const int N, const void *argv, double *result) { double *rfft = NULL; @@ -219,7 +191,7 @@ int xtract_autocorrelation_fft(const float *data, const int N, const void *argv, /* Zero pad the input vector */ rfft = (double *)calloc(M, sizeof(double)); - memcpy(rfft, data, N * sizeof(float)); + memcpy(rfft, data, N * sizeof(double)); rdft(M, 1, rfft, ooura_data_autocorrelation_fft.ooura_ip, ooura_data_autocorrelation_fft.ooura_w); @@ -227,7 +199,7 @@ int xtract_autocorrelation_fft(const float *data, const int N, const void *argv, for(n = 2; n < M; n++) { rfft[n*2] = XTRACT_SQ(rfft[n*2]) + XTRACT_SQ(rfft[n*2+1]); - rfft[n*2+1] = 0.f; + rfft[n*2+1] = 0.0; } rfft[0] = XTRACT_SQ(rfft[0]); @@ -240,14 +212,14 @@ int xtract_autocorrelation_fft(const float *data, const int N, const void *argv, M = M * N; for(n = 0; n < N; n++) - result[n] = rfft[n] / (float)M; + result[n] = rfft[n] / (double)M; free(rfft); return XTRACT_SUCCESS; } -int xtract_mfcc(const float *data, const int N, const void *argv, float *result) +int xtract_mfcc(const double *data, const int N, const void *argv, double *result) { xtract_mel_filter *f; @@ -257,12 +229,12 @@ int xtract_mfcc(const float *data, const int N, const void *argv, float *result) for(filter = 0; filter < f->n_filters; filter++) { - result[filter] = 0.f; + result[filter] = 0.0; for(n = 0; n < N; n++) { result[filter] += data[n] * f->filters[filter][n]; } - result[filter] = logf(result[filter] < XTRACT_LOG_LIMIT ? XTRACT_LOG_LIMIT : result[filter]); + result[filter] = log(result[filter] < XTRACT_LOG_LIMIT ? XTRACT_LOG_LIMIT : result[filter]); } xtract_dct(result, f->n_filters, NULL, result); @@ -270,31 +242,31 @@ int xtract_mfcc(const float *data, const int N, const void *argv, float *result) return XTRACT_SUCCESS; } -int xtract_dct(const float *data, const int N, const void *argv, float *result) +int xtract_dct(const double *data, const int N, const void *argv, double *result) { int n; int m; - float *temp = calloc(N, sizeof(float)); + double *temp = calloc(N, sizeof(double)); for (n = 0; n < N; ++n) { for(m = 1; m <= N; ++m) { - temp[n] += data[m - 1] * cos(M_PI * (n / (float)N) * (m - 0.5)); + temp[n] += data[m - 1] * cos(M_PI * (n / (double)N) * (m - 0.5)); } } return XTRACT_SUCCESS; } -int xtract_autocorrelation(const float *data, const int N, const void *argv, float *result) +int xtract_autocorrelation(const double *data, const int N, const void *argv, double *result) { /* Naive time domain implementation */ int n = N, i; - float corr; + double corr; while(n--) { @@ -309,50 +281,50 @@ int xtract_autocorrelation(const float *data, const int N, const void *argv, flo return XTRACT_SUCCESS; } -int xtract_amdf(const float *data, const int N, const void *argv, float *result) +int xtract_amdf(const double *data, const int N, const void *argv, double *result) { int n = N, i; - float md, temp; + double md, temp; while(n--) { - md = 0.f; + md = 0.0; for(i = 0; i < N - n; i++) { temp = data[i] - data[i + n]; temp = (temp < 0 ? -temp : temp); md += temp; } - result[n] = md / (float)N; + result[n] = md / (double)N; } return XTRACT_SUCCESS; } -int xtract_asdf(const float *data, const int N, const void *argv, float *result) +int xtract_asdf(const double *data, const int N, const void *argv, double *result) { int n = N, i; - float sd; + double sd; while(n--) { - sd = 0.f; + sd = 0.0; for(i = 0; i < N - n; i++) { /*sd = 1;*/ sd += XTRACT_SQ(data[i] - data[i + n]); } - result[n] = sd / (float)N; + result[n] = sd / (double)N; } return XTRACT_SUCCESS; } -int xtract_bark_coefficients(const float *data, const int N, const void *argv, float *result) +int xtract_bark_coefficients(const double *data, const int N, const void *argv, double *result) { int *limits, band, n; @@ -361,7 +333,7 @@ int xtract_bark_coefficients(const float *data, const int N, const void *argv, f for(band = 0; band < XTRACT_BARK_BANDS - 1; band++) { - result[band] = 0.f; + result[band] = 0.0; for(n = limits[band]; n < limits[band + 1]; n++) result[band] += data[n]; } @@ -369,19 +341,19 @@ int xtract_bark_coefficients(const float *data, const int N, const void *argv, f return XTRACT_SUCCESS; } -int xtract_peak_spectrum(const float *data, const int N, const void *argv, float *result) +int xtract_peak_spectrum(const double *data, const int N, const void *argv, double *result) { - float threshold, max, y, y2, y3, p, q, *input = NULL; + double threshold, max, y, y2, y3, p, q, *input = NULL; size_t bytes; int n = N, rv = XTRACT_SUCCESS; - threshold = max = y = y2 = y3 = p = q = 0.f; + threshold = max = y = y2 = y3 = p = q = 0.0; if(argv != NULL) { - q = ((float *)argv)[0]; - threshold = ((float *)argv)[1]; + q = ((double *)argv)[0]; + threshold = ((double *)argv)[1]; } else rv = XTRACT_BAD_ARGV; @@ -394,9 +366,9 @@ int xtract_peak_spectrum(const float *data, const int N, const void *argv, float XTRACT_CHECK_q; - input = (float *)calloc(N, sizeof(float)); + input = (double *)calloc(N, sizeof(double)); - bytes = N * sizeof(float); + bytes = N * sizeof(double); if(input != NULL) input = memcpy(input, data, bytes); @@ -439,30 +411,30 @@ int xtract_peak_spectrum(const float *data, const int N, const void *argv, float return (rv ? rv : XTRACT_SUCCESS); } -int xtract_harmonic_spectrum(const float *data, const int N, const void *argv, float *result) +int xtract_harmonic_spectrum(const double *data, const int N, const void *argv, double *result) { int n = (N >> 1), M = n; - const float *freqs, *amps; - float f0, threshold, ratio, nearest, distance; + const double *freqs, *amps; + double f0, threshold, ratio, nearest, distance; amps = data; freqs = data + n; - f0 = *((float *)argv); - threshold = *((float *)argv+1); + f0 = *((double *)argv); + threshold = *((double *)argv+1); - ratio = nearest = distance = 0.f; + ratio = nearest = distance = 0.0; while(n--) { if(freqs[n]) { ratio = freqs[n] / f0; - nearest = roundf(ratio); + nearest = round(ratio); distance = fabs(nearest - ratio); if(distance > threshold) - result[n] = result[M + n] = 0.f; + result[n] = result[M + n] = 0.0; else { result[n] = amps[n]; @@ -470,19 +442,19 @@ int xtract_harmonic_spectrum(const float *data, const int N, const void *argv, f } } else - result[n] = result[M + n] = 0.f; + result[n] = result[M + n] = 0.0; } return XTRACT_SUCCESS; } -int xtract_lpc(const float *data, const int N, const void *argv, float *result) +int xtract_lpc(const double *data, const int N, const void *argv, double *result) { int i, j, k, M, L; - float r = 0.f, - error = 0.f; + double r = 0.0, + error = 0.0; - float *ref = NULL, + double *ref = NULL, *lpc = NULL ; error = data[0]; @@ -494,11 +466,11 @@ int xtract_lpc(const float *data, const int N, const void *argv, float *result) if(error == 0.0) { - memset(result, 0, M * sizeof(float)); + memset(result, 0, M * sizeof(double)); return XTRACT_NO_RESULT; } - memset(result, 0, M * sizeof(float)); + memset(result, 0, M * sizeof(double)); for (i = 0; i < L; i++) { @@ -513,7 +485,7 @@ int xtract_lpc(const float *data, const int N, const void *argv, float *result) lpc[i] = r; for (j = 0; j < i / 2; j++) { - float tmp = lpc[j]; + double tmp = lpc[j]; lpc[j] = r * lpc[i - 1 - j]; lpc[i - 1 - j] += r * tmp; } @@ -525,14 +497,14 @@ int xtract_lpc(const float *data, const int N, const void *argv, float *result) return XTRACT_SUCCESS; } -int xtract_lpcc(const float *data, const int N, const void *argv, float *result) +int xtract_lpcc(const double *data, const int N, const void *argv, double *result) { /* Given N lpc coefficients extract an LPC cepstrum of size argv[0] */ /* Based on an an algorithm by rabiner and Juang */ int n, k; - float sum; + double sum; int order = N - 1; /* Eventually change this to Q = 3/2 p as suggested in Rabiner */ int cep_length; @@ -540,13 +512,13 @@ int xtract_lpcc(const float *data, const int N, const void *argv, float *result) cep_length = N - 1; /* FIX: if we're going to have default values, they should come from the descriptor */ else cep_length = *(int *)argv; - //cep_length = (int)((float *)argv)[0]; + //cep_length = (int)((double *)argv)[0]; - memset(result, 0, cep_length * sizeof(float)); + memset(result, 0, cep_length * sizeof(double)); for (n = 1; n <= order && n <= cep_length; n++) { - sum = 0.f; + sum = 0.0; for (k = 1; k < n; k++) sum += k * result[k-1] * data[n - k]; result[n-1] = data[n] + sum / n; @@ -555,7 +527,7 @@ int xtract_lpcc(const float *data, const int N, const void *argv, float *result) /* be wary of these interpolated values */ for(n = order + 1; n <= cep_length; n++) { - sum = 0.f; + sum = 0.0; for (k = n - (order - 1); k < n; k++) sum += k * result[k-1] * data[n - k]; result[n-1] = sum / n; @@ -564,11 +536,11 @@ int xtract_lpcc(const float *data, const int N, const void *argv, float *result) return XTRACT_SUCCESS; } -//int xtract_lpcc_s(const float *data, const int N, const void *argv, float *result){ +//int xtract_lpcc_s(const double *data, const int N, const void *argv, double *result){ // return XTRACT_SUCCESS; //} -int xtract_subbands(const float *data, const int N, const void *argv, float *result) +int xtract_subbands(const double *data, const int N, const void *argv, double *result) { int n, bw, xtract_func, nbands, scale, start, lower, *argi, rv; @@ -595,7 +567,7 @@ int xtract_subbands(const float *data, const int N, const void *argv, float *res if(lower >= N || lower + bw >= N) { // printf("n: %d\n", n); - result[n] = 0.f; + result[n] = 0.0; continue; } diff --git a/src/window.c b/src/window.c index bce0770..18a4e15 100644 --- a/src/window.c +++ b/src/window.c @@ -27,108 +27,108 @@ #include "xtract_window_private.h" -void gauss(float *window, const int N, const float sd) +void gauss(double *window, const int N, const double sd) { int n; - const float M = N - 1; - float num, + const double M = N - 1; + double num, den, exponent; for (n = 0; n < N; n++) { - num = n - M / 2.f; - den = sd * M / 2.f; + num = n - M / 2.0; + den = sd * M / 2.0; - exponent = -0.5 * powf(num / den, 2); + exponent = -0.5 * pow(num / den, 2); window[n] = exp(exponent); } } -void hamming(float *window, const int N) +void hamming(double *window, const int N) { int n; - const float M = N - 1; + const double M = N - 1; for (n = 0; n < N; n++) - window[n] = 0.53836 - (0.46164 * cosf(2.0 * PI * (float)n / M)); + window[n] = 0.53836 - (0.46164 * cos(2.0 * PI * (double)n / M)); } -void hann(float *window, const int N) +void hann(double *window, const int N) { int n; - const float M = N - 1; + const double M = N - 1; for (n = 0; n < N; n++) - window[n] = 0.5 * (1.0 - cosf(2.0 * PI * (float)n / M)); + window[n] = 0.5 * (1.0 - cos(2.0 * PI * (double)n / M)); } -void bartlett(float *window, const int N) +void bartlett(double *window, const int N) { int n; - const float M = N - 1; + const double M = N - 1; for (n = 0; n < N; n++) - window[n] = 2.f / M * (M / 2.f - fabsf(n - M / 2.f)); + window[n] = 2.0 / M * (M / 2.0 - fabs(n - M / 2.0)); } -void triangular(float *window, const int N) +void triangular(double *window, const int N) { int n; - const float M = N - 1; + const double M = N - 1; for (n = 0; n < N; n++) - window[n] = 2.f / N * (N / 2.f - fabsf(n - M / 2.f)); + window[n] = 2.0 / N * (N / 2.0 - fabs(n - M / 2.0)); } -void bartlett_hann(float *window, const int N) +void bartlett_hann(double *window, const int N) { int n; - const float M = N - 1, + const double M = N - 1, a0 = 0.62, a1 = 0.5, a2 = 0.38; - float term1 = 0.f, - term2 = 0.f; + double term1 = 0.0, + term2 = 0.0; for (n = 0; n < N; n++) { - term1 = a1 * fabsf(n / M - 0.5); - term2 = a2 * cosf(2.0 * PI * (float)n / M); + term1 = a1 * fabs(n / M - 0.5); + term2 = a2 * cos(2.0 * PI * (double)n / M); window[n] = a0 - term1 - term2; } } -void blackman(float *window, const int N) +void blackman(double *window, const int N) { int n; - const float M = N - 1, + const double M = N - 1, a0 = 0.42, a1 = 0.5, a2 = 0.08; - float term1 = 0.f, - term2 = 0.f; + double term1 = 0.0, + term2 = 0.0; for (n = 0; n < N; n++) { - term1 = a1 * cosf(2.0 * PI * (float)n / M); - term2 = a2 * cosf(4.0 * PI * (float)n / M); + term1 = a1 * cos(2.0 * PI * (double)n / M); + term2 = a2 * cos(4.0 * PI * (double)n / M); window[n] = a0 - term1 + term2; } @@ -137,19 +137,19 @@ void blackman(float *window, const int N) #define BIZ_EPSILON 1E-21 // Max error acceptable /* Based on code from mplayer window.c, and somewhat beyond me */ -float besselI0(float x) +double besselI0(double x) { - float temp; - float sum = 1.0; - float u = 1.0; - float halfx = x/2.0; + double temp; + double sum = 1.0; + double u = 1.0; + double halfx = x/2.0; int n = 1; do { - temp = halfx/(float)n; + temp = halfx/(double)n; u *=temp * temp; sum += u; n++; @@ -161,41 +161,41 @@ float besselI0(float x) } -void kaiser(float *window, const int N, const float alpha) +void kaiser(double *window, const int N, const double alpha) { int n; - const float M = N - 1; - float num; + const double M = N - 1; + double num; for (n = 0; n < N; n++) { - num = besselI0(alpha * sqrtf(1.0 - powf((2.0 * n / M - 1), 2))); + num = besselI0(alpha * sqrt(1.0 - pow((2.0 * n / M - 1), 2))); window[n] = num / besselI0(alpha); } } -void blackman_harris(float *window, const int N) +void blackman_harris(double *window, const int N) { int n; - const float M = N - 1, + const double 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; + double term1 = 0.0, + term2 = 0.0, + term3 = 0.0; 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); + term1 = a1 * cos(2.0 * PI * n / M); + term2 = a2 * cos(4.0 * PI * n / M); + term3 = a3 * cos(6.0 * PI * n / M); window[n] = a0 - term1 + term2 - term3; } diff --git a/src/xtract_window_private.h b/src/xtract_window_private.h index 609010b..5708161 100644 --- a/src/xtract_window_private.h +++ b/src/xtract_window_private.h @@ -32,7 +32,7 @@ * \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); +void gauss(double *window, const int N, const double sd); /** \brief generate a Hamming window * @@ -40,7 +40,7 @@ void gauss(float *window, const int N, const float sd); * \param N the number of elements in the array pointed to by *window * */ -void hamming(float *window, const int N); +void hamming(double *window, const int N); /** \brief generate a Hann window * @@ -48,7 +48,7 @@ void hamming(float *window, const int N); * \param N the number of elements in the array pointed to by *window * */ -void hann(float *window, const int N); +void hann(double *window, const int N); /** \brief generate a Bartlett window * @@ -56,7 +56,7 @@ void hann(float *window, const int N); * \param N the number of elements in the array pointed to by *window * */ -void bartlett(float *window, const int N); +void bartlett(double *window, const int N); /** \brief generate a Triangular window * @@ -64,7 +64,7 @@ void bartlett(float *window, const int N); * \param N the number of elements in the array pointed to by *window * */ -void triangular(float *window, const int N); +void triangular(double *window, const int N); /** \brief generate a Bartlett-Hann window * @@ -72,7 +72,7 @@ void triangular(float *window, const int N); * \param N the number of elements in the array pointed to by *window * */ -void bartlett_hann(float *window, const int N); +void bartlett_hann(double *window, const int N); /** \brief generate a Blackman window * @@ -80,7 +80,7 @@ void bartlett_hann(float *window, const int N); * \param N the number of elements in the array pointed to by *window * */ -void blackman(float *window, const int N); +void blackman(double *window, const int N); /** \brief generate a Kaiser window * @@ -89,7 +89,7 @@ void blackman(float *window, const int N); * \param alpha The larger the value of |alpha|, the narrower the window becomes * */ -void kaiser(float *window, const int N, const float alpha); +void kaiser(double *window, const int N, const double alpha); /** \brief generate a Blackman-Harris window * @@ -97,5 +97,5 @@ void kaiser(float *window, const int N, const float alpha); * \param N the number of elements in the array pointed to by *window * */ -void blackman_harris(float *window, const int N); +void blackman_harris(double *window, const int N); -- cgit v1.2.3