[转载] C-MEX程序编写

原作者，胡荣春 2006-10-11

1  MEX文件简介

在MATLAB中可调用的C或Fortran语言程序称为MEX文件。MATLAB可以直接把MEX文件视为它的内建函数进行调用。MEX文件是动态链接的子例程，MATLAB解释器可以自动载入并执行它。

MEX文件主要有以下用途：

1．  对于大量现有的C或者Fortran程序可以无须改写成MATLAB专用的M文件格式而在MATLAB中执行。

2．  对于那些MATLAB运算速度过慢的算法，可以用C或者Frotran语言编写以提高效率。

2  例子

为格式详解做准备。下为实例1。

#include "mex.h"

/*    timestwo.c   本MEX文件的目的是实现times two的功能*/

void timestwo(double y[], double x[])

{

       y[0] = 2.0*x[0];

}

 

/*下面这个mexFunction的目的是使MATLAB知道如何调用这个timestwo函数

nlhs是MATLAB命令行方式下输出参数的个数；

*plhs[]是MATLAB命令行方式下的输出参数；

nrhs是MATLAB命令行方式下输入参数的个数；

*prhs[]是MATLAB命令行方式下的输入参数； */

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

{

       double *x,*y; //double指针类型不能改变！！

       int mrows,ncols;

       /* Check for proper number of arguments. */

       if(nrhs!=1)

              mexErrMsgTxt("One input required.");

       else if(nlhs>1)

              mexErrMsgTxt("Too many output arguments"); 

 

         /* 在MATLAB命令行方式下，本MEX文件的调用格式是y=timestwo(x)。输入参数（x）个数＝1，输出参数（y）个数＝1，所以在程序一

        开始就检查nrhs是否＝1以及nlhs是否>1（因为MATLAB有一个缺省输出参数ans，所以nlhs可以=0 */

       mrows = mxGetM(prhs[0]); /* 获得输入矩阵的行数 */

       ncols = mxGetN(prhs[0]); /* 获得输入矩阵的列数 */

       if( !mxIsDouble(prhs[0]) || mxIsComplex(prhs[0]) || !(mrows==1 && ncols==1) )

              mexErrMsgTxt("Input must be a noncomplex scalar double."); /* 判断输入矩阵是否是double类，以及它是否只包括单个元素 */

 

       /* 为输出创建一个矩阵，显然这个矩阵也应该是1x1的 */

       plhs[0] = mxCreateDoubleMatrix(mrows,ncols, mxREAL);

       x = mxGetPr(prhs[0]); /* 获得指向输入/输出矩阵数据的指针 */

       y = mxGetPr(plhs[0]);

       timestwo(y,x); /* 调用C 函数timestwo(y,x) */

}

 

把上面这个文件timestwo.c编辑完成后，在matlab命令行里输入：

mex timestwo.c

matlab会提示你选择一个编译器进行编译，如果安装了VC，则选择VC++即可。编译完成后会在同一目录下生成同名的动态链接库文件timestwo.dll。此后再输入“mex ***.c”编译mex文件时将不再提示用户选择编译器，而自动选择默认的编译器编译。若想改变编译器进行编译，可输入“mex timestwo.c –setup”。

编译完成后即可使用此动态链接库了。在MATLAB命令行下输入：

x = 2;

y = timestwo(x)

将会显示：

y =

4

 

3  MEX文件格式详解

首先编制自己的C算法程序，紧跟着定义mexFunction函数，mexFunction的定义法唯一，它只能是如下形式:

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

其名称和参数类型不许有任何改变，在mexFunciton函数中可以调用你刚定义好的C 程序。

3.1 参数定义

以上面的timestwo.c文件为例，当编译完成后在matlab命令行输入“y = timestwo(x)”时，matlab便随即加载timestwo.dll动态链接库文件。timestwo.dll（mexw32、mexw64）加载完成后首先执行mexFunciton函数，并把输入参数x的值“2”传递给prhs[0]，输入参数的个数“1”传递给plhs。由于matlab里的变量都是以double类型存储的，故mex程序里所有输入输出参数对应的C类型均应为double。这也是为什么上例中定义“double *x,*y;”必须为double指针类型的原因。而且还要注意，必须是double指针，不能仅为double。

3.2  输入输出参数检查

接下来可以做输入参数的检查，确定输入参数的个数是否符合定义，否则给出错误提示信息。

3.3 参数挂接

参数检查完成以后就是对应输入输出参数的挂接(一时想不到更好的词，此处暂借用“挂接”一词表达如下所述的意思)：

        /* 为输出创建一个矩阵，显然这个矩阵也应该是1x1的 */

       plhs[0] = mxCreateDoubleMatrix(mrows,ncols, mxREAL);

       x = mxGetPr(prhs[0]); /* 获得指向输入/输出矩阵数据的指针 */

       y = mxGetPr(plhs[0]);

这里首先给输出参数创建了一个矩阵。注意，不论输出的是一个矩阵还是一个标量，都必须为所有的输出参数创建单独的矩阵。如果是标量，则创建的时候把行和列都设为1即可：

mxCreateDoubleMatrix(1,1, mxREAL);

之所以这样做的原因是matlab里面每个变量都是以矩阵形式存在的。

