jamie@141: /*
jamie@141:  * Copyright (C) 2012 Jamie Bullock
jamie@141:  *
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@141:  *
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@141:  *
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@141:  */
jamie@4: 
jamie@4: /* calculates the spectral xtract of one frame, given peak frequency and amplitude to first and second inputs respectively */
jamie@4: 
jamie@41: #include "m_pd.h"
jamie@42: #include <math.h>
jamie@52: #include <string.h>
jamie@41: 
jamie@4: #include "xtract/libxtract.h"
jamie@4: 
jamie@62: #define BLOCKSIZE 1024 /* FIX: this should be dynamic - somehow */
jamie@4: #define NYQUIST 22050.0f
jamie@4: 
jamie@85: #ifndef isnan
jamie@85:     /* FIX: should probably try to handle signalling NaNs */
jamie@85:     int isnan(x){ if(x == x) return 0; else return 1;}
jamie@85: #endif
jamie@85: 
jamie@85: #ifndef isinf
jamie@85:     int isinf(x) {if(x == 1.0 / 0. || x == -1.0 / 0.) return 1; else return 0;}
jamie@85: #endif
jamie@85: 
jamie@4: static t_class *xtract_class;
jamie@4: 
jamie@41: /* Struct for keeping track of memory allocations */
jamie@41: typedef struct _tracked_memory {
jamie@41:     char argv;
jamie@41: } tracked_memory;
jamie@41: 
jamie@4: typedef struct _xtract {
jamie@4:     t_object  x_obj;
jamie@4:     t_float f;
jamie@146:     double *window;
jamie@146:     double *data;
jamie@146:     double *result;
jamie@108:     t_int feature,
jamie@108:           is_scalar,
jamie@108:           is_subframe,
jamie@115:           argv_type,
jamie@108:           init_blocksize,
jamie@108:           done_init;
jamie@98:     t_symbol *feature_name;
jamie@41:     tracked_memory memory;
jamie@4:     void *argv;
jamie@4: } t_xtract_tilde;
jamie@4: 
jamie@4: static t_int *xtract_perform(t_int *w) {
jamie@4:     t_sample *in = (t_sample *)(w[1]);
jamie@4:     t_xtract_tilde *x = (t_xtract_tilde *)(w[2]);
jamie@4:     t_int N = (t_int)(w[3]);
jamie@108:     t_int rv = 0;
jamie@146:     double result = 0.0;
jamie@4: 
jamie@146:     for(n = 0; n < N; ++n) {
jamie@146:         x->data[n] = (double)in[n];
jamie@146:     }
jamie@146: 
jamie@146:     rv = xtract[x->feature](x->data, N, x->argv, &result);
jamie@4: 
jamie@108:     if(rv == XTRACT_FEATURE_NOT_IMPLEMENTED)
jamie@38: 	pd_error(x, "Feature not implemented");
jamie@42: 
jamie@42:     /* set nan, inf or -inf to 0 */
jamie@146:     result = (isinf(result) || isnan(result) ? 0.0 : result);
jamie@38:     
jamie@146:     outlet_float(x->x_obj.ob_outlet, (float)result);
jamie@4:     return (w+4);
jamie@4: }
jamie@4: 
jamie@4: static t_int *xtract_perform_vector(t_int *w) {
jamie@98: 
jamie@4:     t_sample *in = (t_sample *)(w[1]);
jamie@4:     t_sample *out = (t_sample *)(w[2]);
jamie@4:     t_xtract_tilde *x = (t_xtract_tilde *)(w[3]);
jamie@4:     t_int N = (t_int)(w[4]), n;
jamie@108:     t_int rv = 0;
jamie@4: 
jamie@98:     if(N != x->init_blocksize && x->done_init){
jamie@98:         error("xtract~ %s: Blocksize mismatch, try specifying the blocksize as a second argument", x->feature_name->s_name);
jamie@98:         return (w+5);
jamie@98:     }
jamie@98: 
jamie@62:     n = N;
jamie@62: 
jamie@146:     for(n = 0; n < N; ++n) {
jamie@146:         x->data[n] = (double)in[n];
jamie@146:     }
jamie@62: 
jamie@104:     if(x->feature == XTRACT_PEAK_SPECTRUM || x->feature == XTRACT_LPC)
jamie@62: 	N >>= 1;
jamie@4:     
jamie@108:     if(x->is_subframe){
jamie@108: 
jamie@146:         rv = xtract_features_from_subframes(x->data, N, x->feature, 
jamie@146:                 x->argv, x->result);
jamie@108:     }
jamie@108:     else{
jamie@108: 
jamie@146:         rv = xtract[x->feature](x->data, N, x->argv, x->result);
jamie@108:     }
jamie@108: 
jamie@108:     if(rv == XTRACT_FEATURE_NOT_IMPLEMENTED)
jamie@38: 	pd_error(x, "Feature not implemented");
jamie@4: 
jamie@108:     while(n--) 
jamie@146:         out[n] = (float)x->result[n];
jamie@4:     
jamie@4:     return (w+5);
jamie@4: }
jamie@4: 
jamie@4: static void xtract_dsp(t_xtract_tilde *x, t_signal **sp) {
jamie@4: 
jamie@115:     if(!x->is_scalar){
jamie@4:         dsp_add(xtract_perform_vector, 4, 
jamie@4:             sp[0]->s_vec, sp[1]->s_vec, x, sp[0]->s_n);
jamie@115:     }
jamie@4:             
jamie@115:     else  
jamie@115:         dsp_add(xtract_perform, 3, sp[0]->s_vec, x, sp[0]->s_n);
jamie@4: 
jamie@4: }
jamie@4: 
jamie@4: static void *xtract_new(t_symbol *me, t_int argc, t_atom *argv) {
jamie@4:     
jamie@4:     t_xtract_tilde *x = (t_xtract_tilde *)pd_new(xtract_class);
jamie@50:     xtract_mel_filter *mf;
jamie@108:     t_int n, N, M, f, F, 
jamie@108:           n_args, 
jamie@108:           type;
jamie@55:     t_float *argv_max;
jamie@108:     t_symbol *arg1;
jamie@56:     xtract_function_descriptor_t *fd;
jamie@108:     char *p_name, 
jamie@108:          *p_desc, 
jamie@108:          *author;
jamie@62:     int year;
jamie@52: 
jamie@52:     p_name = p_desc = author = NULL;
jamie@4:    
jamie@108:     n_args = type = 0;
jamie@50: 
jamie@50:     f = F = XTRACT_FEATURES;
jamie@50: 
jamie@4:     N = BLOCKSIZE;
jamie@4:     
jamie@4:     x->argv = NULL;
jamie@115:     x->argv_type = 0;
jamie@98:     x->done_init = 0;
jamie@108:     x->is_scalar = 0;
jamie@108:     x->is_subframe = 0;
jamie@108:     x->feature = -1;
jamie@4:     
jamie@146:     /* Allocate data area */
jamie@146:     x->data = (double *)getbytes(N * sizeof(double));
jamie@146:     x->result = (double *)getbytes(N * sizeof(double));
jamie@146: 
jamie@108:     /* Parse arguments */
jamie@108:     if(argc){
jamie@108:         arg1 = atom_getsymbol(argv);
jamie@108:         if(arg1 == gensym("subframe"))
jamie@108:             x->is_subframe = 1;
jamie@108:         else
jamie@108:             x->feature_name = atom_getsymbol(argv);
jamie@108:     }
jamie@113:     else {
jamie@113: 	post("xtract~: No arguments given");
jamie@113:         return (void *)x;
jamie@113:     }
jamie@108:     if(argc > 1){
jamie@108:         if(x->is_subframe)
jamie@108:             x->feature_name = atom_getsymbol(argv+1);
jamie@108:         else
jamie@108:             N = atom_getint(argv+1);
jamie@108:     }
jamie@108:     if(argc > 2)
jamie@108:         N = atom_getint(argv+2);
jamie@98: 
jamie@98:     x->init_blocksize = N;
jamie@108:     M = N >> 1;
jamie@4: 
jamie@50:     /* get function descriptors */
jamie@56:     fd = (xtract_function_descriptor_t *)xtract_make_descriptors();
jamie@50: 
jamie@50:     /* iterate over descriptors */
jamie@50:     while(f--){
jamie@50: 	/* map creation arg to feature */
jamie@98: 	if(x->feature_name == gensym(fd[f].algo.name)){ 
jamie@50: 	    x->feature = f;
jamie@50: 	    break;
jamie@50: 	}
jamie@50:     }
jamie@4: 
jamie@108:     if(x->feature == -1)
jamie@108:         post("xtract~: feature not found: %s", x->feature_name->s_name);
jamie@108: 
jamie@41:     /* allocate memory for feature arguments */
jamie@51:     n_args = fd[f].argc;
jamie@50:     type = fd[f].argv.type;
jamie@50: 
jamie@115:     x->argv_type = type;
jamie@115: 
jamie@50:     if(n_args){
jamie@55: 	for(n = 0; n < n_args; n++){
jamie@55: 		    argv_max = &fd[f].argv.max[n]; 
jamie@85: 		    /*post("Argument %d, max: %.2f", n, *argv_max); */
jamie@55: 	}
jamie@56: 	if(type == XTRACT_MEL_FILTER){
jamie@50: 	    x->memory.argv = (size_t)(n_args * sizeof(xtract_mel_filter));
jamie@50: 	    x->argv = (xtract_mel_filter *)getbytes(x->memory.argv);
jamie@50: 	}
jamie@56: 	else if(type == XTRACT_INT){
jamie@50: 	    x->memory.argv = (size_t)(n_args * sizeof(t_int));
jamie@50: 	    x->argv = (t_int *)getbytes(x->memory.argv);
jamie@50: 	}
jamie@56: 	else if (type == XTRACT_FLOAT){
jamie@50: 	    x->memory.argv = (size_t)(n_args * sizeof(t_float));
jamie@50: 	    x->argv = (t_float *)getbytes(x->memory.argv);
jamie@50: 	}
jamie@50: 	else
jamie@50: 	    x->memory.argv = 0;
jamie@41:     }
jamie@52: 
jamie@52:     p_name = fd[f].algo.p_name;
jamie@52:     p_desc = fd[f].algo.p_desc;
jamie@52:     author = fd[f].algo.author;
jamie@62:     year = fd[f].algo.year; 
jamie@52: 
jamie@52:     if(argc){
jamie@52: 	if(strcmp(p_name, ""))	
jamie@52: 	    post("xtract~: %s", p_name );
jamie@52: 	if(strcmp(p_desc, ""))	
jamie@52: 	    post("xtract~: %s", p_desc );
jamie@52: 	if(strcmp(author, "") && year)	
jamie@62: 	    post("xtract~: %s(%d)", author, year);
jamie@52:     }	
jamie@50:     
jamie@108:     /* Adjust frame size if we are using subframe features */
jamie@108:     if(x->is_subframe)
jamie@108:         N = M;
jamie@108: 
jamie@108:     post("xtract~: window size: %d", N);
jamie@41: 
jamie@41:     /* do init if needed */
jamie@56:     if(x->feature == XTRACT_MFCC){
jamie@41: 
jamie@50:         mf = x->argv;
jamie@4:         
jamie@50:         mf->n_filters = 20;
jamie@4:         
jamie@56:         post("xtract~: mfcc: filters = %d", 
jamie@56: 		((xtract_mel_filter *)x->argv)->n_filters);
jamie@50:         mf->filters = 
jamie@50:             (t_float **)getbytes(mf->n_filters * sizeof(t_float *));
jamie@50:         for(n = 0; n < mf->n_filters; n++)
jamie@50:             mf->filters[n] = (float *)getbytes(N * sizeof(float));
jamie@4:                  
jamie@62:         xtract_init_mfcc(N, NYQUIST, XTRACT_EQUAL_GAIN, 80.0f,
jamie@98:                 18000.0f, mf->n_filters, mf->filters);
jamie@98:         x->done_init = 1;
jamie@4:     }
jamie@98:     else if(x->feature == XTRACT_BARK_COEFFICIENTS){
jamie@4:         xtract_init_bark(N, NYQUIST, x->argv);
jamie@98:         x->done_init = 1;
jamie@98:     }
jamie@108:     else if(x->feature == XTRACT_WINDOWED){
jamie@108:         x->window = xtract_init_window(N, XTRACT_HANN);
jamie@108:         x->argv = x->window;
jamie@108:         x->done_init = 1;
jamie@108:     }
jamie@98: 
jamie@98:     /* Initialise fft_plan if required */
jamie@98:     if(x->feature == XTRACT_AUTOCORRELATION_FFT ||
jamie@98:             x->feature == XTRACT_SPECTRUM ||
jamie@98:             x->feature == XTRACT_DCT){
jamie@98:         xtract_init_fft(N, x->feature);
jamie@98:         x->done_init = 1;
jamie@98:     }
jamie@41:     
jamie@108:     if(fd[f].is_scalar)
jamie@108: 	x->is_scalar = 1;
jamie@108: 
jamie@108: /* 
jamie@62:     if(x->feature == XTRACT_AUTOCORRELATION || 
jamie@62: 	    x->feature == XTRACT_AUTOCORRELATION_FFT || 
jamie@62: 	    x->feature == XTRACT_MFCC || x->feature == XTRACT_AMDF || 
jamie@62: 	    x->feature == XTRACT_ASDF|| x->feature == XTRACT_DCT || 
jamie@62: 	    x->feature == XTRACT_BARK_COEFFICIENTS || 
jamie@62: 	    x->feature == XTRACT_SPECTRUM || 
jamie@62: 	    x->feature == XTRACT_PEAK_SPECTRUM || 
jamie@104: 	    x->feature == XTRACT_HARMONIC_SPECTRUM ||
jamie@104:             x->feature == XTRACT_LPC ||
jamie@108:             x->feature == XTRACT_LPCC ||
jamie@108:             x->feature == XTRACT_WINDOWED) 
jamie@62: 	x->feature_type = XTRACT_VECTOR;
jamie@108:         */              
jamie@108:   /*  else if (x->feature == XTRACT_FLUX || x->feature == XTRACT_ATTACK_TIME || 
jamie@108:             x->feature == XTRACT_DECAY_TIME || x->feature == XTRACT_DIFFERENCE_VECTOR) 
jamie@108:         x->feature_type = XTRACT_DELTA; */
jamie@108: /*       
jamie@56:     else x->feature_type = XTRACT_SCALAR;
jamie@108: */
jamie@4: 
jamie@4:     /* argv through right inlet */
jamie@4:     inlet_new(&x->x_obj, &x->x_obj.ob_pd, gensym("list"), gensym("list"));
jamie@4: 
jamie@4:     /* if feature is vector, create signal out */
jamie@108:     if(!x->is_scalar) 
jamie@108:         outlet_new(&x->x_obj, &s_signal);
jamie@4: 
jamie@4:     /* otherwise: float */
jamie@108:     else 
jamie@108:         outlet_new(&x->x_obj, &s_float);
jamie@108: 
jamie@108:     if(x->is_scalar && x->is_subframe)
jamie@108:         post(
jamie@108:         "xtract~: warning: subframes not yet supported for scalar features");
jamie@4:     
jamie@50:     /* free the function descriptors */
jamie@50:     xtract_free_descriptors(fd);
jamie@50:     
jamie@4:     return (void *)x;
jamie@4: }
jamie@4: 
jamie@4: static void xtract_tilde_get_args(t_xtract_tilde *x, t_symbol *selector, 
jamie@4: t_int argc, t_atom *argv) {
jamie@4: /*
jamie@4:     if(argc > (t_int)sizeof(x->argv) / 
jamie@4:                         (t_int)sizeof(t_float) || x->argv == NULL)
jamie@4:         post("Too many parameters to right inlet");
jamie@4:     else{*/
jamie@4: 
jamie@4:     x->argv = getbytes(argc * sizeof(float));
jamie@4:     
jamie@115:     while(argc--){
jamie@115:         switch(x->argv_type){
jamie@115:             case XTRACT_INT:
jamie@115:                 ((t_int *)x->argv)[argc] = (int)atom_getfloat(&argv[argc]);
jamie@115:                 break;
jamie@115:             case XTRACT_FLOAT:
jamie@115:             default:
jamie@115:                 ((t_float *)x->argv)[argc] = atom_getfloat(&argv[argc]);
jamie@115:                 break;
jamie@115:         }
jamie@115:     }
jamie@4:  /*   }*/
jamie@4: }
jamie@4: 
jamie@26: static void xtract_tilde_show_help(t_xtract_tilde *x, t_symbol *s){
jamie@26:     
jamie@26:     int i;
jamie@26: 
jamie@85:     xtract_function_descriptor_t *fd, *d;
jamie@85: 
jamie@26:     i = XTRACT_FEATURES;
jamie@26: 
jamie@62:     fd = (xtract_function_descriptor_t *)xtract_make_descriptors();
jamie@62:     post("\nxtract~: Feature List\n");
jamie@62:    
jamie@26:     while(i--){
jamie@62: 	d = &fd[i];
jamie@62: 	post("\t%s", d->algo.name);
jamie@62:     }
jamie@62: 
jamie@62:     xtract_free_descriptors(fd);
jamie@26: }
jamie@26: 
jamie@4: static void xtract_tilde_free(t_xtract_tilde *x) {
jamie@41: 
jamie@41:     if(x->argv != NULL && x->memory.argv)
jamie@41:         freebytes(x->argv, x->memory.argv);
jamie@108: 
jamie@108:     if(x->window != NULL)
jamie@108:         xtract_free_window(x->window);
jamie@4: }
jamie@4: 
jamie@4: void xtract_tilde_setup(void) {
jamie@4:     xtract_class = class_new(gensym("xtract~"),
jamie@4:     (t_newmethod)xtract_new,
jamie@4:     (t_method)xtract_tilde_free, 
jamie@4:     sizeof(t_xtract_tilde),
jamie@4:     CLASS_DEFAULT, 
jamie@4:     A_GIMME, 0);
jamie@4: 
jamie@4:     class_addmethod(xtract_class,
jamie@26: 	    (t_method)xtract_dsp, gensym("dsp"), 0);
jamie@4:     class_addmethod(xtract_class,
jamie@26: 	    (t_method)xtract_tilde_get_args, gensym("list"), A_GIMME, 0);
jamie@26:     class_addmethod(xtract_class,
jamie@26: 	    (t_method)xtract_tilde_show_help, gensym("help"), A_DEFSYMBOL, 0);
jamie@4:     CLASS_MAINSIGNALIN(xtract_class, t_xtract_tilde, f);
jamie@33:     class_sethelpsymbol(xtract_class, gensym("xtract-help"));
jamie@4: }