order and matrix pattern is now independent

svn path=/trunk/externals/iemlib/; revision=14654

2 files changed, 167 insertions, 66 deletions

diff --git a/iemlib1/sparse_FIR~-help.pd b/iemlib1/sparse_FIR~-help.pd
index 3eb7ba2..867a9d4 100644
--- a/iemlib1/sparse_FIR~-help.pd
+++ b/iemlib1/sparse_FIR~-help.pd
@@ -1,4 +1,4 @@
-#N canvas 156 22 842 608 10;
+#N canvas 125 49 971 642 10;
 #X floatatom 4 133 5 0 0 0 - - -;
 #X text 32 96 % cpu;
 #X obj 4 93 dsp;
@@ -9,7 +9,7 @@
 #X obj 112 347 gainvu~ 300;
 #X obj 174 369 vu 15 120 empty empty 8 -8 0 10 -66577 -1 1 0;
 #X obj 294 216 vsl 15 128 0 127 0 1 empty empty gain -7 -8 0 10 -225271
--42246 -90881 7500 1;
+-42246 -90881 8400 1;
 #X obj 113 95 noise~;
 #X floatatom 173 512 7 0 0 0 - - -;
 #X floatatom 17 113 5 0 0 0 - - -;
@@ -19,7 +19,6 @@
 #X text 190 555 iem.at;
 #X text 153 577 Graz \, Austria;
 #X floatatom 142 323 5 0 0 0 - - -;
-#X obj 113 179 sparse_FIR~ 256;
 #X msg 177 112 size 256;
 #X msg 367 535 0 1 1 1;
 #X msg 367 556 0 1 40 -1 80 1 120 -1 160 1;
@@ -45,17 +44,22 @@
 0.04 26 0.04 27 0.04 28 0.04 29 0.04;
 #X obj 668 81 pack 0 0.04;
 #X obj 668 59 for++ 0 29;
-#X text 366 228 A list of alternating index- and value- doubles will
+#X text 366 217 A list of alternating index- and value- doubles will
 change the coefficients of convolution. The number of pairs should
 be less than the order of FIR.;
-#X text 367 276 A matrix message will do the same.;
+#X text 365 291 A matrix message will do the same.;
 #X obj 568 93 loadbang;
 #X text 134 16 convolve a signal with a coefficent list of non zero
 index/value pairs. All values between the determined values will be
 set to zero and won't be calculated (sparse).;
-#X text 227 182 1.arg: <float> maximum convolution-length or order
+#X text 227 182 1.arg: <float> convolution-length (FIR order or size)
 ;
-#X text 176 93 will change the max. order of FIR filter.;
+#X text 177 74 will change the max. order of FIR filter (or size).
+;
+#X text 366 263 Each pair or double begins with an index followed by
+its value.;
+#X obj 113 180 sparse_FIR~ 256;
+#X msg 175 90 size 16;
 #X connect 2 0 0 0;
 #X connect 2 1 11 0;
 #X connect 3 0 2 0;
@@ -66,23 +70,24 @@ set to zero and won't be calculated (sparse).;
 #X connect 7 0 10 0;
 #X connect 7 1 4 0;
 #X connect 8 0 6 1;
-#X connect 9 0 18 0;
-#X connect 18 0 6 0;
-#X connect 19 0 18 0;
-#X connect 20 0 18 0;
-#X connect 21 0 18 0;
-#X connect 22 0 18 0;
-#X connect 23 0 36 0;
-#X connect 23 1 24 0;
-#X connect 24 0 34 0;
-#X connect 25 0 34 0;
-#X connect 26 0 23 0;
-#X connect 27 0 18 0;
-#X connect 28 0 18 0;
-#X connect 29 0 18 0;
-#X connect 32 0 18 0;
-#X connect 33 0 18 0;
-#X connect 34 0 18 0;
-#X connect 35 0 25 0;
-#X connect 36 0 35 0;
-#X connect 39 0 34 0;
+#X connect 9 0 43 0;
+#X connect 18 0 43 0;
+#X connect 19 0 43 0;
+#X connect 20 0 43 0;
+#X connect 21 0 43 0;
+#X connect 22 0 35 0;
+#X connect 22 1 23 0;
+#X connect 23 0 33 0;
+#X connect 24 0 33 0;
+#X connect 25 0 22 0;
+#X connect 26 0 43 0;
+#X connect 27 0 43 0;
+#X connect 28 0 43 0;
+#X connect 31 0 43 0;
+#X connect 32 0 43 0;
+#X connect 33 0 43 0;
+#X connect 34 0 24 0;
+#X connect 35 0 34 0;
+#X connect 38 0 33 0;
+#X connect 43 0 6 0;
+#X connect 44 0 43 0;
diff --git a/iemlib1/src/sparse_FIR~.c b/iemlib1/src/sparse_FIR~.c
index 0a1043a..8287179 100644
--- a/iemlib1/src/sparse_FIR~.c
+++ b/iemlib1/src/sparse_FIR~.c
@@ -16,11 +16,13 @@ typedef struct _sparse_FIR_tilde
   t_object  x_obj;
   t_float   *x_coef_beg;
   int       *x_index_beg;
