wolffd@0: /* C mex version of max_mult.m in BPMRF2 directory  */
wolffd@0: /* gcc -Wall -I/mit/matlab_v6.5/distrib/bin/glnx86 -c max_mult.c */
wolffd@0: 
wolffd@0: #include <math.h>
wolffd@0: #include "mex.h"
wolffd@0: 
wolffd@0: void mexFunction(int nlhs, mxArray *plhs[], int nrhs, const mxArray *prhs[])
wolffd@0: {
wolffd@0:   int rows,cols,common,m,n,p;
wolffd@0:   double y1, y2;
wolffd@0:   double *arr1, *arr2, *arr3;
wolffd@0: 
wolffd@0: 
wolffd@0:   if (nrhs!=2 || nlhs>1)
wolffd@0:     mexErrMsgTxt("max_mult requires two inputs and one output");
wolffd@0:   if (mxIsChar(prhs[0]) || mxIsClass(prhs[0], "sparse") || mxIsComplex(prhs[0])
wolffd@0:       || mxIsChar(prhs[1]) || mxIsClass(prhs[1], "sparse") || mxIsComplex(prhs[1]))
wolffd@0:     mexErrMsgTxt("Inputs must be real, full, and nonstring");
wolffd@0:   if (mxGetN(prhs[0])!=mxGetM(prhs[1]))
wolffd@0:     mexErrMsgTxt("The number of columns of A must be the same as the number of rows of x");
wolffd@0:   
wolffd@0: 
wolffd@0:   arr1=mxGetPr(prhs[0]);
wolffd@0:   arr2=mxGetPr(prhs[1]);
wolffd@0:   p=mxGetN(prhs[0]);
wolffd@0:   m=mxGetM(prhs[0]);
wolffd@0:   n=mxGetN(prhs[1]);
wolffd@0:   plhs[0]=mxCreateDoubleMatrix(m, n, mxREAL);
wolffd@0:   arr3=mxMalloc(m*n*sizeof(double));
wolffd@0: 
wolffd@0:   for (rows=0; rows<m ; rows++)
wolffd@0:     for (cols=0; cols<n ; cols++)
wolffd@0:     {
wolffd@0:       y1=arr1[rows]*arr2[cols*p];  
wolffd@0:       for (common=1; common<p; common++)
wolffd@0:       {
wolffd@0: 	y2=arr1[rows+common*m]*arr2[common+cols*p];
wolffd@0:         if (y2>y1)
wolffd@0:           y1=y2;
wolffd@0:       }
wolffd@0:       arr3[rows+cols*m]=y1;
wolffd@0:     }
wolffd@0: 
wolffd@0:   mxSetPr(plhs[0], arr3);
wolffd@0: 	 
wolffd@0: }