如何在Mexfile中的matlab(矩阵,单元格)和c ++(向量或自定义类)之间正确转换变量

Posted

技术标签:

【中文标题】如何在Mexfile中的matlab(矩阵,单元格)和c ++(向量或自定义类)之间正确转换变量【英文标题】:How to convert varibles correctely between matlab (matrix, cell) and c++ (vectors or self-defining class) in Mexfile 【发布时间】:2018-08-03 16:41:19 【问题描述】:

我刚开始使用 mex 文件,在 mexfunction 中转换变量时遇到了一些问题,特别是当我想将矩阵转换为自定义类时(例如“vema.h”中的向量)和MexFunction 中的一个单元格到一个向量。我们提出了这个问题(链接:MATLAB crashes when I run MEX file),这里我们重新整理一下。

下面是计算函数:

#include "mex.h"
#include "matrix.h"
#include <omp.h>
#include "eig3.h"
#include <stdlib.h>
#include <iostream>
#include <cmath>
#include <vector>
#include <fstream>
#include <iomanip>
#include <sys/types.h>
#include <sys/stat.h>
#include "vema.h" //in vema.h, 'Vector' class has been defined inside which is used in the definition of Fb and V in mexfunction.

using namespace std;

Vector closestPointTriangle(Vector&, Vector&, Vector&, Vector&, double&, double&, double&);

// Generates pomwSize-triangle proximity lists using the linked cell algorithm
void createNNLtriangle(vector<int>* NNLt, Vector* Ut, Vector* faces, int* SN, mwSize nsn, mwSize nf, double hs, double bw, double mw) 

    int mx = max(1, (int)(bw/mw)); // ** = 40 cells bw=3.2, mw=0.08
    vector<int> head(mx*mx*mx, -1);
    vector<int> list(nf);
//  std::vector<int> head(mx*mx*mx, -1); //****** mx*mx*mx cells nomber, size mx*mx*mx vector with all values are -1, 40*40*40 = 64000
//  std::vector<int> list(nf); // **** nf = 101882
    int xa, ya, za, xi, yi, zi;
    double ub, vb, wb;
    int pt, tri;
    Vector cog;
    for (int i = 0; i < nf; i++)  // Divide triangle faces mwSizeo cells, i index of face
        //Vector cog = (Ut[faces[i].n1] + Ut[faces[i].n2] + Ut[faces[i].n3])/3.0;
        cog = (Ut[(int)faces[i].x] + Ut[(int)faces[i].y] + Ut[(int)faces[i].z])/3.0;
        int xa = (int)((cog.x + 0.5*bw)/bw*mx);
        int ya = (int)((cog.y + 0.5*bw)/bw*mx);
        int za = (int)((cog.z + 0.5*bw)/bw*mx);
        int tmp =  mx*mx*za + mx*ya + xa; // *** 1641838 > 64000

        list[i]=head[mx*mx*za + mx*ya + xa];
        head[mx*mx*za + mx*ya + xa] = i;
    
    #pragma omp parallel for
    for (int i = 0; i < nsn; i++)  // Search cells around each pomwSize and build proximity list
        int pt = SN[i];
        NNLt[i].clear();
        int xa = (int)((Ut[pt].x + 0.5*bw)/bw*mx);
        int ya = (int)((Ut[pt].y + 0.5*bw)/bw*mx);
        int za = (int)((Ut[pt].z + 0.5*bw)/bw*mx);

        for (int xi = max(0, xa-1); xi <= min(mx-1, xa+1); xi++)// *** Browse head list
        for (int yi = max(0, ya-1); yi <= min(mx-1, ya+1); yi++)
        for (int zi = max(0, za-1); zi <= min(mx-1, za+1); zi++) 
            int tri = head[mx*mx*zi + mx*yi + xi];
            while (tri != -1) 
                if ( pt != (int)faces[tri].x && pt != (int)faces[tri].y && pt != (int)faces[tri].z )               
                    if ( (closestPointTriangle(Ut[pt], Ut[(int)faces[tri].x], Ut[(int)faces[tri].y], Ut[(int)faces[tri].z], ub, vb, wb) - Ut[pt]).length() < hs) 
                        NNLt[i].push_back(tri);
                    
                
                tri = list[tri];
            
        
    