+  int       x_n_coef_resp_order;
   int       x_n_coef;
   int       x_n_coef_malloc;
   t_float   *x_history_beg;
+  int       x_n_order;
+  int       x_n_order_malloc;
   int       x_rw_index;
-  int       x_sparse_FIR_order;
   t_float   x_msi;
 } t_sparse_FIR_tilde;
 
@@ -34,8 +36,8 @@ static t_int *sparse_FIR_tilde_perform(t_int *w)
   int n = (t_int)(w[4]);
   int rw_index = x->x_rw_index;
   int i, j, ix;
-  int order = x->x_sparse_FIR_order;
-  int n_coef = x->x_n_coef;
+  int order = x->x_n_order;
+  int n_coef = x->x_n_coef_resp_order;
   int n_coef8;
   t_float sum=0.0f;
   t_float *coef = x->x_coef_beg;
@@ -105,27 +107,68 @@ sparse_FIR_tilde_perf_zero:
   return(w+5);
 }
 
+static void sparse_FIR_tilde_sort_within(t_sparse_FIR_tilde *x)
+{
+  int cur_order = x->x_n_order;
+  int n_coef = x->x_n_coef;
+  int index, i;
+  int n_coef_resp_order = 0;
+  int *index_pointer_within = x->x_index_beg;
+  t_float *coef_pointer_within = x->x_coef_beg;
+  int *index_pointer = x->x_index_beg + x->x_n_coef_malloc;
+  t_float *coef_pointer = x->x_coef_beg + x->x_n_coef_malloc;
+  t_float coef;
+  
+  for(i=0; i<n_coef; i++)
+  {
+    index = index_pointer[i];
+    coef = coef_pointer[i];
+    if((index >= 0) && (index < cur_order))
+    {
+      index_pointer_within[i] = -index; /* negate index for FIR direction */
+      coef_pointer_within[i] = coef;
+      n_coef_resp_order++;
+    }
+  }
+  x->x_n_coef_resp_order = n_coef_resp_order;
+}
+
 static void sparse_FIR_tilde_list(t_sparse_FIR_tilde *x, t_symbol *s, int argc, t_atom *argv)
 {
-  int order = x->x_sparse_FIR_order;
-  int n_arg2 = argc/2, index, i;
+  int max_order = x->x_n_order_malloc;
+  int n_pair_arg = argc/2, index, i;
   int n_coef = 0;
-  int *index_pointer = x->x_index_beg;
-  t_float *coef_pointer = x->x_coef_beg;
+  int *index_pointer;
+  t_float *coef_pointer;
   t_float coef;
   
-  for(i=0; i<n_arg2; i++)
+  if(n_pair_arg > 0)
   {
-    index = (int)atom_getfloat(argv++);
-    coef = (t_float)atom_getfloat(argv++);
-    if((index >= 0) && (index < order))
+    if(n_pair_arg > x->x_n_coef_malloc) /* resize */
+    {
+      x->x_index_beg =  (int *)resizebytes(x->x_index_beg, 2*x->x_n_coef_malloc*sizeof(int), 2*n_pair_arg*sizeof(int));
+      x->x_coef_beg =  (t_float *)resizebytes(x->x_coef_beg, 2*x->x_n_coef_malloc*sizeof(t_float), 2*n_pair_arg*sizeof(t_float));
+      x->x_n_coef_malloc = n_pair_arg;
+    }
+    
+    index_pointer = x->x_index_beg + x->x_n_coef_malloc;
+    coef_pointer = x->x_coef_beg + x->x_n_coef_malloc;
+    
+    for(i=0; i<n_pair_arg; i++)
     {
-      *index_pointer++ = -index;
-      *coef_pointer++ = coef;
-      n_coef++;
+      index = (int)atom_getfloat(argv++);
+      coef = (t_float)atom_getfloat(argv++);
+      if((index >= 0) && (index < max_order))
+      {
+        *index_pointer++ = index;
+        *coef_pointer++ = coef;
+        n_coef++;
+      }
     }
+    x->x_n_coef = n_coef;
+    
+    sparse_FIR_tilde_sort_within(x);
   }
-  x->x_n_coef = n_coef;
 }
 
 static void sparse_FIR_tilde_matrix(t_sparse_FIR_tilde *x, t_symbol *s, int argc, t_atom *argv)
