How to convert a 4D matrix from Matlab to C++

I am trying to convert a 4D double matrix from Matlab to C++, using the below code however, it looks to me that mxGetPr() does not support linear mapping! are there any mapping indexing formula that can be used with this function (i.e ptr1[linear indexing formula] in the code below)?

MATFile *pmat;
int dt1,dt2,dt3,dt4;
  dt1=3; dt2=4; dt3=5;  dt4=6; % diemntions of the 4D matrix
    vector<vector<vector<vector<double>>>> array4D1(dt1, vector<vector<vector<double>>>
                                           (dt2, vector<vector<double>>
                                           (dt3, vector<double>
 pa2 = matGetVariable(pmat, "G_wp");    % G_wp is a 3x4x5x6 matlab matrix                                         (dt4))));
ptr1 = mxGetPr( pa2 ); 
int idx1=0;
for (int i=0; i<dt1; i++ )
  {for (int j=0; j<dt2; j++)
    {for (int k=0; k<dt3; k++)
      {for (int l=0; l<dt4; l++)
                           {array4D1[i][j][k][l]=ptr1[idx1];idx1++;}
                      }
              }
      }

Answer 1

[James Tursa]

C array storage order is reversed from MATLAB storage order. This holds regardless of the number of dimensions. And you certainly can use "linear indexing" on the MATLAB data memory. E.g., a 2D example:

#include "mex.h"
void mexFunction(int nlhs, mxArray *plhs[], int nrhs, const mxArray *prhs[])
{
/* Assume 2x3 double input, no argument checking */
    double x[3][2]; /* C array dimensions are reverse order from MATLAB */
    double *pr;
    mwSize i, j, m, n;

      pr = mxGetPr(prhs[0]);
      m = mxGetM(prhs[0]); /* 2 */
      n = mxGetN(prhs[0]); /* 3 */
  /* for loops & limits in reverse index order */
      for( j=0; j<n; j++ ) {
          for( i=0; i<m; i++ ) {
              x[j][i] = *pr++;
          }
      }
      for( i=0; i<m; i++ ) {
          for( j=0; j<n; j++ ) {
              mexPrintf("x[%d][%d] = %f\n",i,j,x[j][i]);
          }
      }
      mexCallMATLAB(0,NULL,1,prhs,"display");
  }

And at the command line:

>> mex matrix2d.c
>> x = [1 2 3;4 5 6]
x =
     1     2     3
     4     5     6
>> matrix2d(x)
x[0][0] = 1.000000
x[0][1] = 2.000000
x[0][2] = 3.000000
x[1][0] = 4.000000
x[1][1] = 5.000000
x[1][2] = 6.000000

x =
     1     2     3
     4     5     6

So, everything is as expected as long as you get the reverse order indexing coded correctly.

In your code above, you don't post enough information for me to be able to tell what the dimensions of the mxArray pa2 are, you simply have hard coded the following:

   dt1=3; dt2=4; dt3=5;  dt4=6;

So I can't say for sure how to correct your code. But if you just extend my 2D example to 4D it should work for you. Just do the dimensions on the C++ side in reverse order from MATLAB.

EDIT 7/26/2016:

E.g., here is some C++ code to read a variable from a mat file and use your 4D vector approach:

#include <vector>
using namespace std;
#include "mex.h"
#include "mat.h"
void mexFunction(int nlhs, mxArray *plhs[], int nrhs, const mxArray *prhs[])
{
    MATFile *pmat;
    const mwSize *dims;
    mwSize ndim;
    mwSize dt1, dt2, dt3, dt4;
    mxArray *pa2;
    double *ptr1;
    int a,b,c,d;

      pmat = matOpen("matGwp.mat","r");
      if( !pmat ) {
          mexErrMsgTxt("Cannot open matGwp.mat");
      }
      pa2 = matGetVariable(pmat, "G_wp");
      matClose(pmat);
      if( pa2 == NULL ) {
          mexErrMsgTxt("Cannot get G_wp");
      }
      dims = mxGetDimensions(pa2);
      ndim = mxGetNumberOfDimensions(pa2);
      if( ndim != 4 ) {
          mexErrMsgTxt("G_wp needs to have 4 explicit dimensions");
      }
  // Set up array dimensions in reverse order
      dt4 = dims[0];
      dt3 = dims[1];
      dt2 = dims[2];
      dt1 = dims[3];
      vector<vector<vector<vector<double>>>> array4D1(dt1, vector<vector<vector<double>>>
                                                     (dt2, vector<vector<double>>
                                                     (dt3, vector<double>
                                                     (dt4)               )       )       );
  // Copy elements in reverse order
      ptr1 = mxGetPr( pa2 ); 
      int idx1=0;
      for (int i=0; i<dt1; i++ ) {
          for (int j=0; j<dt2; j++) {
              for (int k=0; k<dt3; k++) {
                  for (int l=0; l<dt4; l++) {
                      array4D1[i][j][k][l] = ptr1[idx1];
                      idx1++;
                  }
              }
          }
      }
  // Some test output
      a = 1; b = 1; c = 1; d = 1;
      mexPrintf("G_wp(%d,%d,%d,%d) = array4D1[%d][%d][%d][%d] = %f\n",a,b,c,d,d-1,c-1,b-1,a-1,array4D1[d-1][c-1][b-1][a-1]);
      a = 2; b = 4; c = 1; d = 6;
      mexPrintf("G_wp(%d,%d,%d,%d) = array4D1[%d][%d][%d][%d] = %f\n",a,b,c,d,d-1,c-1,b-1,a-1,array4D1[d-1][c-1][b-1][a-1]);
      a = 3; b = 2; c = 3; d = 2;
      mexPrintf("G_wp(%d,%d,%d,%d) = array4D1[%d][%d][%d][%d] = %f\n",a,b,c,d,d-1,c-1,b-1,a-1,array4D1[d-1][c-1][b-1][a-1]);
      a = 3; b = 4; c = 5; d = 6;
      mexPrintf("G_wp(%d,%d,%d,%d) = array4D1[%d][%d][%d][%d] = %f\n",a,b,c,d,d-1,c-1,b-1,a-1,array4D1[d-1][c-1][b-1][a-1]);
  }

And here is the example using that code:

>> G_wp = reshape(1:3*4*5*6,3,4,5,6);
>> save matGwp
>> mex mat2vector.cpp
>> mat2vector
G_wp(1,1,1,1) = array4D1[0][0][0][0] = 1.000000
G_wp(2,4,1,6) = array4D1[5][0][3][1] = 311.000000
G_wp(3,2,3,2) = array4D1[1][2][1][2] = 90.000000
G_wp(3,4,5,6) = array4D1[5][4][3][2] = 360.000000
>> G_wp(1,1,1,1)
ans =
     1
>> G_wp(2,4,1,6)
ans =
   311
>> G_wp(3,2,3,2)
ans =
    90
>> G_wp(3,4,5,6)
ans =
   360

So, again, as long as you get the dimensions in reverse order you can make it look like the C example.

---------------------------------

If you want the dimensions to be in the same order in C++ as they are in MATLAB, and are willing to accept that they will be in a different order in memory, you can do this:

#include <vector>
using namespace std;
#include "mex.h"
#include "mat.h"
void mexFunction(int nlhs, mxArray *plhs[], int nrhs, const mxArray *prhs[])
{
    MATFile *pmat;
    const mwSize *dims;
    mwSize ndim;
    mwSize dt1, dt2, dt3, dt4;
    mxArray *pa2;
    double *ptr1;
    int a,b,c,d;

      pmat = matOpen("matGwp.mat","r");
      if( !pmat ) {
          mexErrMsgTxt("Cannot open matGwp.mat");
      }
      pa2 = matGetVariable(pmat, "G_wp");
      matClose(pmat);
      if( pa2 == NULL ) {
          mexErrMsgTxt("Cannot get G_wp");
      }
      dims = mxGetDimensions(pa2);
      ndim = mxGetNumberOfDimensions(pa2);
      if( ndim != 4 ) {
          mexErrMsgTxt("G_wp needs to have 4 explicit dimensions");
      }
  // Set up array dimensions in same order
      dt4 = dims[3];
      dt3 = dims[2];
      dt2 = dims[1];
      dt1 = dims[0];
      vector<vector<vector<vector<double>>>> array4D1(dt1, vector<vector<vector<double>>>
                                                     (dt2, vector<vector<double>>
                                                     (dt3, vector<double>
                                                     (dt4)               )       )       );
  // Copy elements in same order (but will be reversed in memory compared to MATLAB)
      ptr1 = mxGetPr( pa2 ); 
      int idx1=0;
      for (int l=0; l<dt4; l++) {
          for (int k=0; k<dt3; k++) {
              for (int j=0; j<dt2; j++) {
                  for (int i=0; i<dt1; i++ ) {
                      array4D1[i][j][k][l] = ptr1[idx1];
                      idx1++;
                  }
              }
          }
      }
  // Some test output (dimensions in same order)
      a = 1; b = 1; c = 1; d = 1;
      mexPrintf("G_wp(%d,%d,%d,%d) = array4D1[%d][%d][%d][%d] = %f\n",a,b,c,d,a-1,b-1,c-1,d-1,array4D1[a-1][b-1][c-1][d-1]);
      a = 2; b = 4; c = 1; d = 6;
      mexPrintf("G_wp(%d,%d,%d,%d) = array4D1[%d][%d][%d][%d] = %f\n",a,b,c,d,a-1,b-1,c-1,d-1,array4D1[a-1][b-1][c-1][d-1]);
      a = 3; b = 2; c = 3; d = 2;
      mexPrintf("G_wp(%d,%d,%d,%d) = array4D1[%d][%d][%d][%d] = %f\n",a,b,c,d,a-1,b-1,c-1,d-1,array4D1[a-1][b-1][c-1][d-1]);
      a = 3; b = 4; c = 5; d = 6;
      mexPrintf("G_wp(%d,%d,%d,%d) = array4D1[%d][%d][%d][%d] = %f\n",a,b,c,d,a-1,b-1,c-1,d-1,array4D1[a-1][b-1][c-1][d-1]);
  }

