关于袖带 R2C 和 C2R

Question

我已经使用袖带进行研究，但是使用它时遇到了一些问题。我的步骤如下：

1. 使用 R2C 对图像进行正向 FFT
2. 将核系数与复数结果相乘
3. 使用 C2R 对乘法结果进行逆 FFT

但是，当我使用复数结果乘核时，出现了一个严重的问题，cufft复数结果不等于fftw的结果，结果中有很多零。我知道 R2C 的结果大小是 N1(N2/2+1)，但我想得到完整的复杂结果。如何解决这个问题呢？即如何恢复 R2C 结果？以及如何将相乘的结果放入C2R并得到正确答案？

我的实现程序代码如下：

__global__ void MultiplyKernel(cufftComplex *data, float *data1,cufftComplex *data2, unsigned vectorSize) {
    unsigned idx = blockIdx.x*blockDim.x+threadIdx.x;
    if (idx < vectorSize){
        data[idx].x = data2[idx].x*data1[idx];
        data[idx].y = data2[idx].y*data1[idx];
    }
}

__global__ void Scale(cufftReal *data, unsigned vectorSize) {
    unsigned idx = blockIdx.x*blockDim.x+threadIdx.x;
    if (idx < vectorSize){
        data[idx] = data[idx]/vectorSize;
    }
}

void ApplyKernel1(cufftReal *data2, float *ImageBuffer, float *KernelBuffer, unsigned int NX, unsigned int NY,unsigned int NZ)
{
      float *Akernel;
      cufftComplex *data_dev1, *data_dev2;
      cufftReal *data_dev3, *data_dev;
      cudaMalloc((void **)&Akernel, NX * NY * NZ * sizeof(float));
      cudaMalloc((void **)&data_dev3, NX * NY * NZ * sizeof(cufftReal));
      cudaMalloc((void **)&data_dev, NX * NY * NZ * sizeof(cufftComplex));
      cudaMalloc((void **)&data_dev1, NX * NY * NZ * sizeof(cufftComplex));
      cudaMalloc((void **)&data_dev2, NX * NY * NZ * sizeof(cufftComplex));
      cudaMemset(data_dev, 0, NX * NY * NZ * sizeof(cufftReal));
      cudaMemset(data_dev1, 0, NX * NY * NZ * sizeof(cufftComplex));
      cudaMemset(data_dev2, 0, NX * NY * NZ * sizeof(cufftComplex));
      //cufftComplex *resultFFT = (cufftComplex*)malloc(NX * NY * NZ * sizeof(cufftComplex));
      //cufftReal *resultIFFT = (cufftReal*)malloc(NX * NY * NZ * sizeof(cufftReal));

      cudaMemcpy(data_dev, ImageBuffer, NX * NY * NZ * sizeof(cufftReal), cudaMemcpyHostToDevice);

      cufftHandle plan;
      cufftPlan3d(&plan, NZ, NY, NX, CUFFT_R2C);
      cufftExecR2C(plan, data_dev, data_dev1);

      //Multiply kernel
      cudaMemcpy(Akernel, KernelBuffer, NX * NY * NZ * sizeof(float), cudaMemcpyHostToDevice);
      static const int BLOCK_SIZE = 1000;
      const int blockCount = (NX*NY*NZ+BLOCK_SIZE-1)/BLOCK_SIZE;
      MultiplyKernel <<<blockCount, BLOCK_SIZE>>> (data_dev2, Akernel, data_dev1, NX*NY*NZ);


      cufftDestroy(plan);
      //cufftPlan3d(&plan, NZ, NY, NX, CUFFT_C2R);
      cufftPlan3d(&plan, NZ,NY,NX, CUFFT_C2R);
      cufftExecC2R(plan, data_dev2, data_dev3 );
      Scale <<<blockCount, BLOCK_SIZE>>> (data_dev3, NX*NY*NZ);
      cudaMemcpy(data2, data_dev3, NZ * NY * NX * sizeof(cufftReal), cudaMemcpyDeviceToHost);

      cufftDestroy(plan);
      cudaFree(data_dev);
      cudaFree(data_dev1);
      cudaFree(data_dev2);
      cudaFree(data_dev3);
      cudaFree(Akernel);


}

Answer 1

当您将 R2C fft 的结果乘以复数时，结果不再对应于对称数组。