Matlab 連結 C Library

決定寫 C library 之前，建議是先用 profile reporter 確認程式的效能瓶頸。針對造成瓶頸的函式再花時間去打造其 C library 比較有經濟效益。

供 Matlab 連結使用之 C DLL 程式碼可分為兩個部分：第一部分為介面程式，用來處理 Matlab 與 C 之間的變數型態轉換；另一部分為功能主體，用來做函數所需要的運算。該 DLL 必須遵循 mexFunction 架構，並適當轉換輸入及輸出變數型態。

mexFunction 架構

典型的 mexFunction 程式碼的架構為

#include "mex.h"

void mexFunction(int nlhs, mxArray *plhs[], int nrhs, const mxArray *prhs[]) {
    /* parse input arguments */
    /* create output arguments */
    /* base function */
}

首先必須要引入標頭檔 mex.h，這裡面有 Matlab 矩陣型態的定義以及各個輔助函式的 prototype。接著宣告進入點函數 mexFunction(...)，差不多就是 main(...) 的意思，這函數的 prototype 是不可以改變的：

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

mexFunction 包含四個引數，分別為

nlhs：呼叫函式時，等號左邊的變數數目，即 plhs 之個數
plhs：呼叫函式時，等號左邊之變數群
nrhs：呼叫函式時，等號右邊的變數數目，即 prhs 之個數
prhs：呼叫函式時，等號右邊之變數群

由於 Matlab 的函數呼叫允許多個輸出參數，如：

[M, I] = larger(A, B);

呼叫時 nlhs 以及 plhs 即用來存放等號左邊的輸出參數 M 與 I，nrhs 與 prhs 即存放等號右邊的輸入參數 A 及 B。此時，nlhs 及 nrhs 均為 2 。但如果使用者以

M = larger(A, B);

的形式呼叫，這時候 nlhs 會為 1。DLL 裡面可以依照 nlhs 及 nrhs 的不同，判斷使用者呼叫的輸入及輸出型態，以便進行 function overloading。

Matlab 中所有的變數，無論是矩陣或純量型態，在 mexFunction 中都以 mxArray 型態表示，純量就視為一個「大小為1x1之矩陣」。plhs 及 prhs 之變數型態為 mxArray **，亦即一「指向 Matlab 矩陣之指標陣列」。以前述 [M, I]=larger(A, B); 的例子來看，plhs[0] 指向 M 而 plhs[1] 指向 I，prhs[0] 指向 A 而 prhs[1] 則指向 B。

Matlab 矩陣與 C 陣列轉換

mexFunction 中的 mxArray 只是一個指標，實際仍然是 Matlab 之資料型態，並無法直接被 C 程式碼使用。要用 Matlab 提供的一些輔助函式進行轉換，從 mxArray 中取得可供 C 使用的數值。詳細可以參考 Matlab 手冊中關於 MX-Functions 的部分。

要從 mxArray 裡面取得資料，我們至少必須知道兩件事：矩陣的大小和矩陣的內容。矩陣大小可由 mxGetM 和 mxGetN 取得，矩陣內容則由 mxGetPr 以及 mxGetPi 取得。

int mxGetM(const mxArray *array_ptr);
int mxGetN(const mxArray *array_ptr);

mxGetM 及 mxGetN 分別用來取得矩陣之第一維以及第二維尺寸。函數傳回一個 int 整數，即為對應之尺寸大小。若矩陣超過二維，則必須使用 mxGetNumberOfDimensions 及 mxGetDimensions 函數。

double *mxGetPr(const mxArray *array_ptr);
double *mxGetPi(const mxArray *array_ptr);

mxGetPr 及 mxGetPi 則分別用來取出矩陣的實數與虛數部分，並回傳一個 double * 指向所取出資料陣列的開頭。無論原始矩陣的維度(dimension)是多少，所取出之資料均為一個一維陣列，因此必須配合由 mxGetM 及 mxGetN 取得的矩陣大小並計算正確之索引位置。

配合上述輔助函式，可利用下列程式碼，把 [M, I] = larger(A, B); 中矩陣 A 的實數部分變成 C 的變數：

double *A;
int mA, nA;
A = mxGetPr(prhs[0]);
mA = mxGetM(prhs[0]);
nA = mxGetN(prhs[0]);

接著便可以使用 A[i] 來存取矩陣，如果是二維陣列就使用 A[i+j*mA] 之方式。

建立輸出矩陣

