#c #cuda #fft #cufft
#c #cuda #БПФ #cufft
Вопрос:
Я использую cuda версии 7.5 cufft для выполнения некоторого БПФ и обратного БПФ. У меня возникла проблема при выполнении обратного БПФ с использованием cufftExecC2R(.,.) функции.
На самом деле, когда я использую a batch_size = 1 в cufftPlan1d(,) , я получаю правильный результат. Однако, когда я увеличиваю размер пакета, результаты неверны.
Я вставляю пример минимального кода, чтобы проиллюстрировать это. Пожалуйста, не обращайте внимания на грязность кода, поскольку я только что быстро создал это.
#include <cufft.h>
#include <stdlib.h>
#include <stdio.h>
#include <string.h>
#include <math.h>
#include <ctime>
#include <iostream>
typedef float2 Complex;
void iTest(int argc, char** argv);
#define SIGNAL_SIZE 9
#define BATCH_SIZE 2
int main(int argc, char** argv) {
iTest(argc, argv);
return 0;
}
void iProcess(Complex *x, double *y, size_t n) {
cufftComplex *deviceData;
cudaMalloc(reinterpret_cast<void**>(amp;deviceData),
SIGNAL_SIZE * BATCH_SIZE * sizeof(cufftComplex));
cudaMemcpy(deviceData, x, SIGNAL_SIZE * sizeof(cufftComplex) * BATCH_SIZE,
cudaMemcpyHostToDevice);
cufftResult cufftStatus;
cufftHandle handle;
cufftStatus = cufftPlan1d(amp;handle, SIGNAL_SIZE, CUFFT_C2C, BATCH_SIZE);
if (cufftStatus != cudaSuccess) {
printf("cufftPlan1d failed!");
}
cufftComplex *d_complex;
cudaMalloc(reinterpret_cast<void**>(amp;d_complex),
sizeof(cufftComplex) * SIGNAL_SIZE * BATCH_SIZE);
cufftStatus = cufftExecC2C(handle, deviceData, d_complex, CUFFT_FORWARD);
if (cufftStatus != cudaSuccess) {
printf("cufftExecR2C failed!");
}
cufftComplex *hostOutputData = (cufftComplex*)malloc(
(SIGNAL_SIZE) * BATCH_SIZE * sizeof(cufftComplex));
cudaMemcpy(hostOutputData, d_complex,
SIGNAL_SIZE * sizeof(cufftComplex) * BATCH_SIZE,
cudaMemcpyDeviceToHost);
std::cout << "nPrinting COMPLEX" << "n";
for (int j = 0; j < (SIGNAL_SIZE) * BATCH_SIZE; j )
printf("%i t %f t %fn", j, hostOutputData[j].x, hostOutputData[j].y);
//! convert complex to real
cufftHandle c2r_handle;
cufftStatus = cufftPlan1d(amp;c2r_handle, SIGNAL_SIZE, CUFFT_C2R, BATCH_SIZE);
if (cufftStatus != cudaSuccess) {
printf("cufftPlan1d failed!");
}
cufftReal *d_odata;
cudaMalloc(reinterpret_cast<void**>(amp;d_odata),
sizeof(cufftReal) * SIGNAL_SIZE * BATCH_SIZE);
cufftStatus = cufftExecC2R(c2r_handle, d_complex, d_odata);
cufftReal odata[SIGNAL_SIZE * BATCH_SIZE];
cudaMemcpy(odata, d_odata, sizeof(cufftReal) * SIGNAL_SIZE * BATCH_SIZE,
cudaMemcpyDeviceToHost);
std::cout << "nPrinting REAL" << "n";
for (int i = 0; i < SIGNAL_SIZE * BATCH_SIZE; i ) {
std::cout << i << " t" << odata[i]/(SIGNAL_SIZE) << "n";
}
cufftDestroy(handle);
cudaFree(deviceData);
}
void iTest(int argc, char** argv) {
Complex* h_signal = reinterpret_cast<Complex*>(
malloc(sizeof(Complex) * SIGNAL_SIZE * BATCH_SIZE));
std::cout << "nPrinting INPUT" << "n";
for (unsigned int i = 0; i < SIGNAL_SIZE * BATCH_SIZE; i) {
h_signal[i].x = rand() / static_cast<float>(RAND_MAX);
h_signal[i].y = 0;
std::cout << i << "t" << h_signal[i].x << "n";
}
std::cout << "n";
double y[SIGNAL_SIZE * BATCH_SIZE];
iProcess(h_signal, y, 1);
}
Я не могу выяснить, где ошибка в моем коде и какую информацию я упускаю.
Пример вывода при использовании BATCH_SIZE = 1
Пример вывода при использовании BATCH_SIZE = 2
Ответ №1:
Информация, которую вам не хватает, заключается в том, что вы не понимаете, что существуют различия в формате входных данных, ожидаемых для преобразования C2C, по сравнению с C2R (или R2C).
Вам следует начать с чтения этого раздела и этого раздела документации CUFFT.
Обратите внимание, что в нем говорится:
Каждая из этих функций требует другого расположения входных данных
Но вы передаете входные данные, которые были правильными для преобразования C2C, непосредственно в преобразование C2R. Это не сработает.
Самое прямое решение IMO — преобразовать всю вашу работу в типы преобразования C2C. Преобразование C2C может поддерживать как прямое (например, «реальное в комплексное»), так и обратное (например, «комплексное в вещественное»). Используемый вами тип преобразования C2R также может поддерживать преобразование «от сложного к реальному», но расположение данных, которое вы использовали бы для C2R, отличается от расположения данных, которое вы использовали бы для C2C с указанным обратным путем, для того, что в остальном является тем же самым преобразованием. Вы этого не учли.
Вот отработанный пример, показывающий модифицированную версию вашего кода, которая использует C2C как для прямого, так и для обратного путей и корректно воспроизводит входные данные для пакета размером 2:
$ cat t19.cu
#include <cufft.h>
#include <stdlib.h>
#include <stdio.h>
#include <string.h>
#include <math.h>
#include <ctime>
#include <iostream>
typedef float2 Complex;
void iTest(int argc, char** argv);
#define SIGNAL_SIZE 9
#define BATCH_SIZE 2
int main(int argc, char** argv) {
iTest(argc, argv);
return 0;
}
void iProcess(Complex *x, double *y, size_t n) {
cufftComplex *deviceData;
cudaMalloc(reinterpret_cast<void**>(amp;deviceData),
SIGNAL_SIZE * BATCH_SIZE * sizeof(cufftComplex));
cudaMemcpy(deviceData, x, SIGNAL_SIZE * sizeof(cufftComplex) * BATCH_SIZE,
cudaMemcpyHostToDevice);
cufftResult cufftStatus;
cufftHandle handle;
cufftStatus = cufftPlan1d(amp;handle, SIGNAL_SIZE, CUFFT_C2C, BATCH_SIZE);
if (cufftStatus != cudaSuccess) {
printf("cufftPlan1d failed!");
}
cufftComplex *d_complex;
cudaMalloc(reinterpret_cast<void**>(amp;d_complex),
sizeof(cufftComplex) * SIGNAL_SIZE * BATCH_SIZE);
cufftStatus = cufftExecC2C(handle, deviceData, d_complex, CUFFT_FORWARD);
if (cufftStatus != cudaSuccess) {
printf("cufftExecR2C failed!");
}
cufftComplex *hostOutputData = (cufftComplex*)malloc(
(SIGNAL_SIZE) * BATCH_SIZE * sizeof(cufftComplex));
cudaMemcpy(hostOutputData, d_complex,
SIGNAL_SIZE * sizeof(cufftComplex) * BATCH_SIZE,
cudaMemcpyDeviceToHost);
std::cout << "nPrinting COMPLEX" << "n";
for (int j = 0; j < (SIGNAL_SIZE) * BATCH_SIZE; j )
printf("%i t %f t %fn", j, hostOutputData[j].x, hostOutputData[j].y);
//! convert complex to real
/* cufftHandle c2r_handle;
cufftStatus = cufftPlan1d(amp;c2r_handle, SIGNAL_SIZE, CUFFT_C2R, BATCH_SIZE);
if (cufftStatus != cudaSuccess) {
printf("cufftPlan1d failed!");
}
*/
cufftComplex *d_odata;
cudaMalloc(reinterpret_cast<void**>(amp;d_odata),
sizeof(cufftComplex) * SIGNAL_SIZE * BATCH_SIZE);
cufftStatus = cufftExecC2C(handle, d_complex, d_odata, CUFFT_INVERSE);
cufftComplex odata[SIGNAL_SIZE * BATCH_SIZE];
cudaMemcpy(odata, d_odata, sizeof(cufftComplex) * SIGNAL_SIZE * BATCH_SIZE,
cudaMemcpyDeviceToHost);
std::cout << "nPrinting REAL" << "n";
for (int i = 0; i < SIGNAL_SIZE * BATCH_SIZE; i ) {
std::cout << i << " t" << odata[i].x/(SIGNAL_SIZE) << "n";
}
cufftDestroy(handle);
cudaFree(deviceData);
}
void iTest(int argc, char** argv) {
Complex* h_signal = reinterpret_cast<Complex*>(
malloc(sizeof(Complex) * SIGNAL_SIZE * BATCH_SIZE));
std::cout << "nPrinting INPUT" << "n";
for (unsigned int i = 0; i < SIGNAL_SIZE * BATCH_SIZE; i) {
h_signal[i].x = rand() / static_cast<float>(RAND_MAX);
h_signal[i].y = 0;
std::cout << i << "t" << h_signal[i].x << "n";
}
std::cout << "n";
double y[SIGNAL_SIZE * BATCH_SIZE];
iProcess(h_signal, y, 1);
}
$ nvcc -arch=sm_61 -o t19 t19.cu -lcufft
t19.cu: In function ‘void iProcess(Complex*, double*, size_t)’:
t19.cu:34:32: warning: comparison between ‘cufftResult {aka enum cufftResult_t}’ and ‘enum cudaError’ [-Wenum-compare]
if (cufftStatus != cudaSuccess) {
^
t19.cu:43:32: warning: comparison between ‘cufftResult {aka enum cufftResult_t}’ and ‘enum cudaError’ [-Wenum-compare]
if (cufftStatus != cudaSuccess) {
^
$ cuda-memcheck ./t19
========= CUDA-MEMCHECK
Printing INPUT
0 0.840188
1 0.394383
2 0.783099
3 0.79844
4 0.911647
5 0.197551
6 0.335223
7 0.76823
8 0.277775
9 0.55397
10 0.477397
11 0.628871
12 0.364784
13 0.513401
14 0.95223
15 0.916195
16 0.635712
17 0.717297
Printing COMPLEX
0 5.306536 0.000000
1 0.015338 -0.734991
2 -0.218001 0.740248
3 0.307508 -0.706533
4 1.022732 0.271765
5 1.022732 -0.271765
6 0.307508 0.706533
7 -0.218001 -0.740248
8 0.015338 0.734991
9 5.759857 0.000000
10 -0.328981 0.788566
11 0.055356 -0.521014
12 -0.127504 0.581872
13 0.014066 0.123027
14 0.014066 -0.123027
15 -0.127504 -0.581872
16 0.055356 0.521014
17 -0.328981 -0.788566
Printing REAL
0 0.840188
1 0.394383
2 0.783099
3 0.79844
4 0.911647
5 0.197551
6 0.335223
7 0.76823
8 0.277775
9 0.55397
10 0.477397
11 0.628871
12 0.364784
13 0.513401
14 0.95223
15 0.916195
16 0.635712
17 0.717297
========= ERROR SUMMARY: 0 errors
$
Комментарии:
1. Я протестировал, и это работает отлично. Спасибо, и я думаю, что тщательное прочтение предоставленных вами ссылок на документацию очень полезно. Я рекомендую всем, у кого возникли проблемы, сначала прочитать их.