MATLAB Language
Fourier-Transformationen und Inverse Fourier-Transformationen

Syntax

1. Y = fft (X)% berechnet die FFT des Vektors oder der Matrix X mit einer Standardtransformationslänge von 256 (für Version zu bestätigen)

2. Y = fft (X, n)% berechnet die FFT von X unter Verwendung von n als Transformationslänge, n muss eine Zahl mit 2 Potenz sein. Ist die Länge von X kleiner als n, füllt Matlab X automatisch mit Nullen auf, so dass Länge (X) = n ist

3. Y = fft (X, n, dim)% berechnet die FFT von X unter Verwendung von n als Transformationslänge entlang der Dimensionsdimension (kann 1 oder 2 für horizontal bzw. vertikal sein)

4. Y = fft2 (X)% Berechnet die 2D-FFT von X

5. Y = fftn (X, dim)% Berechnet die dim-dimensionale FFT von X in Bezug auf den Dimensionsvektor dim.

6. y = ifft (X)% berechnet die Inverse der FFT von X (die eine Matrix / einen Vektor von Zahlen ist) unter Verwendung der Standard-256-Transformationslänge

7. y = ifft (X, n)% berechnet die IFFT von X unter Verwendung von n als Transformationslänge

8. y = ifft (X, n, dim)% berechnet die IFFT von X unter Verwendung von n als Transformationslänge über Maß dim (kann 1 oder 2 für horizontal bzw. vertikal sein)

9. y = ifft (X, n, dim, 'symmetrisch')% Die Option Symmetric bewirkt, dass ifft X als konjugiert symmetrisch entlang der aktiven Dimension behandelt. Diese Option ist nützlich, wenn X nur aufgrund von Rundungsfehlern nicht exakt konjugiert ist.

10. y = ifft2 (X)% Berechnet die inverse 2D-ft von X

11. y = ifftn (X, dim)% Berechnet die inverse dim-dimensionale fft von X.

Parameter

Bemerkungen

ifft(fft(X)) = X

Implementieren Sie eine einfache Fourier-Transformation in Matlab

Die Fourier-Transformation ist wahrscheinlich die erste Lektion in der digitalen Signalverarbeitung. Sie ist überall einsetzbar und ist ein leistungsfähiges Werkzeug für die Analyse von Daten (in allen Bereichen) oder Signalen. Matlab verfügt über eine Reihe leistungsstarker Toolboxen für die Fourier-Transformation. In diesem Beispiel verwenden wir die Fourier-Transformation, um ein grundlegendes Sinuswellensignal zu analysieren und mithilfe von FFT ein so genanntes Periodogramm zu erzeugen:

