From cdf8ded57019d5c905f16422d40be7b1a18ab3bc Mon Sep 17 00:00:00 2001 From: Tom Schouten Date: Tue, 21 Jan 2003 10:18:19 +0000 Subject: This commit was generated by cvs2svn to compensate for changes in r350, which included commits to RCS files with non-trunk default branches. svn path=/trunk/externals/creb/; revision=351 --- modules/bdiag.c | 272 ++++++++++++++++++++++++++++++++++++++++++++++++++++++++ 1 file changed, 272 insertions(+) create mode 100644 modules/bdiag.c (limited to 'modules/bdiag.c') diff --git a/modules/bdiag.c b/modules/bdiag.c new file mode 100644 index 0000000..37a0349 --- /dev/null +++ b/modules/bdiag.c @@ -0,0 +1,272 @@ +/* + * bdiag.c - block diagonal state space system + * treats input dsp block as n parallel signals + * + * s1 = (a * s1) + (b * s2) + u1; + * s2 = (a * s2) - (b * s1) + u2; + * + * Copyright (c) 2000-2003 by Tom Schouten + * + * This program is free software; you can redistribute it and/or modify + * it under the terms of the GNU General Public License as published by + * the Free Software Foundation; either version 2 of the License, or + * (at your option) any later version. + * + * This program is distributed in the hope that it will be useful, + * but WITHOUT ANY WARRANTY; without even the implied warranty of + * MERCHANTABILITY or FITNESS FOR A PARTICULAR PURPOSE. See the + * GNU General Public License for more details. + * + * You should have received a copy of the GNU General Public License + * along with this program; if not, write to the Free Software + * Foundation, Inc., 675 Mass Ave, Cambridge, MA 02139, USA. + */ + + +#include "m_pd.h" +#include +#include +#include + +#define MAXORDER 64 + +typedef struct bdiagctl +{ + t_float *c_state; + t_float *c_eigen; + t_int c_order; +} t_bdiagctl; + +typedef struct bdiag +{ + t_object x_obj; + t_float x_f; + t_bdiagctl x_ctl; +} t_bdiag; + + +static float randfloat(void){ + float r = rand (); + r /= (RAND_MAX/2); + r -= 1; + return r; + +} + +static void bdiag_random(t_bdiag *x) +{ + int i; + + for (i=0; ix_ctl.c_order; i++) + { + x->x_ctl.c_state[i] = randfloat(); + } + +} + + +static void bdiag_reset(t_bdiag *x) +{ + int i; + + for (i=0; ix_ctl.c_order; i++) + { + x->x_ctl.c_state[i] = 0; + } + +} + + + + + + +static void bdiag_eigen(t_bdiag *x, t_floatarg index, t_floatarg aval, t_floatarg bval) +{ + int i = (int)index; + if (i<0) return; + if (i>=x->x_ctl.c_order/2) return; + x->x_ctl.c_eigen[2*i+0] = aval; + x->x_ctl.c_eigen[2*i+1] = bval; + +} + +/* set decay time and frequency of pole at index */ +static void bdiag_timefreq(t_bdiag *x, t_floatarg index, t_floatarg time, t_floatarg freq) +{ + float r,a,b,n; + float sr = sys_getsr() / (float)x->x_ctl.c_order; + + /* time in ms */ + time *= 0.001; + + if (time < 0.0f) time = 0.0f; + r = pow(0.001f, 1.0f / (time * sr)); + if (r < 0.0f) r = 0.0f; + if (r > 1.0f) r = 1.0f; + + a = cos(2*M_PI*freq/sr); + b = sin(2*M_PI*freq/sr); + + /* normalize to be sure */ + n = 1.0f / sqrt(a*a + b*b); + a *= n; + b *= n; + + bdiag_eigen(x, index, r*a, r*b); +} + +static void bdiag_preset(t_bdiag *x, t_floatarg preset) +{ + int p = preset; + int i; + float a, b, w, r; + + switch(p){ + case 0: + post("preset 0"); + for (i=0; ix_ctl.c_order/2; i++){ + w = randfloat() * .001; + r = 1. - (((float)i + 1.)/1000.); + a = cos(w) * r; + b = sin(w) * r; + post("%f %f %f %f", w, r, a, b); + bdiag_eigen(x,i,a,b); + } + break; + case 1: + default: + break; + + } +} + +static t_int *bdiag_perform(t_int *w) +{ + + + t_float *in = (float *)(w[3]); + t_float *out = (float *)(w[4]); + t_bdiagctl *ctl = (t_bdiagctl *)(w[1]); + + t_float *eigen = ctl->c_eigen; + t_float *state = ctl->c_state; + t_int n = (t_int)(w[2]); + + t_float u1,u2,a,b,s1,s2,s1new,s2new; + + int i; + + for (i=0; is_n; + int i; + + if (n == 1) + { + post("bdiag: doesnt work with blocksize == 1"); + dsp_add_copy(sp[0]->s_vec, sp[1]->s_vec, sp[0]->s_n); + } + else + { + if (x->x_ctl.c_order != n) + { + if (x->x_ctl.c_state) free(x->x_ctl.c_state); + if (x->x_ctl.c_eigen) free(x->x_ctl.c_eigen); + + x->x_ctl.c_state = (t_float *)malloc(n*sizeof(t_float)); + x->x_ctl.c_eigen = (t_float *)malloc(n*sizeof(t_float)); + + for(i=0;ix_ctl.c_state[i] = 0; + x->x_ctl.c_eigen[i] = 0; + } + + x->x_ctl.c_order = n; + } + + + dsp_add(bdiag_perform, 4, &x->x_ctl, sp[0]->s_n, sp[0]->s_vec, sp[1]->s_vec); + } + +} + +static void bdiag_free(t_bdiag *x) +{ + + if (x->x_ctl.c_state) free(x->x_ctl.c_state); + if (x->x_ctl.c_eigen) free(x->x_ctl.c_eigen); + + +} + +t_class *bdiag_class; + + +static void *bdiag_new(t_floatarg permute) +{ + t_bdiag *x = (t_bdiag *)pd_new(bdiag_class); + int i, n=64; + + outlet_new(&x->x_obj, gensym("signal")); + + x->x_ctl.c_state = (t_float *)malloc(n*sizeof(t_float)); + x->x_ctl.c_eigen = (t_float *)malloc(n*sizeof(t_float)); + + for(i=0;ix_ctl.c_state[i] = 0; + x->x_ctl.c_eigen[i] = 0; + } + + x->x_ctl.c_order = n; + + + return (void *)x; + + +} + + +void bdiag_tilde_setup(void) +{ + //post("bdiag~ v0.1"); + bdiag_class = class_new(gensym("bdiag~"), (t_newmethod)bdiag_new, + (t_method)bdiag_free, sizeof(t_bdiag), 0, A_DEFFLOAT, 0); + CLASS_MAINSIGNALIN(bdiag_class, t_bdiag, x_f); + class_addmethod(bdiag_class, (t_method)bdiag_random, gensym("random"), 0); + class_addmethod(bdiag_class, (t_method)bdiag_random, gensym("bang"), 0); + class_addmethod(bdiag_class, (t_method)bdiag_reset, gensym("reset"), 0); + class_addmethod(bdiag_class, (t_method)bdiag_dsp, gensym("dsp"), 0); + + class_addmethod(bdiag_class, (t_method)bdiag_eigen, gensym("eigen"), A_DEFFLOAT, A_DEFFLOAT, A_DEFFLOAT, 0); + class_addmethod(bdiag_class, (t_method)bdiag_timefreq, gensym("timefreq"), A_DEFFLOAT, A_DEFFLOAT, A_DEFFLOAT, 0); + class_addmethod(bdiag_class, (t_method)bdiag_preset, gensym("preset"), A_DEFFLOAT, 0); +} + -- cgit v1.2.1