This commit was generated by cvs2svn to compensate for changes in r2637,svn2git-root

which included commits to RCS files with non-trunk default branches.

svn path=/trunk/externals/iem/iemmatrix/; revision=2638

1 files changed, 171 insertions, 0 deletions

diff --git a/src/mtx_pivot.c b/src/mtx_pivot.c
new file mode 100644
index 0000000..1ac4804
--- /dev/null
+++ b/src/mtx_pivot.c
@@ -0,0 +1,171 @@
+/*
+ *  iemmatrix
+ *
+ *  objects for manipulating simple matrices
+ *  mostly refering to matlab/octave matrix functions
+ *
+ * (c) IOhannes m zmölnig, forum::für::umläute
+ * 
+ * IEM, Graz
+ *
+ * this code is published under the LGPL
+ *
+ */
+#include "iemmatrix.h"
+
+
+/* mtx_pivot */
+static t_class *mtx_pivot_class;
+
+typedef struct _mtx_pivot
+{
+  t_object x_obj;
+
+  t_matrix m;  // the output matrix
+  t_matrix m_pre;  // the pre -multiply matrix
+  t_matrix m_post; // the post-multiply matrix
+
+  t_outlet *pivo, *pre, *post;
+  
+} t_mtx_pivot;
+
+static void mtx_pivot_matrix(t_mtx_pivot *x, t_symbol *s, int argc, t_atom *argv)
+{
+  int row=atom_getfloat(argv);
+  int col=atom_getfloat(argv+1);
+  t_atom *m_pre, *m_post;
+  int i, j, k;
+  int min_rowcol = (row<col)?row:col;
+  t_matrixfloat *buffer, *buf;
+  int *i_pre, *i_post, *i_buf;
+
+  int pivot_row, pivot_col;
+
+  if (argc<2){    post("mtx_pivot: crippled matrix");    return;  }
+  if ((col<1)||(row<1)) {    post("mtx_pivot: invalid dimensions");    return;  }
+  if (col*row>argc-2){    post("sparse matrix not yet supported : use \"mtx_check\"");    return;  }
+
+  adjustsize(&x->m, row, col);
+  adjustsize(&x->m_pre, row, row);
+  adjustsize(&x->m_post,col, col);
+  matrix_set(&x->m_pre, 0);
+  matrix_set(&x->m_post, 0);
+
+  buffer = matrix2float(argv);
+  i_pre    = (int *)getbytes(sizeof(int)*row);
+  i_post   = (int *)getbytes(sizeof(int)*col);
+
+  /* clear pre&post matrices */
+  i=row;
+  i_buf=i_pre;
+  while(i--)*i_buf++=row-i-1;
+  i=col;
+  i_buf=i_post;
+  while(i--)*i_buf++=col-i-1;
+
+  /* do the pivot thing */
+
+  for (k=0; k<min_rowcol-1; k++){
+    // 1. find max_element
+    t_matrixfloat max = 0;
+    pivot_row = pivot_col = k;
+
+    for(i=k; i<row; i++){
+      buf=buffer+col*i+k;
+
+      j=col-k;
+      while(j--){
+	t_matrixfloat f = fabsf(*buf++);
+	if (f>max) {
+	  max=f;
+	  pivot_row = i;
+	  pivot_col = col-j-1;
+	}
+      }
+    }
+    // 2. move max el to [k,k]
+    // 2a swap rows
+    if (k-pivot_row) {
+      t_matrixfloat *oldrow=buffer+col*k;
+      t_matrixfloat *newrow=buffer+col*pivot_row;
+
+      i=col;
+      while(i--){
+	t_matrixfloat f=*oldrow;
+	*oldrow++=*newrow;
+	*newrow++=f;
+      }
+      i=i_pre[k];
+      i_pre[k]=i_pre[pivot_row];
+      i_pre[pivot_row]=i;
+    }
+    // 2b swap columns
+    if (k-pivot_col) {
+      t_matrixfloat *oldcol=buffer+k;
+      t_matrixfloat *newcol=buffer+pivot_col;
+
+      i=row;
+      while(i--){
+	t_matrixfloat f=*oldcol;
+	*oldcol=*newcol;
+	*newcol=f;
+	oldcol+=col;
+	newcol+=col;
+      }
+      i=i_post[k];
+      i_post[k]=i_post[pivot_col];
+      i_post[pivot_col]=i;
+    }
+  }
+
+  float2matrix(x->m.atombuffer, buffer);
+
+  i=col;
+  m_post = x->m_post.atombuffer+2;
+  while(i--){
+    SETFLOAT(m_post+i_post[i]*col+i, 1);
+  }
+  i=row;
+  m_pre = x->m_pre.atombuffer+2;
+  while(i--)SETFLOAT(m_pre+i_pre[i]+i*col, 1);
+
+  
+  outlet_anything(x->post, gensym("matrix"), 2+col*col, x->m_post.atombuffer);
+  outlet_anything(x->pre, gensym("matrix"), 2+row*row, x->m_pre.atombuffer);
+  outlet_anything(x->pivo , gensym("matrix"), 2+row*col, x->m.atombuffer );
+}
+
+static void mtx_pivot_free(t_mtx_pivot *x)
+{
+  matrix_free(&x->m);
+  matrix_free(&x->m_pre);
+  matrix_free(&x->m_post);
+}
+
+static void *mtx_pivot_new(void)
+{
+  t_mtx_pivot *x = (t_mtx_pivot *)pd_new(mtx_pivot_class);
+
+  x->pivo = outlet_new(&x->x_obj, 0);
+  x->pre  = outlet_new(&x->x_obj, 0);
+  x->post = outlet_new(&x->x_obj, 0);
+
+  x->m.atombuffer = x->m_pre.atombuffer = x->m_post.atombuffer = 0;
+  x->m.row = x->m.col = x->m_pre.row = x->m_pre.col = x->m_post.row = x->m_post.col = 0;
+
+  return(x);
+}
+
+void mtx_pivot_setup(void)
+{
+  mtx_pivot_class = class_new(gensym("mtx_pivot"), (t_newmethod)mtx_pivot_new, (t_method)mtx_pivot_free,
+			      sizeof(t_mtx_pivot), 0, 0, 0);
+  class_addmethod(mtx_pivot_class, (t_method)mtx_pivot_matrix, gensym("matrix"), A_GIMME, 0);
+
+  class_sethelpsymbol(mtx_pivot_class, gensym("iemmatrix/mtx_transpose"));
+}
+
+void iemtx_pivot_setup(void)
+{
+  mtx_pivot_setup();
+}