2

I implemented a fftw (fftw.org) example to use Fast Fourier transforms... This is the code....

I load an image that I convert from uint8_t to double (this code works fine...).

string bmpFileNameImage = "files/testDummyFFTWWithWisdom/onechannel_image.bmp"; 
BMPImage bmpImage(bmpFileNameImage);
vector<double>pixelColors; 
vector<uint8_t> image = bmpImage.copyBits();

toDouble(image,pixelColors,256,256, 1); 
int width = bmpImage.width();
int height = bmpImage.height();

I use wisdom files to improve the performance

FILE * file = fopen("wisdom.fftw", "r");
if (file) {
    fftw_import_wisdom_from_file(file);
    fclose(file);
} 

///*  fftw variables  */
fftw_complex *out;

double *wisdomInput = (double *) fftw_malloc(sizeof(double)*width*2*(height/2 +1 ));

const fftw_plan forward =fftw_plan_dft_r2c_2d(width,height, wisdomInput,reinterpret_cast<fftw_complex *>(wisdomInput),FFTW_PATIENT);
const fftw_plan inverse = fftw_plan_dft_c2r_2d(width, height,reinterpret_cast<fftw_complex *>(wisdomInput),wisdomInput, FFTW_PATIENT);

file = fopen("wisdom.fftw", "w");

if (file) {
    fftw_export_wisdom_to_file(file);
    fclose(file);
}

Finally, I execute the fftw library.... I receive an Access violation error with the first function (fftw_execute_dft_r2c) and I don't know why... I read this tutorial: http://www.fftw.org/fftw3_doc/Multi_002dDimensional-DFTs-of-Real-Data.html#Multi_002dDimensional-DFTs-of-Real-Data. I do a malloc with (ny/2+1) how it is explained.... . I don't understand why it is not working.... I am testing different sizes...

out = (fftw_complex*) fftw_malloc(sizeof(fftw_complex) * width *(height / 2 + 1));
double *result =(double *)fftw_malloc(width * (height+2)  * sizeof(double));

fftw_execute_dft_r2c(forward,&pixelColors[0],out);
fftw_execute_dft_c2r(inverse,out,result);

Regards.

carlos.baez
  • 1,063
  • 2
  • 11
  • 31

3 Answers3

1

This is the corrected code. It had a few mistakes:

  • It was reading a wrong wisdom.fftw file (from some old test...). Now, It always creates a new fftw_plan and a new file.
  • I misunderstood how it works the fftw library with in-place and out-of-place parameters. I had to change mallocs for the correct padding for "in-place" (I added +2 in malloc functions).
  • In order to restore the image, I had to divide by its size ((width+2) * height) how it is explained in this link.

`

/* load image */
string bmpFileNameImage = "files/polyp.bmp";

BMPImage bmpImage(bmpFileNameImage);

int width = bmpImage.width();
int height = bmpImage.height();

vector<double> pixelColors;
vector<uint8_t> image = bmpImage.copyBits();
//get one channel from the image
Uint8ToDouble(image,pixelColors,bmpImage.width(),bmpImage.height(),1);

//We don't reuse old wisdom.fftw... It can be corrupt 
/*
FILE * file = fopen("wisdom.fftw", "r");
if (file) {
    fftw_import_wisdom_from_file(file);
    fclose(file);
} */

double *wisdomInput = (double *) fftw_malloc(sizeof(double)*height*(width+2));

const fftw_plan forward =fftw_plan_dft_r2c_2d(width,height,wisdomInput,reinterpret_cast<fftw_complex *>(wisdomInput),FFTW_PATIENT);
const fftw_plan inverse = fftw_plan_dft_c2r_2d(width,height,reinterpret_cast<fftw_complex *>(wisdomInput),wisdomInput, FFTW_PATIENT);
double *bitsColors =(double *)fftw_malloc((width) * height * sizeof(double));

for (int y = 0; y < height; y++) {
    for (int x = 0; x < width+2; x++) {
        if (x < width) {
            int currentIndex = ((y * width) + (x));
            bitsColors[currentIndex] = (static_cast<double>(result[y * (width+2) + x])) / (height*width);
        }
    }
}
fftw_free (wisdomInput);
fftw_free (out);
fftw_free (result);
fftw_free (bitsColors);

fftw_destroy_plan(forward);
fftw_destroy_plan(inverse);
fftw_cleanup();
}

`

Community
  • 1
  • 1
carlos.baez
  • 1,063
  • 2
  • 11
  • 31
0
fftw_execute_dft_r2c(forward,&pixelColors[0],out);

What are you doing here ? The array has already a pointer. Change it to fftw_execute_dft_r2c(forward,pixelColors[0],out); it should work now.

Blacktempel
  • 3,935
  • 3
  • 29
  • 53
  • It is not an array, It is a vector. I am getting the pointer for the first element... [related](http://stackoverflow.com/questions/2923272/how-to-convert-vector-to-array-c) – carlos.baez Jun 25 '13 at 10:00
  • `fftw_execute_dft_r2c(forward,&pixelColors[0],*out);` As I'm not sure what fttw_complex is used for, in this code. – Blacktempel Jun 25 '13 at 10:15
  • This is the method definition `void fftw_execute_dft_r2c( const fftw_plan p, double *in, fftw_complex *out);` My idea is to use it with the inverse method (fftw_execute_dft_c2r) in order to restore the image and compare with the original image. – carlos.baez Jun 25 '13 at 12:02
0

Maybe the problem is here (http://www.fftw.org/doc/New_002darray-Execute-Functions.html):

[...] that the following conditions are met:

  • The input and output arrays are the same (in-place) or different (out-of-place) if the plan was originally created to be in-place or out-of-place, respectively.

In the plan you are using in-place transformation parameters (with bad allocation, BTW, since:

double *wisdomInput = (double *) fftw_malloc(sizeof(double)*width*2*(height/2 +1 ));

should be:

double *wisdomInput = (double *) fftw_malloc(sizeof(fftw_complex)*width*2*(height/2 +1 ));

to be suitable for output too).

But you're calling fftw_execute_dft_r2c function with out-of-place parameters.

Gonmator
  • 760
  • 6
  • 15