%Signal Generation
A1=10;                % Amplitude 1
A2=10;                % Amplitude 2
w1=2*pi*0.2;          % Angular frequency 1
w2=2*pi*0.225;        % Angular frequency 2
Ts=1;                 % Sampling time
N=64;                 % Number of process samples to be generated
K=5;                  % Number of independent process realizations
sgm=1;                % Standard deviation of the noise
n=repmat([0:N-1].',1,K);             % Generate resolution
phi1=repmat(rand(1,K)*2*pi,N,1);     % Random phase matrix 1
phi2=repmat(rand(1,K)*2*pi,N,1);     % Random phase matrix 2
x=A1*sin(w1*n*Ts+phi1)+A2*sin(w2*n*Ts+phi2)+sgm*randn(N,K);   % Resulting Signal

NFFT=256;            % FFT length
F=fft(x,NFFT);       % Fast Fourier Transform Result
Z=1/N*abs(F).^2;     % Convert FFT result into a Periodogram

Beachten Sie, dass die diskrete Fourier-Transformation in Matlab durch die schnelle Fourier-Transformation (FFT) implementiert wird. Beide Ergebnisse liefern das gleiche Ergebnis, aber FFT ist eine schnelle Implementierung von DFT.

figure
w=linspace(0,2,NFFT);
plot(w,10*log10(Z)),grid;
xlabel('w [\pi rad/s]')
ylabel('Z(f) [dB]')
title('Frequency Range: [ 0 ,  \omega_s ]')

Inverse Fourier-Transformationen

Einer der größten Vorteile der Fourier-Transformation besteht darin, dass sie ohne Zeitverlust in die Zeitdomäne zurückkehren können. Betrachten wir das gleiche Signal, das wir im vorherigen Beispiel verwendet haben:

A1=10;                % Amplitude 1
A2=10;                % Amplitude 2
w1=2*pi*0.2;          % Angular frequency 1
w2=2*pi*0.225;        % Angular frequency 2
Ts=1;                 % Sampling time
N=64;                 % Number of process samples to be generated
K=1;                  % Number of independent process realizations
sgm=1;                % Standard deviation of the noise
n=repmat([0:N-1].',1,K);             % Generate resolution
phi1=repmat(rand(1,K)*2*pi,N,1);     % Random phase matrix 1
phi2=repmat(rand(1,K)*2*pi,N,1);     % Random phase matrix 2
x=A1*sin(w1*n*Ts+phi1)+A2*sin(w2*n*Ts+phi2)+sgm*randn(N,K);   % Resulting Signal

NFFT=256;            % FFT length
F=fft(x,NFFT);       % FFT result of time domain signal

Wenn wir F in Matlab öffnen, werden wir feststellen, dass es sich um eine Matrix aus komplexen Zahlen handelt, einen Realteil und einen Imaginärteil. Um das ursprüngliche Zeitbereichssignal wiederzugewinnen, benötigen wir definitionsgemäß sowohl den reellen Wert (der die Größenänderung darstellt) als auch den Imaginärbereich (der die Phasenänderung darstellt). Um zum Zeitbereich zurückzukehren, möchte man einfach:

TD = ifft(F,NFFT);   %Returns the Inverse of F in Time Domain

Beachten Sie hier, dass das zurückgegebene TD die Länge 256 haben würde, da wir NFFT auf 256 setzen. Die Länge von x beträgt jedoch nur 64, sodass Matlab Nullen bis zum Ende der TD-Transformation füllt. Wenn beispielsweise NFFT 1024 ist und die Länge 64 beträgt, dann ergibt das zurückgegebene TD 64 + 960 Nullen. Beachten Sie auch, dass Sie aufgrund von Rundungen mit Fließkommazahlen möglicherweise etwas wie 3,1 * 10e-20 erhalten. Dies ist jedoch allgemein gedacht: Bei jedem X entspricht ifft (fft (X)) dem X innerhalb des Rundungsfehlers.

Lassen Sie uns für einen Moment sagen, dass wir nach der Transformation etwas getan haben und nur noch den REALEN Teil der FFT haben:

R = real(F);         %Give the Real Part of the FFT
TDR = ifft(R,NFFT);  %Give the Time Domain of the Real Part of the FFT

Dies bedeutet, dass wir den Imaginärteil unserer FFT verlieren und daher Informationen in diesem umgekehrten Prozess verlieren. Um das Original zu erhalten, ohne Informationen zu verlieren, sollten Sie immer den Imaginärteil der FFT mit imag und Ihre Funktionen entweder auf beide oder auf den Realteil anwenden.

figure
subplot(3,1,1)
plot(x);xlabel('time samples');ylabel('magnitude');title('Original Time Domain Signal')
subplot(3,1,2)
plot(TD(1:64));xlabel('time samples');ylabel('magnitude');title('Inverse Fourier Transformed - Time Domain Signal')
subplot(3,1,3)
plot(TDR(1:64));xlabel('time samples');ylabel('magnitude');title('Real part of IFFT transformed Time Domain Signal')

Bilder und multidimensionale FTs

In der medizinischen Bildgebung, Spektroskopie, Bildverarbeitung, Kryptographie und anderen Bereichen der Wissenschaft und Technik ist es oft der Fall, dass man mehrdimensionale Fourier-Transformationen von Bildern berechnen möchte. Dies ist in Matlab recht einfach: (mehrdimensionale) Bilder sind immerhin nur n-dimensionale Matrizen, und Fourier-Transformationen sind lineare Operatoren: Eine iterativ Fourier-Transformation entlang anderer Dimensionen. Matlab bietet fft2 und ifft2 , um dies in 2-d oder fftn in n-Dimensionen fftn .

Ein möglicher Fallstrick ist, dass die Fourier-Transformation von Bildern normalerweise "zentrisch geordnet" dargestellt wird, dh mit dem Ursprung des k-Raums in der Bildmitte. Matlab bietet den Befehl fftshift , um die Position der DC-Komponenten der Fourier-Transformation entsprechend fftshift . Diese Ordnungsschreibweise macht es wesentlich einfacher, übliche Bildverarbeitungstechniken auszuführen, von denen eine nachstehend veranschaulicht wird.

Nullfüllung

Eine "schnelle und schmutzige" Möglichkeit, ein kleines Bild auf eine größere Größe zu interpolieren, besteht darin, es Fourier-Transformationen, die Fourier-Transformation mit Nullen auffüllen und dann die inverse Transformation durchführen. Dies interpoliert effektiv zwischen jedem Pixel mit einer ähnlich geformten Basisfunktion und wird üblicherweise verwendet, um medizinische Bildgebungsdaten mit niedriger Auflösung zu vergrößern. Beginnen wir mit dem Laden eines integrierten Beispiels

%Load example image 
I=imread('coins.png'); %Load example data -- coins.png is builtin to Matlab 
I=double(I); %Convert to double precision -- imread returns integers 
imageSize = size(I); % I is a 246 x 300 2D image

%Display it
imagesc(I); colormap gray; axis equal; 
%imagesc displays images scaled to maximum intensity

Wir können jetzt die Fourier-Transformation von I erhalten. Um zu veranschaulichen, was fftshift bewirkt, wollen wir die beiden Methoden miteinander vergleichen:

% Fourier transform 
%Obtain the centric- and non-centric ordered Fourier transform of I
k=fftshift(fft2(fftshift(I))); 
kwrong=fft2(I); 

%Just for the sake of comparison, show the magnitude of both transforms: 
figure; subplot(2,1,1); 
imagesc(abs(k),[0 1e4]); colormap gray; axis equal; 
subplot(2,1,2); 
imagesc(abs(kwrong),[0 1e4]); colormap gray; axis equal;
%(The second argument to imagesc sets the colour axis to make the difference clear). 

Wir haben jetzt die 2D-FT eines Beispielbildes erhalten. Um es mit Nullen zu füllen, möchten wir jeden k-Raum nehmen, die Kanten mit Nullen füllen und dann die Rücktransformation nehmen:

%Zero fill
kzf = zeros(imageSize .* 2); %Generate a  492x600 empty array to put the result in 
kzf(end/4:3*end/4-1,end/4:3*end/4-1) = k; %Put k in the middle 
kzfwrong = zeros(imageSize .* 2); %Generate a  492x600 empty array to put the result in 
kzfwrong(end/4:3*end/4-1,end/4:3*end/4-1) = kwrong; %Put k in the middle 

%Show the differences again 
%Just for the sake of comparison, show the magnitude of both transforms: 
figure; subplot(2,1,1); 
imagesc(abs(kzf),[0 1e4]); colormap gray; axis equal; 
subplot(2,1,2); 
imagesc(abs(kzfwrong),[0 1e4]); colormap gray; axis equal;
%(The second argument to imagesc sets the colour axis to make the difference clear). 

An diesem Punkt ist das Ergebnis ziemlich unauffällig:

Wenn wir dann die Rücktransformationen nehmen, können wir sehen, dass (richtig!) Daten mit Nullfüllung eine sinnvolle Methode für die Interpolation darstellen:

% Take the back transform and view 
Izf = fftshift(ifft2(ifftshift(kzf))); 
Izfwrong = ifft2(kzfwrong);

figure; subplot(1,3,1); 
imagesc(abs(Izf)); colormap gray; axis equal; 
title('Zero-filled image'); 
subplot(1,3,2); 
imagesc(abs(Izfwrong)); colormap gray; axis equal;
title('Incorrectly zero-filled image'); 
subplot(1,3,3); 
imagesc(I); colormap gray; axis equal; 
title('Original image'); 
set(gcf,'color','w'); 

Beachten Sie, dass die mit Nullen gefüllte Bildgröße doppelt so groß ist wie die des Originals. Man kann in jeder Dimension mehr als einen Faktor von Null auffüllen, obwohl dies offensichtlich die Größe eines Bildes nicht willkürlich erhöht.

Hinweise, Tipps, 3D und mehr

Das obige Beispiel gilt für 3D-Bilder (wie sie beispielsweise häufig durch medizinische Bildgebungstechniken oder durch konfokale Mikroskopie erzeugt werden), erfordert jedoch beispielsweise, dass fft2 durch fftn(I, 3) wird. Aufgrund der etwas umständlichen Art des Schreibens von fftshift(fft(fftshift(... mehrmals)) ist es üblich, Funktionen wie fft2c lokal zu definieren, um lokal eine einfachere Syntax bereitzustellen.

function y = fft2c(x)

y = fftshift(fft2(fftshift(x)));

Beachten Sie, dass die FFT schnell ist, große multidimensionale Fourier-Transformationen jedoch auf einem modernen Computer immer noch Zeit benötigen. Es ist außerdem inhärent komplex: Die Größe des k-Raums wurde oben gezeigt, aber die Phase ist absolut lebenswichtig; Übersetzungen im Bildbereich entsprechen einer Phasenrampe im Fourierbereich. Es gibt mehrere weitaus komplexere Operationen, die man im Fourier-Bereich ausführen möchte, beispielsweise das Filtern hoher oder niedriger Ortsfrequenzen (durch Multiplizieren mit einem Filter) oder das Ausblenden diskreter Punkte, die dem Rauschen entsprechen. Dementsprechend gibt es eine große Menge von Community-generiertem Code für die Abwicklung gängiger Fourier-Operationen auf Matlabs Haupt-Community-Repository-Site File Exchange .