1. 程式人生 > >Matlab以MEX方式呼叫C原始碼

Matlab以MEX方式呼叫C原始碼

如果我有一個用C語言寫的函式,實現了一個功能,如一個簡單的函式:

double add(double x, double y)

{

return x + y;

}

 
現在我想要在Matlab中使用它,比如輸入:
 
>> a = add(1.1, 2.2)

    3.3000

要得出以上的結果,那應該怎樣做呢?

解決方法之一是要通過使用MEX檔案,MEX檔案使得呼叫C函式和呼叫Matlab的內建函式一樣方便。MEX檔案是由原C程式碼加上MEX檔案專用的介面函式後編譯而成的。可以這樣理解,MEX檔案實現了一種介面,它把在Matlab中呼叫函式時輸入的自變數通過特定的介面調入了C函式,得出的結果再通過該介面調回Matlab。該特定介面的操作,包含在mexFunction這個函式中,由使用者具體設定。

所以現在我們要寫一個包含add和mexFunction的C檔案,Matlab呼叫函式,把函式中的自變數(如上例中的1.1和2.2)傳給 mexFunction的一個引數,mexFunction把該值傳給add,把得出的結果傳回給mexFunction的另一個引數,Matlab通過該引數來給出在Matlab語句中呼叫函式時的輸出值(如上例中的a)。

值得注意的是,mex檔案是與平臺有關的,以我的理解,mex檔案就是另類的動態連結庫。在matlab6.5中使用mex -v  選項,你可以看到最後mex階段有類似如下的資訊:

--> "del _lib94902.obj"   
--> "del "test.exp""   
--> "del "test.lib""

也就是說,雖然在matlab6.5生成的是dll檔案,但是中間確實有過lib檔案生成。

比如該C檔案已寫好,名為add.c。那麼在Matlab中,輸入:

>> mex add.c

就能把add.c編譯為MEX檔案(編譯器的設定使用指令mex -setup),在Windows中,MEX檔案型別為mexw32,即現在我們得出add.mexw32檔案。現在,我們就可以像呼叫M函式那樣呼叫 MEX檔案,如上面說到的例子。所以,通過MEX檔案,使用C函式就和使用M函式是一樣的了。

我們現在來說mexFunction怎樣寫。

mexFunction的定義為:

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

{

}

可以看到,mexFunction是沒返回值的,它不是通過返回值把結果傳回Matlab的,而是通過對引數plhs的賦值。mexFunction的四個引數皆是說明Matlab呼叫MEX檔案時的具體資訊,如這樣呼叫函式時:

>> b = 1.1; c = 2.2;

>> a = add(b, c)

mexFunction四個引數的意思為:

nlhs = 1,說明呼叫語句左手面(lhs-left hand side)有一個變數,即a。

nrhs = 2,說明呼叫語句右手面(rhs-right hand side)有兩個自變數,即b和c。

plhs是一個數組,其內容為指標,該指標指向資料型別mxArray。因為現在左手面只有一個變數,即該陣列只有一個指標,plhs[0]指向的結果會賦值給a。

prhs和plhs類似,因為右手面有兩個自變數,即該陣列有兩個指標,prhs[0]指向了b,prhs[1]指向了c。要注意prhs是const的指標陣列,即不能改變其指向內容。

因為Matlab最基本的單元為array,無論是什麼型別也好,如有double array、 cell array、 struct array……所以a,b,c都是array,b = 1.1便是一個1x1的double array。而在C語言中,Matlab的array使用mxArray型別來表示。所以就不難明白為什麼plhs和prhs都是指向mxArray型別的指標陣列。

完整的add.c如下:

#include "mex.h" // 使用MEX檔案必須包含的標頭檔案

// 執行具體工作的C函式

double add(double x, double y)

{

    return x + y;

}

// MEX檔案介面函式

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

{

    double *a;

    double b, c;

    plhs[0] = mxCreateDoubleMatrix(1, 1, mxREAL);

    a = mxGetPr(plhs[0]);

    b = *(mxGetPr(prhs[0]));

    c = *(mxGetPr(prhs[1]));

    *a = add(b, c);

}

mexFunction的內容是什麼意思呢?我們知道,如果這樣呼叫函式時:

