|
wolffd@0
|
1 /****************************************************
|
|
wolffd@0
|
2 A = mult_by_array(big, small)
|
|
wolffd@0
|
3 implicitely copies small |big|/|small| times
|
|
wolffd@0
|
4 and then does element-wise multiplication.
|
|
wolffd@0
|
5
|
|
wolffd@0
|
6 i.e.,
|
|
wolffd@0
|
7 C = repmat(small(:), 1, length(big(:))/length(small(:)))
|
|
wolffd@0
|
8 A = reshape(big(:) .* C(:), size(big))
|
|
wolffd@0
|
9
|
|
wolffd@0
|
10 However, this C version avoids the expense of the repmat.
|
|
wolffd@0
|
11
|
|
wolffd@0
|
12 Written by wei.hu@intel.com, 28 Jan 2002.
|
|
wolffd@0
|
13 /****************************************************/
|
|
wolffd@0
|
14
|
|
wolffd@0
|
15
|
|
wolffd@0
|
16 #include "mex.h"
|
|
wolffd@0
|
17
|
|
wolffd@0
|
18 void mexFunction(int nlhs, mxArray *plhs[], int nrhs, const mxArray *prhs[])
|
|
wolffd@0
|
19 {
|
|
wolffd@0
|
20 double *sp, *zp;
|
|
wolffd@0
|
21 int i, j, NB, NS, xnd, ynd, ndim;
|
|
wolffd@0
|
22 const int *xdim, *ydim;
|
|
wolffd@0
|
23 int *s, *sx, *sy, *cpsy, *subs, *cpsy2;
|
|
wolffd@0
|
24
|
|
wolffd@0
|
25 if (nrhs != 2)
|
|
wolffd@0
|
26 mexErrMsgTxt("Incorrect number of inputs.");
|
|
wolffd@0
|
27
|
|
wolffd@0
|
28 if (nlhs > 1)
|
|
wolffd@0
|
29 mexErrMsgTxt("Too many output arguments.");
|
|
wolffd@0
|
30
|
|
wolffd@0
|
31 plhs[0] = mxDuplicateArray(prhs[0]);
|
|
wolffd@0
|
32 zp = mxGetPr(plhs[0]);
|
|
wolffd@0
|
33 sp = mxGetPr(prhs[1]);
|
|
wolffd@0
|
34
|
|
wolffd@0
|
35 xnd = mxGetNumberOfDimensions(prhs[0]);
|
|
wolffd@0
|
36 ynd = mxGetNumberOfDimensions(prhs[1]);
|
|
wolffd@0
|
37 xdim = mxGetDimensions(prhs[0]);
|
|
wolffd@0
|
38 ydim = mxGetDimensions(prhs[1]);
|
|
wolffd@0
|
39 ndim = xnd;
|
|
wolffd@0
|
40
|
|
wolffd@0
|
41 NB = mxGetNumberOfElements(prhs[0]);
|
|
wolffd@0
|
42 NS = mxGetNumberOfElements(prhs[1]);
|
|
wolffd@0
|
43
|
|
wolffd@0
|
44 if(NS == 1){
|
|
wolffd@0
|
45 for(i=0; i<NB; i++){
|
|
wolffd@0
|
46 *zp++ *= *sp;
|
|
wolffd@0
|
47 }
|
|
wolffd@0
|
48 return;
|
|
wolffd@0
|
49 }
|
|
wolffd@0
|
50
|
|
wolffd@0
|
51 if(NS == NB){
|
|
wolffd@0
|
52 for(i=0; i<NB; i++){
|
|
wolffd@0
|
53 *zp++ *= *sp++;
|
|
wolffd@0
|
54 }
|
|
wolffd@0
|
55 return;
|
|
wolffd@0
|
56 }
|
|
wolffd@0
|
57
|
|
wolffd@0
|
58 sx = (int *)malloc(sizeof(int)*ndim);
|
|
wolffd@0
|
59 sy = (int *)malloc(sizeof(int)*ndim);
|
|
wolffd@0
|
60 s = (int *)malloc(sizeof(int)*ndim);
|
|
wolffd@0
|
61 *(cpsy = (int *)malloc(sizeof(int)*ndim)) = 1;
|
|
wolffd@0
|
62 subs = (int *)malloc(sizeof(int)*ndim);
|
|
wolffd@0
|
63 cpsy2 = (int *)malloc(sizeof(int)*ndim);
|
|
wolffd@0
|
64 for(i=0; i<ndim; i++){
|
|
wolffd@0
|
65 subs[i] = 0;
|
|
wolffd@0
|
66 sx[i] = xdim[i];
|
|
wolffd@0
|
67 sy[i] = (i < ynd) ? ydim[i] : 1;
|
|
wolffd@0
|
68 s[i] = sx[i] - 1;
|
|
wolffd@0
|
69 }
|
|
wolffd@0
|
70
|
|
wolffd@0
|
71 for (i = 0; i < ndim-1; i++){
|
|
wolffd@0
|
72 cpsy[i+1] = cpsy[i]*sy[i]--;
|
|
wolffd@0
|
73 cpsy2[i] = cpsy[i]*sy[i];
|
|
wolffd@0
|
74 }
|
|
wolffd@0
|
75 cpsy2[ndim-1] = cpsy[ndim-1]*(--sy[ndim-1]);
|
|
wolffd@0
|
76
|
|
wolffd@0
|
77 for(j=0; j<NB; j++){
|
|
wolffd@0
|
78 *zp++ *= *sp;
|
|
wolffd@0
|
79 for(i=0; i<ndim; i++){
|
|
wolffd@0
|
80 if(subs[i] == s[i]){
|
|
wolffd@0
|
81 subs[i] = 0;
|
|
wolffd@0
|
82 if(sy[i]) sp -= cpsy2[i];
|
|
wolffd@0
|
83 }
|
|
wolffd@0
|
84 else{
|
|
wolffd@0
|
85 subs[i]++;
|
|
wolffd@0
|
86 if(sy[i]) sp += cpsy[i];
|
|
wolffd@0
|
87 break;
|
|
wolffd@0
|
88 }
|
|
wolffd@0
|
89 }
|
|
wolffd@0
|
90 }
|
|
wolffd@0
|
91 free(sx);
|
|
wolffd@0
|
92 free(sy);
|
|
wolffd@0
|
93 free(s);
|
|
wolffd@0
|
94 free(cpsy);
|
|
wolffd@0
|
95 free(subs);
|
|
wolffd@0
|
96 free(cpsy2);
|
|
wolffd@0
|
97 }
|
