CUFFT輸出不正確

Question

我有這個程序有問題：

#include <stdlib.h> 
#include <stdio.h> 
#include <string.h> 
#include <math.h> 
#include <cufft.h> 
#include <cuComplex.h> 

#define SIGNAL_SIZE  1024 

int main(int argc, char **argv) { 
    cudaEvent_t start, stop; 
    cudaEventCreate(&start); 
    cudaEventCreate(&stop); 

    // Allocate host memory for the signal 
    cuDoubleComplex *h_signal = (cuDoubleComplex *) malloc(sizeof(cuDoubleComplex) * SIGNAL_SIZE); 

    // Initalize the memory for the signal 
    for (unsigned int i = 0; i < SIGNAL_SIZE; ++i) { 
     if((double)i/SIGNAL_SIZE>=0 && (double)i/SIGNAL_SIZE<0.5) h_signal[i].x = (double)i/SIGNAL_SIZE; 
     else if((double)i/SIGNAL_SIZE>=0.5 && (double)i/SIGNAL_SIZE<1) h_signal[i].x = (double)i/SIGNAL_SIZE-1; 
     h_signal[i].y = 0; 
    } 

// Allocate device memory for signal 
    cuDoubleComplex *d_signal; 

    cudaMalloc((void **) &d_signal, SIGNAL_SIZE*sizeof(cuDoubleComplex)); 
    // Copy host memory to device 
    cudaMemcpy(d_signal, h_signal, SIGNAL_SIZE*sizeof(cuDoubleComplex), cudaMemcpyHostToDevice); 


cudaEventRecord(start, 0); 
    cufftHandle plan; 
    cufftPlan1d(&plan, SIGNAL_SIZE , CUFFT_C2C, 1); 

    // FFT computation 
    cufftExecC2C(plan, (cufftComplex *) d_signal, (cufftComplex *) d_signal, 
     CUFFT_FORWARD); 

    cufftExecC2C(plan, (cufftComplex *) d_signal, (cufftComplex *) d_signal, CUFFT_INVERSE); 

    cuDoubleComplex *h_signal_inv =(cuDoubleComplex *) malloc(sizeof(cuDoubleComplex) * SIGNAL_SIZE); 
    cudaMemcpy(h_signal_inv, d_signal, sizeof(cuDoubleComplex) * SIGNAL_SIZE, cudaMemcpyDeviceToHost); 


    cudaEventRecord(stop, 0); 
    cudaEventSynchronize(stop); 

    float elapsedTime; 
    cudaEventElapsedTime(&elapsedTime, start, stop); 
    printf("Elapsed Time: %3.1f ms\n", elapsedTime); 


    for(int i=0;i<SIGNAL_SIZE;i++) printf("\n%f %f", h_signal[i].x, h_signal_inv[i].x); 

    cufftDestroy(plan); 

    free(h_signal); 
    free(h_signal_inv); 

    cudaFree(d_signal); 

    cudaDeviceReset(); 
    return 0; 
} 

我想變換的信號，然後再回來與逆，但輸出是錯誤的上半年。

你能幫我找到錯誤嗎？

非常感謝！

Answer 1

你正在讓你的數據類型感到困惑。

cufftDoubleComplex與cufftComplex不一樣。當使用cufftDoubleComplex時，your transform type should be Z2Z, not C2C。

另外，爲了看到數據奇偶性做時通過逆一個正向變換，然後使用變換CUFFT，它necessary to divide the result by the signal size：

CUFFT執行未歸一化的FFT;也就是說，對輸入數據組執行正向FFT，然後對所得到的組進行逆FFT，得到等於輸入的數據，按照元素的數量縮放。通過數據集尺寸的倒數來縮放變換，留給用戶以適合的方式執行。

下面的代碼具有上述問題解決，應該給你更好的結果：

#include <stdlib.h> 
#include <stdio.h> 
#include <string.h> 
#include <math.h> 
#include <cufft.h> 
#include <cuComplex.h> 

#define SIGNAL_SIZE  1024 

int main(int argc, char **argv) { 
    cudaEvent_t start, stop; 
    cudaEventCreate(&start); 
    cudaEventCreate(&stop); 

    // Allocate host memory for the signal 
    cuDoubleComplex *h_signal = (cuDoubleComplex *) malloc(sizeof(cuDoubleComplex) * SIGNAL_SIZE); 

    // Initalize the memory for the signal 
    for (unsigned int i = 0; i < SIGNAL_SIZE; ++i) { 
     if((double)i/SIGNAL_SIZE>=0 && (double)i/SIGNAL_SIZE<0.5) h_signal[i].x = (double)i/SIGNAL_SIZE; 
     else if((double)i/SIGNAL_SIZE>=0.5 && (double)i/SIGNAL_SIZE<1) h_signal[i].x = (double)i/SIGNAL_SIZE-1; 
     h_signal[i].y = 0; 
    } 

// Allocate device memory for signal 
    cuDoubleComplex *d_signal; 

    cudaMalloc((void **) &d_signal, SIGNAL_SIZE*sizeof(cuDoubleComplex)); 
    // Copy host memory to device 
    cudaMemcpy(d_signal, h_signal, SIGNAL_SIZE*sizeof(cuDoubleComplex), cudaMemcpyHostToDevice); 


cudaEventRecord(start, 0); 
    cufftHandle plan; 
    cufftPlan1d(&plan, SIGNAL_SIZE , CUFFT_Z2Z, 1); 

    // FFT computation 
    cufftExecZ2Z(plan, d_signal, d_signal, CUFFT_FORWARD); 

    cufftExecZ2Z(plan, d_signal, d_signal, CUFFT_INVERSE); 

    cuDoubleComplex *h_signal_inv =(cuDoubleComplex *) malloc(sizeof(cuDoubleComplex) * SIGNAL_SIZE); 
    cudaMemcpy(h_signal_inv, d_signal, sizeof(cuDoubleComplex) * SIGNAL_SIZE, cudaMemcpyDeviceToHost); 


    cudaEventRecord(stop, 0); 
    cudaEventSynchronize(stop); 

    float elapsedTime; 
    cudaEventElapsedTime(&elapsedTime, start, stop); 
    printf("Elapsed Time: %3.1f ms\n", elapsedTime); 


    for(int i=0;i<SIGNAL_SIZE;i++) printf("\n%f %f", h_signal[i].x, h_signal_inv[i].x/SIGNAL_SIZE); 

    cufftDestroy(plan); 

    free(h_signal); 
    free(h_signal_inv); 

    cudaFree(d_signal); 

    cudaDeviceReset(); 
    return 0; 
}