: \$Id: vecst.mod,v 1.66 2001/05/15 14:54:49 billl Exp \$ :* COMMENT COMMENT randwd randomly chooses n bits to set to 1 hamming v.hamming(v1) is hamming distance between 2 vecs flipbits v.flipbits(scratch,num) flips num rand chosen bits flipbalbits v.flipbalbits(scratch,num) balanced flipping vpr v.vpr prints out vector as 1 (x[i]>0) or 0 (x[i]<=0) thresh turns analog vec to BVBASE,1 vec separating at thresh (scalar or vec) triplet return location of a triplet in a vector onoff turns state vec on or off depending on values in other vecs bpeval used by backprop: vo=outp*(1-outp)*del w like where but sets chosen nums // modified from .wh in 1.23 xing determines where vector crosses in a positive direction (destructive) xzero count 0 crossings by putting in a 1 or -1 indset(ind,val) sets spots indexed by ind to val ismono([arg]) checks if a vector is monotonically increaseing (-1 decreasing) count(num) count the number of nums in vec scr.fewind(ind,vec1,vec2[,vec3,vec4]) // uses ind as index into other vecs scr.nind(ind,vec1,vec2[,vec3,vec4]) // uses ind to elim elems in other vecs vdest.intrp(flag) // interpolate numbers replacing numbers given as flag v.insct(v1,v2) // return v1 intersect v2 v.cvlv(v1,v2) // convolve v1 with v2 fac not vec related - returns factorial logfac not vec related - returns log factorial v2d copy into a double block d2v copy out of a double block NB: same as vec.from_double() frd perform an fread ENDCOMMENT NEURON { SUFFIX nothing GLOBAL BVBASE : bit vector base number (typically 0 or -1) } PARAMETER { BVBASE = -1. } VERBATIM #include #include #include /* contains MAXLONG */ #ifndef NRN_VERSION_GTEQ_8_2_0 extern double* hoc_pgetarg(); extern double hoc_call_func(Symbol*, int narg); extern FILE* hoc_obj_file_arg(int narg); extern void vector_resize(); extern int vector_instance_px(); extern void* vector_arg(); #endif /* some machines do not have drand48 and srand48 so use the implementation at the end of this file */ double my_drand48(); void my_srand48(long seedval); #undef drand48 #undef srand48 #define drand48 my_drand48 #define srand48 my_srand48 ENDVERBATIM :* v1.slope(num) does a linear regression to find the slope, assuming num=timestep of vector VERBATIM static double slope(void* vv) { int i, n; double *x, *y; double timestep, sigxy, sigx, sigy, sigx2; /* how to get the instance data */ n = vector_instance_px(vv, &y); if(ifarg(1)) { timestep = *getarg(1); } else { printf("You must supply a timestep\n"); return 0; } sigxy= sigx= sigy= sigx2=0; // initialize these x = (double *) malloc(sizeof(double)*n); for(i=0; i=nx) { hoc_execerror("# of flips exceeds (or ==) vector size", 0); } for (i=0; i < nx; i++) { x[i] = BVBASE; } for (i=0,jj=0; i < flip; i++) { /* flip these bits */ ii = (int) ((nx+1)*drand48()); if (x[ii]==BVBASE) { x[ii] = 1.; if (flag) { y[jj] = ii; jj++; } } else { i--; } } return flip; } ENDVERBATIM :* v1.hamming(v2[,v3]) compares v1 and v2 for matches, v3 gives diff vector VERBATIM static double hamming(void* vv) { int i, nx, ny, nz, prflag; double* x, *y, *z,sum; sum = 0.; nx = vector_instance_px(vv, &x); ny = vector_arg_px(1, &y); if (ifarg(2)) { /* write a diff vector to z */ prflag = 1; nz = vector_arg_px(2, &z); } else { prflag = 0; } if (nx!=ny || (prflag && nx!=nz)) { hoc_execerror("Vectors must be same size", 0); } for (i=0; i < nx; ++i) { if (x[i] != y[i]) { sum++; if (prflag) { z[i] = 1.; } } else if (prflag) { z[i] = 0.; } } return sum; } ENDVERBATIM :* v1.indset(ind,x[,y]) sets indexed values to x and other values to optional y VERBATIM static double indset(void* vv) { int i, nx, ny, nz, flag; double* x, *y, *z, val, val2 ; nx = vector_instance_px(vv, &x); ny = vector_arg_px(1, &y); if (hoc_is_object_arg(2)) { flag=1; nz = vector_arg_px(2, &z); if (ny!=nz) { hoc_execerror("v.indset: Vector sizes don't match.", 0); } } else { flag=0; val=*getarg(2); } if (ifarg(3)) { val2 = *getarg(3); for (i=0; i nx) { hoc_execerror("v.indset: Index exceeds vector size", 0); } if (flag) x[(int)y[i]]=z[i]; else x[(int)y[i]]=val; } return i; } ENDVERBATIM VERBATIM /* Maintain parallel int vector to avoid slowness of repeated casts */ static int scrsz=0; static int *scr; ENDVERBATIM :* tmp.fewind(ind,vec1,vec2[,vec3,vec4,...]) : picks out numbers from multiple vectors using index ind VERBATIM static double fewind(void* vv) { int i, j, nx, ni, nv[10], num; double *x, *ind, *vvo[10]; nx = vector_instance_px(vv, &x); for (i=0;ifarg(i);i++); if (i>9) hoc_execerror("ERR: fewind can only handle 9 vectors", 0); num = i-2; /* number of vectors to be picked apart */ for (i=0;iscrsz) { if (scrsz>0) { free(scr); scr=(int *)NULL; } scrsz=ni+10000; scr=(int *)ecalloc(scrsz, sizeof(int)); } for (i=0;insrc) { hoc_execerror("Filter bigger than source ", 0); } for (i=0;i0 && k9) hoc_execerror("ERR: nind can only handle 9 vectors", 0); num = i-2; /* number of vectors to be picked apart */ for (i=0;iscrsz) { if (scrsz>0) { free(scr); scr=(int *)NULL; } scrsz=ni+10000; scr=(int *)ecalloc(scrsz, sizeof(int)); } for (i=0,last=-1;i=nx) hoc_execerror("nind(): Index out of bounds", 0); if (scr[i]<=last) hoc_execerror("nind(): indices should mono increase", 0); last=scr[i]; } for (j=0;jBVBASE) { fprintf(f,"%d",1); } else { fprintf(f,"%d",0); } } fprintf(f,"\n"); } else { for (i=0; iBVBASE) { printf("%d",1); } else { printf("%d",0); } } printf("\n"); } return 1.; } ENDVERBATIM :* v1.thresh() threshold above and below thresh VERBATIM static double thresh(void* vv) { int i, nx, ny, cnt; double *x, *y, th; nx = vector_instance_px(vv, &x); cnt=0; if (hoc_is_object_arg(1)) { ny = vector_arg_px(1, &y); th=0; if (nx!=ny) { hoc_execerror("Vector sizes don't match in thresh.", 0); } for (i=0; i=y[i]) { x[i]= 1.; cnt++;} else { x[i]= BVBASE; } } } else { th = *getarg(1); for (i=0; i= th) { x[i]= 1.; cnt++;} else { x[i]= BVBASE; } } } return cnt; } ENDVERBATIM :* v1.triplet() return location of a triplet VERBATIM static double triplet(void* vv) { int i, nx; double *x, *y, a, b; nx = vector_instance_px(vv, &x); a = *getarg(1); b = *getarg(2); for (i=0; i0.) { st[i]=1.; continue; } /* cell must fire */ if (vol[i]>=thr[i] && obon[i]<= -refr[i]) { /* past refractory period */ st[i]=1.; obon[i]=dur[i]; num++; } else { st[i]= BVBASE; } } return num; } ENDVERBATIM :* vo.bpeval(outp,del) : service routine for back-prop: vo=outp*(1-outp)*del VERBATIM static double bpeval(void* vv) { int i, n, no, nd, flag=0; double add,div; double *vo, *outp, *del; n = vector_instance_px(vv, &vo); no = vector_arg_px(1, &outp); nd = vector_arg_px(2, &del); if (ifarg(3) && ifarg(4)) { add= *getarg(3); div= *getarg(4); flag=1;} if (n!=no||n!=nd) hoc_execerror("v.bpeval: vectors not all same size", 0); if (flag) { for (i=0;i only look at these elements\n"); printf(" 'op'=='function name' is a .apply targeted by v2 called as func(x[i],thresh,val)\n"); return -1.; } op = gargstr(1); n = vector_instance_px(vv, &x); th = *getarg(2); if (ifarg(3)) { val = *getarg(3); } else { val = 0.0; } if (ifarg(4)) {ni = vector_arg_px(4, &ind); f=1;} // just look at the spots indexed if (!strcmp(op,"==")) { if (f==1) {for (i=0; i")) { if (f==1) {for (i=0; ith) { x[(int)ind[i]] = val;}} } else {for (i=0; ith) { x[i] = val;}}} } else if (!strcmp(op,"<")) { if (f==1) {for (i=0; i=")) { if (f==1) {for (i=0; i=th) { x[(int)ind[i]] = val;}} } else {for (i=0; i=th) { x[i] = val;}}} } else if (!strcmp(op,"<=")) { if (f==1) {for (i=0; ith) { /* ? passing thresh */ if (f==0) { x[i]=1.; f=1; } else { x[i]=0.; } } else { /* below thresh */ x[i]=0.; if (f==1) { f=0; } /* just passed going down */ } } return (double)i; } ENDVERBATIM :* v1.xzero() a .where that looks for zero crossings VERBATIM static double xzero(void* vv) { int i, n, nv, up; double *x, *vc, th, cnt=0.; n = vector_instance_px(vv, &x); nv = vector_arg_px(1, &vc); if (ifarg(2)) { th = *getarg(2); } else { th=0.0;} if (vc[0]th) { up=x[i]=1; cnt++; } } return cnt; } ENDVERBATIM :* v1.sw(FROM,TO) switchs all FROMs to TO VERBATIM static double sw(void* vv) { int i, n; double *x, fr, to; n = vector_instance_px(vv, &x); fr = *getarg(1); to = *getarg(2); for (i=0; i-1 && num!=n) { hoc_execerror("Vector size doesn't match.", 0); } for (i=0; i-1 && num!=n) { hoc_execerror("Vector size doesn't match.", 0); } for (i=0; i=last; i++) last=x[i]; } else if (flag==-1) { for (i=1; i100) { /* gamma function */ double x,tmp,ser; x = _ln; tmp=x+5.5; tmp -= (x+0.5)*log(tmp); ser=1.0; for (j=0;j<=5;j++) { x += 1.0; ser += cof[j]/x; } return exp(-tmp+log(2.50662827465*ser)); } else { while (ntop100) { /* gamma function */ double x,tmp,ser; x = _ln; tmp=x+5.5; tmp -= (x+0.5)*log(tmp); ser=1.0; for (j=0;j<=5;j++) { x += 1.0; ser += cof[j]/x; } return (-tmp+log(2.50662827465*ser)); } else { while (ntop> N) #define MUL(x, y, z) { long l = (long)(x) * (long)(y); \ (z)[0] = LOW(l); (z)[1] = HIGH(l); } #define CARRY(x, y) ((long)(x) + (long)(y) > MASK) #define ADDEQU(x, y, z) (z = CARRY(x, (y)), x = LOW(x + (y))) #define X0 0x330E #define X1 0xABCD #define X2 0x1234 #define A0 0xE66D #define A1 0xDEEC #define A2 0x5 #define C 0xB #define SET3(x, x0, x1, x2) ((x)[0] = (x0), (x)[1] = (x1), (x)[2] = (x2)) #define SEED(x0, x1, x2) (SET3(x, x0, x1, x2), SET3(a, A0, A1, A2), c = C) static unsigned x[3] = { X0, X1, X2 }, a[3] = { A0, A1, A2 }, c = C; static unsigned short lastx[3]; static void next(); double my_drand48() { static double two16m = 1.0 / (1L << N); next(); return (two16m * (two16m * (two16m * x[0] + x[1]) + x[2])); } static void next() { unsigned p[2], q[2], r[2], carry0, carry1; MUL(a[0], x[0], p); ADDEQU(p[0], c, carry0); ADDEQU(p[1], carry0, carry1); MUL(a[0], x[1], q); ADDEQU(p[1], q[0], carry0); MUL(a[1], x[0], r); x[2] = LOW(carry0 + carry1 + CARRY(p[1], r[0]) + q[1] + r[1] + a[0] * x[2] + a[1] * x[1] + a[2] * x[0]); x[1] = LOW(p[1] + r[0]); x[0] = LOW(p[0]); } void my_srand48(long seedval) { SEED(X0, LOW(seedval), HIGH(seedval)); } #if 0 #ifdef DRIVER /* This should print the sequences of integers in Tables 2 and 1 of the TM: 1623, 3442, 1447, 1829, 1305, ... 657EB7255101, D72A0C966378, 5A743C062A23, ... */ #include main() { int i; for (i = 0; i < 80; i++) { printf("%4d ", (int)(4096 * my_drand48())); printf("%.4X%.4X%.4X\n", x[2], x[1], x[0]); } } #endif #endif ENDVERBATIM