# Cufft\_R2C and Cufft\_C2R are inaccurate 

**URL:** <https://forums.developer.nvidia.com/t/cufft-r2c-and-cufft-c2r-are-inaccurate/32810>\
**Category:** GPU-Accelerated Libraries\
**Created:** [April 7, 2014, 4:03am UTC](https://forums.developer.nvidia.com/t/cufft-r2c-and-cufft-c2r-are-inaccurate/32810 "2014-04-07T04:03:09Z")\
**Posts on this page:** 3\
**Page:** 1

<div class="post-metadata">

**Author:** ![kwang16](https://developer.download.nvidia.com/images/forums/profile-default-devtalk-84.png) [@kwang16](https://forums.developer.nvidia.com/u/kwang16)\
**Post date:** [April 7, 2014, 4:03am UTC](https://forums.developer.nvidia.com/t/cufft-r2c-and-cufft-c2r-are-inaccurate/32810/1 "2014-04-07T04:03:09Z")

</div>

I described my problem here: [Instability of CUFFT\_R2C and CUFFT\_C2R | Medical Imaging Solution](http://imagingsci.wordpress.com/2014/04/07/instability-of-cufft-r2c-and-c2r/)

My testing codes for ifft (C2R) are attached.

using namespace std;

#include \<stdio.h\>  
#include   
#include \<math.h\>  
#include \<stdlib.h\>  
#include \<string.h\>

#include \<cufft.h\>

#define INFILE “x.DAT”  
#define OUTFILE1 “X.DAT”  
#define OUTFILE2 “xx.DAT”

#define NO\_x1 (1024)  
#define NO\_x2 (1024)  
#define NO\_X1 (1024)// Run faster  
#define NO\_X2 (1024)

int main(int argc, char\*\*argv){  
int device,nx1,nx2;  
FILE \*fp;  
if (cudaSetDevice(0)!=cudaSuccess){  
cout\<\<“Error when initializing device2!”\<\<endl;  
exit(-1);  
}  
cudaGetDevice(&device);  
cout\<\<“I am occupying device-”\<\<device\<\<"; "\<\<endl;

float _x = (cufftReal_)malloc(NO\_x1_NO\_x2_sizeof(float));  
cufftReal _d\_x;  
cudaMalloc((void_\*)&d\_x,sizeof(cufftReal)_NO\_x1_NO\_x2);  
cufftComplex _d\_X;  
cudaMalloc((void_\*)&d\_X,sizeof(cufftComplex)_NO\_X1_NO\_X2);  
cufftComplex _X = (cufftComplex )malloc(NO\_X1NO\_X2_sizeof(cufftComplex));  
cufftHandle plan\_cr;  
if (cufftPlan2d(&plan\_cr, NO\_x1, NO\_x2, CUFFT\_C2R) != CUFFT\_SUCCESS){  
cout\<\<“CUFFT error: Plan creation failed???”\<\<endl;  
return (1);  
}  
if((fp=fopen(“X.DAT”,“rb”))==NULL){  
printf(“Can not open file %s.\n”,INFILE);  
exit(1);  
}  
if (fread(X,sizeof(cufftComplex),NO\_X1_NO\_X2,fp)!=NO\_X1_NO\_X2){  
printf(“Can not read file &s.\n”,INFILE);  
exit(1);  
}  
fclose(fp);

cudaMemset(d\_X,0,sizeof(cufftComplex)_NO\_X1_NO\_X2); // This line is for testing only  
cudaMemcpy(d\_X,X,sizeof(cufftComplex)_513_512,cudaMemcpyHostToDevice);

cudaMemset(d\_x,0,sizeof(float)_NO\_x1_NO\_x2);  
if (cufftExecC2R(plan\_cr, d\_X, d\_x) != CUFFT\_SUCCESS){  
cout\<\<“CUFFT error: ExecR2C forward failed”\<\<endl;  
return;  
}

cudaMemcpy(x,d\_x,sizeof(float)_NO\_x1_NO\_x2,cudaMemcpyDeviceToHost);

if ((fp=fopen(“xestimate.DAT”,“wb”))==NULL){  
cout\<\<"Can not open file "\<\<OUTFILE1\<\<endl;  
exit(0);  
}  
fwrite(x,sizeof(float),NO\_x1\*NO\_x2,fp);  
fclose(fp);  
free(x);  
free(X);  
cudaFree(d\_X);cudaFree(d\_x);  
cufftDestroy(plan\_cr);  
cudaThreadExit();  
return(0);  
}

---

<div class="post-metadata">

**Author:** ![Moa53](https://developer.download.nvidia.com/images/forums/profile-default-devtalk-84.png) [@Moa53](https://forums.developer.nvidia.com/u/Moa53)\
**Post date:** [April 11, 2014, 9:11am UTC](https://forums.developer.nvidia.com/t/cufft-r2c-and-cufft-c2r-are-inaccurate/32810/2 "2014-04-11T09:11:07Z")

</div>

Hi,

here is a small code I made (based on yours, so if some stuffs look similar it’s normal ;) ). It’s a R2C transform, followed by a C2R one. Obviously, the final output is the same as the initial input. So you can see how to use that properly, because I think there’re few things that are not really clear. Or maybe I just didn’t understand the purpose of your code, I must admit that I’m not on arrow these days (in that case just ignore my message please ^^)

```auto
using namespace std;

#include <iostream>
#include <cufft.h>
#include <stdio.h>

#define INFILE "x.in"
#define OUTFILE1 "complex_x.out"
#define OUTFILE2 "x.out"
#define N1 1024
#define N2 1024

int main(int argc, char **argv){

int real_size = N1*N2;
int complex_size= N1*(N2/2+1); //Output memory is N1(N/2+1), not (N1/2)(N2/2+1)

cufftReal *real_x,*d_real_x;
cufftComplex *complex_x,*d_complex_x;

cufftHandle plan_r2c,plan_c2r;

FILE *fp;

//Creation of the plans
if(cufftPlan2d(&plan_r2c,N1,N2,CUFFT_R2C) != CUFFT_SUCCESS){
        cout<<"CUFFT ERROR : R2C plan creation failed"<<endl;
        exit(-1);
}

if(cufftPlan2d(&plan_c2r,N1,N2,CUFFT_C2R) != CUFFT_SUCCESS){
        cout<<"CUFFT ERROR : C2R plan creation failed"<<endl;
        exit(-1);
}

//Memory allocation
real_x = (cufftReal*)malloc(real_size*sizeof(cufftReal));
cudaMalloc((void**)&d_real_x,real_size*sizeof(cufftReal));

complex_x = (cufftComplex*)malloc(complex_size*sizeof(cufftComplex));
cudaMalloc((void**)&d_complex_x,complex_size*sizeof(cufftComplex));

//Data intialisation
for(int ii=0; ii<real_size; ii++) real_x[ii]=(double)ii/2.0; //Not the courage to put random stuff :p

if((fp=fopen(INFILE,"w"))==NULL){
        cout<<"Can't open file "<<INFILE<<endl;
        exit(-1);
}

for (int ii=0;ii<real_size;ii++) fprintf(fp,"%.1f\n",real_x[ii]); //Just in order to have a file for comparison

fclose(fp);

//FFT real to complex *************************************************
cudaMemcpy(d_real_x,real_x,real_size*sizeof(cufftReal),cudaMemcpyHostToDevice);

if(cufftExecR2C(plan_r2c,d_real_x,d_complex_x) != CUFFT_SUCCESS){
        cout<<"CUFFT error : ExecR2C failed"<<endl;
        exit(-1);
}

cudaMemcpy(complex_x,d_complex_x,complex_size*sizeof(cufftComplex),cudaMemcpyDeviceToHost);
// ********************************************************************

if((fp=fopen(OUTFILE1,"w"))==NULL){
        cout<<"Can't open file "<<OUTFILE1<<endl;
        exit(-1);
}

for (int ii=0;ii<complex_size;ii++) fprintf(fp,"%f %f\n",complex_x[ii].x,complex_x[ii].y);

fclose(fp);

//FFT complex to real+++++++++++++++++++++++++++++++++++++++++++++++++
cudaMemcpy(d_complex_x,complex_x,complex_size*sizeof(cufftComplex),cudaMemcpyHostToDevice);

if(cufftExecC2R(plan_c2r,d_complex_x,d_real_x) != CUFFT_SUCCESS){
        cout<<"CUFFT error : ExecC2R failed"<<endl;
        exit(-1);
}

cudaMemcpy(real_x,d_real_x,real_size*sizeof(cufftReal),cudaMemcpyDeviceToHost);
//++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++

if((fp=fopen(OUTFILE2,"w"))==NULL){
        cout<<"Can't open file "<<OUTFILE2<<endl;
        exit(-1);
}

for (int ii=0;ii<real_size;ii++) fprintf(fp,"%.1f\n",real_x[ii]/real_size); //cuFFT is unnormalized, you have to divise by N1*N2 to get the inital values

fclose(fp);

free(real_x);
free(complex_x);
cudaFree(d_real_x);
cudaFree(d_complex_x);

cufftDestroy(plan_r2c);
cufftDestroy(plan_c2r);

return(0);

}

```

PS : I’m french, sorry for my english .\_.

---

<div class="post-metadata">

**Author:** ![pasoleatis](https://developer.download.nvidia.com/images/forums/profile-default-devtalk-84.png) [@pasoleatis](https://forums.developer.nvidia.com/u/pasoleatis)\
**Post date:** [April 11, 2014, 3:36pm UTC](https://forums.developer.nvidia.com/t/cufft-r2c-and-cufft-c2r-are-inaccurate/32810/3 "2014-04-11T15:36:40Z")

</div>

In my code I have to define the real matrix of size [Nx][Ny][Nz+2], the last 2 elements coresponding to Nz and Nz+1 are just padding so that one can do a inplace transform.  
The complex array is of size [Nx][Ny][Nz/2+1].  
You have to use this or set the compatibility mode in a such way that the padding is disabled. The same applies for 2D arrays.

I posted this message without checking if there was a reply meanwhile. :)