Pergunta

Estou confuso com os resultados que estou obtendo da FFT e apreciaria qualquer ajuda.

Estou usando o FFTW 3.2.2, mas obtive resultados semelhantes com outras implementações da FFT (em Java). Quando tomo a FFT de uma onda senoidal, a escala do resultado depende da frequência (Hz) da onda-especificamente, seja perto de um número inteiro ou não. Os valores resultantes são muito pequenos quando a frequência está próxima de um número inteiro e são ordens de magnitude maiores quando a frequência está entre os números inteiros. Este gráfico mostra a magnitude do pico no resultado da FFT correspondente à frequência da onda, para diferentes frequências. Isto está certo??

Eu verifiquei se a FFT inversa da FFT é igual à onda seno -seno original vezes o número de amostras, e é. A forma da FFT também parece estar correta.

Não seria tão ruim se eu estivesse analisando ondas senoidais individuais, porque eu poderia apenas procurar o pico na FFT, independentemente de sua altura. O problema é que eu quero analisar quantias de ondas senoidais. Se estou analisando uma soma de ondas senoidais em, digamos, 440 Hz e 523,25 Hz, apenas o pico para o de 523,25 Hz aparece. O pico para o outro é tão pequeno que parece barulho. Deve haver alguma maneira de fazer isso funcionar, porque no Matlab funciona-eu recebo picos de tamanho semelhante nas duas frequências. Como posso alterar o código abaixo para equalizar a escala para diferentes frequências?

#include <cstdlib>
#include <cstring>
#include <cmath> 
#include <fftw3.h>
#include <cstdio>
using namespace std; 

const double PI = 3.141592;

/* Samples from 1-second sine wave with given frequency (Hz) */
void sineWave(double a[], double frequency, int samplesPerSecond, double ampFactor); 

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

 /* Args: frequency (Hz), samplesPerSecond, ampFactor */
 if (argc != 4)  return -1; 
 double frequency  = atof(argv[1]); 
 int samplesPerSecond = atoi(argv[2]); 
 double ampFactor  = atof(argv[3]); 

 /* Init FFT input and output arrays. */
 double * wave = new double[samplesPerSecond]; 
 sineWave(wave, frequency, samplesPerSecond, ampFactor); 
 double * fftHalfComplex = new double[samplesPerSecond]; 
 int fftLen = samplesPerSecond/2 + 1; 
 double * fft = new double[fftLen]; 
 double * ifft = new double[samplesPerSecond]; 

 /* Do the FFT. */
 fftw_plan plan = fftw_plan_r2r_1d(samplesPerSecond, wave, fftHalfComplex, FFTW_R2HC, FFTW_ESTIMATE);
 fftw_execute(plan); 
 memcpy(fft, fftHalfComplex, sizeof(double) * fftLen); 
 fftw_destroy_plan(plan);

 /* Do the IFFT. */
 fftw_plan iplan = fftw_plan_r2r_1d(samplesPerSecond, fftHalfComplex, ifft, FFTW_HC2R, FFTW_ESTIMATE); 
 fftw_execute(iplan); 
 fftw_destroy_plan(iplan);

 printf("%s,%s,%s", argv[1], argv[2], argv[3]); 
 for (int i = 0; i < samplesPerSecond; i++) {
  printf("\t%.6f", wave[i]); 
 }
 printf("\n"); 
 printf("%s,%s,%s", argv[1], argv[2], argv[3]); 
 for (int i = 0; i < fftLen; i++) {
  printf("\t%.9f", fft[i]); 
 }
 printf("\n"); 
 printf("\n"); 
 printf("%s,%s,%s", argv[1], argv[2], argv[3]); 
 for (int i = 0; i < samplesPerSecond; i++) {
  printf("\t%.6f (%.6f)", ifft[i], samplesPerSecond * wave[i]);  // actual and expected result
 }

 delete[] wave; 
 delete[] fftHalfComplex; 
 delete[] fft; 
 delete[] ifft; 
}

void sineWave(double a[], double frequency, int samplesPerSecond, double ampFactor) {
 for (int i = 0; i < samplesPerSecond; i++) {
  double time = i / (double) samplesPerSecond; 
  a[i] = ampFactor * sin(2 * PI * frequency * time); 
 }
}
Foi útil?

Solução

Os valores resultantes são muito pequenos quando a frequência está próxima de um número inteiro e são ordens de magnitude maiores quando a frequência está entre os números inteiros.

Isso ocorre porque uma transformação rápida de Fourier assume que a entrada é periódica e é repetida infinitamente. Se você tem um número não integral de ondas senoidais e repetir essa forma de onda, não é uma onda senoidal perfeita. Isso causa um resultado de FFT que sofre de "vazamento espectral"

Investigar funções de janela. Estes atenuam a entrada no início e no final, para que o vazamento espectral seja diminuído.

PS: Se você deseja obter conteúdo preciso de frequência em torno do fundamental, capture muitos ciclos de ondas e não precisa capturar muitos pontos por ciclo (32 ou 64 pontos por ciclo provavelmente é suficiente). Se você deseja obter conteúdo preciso de frequência em harmônicos mais altos, capture um número menor de ciclos e mais pontos por ciclo.

Outras dicas

Eu só posso recomendar que você olhe para o código de rádio GNU. O arquivo que pode ser de particular interesse para você é usrp_fft.py.

Licenciado em: CC-BY-SA com atribuição
Não afiliado a StackOverflow
scroll top