// Returns the closest pomwSize of triangle abc to pomwSize p ***** a or b or c, if not, pt projection through the barycenter inside the triangle 
Vector closestPointTriangle(Vector& p, Vector& a, Vector& b, Vector& c, double& u, double& v, double& w) 

    Vector ab = b - a;
    Vector ac = c - a;
    Vector ap = p - a;
    double d1 = ab.dot(ap);
    double d2 = ac.dot(ap);
    if (d1 <= 0.0 && d2 <= 0.0) 
        u = 1.0;
        v = 0.0;
        w = 0.0;
        return a;
    
    Vector bp = p - b;
    double d3 = ab.dot(bp);
    double d4 = ac.dot(bp);
    if (d3 >= 0.0 && d4 <= d3) 
        u = 0.0;
        v = 1.0;
        w = 0.0;
        return b;
    
    double vc = d1*d4 - d3*d2;
    if (vc <= 0.0 && d1 >= 0.0 && d3 <= 0.0) 
        v = d1 / (d1 - d3);
        u = 1.0 - v;
        w = 0.0;
        return a + ab * v;
    
    Vector cp = p - c;
    double d5 = ab.dot(cp);
    double d6 = ac.dot(cp);
    if (d6 >= 0.0 && d5 <= d6) 
        u = 0.0;
        v = 0.0;
        w = 1.0;
        return c;
    
    double vb = d5*d2 - d1*d6;
    if (vb <= 0.0 && d2 >= 0.0 && d6 <= 0.0) 
        w = d2 / (d2 - d6);
        u = 1.0 - w;
        v = 0.0;    
        return a + ac * w;
    
    double va = d3*d6 - d5*d4;
    if (va <= 0.0 && (d4 - d3) >= 0.0 && (d5 - d6) >= 0.0) 
        w = (d4 - d3) / ((d4 - d3) + (d5 - d6));
        u = 0.0;
        v = 1.0 - w;
        return b + (c - b) * w;
    
    double denom = 1.0 / (va + vb + vc);
    v = vb * denom;
    w = vc * denom;
    u = 1.0 - v - w;
    return a + ab * v + ac * w;

mexfunction:

 /* The gateway function */
    void mexFunction(mwSize nlhs, mxArray *plhs[],
                     mwSize nrhs, const mxArray *prhs[])
          
    vector<int> *NNLt;
    Vector *V;
    Vector *Fb;
    int *sn;
    mwSize nsn; 
    mwSize nf; 
    double hs;
    double bw;
    double mw; 
    mwSize i;

    /* check for proper number of arguments */
    if(nrhs!=9) 
        mexErrMsgIdAndTxt("MyToolbox:arrayProduct:nrhs","Nine inputs required.");
    
    if(nlhs!=1) 
        mexErrMsgIdAndTxt("MyToolbox:arrayProduct:nlhs","One output required.");
        

    /* get the value of the scalar input  */
    nsn = mxGetScalar(prhs[4]);
    nf = mxGetScalar(prhs[5]);
    hs = mxGetScalar(prhs[6]);  
    bw = mxGetScalar(prhs[7]);
    mw = mxGetScalar(prhs[8]);
    /* create a pomwSizeer to the real data in the input matrix  */
    V = (Vector *)mxGetData(prhs[1]);
    Fb = (Vector *)mxGetData(prhs[2]);
    sn = (int *)mxGetData(prhs[3]);

    for(i=0; i<nsn; i++) 
        mxArray *tmp = mxGetCell(prhs[0],i);
        std::size_t N = mxGetN(tmp); // number of columns
        vector<std::vector<int>> NNLt(nsn, std::vector<int>(N));
        double* ptr = mxGetPr(tmp);
        copy(ptr, ptr+N, NNLt[i].begin());  
    ;
    /* call the computational routine */    
    createNNLtriangle(NNLt, V, Fb, sn, nsn, nf, hs, bw, mw);    

    /* create the output matrix */
    plhs[0] = mxCreateCellMatrix(nsn,50);
    /* get a pomwSizeer to the real data in the output matrix */
    for(i=0;i<nsn;i++)      
         mxArray* tmp = mxCreateDoubleMatrix(1, NNLt[i].size(), mxREAL);
         copy(NNLt[i].begin(), NNLt[i].end(), mxGetPr(tmp));     
         mxSetCell(plhs[0], i, tmp);
    

