Mex dll and c file in Win7 64bit (Matlab 2010b)

I just failed to mex the dll and c file correctely. I spent days in the inverstigation from a very beginner to "mex". The miscrosoft visual 2010 professional was used, and the file was downloaded at http://nislab.ee.duke.edu/NISLab/Platelets_%28with_software%29.html. I just follow the instruction and run "compile.m" where there are only two lines: "mex -output plateletApprox.dll plateletApprox.c C:\Progra~1\MATLAB71\extern\lib\win32\lcc\libmwlapack.lib

mex plateletApprox.c C:\Progra~1\MATLAB71\extern\lib\win32\lcc\libmwlapack.lib". But error always occure. Is there any one who can kindly help me for this issue. you just need to download the file and run the compile.m file. Any specified knowledge is NOT necessary. Thanks so much!

[Kaustubha Govind]

What is the exact error that you receive? It seems that the path used "C:\Progra~1\MATLAB71\extern\lib\win32\lcc" should be adjusted according your MATLAB location. Also, it assumes that you are using 32-bit MATLAB.

[Walter Roberson]

Also wasn't lapack split into two libraries?

[MENG]

Thanks for all your reply. I tried to adjust "C:\Progra~1\MATLAB71\extern\lib\win32\lcc" to my MaTLAB location, but I can't find "lcc" under "win64" (since it's 64bit version). Instead of "lcc", I used Microsoft VB 2010 Prof. as the compiler. Is it the problem? Can I use 64bit here?

[MENG]

Below is the error message:

plateletApprox.c

plateletApprox.c(9) : warning C4005: 'max' : macro redefinition

c:\Program Files (x86)\Microsoft Visual Studio 10.0\VC\INCLUDE\stdlib.h(854) : see previous definition of 'max'

plateletApprox.c(10) : warning C4005: 'min' : macro redefinition

c:\Program Files (x86)\Microsoft Visual Studio 10.0\VC\INCLUDE\stdlib.h(855) : see previous definition of 'min'

plateletApprox.c(86) : warning C4013: 'dgels_' undefined; assuming extern returning int

plateletApprox.c(87) : warning C4013: 'dgemv_' undefined; assuming extern returning int

plateletApprox.c(54) : warning C4101: 'Tn' : unreferenced local variable

plateletApprox.c(54) : warning C4101: 'Tm' : unreferenced local variable

plateletApprox.c(109) : error C2065: 'M_PI' : undeclared identifier

plateletApprox.c(114) : error C2065: 'M_PI' : undeclared identifier

plateletApprox.c(336) : warning C4244: 'function' : conversion from 'double' to 'int', possible loss of data

plateletApprox.c(407) : warning C4267: '=' : conversion from 'size_t' to 'int', possible loss of data

plateletApprox.c(408) : warning C4267: '=' : conversion from 'size_t' to 'int', possible loss of data

plateletApprox.c(476) : warning C4267: '=' : conversion from 'size_t' to 'int', possible loss of data

C:\PROGRA~1\MATLAB\R2010B\BIN\MEX.PL: Error: Compile of 'plateletApprox.c' failed.

Answer 1

[James Tursa]

As Walter has already pointed out, versions of MATLAB through R2006b had the BLAS and LAPACK routines in one library, the libmwlapack.lib file. Versions of MATLAB R2007a and later split them into two separate libraries, libmwblas.lib and libmwlapack.lib. The routines you are missing from your link, dgels and dgemv, are BLAS routines. You need to include the libmwblas.lib file in your mex command (in addition to keeping the libmwlapack.lib for the LAPACK routines you use if applicable).

For the M_PI undeclared identifier error, simply add this to the top of the source code file plateletApprox.c:

#ifndef M_PI
#  define M_PI 3.141592653589793238
#endif

After you do these two things you will probably be able to compile and run the code ... and very quickly crash MATLAB in the process. I scanned through the source code and there has been a hard-coded assumption made in the code that the integer types in the BLAS and LAPACK library interfaces are 32-bit int's. For your 64-bit platform, however, this is likely not true ... they are probably 64-bit integers (size_t). This will cause MATLAB to crash. The source code will need to be scrubbed of all the int definitions for the BLAS interface stuff and replaced with mwSignedIndex instead. If you don't know C well enough to do this then you must ask someone for help (either on the website that has the source code or perhaps here in this forum). The simple thing to try would be to replace every instance of " int " with " mwSignedIndex " except for the arguments in mexFunction itself. I.e., leave the following line alone:

void mexFunction(int nlhs,mxArray *plhs[],int nrhs,const mxArray *prhs[])

[James Tursa]

Looking at netlib.org it looks like dgels is an LAPACK routine, not a BLAS routine. In any event, just include both libraries in your mex command and you should have all the bases covered.

[MENG]

Hi James, thanks very much for your answer. I (1) include both libmwblas.lib and libmwlapack.lib; (2) add new definition of M_PI by adding your cod eto the top of the c file. When I run compile.m again, it gives me:

Creating library C:\Users\Meng\AppData\Local\Temp\mex_HtjpaE\templib.x and object C:\Users\Meng\AppData\Local\Temp\mex_HtjpaE\templib.exp

plateletApprox.obj : error LNK2019: unresolved external symbol dgemv_ referenced in function fitPlatelet

plateletApprox.obj : error LNK2019: unresolved external symbol dgels_ referenced in function fitPlatelet

plateletApprox.mexw64 : fatal error LNK1120: 2 unresolved externals

These two errors are what you saw and expected? But I tried replace all "int" with "mwSignedIndex" only expcet the line mexFunction, and still errors(many more than the situation before replacing all "int")? Any suggestions? Thanks so much, and you just make me see the hope of solving this problem!

Answer 2

[James Tursa]

I don't have a 64-bit system to test with, but the following does produce results on a 32-bit system. I have no idea if they are correct. Here is the compile file I used (name it plateletApprox.m and place it in the same directory as the C source file plateletApprox.c):

% plateletApprox.m file: Compile file for Windows & Visual Studio
function varargout = plateletApprox(varargin)
mfile = mfilename;
try
    mname = mfilename('fullpath');
catch
    mname = mfilename;
end
cname = [mname '.c'];
if( isempty(dir(cname)) )
    disp(['Cannot find the file ' cname ' in the same directory as the']);
    disp(['file ' mname '.m. Please ensure that they are in the same']);
    disp('directory and try again. The following file was not found:');
    disp(' ');
    disp(cname);
    disp(' ');
    error(['Unable to compile ' cname]);
end

MEXOPTS = {cname};
comp = computer;
mext = mexext;
lc = length(comp);
lm = length(mext);
cbits = comp(max(1:lc-1):lc);
mbits = mext(max(1:lm-1):lm);
if( isequal(cbits,'64') | isequal(mbits,'64') )
    compdir = 'win64';
    MEXOPTS(end+1) = {'-largeArrayDims'};
else
    compdir = 'win32';
end
libdir = 'microsoft';

