aboutsummaryrefslogtreecommitdiff
path: root/src/z_mtx.c
diff options
context:
space:
mode:
Diffstat (limited to 'src/z_mtx.c')
-rw-r--r--src/z_mtx.c669
1 files changed, 669 insertions, 0 deletions
diff --git a/src/z_mtx.c b/src/z_mtx.c
new file mode 100644
index 0000000..785c7e9
--- /dev/null
+++ b/src/z_mtx.c
@@ -0,0 +1,669 @@
+/* 1605:forum::für::umläute:2001 */
+
+/* objects for manipulating matrices */
+/* mostly i refer to matlab/octave matrix functions */
+
+/*
+ matrix : basic object : create and store matrices
+ mtx : alias for matrix
+
+ mtx_resize
+ mtx_row
+ mtx_col
+ mtx_element
+
+ mtx_ones
+ mtx_zeros
+ mtx_eye
+ mtx_egg
+
+ mtx_diag
+ mtx_diegg
+ mtx_trace
+
+ mtx_mean
+ mtx_rand
+
+ mtx_transpose
+ mtx_scroll
+ mtx_roll
+
+ mtx_add
+ mtx_+
+ mtx_mul
+ mtx_*
+ mtx_.*
+ mtx_./
+
+ mtx_inverse
+ mtx_pivot
+
+ mtx_size
+
+ mtx_check
+ mtx_print
+*/
+
+#include "zexy.h"
+#include <math.h>
+
+#ifdef NT
+#include <memory.h>
+#define fabsf fabs
+#pragma warning( disable : 4244 )
+#pragma warning( disable : 4305 )
+#endif
+
+#define VALID_MTX 12450
+
+typedef struct _matrix
+{
+ int row;
+ int col;
+
+ int valid;
+
+ t_float *buffer;
+
+} t_matrix;
+
+/* intern utility functions */
+static void MTX_setdimen(t_matrix *x, int row, int col)
+{
+ x->col = col;
+ x->row = row;
+}
+
+static void MTX_adjustsize(t_matrix *m, int desiredRow, int desiredCol)
+{
+ int col=m->col, row=m->row;
+
+ if (desiredRow<1){
+ post("cannot make less than 1 rows");
+ desiredRow=1;
+ }
+ if (desiredCol<1){
+ post("cannot make less than 1 columns");
+ desiredCol=1;
+ }
+
+ if (col*row!=desiredRow*desiredCol){
+ if(m->buffer)freebytes(m->buffer, (col*row)*sizeof(t_float));
+ m->buffer=(t_float *)getbytes((desiredCol*desiredRow)*sizeof(t_float));
+ }
+
+ MTX_setdimen(m, desiredRow, desiredCol);
+ return;
+}
+
+static void MTX_debugmtx(t_matrix *m, int id)
+{
+ int i=m->col;
+ t_float *buf = m->buffer;
+ while(i--){
+ int j=m->row;
+ startpost("debug%d: ", id);
+ while(j--)
+ startpost("%f ", *buf++);
+ endpost();
+ }
+}
+
+static void MTX_freematrix(t_matrix *x)
+{
+ freebytes(x->buffer, x->row*x->col*sizeof(t_float));
+ x->row = x->col = x->valid = 0;
+}
+
+static void outlet_matrix(t_outlet *x, t_matrix *m)
+{
+ t_gpointer *gp=(t_gpointer *)m;
+ outlet_pointer(x, gp);
+}
+static t_matrix *get_inmatrix(t_gpointer *gp)
+{
+ t_matrix *m = (t_matrix *)gp;
+ if (m->valid != VALID_MTX){
+ error("matrix: no valid matrix passed !");
+ return 0;
+ }
+ if ((m->col <= 0) || (m->row <= 0)){
+ error("matrix: strange dimensions");
+ return 0;
+ }
+ return m;
+}
+static void MTX_set(t_matrix *m, t_float f)
+{
+ t_float *buf=m->buffer;
+ int size = m->col*m->row;
+ if(buf)while(size--)*buf++=f;
+}
+
+/* -------------------- matrix ------------------------------ */
+
+static t_class *matrix_class;
+
+typedef struct _mtx
+{
+ t_object x_obj;
+
+ t_matrix *m;
+
+ int current_row, current_col; /* this makes things easy for the mtx_row & mtx_col...*/
+
+ t_canvas *x_canvas; /* and for reading / writing */
+} t_mtx;
+
+/* core functions */
+
+static void matrix_bang(t_mtx *x)
+{ /* output the matrix */
+ post("matrix %dx%d @ %x", x->m->row, x->m->col, x->m->buffer);
+ outlet_matrix(x->x_obj.ob_outlet, x->m);
+}
+static void matrix_matrix2(t_mtx *x, t_gpointer *gp)
+{
+ t_matrix *m = get_inmatrix(gp);
+ if(!m)return;
+
+ if (x->m->row*x->m->col != m->row*m->col) {
+ freebytes(x->m->buffer, x->m->row*x->m->col*sizeof(t_float));
+ x->m->buffer = copybytes(m->buffer, m->row*m->col*sizeof(t_float));
+ } else memcpy(x->m->buffer, m->buffer, m->row*m->col*sizeof(t_float));
+
+ MTX_setdimen(x->m, m->row, m->col);
+}
+static void matrix_matrix(t_mtx *x, t_gpointer *gp)
+{
+ t_matrix *m = get_inmatrix(gp);
+ if(!m)return;
+
+ matrix_matrix2(x, gp);
+ outlet_matrix(x->x_obj.ob_outlet, x->m);
+}
+
+
+/* basic functions */
+
+static void matrix_zeros(t_mtx *x, t_symbol *s, int argc, t_atom *argv)
+{
+ int col, row;
+ switch(argc) {
+ case 0: /* zero out the actual matrix */
+ break;
+ case 1:
+ row=atom_getfloat(argv);
+ MTX_adjustsize(x->m, row, row);
+ break;
+ default:
+ row=atom_getfloat(argv++);
+ col=atom_getfloat(argv);
+ MTX_adjustsize(x->m, row, col);
+ }
+ MTX_set(x->m, 0);
+ outlet_matrix(x->x_obj.ob_outlet, x->m);
+}
+static void matrix_ones(t_mtx *x, t_symbol *s, int argc, t_atom *argv)
+{
+ int col, row;
+ switch(argc) {
+ case 0: /* zero out the actual matrix */
+ break;
+ case 1:
+ row=atom_getfloat(argv);
+ MTX_adjustsize(x->m, row, row);
+ break;
+ default:
+ row=atom_getfloat(argv++);
+ col=atom_getfloat(argv);
+ MTX_adjustsize(x->m, row, col);
+ }
+ MTX_set(x->m, 1);
+ outlet_matrix(x->x_obj.ob_outlet, x->m);
+}
+static void matrix_eye(t_mtx *x, t_symbol *s, int argc, t_atom *argv)
+{
+ int col, row, n;
+ switch(argc) {
+ case 0: /* zero out the actual matrix */
+ break;
+ case 1:
+ row=atom_getfloat(argv);
+ MTX_adjustsize(x->m, row, row);
+ break;
+ default:
+ row=atom_getfloat(argv++);
+ col=atom_getfloat(argv);
+ MTX_adjustsize(x->m, row, col);
+ }
+ MTX_set(x->m, 0);
+
+ col=x->m->col;
+ row=x->m->row;
+ n = (col<row)?col:row;
+ while(n--)x->m->buffer[n*(1+col)]=1;
+ outlet_matrix(x->x_obj.ob_outlet, x->m);
+}
+static void matrix_egg(t_mtx *x, t_symbol *s, int argc, t_atom *argv)
+{
+ int col, row, n;
+ switch(argc) {
+ case 0: /* zero out the actual matrix */
+ break;
+ case 1:
+ row=atom_getfloat(argv);
+ MTX_adjustsize(x->m, row, row);
+ break;
+ default:
+ row=atom_getfloat(argv++);
+ col=atom_getfloat(argv);
+ MTX_adjustsize(x->m, row, col);
+ }
+ MTX_set(x->m, 0);
+
+ col=x->m->col;
+ row=x->m->row;
+ n = (col<row)?col:row;
+ while(n--)x->m->buffer[(n+1)*(col-1)]=1;
+ outlet_matrix(x->x_obj.ob_outlet, x->m);
+}
+static void matrix_diag(t_mtx *x, t_symbol *s, int argc, t_atom *argv)
+{
+ int col=argc;
+ argv+=argc-1;
+ if (argc<1) {
+ post("matrix: no diagonale present");
+ return;
+ }
+ MTX_adjustsize(x->m, argc, argc);
+ MTX_set(x->m, 0);
+
+ while(argc--)x->m->buffer[argc*(1+col)]=atom_getfloat(argv--);
+ outlet_matrix(x->x_obj.ob_outlet, x->m);
+}
+static void matrix_diegg(t_mtx *x, t_symbol *s, int argc, t_atom *argv)
+{
+ int col=argc;
+ argv+=argc-1;
+ if (argc<1) {
+ post("matrix: no dieggonale present");
+ return;
+ }
+ MTX_adjustsize(x->m, argc, argc);
+ MTX_set(x->m, 0);
+ while(argc--)x->m->buffer[(argc+1)*(col-1)]=atom_getfloat(argv--);
+ outlet_matrix(x->x_obj.ob_outlet, x->m);
+}
+static void matrix_float(t_mtx *x, t_float f)
+{
+ MTX_set(x->m, f);
+ outlet_matrix(x->x_obj.ob_outlet, x->m);
+}
+static void matrix_row(t_mtx *x, t_symbol *s, int argc, t_atom *argv)
+{
+ t_atom *ap=(t_atom *)getbytes(sizeof(t_atom)*x->m->col);
+ t_float *mtx_flt;
+ int row=x->m->row, col=x->m->col;
+ int r, c;
+ t_float f;
+
+ switch (argc){
+ case 0: /* show all rows as a list of floats */
+ for (r=0;r<row;r++){
+ mtx_flt = x->m->buffer+r*col;
+ argv = ap;
+ c=col;
+ while(c--){
+ SETFLOAT(argv, *mtx_flt);
+ argv++;
+ }
+ outlet_list(x->x_obj.ob_outlet, gensym("row"), col, ap);
+ }
+ break;
+ case 1: /* show this row as a list of floats */
+ r=atom_getfloat(argv)-1;
+ if ((r<0)||(r>=row)){
+ post("matrix: row index %d is out of range", r+1);
+ return;
+ }
+ argv=ap;
+ mtx_flt=x->m->buffer+r*col;
+ c=col;
+ while(c--){
+ SETFLOAT(argv, *mtx_flt);
+ argv++;
+ }
+ outlet_list(x->x_obj.ob_outlet, gensym("row"), col, ap);
+ break;
+ case 2: /* set this row to a constant value */
+ r=atom_getfloat(argv)-1;
+ f=atom_getfloat(argv+1);
+ if ((r<0)||(r>=row)){
+ post("matrix: row index %d is out of range", r+1);
+ return;
+ }
+ mtx_flt = x->m->buffer+r*col;
+ c=col;
+ while(c--)*mtx_flt++=f;
+ outlet_matrix(x->x_obj.ob_outlet, x->m);
+ default: /* set this row to new values */
+ r=atom_getfloat(argv++)-1;
+ if (argc--<col){
+ post("matrix: sparse rows not yet supported : use \"mtx_check\"");
+ return;
+ }
+ if ((r<0)||(r>=row)){
+ post("matrix: row index %d is out of range", r+1);
+ return;
+ }
+ mtx_flt=x->m->buffer+col*r;
+ c=col;
+ while(c--)*mtx_flt++=atom_getfloat(argv++);
+ outlet_matrix(x->x_obj.ob_outlet, x->m);
+ }
+ freebytes(ap, x->m->col*sizeof(t_atom));
+}
+static void matrix_col(t_mtx *x, t_symbol *s, int argc, t_atom *argv)
+{
+ t_float *mtx_flt;
+ int row=x->m->row, col=x->m->col;
+ int c=0, r=0;
+ t_float f;
+ t_atom *ap=(t_atom *)getbytes(row*sizeof(t_atom));
+
+ switch (argc){
+ case 0: /* show all columns as a list of floats */
+ for (c=0;c<col;c++) {
+ mtx_flt = x->m->buffer+c;
+
+ for (r=0;r<row;r++)SETFLOAT(&ap[r], mtx_flt[col*r]);
+ outlet_list(x->x_obj.ob_outlet, gensym("col"), row, ap);
+ }
+ break;
+ case 1:/* show this column as a list of floats */
+ c=atom_getfloat(argv)-1;
+ if ((c<0)||(c>=col)){
+ post("matrix: col index %d is out of range", c+1);
+ return;
+ }
+ mtx_flt = x->m->buffer+c;
+
+ for (r=0;r<row;r++)SETFLOAT(&ap[r],x->m->buffer[col*r]);
+ outlet_list(x->x_obj.ob_outlet, gensym("col"), row, ap);
+ break;
+ case 2: /* set this column to a constant value */
+ c=atom_getfloat(argv)-1;
+ f=atom_getfloat(argv+1);
+ if ((c<0)||(c>=row)){
+ post("matrix: col index %d is out of range", c+1);
+ return;
+ }
+ mtx_flt = x->m->buffer+r*col;
+ c=col;
+ while(c--){
+ *mtx_flt=f;
+ mtx_flt+=col;
+ }
+ outlet_matrix(x->x_obj.ob_outlet, x->m);
+ default:/* set this column to new values */
+ c=atom_getfloat(argv++)-1;
+ if (argc--<row){
+ post("matrix: sparse cols not yet supported : use \"mtx_check\"");
+ return;
+ }
+ if ((c<0)||(c>=col)){
+ post("matrix: col index %d is out of range", c+1);
+ return;
+ }
+ mtx_flt=x->m->buffer+c;
+ r=row;
+ while(r--){
+ *mtx_flt=atom_getfloat(argv++);
+ mtx_flt+=col;
+ }
+ outlet_matrix(x->x_obj.ob_outlet, x->m);
+ }
+ freebytes(ap, row*sizeof(t_atom));
+}
+static void matrix_element(t_mtx *x, t_symbol *s, int argc, t_atom *argv)
+{
+ t_float *mtx_flt=x->m->buffer;
+ int row=x->m->row, col=x->m->col;
+ int r, c, i=row*col;
+
+ switch (argc){
+ case 0:
+ while(i--)outlet_float(x->x_obj.ob_outlet, *mtx_flt++);
+ break;
+ case 1:
+ r=c=atom_getfloat(argv)-1;
+ if ((r<0)||(r>=row)){
+ post("matrix: row index %d is out of range", r+1);
+ return;
+ }
+ if ((c<0)||(c>=col)){
+ post("matrix: col index %d is out of range", c+1);
+ return;
+ }
+ outlet_float(x->x_obj.ob_outlet, mtx_flt[c+r*col]);
+ break;
+ case 2:
+ r=atom_getfloat(argv++)-1;
+ c=atom_getfloat(argv++)-1;
+ if ((r<0)||(r>=row)){ post("matrix: row index %d is out of range", r+1); return; }
+ if ((c<0)||(c>=col)){ post("matrix: col index %d is out of range", c+1); return; }
+ outlet_float(x->x_obj.ob_outlet, mtx_flt[c+r*col]);
+ break;
+ default:
+ r=atom_getfloat(argv++)-1;
+ c=atom_getfloat(argv++)-1;
+ if ((r<0)||(r>=row)){ post("matrix: row index %d is out of range", r+1); return; }
+ if ((c<0)||(c>=col)){ post("matrix: col index %d is out of range", c+1); return; }
+ mtx_flt[c+r*col]=atom_getfloat(argv);
+ outlet_matrix(x->x_obj.ob_outlet, x->m);
+ }
+}
+
+/* ------------- file I/O ------------------ */
+static void matrix_read(t_mtx *x, t_symbol *filename)
+{
+ t_binbuf *bbuf = binbuf_new();
+ t_atom *ap;
+ int n;
+
+ if (binbuf_read_via_path(bbuf, filename->s_name, canvas_getdir(x->x_canvas)->s_name, 0))
+ error("matrix: failed to read %s", filename->s_name);
+
+ ap=binbuf_getvec(bbuf);
+ n =binbuf_getnatom(bbuf)-1;
+
+ if ((ap->a_type == A_SYMBOL) && !strcmp(ap->a_w.w_symbol->s_name,"matrix") ){
+ /* ok, this looks like a matrix */
+ int row = atom_getfloat(ap+1);
+ int col = atom_getfloat(ap+2);
+ t_float *mtx_flt;
+ if (n-2<row*col){
+ error("matrix: sparse matrices not supported");
+ binbuf_free(bbuf);
+ return;
+ }
+ if (row<1 || col<1){
+ error("matrix: illegal matrix dimensions");
+ binbuf_free(bbuf);
+ return;
+ }
+
+ MTX_adjustsize(x->m, row, col);
+ ap+=2;
+ mtx_flt=x->m->buffer;
+ n=row*col;
+ while(n--)*mtx_flt++=atom_getfloat(ap++);
+ }
+
+ binbuf_free(bbuf);
+}
+static void matrix_write(t_mtx *x, t_symbol *filename)
+{
+ t_binbuf *bbuf = binbuf_new();
+ t_atom atom;
+ t_float *mtx_flt=x->m->buffer;
+ char buf[MAXPDSTRING];
+ int n;
+ int r = x->m->row;
+ int c = x->m->col;
+ canvas_makefilename(x->x_canvas, filename->s_name, buf, MAXPDSTRING);
+
+ SETSYMBOL(&atom, gensym("matrix"));
+ binbuf_add(bbuf, 1, &atom);
+ SETFLOAT (&atom, x->m->row);
+ binbuf_add(bbuf, 1, &atom);
+ SETFLOAT (&atom, x->m->col);
+ binbuf_add(bbuf, 1, &atom);
+ binbuf_addsemi(bbuf);
+ n=x->m->row*x->m->col;
+ while(r--){
+ c=x->m->col;
+ while(c--){
+ SETFLOAT(&atom, *mtx_flt++);
+ binbuf_add(bbuf, 1, &atom);
+ }
+ binbuf_addsemi(bbuf);
+ }
+
+ if (binbuf_write(bbuf, buf, "", 1)){
+ error("matrix: failed to write %s", filename->s_name);
+ }
+
+ binbuf_free(bbuf);
+}
+
+/* ----------------- setup ------------------- */
+
+static void matrix_free(t_mtx *x)
+{
+ MTX_freematrix(x->m);
+ freebytes(x->m, sizeof(t_matrix));
+}
+static void *matrix_new(t_symbol *s, int argc, t_atom *argv)
+{
+ t_mtx *x = (t_mtx *)pd_new(matrix_class);
+ int row, col;
+
+ x->m = (t_matrix *)getbytes(sizeof(t_matrix));
+
+ inlet_new(&x->x_obj, &x->x_obj.ob_pd, gensym("pointer"), gensym(""));
+ outlet_new(&x->x_obj, 0);
+
+ x->m->buffer= 0;
+ x->m->row = x->m->col = 0;
+ x->m->valid = VALID_MTX;
+ x->x_canvas = canvas_getcurrent();
+
+ switch (argc) {
+ case 0:
+ row = col = 0;
+ break;
+ case 1:
+ if (argv->a_type == A_SYMBOL) {
+ matrix_read(x, argv->a_w.w_symbol);
+ return(x);
+ }
+ row = col = atom_getfloat(argv);
+ break;
+ default:
+ row = atom_getfloat(argv++);
+ col = atom_getfloat(argv++);
+ }
+
+ if(row*col){
+ MTX_adjustsize(x->m, row, col);
+ MTX_set(x->m, 0);
+ }
+
+ return (x);
+}
+
+static void matrix_setup(void)
+{
+ matrix_class = class_new(gensym("mtx"), (t_newmethod)matrix_new,
+ (t_method)matrix_free, sizeof(t_mtx), 0, A_GIMME, 0);
+
+ /* the core : functions for matrices */
+ class_addmethod (matrix_class, (t_method)matrix_matrix, gensym("matrix"), A_GIMME, 0);
+ class_addmethod (matrix_class, (t_method)matrix_matrix2, gensym(""), A_GIMME, 0);
+
+ /* the basics : functions for creation */
+ class_addmethod (matrix_class, (t_method)matrix_eye, gensym("eye"), A_GIMME, 0);
+ class_addmethod (matrix_class, (t_method)matrix_diag, gensym("diag"), A_GIMME, 0);
+ class_addmethod (matrix_class, (t_method)matrix_ones, gensym("ones"), A_GIMME, 0);
+ class_addmethod (matrix_class, (t_method)matrix_zeros, gensym("zeros"), A_GIMME, 0);
+ class_addmethod (matrix_class, (t_method)matrix_egg, gensym("egg"), A_GIMME, 0);
+ class_addmethod (matrix_class, (t_method)matrix_diegg, gensym("diegg"), A_GIMME, 0);
+
+ /* the rest : functions for anything */
+ class_addbang (matrix_class, matrix_bang);
+ class_addfloat (matrix_class, matrix_float);
+ // class_addlist (matrix_class, matrix_list);
+ class_addmethod (matrix_class, (t_method)matrix_row, gensym("row"), A_GIMME, 0);
+ class_addmethod (matrix_class, (t_method)matrix_col, gensym("column"), A_GIMME, 0);
+ class_addmethod (matrix_class, (t_method)matrix_col, gensym("col"), A_GIMME, 0);
+ class_addmethod (matrix_class, (t_method)matrix_element, gensym("element"), A_GIMME, 0);
+
+ /* the file functions */
+ class_addmethod (matrix_class, (t_method)matrix_write, gensym("write"), A_SYMBOL, 0);
+ class_addmethod (matrix_class, (t_method)matrix_read , gensym("read") , A_SYMBOL, 0);
+
+
+ class_sethelpsymbol(matrix_class, gensym("zexy/mtx"));
+}
+
+
+/* mtx_print */
+
+static t_class *mtx_print_class;
+
+typedef struct _m_print
+{
+ t_object x_obj;
+
+} t_m_print;
+static void mtx_print(t_m_print *x, t_gpointer *gp)
+{
+ t_matrix *m = get_inmatrix(gp);
+ int c, r;
+ t_float *f;
+ if (!m) return;
+
+ c=m->col;
+ r=m->row;
+ f=m->buffer;
+
+ post("matrix %dx%d @ %x", r,c, f);
+
+ while(r--){
+ c=m->col;
+ while(c--)startpost("%f\t", *f++);
+ endpost();
+ }
+
+}
+static void *mtx_print_new(void)
+{
+ t_m_print *x = (t_m_print *)pd_new(mtx_print_class);
+
+ return (x);
+}
+static void mtx_print_setup(void)
+{
+ mtx_print_class = class_new(gensym("m_print"), (t_newmethod)mtx_print_new,
+ 0, sizeof(t_m_print), 0, 0, 0);
+ class_addpointer(mtx_print_class, mtx_print);
+}
+
+
+void z_mtx_setup(void)
+{
+ matrix_setup();
+ mtx_print_setup();
+}