>> output = add(1.1, 2.2);

在未涉及具體的計算時,output的值是未知的,是未賦值的。所以在具體的程式中,我們建立一個1x1的實double矩陣(使用 mxCreateDoubleMatrix函式,其返回指向剛建立的mxArray的指標),然後令plhs[0]指向它。接著令指標a指向plhs [0]所指向的mxArray的第一個元素(使用mxGetPr函式,返回指向mxArray的首元素的指標)。同樣地,我們把prhs[0]和prhs [1]所指向的元素(即1.1和2.2)取出來賦給b和c。於是我們可以把b和c作自變數傳給函式add,得出給果賦給指標a所指向的mxArray中的元素。因為a是指向plhs[0]所指向的mxArray的元素,所以最後作輸出時,plhs[0]所指向的mxArray賦值給output,則 output便是已計算好的結果了。

上面說的一大堆指向這指向那,什麼mxArray,初學者肯定都會被弄到頭暈眼花了。很抱歉,要搞清楚這些亂糟糟的關係,只有多看多練。

實際上mexFunction是沒有這麼簡單的,我們要對使用者的輸入自變數的個數和型別進行測試,以確保輸入正確。如在add函式的例子中,使用者輸入char array便是一種錯誤了。

從上面的講述中我們總結出,MEX檔案實現了一種介面,把C語言中的計算結果適當地返回給Matlab罷了。當我們已經有用C編寫的大型程式時,大可不必在 Matlab裡重寫,只寫個介面,做成MEX檔案就成了。另外,在Matlab程式中的部份計算瓶頸(如迴圈),可通過MEX檔案用C語言實現,以提高計算速度。

以上是對mex檔案的初步認識,下面詳細介紹如何用c語言編寫mex檔案:

1 為什麼要用C語言編寫MEX檔案

 MATLAB是矩陣語言,是為向量和矩陣操作設計的,一般來說,如果運算可以用向量或矩陣實現,其運算速度是非常快的。但若運算中涉及到大量的迴圈處理,MATLAB的速度的令人難以忍受的。解決方法之一為,當必須使用for迴圈時,把它寫為MEX檔案,這樣不必在每次執行迴圈中的語句時MATLAB都對它們進行解釋。

 2 編譯器的安裝與配置

 要使用MATLAB編譯器,使用者計算機上應用事先安裝與MATLAB適配的以下任何一種ANSI C/C++編譯器:

5.0、6.0版的MicroSoft Visual C++(MSVC)

5.0、5.2、5.3、5.4、5.5版的Borland C++

LCC(由MATLAB自帶,只能用來產生MEX檔案)

下面是安裝與配置MATLAB編譯器應用程式MEX的設定的步驟:

(1)在MATLAB命令視窗中執行mex –setup,出現下列提示:

Please choose your compiler for building external interface (MEX) files:

Would you like mex to locate installed compilers [y]/n?

 
(2)選擇y,MATLAB將自動搜尋計算機上已安裝的外部編譯器的型別、版本及所在路徑,並列出來讓使用者選擇:

Select a compiler:

[1] Borland C++Builder version 6.0 in C:\Program Files\Borland

[2] Digital Visual Fortran version 6.0 in C:\Program Files\Microsoft Visual Studio

[3] Lcc C version 2.4 in D:\MATLAB6P5P1\sys\lcc

[4] Microsoft Visual C/C++ version 6.0 in C:\Program Files\Microsoft Visual Studio

[0] None

Compiler:

 (3)選擇其中一種(在這裡選擇了3),MATLAB讓使用者進行確認:

Please verify your choices:

Compiler: Lcc C 2.4

Location: D:\MATLAB6P5P1\sys\lcc

Are these correct?([y]/n):

 
(4)選擇y,結束MATLAB編譯器的配置。

3 一個簡單的MEX檔案例子

【例1】用m檔案建立一個1000×1000的Hilbert矩陣。

tic

m=1000;

n=1000;

a=zeros(m,n);

for i=1:1000

     for j=1:1000

         a(i,j)=1/(i+j);

     end

end

toc

 在matlab中新建一個Matlab_1.cpp 檔案並輸入以下程式:

#include "mex.h"

//計算過程

