Dlaczego transformacja Fouriera pojedynczego cyklu fali sinusoidalnej nie jest pojedynczym słupkiem?

12

Próbowałem różnych kodów transformacji Fouriera na pojedynczych falach sinusoidalnych i wszystkie wytwarzają rozproszone widmo z rezonansem na częstotliwości sygnału, kiedy teoretycznie powinny wyświetlać pojedynczy słupek.

Częstotliwość próbkowania ma niewielki wpływ (tutaj 10 kHz), jednak liczba cykli:

Jeden cykl:

wprowadź opis zdjęcia tutaj

100 cykli:

wprowadź opis zdjęcia tutaj

100000 cykli:

wprowadź opis zdjęcia tutaj

Wygląda na to, że transformata Fouriera zbiega się tylko dla nieskończonej liczby cykli, dlaczego tak jest? Czy okno czasowe dokładnie jednego cyklu nie powinno dawać takich samych wyników jak w przypadku N cykli?

Zastosowanie: Wynika to zarówno z ciekawości, jak i dlatego, że chcę dowiedzieć się, jak bardzo reakcja krokowa systemu pierwszego rzędu będzie ekscytująca rezonans zespołu mechanicznego. Dlatego potrzebuję dokładnej transformaty Fouriera odpowiedzi ... której już nie ufam. Co mogę zrobić, aby poprawić dokładność w oparciu o przypadek „fali sinusoidalnej”?

wprowadź opis zdjęcia tutaj

PS: Te konkretne zrzuty ekranu są oparte na kodzie tutaj .

Pan Mystère
źródło
6
Oprócz zaakceptowanej odpowiedzi, zauważ, że nie ma powodu, aby sądzić, że dyskretna transformata Fouriera (którą obliczasz za pomocą DFT) byłaby impulsem dla sygnału wejściowego, który jest jednym okresem sinusoidy. Ciągła transformata Fouriera sinusoidy jest impulsem, tak, ale sinusoida jest nieskończona. Ograniczenie sygnału w czasie jest równoznaczne z pomnożeniem przez prostokątną funkcję okna. Rezultatem w dziedzinie częstotliwości jest splot impulsu i transformata Fouriera okna, które zasadniczo obserwujesz.
Jason R
Dziękuję za uwagę. Jak więc wytłumaczyć, że jeśli zmienię liczbę NFFT na długość wektora, wynikiem będzie pojedynczy słupek?
Mister Mystère,
1
Dobre pytanie. Dzieje się tak z powodu nieodłącznego założenia DFT. Przyjmuje się, że sygnał o skończonej długości podawany na wejściu DFT okresowo rozciąga się w obu kierunkach z nieskończonym czasem trwania. Dlatego, gdy masz całkowitą liczbę cykli w „otworze” DFT, kończysz się transformacją sinusoidy o nieskończonym czasie trwania: pojedynczego impulsu. Odpowiada to przypadkowi dokładnie zerowego wycieku widmowego i rzadko występuje w praktyce.
Jason R

Odpowiedzi:

30

To artefakt okienkowania.

Połączony kod uzupełnia sygnał próbki o wartości 10.000 zerami, dzięki czemu długość wynosi potęgę dwóch.

%% Author :- Embedded Laboratory

%%This Project shows how to apply FFT on a signal and its physical 
% significance.

fSampling = 10000;          %Sampling Frequency
tSampling = 1/fSampling;    %Sampling Time
L = 10000;                  %Length of Signal
t = (0:L-1)*tSampling;      %Time Vector
F = 100;                    %Frequency of Signal

%% Signal Without Noise
xsig = sin(2*pi*F*t);
...

%%Frequency Transform of above Signal
subplot(2,1,2)
NFFT = 2^nextpow2(L);
Xsig = fft(xsig,NFFT)/L;
...

Zauważ, że w powyższym kodzie FFT jest pobierany z rozmiarem FFT, NFFTktóry jest następną potęgą 2 większą niż długość sygnału (w tym przypadku 16 384). Z dokumentacji Mathworksfft() :

Y = fft(X,n)zwraca n-punkt DFT. fft(X)jest równoważne fft(X, n)gdzie njest rozmiar Xpierwszego wymiaru noningingleton. Jeśli długość Xjest mniejsza niż n, Xjest dopełniana zerami na końcu do długości n. Jeśli długość Xjest większa niż n, sekwencja Xjest obcinana. Gdy Xjest macierzą, długość kolumn jest dostosowywana w ten sam sposób.

Oznacza to, że tak naprawdę nie bierzesz FFT „czystej fali sinusoidalnej” - bierzesz FFT fali sinusoidalnej z płaskim sygnałem po niej.

Jest to równoważne z wzięciem FFT fali sinusoidalnej pomnożonej przez funkcję kwadratowego okna. Widmo FFT jest wówczas splotem widma częstotliwości fali sinusoidalnej (funkcja impulsu) z widmem częstotliwości fali prostokątnej (sinc (f).)

Jeśli zmienisz się L = 16,384tak, że nie będzie padania zera sygnału, zaobserwujesz perfectFFT.

Dalsze słowa kluczowe wyszukiwania: „Wyciek spektralny”, „Funkcja okna”, „Okno Hamminga”.


Edycja: Oczyściłem trochę materiału, który napisałem na ten temat na uniwersytecie, który jest znacznie bardziej szczegółowy. Opublikowałem to na swoim blogu .

Li-aung Yip
źródło
Przez cały czas znajdowałem się tuż przed moją twarzą. Dobra robota, właśnie zmieniłem numer NFFT na długość wektora i to zrobiło.
Mister Mystère,
1
@ MisterMystère: Zobacz link do odpowiedniego materiału, który napisałem na uniwersytecie. Znacznie dokładniejsze wyjaśnienie, w tym zdjęcia.
Li-aung Yip,
(Chociaż zapomniałem ponownie wpisać formuły matematyczne - już naprawione.)
Li-aung Yip
Przypomnienie, że nie ma żadnej przewagi, by uzupełnić do nextpow2 przy użyciu algorytmów Matlab FFT, które moim zdaniem to fftw (najszybsza transformacja Fouriera na zachodzie)
Scott Seidman