Jak sprawdzić wyniki FFT fali sinusoidalnej?

Dostałem plik audio (fala sinusoidalna) 1000 Hz jako dane wejściowe do mojego algorytmu FFT. Mam 8192 próbek widma mocy w jednym zestawie.

Jaki jest najlepszy i najprostszy sposób sprawdzenia, czy moje wyniki są prawidłowe, czy złe?

Jeśli podam cichy plik audio, wówczas wynik wyjściowy wynosi zero dla wszystkich próbek. W fali sinusoidalnej o / p wzrasta od 20 (0-ta próbka) do 26059811 (743-ta próbka) i stopniowo maleje do 40.

Jeśli mam pojęcie o zakresie wyjściowym, mogę technicznie udowodnić, czy FFT działa.

Wszelkie pomysły byłyby pomocne.

Odnoszą się do tego linku do jakichkolwiek wątpliwości technicznych.

Wygląda na to, że obliczasz widmo przez uśrednienie 10 okien (nie nakładających się?), Aby uzyskać kwadrat podniesiony do kwadratu przy częstotliwościach 8192 lub 8193 (od 0 do Nyquist, ale niektóre algorytmy mogą obniżyć częstotliwość Nyquist do bin 8192).

Pierwszą rzeczą do sprawdzenia jest to, że szczyt znajduje się w odpowiednim koszu. Nie powiedziałeś, jaka jest częstotliwość próbkowania, ale bin 743 byłby 743/16384 razy szybszy. Jeśli sygnał rzeczywiście ma częstotliwość 800 Hz, oznacza to, że Fs wynosi około 17640 próbek na sekundę. To wydaje się złe. Twój sygnał testowy prawdopodobnie miałby standardową szybkość, taką jak 8000, 16000, 22050, 32000, 44100 lub 48000. Dla Fs = 22050 szczyt byłby ostro w bin 800/22050 * 16384 = 594.

Kolejnym kryterium do sprawdzenia jest to, że całkowita energia w sygnale jest w przybliżeniu taka sama zarówno w dziedzinie czasu, jak i częstotliwości. Oto przykład w Pythonie:

In [1]: NFFT = 2048; N = 10*NFFT; n = arange(N); Fs = 22050
In [2]: x = 0.4*cos(2*pi*400/Fs*n) + 0.6*cos(2*pi*800/Fs*n)

In [3]: y,freqs = psd(x, NFFT=NFFT, Fs=Fs, pad_to=16384)  # PSD by Welch's Method

In [4]: sum(x**2)/Fs           # time-domain energy
Out[4]: 0.24149869319296949
In [5]: sum(y) * N/16384       # frequency-domain energy
Out[5]: 0.24148752834391252

Sygnał wejściowy x, który składa się z dwóch sinusoid próbkowanych przy Fs = 22050 próbek na sekundę, jest podzielony na 10 nie nakładających się okien o wielkości NFFT = 2048 próbek. Wezwanie do psd (gęstość widmowa mocy) oblicza widmo y jako średnią wielkości do kwadratu z dziesięciu 16384-punktowych DFT (w rzeczywistości jest to 8193 punktów, ponieważ x ma wartość rzeczywistą).

Obliczona energia w dziedzinie częstotliwości ma współczynnik skalowania N / 16384, ponieważ funkcja psd skalowała y do rozmiaru DFT zamiast do całkowitej długości sygnału. To, czy jest to problem, zależy od tego, jak system obsługuje normalizację PSD. Inną opcjonalną normalizacją jest skalowanie o 1 / Fs. To dopasowuje energię do oryginalnego sygnału analogowego. Domyślne normalizacje powinny być dobrze udokumentowane w bibliotece.

Musisz wykreślić wielkość wyniku FFT. Nie znam twojego języka programowania, ale w Pythonie użyłbyś czegoś takiego plot(abs(fft(a))). W przypadku cichego wprowadzania dane wyjściowe powinny składać się z samych zer. Dla wejścia fali sinusoidalnej powinieneś zobaczyć dwa szczyty:

Dla prawdziwego sygnału, kolce będą symetryczne od lewej do prawej. Jeśli robisz prawdziwy FFT (który jest bardziej wydajny obliczeniowo), otrzymasz tylko lewą połowę wykresu jako dane wyjściowe, ponieważ ignoruje nadmiarowe odbicie lustrzane.

Jeśli częstotliwość jest wyższa, kolce będą bliżej środka. Jeśli częstotliwość jest idealnie zsynchronizowana z rozmiarem fragmentu , skok będzie miał tylko jeden punkt szerokości, a wszystko inne będzie dokładnie równe 0. W przeciwnym razie będzie miał zwężającą się „spódnicę” jak powyżej.

Użyłem narzędzia analizy Fouriera w pakiecie Excel Analysis Toolpak, aby szybko sprawdzić dane i wyniki.