void hilb(double *y,int n)

{

    int i,j;

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

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

            *(y+j+i*n)=1/((double)i+(double)j+1);

}

//介面過程

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

{

    double x,*y;

    int n;

    if (nrhs!=1)

        mexErrMsgTxt("One inputs required.");

    if (nlhs != 1)

        mexErrMsgTxt("One output required.");

    if (!mxIsDouble(prhs[0])||mxGetN(prhs[0])*mxGetM(prhs[0])!=1)

        mexErrMsgTxt("Input must be scalars.");

    x=mxGetScalar(prhs[0]);

    plhs[0]=mxCreateDoubleMatrix(x,x,mxREAL);

    n=mxGetM(plhs[0]);

    y=mxGetPr(plhs[0]);

    hilb(y,n);

}

 該程式是一個C語言程式,它也實現了建立Hilbert矩陣的功能。在MATLAB命令視窗輸入以下命令:mex Matlab_1.cpp,即可編譯成功。進入該資料夾,會發現多了兩個檔案:Matlab_1.asv和Matlab_1.dll,其中Matlab_1.dll即是MEX檔案。執行下面程式:

tic

a=Matlab_1(1000);

toc

 elapsed_time =

     0.0470

 由上面看出,同樣功能的MEX檔案比m檔案快得多。

4 MEX檔案的組成與引數

MEX檔案的原始碼一般由兩部分組成:

(1)計算過程。該過程包含了MEX檔案實現計算功能的程式碼,是標準的C語言子程式。

(2)入口過程。該過程提供計算過程與MATLAB之間的介面,以入口函式mxFunction實現。在該過程中,通常所做的工作是檢測輸入、輸出引數個數和型別的正確性,然後利用mx-函式得到MATLAB傳遞過來的變數(比如矩陣的維數、向量的地址等),傳遞給計算過程。

MEX檔案的計算過程和入口過程也可以合併在一起。但不管那種情況,都要包含#include "mex.h",以保證入口點和介面過程的正確宣告。注意,入口過程的名稱必須是mexFunction,並且包含四個引數,即:

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

其中,引數nlhs和nrhs表示MATLAB在呼叫該MEX檔案時等式左端和右端變數的個數,例如在MATLAB命令視窗中輸入以下命令:

[a,b,c]=Matlab_1(d,e,f,g)

則nlhs為3,nrhs為4。

MATLAB在呼叫MEX檔案時,輸入和輸出引數儲存在兩個mxArray*型別的指標陣列中,分別為prhs[]和plhs[]。prhs[0]表示第一個輸入引數,prhs[1]表示第二個輸入引數,…,以此類推。如上例中,d→prhs[0],e→prhs[1],f→prhs[2],f→prhs[3]。同時注意,這些引數的型別都是mxArray *。

介面過程要把引數傳遞給計算過程,還需要從prhs中讀出矩陣的資訊,這就要用到下面的mx-函式和mex-函式。

 5 常用的mex-函式和mx-函式

在MATLAB6.5版本中,提供的mx-函式有106個,mex-函式有38個,下面我們僅介紹常用的函式。

5.1入口函式mexFunction

該函式是C MEX檔案的入口函式,它的格式是固定的:

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

說明:MATLAB函式的呼叫方式一般為:[a,b,c,…]=被呼叫函式名稱(d,e,f,…),nlhs儲存了等號左端輸出引數的個數,指標陣列plhs具體儲存了等號左端各引數的地址,注意在plhs各元素針向的mxArray記憶體未分配,需在介面過程中分配記憶體;prhs儲存了等號右端輸入引數的個數,指標陣列prhs具體儲存了等號右端各引數的地址,注意MATLAB在呼叫該MEX檔案時,各輸入引數已存在,所以在介面過程中不需要再為這些引數分配記憶體。

5.2出錯資訊釋出函式mexErrMsgTxt,mexWarnMsgTxt

兩函式的具體格式如下:

#include "mex.h"

void mexErrMsgTxt(const char *error_msg);

void mexWarnMsgTxt(const char *warning_msg);

其中error_msg包含了要顯示錯誤資訊,warning_msg包含要顯示的警告資訊。兩函式的區別在於mexErrMsgTxt顯示出錯資訊後即返回到MATLAB,而mexWarnMsgTxt顯示警告資訊後繼續執行。