函式既然是給 Matlab 用的，輸出當然也必須是 Matlab 矩陣。輸出矩陣的大小及內容是由 C DLL 決定，所以輸出必須由 DLL 負責建立。利用 mxCreateDoubleMatrix 我們可以建立 double 型態的輸出矩陣。

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

這函數可創造一個任意大小、實數或複數型態的二維 Matlab 矩陣。使用時給定矩陣的大小(m-by-n)以及矩陣的型態(mxREAL或mxCOMPLEX擇一)。函式即會創造一個新的 Matlab 矩陣，配置足夠的空間，並回傳一個 mxArray 指標指向所創造的新矩陣。若要創造多維矩陣或非實數/虛數型態之矩陣，可參考 mxCreateNumericArray 以及 mxCreateNumericMatrix。

此處用 mxCreateDoubleMatrix 建立輸出矩陣 M 及 I，並利用前述 mxGetPr 取得輸出矩陣對應之 C 陣列：

double *M, *I;
plhs[0] = mxCreateDoubleMatrix(mA, nA, mxREAL);
M = mxGetPr(plhs[0]);
plhs[1] = mxCreateDoubleMatrix(mA, nA, mxREAL);
I = mxGetPr(plhs[1]);

上面程式碼將 M 與 I 設定為與 A 相同大小之實數矩陣，並取出其資料指標 M 及 I 以備寫入輸出資料用。做好輸入及輸出參數的設置，接下來就依照一般 C 語法撰寫函式的主體功能部分，寫完後存檔為 larger.c：

#include "mex.h"
 
void mexFunction(int nlhs, mxArray *plhs[], int nrhs, const mxArray *prhs[])
{
    double *A, *B, *M, *I;
    int m, n, s, i;

    /* parse input arguments */ 
    A = mxGetPr(prhs[0]);
    B = mxGetPr(prhs[1]);
    m = mxGetM(prhs[0]);
    n = mxGetN(prhs[0]);

    /* create output arguments */ 
    plhs[0] = mxCreateDoubleMatrix(m, n, mxREAL);
    M = mxGetPr(plhs[0]);
    if(nlhs>1) {
      plhs[1] = mxCreateDoubleMatrix(m, n, mxREAL);
      I = mxGetPr(plhs[1]);
    }

    /* base function */
    s = m*n;
    if(nlhs==1) {
      for(i=0; i<s; i++) {
        if(A[i]>=B[i]) M[i] = A[i];
        else M[i] = B[i];
      }
    } else {
      for(i=0; i<s; i++) {
        if(A[i]>=B[i]) {
          M[i] = A[i];
          I[i] = 0;
        } else {
          M[i] = B[i];
          I[i] = 1;
        }
      }
    }
}

設定編譯器

在 Matlab 下編譯 C DLL 的指令為 mex。使用 mex 之前必須設定編譯使用的 C 編譯器。在 Matlab 下輸入 mex -setup 指令進行設定：

>> mex -setup
Select a compiler:
[1] Lcc C version 2.4 in C:\\PROGRAM FILES\\MATLAB704\\syslcc 

[0] None

Compiler:

mex 會偵測電腦中已安裝適合的C編譯器，Lcc 是 Matlab 內建的，可以用不用怕。

編譯

選好編譯器後就可以開始編譯了，假設剛剛寫好的檔案存檔為 larger.c，將 Matlab 工作目錄切換至剛剛存檔 larger.c 的目錄，執行：

mex larger.c

執行完畢後就可以在該目錄下看見編譯好的函式庫 larger.dll，如果失敗就要依照錯誤訊息檢查 DLL 程式碼以及 mex 設定。

測試

現在可以來測試看看編譯好的 DLL：

>> A=[1 2 3; 4 5 6; 7 8 9];
>> B=[9 8 7; 6 5 4; 3 2 1];
>> [M,I]=larger(A, B);
>> M

M =

     9     8     7
     6     5     6
     7     8     9

>> I

I =

     1     1     1
     1     0     0
     0     0     0

>>

也可以測試一下 M=larger(A,B); 是否能正常工作。接著試試看

>> larger(A,B)

------------------------------------------------------------------------
       Segmentation violation detected at Sat Dec 10 09:52:22 2005
------------------------------------------------------------------------

(以下刪除)

出現 segmentation fault，這是因為程式碼並沒有考慮「未指定任何輸出」，即 nlhs 為 0 時的情形。一般函式在沒有指定輸出的情況下，應該要傳回一個輸出結果。把程式碼中，

if(nlhs==1) {

這一行，修正為

if(nlhs<=1)

就可以解決這個問題。