prhs[]和plhs[]是不能直接在程序中使用的，必须要定义对应的指针变量。故而要对接口参数的指针和程序中使用参数的指针进行挂接，也就是“x = mxGetPr(prhs[0]);”语句所做的工作。看到这里就应该明白为什么void mexFunction函数里面所有的输入输出参数都必须定义为double指针的原因了吧。因为所有定义的参数必须要和prhs[]、plhs[]进行挂接。

3.4  调用功能函数

参数挂接完成以后就可以调用自己编写的功能函数了：

timestwo(y,x); /* 调用C 函数timestwo(y,x) */

这里需要注意，调用的输入参数与定义的参数要相符。如timestwo()函数的定义如下：

void timestwo(double y[], double x[])

两个参数均为指针类型，故调用的时候“timestwo(y,x);”输入的是“y,x”两个指针。而如果timestwo()函数的定义为：

void timestwo(double y[], double x)

则说明输入参数x应该为double类型的变量，这时的输入参数可写为：

timestwo(y,x[0]);

3.5  功能函数的编写

功能函数的编写取决于用户要完成的功能。这里值得提出来说明的有以下2点：

(1)   功能函数名和最终在matlab里调用的函数名没有直接联系，而只于mexFunction里调用的函数名有关系。Matlab调用的函数名取决于最后生成的dll文件名，而这个文件名可以任意更改(但首字母不能为数字)。

(2)   功能函数没有返回值，是void类型，而整个dll模块最终返回给matlab的输出参数是包含在功能函数的输入参数中的。由前面介绍的输入输出参数需要挂接的原因，故功能函数的输入参数中对应的最终输出参数必须定义为double指针类型，如：“void timestwo(double y[], double x)”。其中y是最终的输出参数，必须为double指针类型，而从matlab里传来的输入参数，则视情况而定是否需要定义为指针类型。

3.6   实例2

这个例子是完成dijkstra最短路径算法的程序。

#include "mex.h"

/* 求网络中所有节点到指定节点的最短路，Dijkstra算法*/

 

#define INFN    1000

#define INFI    10000

 

void f_dijkstra_all(double *w, double *D, unsigned long n, unsigned long p)

{

    unsigned long i,j,wtmp,pnew,nu=n-1,np=1;//np为已知距离节点数，nu为未知距离节点数,IDp为已知节点ID

    unsigned long IDu[INFN], IDp[INFN];       //未知节点IDu, 已知节点IDp

 

    p--;

    IDp[0]=p;

   

    for (i=0;i<n;i++)

    {

        w[i]=D[n*i+p];

        if ( i<p )

            IDu[i]=i;

        else

            IDu[i]=i+1;       

    }

   

    while (nu)

    {   wtmp=INFI;

        for (i=0;i<nu;i++)

        {   for (j=0;j<np;j++)

                if ( w[IDp[j]]+D[ n*IDp[j]+IDu[i] ] < w[IDu[i]] )

                    w[IDu[i]]=w[IDp[j]]+D[ n*IDp[j]+IDu[i] ];

           

            if ( w[IDu[i]]<INFI && w[IDu[i]] < wtmp )

            {   wtmp=w[IDu[i]];

                pnew=i;

            }

        }

        if ( wtmp<INFI )

        {   IDp[j]=IDu[pnew];

            np++;

           

            for (i=0;i<nu;i++)               //IDu(pnew)=[];

               if ( i>=pnew )

                    IDu[i]=IDu[i+1];

           

            nu--;

        }

        else

            nu=0;

    }

}

 

/*下面这个mexFunction的目的是使MATLAB知道如何调用这个函数*/

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

 { 

　　if(nrhs!=3)

           mexErrMsgTxt("Three input required.");

　　else if(nlhs>1)

           mexErrMsgTxt("Too many output arguments");

 

　　//mrows = mxGetM(prhs[0]); 　　　　　　　　　　　　　　 /* 获得输入矩阵的行数 */

　　ncols = mxGetN(prhs[0]); 　　　　　　　　　　　　　　　　/* 获得输入矩阵的列数 */

　　if ( mxIsComplex(prhs[0]) )

           mexErrMsgTxt("Input must be a noncomplex.");    　　/*判断输入矩阵是否为复数*/

   

　　plhs[0] = mxCreateDoubleMatrix(1,ncols, mxREAL);    　　/* 为输出创建一个矩阵*/

　　D = mxGetPr(prhs[0]);                              　　　　 　　  /* 获得指向输入/输出矩阵数据的指针 */

　　n = mxGetPr(prhs[1]);

　　p = mxGetPr(prhs[2]);

　　w = mxGetPr(plhs[0]);

 

　　f_dijkstra_all(w, D, n[0], p[0]);                          /* 调用C 函数*/

}

4  补充说明

在编写功能函数时，需要特别注意得是：matlab中的矩阵传递到C中变成数组时，元素的排列顺序是先列后行！

举个例子：把某个矩阵D作为输入参数，此矩阵为：

D＝[　　1   2     3

　　　　4   5     6]

而把此输入参数挂载到某个定义好的double指针变量Dp后，依次写出Dp所指向的元素，你会发现：

Dp[0]=1

Dp[1]=4

Dp[2]=2

Dp[3]=5

Dp[4]=3

Dp[5]=6

所以要正确无误的引用原来矩阵中第i行第j列的元素，应使用Dp[i-1+(j-1)*n]。其中n为矩阵的行数。

 

参考文献

MATLAB与C 的接口问题分类:我的技术我做主, 北京理工大学BBS