And the command line results are again as expected:

>> mex mat2vector2.cpp
>> mat2vector2
G_wp(1,1,1,1) = array4D1[0][0][0][0] = 1.000000
G_wp(2,4,1,6) = array4D1[1][3][0][5] = 311.000000
G_wp(3,2,3,2) = array4D1[2][1][2][1] = 90.000000
G_wp(3,4,5,6) = array4D1[2][3][4][5] = 360.000000

------------------------------------------------

CAUTION: All of this assumes that you want to go through all of this vector stuff simply in order to preserve the 4D coding [][][][] syntax. You do gain the syntax, but you lose other capabilities because of this. E.g., you can no longer simply pass the pointer to a 2D slice of a double matrix to another routine for processing (e.g., BLAS or LAPACK).

[alaa Afeef]

Thank you James Tursa for your answer, however, It did not work for 4D matrix. I managed to do a vectorised solution for the 2D matrix problem as below, but it still not working for the 4D case. can you please provide your code for the 4D case?

 pa2 = matGetVariable(pmat, "Q"); % Q : 90x90 matlab matrix
  Col = mxGetN(pa2);
  Row = mxGetM(pa2);
    printf("Col= %d, Row= %d\n", Col,Row);

vector<vector<double> > array2D;
// Set up sizes. (HEIGHT x WIDTH)
array2D.resize(HEIGHT);
for (int i = 0; i < HEIGHT; ++i)
  array2D[i].resize(WIDTH);
double* ptr1 = mxGetPr( pa2 );
for ( int i=0; i<90; i++ )
for (int j=0; j<90; j++)
          array2D[i][j]=ptr1[i*90+j];

[James Tursa]

First, can you tell me what you are trying to do? Simply copy the MATLAB variable values into a C++ 4D vector object? What are the dimensions of the MATLAB variable?

[alaa Afeef]

Thank you James Tursa, that is brilliant. Actually, what I am trying to do is to copy a large 4D matrix from a mat file into C++ in order to be passed to another routine for processing (cuFFT,LAPACK). I can see that passing a vector to such routines will not work, is there a way to convert the vector into pointer based matrix in C++? actually, I was trying to avoid the array or arrays approach in C++ to copy the Mat data. and looking for a computationally efficient way to handle it.

[James Tursa]

What slices of the 4D matrix are you planning to pass to these other routines (cuFFT and/or LAPACK)? How you intend to do the slicing will affect how you want to organize the data in memory at the C++ level. The easy answer is to use the 1st dimension of the MATLAB variable as your 1D vector slice since that data is already contiguous in memory. The harder answer is you want to use a dimension other than the 1st, in which case you will need to reorder the data in memory before passing pointers to cuFFT or LAPACK. So, the big question is then, what dimension of the MATLAB variable do you want to use for the 1D vector? Maybe you could show, for your small 3x4x5x6 example, which elements would actually form your 1D vector that you want to do further processing with?