jamie@141: /* jamie@141: * Copyright (C) 2012 Jamie Bullock jamie@86: * jamie@141: * Permission is hereby granted, free of charge, to any person obtaining a copy jamie@141: * of this software and associated documentation files (the "Software"), to jamie@141: * deal in the Software without restriction, including without limitation the jamie@141: * rights to use, copy, modify, merge, publish, distribute, sublicense, and/or jamie@141: * sell copies of the Software, and to permit persons to whom the Software is jamie@141: * furnished to do so, subject to the following conditions: jamie@86: * jamie@141: * The above copyright notice and this permission notice shall be included in jamie@141: * all copies or substantial portions of the Software. jamie@86: * jamie@141: * THE SOFTWARE IS PROVIDED "AS IS", WITHOUT WARRANTY OF ANY KIND, EXPRESS OR jamie@141: * IMPLIED, INCLUDING BUT NOT LIMITED TO THE WARRANTIES OF MERCHANTABILITY, jamie@141: * FITNESS FOR A PARTICULAR PURPOSE AND NONINFRINGEMENT. IN NO EVENT SHALL THE jamie@141: * AUTHORS OR COPYRIGHT HOLDERS BE LIABLE FOR ANY CLAIM, DAMAGES OR OTHER jamie@141: * LIABILITY, WHETHER IN AN ACTION OF CONTRACT, TORT OR OTHERWISE, ARISING jamie@141: * FROM, OUT OF OR IN CONNECTION WITH THE SOFTWARE OR THE USE OR OTHER DEALINGS jamie@141: * IN THE SOFTWARE. jamie@141: * jamie@86: */ jamie@29: jamie@29: #include "xtract/libxtract.h" jamie@29: #include jamie@135: #include jamie@165: #include jamie@165: jamie@166: #ifndef M_PI jamie@166: #define M_PI 3.14159265358979323846264338327 jamie@166: #endif jamie@166: jamie@165: typedef enum waveform_type_ jamie@165: { jamie@165: SINE, jamie@165: SAWTOOTH, jamie@185: SQUARE, jamie@185: NOISE jamie@165: } jamie@165: waveform_type; jamie@29: jamie@135: #define BLOCKSIZE 1024 jamie@185: #define HALF_BLOCKSIZE BLOCKSIZE >> 1 jamie@135: #define SAMPLERATE 44100 jamie@165: #define PERIOD 102 jamie@135: #define MFCC_FREQ_BANDS 13 jamie@135: #define MFCC_FREQ_MIN 20 jamie@135: #define MFCC_FREQ_MAX 20000 jamie@29: jamie@165: jamie@165: double wavetable[BLOCKSIZE]; jamie@165: jamie@165: void fill_wavetable(const float frequency, waveform_type type) jamie@165: { jamie@165: jamie@165: int samples_per_period = SAMPLERATE / frequency; jamie@165: jamie@165: for (int i = 0; i < BLOCKSIZE; ++i) jamie@165: { jamie@165: int phase = i % samples_per_period; jamie@165: jamie@165: switch (type) jamie@165: { jamie@165: case SINE: jamie@165: wavetable[i] = sin((phase / (double)PERIOD) * 2 * M_PI); jamie@165: break; jamie@165: case SQUARE: jamie@165: if (phase < (samples_per_period / 2.f)) jamie@165: { jamie@165: wavetable[i] = -1.0; jamie@165: } jamie@165: else jamie@165: { jamie@165: wavetable[i] = 1.0; jamie@165: } jamie@165: break; jamie@165: case SAWTOOTH: jamie@165: wavetable[i] = ((phase / (double)PERIOD) * 2) - 1.; jamie@165: break; jamie@185: case NOISE: jamie@185: wavetable[i] = (arc4random_uniform(1000) / 500.0) - 1; jamie@185: break; jamie@165: } jamie@165: } jamie@165: } jamie@165: jamie@165: void print_wavetable(void) jamie@165: { jamie@165: for (int i = 0; i < BLOCKSIZE; ++i) jamie@165: { jamie@165: printf("%f\n", wavetable[i]); jamie@165: } jamie@165: } jamie@165: jamie@135: int main(void) jamie@135: { jamie@165: double mean = 0.0; jamie@165: double f0 = 0.0; jamie@185: double flux = 0.0; jamie@165: double centroid = 0.0; jamie@174: double spectrum[BLOCKSIZE] = {0}; jamie@174: double windowed[BLOCKSIZE] = {0}; jamie@174: double peaks[BLOCKSIZE] = {0}; jamie@174: double harmonics[BLOCKSIZE] = {0}; jamie@185: double subframes[BLOCKSIZE] = {0}; jamie@185: double difference[HALF_BLOCKSIZE] = {0}; jamie@165: double *window = NULL; jamie@174: double mfccs[MFCC_FREQ_BANDS] = {0}; jamie@172: double argd[4] = {0}; jamie@165: double samplerate = 44100.0; jamie@135: int n; jamie@135: xtract_mel_filter mel_filters; jamie@120: jamie@185: fill_wavetable(344.53125f, NOISE); // 344.53125f = 128 samples @ 44100 Hz jamie@185: print_wavetable(); jamie@165: jamie@165: /* get the F0 */ jamie@165: xtract[XTRACT_WAVELET_F0](wavetable, BLOCKSIZE, &samplerate, &f0); jamie@165: printf("\nF0: %f\n", f0); jamie@135: jamie@135: /* get the mean of the input */ jamie@165: xtract[XTRACT_MEAN](wavetable, BLOCKSIZE, NULL, &mean); jamie@165: printf("\nInput mean = %.2f\n\n", mean); /* We expect this to be zero for a square wave */ jamie@165: jamie@165: /* create the window function */ jamie@165: window = xtract_init_window(BLOCKSIZE, XTRACT_HANN); jamie@165: xtract_windowed(wavetable, BLOCKSIZE, window, windowed); jamie@174: xtract_free_window(window); jamie@135: jamie@135: /* get the spectrum */ jamie@146: argd[0] = SAMPLERATE / (double)BLOCKSIZE; jamie@146: argd[1] = XTRACT_MAGNITUDE_SPECTRUM; jamie@177: argd[2] = 0.f; /* DC component - we expect this to zero for square wave */ jamie@146: argd[3] = 0.f; /* No Normalisation */ jamie@29: jamie@135: xtract_init_fft(BLOCKSIZE, XTRACT_SPECTRUM); jamie@165: xtract[XTRACT_SPECTRUM](windowed, BLOCKSIZE, &argd[0], spectrum); jamie@185: xtract_free_fft(); jamie@165: jamie@165: xtract[XTRACT_SPECTRAL_CENTROID](spectrum, BLOCKSIZE, NULL, ¢roid); jamie@165: printf("\nSpectral Centroid: %f\n", centroid); jamie@165: jamie@165: argd[1] = 10.0; /* peak threshold as % of maximum peak */ jamie@165: xtract[XTRACT_PEAK_SPECTRUM](spectrum, BLOCKSIZE / 2, argd, peaks); jamie@165: jamie@165: argd[0] = f0; jamie@165: argd[1] = .3; /* harmonic threshold */ jamie@165: xtract[XTRACT_HARMONIC_SPECTRUM](peaks, BLOCKSIZE, argd, harmonics); jamie@120: jamie@135: /* print the spectral bins */ jamie@165: printf("\nSpectrum:\n"); jamie@165: for(n = 0; n < (BLOCKSIZE >> 1); ++n) jamie@165: { jamie@165: printf("freq: %.1f\tamp: %.6f", spectrum[n + (BLOCKSIZE >> 1)], spectrum[n]); jamie@165: if (peaks[n + (BLOCKSIZE >> 1)] != 0.f) jamie@165: { jamie@165: printf("\tpeak:: freq: %.1f\tamp: %.6f\n", peaks[n + (BLOCKSIZE >> 1)], peaks[n]); jamie@165: } jamie@165: else jamie@165: { jamie@165: printf("\n"); jamie@165: } jamie@120: } jamie@120: printf("\n"); jamie@120: jamie@135: /* compute the MFCCs */ jamie@135: mel_filters.n_filters = MFCC_FREQ_BANDS; jamie@155: mel_filters.filters = (double **)malloc(MFCC_FREQ_BANDS * sizeof(double *)); jamie@135: for(n = 0; n < MFCC_FREQ_BANDS; ++n) jamie@135: { jamie@155: mel_filters.filters[n] = (double *)malloc(BLOCKSIZE * sizeof(double)); jamie@135: } jamie@135: jamie@135: xtract_init_mfcc(BLOCKSIZE >> 1, SAMPLERATE >> 1, XTRACT_EQUAL_GAIN, MFCC_FREQ_MIN, MFCC_FREQ_MAX, mel_filters.n_filters, mel_filters.filters); jamie@135: xtract_mfcc(spectrum, BLOCKSIZE >> 1, &mel_filters, mfccs); jamie@135: jamie@135: /* print the MFCCs */ jamie@135: printf("MFCCs:\n"); jamie@135: for(n = 0; n < MFCC_FREQ_BANDS; ++n) jamie@135: { jamie@135: printf("band: %d\t", n); jamie@135: if(n < 10) { jamie@135: printf("\t"); jamie@135: } jamie@135: printf("coeff: %f\n", mfccs[n]); jamie@135: } jamie@135: jamie@185: /* compute Spectral Flux */ jamie@185: argd[0] = SAMPLERATE / HALF_BLOCKSIZE; jamie@185: argd[1] = XTRACT_MAGNITUDE_SPECTRUM; jamie@185: argd[2] = 0.f; /* DC component */ jamie@185: argd[3] = 0.f; /* No Normalisation */ jamie@185: jamie@185: xtract_init_fft(HALF_BLOCKSIZE, XTRACT_SPECTRUM); jamie@185: xtract_features_from_subframes(wavetable, BLOCKSIZE, XTRACT_SPECTRUM, argd, subframes); jamie@185: xtract_difference_vector(subframes, BLOCKSIZE, NULL, difference); jamie@185: jamie@185: argd[0] = 1.0; /* norm order */ jamie@185: argd[1] = XTRACT_POSITIVE_SLOPE; /* positive slope */ jamie@185: jamie@185: xtract_flux(difference, HALF_BLOCKSIZE, argd, &flux); jamie@185: jamie@185: printf("Flux: %f\n", flux); jamie@185: jamie@135: /* cleanup */ jamie@135: for(n = 0; n < MFCC_FREQ_BANDS; ++n) jamie@135: { jamie@135: free(mel_filters.filters[n]); jamie@135: } jamie@135: free(mel_filters.filters); jamie@135: jamie@29: return 0; jamie@135: jamie@29: }