Каффт установить частоту?

Я использую CUDA Cufft для обработки данных, которые я получаю от гидрофона (500 000 целых чисел в секунду при 250 Гц, верхний и нижний каналы). Теперь основной пример того, как работает Cufft, здесь …

 void runTest(int argc, char** argv)

{

printf("[1DCUFFT] is starting...\n");cufftComplex* h_signal = (cufftComplex*)malloc(sizeof(cufftComplex)* SIGNAL_SIZE);
// Allocate host memory for the signal
//Complex* h_signal = (Complex*)malloc(sizeof(Complex) * SIGNAL_SIZE);
// Initalize the memory for the signal
for (unsigned int i = 0; i < SIGNAL_SIZE; ++i) {
h_signal[i].x = rand() / (float)RAND_MAX;
h_signal[i].y = 0;
}int mem_size = sizeof(cufftComplex)* SIGNAL_SIZE;

// Allocate device memory for signal
cufftComplex* d_signal;
cudaMalloc((void**)&d_signal, mem_size);

// Copy host memory to device
cudaMemcpy(d_signal, h_signal, mem_size,
cudaMemcpyHostToDevice);// CUFFT plan
cufftHandle plan;
cufftPlan1d(&plan, mem_size, CUFFT_C2C, 1);

// Transform signal
printf("Transforming signal cufftExecC2C\n");
cufftExecC2C(plan, (cufftComplex *)d_signal, (cufftComplex *)d_signal, CUFFT_FORWARD);// Transform signal back
printf("Transforming signal back cufftExecC2C\n");
cufftExecC2C(plan, (cufftComplex *)d_signal, (cufftComplex *)d_signal, CUFFT_INVERSE);

// Copy device memory to host
cufftComplex* h_inverse_signal = (cufftComplex*)malloc(sizeof(cufftComplex)* SIGNAL_SIZE);;
cudaMemcpy(h_inverse_signal, d_signal, mem_size,
cudaMemcpyDeviceToHost);

for (int i = 0; i < SIGNAL_SIZE; i++){
h_inverse_signal[i].x = h_inverse_signal[i].x / (float)SIGNAL_SIZE;
h_inverse_signal[i].y = h_inverse_signal[i].y / (float)SIGNAL_SIZE;

printf("first : %f %f  after %f %f \n", h_signal[i].x, h_signal[i].y, h_inverse_signal[i].x, h_inverse_signal[i].y);
}
//Destroy CUFFT context
cufftDestroy(plan);

// cleanup memory
free(h_signal);

free(h_inverse_signal);
cudaFree(d_signal);
cudaDeviceReset();
}

Теперь все, что я хочу знать, это как установить частоту БПФ (cufft) на 250 Гц?

Спасибо

Джеймс

0

Решение

Вы не БПФ из N точек одинаково, независимо от частоты, с которой эти N точек были отобраны.

Кроме того, 500.000 целых чисел в секунду — это частота дискретизации 500.000 Гц, то есть 500 кГц. Это дает вам предел Найквиста в 250 кГц.

2

Другие решения

Если я вас правильно понял, вам просто нужно знать, какой элемент в выходном векторе равен 250 Гц.

БПФ дает вам все частоты, которые могут быть рассчитаны на основе длины и временного разрешения вашего вектора времени.
Простое правило для расчета:
— диапазон частот = 1 / разрешение по времени.
— разрешение по частоте = 1 / длительность.

Кроме того, нужно знать, что БПФ реальной функции (без мнимой части данных вектора времени) дает симметричный спектр с избыточностью. Спектр достигает от (- 1/2 частотного диапазона до +1/2 частотного диапазона). Данные отрицательной частоты могут быть отброшены в случае вектора реального времени. Хотя это немного сложнее. Стандартная реализация FFT (которая является операцией на месте) сначала дает вам положительные частоты, а затем отрицательные частоты. Поскольку вас интересуют только положительные частоты, вторая половина вектора БПФ может быть отброшена. В вашем случае просто игнорируйте данные выше индекса 250к.

В вашем случае частоты колеблются от -250 кГц до 250 кГц с разрешением 1 Гц, но из-за вышеизложенного первые 250 тыс. Точек на самом деле являются положительными частотами на расстоянии 1 Гц.

Поэтому возьмите 250-ю точку в БПФ (без сдвига, т.е. в сыром виде), и вы получите сигнал с частотой 250 Гц. Я бы нанес данные на график от 0 до 500, чтобы увидеть, насколько широк этот пик составляет около 250 Гц. Сила сигнала является интегралом этих ненулевых частот (ненулевое значение применяется здесь для обозначения всего, что выше шума). Ширина сигнала указывает модуляцию, которая применяется к сигналу (который может включать в себя другие артефакты измерения). Если сигнал сдвинут от 250 Гц, у вас может быть доплеровский сдвиг (либо ваш источник, либо вы движетесь).

Если вас интересует только конечный частотный диапазон, может быть быстрее вычислить интеграл Фурье (O (n ^ 2)) только для этих нескольких частотных точек. Обычно люди используют FFT, потому что это O (n * log (n)), но если вам нужно только сказать 10 частотных точек, то O (10 * n) не сильно отличается.

2