Segmentation fault when recoving data from a C code through Mex programming

Hi everyone,

I'm working on a C code that I compile on Matlab by the mex engine and then I want to recover the result on Matlab, here are the lines that I type on Matlab :

mex -g -largeArrayDims MyFile.c
SparseMatrix=MyFile();

After the execution the variable SparseMatrix is created in the matlab window reserved for variables but when I try to access it I get a segmentation fault and Matlab shuts down and I don't know where the error comes from because the compilation and the execution of the C file goes well.

Here's my C code:

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

    MATFile *pmat_in;
    mxArray *pa_in;
    const char **name_in;

    /* Open file to get directory*/
    pmat_in=matOpen("MatFile.mat", "r");

    if (pmat_in == NULL)
    {
        printf("Error creating file : MatFile \n");
        return(EXIT_FAILURE);
    }

    pa_in = matGetNextVariable(pmat_in, &name_in);
    if (pa_in == NULL)
    {
        printf("Error reading in file \n");
        return(EXIT_FAILURE);
    }

    if (matClose(pmat_in) != 0)
    {
        printf("Error closing file \n");
        return(EXIT_FAILURE);
    }

    mxArray* filled_sparse_matrix;
    filled_sparse_matrix=MyFunction(arguments);

    plhs[0]=filled_sparse_matrix;
    return;   
}

Does anyone have any idea of where the segmentation fault might come from please? I would appreciate your help very much.

Thanks in advance.

[James Tursa]

These lines probably aren't really doing what you think they are, since the mexFunction return is void.

        return(EXIT_FAILURE);

Instead, you should call the mex error function to return to the caller gracefully. E.g.,

        mexErrMsgTxt("Error creating file : MatFile \n");

That being said, the coding error causing the crash is probably in MyFunction, and you haven't shown us that code. Likely the sparse matrix is not being built properly, but we can't tell what exactly is wrong without seeing the code.

[JM]

Thank you James for your comment, I replaced the return statement with the error message as you suggested and here's what MyFunction does:

mxArray *MyFunction(mxArray *pa_in, mwSize *array_rad, double alpha, double gamma, double seuil)
{

    const mwSize *dim_image;

    if (pa_in == NULL)
    {
        mexErrMsgTxt("pa_in is NULL \n");
    }

    double *pr_in=mxGetPr(pa_in);

    if (pr_in == NULL)
    {
        mexErrMsgTxt("pr_in is NULL \n");
    }

    mwSize array_rad_at=array_rad[0];
    mwSize array_rad_ec=array_rad[1];

    /* The dimensions of the image before transformation*/
    dim_image=mxGetDimensions(pa_in);

    mwSize nr_new=dim_image[0];
    mwSize nc_new=dim_image[1];
    mwSize b_new=dim_image[2];

    mwSize nr=nr_new-(2*array_rad_at_);
    mwSize nc=nc_new-(2*array_rad_at_);
    mwSize b=b_new-(2*array_rad_ec);

    mwSize nzmax=nr*nr*nc*nc*b*b;

    /*creatino of the sparse matrix*/
    mxArray *sparse_matrix=mxCreateSparse(nr*nc*b,nr*nc*b,nzmax,mxREAL);

    /*testing if the matrix was created successfully*/
    if (sparse_matrix== NULL)
    {
        mexErrMsgTxt("The sparse matrix was not created successfully \n");
    }

    double *pr_out, *ir_out, *jc_out, idx_ir=0,idx_jc=1, val_jc=0,cpt_jc=0, idx_col,idx_row, cpt_col=-1;

    pr_out = mxGetPr (sparse_matrix);
    ir_out = mxGetIr (sparse_matrix);
    jc_out = mxGetJc (sparse_matrix);
    jc_out[0]=0, jc_out++;

    /*computing the radiuses */
    double r=-gamma*log(seuil);
    double r_x=floor(nr*sqrt(r));
    double r_y=floor(nc*sqrt(r));
    double r_z=floor(b*sqrt(r/alpha));

    int d=(2*array_rad_at_+1)*(2*array_rad_at_+1)*(2*array_rad_ec+1);

    double i,j,f,i_prime,j_prime,f_prime,bound_start_f, bound_end_f,    bound_start_i, bound_end_i,bound_start_j,bound_end_j;
    mwIndex *idx_start_1, *idx_end_1, *idx_start_2, *idx_end_2, *tmp_pix_1,     *tmp_pix_2, cpt=0;

    for (f=array_rad_ec;f<b_new-array_rad_ec;f++)
    {

        bound_start_f=f-r_z;
        bound_end_f=f+r_z;

        if(bound_start_f < array_rad_ec)
        {
            bound_start_f =array_rad_ec;
        }

        if(bound_end_f > b_new-(array_rad_ec+1))
        {
            bound_end_f = b_new-(array_rad_ec+1);
        }

        for(i=array_rad_at_;i<nr_new-array_rad_at_;i++)
        {
            bound_start_i=i-r_x;
            bound_end_i=i+r_x;

            if (bound_start_i<array_rad_at_)
            {
                bound_start_i=array_rad_ec;
            }

            if(bound_end_i > nr_new-(array_rad_at_+1))
            {
                bound_end_i=nr_new-(array_rad_at_+1);
            }

            for(j=array_rad_at_;j<nc_new-array_rad_at_;j++)
            {
                bound_start_j=j-r_y;
                bound_end_j=j+r_y;

                if(bound_start_j < array_rad_at_)
                {
                    bound_start_j=array_rad_at_;
                }

                if(bound_end_j > nc_new-(array_rad_at_+1))
                {
                    bound_end_j= nc_new-(array_rad_at_+1);
                }

                for(f_prime=bound_start_f; f_prime<=bound_end_f;f_prime++)
                {
                    for(i_prime=bound_start_i; i_prime<=bound_end_i;i_prime++)
                    {
                        for(j_prime=bound_start_j; j_prime<=bound_end_j;j_prime++)
                        {

                            idx_row=(j_prime-array_rad_at_)+(i_prime-      array_rad_at_)*nc + (f_prime-array_rad_ec)*nr*nc;

                            *ir_out=idx_row;
                            *pr_out=idx_row;

                            ir_out++;
                            pr_out++;
                            cpt_jc++;

                        }
                    }
                }

                val_jc=val_jc+cpt_jc;
                *jc_out=val_jc;
                * jc_out++;

                cpt_jc=0;
                cpt++;

            }

        }
    }

    return sparse_matrix;  
}