@@ -160,19 +203,19 @@ static void sparse_FIR_tilde_matrix(t_sparse_FIR_tilde *x, t_symbol *s, int argc
 
 static void sparse_FIR_tilde_order(t_sparse_FIR_tilde *x, t_floatarg fn)
 {
-  int n = (int)fn;
+  int n_order = (int)fn;
   
-  if(n > 0)
+  if(n_order > 0)
   {
-    if(n > x->x_sparse_FIR_order)/* resize */
+    if(n_order > x->x_n_order_malloc) /* resize */
     {
-      x->x_history_beg =  (t_float *)resizebytes(x->x_history_beg, 2*x->x_n_coef_malloc*sizeof(t_float), 2*n*sizeof(t_float));
-      x->x_index_beg =  (int *)resizebytes(x->x_index_beg, x->x_n_coef_malloc*sizeof(int), n*sizeof(int));
-      x->x_coef_beg =  (t_float *)resizebytes(x->x_coef_beg, x->x_n_coef_malloc*sizeof(t_float), n*sizeof(t_float));
-      x->x_n_coef_malloc = n;
+      x->x_history_beg =  (t_float *)resizebytes(x->x_history_beg, 2*x->x_n_order_malloc*sizeof(t_float), 2*n_order*sizeof(t_float));
+      x->x_n_order_malloc = n_order;
     }
-    x->x_sparse_FIR_order = n;
+    x->x_n_order = n_order;
     x->x_rw_index = 0;
+    
+    sparse_FIR_tilde_sort_within(x);
   }
 }
 
@@ -184,32 +227,84 @@ static void sparse_FIR_tilde_dsp(t_sparse_FIR_tilde *x, t_signal **sp)
 static void *sparse_FIR_tilde_new(t_floatarg fn)
 {
   t_sparse_FIR_tilde *x = (t_sparse_FIR_tilde *)pd_new(sparse_FIR_tilde_class);
-  int n=(int)fn;
+  int n_order=(int)fn;
+  int i;
   
   outlet_new(&x->x_obj, &s_signal);
-  x->x_msi = 0;
-  x->x_n_coef = 0;
-  if(n < 1)
-    n = 1;
-  x->x_sparse_FIR_order = n;
-  x->x_n_coef_malloc = n;
-  x->x_history_beg = (t_float *)getbytes((2*x->x_n_coef_malloc)*sizeof(t_float));
-  x->x_index_beg = (int *)getbytes(x->x_n_coef_malloc*sizeof(int));
-  x->x_coef_beg = (t_float *)getbytes(x->x_n_coef_malloc*sizeof(t_float));
+  
+  x->x_n_coef = 1;
+  x->x_n_coef_resp_order = 1;
+  x->x_n_coef_malloc = 1;
+  x->x_index_beg = (int *)getbytes(2*x->x_n_coef_malloc*sizeof(int));
+  x->x_coef_beg = (t_float *)getbytes(2*x->x_n_coef_malloc*sizeof(t_float));
+  x->x_index_beg[0] = 0;
+  x->x_index_beg[1] = 0;
+  x->x_coef_beg[0] = 0.0f;
+  x->x_coef_beg[1] = 0.0f;
+  if(n_order < 1)
+    n_order = 1;
+  x->x_n_order = n_order;
+  x->x_n_order_malloc = n_order;
+  x->x_history_beg = (t_float *)getbytes((2*x->x_n_order_malloc)*sizeof(t_float));
   x->x_rw_index = 0;
+  n_order = 2*x->x_n_order_malloc;
+  for(i=0; i<n_order; i++)
+    x->x_history_beg[i] = 0.0f;
+  
+  x->x_msi = 0;
+  
+  post("NEW: n_coef_resp_order = %d, n_coef = %d, n_coef_malloc = %d, n_order = %d, n_order_malloc = %d", x->x_n_coef_resp_order, x->x_n_coef, x->x_n_coef_malloc, x->x_n_order, x->x_n_order_malloc);
+  
   return(x);
 }
 
 static void sparse_FIR_tilde_free(t_sparse_FIR_tilde *x)
 {
-  if(x->x_history_beg)
-    freebytes(x->x_history_beg, (2*x->x_n_coef_malloc)*sizeof(t_float));
-  if(x->x_index_beg)
-    freebytes(x->x_index_beg, x->x_n_coef_malloc*sizeof(int));
-  if(x->x_coef_beg)
-    freebytes(x->x_coef_beg, x->x_n_coef_malloc*sizeof(t_float));
+  freebytes(x->x_history_beg, (2*x->x_n_order_malloc)*sizeof(t_float)); /* twice, because of my simple circle-buffer */
+  freebytes(x->x_index_beg, 2*x->x_n_coef_malloc*sizeof(int)); /* twice, because of buffering both, all coefficients and only the relevant for current order */
+  freebytes(x->x_coef_beg, 2*x->x_n_coef_malloc*sizeof(t_float)); /* twice, because of buffering both, all coefficients and only the relevant for current order */
 }
 
+/*static void sparse_FIR_tilde_dump(t_sparse_FIR_tilde *x)
+{
+  t_float *hist=x->x_history_beg;
+  int *ix=x->x_index_beg;
+  int n=x->x_n_order;
+  
+  post("n_coef_resp_order = %d, n_coef = %d, n_coef_malloc = %d, n_order = %d, n_order_malloc = %d", x->x_n_coef_resp_order, x->x_n_coef, x->x_n_coef_malloc, x->x_n_order, x->x_n_order_malloc);
+  post("HIST:");
+  
+  while(n > 8)
+  {
+    post("hist = %g, %g, %g, %g, %g, %g, %g, %g", hist[n-1], hist[n-2], hist[n-3], hist[n-4], hist[n-5], hist[n-6], hist[n-7], hist[n-8]);
+    n -= 8;
+    hist -= 8;
+  }
+  while(n > 0)
+  {
+    post("hist = %g", hist[n-1]);
+    n--;
+    hist--;
+  }
+  post("COEF:");
+  
+  hist = x->x_coef_beg;
+  n = x->x_n_coef_resp_order;
+  while(n > 8)
+  {
+    post("coef = %d@%g, %d@%g, %d@%g, %d@%g, %d@%g, %d@%g, %d@%g, %d@%g", ix[n-1],hist[n-1], ix[n-2],hist[n-2], ix[n-3],hist[n-3], ix[n-4],hist[n-4], ix[n-5],hist[n-5], ix[n-6],hist[n-6], ix[n-7],hist[n-7], ix[n-8],hist[n-8]);
+    n -= 8;
+    hist -= 8;
+  }
+  while(n > 0)
+  {
+    post("coef = %d@%g", ix[n-1],hist[n-1]);
+    n--;
+    hist--;
+  }
+  post("***********************");
+}*/
+
 void sparse_FIR_tilde_setup(void)
 {
   sparse_FIR_tilde_class = class_new(gensym("sparse_FIR~"), (t_newmethod)sparse_FIR_tilde_new,
@@ -220,4 +315,5 @@ void sparse_FIR_tilde_setup(void)
   class_addmethod(sparse_FIR_tilde_class, (t_method)sparse_FIR_tilde_matrix, gensym("matrix"), A_GIMME, 0);
   class_addmethod(sparse_FIR_tilde_class, (t_method)sparse_FIR_tilde_order, gensym("order"), A_FLOAT, 0);
   class_addmethod(sparse_FIR_tilde_class, (t_method)sparse_FIR_tilde_order, gensym("size"), A_FLOAT, 0);
+  //class_addmethod(sparse_FIR_tilde_class, (t_method)sparse_FIR_tilde_dump, gensym("dump"), 0);
 }