Displaying a cvMatrix containing complex numbers (CV_64FC2)

Question

Hi,

I'm new to OpenCV, and I would like to compare the results of a python program with my calculations in OpenCV. My matrix contains complex numbers since its the result of a cvDFT. Python handles complex numbers well and displays it with scientific notation. My C++ program is not effective when trying to use std::cout.

I tried to store my numbers array in a std::complex[] instead of a double[] but it is not compiling.

Here is my code, and its result :

    CvMat *dft_A;

    dft_A = cvCreateMat(5, 5, CV_64FC2); // complex matrix
    double a[] = {
        0, 0, 0, 0, 0,
            1, 1, 1, 1, 1,
            2, 2, 2, 2, 2,
            3, 3, 3, 3, 3,
            4, 4, 4, 4, 4
           };
    dft_A->data.db = a;
    std::cout << "before : " << a[0] << std::endl;
    cvDFT( dft_A, dft_A, CV_DXT_FORWARD);  // DFT !
    std::cout << "after : " << a[0] << std::endl;


        >> before : 0

Here is the same in python, with the output :

>>> a = np.mgrid[:5, :5][0]
>>> a
array([[0, 0, 0, 0, 0],
       [1, 1, 1, 1, 1],
       [2, 2, 2, 2, 2],
       [3, 3, 3, 3, 3],
       [4, 4, 4, 4, 4]])
>>> np.fft.fft2(a)
array([[ 50.0 +0.j        ,   0.0 +0.j        ,   0.0 +0.j        ,
          0.0 +0.j        ,   0.0 +0.j        ],
       [-12.5+17.20477401j,   0.0 +0.j        ,   0.0 +0.j        ,
          0.0 +0.j        ,   0.0 +0.j        ],
       [-12.5 +4.0614962j ,   0.0 +0.j        ,   0.0 +0.j        ,
          0.0 +0.j        ,   0.0 +0.j        ],
       [-12.5 -4.0614962j ,   0.0 +0.j        ,   0.0 +0.j        ,
          0.0 +0.j        ,   0.0 +0.j        ],
       [-12.5-17.20477401j,   0.0 +0.j        ,   0.0 +0.j        ,
          0.0 +0.j        ,   0.0 +0.j        ]])
>>>

The problem is obviously coming from the second cout which is inefficient with the type of date (CV_64FC2 for complex number).

My question is : how can I dump the result so I can check that my python code is doing the same as my cpp/opencv code ?

Thanks !

Answer 1

Have you tried the python bindings for OpenCV? http://www.exothermia.net/monkeys_and_robots/2009/12/11/working-opencv-python-bindings/

With the bindings you can call OpenCV functions from python and have the results as numpy arrays and then compare them with those from your pure python code. With some tinkering you can wrap your own C code and have it available in python as well.

But if you only want to dump the data, you can probably save the real and imaginary parts as images and read them in python (i'm not very familiar with OpenCV, you have to check its support - and python's - for float images).

Answer 2

There is a dft example in OpenCV 2.0 code, which I am also studying right now. Here is a copy paste for you that might give you an idea. As you can see, it uses cvSplit to spilit to real and imaginary components. Hope that helps:

im = cvLoadImage( filename, CV_LOAD_IMAGE_GRAYSCALE );
if( !im )
    return -1;

realInput = cvCreateImage( cvGetSize(im), IPL_DEPTH_64F, 1);
imaginaryInput = cvCreateImage( cvGetSize(im), IPL_DEPTH_64F, 1);
complexInput = cvCreateImage( cvGetSize(im), IPL_DEPTH_64F, 2);

cvScale(im, realInput, 1.0, 0.0);
cvZero(imaginaryInput);
cvMerge(realInput, imaginaryInput, NULL, NULL, complexInput);

dft_M = cvGetOptimalDFTSize( im->height - 1 );
dft_N = cvGetOptimalDFTSize( im->width - 1 );

dft_A = cvCreateMat( dft_M, dft_N, CV_64FC2 );
image_Re = cvCreateImage( cvSize(dft_N, dft_M), IPL_DEPTH_64F, 1);
image_Im = cvCreateImage( cvSize(dft_N, dft_M), IPL_DEPTH_64F, 1);

// copy A to dft_A and pad dft_A with zeros
cvGetSubRect( dft_A, &tmp, cvRect(0,0, im->width, im->height));
cvCopy( complexInput, &tmp, NULL );
if( dft_A->cols > im->width )
{
    cvGetSubRect( dft_A, &tmp, cvRect(im->width,0, dft_A->cols - im->width, im->height));
    cvZero( &tmp );
}

// no need to pad bottom part of dft_A with zeros because of
// use nonzero_rows parameter in cvDFT() call below

cvDFT( dft_A, dft_A, CV_DXT_FORWARD, complexInput->height );

cvNamedWindow("win", 0);
cvNamedWindow("magnitude", 0);
cvShowImage("win", im);

// Split Fourier in real and imaginary parts
cvSplit( dft_A, image_Re, image_Im, 0, 0 );

// Compute the magnitude of the spectrum Mag = sqrt(Re^2 + Im^2)
cvPow( image_Re, image_Re, 2.0);
cvPow( image_Im, image_Im, 2.0);
cvAdd( image_Re, image_Im, image_Re, NULL);
cvPow( image_Re, image_Re, 0.5 );

// Compute log(1 + Mag)
cvAddS( image_Re, cvScalarAll(1.0), image_Re, NULL ); // 1 + Mag
cvLog( image_Re, image_Re ); // log(1 + Mag)