Mercurial > hg > camir-aes2014
diff toolboxes/FullBNT-1.0.7/bnt/potentials/Tables/rep_mult.c @ 0:e9a9cd732c1e tip
first hg version after svn
author | wolffd |
---|---|
date | Tue, 10 Feb 2015 15:05:51 +0000 |
parents | |
children |
line wrap: on
line diff
--- /dev/null Thu Jan 01 00:00:00 1970 +0000 +++ b/toolboxes/FullBNT-1.0.7/bnt/potentials/Tables/rep_mult.c Tue Feb 10 15:05:51 2015 +0000 @@ -0,0 +1,92 @@ +/* rep_mult.c repmat first two operands to the size provided by */ +/* the third operand, then perform point multiply */ +/* 3 input, 1 output */ +/* C = rep_mult(A, B, sizes) */ + +#include "mex.h" + +void mexFunction(int nlhs, mxArray *plhs[], int nrhs, const mxArray *prhs[]) +{ + double *xp, *yp, *zp, *pSizes; + int xnd, ynd, numElements = 1; + const int *xdim, *ydim; + int i, j, ndim; + int *s, *sx, *sy, *cpsx, *cpsy; + int *subs, *s1, *cpsx2, *cpsy2; + + if (nrhs != 3) + mexErrMsgTxt("Incorrect number of inputs."); + + if (nlhs > 1) + mexErrMsgTxt("Too many output arguments."); + + xnd = mxGetNumberOfDimensions(prhs[0]); + ynd = mxGetNumberOfDimensions(prhs[1]); + xdim = mxGetDimensions(prhs[0]); + ydim = mxGetDimensions(prhs[1]); + ndim = mxGetNumberOfElements(prhs[2]); + + pSizes = mxGetPr(prhs[2]); + + sx = (int *)malloc(sizeof(int)*ndim); + sy = (int *)malloc(sizeof(int)*ndim); + s = (int *)malloc(sizeof(int)*ndim); + s1 = (int *)malloc(sizeof(int)*ndim); + *(cpsx = (int *)malloc(sizeof(int)*ndim)) = 1; + *(cpsy = (int *)malloc(sizeof(int)*ndim)) = 1; + subs = (int *)malloc(sizeof(int)*ndim); + cpsx2 = (int *)malloc(sizeof(int)*ndim); + cpsy2 = (int *)malloc(sizeof(int)*ndim); + for(i=0; i<ndim; i++){ + subs[i] = 0; + sx[i] = (i < xnd) ? xdim[i] : 1; + sy[i] = (i < ynd) ? ydim[i] : 1; + s[i] = (int)pSizes[i]; + s1[i] = s[i] - 1; + numElements *= s[i]; + } + + for(i=0; i<ndim-1; i++){ + cpsx[i+1] = cpsx[i]*sx[i]--; + cpsy[i+1] = cpsy[i]*sy[i]--; + cpsx2[i] = cpsx[i]*sx[i]; + cpsy2[i] = cpsy[i]*sy[i]; + } + cpsx2[ndim-1] = cpsx[ndim-1]*(--sx[ndim-1]); + cpsy2[ndim-1] = cpsy[ndim-1]*(--sy[ndim-1]); + + plhs[0] = mxCreateNumericArray(ndim, s, mxDOUBLE_CLASS, mxREAL); + zp = mxGetPr(plhs[0]); + xp = mxGetPr(prhs[0]); + yp = mxGetPr(prhs[1]); + + for(j=0; j<numElements; j++){ + *zp++ = *xp * *yp; + for(i=0; i<ndim; i++){ + if(subs[i] == s1[i]){ + subs[i] = 0; + if(sx[i]) + xp -= cpsx2[i]; + if(sy[i]) + yp -= cpsy2[i]; + } + else{ + subs[i]++; + if(sx[i]) + xp += cpsx[i]; + if(sy[i]) + yp += cpsy[i]; + break; + } + } + } + free(sx); + free(sy); + free(s); + free(s1); + free(cpsx); + free(cpsy); + free(subs); + free(cpsx2); + free(cpsy2); +}