Fourier-Transformationen
Die Fourier-Transformation ist eine mathematische Formel, die ein zeitlich oder räumlich abgetastetes Signal in dasselbe in zeitlicher oder räumlicher Frequenz abgetastete Signal transformiert. In der Signalverarbeitung kann die Fourier-Transformation wichtige Eigenschaften eines Signals aufzeigen, nämlich seine Frequenzkomponenten.
Die Fourier-Transformation für einen Vektor mit gleichmäßig abgetasteten Punkten ist wie folgt definiert:
ist eine der komplexen Wurzeln der Einheit, wobei die imaginäre Einheit ist. Für und liegen die Indizes und im Bereich von bis .
Die fft
-Funktion in MATLAB® verwendet einen schnellen Fourier-Transformationsalgorithmus, um die Fourier-Transformation von Daten zu berechnen. Betrachten Sie ein sinusförmiges Signal x
, das eine Funktion der Zeit t
ist und Frequenzkomponenten von 15 Hz und 20 Hz besitzt. Verwenden Sie einen Zeitvektor, der in Schritten von 1/50 Sekunden über einen Zeitraum von 10 Sekunden abgetastet wird.
Ts = 1/50; t = 0:Ts:10-Ts; x = sin(2*pi*15*t) + sin(2*pi*20*t); plot(t,x) xlabel('Time (seconds)') ylabel('Amplitude')
Berechnen Sie die Fourier-Transformation des Signals und erstellen Sie den Vektor f
, der der Abtastung des Signals im Frequenzraum entspricht.
y = fft(x); fs = 1/Ts; f = (0:length(y)-1)*fs/length(y);
Stellen Sie den Betrag des transformierten Signals als Funktion der Frequenz dar.
plot(f,abs(y)) xlabel('Frequency (Hz)') ylabel('Magnitude') title('Magnitude')
Die Grafik zeigt vier Frequenzspitzen, obwohl das Signal zwei Frequenzspitzen bei 15 Hz und 20 Hz haben sollte. Hier ist die zweite Hälfte des Diagramms das Spiegelbild der ersten Hälfte. Die diskrete Fourier-Transformation eines Signals im Zeitbereich hat eine periodische Natur, bei der die erste Hälfte des Spektrums in den positiven Frequenzen und die zweite Hälfte in den negativen Frequenzen liegt. Die Frequenzkomponenten 30 Hz und 35 Hz in der Grafik entsprechen den Frequenzkomponenten -20 Hz und -15 Hz. Um diese Periodizität besser zu visualisieren, können Sie die Funktion fftshift
verwenden, die eine kreisförmige Verschiebung der Transformation um den Nullpunkt herum vornimmt.
n = length(x); fshift = (-n/2:n/2-1)*(fs/n); yshift = fftshift(y); plot(fshift,abs(yshift)) xlabel('Frequency (Hz)') ylabel('Magnitude')
Verrauschte Signale
Bei wissenschaftlichen Anwendungen werden Signale oft durch zufälliges Rauschen verfälscht, wodurch ihre Frequenzkomponenten verschleiert werden. Die Fourier-Transformation kann zufälliges Rauschen herausrechnen und die Frequenzen sichtbar machen. Erzeugen Sie zum Beispiel ein neues Signal xnoise
, indem Sie Gaußsches Rauschen in das ursprüngliche Signal x
einspeisen.
rng('default')
xnoise = x + 2.5*randn(size(t));
Die Signalleistung als Funktion der Frequenz ist eine in der Signalverarbeitung häufig verwendete Metrik. Die Leistung ist der quadrierte Betrag der Fourier-Transformation eines Signals, normiert auf die Anzahl der Frequenzabtastungen. Berechnen Sie das Leistungsspektrum des verrauschten Signals, das auf der Nullfrequenz zentriert ist, und plotten Sie es. Trotz des Rauschens können Sie die Frequenzen des Signals aufgrund der Leistungsspitzen noch erkennen.
ynoise = fft(xnoise); ynoiseshift = fftshift(ynoise); power = abs(ynoiseshift).^2/n; plot(fshift,power) title('Power') xlabel('Frequency (Hz)') ylabel('Power')
Rechnerische Effizienz
Wenn Sie die Fourier-Transformationsformel direkt zur Berechnung der einzelnen Elemente von verwenden, erfordert dies eine Größenordnung von Gleitkomma-Operationen. Der Algorithmus der schnellen Fourier-Transformation benötigt nur eine Größenordnung von Operationen zur Berechnung. Diese Berechnungseffizienz ist ein großer Vorteil bei der Verarbeitung von Daten mit Millionen von Datenpunkten. Viele spezialisierte Implementierungen des Algorithmus der schnellen Fourier-Transformation sind sogar noch effizienter, wenn kleine Primfaktoren besitzt, z. B. wenn eine Potenz von 2 ist.
Betrachten Sie beispielsweise Audiodaten, die von Unterwassermikrofonen vor der Küste Kaliforniens gesammelt wurden. Diese Daten finden Sie in einer vom Bioacoustics Research Program der Cornell University verwalteten Bibliothek. Laden und formatieren Sie eine Teilmenge der Daten in bluewhale.au
, die Laute des Pazifischen Blauwals enthält. Da die Gesänge der Blauwale aus niederfrequenten Töne bestehen, sind sie für den Menschen kaum hörbar. Die Zeitskala in den Daten ist um den Faktor 10 komprimiert, um die Tonhöhe anzuheben und die Vokalisationen deutlicher hörbar zu machen. Mit dem Befehl sound(x,fs)
können Sie sich die gesamte Audiodatei anhören.
whaleFile = 'bluewhale.au'; [x,fs] = audioread(whaleFile); whaleMoan = x(2.45e4:3.10e4); t = 10*(0:1/fs:(length(whaleMoan)-1)/fs); plot(t,whaleMoan) xlabel('Time (seconds)') ylabel('Amplitude') xlim([0 t(end)])
Geben Sie eine neue Signallänge an, die um die nächste 2er-Potenz größer ist als die ursprüngliche Länge. Berechnen Sie dann mit fft
die Fourier-Transformation unter Verwendung der neuen Signallänge. fft
füllt die Daten automatisch mit Nullen auf, um die Stichprobengröße zu erhöhen. Diese Auffüllung kann die Berechnung der Transformation erheblich beschleunigen, insbesondere bei Stichprobengrößen mit großen Primfaktoren.
m = length(whaleMoan); n = pow2(nextpow2(m)); y = fft(whaleMoan,n);
Plotten Sie das Leistungsspektrum des Signals. Die Darstellung zeigt, dass der Ton aus einer Grundfrequenz um 17 Hz und einer Folge von Obertönen besteht, wobei die Betonung auf dem zweiten Oberton liegt.
f = (0:n-1)*(fs/n)/10; % frequency vector power = abs(y).^2/n; % power spectrum plot(f(1:floor(n/2)),power(1:floor(n/2))) xlabel('Frequency (Hz)') ylabel('Power')
Phase der Sinusoiden
Mit der Fourier-Transformation können Sie auch das Phasenspektrum des ursprünglichen Signals extrahieren. Erstellen Sie zum Beispiel ein Signal, das aus zwei Sinuskurven mit den Frequenzen 15 Hz und 40 Hz besteht. Die erste Sinuskurve ist eine Kosinuswelle mit der Phase und die zweite ist eine Kosinuswelle mit der Phase . Tasten Sie das Signal 1 Sekunde lang bei 100 Hz ab.
fs = 100; t = 0:1/fs:1-1/fs; x = cos(2*pi*15*t - pi/4) - sin(2*pi*40*t);
Berechnen Sie die Fourier-Transformation des Signals. Plotten Sie den Betrag der Transformation als Funktion der Frequenz.
y = fft(x); z = fftshift(y); ly = length(y); f = (-ly/2:ly/2-1)/ly*fs; stem(f,abs(z)) xlabel("Frequency (Hz)") ylabel("|y|") grid
Berechnen Sie die Phase der Transformation und entfernen Sie dabei Transformationswerte mit kleinem Betrag. Plotten Sie die Phase als Funktion der Frequenz.
tol = 1e-6; z(abs(z) < tol) = 0; theta = angle(z); stem(f,theta/pi) xlabel("Frequency (Hz)") ylabel("Phase / \pi") grid
Siehe auch
fft
| fftshift
| nextpow2
| ifft
| fft2
| fftn
| fftw