For the moment, the ir_out and the sparse matrix are filled with the same exact values. This code was made as a first step to test the recovery of the created sparse matrix from Matlab.

Thank you in advance.

Answer 1

[James Tursa]

I haven't tried to run, or even compile, your code. But based on a very quick visual scan it appears to me that this line:

                * jc_out++;

should be this instead:

                jc_out++;

It is the jc analogous line to the ir and pr lines:

                            ir_out++;
                            pr_out++;

There should be no dereferencing in front of the jc_out++ just as there is no dereferencing in front of the ir_out++ and pr_out++ lines. (I don't think that initial * really hurts you since the jc_out++ will happen regardless, but neverthless the * really doesn't belong on that line).

I would put in a row index check as well. E.g.,

                            if( idx_row < 0 || idx_row >= nr*nc*b ) {
                                mexErrMsgTxt("Row index invalid in sparse matrix build.");
                            }
                            *ir_out=idx_row;

What is cpt used for?

Then there are the bounds checking in the j loop. I haven't had the time yet, but since the number of columns of the sparse matrix is nr*nc*b, I would expect the j loop to have this as the bounds limit but I see this:

            for(j=array_rad_at_;j<nc_new-array_rad_at_;j++)

So, this may be perfectly legit ... I am just pointing out that at first glance it is not clear if the j loop will iterate exactly nr*nc*b times. Something to check.

Then there is the nzmax calculation and the sparse matrix size:

    mwSize nzmax=nr*nr*nc*nc*b*b;
    mxArray *sparse_matrix=mxCreateSparse(nr*nc*b,nr*nc*b,nzmax,mxREAL);

If notionally the sparse matrix has size MxN, you have created the data area to hold MxN elements. This will result in absolutely no memory savings at all, and in fact will take up significantly more memory than a non-sparse full MxN matrix would. Is this just test code to make sure you didn't overrun the data area downstream, or is this really how your final code will be? Because you will take a memory hit with your current code.

The following code, while it doesn't hurt anything, really isn't needed at all:

    if (sparse_matrix== NULL)
    {
        mexErrMsgTxt("The sparse matrix was not created successfully \n");
    }

In a mex routine, if there isn't enough memory to create the sparse matrix then the mex routine will return to the caller immediately with an error. I.e., you will not even return from the mxCreateSparse function call back into your code in this event. The sparse_matrix variable will never be NULL in your downstream code.

Finally, you should put in more argument checks to make your code robust in the presence of unexpected inputs. E.g., for pa_in you need to check several things:

if( !mxIsDouble(pa_in) || mxIsComplex(pa_in) || mxIsSparse(pa_in) || mxGetNumberOfDimensions(pa_in) != 3 ) {
    mexErrMsgTxt("Input variable is not as expected.");
}

Failing any one of these checks could easily bomb the code downstream if you didn't check for them.

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

In addition, you can use the SPOK utility written by Tim Davis from the FEX:

https://www.mathworks.com/matlabcentral/fileexchange/20186-spok--checks-if-a-matlab-sparse-matrix-is-ok

*** EDIT ***

spok has been removed from the FEX. You can use my replacement function sarek instead:

https://www.mathworks.com/matlabcentral/fileexchange/158161-sarek-sparse-analyzer-real-et-komplex?s_tid=srchtitle