5.3 mexCallMATLAB和mexString

兩函式具體格式如下:

#include "mex.h"

int mexCallMATLAB(int nlhs, mxArray *plhs[],

int nrhs, mxArray *prhs[], const char *command_name);

int mexString(const char *command);

mexCallMATLAB前四個引數的含義與mexFunction的引數相同,command_name可以MATLAB內建函式名、使用者自定義函式、M檔案或MEX檔名構成的字串,也可以MATLAB合法的運算子。

mexString用來操作MATLAB空間已存在的變數,它不返回任何引數。

mexCallMATLAB與mexString差異較大,請看下面的例子。

【例2】試用MEX檔案求5階完全圖鄰接矩陣 的特徵值及對應的特徵向量。
5階完全圖的鄰接矩陣為:(這裡找不到圖片了,抱歉。不過不會影響您對本文的理解。)

下面是求該矩陣的MEX檔案。

[Matlab_2.cpp]

#include "mex.h"

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

{

    double x;

    mxArray *y,*z,*w;

    int n;

    if (nrhs!=1)

        mexErrMsgTxt("One inputs required.");

    if (nlhs != 3)

        mexErrMsgTxt("Three output required.");

    if (!mxIsDouble(prhs[0])||mxGetN(prhs[0])*mxGetM(prhs[0])!=1)

        mexErrMsgTxt("Input must be a scalar.");

    x=mxGetScalar(prhs[0]);

    plhs[0]=mxCreateDoubleMatrix(x,x,mxREAL);

    plhs[1]=mxCreateDoubleMatrix(x,x,mxREAL);

    plhs[2]=mxCreateDoubleMatrix(x,x,mxREAL);

    n=mxGetM(plhs[0]);

    y=plhs[0];

    z=plhs[1];

    w=plhs[2];

    //利用mexCallMATLAB計算特徵值

    mexCallMATLAB(1,&plhs[1],1,prhs,"ones");

    mexCallMATLAB(1,&plhs[2],1,prhs,"eye");

    mexCallMATLAB(1,&plhs[0],2,&plhs[1],"-");

    mexCallMATLAB(2,&plhs[1],1,&plhs[0],"eig");

    //演示mexString的功能

    mexString("y=y*2");

    mexString("a=a*2");

}

在MATLAB命令視窗輸入以下命令:

>> mex Matlab_2.cpp

>> clear

>> a=magic(5)

a =

     17     24      1      8     15

     23      5      7     14     16

      4      6     13     20     22

     10     12     19     21      3

     11     18     25      2      9

>> [y,z,w]=Matlab_2(5)

??? Undefined function or variable 'y'.

a =

     34     48      2     16     30

     46     10     14     28     32

      8     12     26     40     44

     20     24     38     42      6

     22     36     50      4     18

y =

      0      1      1      1      1

      1      0      1      1      1

      1      1      0      1      1

      1      1      1      0      1

      1      1      1      1      0

z =

     0.8333    -0.1667    -0.1667     0.2236     0.4472

    -0.1667     0.8333    -0.1667     0.2236     0.4472

    -0.1667    -0.1667     0.8333     0.2236     0.4472

    -0.5000    -0.5000    -0.5000     0.2236     0.4472

          0          0          0    -0.8944     0.4472

w =

     -1      0      0      0      0

      0     -1      0      0      0

      0      0     -1      0      0

      0      0      0     -1      0

      0      0      0      0      4

 由上面可以看出,K5的特徵值為–1和4,其中–1是四重根。MATLAB提供了mexGetVariable、mexPutVariable函式,以實現MEX空間與其它空間交換資料的任務,具體可以參看MATLAB幫助文件。

5.4建立二維雙精度矩陣函式mxCreateDoubleMatrix

其格式具體如下:

#include "matrix.h"

mxArray *mxCreateDoubleMatrix(int m, int n, mxComplexity ComplexFlag);

其中m代表行數,n代表列數,ComplexFlag可取值mxREAL 或mxCOMPLEX。如果建立的矩陣需要虛部,選擇mxCOMPLEX,否則選用mxREAL。

 類似的函式有:

mxCreateCellArray

建立n維元胞mxArray

mxCreateCellMatrix

建立二維元胞mxArray

mxCreateCharArray

建立n維字串mxArray

mxCreateCharMatrixFromStrings

建立二維字串mxArray

mxCreateDoubleMatrix

建立二維雙精度浮點mxArray

mxCreateDoubleScalar

建立指定值的二維精度浮點mxArray

mxCreateLogicalArray

建立n維邏輯mxArray,初值為false

mxCreateLogicalMatrix

建立二維邏輯mxArray,初值為false

mxCreateLogicalScalar

建立指定值的二維邏輯mxArray

mxCreateNumericArray

建立n維數值mxArray

mxCreateNumericMatrix

建立二維數值mxArray,初值為0

mxCreateScalarDouble

建立指定值的雙精度mxArray

MxCreateSparse

建立二維稀疏mxArray

mxCreateSparseLogicalMatrix

建立二維稀疏邏輯mxArray

MxCreateString

建立指定字串的1 n的串mxArray

mxCreateStructArray

建立n維架構mxArray

mxCreateStructMatrix

建立二維架構mxArray

5.5 獲取行維和列維函式mxGetM、mxGetN

其格式如下:

#include "matrix.h"

int mxGetM(const mxArray *array_ptr);

int mxGetN(const mxArray *array_ptr);

與之相關的還有:

mxSetM:設定矩陣的行維

mxSetN:設定矩陣的列維

5.6 獲取矩陣實部和虛部函式mxGetPr、mxGetPi

其格式如下:

#include "matrix.h"

double *mxGetPr(const mxArray *array_ptr);

double *mxGetPi(const mxArray *array_ptr);

與之相關的函式還有:

mxSetPr:設定矩陣的實部

mxSetPi:設定矩陣的虛部

【例3】實現字串的倒序輸出。

#include "mex.h"

void revord(char *input_buf,int buflen,char *output_buf)

{

    int i;

    //實現字串倒序

    for(i=0;i<buflen-1;i++)

        *(output_buf+i)=*(input_buf+buflen-i-2);

}

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

{

    //定義輸入和輸出參量的指標

    char *input_buf,*output_buf;

    int buflen,status;

    //檢查輸入引數個數

    if(nrhs!=1)

        mexErrMsgTxt("One input required.");

    else if(nlhs>1)

        mexErrMsgTxt("Too many output arguments.");

    //檢查輸入引數是否是一個字串

    if(mxIsChar(prhs[0])!=1)

        mexErrMsgTxt("Input must be a string.");

    //檢查輸入引數是否是一個行變數

    if(mxGetM(prhs[0])!=1)

        mexErrMsgTxt("Input must a row vector.");

    //得到輸入字串的長度

    buflen=(mxGetM(prhs[0])*mxGetN(prhs[0]))+1;

    //為輸入和輸出字串分配記憶體

    input_buf=mxCalloc(buflen,sizeof(char));

    output_buf=mxCalloc(buflen,sizeof(char));

    //將輸入參量的mxArray結構中的數值拷貝到C型別字串指標

    status=mxGetString(prhs[0],input_buf,buflen);

    if(status!=0)

        mexWarnMsgTxt("Not enough space. String is truncated.");

    //呼叫C程式

    revord(input_buf,buflen,output_buf);

    plhs[0]=mxCreateString(output_buf);

}

 這個程式中需要注意的地方是mxCalloc函式,它代替了標準C程式中的calloc函式用於動態分配記憶體,而mxCalloc函式採用的是MATLAB的記憶體管理機制,並將所有申請的記憶體初始化為0,因此凡是C程式碼需要使用calloc函式的地方,對應的Mex檔案應該使用mxCalloc函式。同樣,凡是C程式碼需要使用realloc函式的地方,對應的Mex檔案應該使用mxRealloc函式。

在MATLAB命令視窗中對revord.cpp程式程式碼編譯連結:

>> mex revord.cpp

在MATLAB命令視窗中對C-MEX檔案revord.dll進行測試:

>> x='I am student.';

>> revord(x)

ans =

.tneduts ma I

終於寫完了,相信大家對mex檔案應該有點熟悉了,具體還要到實際應用中慢慢體會。