aboutsummaryrefslogtreecommitdiff
path: root/externals/gridflow/base/flow_objects.c
diff options
context:
space:
mode:
Diffstat (limited to 'externals/gridflow/base/flow_objects.c')
-rw-r--r--externals/gridflow/base/flow_objects.c1195
1 files changed, 0 insertions, 1195 deletions
diff --git a/externals/gridflow/base/flow_objects.c b/externals/gridflow/base/flow_objects.c
deleted file mode 100644
index f971553d..00000000
--- a/externals/gridflow/base/flow_objects.c
+++ /dev/null
@@ -1,1195 +0,0 @@
-/*
- $Id: flow_objects.c,v 1.2 2006-03-15 04:37:08 matju Exp $
-
- GridFlow
- Copyright (c) 2001,2002,2003,2004,2005 by Mathieu Bouchard
-
- 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.
-
- See file ../COPYING for further informations on licensing terms.
-
- 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., 59 Temple Place - Suite 330, Boston, MA 02111-1307, USA.
-*/
-
-#include <sys/time.h>
-#include <stdlib.h>
-#include <math.h>
-#include "grid.h.fcs"
-
-// BAD HACK: GCC complains: unimplemented (--debug|--debug-harder mode only)
-#ifdef HAVE_DEBUG
-#define SCOPY(a,b,n) COPY(a,b,n)
-#else
-#define SCOPY(a,b,n) SCopy<n>::f(a,b)
-#endif
-
-template <int n> class SCopy {
-public: template <class T> static inline void __attribute__((always_inline)) f(Pt<T> a, Pt<T> b) {
- *a=*b; SCopy<n-1>::f(a+1,b+1);}};
-template <> class SCopy<0> {
-public: template <class T> static inline void __attribute__((always_inline)) f(Pt<T> a, Pt<T> b) {}};
-
-/*template <> class SCopy<4> {
-public: template <class T>
- static inline void __attribute__((always_inline)) f(Pt<T> a, Pt<T> b) {
- *a=*b; SCopy<3>::f(a+1,b+1);}
- // wouldn't gcc 2.95 complain here?
- static inline void __attribute__((always_inline)) f(Pt<uint8> a, Pt<uint8> b)
- { *(int32 *)a=*(int32 *)b; }
-};*/
-
-Numop *op_add, *op_sub, *op_mul, *op_div, *op_mod, *op_shl, *op_and, *op_put;
-
-static void expect_dim_dim_list (P<Dim> d) {
- if (d->n!=1) RAISE("dimension list should be Dim[n], not %s",d->to_s());}
-//static void expect_min_one_dim (P<Dim> d) {
-// if (d->n<1) RAISE("minimum 1 dimension");}
-static void expect_max_one_dim (P<Dim> d) {
- if (d->n>1) { RAISE("expecting Dim[] or Dim[n], got %s",d->to_s()); }}
-//static void expect_exactly_one_dim (P<Dim> d) {
-// if (d->n!=1) { RAISE("expecting Dim[n], got %s",d->to_s()); }}
-
-//****************************************************************
-\class GridCast < GridObject
-struct GridCast : GridObject {
- \attr NumberTypeE nt;
- \decl void initialize (NumberTypeE nt);
- \grin 0
-};
-
-GRID_INLET(GridCast,0) {
- out = new GridOutlet(this,0,in->dim,nt);
-} GRID_FLOW {
- out->send(n,data);
-} GRID_END
-
-\def void initialize (NumberTypeE nt) {
- rb_call_super(argc,argv);
- this->nt = nt;
-}
-
-\classinfo { IEVAL(rself,"install '#cast',1,1"); }
-\end class GridCast
-
-//****************************************************************
-//{ ?,Dim[B] -> Dim[*Cs] }
-// out0 nt to be specified explicitly
-\class GridImport < GridObject
-struct GridImport : GridObject {
- \attr NumberTypeE cast;
- \attr P<Dim> dim; // size of grids to send
- PtrGrid dim_grid;
- GridImport() { dim_grid.constrain(expect_dim_dim_list); }
- ~GridImport() {}
- \decl void initialize(Ruby x, NumberTypeE cast=int32_e);
- \decl void _0_cast(NumberTypeE cast);
- \decl void _0_reset();
- \decl void _0_symbol(Symbol x);
- \decl void _0_list(...);
- \decl void _1_per_message();
- \grin 0
- \grin 1 int32
- template <class T> void process (int n, Pt<T> data) {
- while (n) {
- if (!out || !out->dim) out = new GridOutlet(this,0,dim?dim:in[0]->dim,cast);
- int32 n2 = min((int32)n,out->dim->prod()-out->dex);
- out->send(n2,data);
- n-=n2; data+=n2;
- }
- }
-};
-
-GRID_INLET(GridImport,0) {} GRID_FLOW { process(n,data); } GRID_END
-GRID_INPUT(GridImport,1,dim_grid) { dim = dim_grid->to_dim(); } GRID_END
-
-\def void _0_symbol(Symbol x) {
- const char *name = rb_sym_name(argv[0]);
- int n = strlen(name);
- if (!dim) out=new GridOutlet(this,0,new Dim(n));
- process(n,Pt<uint8>((uint8 *)name,n));
-}
-
-\def void _0_cast(NumberTypeE cast) { this->cast = cast; }
-
-\def void _0_list(...) {
- if (in.size()<1 || !in[0]) _0_grid(0,0); //HACK: enable grid inlet...
- in[0]->from_ruby_list(argc,argv,cast);
-}
-
-\def void _1_per_message() { dim=0; dim_grid=0; }
-
-\def void initialize(Ruby x, NumberTypeE cast) {
- rb_call_super(argc,argv);
- this->cast = cast;
- if (argv[0]!=SYM(per_message)) {
- dim_grid=new Grid(argv[0]);
- dim = dim_grid->to_dim();
- }
-}
-
-\def void _0_reset() {
- STACK_ARRAY(int32,foo,1); *foo=0;
- while (out->dim) out->send(1,foo);
-}
-
-\classinfo { IEVAL(rself,"install '#import',2,1"); }
-\end class GridImport
-
-//****************************************************************
-/*{ Dim[*As] -> ? }*/
-/* in0: integer nt */
-\class GridExport < GridObject
-struct GridExport : GridObject {
- \grin 0
-};
-
-template <class T>
-static Ruby INTORFLOAT2NUM(T value) {return INT2NUM(value);}
-static Ruby INTORFLOAT2NUM(int64 value) {return gf_ll2num(value);}
-static Ruby INTORFLOAT2NUM(float32 value) {return rb_float_new(value);}
-static Ruby INTORFLOAT2NUM(float64 value) {return rb_float_new(value);}
-static Ruby INTORFLOAT2NUM(ruby value) {return value.r;}
-
-GRID_INLET(GridExport,0) {
-} GRID_FLOW {
- for (int i=0; i<n; i++) {
- Ruby a[] = { INT2NUM(0), INTORFLOAT2NUM(data[i]) };
- send_out(COUNT(a),a);
- }
-} GRID_END
-\classinfo { IEVAL(rself,"install '#export',1,1"); }
-\end class GridExport
-
-/* **************************************************************** */
-/*{ Dim[*As] -> ? }*/
-/* in0: integer nt */
-\class GridExportList < GridObject
-struct GridExportList : GridObject {
- Ruby /*Array*/ list;
- int n;
- \grin 0
-};
-
-GRID_INLET(GridExportList,0) {
- int n = in->dim->prod();
- if (n>250000) RAISE("list too big (%d elements)", n);
- list = rb_ary_new2(n+2);
- this->n = n;
- rb_ivar_set(rself,SI(@list),list); // keep
- rb_ary_store(list,0,INT2NUM(0));
- rb_ary_store(list,1,bsym._list);
-} GRID_FLOW {
- for (int i=0; i<n; i++, data++)
- rb_ary_store(list,in->dex+i+2,INTORFLOAT2NUM(*data));
-} GRID_FINISH {
- send_out(rb_ary_len(list),rb_ary_ptr(list));
- list = 0;
- rb_ivar_set(rself,SI(@list),Qnil); // unkeep
-} GRID_END
-
-\classinfo { IEVAL(rself,"install '#export_list',1,1"); }
-\end class GridExportList
-
-/* **************************************************************** */
-// GridStore ("@store") is the class for storing a grid and restituting
-// it on demand. The right inlet receives the grid. The left inlet receives
-// either a bang (which forwards the whole image) or a grid describing what
-// to send.
-//{ Dim[*As,B],Dim[*Cs,*Ds] -> Dim[*As,*Ds] }
-// in0: integer nt
-// in1: whatever nt
-// out0: same nt as in1
-\class GridStore < GridObject
-struct GridStore : GridObject {
- PtrGrid r; // can't be \attr
- PtrGrid put_at; // can't be //\attr
- \attr Numop *op;
- int32 wdex [Dim::MAX_DIM]; // temporary buffer, copy of put_at
- int32 fromb[Dim::MAX_DIM];
- int32 to2 [Dim::MAX_DIM];
- int lsd; // lsd = Last Same Dimension (for put_at)
- int d; // goes with wdex
- \decl void initialize (Grid *r=0);
- \decl void _0_bang ();
- \decl void _0_op (Numop *op);
- \decl void _1_reassign ();
- \decl void _1_put_at (Grid *index);
- \grin 0 int
- \grin 1
- GridStore() { put_at.constrain(expect_max_one_dim); }
- template <class T> void compute_indices(Pt<T> v, int nc, int nd);
-};
-
-// takes the backstore of a grid and puts it back into place. a backstore
-// is a grid that is filled while the grid it would replace has not
-// finished being used.
-static void snap_backstore (PtrGrid &r) {
- if (r.next) {r=r.next.p; r.next=0;}
-}
-
-template <class T> void GridStore::compute_indices(Pt<T> v, int nc, int nd) {
- for (int i=0; i<nc; i++) {
- uint32 wrap = r->dim->v[i];
- bool fast = lowest_bit(wrap)==highest_bit(wrap); // is power of two?
- if (i) {
- if (fast) op_shl->map(nd,v,(T)highest_bit(wrap));
- else op_mul->map(nd,v,(T)wrap);
- }
- if (fast) op_and->map(nd,v+nd*i,(T)(wrap-1));
- else op_mod->map(nd,v+nd*i,(T)(wrap));
- if (i) op_add->zip(nd,v,v+nd*i);
- }
-}
-
-// !@#$ i should ensure that n is not exceedingly large
-// !@#$ worse: the size of the foo buffer may still be too large
-GRID_INLET(GridStore,0) {
- // snap_backstore must be done before *anything* else
- snap_backstore(r);
- int na = in->dim->n;
- int nb = r->dim->n;
- int nc = in->dim->get(na-1);
- STACK_ARRAY(int32,v,Dim::MAX_DIM);
- if (na<1) RAISE("must have at least 1 dimension.",na,1,1+nb);
- int lastindexable = r->dim->prod()/r->dim->prod(nc) - 1;
- int ngreatest = nt_greatest((T *)0);
- if (lastindexable > ngreatest) {
- RAISE("lastindexable=%d > ngreatest=%d (ask matju)",lastindexable,ngreatest);
- }
- if (nc > nb)
- RAISE("wrong number of elements in last dimension: "
- "got %d, expecting <= %d", nc, nb);
- int nd = nb - nc + na - 1;
- COPY(v,in->dim->v,na-1);
- COPY(v+na-1,r->dim->v+nc,nb-nc);
- out=new GridOutlet(this,0,new Dim(nd,v),r->nt);
- if (nc>0) in->set_factor(nc);
-} GRID_FLOW {
- int na = in->dim->n;
- int nc = in->dim->get(na-1);
- int size = r->dim->prod(nc);
- assert((n % nc) == 0);
- int nd = n/nc;
- STACK_ARRAY(T,w,n);
- Pt<T> v=w;
- if (sizeof(T)==1 && nc==1 && r->dim->v[0]<=256) {
- // bug? shouldn't modulo be done here?
- v=data;
- } else {
- COPY(v,data,n);
- for (int k=0,i=0; i<nc; i++) for (int j=0; j<n; j+=nc) v[k++] = data[i+j];
- compute_indices(v,nc,nd);
- }
-#define FOO(type) { \
- Pt<type> p = (Pt<type>)*r; \
- if (size<=16) { \
- Pt<type> foo = ARRAY_NEW(type,nd*size); \
- int i=0; \
- switch (size) { \
- case 1: for (; i<nd&-4; i+=4, foo+=4) { \
- foo[0] = p[v[i+0]]; \
- foo[1] = p[v[i+1]]; \
- foo[2] = p[v[i+2]]; \
- foo[3] = p[v[i+3]]; \
- } break; \
- case 2: for (; i<nd; i++, foo+=2) SCOPY(foo,p+2*v[i],2); break; \
- case 3: for (; i<nd; i++, foo+=3) SCOPY(foo,p+3*v[i],3); break; \
- case 4: for (; i<nd; i++, foo+=4) SCOPY(foo,p+4*v[i],4); break; \
- default:; }; \
- for (; i<nd; i++, foo+=size) COPY(foo,p+size*v[i],size); \
- out->give(size*nd,foo-size*nd); \
- } else { \
- for (int i=0; i<nd; i++) out->send(size,p+size*v[i]); \
- } \
-}
- TYPESWITCH(r->nt,FOO,)
-#undef FOO
-} GRID_FINISH {
- if (in->dim->prod()==0) {
- int n = in->dim->prod(0,-2);
- int size = r->dim->prod();
-#define FOO(T) while (n--) out->send(size,(Pt<T>)*r);
- TYPESWITCH(r->nt,FOO,)
-#undef FOO
- }
-} GRID_END
-
-GRID_INLET(GridStore,1) {
- NumberTypeE nt = NumberTypeE_type_of(*data);
- if (!put_at) { // reassign
- if (in[0].dim)
- r.next = new Grid(in->dim,nt);
- else
- r = new Grid(in->dim,nt);
- return;
- }
- // put_at ( ... )
- //!@#$ should check types. if (r->nt!=in->nt) RAISE("shoo");
- int nn=r->dim->n, na=put_at->dim->v[0], nb=in->dim->n;
- STACK_ARRAY(int32,sizeb,nn);
- for (int i=0; i<nn; i++) { fromb[i]=0; sizeb[i]=1; }
- COPY(Pt<int32>(wdex,nn) ,(Pt<int32>)*put_at ,put_at->dim->prod());
- COPY(Pt<int32>(fromb,nn)+nn-na,(Pt<int32>)*put_at ,na);
- COPY(Pt<int32>(sizeb,nn)+nn-nb,(Pt<int32>)in->dim->v,nb);
- for (int i=0; i<nn; i++) to2[i] = fromb[i]+sizeb[i];
- d=0;
- // find out when we can skip computing indices
- //!@#$ should actually also stop before blowing up packet size
- lsd=nn;
- while (lsd>=nn-in->dim->n) {
- lsd--;
- int cs = in->dim->prod(lsd-nn+in->dim->n);
- if (cs>GridOutlet::MAX_PACKET_SIZE || fromb[lsd]!=0 || sizeb[lsd]!=r->dim->v[lsd]) break;
- }
- lsd++;
- int cs = in->dim->prod(lsd-nn+in->dim->n);
- in->set_factor(cs);
-} GRID_FLOW {
- if (!put_at) { // reassign
- COPY(((Pt<T>)*(r.next ? r.next.p : &*r.p))+in->dex, data, n);
- return;
- }
- // put_at ( ... )
- int nn=r->dim->n;
- int cs = in->factor(); // chunksize
- STACK_ARRAY(int32,v,lsd);
- Pt<int32> x = Pt<int32>(wdex,nn);
- while (n) {
- // here d is the dim# to reset; d=n for none
- for(;d<lsd;d++) x[d]=fromb[d];
- COPY(v,x,lsd);
- compute_indices(v,lsd,1);
- op->zip(cs,(Pt<T>)*r+v[0]*cs,data);
- data+=cs;
- n-=cs;
- // find next set of indices; here d is the dim# to increment
- for(;;) {
- d--;
- if (d<0) goto end;
- x[d]++;
- if (x[d]<to2[d]) break;
- }
- end:; // why here ??? or why at all?
- d++;
- }
- //end:; // why not here ???
-} GRID_END
-\def void _0_op(Numop *op) { this->op=op; }
-\def void _0_bang () { rb_funcall(rself,SI(_0_list),3,INT2NUM(0),SYM(#),INT2NUM(0)); }
-\def void _1_reassign () { put_at=0; }
-\def void _1_put_at (Grid *index) { put_at=index; }
-\def void initialize (Grid *r) {
- rb_call_super(argc,argv);
- this->r = r?r:new Grid(new Dim(),int32_e,true);
- op = op_put;
-}
-\classinfo { IEVAL(rself,"install '#store',2,1"); }
-\end class GridStore
-
-//****************************************************************
-//{ Dim[*As]<T> -> Dim[*As]<T> } or
-//{ Dim[*As]<T>,Dim[*Bs]<T> -> Dim[*As]<T> }
-\class GridOp < GridObject
-struct GridOp : GridObject {
- \attr Numop *op;
- PtrGrid r;
- \decl void initialize(Numop *op, Grid *r=0);
- \grin 0
- \grin 1
- \decl void _0_op(Numop *op);
-};
-
-GRID_INLET(GridOp,0) {
- snap_backstore(r);
- SAME_TYPE(in,r);
- out=new GridOutlet(this,0,in->dim,in->nt);
- in->set_mode(6);
-} GRID_ALLOC {
- //out->ask(in->allocn,(Pt<T> &)in->alloc,in->allocfactor,in->allocmin,in->allocmax);
-} GRID_FLOW {
- Pt<T> rdata = (Pt<T>)*r;
- int loop = r->dim->prod();
- if (sizeof(T)==8) {
- fprintf(stderr,"1: data=%p rdata=%p\n",data.p,rdata.p);
- WATCH(n,data);
- }
- if (loop>1) {
- if (in->dex+n <= loop) {
- op->zip(n,data,rdata+in->dex);
- } else {
- // !@#$ should prebuild and reuse this array when "loop" is small
- STACK_ARRAY(T,data2,n);
- int ii = mod(in->dex,loop);
- int m = min(loop-ii,n);
- COPY(data2,rdata+ii,m);
- int nn = m+((n-m)/loop)*loop;
- for (int i=m; i<nn; i+=loop) COPY(data2+i,rdata,loop);
- if (n>nn) COPY(data2+nn,rdata,n-nn);
- if (sizeof(T)==8) {
- fprintf(stderr,"2: data=%p data2=%p\n",data.p,data2.p);
- WATCH(n,data); WATCH(n,data2);
- }
- op->zip(n,data,data2);
- if (sizeof(T)==8) {WATCH(n,data); WATCH(n,data2);}
- }
- } else {
- op->map(n,data,*rdata);
- }
- out->give(n,data);
-} GRID_END
-
-GRID_INPUT2(GridOp,1,r) {} GRID_END
-\def void _0_op(Numop *op) { this->op=op; }
-
-\def void initialize(Numop *op, Grid *r=0) {
- rb_call_super(argc,argv);
- this->op=op;
- this->r = r?r:new Grid(new Dim(),int32_e,true);
-}
-
-\classinfo { IEVAL(rself,"install '#',2,1"); }
-\end class GridOp
-
-//****************************************************************
-\class GridFold < GridObject
-struct GridFold : GridObject {
- \attr Numop *op;
- \attr PtrGrid seed;
- \decl void initialize (Numop *op);
- \decl void _0_op (Numop *op);
- \decl void _0_seed (Grid *seed);
- \grin 0
-};
-
-GRID_INLET(GridFold,0) {
- //{ Dim[*As,B,*Cs]<T>,Dim[*Cs]<T> -> Dim[*As,*Cs]<T> }
- if (seed) SAME_TYPE(in,seed);
- int an = in->dim->n;
- int bn = seed?seed->dim->n:0;
- if (an<=bn) RAISE("minimum 1 more dimension than the seed (%d vs %d)",an,bn);
- STACK_ARRAY(int32,v,an-1);
- int yi = an-bn-1;
- COPY(v,in->dim->v,yi);
- COPY(v+yi,in->dim->v+an-bn,bn);
- if (seed) SAME_DIM(an-(yi+1),in->dim,(yi+1),seed->dim,0);
- out=new GridOutlet(this,0,new Dim(an-1,v),in->nt);
- int k = seed ? seed->dim->prod() : 1;
- in->set_factor(in->dim->get(yi)*k);
-} GRID_FLOW {
- int an = in->dim->n;
- int bn = seed?seed->dim->n:0;
- int yn = in->dim->v[an-bn-1];
- int zn = in->dim->prod(an-bn);
- STACK_ARRAY(T,buf,n/yn);
- int nn=n;
- int yzn=yn*zn;
- for (int i=0; n; i+=zn, data+=yzn, n-=yzn) {
- if (seed) COPY(buf+i,((Pt<T>)*seed),zn);
- else CLEAR(buf+i,zn);
- op->fold(zn,yn,buf+i,data);
- }
- out->send(nn/yn,buf);
-} GRID_END
-
-\def void _0_op (Numop *op ) { this->op =op; }
-\def void _0_seed (Grid *seed) { this->seed=seed; }
-\def void initialize (Numop *op) { rb_call_super(argc,argv); this->op=op; }
-\classinfo { IEVAL(rself,"install '#fold',1,1"); }
-\end class GridFold
-
-\class GridScan < GridObject
-struct GridScan : GridObject {
- \attr Numop *op;
- \attr PtrGrid seed;
- \decl void initialize (Numop *op);
- \decl void _0_op (Numop *op);
- \decl void _0_seed (Grid *seed);
- \grin 0
-};
-
-GRID_INLET(GridScan,0) {
- //{ Dim[*As,B,*Cs]<T>,Dim[*Cs]<T> -> Dim[*As,B,*Cs]<T> }
- if (seed) SAME_TYPE(in,seed);
- int an = in->dim->n;
- int bn = seed?seed->dim->n:0;
- if (an<=bn) RAISE("minimum 1 more dimension than the right hand");
- if (seed) SAME_DIM(bn,in->dim,an-bn,seed->dim,0);
- out=new GridOutlet(this,0,in->dim,in->nt);
- in->set_factor(in->dim->prod(an-bn-1));
-} GRID_FLOW {
- int an = in->dim->n;
- int bn = seed?seed->dim->n:0;
- int yn = in->dim->v[an-bn-1];
- int zn = in->dim->prod(an-bn);
- int factor = in->factor();
- STACK_ARRAY(T,buf,n);
- COPY(buf,data,n);
- if (seed) {
- for (int i=0; i<n; i+=factor) op->scan(zn,yn,(Pt<T>)*seed,buf+i);
- } else {
- STACK_ARRAY(T,seed,zn);
- CLEAR(seed,zn);
- for (int i=0; i<n; i+=factor) op->scan(zn,yn,seed,buf+i);
- }
- out->send(n,buf);
-} GRID_END
-
-\def void _0_op (Numop *op ) { this->op =op; }
-\def void _0_seed (Grid *seed) { this->seed=seed; }
-\def void initialize (Numop *op) { rb_call_super(argc,argv); this->op = op; }
-\classinfo { IEVAL(rself,"install '#scan',1,1"); }
-\end class GridScan
-
-//****************************************************************
-//{ Dim[*As,C]<T>,Dim[C,*Bs]<T> -> Dim[*As,*Bs]<T> }
-\class GridInner < GridObject
-struct GridInner : GridObject {
- \attr Numop *op_para;
- \attr Numop *op_fold;
- \attr PtrGrid seed;
- PtrGrid r;
- PtrGrid r2;
- GridInner() {}
- \decl void initialize (Grid *r=0);
- \decl void _0_op (Numop *op);
- \decl void _0_fold (Numop *op);
- \decl void _0_seed (Grid *seed);
- \grin 0
- \grin 1
-};
-
-template <class T> void inner_child_a (Pt<T> buf, Pt<T> data, int rrows, int rcols, int chunk) {
- Pt<T> bt = buf, dt = data;
- for (int j=0; j<chunk; j++, bt+=rcols, dt+=rrows) op_put->map(rcols,bt,*dt);
-}
-template <class T, int rcols> void inner_child_b (Pt<T> buf, Pt<T> data, int rrows, int chunk) {
- Pt<T> bt = buf, dt = data;
- for (int j=0; j<chunk; j++, bt+=rcols, dt+=rrows) {
- for (int k=0; k<rcols; k++) bt[k] = *dt;
- }
-}
-GRID_INLET(GridInner,0) {
- SAME_TYPE(in,r);
- SAME_TYPE(in,seed);
- P<Dim> a = in->dim;
- P<Dim> b = r->dim;
- if (a->n<1) RAISE("a: minimum 1 dimension");
- if (b->n<1) RAISE("b: minimum 1 dimension");
- if (seed->dim->n != 0) RAISE("seed must be a scalar");
- int a_last = a->get(a->n-1);
- int n = a->n+b->n-2;
- SAME_DIM(1,a,a->n-1,b,0);
- STACK_ARRAY(int32,v,n);
- COPY(v,a->v,a->n-1);
- COPY(v+a->n-1,b->v+1,b->n-1);
- out=new GridOutlet(this,0,new Dim(n,v),in->nt);
- in->set_factor(a_last);
-
- int rrows = in->factor();
- int rsize = r->dim->prod();
- int rcols = rsize/rrows;
- Pt<T> rdata = (Pt<T>)*r;
- int chunk = GridOutlet::MAX_PACKET_SIZE/rsize;
- r2=new Grid(new Dim(chunk*rsize),r->nt);
- Pt<T> buf3 = (Pt<T>)*r2;
- for (int i=0; i<rrows; i++)
- for (int j=0; j<chunk; j++)
- COPY(buf3+(j+i*chunk)*rcols,rdata+i*rcols,rcols);
-} GRID_FLOW {
- int rrows = in->factor();
- int rsize = r->dim->prod();
- int rcols = rsize/rrows;
- int chunk = GridOutlet::MAX_PACKET_SIZE/rsize;
- STACK_ARRAY(T,buf ,chunk*rcols);
- STACK_ARRAY(T,buf2,chunk*rcols);
- int off = chunk;
- while (n) {
- if (chunk*rrows>n) chunk=n/rrows;
- op_put->map(chunk*rcols,buf2,*(T *)*seed);
- for (int i=0; i<rrows; i++) {
- switch (rcols) {
- case 1: inner_child_b<T,1>(buf,data+i,rrows,chunk); break;
- case 2: inner_child_b<T,2>(buf,data+i,rrows,chunk); break;
- case 3: inner_child_b<T,3>(buf,data+i,rrows,chunk); break;
- case 4: inner_child_b<T,4>(buf,data+i,rrows,chunk); break;
- default: inner_child_a(buf,data+i,rrows,rcols,chunk);
- }
- op_para->zip(chunk*rcols,buf,(Pt<T>)*r2+i*off*rcols);
- op_fold->zip(chunk*rcols,buf2,buf);
- }
- out->send(chunk*rcols,buf2);
- n-=chunk*rrows;
- data+=chunk*rrows;
- }
-} GRID_FINISH {
- r2=0;
-} GRID_END
-
-GRID_INPUT(GridInner,1,r) {} GRID_END
-
-\def void initialize (Grid *r) {
- rb_call_super(argc,argv);
- this->op_para = op_mul;
- this->op_fold = op_add;
- this->seed = new Grid(new Dim(),int32_e,true);
- this->r = r ? r : new Grid(new Dim(),int32_e,true);
-}
-
-\def void _0_op (Numop *op ) { this->op_para=op; }
-\def void _0_fold (Numop *op ) { this->op_fold=op; }
-\def void _0_seed (Grid *seed) { this->seed=seed; }
-\classinfo { IEVAL(rself,"install '#inner',2,1"); }
-\end class GridInner
-
-/* **************************************************************** */
-/*{ Dim[*As]<T>,Dim[*Bs]<T> -> Dim[*As,*Bs]<T> }*/
-\class GridOuter < GridObject
-struct GridOuter : GridObject {
- \attr Numop *op;
- PtrGrid r;
- \decl void initialize (Numop *op, Grid *r=0);
- \grin 0
- \grin 1
-};
-
-GRID_INLET(GridOuter,0) {
- SAME_TYPE(in,r);
- P<Dim> a = in->dim;
- P<Dim> b = r->dim;
- int n = a->n+b->n;
- STACK_ARRAY(int32,v,n);
- COPY(v,a->v,a->n);
- COPY(v+a->n,b->v,b->n);
- out=new GridOutlet(this,0,new Dim(n,v),in->nt);
-} GRID_FLOW {
- int b_prod = r->dim->prod();
- if (b_prod > 4) {
- STACK_ARRAY(T,buf,b_prod);
- while (n) {
- for (int j=0; j<b_prod; j++) buf[j] = *data;
- op->zip(b_prod,buf,(Pt<T>)*r);
- out->send(b_prod,buf);
- data++; n--;
- }
- return;
- }
- n*=b_prod;
- Pt<T> buf = ARRAY_NEW(T,n);
- STACK_ARRAY(T,buf2,b_prod*64);
- for (int i=0; i<64; i++) COPY(buf2+i*b_prod,(Pt<T>)*r,b_prod);
- switch (b_prod) {
- #define Z buf[k++]=data[i]
- case 1: for (int i=0,k=0; k<n; i++) {Z;} break;
- case 2: for (int i=0,k=0; k<n; i++) {Z;Z;} break;
- case 3: for (int i=0,k=0; k<n; i++) {Z;Z;Z;} break;
- case 4: for (int i=0,k=0; k<n; i++) {Z;Z;Z;Z;} break;
- default:for (int i=0,k=0; k<n; i++) for (int j=0; j<b_prod; j++, k++) Z;
- }
- #undef Z
- int ch=64*b_prod;
- int nn=(n/ch)*ch;
- for (int j=0; j<nn; j+=ch) op->zip(ch,buf+j,buf2);
- op->zip(n-nn,buf+nn,buf2);
- out->give(n,buf);
-} GRID_END
-
-GRID_INPUT(GridOuter,1,r) {} GRID_END
-
-\def void initialize (Numop *op, Grid *r) {
- rb_call_super(argc,argv);
- this->op = op;
- this->r = r ? r : new Grid(new Dim(),int32_e,true);
-}
-
-\classinfo { IEVAL(rself,"install '#outer',2,1"); }
-\end class GridOuter
-
-//****************************************************************
-//{ Dim[]<T>,Dim[]<T>,Dim[]<T> -> Dim[A]<T> } or
-//{ Dim[B]<T>,Dim[B]<T>,Dim[B]<T> -> Dim[*As,B]<T> }
-\class GridFor < GridObject
-struct GridFor : GridObject {
- \attr PtrGrid from;
- \attr PtrGrid to;
- \attr PtrGrid step;
- GridFor () {
- from.constrain(expect_max_one_dim);
- to .constrain(expect_max_one_dim);
- step.constrain(expect_max_one_dim);
- }
- \decl void initialize (Grid *from, Grid *to, Grid *step);
- \decl void _0_set (Grid *r=0);
- \decl void _0_bang ();
- \grin 0 int
- \grin 1 int
- \grin 2 int
- template <class T> void trigger (T bogus);
-};
-
-\def void initialize (Grid *from, Grid *to, Grid *step) {
- rb_call_super(argc,argv);
- this->from=from;
- this->to =to;
- this->step=step;
-}
-
-template <class T>
-void GridFor::trigger (T bogus) {
- int n = from->dim->prod();
- int32 nn[n+1];
- STACK_ARRAY(T,x,64*n);
- Pt<T> fromb = (Pt<T>)*from;
- Pt<T> tob = (Pt<T>)*to ;
- Pt<T> stepb = (Pt<T>)*step;
- STACK_ARRAY(T,to2,n);
-
- for (int i=step->dim->prod()-1; i>=0; i--)
- if (!stepb[i]) RAISE("step must not contain zeroes");
- for (int i=0; i<n; i++) {
- nn[i] = (tob[i] - fromb[i] + stepb[i] - cmp(stepb[i],(T)0)) / stepb[i];
- if (nn[i]<0) nn[i]=0;
- to2[i] = fromb[i]+stepb[i]*nn[i];
- }
- P<Dim> d;
- if (from->dim->n==0) { d = new Dim(*nn); }
- else { nn[n]=n; d = new Dim(n+1,nn); }
- int total = d->prod();
- out=new GridOutlet(this,0,d,from->nt);
- if (total==0) return;
- int k=0;
- for(int d=0;;d++) {
- // here d is the dim# to reset; d=n for none
- for(;d<n;d++) x[k+d]=fromb[d];
- k+=n;
- if (k==64*n) {out->send(k,x); k=0; COPY(x,x+63*n,n);}
- else { COPY(x+k,x+k-n,n);}
- d--;
- // here d is the dim# to increment
- for(;;d--) {
- if (d<0) goto end;
- x[k+d]+=stepb[d];
- if (x[k+d]!=to2[d]) break;
- }
- }
- end: if (k) out->send(k,x);
-}
-
-\def void _0_bang () {
- SAME_TYPE(from,to);
- SAME_TYPE(from,step);
- if (!from->dim->equal(to->dim) || !to->dim->equal(step->dim))
- RAISE("dimension mismatch");
-#define FOO(T) trigger((T)0);
- TYPESWITCH_JUSTINT(from->nt,FOO,);
-#undef FOO
-}
-
-\def void _0_set (Grid *r) { from=new Grid(argv[0]); }
-GRID_INPUT(GridFor,2,step) {} GRID_END
-GRID_INPUT(GridFor,1,to) {} GRID_END
-GRID_INPUT(GridFor,0,from) {_0_bang(0,0);} GRID_END
-\classinfo { IEVAL(rself,"install '#for',3,1"); }
-\end class GridFor
-
-//****************************************************************
-\class GridFinished < GridObject
-struct GridFinished : GridObject {
- \grin 0
-};
-GRID_INLET(GridFinished,0) {
- in->set_mode(0);
-} GRID_FINISH {
- Ruby a[] = { INT2NUM(0), bsym._bang };
- send_out(COUNT(a),a);
-} GRID_END
-\classinfo { IEVAL(rself,"install '#finished',1,1"); }
-\end class GridFinished
-
-\class GridDim < GridObject
-struct GridDim : GridObject {
- \grin 0
-};
-GRID_INLET(GridDim,0) {
- GridOutlet out(this,0,new Dim(in->dim->n));
- out.send(in->dim->n,Pt<int32>(in->dim->v,in->dim->n));
- in->set_mode(0);
-} GRID_END
-\classinfo { IEVAL(rself,"install '#dim',1,1"); }
-\end class GridDim
-
-\class GridType < GridObject
-struct GridType : GridObject {
- \grin 0
-};
-GRID_INLET(GridType,0) {
- Ruby a[] = { INT2NUM(0), SYM(symbol), number_type_table[in->nt].sym };
- send_out(COUNT(a),a);
- in->set_mode(0);
-} GRID_END
-\classinfo { IEVAL(rself,"install '#type',1,1"); }
-\end class GridType
-
-//****************************************************************
-//{ Dim[*As]<T>,Dim[B] -> Dim[*Cs]<T> }
-\class GridRedim < GridObject
-struct GridRedim : GridObject {
- \attr P<Dim> dim;
- PtrGrid dim_grid;
- PtrGrid temp; // temp->dim is not of the same shape as dim
- GridRedim() { dim_grid.constrain(expect_dim_dim_list); }
- ~GridRedim() {}
- \decl void initialize (Grid *d);
- \grin 0
- \grin 1 int32
-};
-
-GRID_INLET(GridRedim,0) {
- int a = in->dim->prod(), b = dim->prod();
- if (a<b) temp=new Grid(new Dim(a),in->nt);
- out=new GridOutlet(this,0,dim,in->nt);
-} GRID_FLOW {
- int i = in->dex;
- if (!temp) {
- int b = dim->prod();
- int n2 = min(n,b-i);
- if (n2>0) out->send(n2,data);
- // discard other values if any
- } else {
- int a = in->dim->prod();
- int n2 = min(n,a-i);
- COPY((Pt<T>)*temp+i,data,n2);
- if (n2>0) out->send(n2,data);
- }
-} GRID_FINISH {
- if (!!temp) {
- int a = in->dim->prod(), b = dim->prod();
- if (a) {
- for (int i=a; i<b; i+=a) out->send(min(a,b-i),(Pt<T>)*temp);
- } else {
- STACK_ARRAY(T,foo,1);
- foo[0]=0;
- for (int i=0; i<b; i++) out->send(1,foo);
- }
- }
- temp=0;
-} GRID_END
-
-GRID_INPUT(GridRedim,1,dim_grid) { dim = dim_grid->to_dim(); } GRID_END
-
-\def void initialize (Grid *d) {
- rb_call_super(argc,argv);
- dim_grid=d;
- dim = dim_grid->to_dim();
-}
-
-\classinfo { IEVAL(rself,"install '#redim',2,1"); }
-\end class GridRedim
-
-//****************************************************************
-\class GridJoin < GridObject
-struct GridJoin : GridObject {
- \attr int which_dim;
- PtrGrid r;
- \grin 0
- \grin 1
- \decl void initialize (int which_dim=-1, Grid *r=0);
-};
-
-GRID_INLET(GridJoin,0) {
- NOTEMPTY(r);
- SAME_TYPE(in,r);
- P<Dim> d = in->dim;
- if (d->n != r->dim->n) RAISE("wrong number of dimensions");
- int w = which_dim;
- if (w<0) w+=d->n;
- if (w<0 || w>=d->n)
- RAISE("can't join on dim number %d on %d-dimensional grids",
- which_dim,d->n);
- STACK_ARRAY(int32,v,d->n);
- for (int i=0; i<d->n; i++) {
- v[i] = d->get(i);
- if (i==w) {
- v[i]+=r->dim->v[i];
- } else {
- if (v[i]!=r->dim->v[i]) RAISE("dimensions mismatch: dim #%i, left is %d, right is %d",i,v[i],r->dim->v[i]);
- }
- }
- out=new GridOutlet(this,0,new Dim(d->n,v),in->nt);
- if (d->prod(w)) in->set_factor(d->prod(w));
-} GRID_FLOW {
- int w = which_dim;
- if (w<0) w+=in->dim->n;
- int a = in->factor();
- int b = r->dim->prod(w);
- Pt<T> data2 = (Pt<T>)*r + in->dex*b/a;
- if (a==3 && b==1) {
- int m = n+n*b/a;
- STACK_ARRAY(T,data3,m);
- Pt<T> data4 = data3;
- while (n) {
- SCOPY(data4,data,3); SCOPY(data4+3,data2,1);
- n-=3; data+=3; data2+=1; data4+=4;
- }
- out->send(m,data3);
- } else if (a+b<=16) {
- int m = n+n*b/a;
- STACK_ARRAY(T,data3,m);
- int i=0;
- while (n) {
- COPY(data3+i,data,a); data+=a; i+=a; n-=a;
- COPY(data3+i,data2,b); data2+=b; i+=b;
- }
- out->send(m,data3);
- } else {
- while (n) {
- out->send(a,data);
- out->send(b,data2);
- data+=a; data2+=b; n-=a;
- }
- }
-} GRID_FINISH {
- if (in->dim->prod()==0) out->send(r->dim->prod(),(Pt<T>)*r);
-} GRID_END
-
-GRID_INPUT(GridJoin,1,r) {} GRID_END
-
-\def void initialize (int which_dim, Grid *r) {
- rb_call_super(argc,argv);
- this->which_dim = which_dim;
- if (r) this->r=r;
-}
-
-\classinfo { IEVAL(rself,"install '@join',2,1"); }
-\end class GridJoin
-
-//****************************************************************
-\class GridGrade < GridObject
-struct GridGrade : GridObject {
- \grin 0
-};
-
-template <class T> struct GradeFunction {
- static int comparator (const void *a, const void *b) {
- return **(T**)a - **(T**)b;}};
-#define FOO(S) \
-template <> struct GradeFunction<S> { \
- static int comparator (const void *a, const void *b) { \
- S x = **(S**)a - **(S**)b; \
- return x<0 ? -1 : x>0;}};
-FOO(int64)
-FOO(float32)
-FOO(float64)
-#undef FOO
-
-GRID_INLET(GridGrade,0) {
- out=new GridOutlet(this,0,in->dim,in->nt);
- in->set_factor(in->dim->get(in->dim->n-1));
-} GRID_FLOW {
- int m = in->factor();
- STACK_ARRAY(T*,foo,m);
- STACK_ARRAY(T,bar,m);
- for (; n; n-=m,data+=m) {
- for (int i=0; i<m; i++) foo[i] = &data[i];
- qsort(foo,m,sizeof(T),GradeFunction<T>::comparator);
- for (int i=0; i<m; i++) bar[i] = foo[i]-(T *)data;
- out->send(m,bar);
- }
-} GRID_END
-
-\classinfo { IEVAL(rself,"install '#grade',1,1"); }
-\end class GridGrade
-
-//****************************************************************
-//\class GridMedian < GridObject
-//****************************************************************
-
-\class GridTranspose < GridObject
-struct GridTranspose : GridObject {
- \attr int dim1;
- \attr int dim2;
- int d1,d2,na,nb,nc,nd; // temporaries
- \decl void initialize (int dim1=0, int dim2=1);
- \decl void _1_float (int dim1);
- \decl void _2_float (int dim2);
- \grin 0
-};
-
-\def void _1_float (int dim1) { this->dim1=dim1; }
-\def void _2_float (int dim2) { this->dim2=dim2; }
-
-GRID_INLET(GridTranspose,0) {
- STACK_ARRAY(int32,v,in->dim->n);
- COPY(v,in->dim->v,in->dim->n);
- d1=dim1; d2=dim2;
- if (d1<0) d1+=in->dim->n;
- if (d2<0) d2+=in->dim->n;
- if (d1>=in->dim->n || d2>=in->dim->n || d1<0 || d2<0)
- RAISE("would swap dimensions %d and %d but this grid has only %d dimensions",
- dim1,dim2,in->dim->n);
- memswap(v+d1,v+d2,1);
- if (d1==d2) {
- out=new GridOutlet(this,0,new Dim(in->dim->n,v), in->nt);
- } else {
- nd = in->dim->prod(1+max(d1,d2));
- nc = in->dim->v[max(d1,d2)];
- nb = in->dim->prod(1+min(d1,d2))/nc/nd;
- na = in->dim->v[min(d1,d2)];
- out=new GridOutlet(this,0,new Dim(in->dim->n,v), in->nt);
- in->set_factor(na*nb*nc*nd);
- }
- // Turns a Grid[*,na,*nb,nc,*nd] into a Grid[*,nc,*nb,na,*nd].
-} GRID_FLOW {
- STACK_ARRAY(T,res,na*nb*nc*nd);
- if (dim1==dim2) { out->send(n,data); return; }
- int prod = na*nb*nc*nd;
- for (; n; n-=prod, data+=prod) {
- for (int a=0; a<na; a++)
- for (int b=0; b<nb; b++)
- for (int c=0; c<nc; c++)
- COPY(res +((c*nb+b)*na+a)*nd,
- data+((a*nb+b)*nc+c)*nd,nd);
- out->send(na*nb*nc*nd,res);
- }
-} GRID_END
-
-\def void initialize (int dim1=0, int dim2=1) {
- rb_call_super(argc,argv);
- this->dim1 = dim1;
- this->dim2 = dim2;
-}
-
-\classinfo { IEVAL(rself,"install '#transpose',3,1"); }
-\end class GridTranspose
-
-//****************************************************************
-\class GridReverse < GridObject
-struct GridReverse : GridObject {
- \attr int dim1; // dimension to act upon
- int d; // temporaries
- \decl void initialize (int dim1=0);
- \decl void _1_float (int dim1);
- \grin 0
-};
-
-\def void _1_float (int dim1) { this->dim1=dim1; }
-
-GRID_INLET(GridReverse,0) {
- d=dim1;
- if (d<0) d+=in->dim->n;
- if (d>=in->dim->n || d<0)
- RAISE("would reverse dimension %d but this grid has only %d dimensions",
- dim1,in->dim->n);
- out=new GridOutlet(this,0,new Dim(in->dim->n,in->dim->v), in->nt);
- in->set_factor(in->dim->prod(d));
-} GRID_FLOW {
- int f1=in->factor(), f2=in->dim->prod(d+1);
- while (n) {
- int hf1=f1/2;
- Pt<T> data2 = data+f1-f2;
- for (int i=0; i<hf1; i+=f2) memswap(data+i,data2-i,f2);
- out->send(f1,data);
- data+=f1; n-=f1;
- }
-} GRID_END
-
-\def void initialize (int dim1=0) {
- rb_call_super(argc,argv);
- this->dim1 = dim1;
-}
-
-\classinfo { IEVAL(rself,"install '#reverse',2,1"); }
-\end class GridReverse
-
-//****************************************************************
-\class GridCentroid < GridObject
-struct GridCentroid : GridObject {
- \decl void initialize ();
- \grin 0 int
- int sumx,sumy,sum,y; // temporaries
-};
-
-GRID_INLET(GridCentroid,0) {
- if (in->dim->n != 3) RAISE("expecting 3 dims");
- if (in->dim->v[2] != 1) RAISE("expecting 1 channel");
- in->set_factor(in->dim->prod(1));
- out=new GridOutlet(this,0,new Dim(2), in->nt);
- sumx=0; sumy=0; sum=0; y=0;
-} GRID_FLOW {
- int sx = in->dim->v[1];
- while (n) {
- for (int x=0; x<sx; x++) {
- sumx+=x*data[x];
- sumy+=y*data[x];
- sum += data[x];
- }
- n-=sx;
- data+=sx;
- y++;
- }
-} GRID_FINISH {
- STACK_ARRAY(int32,blah,2);
- blah[0] = sum ? sumy/sum : 0;
- blah[1] = sum ? sumx/sum : 0;
- out->send(2,blah);
- rb_funcall(rself,SI(send_out),2,INT2NUM(1),INT2NUM(blah[0]));
- rb_funcall(rself,SI(send_out),2,INT2NUM(2),INT2NUM(blah[1]));
-} GRID_END
-
-\def void initialize () {
- rb_call_super(argc,argv);
-}
-
-\classinfo { IEVAL(rself,"install '#centroid',1,3"); }
-\end class GridCentroid
-
-//****************************************************************
-\class GridPerspective < GridObject
-struct GridPerspective : GridObject {
- \attr int32 z;
- \grin 0
- \decl void initialize (int32 z=256);
-};
-
-GRID_INLET(GridPerspective,0) {
- int n = in->dim->n;
- STACK_ARRAY(int32,v,n);
- COPY(v,in->dim->v,n);
- v[n-1]--;
- in->set_factor(in->dim->get(in->dim->n-1));
- out=new GridOutlet(this,0,new Dim(n,v),in->nt);
-} GRID_FLOW {
- int m = in->factor();
- for (; n; n-=m,data+=m) {
- op_mul->map(m-1,data,(T)z);
- op_div->map(m-1,data,data[m-1]);
- out->send(m-1,data);
- }
-} GRID_END
-
-\def void initialize (int32 z) {rb_call_super(argc,argv); this->z=z; }
-\classinfo { IEVAL(rself,"install '#perspective',1,1"); }
-\end class GridPerspective
-
-static Numop *OP(Ruby x) { return FIX2PTR(Numop,rb_hash_aref(op_dict,x)); }
-
-void startup_flow_objects () {
- op_add = OP(SYM(+));
- op_sub = OP(SYM(-));
- op_mul = OP(SYM(*));
- op_shl = OP(SYM(<<));
- op_mod = OP(SYM(%));
- op_and = OP(SYM(&));
- op_div = OP(SYM(/));
- op_put = OP(SYM(put));
- \startall
-}