下面是Matlab中调用Mexfile的部分:

NNLt = cell(1,nsn);
V = A(1:n,1:3); //A: double Matrix
Fb = A(n:2*n,1:3);
nf = size(Fb,1);
sn = B(1,:); //B: int Matrix
mw = 8.0*a;
hs = 0.6*a; 
hc = 0.2*a; 
nsn=length(sn);
parfor i = 1:nsn
    maxDist = max(maxDist, length(V(sn(i),:) - Vtold(i,:)));
end;
if maxDist > 0.5*(hs-hc)
    [NNLt] = createNNLtriangle(NNLt, V, Fb, sn, nsn, nf, hs, bw, mw);
    for i = 1:nsn
        Vtold(i,:) = V(sn(i),:);
    end;
end;

还有 vema.h 的一部分:

#include <cmath>

class Vector
public:
    double x, y, z;
    Vector(): x(0.0), y(0.0), z(0.0) ;
    Vector(double ax, double ay, double az): x(ax), y(ay), z(az) ;
...

谢谢。

小雨。

【问题讨论】:

您好。请注意,*** 的要求之一是为您可能遇到的任何错误或问题提供minimal reproducible example。请注意,Minimal 是必需的。否则任何人都很难提供帮助,因为代码太长了。尝试将所有发布的代码减少到最低限度,以重现您遇到的错误。 请注意,您还没有在这个“问题”中提出任何问题 另请注意,您似乎忽略了上一篇文章中给出的任何建议。 V = (Vector *)mxGetData(prhs[1]); 是不允许的。 mxGetData 返回指向 prhs[1] 具有的相同类型数据的指针。如您所知,MATLAB 没有类型 Vector*。例如 double V=static_cast&lt;double *&gt;(mexGetData(prhs[1])) 会起作用,因为 V 在 C++ 和 MATLAB 中都是双精度的。 【参考方案1】:

让我向您展示如何使用mex 将数据从 MATLAB 获取到 C++ 的简单示例:

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

    float * myvar= static_cast<float *>(mxGetData(prhs[0]));
    mexPrintf("%f \n",myvar[0]);

但是,请注意,这仅在您的 mexfucntion 的输入为 float 或在 MATLAB 中为 single 时才有效。

在 MATLAB 中:

A=rand(10,1);
mymexfucntion(A);

将导致崩溃,因为Adouble 类型。而是

A=rand(10,1);
mymexfucntion(single(A));

会起作用。因为它与 C++ 变量的数据类型相同。如果 C++ 变量是类、向量或简单类型,则无关紧要,如果要转换和mxGetData 它,MATLAB 变量需要是相同的类型。你不能用Vector 做到这一点,因为这是你编造的。

这一切对您的问题意味着什么?简单:您需要编写自己的 MATLAB-Vector 转换器。


附带说明:mexGetData不会为您提供数据的副本。它会给你一个指向数据的指针。如果修改了,MATLAB中的变量也会被修改,注意!

【讨论】:

以上是关于如何在Mexfile中的matlab(矩阵,单元格)和c ++(向量或自定义类)之间正确转换变量的主要内容,如果未能解决你的问题,请参考以下文章

在 Matlab 中访问单元格中的向量

Matlab:在for循环中排除零值矩阵单元的if条件

对具有相同 zoneID 的相邻单元格进行分组,MATLAB

Matlab - 在邻域中应用函数

定期从 mexFile 向 MATLAB 发送数据

在matlab中连接单元格数组的向量