4

This question shows how to return a double matrix to MATLAB using a map object.

The following works for non-complex data:

double *outputPtr;
plhs[0] = mxCreateDoubleMatrix((mwSize)n, (mwSize)m, mxREAL);
outputPtr = mxGetPr(plhs[0]);
Map<MatrixXd> output(outputPtr, n, m);

Since MATLAB stores real and complex elements of a matrix separately rather than interleaved, I don't think you can directly map to MatrixXcd. This is what I'm trying but it doesn't write to the data at plhs[0]:

double *outputReal, *outputImag;
plhs[0] = mxCreateDoubleMatrix((mwSize)n, (mwSize)m, mxCOMPLEX);
outputReal = mxGetPr(plhs[0]);
outputImag = mxGetPi(plhs[0]);
MatrixXcd outputMat(n,m);
outputMat.real() = Map<MatrixXd>(outputReal,n,m);
outputMat.imag() = Map<MatrixXd>(outputImag,n,m);

I think The problem with that is that it's copying the data in those maps to the real and imaginary parts of outputMat. So not only does this do an unnecessary copy on initializing but changes to a matrix assigned this way won't touch the data at the output pointer.

Is there a way to initialize the output matrix in such a way that its real and imaginary data will be stored at mxGetPr(plhs[0]) and mxGetPi(plhs[0])?

Cris Luengo
  • 55,762
  • 10
  • 62
  • 120
thorbjorn88
  • 71
  • 4
  • Eigen currently only supports interleaved storage for complex values, so you will need to set up two separate output Maps. – chtz Apr 14 '21 at 18:17
  • 2
    Since MATLAB R2018a complex numbers are stored interleaved. You need to compile your MEX-files with the `-R2018a` option, otherwise complex data will be copied into separate real and imaginary arrays for backwards compatibility. https://www.mathworks.com/help/matlab/matlab_external/matlab-support-for-interleaved-complex.html – Cris Luengo Apr 15 '21 at 04:49

1 Answers1

3

Thanks to Cris Luengo for pointing me in the right direction. Here's what I got to work.

plhs[0] = mxCreateDoubleMatrix((mwSize)n, (mwSize)m, mxCOMPLEX);
auto* Poutput = reinterpret_cast<complex<double>*>(mxGetComplexDoubles(plhs[0]));
Map<MatrixXcd> output(Poutput, n, m);

Then when compiling make sure to give mex the -R2018a flag which will build for interleaved complex data.

RHertel
  • 23,412
  • 5
  • 38
  • 64
thorbjorn88
  • 71
  • 4