lib_blas = [matlabroot '\extern\lib\' compdir '\' libdir '\libmwblas.lib'];
lib_lapack = [matlabroot '\extern\lib\' compdir '\' libdir '\libmwlapack.lib'];
dir_blas = dir(lib_blas);
dir_lapack = dir(lib_lapack);
if( ~isempty(dir_blas) )
    MEXOPTS(end+1) = {lib_blas};
end
if( isempty(dir_lapack) )
    error('LAPACK library not found');
end
MEXOPTS(end+1) = {lib_lapack};

cdold = cd;
cdnew = mname(1:end-length(mfile));
if( ~isempty(cdnew) )
    cd(cdnew);
end

mex(MEXOPTS{:});
cd(cdold);
[varargout{1:nargout}] = plateletApprox(varargin{:});
return
end

And here is the modified source file to hopefully account for the mwSignedIndex stuff:

/* plateletApprox.c */
  /*
 [y,previewEstimate,previewPartition] = CTFPlateletApprox(x[,pen,pen_w,previewBoxSidelength,noiseType]);
 */
#include <math.h>
#include <stdio.h>
#include "mex.h"
#include "matrix.h"
#include <stdlib.h>
#define max(A,B) (A > B ? A : B)
#define min(A,B) (A < B ? A : B)
#define even(x)  ((x & 1) ? 0 : 1)
#define isint(x) ((x - floor(x)) > 0.0 ? 0 : 1)
#define mat(a, i, j) (*(a + (m*(j)+i)))

#ifndef M_PI
#  define M_PI 3.1415926535897932384
#endif

#ifndef mwSignedIndex
#  define mwSignedIndex int
#endif

mwSignedIndex mod(mwSignedIndex A, mwSignedIndex B){
  mwSignedIndex y;
  double a = (double)A;
  double b = (double)B;
  y = (mwSignedIndex)((a/b - floor(a/b))*b);
  return y;
}

double err(double **data, double **est,mwSignedIndex m, mwSignedIndex n,mwSignedIndex oi, mwSignedIndex oj,char *noiseType)
{
  double error = 0;
  mwSignedIndex i,j;
  if (strcmp(noiseType,"Poisson")==0) {
    for (i = 0; i < m; i++) {
      for (j = 0; j < n; j++) {
        if ((est[i+oi][j+oj]<=0.0) && (data[i+oi][j+oj]<=0.0)) {
          error += 0.0;
        }
        else {
          error =  error+est[i+oi][j+oj] - data[i+oi][j+oj]*log(est[i+oi][j+oj]);
        }
      }
    }
  }
  else if (strcmp(noiseType,"Gaussian")==0)
  {
    for (i = 0; i < m; i++) {
      for (j = 0; j < n; j++) {
        error += pow(est[i+oi][j+oj] - data[i+oi][j+oj], 2.0);
      }
    }
  }
  else
    mexErrMsgTxt("noiseType must be either 'Poisson' or 'Gaussian'.");
  return error;
}

void fitPlatelet(double **obs, double **plateFit, mwSignedIndex oi, mwSignedIndex oj, mwSignedIndex m, mwSignedIndex n) {
  char *trans = "N";
  double *T, *details, *obs_1D, *plateFit_1D, *work, Tm, Tn;
  mwSignedIndex i, j, lwork, k = 0;
  mwSignedIndex M, N, NRHS, LDA, LDB, INFO, INCX, INCY;
  double ALPHA, BETA;
  T = mxCalloc(m*n*3,sizeof(double));
  details = mxCalloc(m*n*3,sizeof(double));
  obs_1D = mxCalloc(m*n,sizeof(double));
  plateFit_1D = mxCalloc(m*n,sizeof(double));
  lwork = m*n*2;
  work = mxCalloc(lwork,sizeof(double));
  k =0;
  for (i = 0; i < m; i++)  {
    for (j = 0; j < n; j++)    {
      T[k] = 1.0;
      T[k+m*n] = (double)i;
      T[k+m*n*2] = (double)j;
      details[k] = 1.0;
      details[k+m*n] = (double)i;
      details[k+m*n*2] = (double)j;
      obs_1D[k] = obs[i+oi][j+oj];
      k++;
    }
  }
  M = k;
  N = 3;
  NRHS = 1;
  LDA = M;
  LDB = M;
  ALPHA = 1.0;
  INCX = 1;
  BETA = 0.0;
  INCY = 1;
  dgels(trans,&M,&N,&NRHS,details,&LDA,obs_1D,&LDB,work,&lwork,&INFO);
  dgemv(trans,&M,&N,&ALPHA,T,&M,obs_1D,&INCX,&BETA,plateFit_1D,&INCY);
  k = 0;
  for (i = 0; i < m; i++)  {
    for (j = 0; j < n; j++)    {
      plateFit[i+oi][j+oj] = plateFit_1D[k];
      k++;
    }
  }
  mxFree(T);
  mxFree(details);
  mxFree(obs_1D);
  mxFree(plateFit_1D);
  mxFree(work);
}

void getMask(double** mask, double r, double th, mwSignedIndex N){
  double thr; 
  double th2; 
  double th2r;
  double x,y,R;
  mwSignedIndex i,j;
  double b1, b2,b3,b4;
  thr = th*M_PI/180;
  th2 = fmod(th, 90.0);
  if (th2 >= 45){
    th2 = 90 - th2;
  }
  th2r = th2*M_PI/180;
  b1 = -.5*(cos(th2r)+sin(th2r));
  b2 = -.5*(cos(th2r)-sin(th2r));
  b3 = .5*(cos(th2r)-sin(th2r));
  b4 = .5*(cos(th2r)+sin(th2r));
  x = -N/2 - .5;
  for (i=0;i<N;i++){
    x = x + 1;
    y = N/2 + .5;
    for (j=0;j<N;j++){
      y = y - 1;
      R = r - x*cos(thr)-y*sin(thr);
      if (R < b1){
        mask[i][j] = 0;
      }
      else if ((R >= b1) && (R < b2)){
        mask[i][j] = .25 + ((R*R + R*(cos(th2r) + sin(th2r)) + .25)/(2*cos(th2r)*sin(th2r)));
      }
      else if ((R >= b2) && (R < b3)){
        mask[i][j] = .5 + R/cos(th2r);
      }
      else if ((R >= b3) && (R < b4)){
        mask[i][j] = .75 - tan(th2r)/4 - 1/(4*tan(th2r)) + ((-1*R*R + R*(cos(th2r) + sin(th2r)) + .25)/(2*cos(th2r)*sin(th2r)));
      }
      else{
          mask[i][j] = 1;
        }
        mask[i][j] = (mask[i][j]>0.5);
      }
    }
}

void fitMaskedPlatelet(double **obs, double **plateFit, double **mask, mwSignedIndex oi, mwSignedIndex oj, mwSignedIndex m, mwSignedIndex n) {
  char *trans = "N";
  double *T, *details, *obs_1D, *plateFit_1D, *work, Tm, Tn;
  mwSignedIndex i, j, lwork, k = 0;
  mwSignedIndex M, N, NRHS, LDA, LDB, INFO, INCX, INCY;
  double ALPHA, BETA;
  T = mxCalloc(m*n*6,sizeof(double));
  details = mxCalloc(m*n*6,sizeof(double));
  obs_1D = mxCalloc(m*n,sizeof(double));
  plateFit_1D = mxCalloc(m*n,sizeof(double));
  lwork = m*n*2;
  work = mxCalloc(lwork,sizeof(double));
  k =0;
  for (i = 0; i < m; i++)  {
    for (j = 0; j < n; j++)    {
      T[k] = 1.0*mask[i][j];
      T[k+m*n] = (double)i*mask[i][j];
      T[k+m*n*2] = (double)j*mask[i][j];
      T[k+m*n*3] = 1.0*(1-mask[i][j]);
      T[k+m*n*4] = (double)i*(1-mask[i][j]);
      T[k+m*n*5] = (double)j*(1-mask[i][j]);
      details[k] = 1.0*mask[i][j];
      details[k+m*n] = (double)i*mask[i][j];
      details[k+m*n*2] = (double)j*mask[i][j];
      details[k+m*n*3] = 1.0*(1-mask[i][j]);
      details[k+m*n*4] = (double)i*(1-mask[i][j]);
      details[k+m*n*5] = (double)j*(1-mask[i][j]);
      obs_1D[k] = obs[i+oi][j+oj];
      k++;
    }
  }
  M = k;
  N = 6;
  NRHS = 1;
  LDA = M;
  LDB = M;
  ALPHA = 1.0;
  INCX = 1;
  BETA = 0.0;
  INCY = 1;
  dgels(trans,&M,&N,&NRHS,details,&LDA,obs_1D,&LDB,work,&lwork,&INFO);
  dgemv(trans,&M,&N,&ALPHA,T,&M,obs_1D,&INCX,&BETA,plateFit_1D,&INCY);
  k = 0;
  for (i = 0; i < m; i++)  {
    for (j = 0; j < n; j++)    {
      plateFit[i+oi][j+oj] = plateFit_1D[k];
      k++;
    }
  }
  mxFree(T);
  mxFree(details);
  mxFree(obs_1D);
  mxFree(plateFit_1D);
  mxFree(work);
}

void fitWedgePlate(double** I, mwSignedIndex DIM, double** OUTIMAGE, double** OUTMASK, mwSignedIndex num_r, char *noiseType){
  mwSignedIndex i,j,s,s2 = 0;
  double **mask, **wedgePlateEst;
  double r, th, WC1, WC2, error;
  double bestError, bestR, bestTH = 0;
  mask = (double**)mxCalloc(DIM,sizeof(double*));
  wedgePlateEst = (double**)mxCalloc(DIM,sizeof(double*));
  for (i = 0; i < DIM; i++){
    mask[i] = (double*)mxCalloc(DIM,sizeof(double));
    wedgePlateEst[i] = (double*)mxCalloc(DIM,sizeof(double));
  }
  bestError = 100000000;
  for (s=0;s<num_r;s++){
    for (s2=0;s2<num_r;s2++){
      r = (double)DIM/sqrt(2.0)*(double)s/(double)(num_r+1.0);
      th = 360.0*(double)s2/(double)(num_r+1.0);
      getMask(mask,r,th,DIM);

        fitMaskedPlatelet(I,wedgePlateEst,mask,0,0,DIM,DIM);

        error = err(I,wedgePlateEst,DIM,DIM,0,0,noiseType);
        if (error < bestError) {
          bestError = error;
          bestR = r;
          bestTH = th;
          for (i = 0; i < DIM; i++){
            for (j = 0; j < DIM; j++){
              OUTIMAGE[i][j] = wedgePlateEst[i][j];
              OUTMASK[i][j] = mask[i][j];
            }
          }
        }
      }
    }
    for (i = 0; i < DIM; i++){
    mxFree(mask[i]);
    mxFree(wedgePlateEst[i]);
    }
    mxFree(mask);
    mxFree(wedgePlateEst);
  }  

void wpEst(double **est, double **data,
     double pen, double pen_w, mwSignedIndex m,
     mwSignedIndex n, mwSignedIndex oi, mwSignedIndex oj, 
     char *noiseType, mwSignedIndex numWedges)
{
  mwSignedIndex i,j,k,dyadLen;
  double **dataBlock, **plateEst, **wedgePlateEst, **wedgePlateMask, **data2, **est2;

    dyadLen = m;
    dataBlock = (double**)mxCalloc(m,sizeof(double*));
    wedgePlateEst = (double**)mxCalloc(m,sizeof(double*));
    wedgePlateMask = (double**)mxCalloc(m,sizeof(double*));
    plateEst = (double**)mxCalloc(m,sizeof(double*));
    for (i = 0; i < m; i++) {
      dataBlock[i]= (double *)mxCalloc(n,sizeof(double));
      wedgePlateEst[i] = (double *)mxCalloc(n,sizeof(double));
      wedgePlateMask[i] = (double *)mxCalloc(n,sizeof(double));
      plateEst[i] = (double *)mxCalloc(n,sizeof(double));
    }
    for (i = 0;i<m;i++) {
      for (j = 0;j<n;j++) {
        dataBlock[i][j] = data[i+oi][j+oj];
      }
    }

    fitPlatelet(dataBlock,plateEst,0,0,m,n);

    fitWedgePlate(dataBlock, dyadLen, wedgePlateEst, wedgePlateMask, numWedges ,noiseType);

    if ((err(dataBlock,plateEst,dyadLen,dyadLen,0,0,noiseType)+pen) <
        (err(dataBlock,wedgePlateEst,dyadLen,dyadLen,0,0,noiseType)
         + 2.0*pen_w)) {
      for (j = 0; j < m; j++) {
        for (k = 0; k < n; k++) {
    est[j+oi][k+oj] = plateEst[j][k];
        }
      }
    }
    else {
      for (j = 0; j < m; j++) {
        for (k = 0; k < n; k++) {
    est[j+oi][k+oj] = wedgePlateEst[j][k];
        }
      }
    }      

    for (i = 0; i < dyadLen; i++) {
    mxFree(dataBlock[i]);
    mxFree(wedgePlateEst[i]);
    mxFree(wedgePlateMask[i]);
    mxFree(plateEst[i]);
    }
    mxFree(dataBlock);
    mxFree(wedgePlateEst);
    mxFree(wedgePlateMask);
    mxFree(plateEst);
  }

void plateletApprox(double **x, mwSignedIndex m, mwSignedIndex n, double pen, double pen_w, double **y,
                    double **dispPartition, mwSignedIndex numWedges, char *noiseType)
{
  mwSignedIndex L, im, in, iL, j, k, dyadLen;
  double pMerge,ptmp,pSplit;
  double **merge,**partition, **corners, **cx, **cy;
    if (n==1){
      n = m;
      m = 1;
    }
    L = (int) floor(log(n)/log((double)2));
    merge = (double **)mxCalloc(m,sizeof(double *));
    partition = (double **)mxCalloc(m,sizeof(double *));
    for (im = 0; im < m; im++) {
      merge[im] = (double *)mxCalloc(n,sizeof(double));
      partition[im] = (double *)mxCalloc(n,sizeof(double));
      for (in=0;in<n;in++) {
        y[im][in] = x[im][in];
        partition[im][in] = 1.0;
        dispPartition[im][in] = rand();
      }
    }
    for (iL=0;iL<L;iL++) {
      dyadLen = (int)pow(2.0,iL+1);
      for (im=0;im<m;im+=dyadLen){
        for (in=0;in<n;in+=dyadLen) {
          ptmp = 0;
          for (j=0;j<dyadLen;j++) {
            for (k=0;k<dyadLen;k++){
              /*merge[im][in] += x[im+j][in+k]/dyadLen/dyadLen;*/
              ptmp += dispPartition[im+j][in+k]/dyadLen/dyadLen;
            }
          }
    if ((dyadLen>2)&&(dyadLen<64)) {
      wpEst(merge,x,pen,pen_w,dyadLen,dyadLen,im,in,noiseType,floor(numWedges/sqrt((double)m/(double)dyadLen)));
    }
    else {
      fitPlatelet(x,merge,im,in,dyadLen,dyadLen);
    }

            pMerge = -err(x,merge,dyadLen,dyadLen,im,in,noiseType)-pen;
            pSplit = -err(x,y,dyadLen,dyadLen,im,in,noiseType);
            for (j=0;j<dyadLen;j++) {
              for (k=0;k<dyadLen;k++) {
          pSplit -= partition[im+j][in+k]*pen;
        }
            }
            if (pMerge >= pSplit) {
              for (j=0;j<dyadLen;j++){
                for (k=0;k<dyadLen;k++) {
                  y[im+j][in+k] = merge[im+j][in+k];
                  partition[im+j][in+k] = 1.0/dyadLen/dyadLen;
                  dispPartition[im+j][in+k] = ptmp;
                }
              }
            }
          }
        }
      }
      for (im = 0; im < m; im++) {
        mxFree(merge[im]);
        mxFree(partition[im]);
      }
      mxFree(merge);
      mxFree(partition);
  }

/******************************************/
/* main                                    */
/******************************************/
void mexFunction(int nlhs,mxArray *plhs[],int nrhs,const mxArray *prhs[])
{
  /******************************************/
  /* call in MATLAB:
  [refinedEst, previewEst, previewPartition, corners, shimmeyPreviewEst, shimmeyPreviewPartition] =
   CTFWedgeletApprox(x,previewBoxSidelength,noiseType,pen,pen_w)
   /******************************************/

     double *xCol, *yCol, *y2Col, *yPartitionCol, *y2PartitionCol,*zCol;
     double **x, **xshimmey, **y, **y2, **z;
     double **yPartition, **y2Partition, *y12Corners;
     double **yshimmey, **yshimmeyPartition, **y1Corners, **y2Corners;
     double **corners, *cornersCol, **coarsePlatelets;
     mwSignedIndex m, n, i, j, L, numWedges;
     double mtest, ntest, pen, pen_w, pen_preview, sum = 0;
     mwSignedIndex previewBoxSidelength, numSmallBoxes;
     char *noiseType;
     mxArray *uniqueArgs[2], *uniqueCorners[1];
       /******************************************/
     /* check for correct # of input variables */
     /******************************************/
     if (nrhs>6){
       mexErrMsgTxt("There are at most 6 input parameters allowed!");
       return;
     }
     if (nrhs<1){
       mexErrMsgTxt("There is at least 1 input parameter required!");
       return;
     }
       /**************************************/
     /* get data */
     /**************************************/
       xCol = mxGetPr(prhs[0]);
     n = mxGetN(prhs[0]);
     m = mxGetM(prhs[0]);
     /**************************************/
     /* Check the ROW dimension of input */
     /**************************************/
     if ((m < 2)||(n<2))
     mexErrMsgTxt("This is an image processing method; it cannot process a vector.");
     if(m > 1){
       mtest = (double) log(m)/log((double)2);
       if (!isint(mtest))
         mexErrMsgTxt("The matrix row dimension must be of size m*2^(L)");
     }
     /**************************************/
     /* Check the COLUMN dimension of input  */
     /**************************************/
     if(n > 1){
       ntest = (double) log(n)/log(2.0);
       if (!isint(ntest))
         mexErrMsgTxt("The matrix column dimension must be of size n*2^(L)");
     }  
     /**************************************/
     /* allocate memory */
     /**************************************/
     x = (double **)mxCalloc(m,sizeof(double *));
     xshimmey = (double **)mxCalloc(m,sizeof(double *));
     y = (double **)mxCalloc(m,sizeof(double *));
     yPartition = (double **)mxCalloc(m,sizeof(double *));
     y2 = (double **)mxCalloc(m,sizeof(double *));
     y2Partition = (double **)mxCalloc(m,sizeof(double *));
     yshimmey = (double **)mxCalloc(m,sizeof(double *));
     yshimmeyPartition = (double **)mxCalloc(m,sizeof(double *));
     z = (double **)mxCalloc(m,sizeof(double *));
     coarsePlatelets = (double **)mxCalloc(m,sizeof(double*));
     for (i = 0; i < m; i++) {
       x[i] = (double *)mxCalloc(n,sizeof(double));
       xshimmey[i] = (double *)mxCalloc(n,sizeof(double));
       y[i] = (double *)mxCalloc(n,sizeof(double));
       yPartition[i] = (double *)mxCalloc(n,sizeof(double));
       yshimmey[i] = (double *)mxCalloc(n,sizeof(double));
       yshimmeyPartition[i] = (double *)mxCalloc(n,sizeof(double));
       y2[i] = (double *)mxCalloc(n,sizeof(double));
       y2Partition[i] = (double *)mxCalloc(n,sizeof(double));
       z[i] = (double *)mxCalloc(n,sizeof(double));
       coarsePlatelets[i] = (double *)mxCalloc(n,sizeof(double));
     }
     y1Corners = (double **)mxCalloc(m*n,sizeof(double *));
     y2Corners = (double **)mxCalloc(m*n,sizeof(double *));
     for (i = 0; i < m*n; i++) {
       y1Corners[i] = (double *)mxCalloc(3,sizeof(double));
       y2Corners[i] = (double *)mxCalloc(3,sizeof(double));
     }
     y12Corners = (double *)mxCalloc(m*n*6,sizeof(double));

     for (i = 0; i < m; i++) {
       for (j = 0; j < n; j++) {
         x[i][j] = xCol[i*n+j];
       }
     }

     /**************************************/
     /* set noiseType */
     /**************************************/
     if (nrhs < 2) {
       noiseType = "Gaussian";
     }
     else {
       if (mxIsChar(prhs[1]) != 1) {
         mexErrMsgTxt("noiseType must be a string");
       }
       i = (mxGetM(prhs[1])*mxGetN(prhs[1])) + 1;
       noiseType = mxCalloc(i,sizeof(char));
       mxGetString(prhs[1],noiseType,i);
     }
     if ((strcmp(noiseType,"Poisson") != 0) && (strcmp(noiseType,"Gaussian")!=0)) {
       mexErrMsgTxt("noiseType must be either 'Poisson' or 'Gaussian'.");
     }
     /**************************************/
     /*
      set penalty
      Gaussian
      pen = 2*noise_var*((1/3)*log(n));
        Poisson
      pen = (1/3)log(n)
      */
     /**************************************/
     if ((nrhs < 3)||(mxGetN(prhs[2])==0)) {
       if (strcmp(noiseType,"Poisson")==0) {
         for (i = 0; i < n*m; i++) {
           sum += xCol[i];
         }
         if (sum < 0) {
           mexErrMsgTxt("Poisson data cannot be negative.");
         }
         pen = (double) log(sum)*3.0/3.0;
         pen_w = (double) log(sum)*4.0/3.0;
       }
       else if (strcmp(noiseType,"Gaussian")==0) {
         pen = (double) log(m*n)*3.0/3.0;
         pen_w = (double) log(m*n)*4.0/3.0;
       }
       else {
         printf("noiseType = %s\n",noiseType);
         mexErrMsgTxt("noiseType must be either 'Poisson' or 'Gaussian'.");
       }
     }
     else {
       pen = (double) mxGetScalar(prhs[2]);
       pen_w = (double) mxGetScalar(prhs[3]);
     }
     if ((pen < 0)||(pen_w < 0)) {
       mexErrMsgTxt("The penalty must be positive.");
     }

     /**************************************/
     /* set numWedges */
     /**************************************/
     if ((nrhs < 5) || mxGetN(prhs[4])==0) {
       numWedges = (mwSignedIndex) floor(pow(n*m,1.0/3.0));
     }
     else {
       numWedges = (mwSignedIndex) floor(sqrt(mxGetScalar(prhs[4])));
     }
     /**************************************/
     /* set up outputs */
     /**************************************/

     /* estimate */
     plhs[0] = mxCreateDoubleMatrix(m,n,mxREAL);
     yCol = mxGetPr(plhs[0]);

     /* preview partition */
     plhs[1] = mxCreateDoubleMatrix(m,n,mxREAL);
     yPartitionCol = mxGetPr(plhs[1]);

     /**************************************/
     /* get the preview estimate */
     /**************************************/

     plateletApprox(x, m, n, pen,pen_w, y, yPartition,
                    numWedges, noiseType);

     /**************************************/
     /* refine the estimate */
     /**************************************/
     /*refineEst(z, x, y, pen, pen_w, m, n, previewBoxSidelength, corners,
       numSmallBoxes, noiseType, numWedges);*/
     /**************************************/
     /* proper output format */
     /**************************************/
     for (i = 0; i < m; i++) {
       for (j = 0; j < n; j++) {
         yCol[i*n+j] = y[i][j];
         yPartitionCol[i*n+j] = yPartition[i][j];
       }
     }

     /**************************************/
     /* return to matlab */
     /**************************************/
  }

To compile the mex routine, simply run the example. When plateletApprox is first called it will invoke the m-file and cause the compilation to take place. Once the compilation is done the mex routine will be called from that point on.

[MENG]

Thanks so much, James! Finally I can mex the c file by your code without no error! However, the matlab crashes when I run the "Example.m":

Intel(r) Visual Fortran run-time error

forrtl:severe(38): error during write, unit -1, file CONOUT$

Stack trace terminated abnormally.

The error diaglogue just close Matlab. Any idea for that?

Or I may ask it as another question?

You are so great!

Answer 3

[Ken Atwell]

This mex script presumes the use of the bundled LCC compiler and hard-codes the installation folder of MATLAB. Try modifying compile.m to be just:

 mex plateletApprox.c libmwlapack.lib

Does that work?

[MENG]

Dear Ken, thanks for your reply. Unfortunately it doesn't work. I tried to adjust "C:\Progra~1\MATLAB71\extern\lib\win32\lcc" to my MaTLAB location, but I can't find "lcc" under "win64" (since it's 64bit version). Instead of "lcc", I used Microsoft VB 2010 Prof. as the compiler. You mean I just disregard the first mex: "mex -output plateletApprox.dll plateletApprox.c C:\Progra~1\MATLAB71\extern\lib\win32\lcc\libmwlapack.lib"? I adjusted it as mex -output plateletApprox.dll plateletApprox.c libmwlapack.lib", it still doesn't work.

[Ken Atwell]

Try exactly what I specified (no -output, no path to the .lib file) -- does that work? If not, paste the exact error messages you get.

[MENG]

I first select the compiler: Microsoft Visual C++ 2010

Then run your code exactly to get the error:

plateletApprox.c

plateletApprox.c(9) : warning C4005: 'max' : macro redefinition

c:\Program Files (x86)\Microsoft Visual Studio 10.0\VC\INCLUDE\stdlib.h(854) : see previous definition of 'max'

plateletApprox.c(10) : warning C4005: 'min' : macro redefinition

c:\Program Files (x86)\Microsoft Visual Studio 10.0\VC\INCLUDE\stdlib.h(855) : see previous definition of 'min'

plateletApprox.c(86) : warning C4013: 'dgels_' undefined; assuming extern returning int

plateletApprox.c(87) : warning C4013: 'dgemv_' undefined; assuming extern returning int

plateletApprox.c(54) : warning C4101: 'Tn' : unreferenced local variable

plateletApprox.c(54) : warning C4101: 'Tm' : unreferenced local variable

plateletApprox.c(109) : error C2065: 'M_PI' : undeclared identifier

plateletApprox.c(114) : error C2065: 'M_PI' : undeclared identifier

plateletApprox.c(336) : warning C4244: 'function' : conversion from 'double' to 'int', possible loss of data

plateletApprox.c(407) : warning C4267: '=' : conversion from 'size_t' to 'int', possible loss of data

plateletApprox.c(408) : warning C4267: '=' : conversion from 'size_t' to 'int', possible loss of data

plateletApprox.c(476) : warning C4267: '=' : conversion from 'size_t' to 'int', possible loss of data

C:\PROGRA~1\MATLAB\R2010B\BIN\MEX.PL: Error: Compile of 'plateletApprox.c' failed.

??? Error using ==> mex at 208

Unable to complete successfully.