Transformationen

Schnelle Fouriertransformation (FFT)

Die Fouriertransformation ist ein fundamentales Verfahren in der Signalverarbeitung. Durch die Fouriertransformation lassen sich Signale von der Darstellung {(Zeitpunkt, Abtastwert)} in die Darstellung {(Frequenzanteil, Amplitude, Phase)} überführen. Viele Operationen, z.B. Filter, lassen sich im Frequenzraum leichter durchführen. Anschließend wird das Signal mit der inversen Fouriertransformation wieder zurück transformiert. [Beispiel]

Außerdem hat die Fouriertransformation Anwendungen in der Numerik; z.B. lässt sich eine Polynommultiplikation schneller durchführen, wenn die Polynome in Stützstellendarstellung statt in Koeffizientendarstellung vorliegen.

Das Verfahren der schnellen Fouriertransformation (engl.: Fast Fourier Transform – FFT) hat eine Zeitkomplexität von O(n log(n)). Dadurch ist die Polynommultiplikation sogar einschließlich Transformation in Stützstellendarstellung und Rücktransformation noch schneller als die direkte Multiplikation in Koeffizientendarstellung.

Definition: Sei komplexe Zahlen der Körper der komplexen Zahlen. Ein Element w Element heißt n-te Einheitswurzel, wenn wⁿ = 1 ist, und w heißt primitive n-te Einheitswurzel, wenn wⁿ = 1 ist, aber w^k ungleich 1 für alle k Element {1, ..., n-1}.

Beispiel: Sei n = 4. Dann ist i primitive 4-te Einheitswurzel¹). Alle 4-ten Einheitswurzeln sind:

i⁰ = 1, i¹ = i, i² = -1, i³ = -i.

Für n = 6 ist cos(2π/6) + i sin(2π/6) primitive 6-te Einheitswurzel (Bild 1).

Allgemein gilt in komplexe Zahlen :

w = cos(k·2π/n) + i sin(k·2π/n) mit k Element {0, ..., n-1}

ist n-te Einheitswurzel; ist k teilerfremd zu n, so ist w primitiv.

In der eulerschen Schreibweise lässt sich eine n-te Einheitswurzeln w auch ausdrücken als

w = e^ik2π/n mit k Element {0, ..., n-1}.

Bild 1: Einheitskreis in der gaußschen Zahlenebene mit 4-ten und 6-ten Einheitswurzeln

Eigenschaften von `n`-ten Einheitswurzeln

Es gibt in genau n verschiedene n-te Einheitswurzeln, diese sind darstellbar als die Potenzen einer primitiven n-ten Einheitswurzel w:
w⁰, w¹, w², ..., w^n-1.
Jede ganzzahlige Potenz w^k einer n-ten Einheitswurzel w ist wieder n-te Einheitswurzel, denn
(w^k)ⁿ = w^k·n = (wⁿ)^k = 1^k = 1.

Dies gilt auch für negative k.
Ist n gerade, so gilt für jede primitive n-te Einheitswurzel w:
w^n/2 = -1,

denn (w^n/2)² = wⁿ = 1, d.h. w^n/2 ist 2-te Einheitswurzel, also 1 oder -1. Da aber w^n/21 ist, da w primitiv ist, gilt w^n/2 = -1.
Das Quadrat w² einer primitiven n-ten Einheitswurzel (n gerade) ist primitive n/2-te Einheitswurzel, denn
1. (w²)^n/2 = wⁿ = 1.
2. Angenommen, w² sei nicht primitiv, dann existiert ein k {1, ..., n/2-1} mit (w²)^k = 1. Dann ist aber w^2k = 1 mit 2k < n, im Widerspruch dazu, dass w primitiv ist.
Ist w primitive n-te Einheitswurzel, so ist w^-1 ebenfalls primitive n-te Einheitswurzel, denn
1. (w^-1)ⁿ = w^-n = 1/wⁿ = 1/1 = 1.
2. Angenommen, w^-1 sei nicht primitiv, dann existiert ein k {1, ..., n-1} mit (w^-1)^k = w^-k = 1. Dann ist aber w^k = 1/w^-k = 1/1 = 1, im Widerspruch dazu, dass w primitiv ist.
In gilt für die zu einer n-ten Einheitswurzel w konjugierte Zahl w
w = w^-1

denn es ist

w · w = (cos(k·2π/n) + i sin(k·2π/n)) · (cos(k·2π/n) – i sin(k·2π/n))

= cos(k·2π/n)² + sin(k·2π/n)²

= 1.

Definition: Sei n Element natürliche Zahlen und w primitive n-te Einheitswurzel in komplexe Zahlen . Eine n kreuz n-Matrix F mit

F_i,j = w^i·j

für alle i, j Element {0, ..., n-1} heißt Fouriermatrix.²)

Die lineare Abbildung f : komplexe Zahlen ⁿ ⁿ mit

f(a) = a·F

für alle (Zeilen-)vektoren a Element komplexe Zahlen ⁿ heißt diskrete Fouriertransformation (DFT).³)

Beispiel: Sei n = 4. Dann ist i primitive n-te Einheitswurzel. Die zugehörige Fouriermatrix ist

F =

1	1	1	1
1	i	-1	-i
1	-1	1	-1
1	-i	-1	i

So ist etwa die -1 in der letzten Zeile der Matrix das Element

F_3,2 = w^3·2 = w⁶ = (-1)³ = -1

(die Elemente der Matrix werden von 0 bis n-1 indiziert).

Die Fouriertransformation des Vektors a = [1 1 1 0] ergibt

y = a·F = [3 i 1 -i]

Offenbar ist die Fouriermatrix F die Vandermonde-Matrix des Vektors w⁰, ..., w^n-1. Die Matrix-Vektor-Multiplikation y = a·F lässt sich somit als Polynomauswertung an den Stellen w⁰, ..., w^n-1 auffassen, wobei der Vektor a die Koeffizienten des Polynoms enthält. Das Ergebnis y₀, ..., y_n-1 ist die Stützstellendarstellung des Polynoms.

Satz: Die inverse Fouriermatrix F^-1 existiert und ist gleich

F^-1_i,j = 1/n · w^-i·j

für alle i, j Element {0, ..., n-1}. Die inverse Fouriermatrix enthält also die zu den Elementen der Fouriermatrix inversen Elemente, dividiert durch n.

Definition: Die lineare Abbildung f^-1 : komplexe Zahlen ⁿ ⁿ mit

f^-1(a) = a·F^-1

für alle a Element komplexe Zahlen ⁿ heißt inverse Fouriertransformation.

Beispiel: Sei n = 4. Die inverse Fouriermatrix ist

F^-1 = 1/4 ·

1	1	1	1
1	-i	-1	i
1	-1	1	-1
1	i	-1	-i

Die inverse Fouriertransformation des Vektors y = [3 i 1 -i] ergibt

a = y·F^-1 = [1 1 1 0].

In dieser Form ist die Fouriertransformation eine Matrix-Vektor-Multiplikation mit der Komplexität O(n²). Durch Ausnutzung der Symmetrie der n-ten Einheitswurzeln lässt sich die Berechnung auf O(n log(n)) beschleunigen. Dieses Verfahren heißt schnelle Fouriertransformation (Fast Fourier Transform – FFT) [CT 65].

Die Idee des Verfahrens der schnellen Fouriertransformation ist, die einzelnen Berechnungen der Matrix-Vektor-Multiplikation y = a·F in einer speziellen Reihenfolge auszuführen, so dass jeweils auf schon berechnete Zwischenergebnisse zurückgegriffen werden kann. Dabei werden die o.a. Eigenschaften (3) und (4) ausgenutzt, nämlich dass w^n/2 = -1 ist und dass das Quadrat w² der primitiven n-ten Einheitswurzel w primitive n/2-te Einheitswurzel ist. Das Verfahren setzt voraus, dass n eine Zweierpotenz ist.

Gegeben ist ein Vektor a = [a₀ ... a_n-1]. Berechnet werden soll y = a·F, d.h. gesucht ist der Vektor

y = [y₀ ... y_n-1]

Der Vektor y wird nun durch zwei ineinander verschränkte Vektoren y' und y'' halber Länge dargestellt:

y = [y'₀ y''₀ ... y'_n/2-1 y''_n/2-1].

D.h. die Komponenten von y an geraden Indexpositionen bilden den Vektor y' und diejenigen an ungeraden Indexpositionen den Vektor y''.

Die Vektoren y' und y'' werden getrennt voneinander berechnet. Zunächst werden die Komponenten von y' berechnet, indem die entsprechenden Komponenten des Vektors a mit den entsprechenden Spalten der Fouriermatrix multipliziert werden.

Im Folgenden sei m = n/2.

Zunächst werden die Komponenten von y mit geradem Index berechnet. Es gilt für alle k Element {0, ..., m-1}:

y_k' = y_2k = Summe _{i = 0, ..., n-1} a_i w^i·2k.

Die Summe wird in zwei Hälften zerlegt:

y_k' = Summe _{i = 0, ..., m-1} a_i w^i·2k + Summe _{i = 0, ..., m-1} a_i+m w^(i+m)·2k.

Nun ist aber

w^(i+m)·2k = w^{i·2k + m·2k} = w^i·2k·w^m·2k = w^i·2k·w^nk = w^i·2k,

da w^nk = 1 ist, und damit

y_k' = Summe _{i = 0, ..., m-1} (a_i + a_i+m) w^i·2k.

Mit v = w², v primitive m-te Einheitswurzel, gilt:

y_k' = Summe _{i = 0, ..., m-1} (a_i + a_i+m) v^i·k,

d.h. y_k' ist nichts anderes als die k-te Komponente der Fouriertransformation des Vektors

[a_i + a_i+m]_{i = 0, ..., m-1}

der Länge m.

Ähnlich werden die Komponenten von y mit ungeradem Index berechnet. Es gilt für alle k Element {0, ..., m-1} :

y_k'' = y_2k+1 = Summe _{i = 0, ..., n-1} a_i w^i·(2k+1).

Wiederum wird die Summe in zwei Hälften zerlegt:

y_k'' = Summe _{i = 0, ..., m-1} a_i w^i·(2k+1) + Summe _{i = 0, ..., n/2-1} a_i+m w^{(i+m)·(2k+1)}.

Nun ist aber

w^{(i+m)·(2k+1)} = w^{i·2k + m·2k + i + m} = w^i·2k · w^nk · wⁱ · w^n/2 = – wⁱ·w^i·2k,

da w^nk = 1 und w^n/2 = -1 ist.

Somit gilt:

y_k'' = Summe _{i = 0, ..., m-1} a_i wⁱ·w^i·2k + Summe _{i = 0, ..., m-1} – a_i+m wⁱ·w^i·2k

= Summe _{i = 0, ..., m-1} wⁱ·(a_i – a_i+m) w^i·2k.

Mit m = n/2 sowie v = w² gilt wiederum:

y_k'' = Summe _{i = 0, ..., m-1} wⁱ·(a_i – a_i+m) v^i·k,

d.h. y_k'' ist nichts anderes als die k-te Komponente der Fouriertransformation des Vektors

w_i·[a_i – a_i+m]_{i = 0, ..., m-1}.

Durch rekursive Anwendung dieses Verfahrens auf Vektoren jeweils halber Länge wird im Ergebnis die Fouriertransformation berechnet. Bild 2 zeigt schematisch den Ablauf der Berechnung für n = 8.

Bild 2: Datenfluss der schnellen Fouriertransformation

Um die Fouriertransformation eines Vektors a zu berechnen, ist zunächst ein Vektor a' zu berechnen, dessen Komponenten a_i' für i = 0, ..., m-1 wie oben gesehen folgende Werte haben:

a_i' = a_i + a_i+m sowie

a_i+m' = wⁱ·(a_i – a_i+m).

Hierfür sind m Additionen, m Subtraktionen und m Multiplikationen erforderlich sowie nochmals m Multiplikationen, um jeweils wⁱ aus w^i-1 zu berechnen. Insgesamt ergibt dies also 2m = n Additionen und 2m = n Multiplikationen.

Anschließend ist rekursiv auf die beiden Hälften von a' die Fouriertransformation anzuwenden.

Die Ergebnisvektoren y' und y'' stellen die geraden und die ungeraden Komponenten von y dar; sie sind noch ineinander zu verschränken (perfect shuffle), um den gewünschten Vektor y zu erhalten. Die Permutation perfect shuffle lässt sich, etwa mit der unten angegebenen Prozedur, in 3/2 n Schritten durchführen.

Die Zeitkomplexität T(n) der FFT ergibt sich somit als

T(n) = 3.5n + 2·T(n/2) sowie

T(1) = 0.

Die Auflösung dieser Rekursion ergibt

T(n) = 3.5n·log(n), d.h.

T(n) Element O(n log(n)).

Die folgende Prozedur berechnet die Fouriertransformation eines komplexen Vektors a, beginnend beim Index lo und der Länge n. Der Parameter w steht für die primitive n-te Einheitswurzel. Rechenoperationen mit komplexen Zahlen sind der Übersichtlichkeit halber mit normalen Rechenzeichen (+, -, *) dargestellt, obwohl diese in Java eigentlich nicht zur Verfügung stehen.

void fft(Complex[] a, int n, int lo, Complex w)
{
    int i, m;
    Complex z, v, h;

    if (n>1)
    {
        m=n/2;
        z=1;
        for (i=lo; i<lo+m; i++)
        {
            h=a[i]-a[i+m];
            a[i]=a[i]+a[i+m];
            a[i+m]=h*z;
            z=z*w;
        }
        v=w*w;
        fft(a, m, lo, v);
        fft(a, m, lo+m, v);
        shuffle (a, n, lo);
    }
}

Die Prozedur shuffle verschränkt die beiden, von den rekursiven Aufrufen von fft erzeugten Hälften reißverschlussartig ineinander. Die entsprechende Permutation für n = 8 lautet

0 1 2 3 4 5 6 7

0 4 1 5 2 6 3 7

Zur Ausführung der Permutation wird ein Hilfsarray b verwendet, in das zunächst die eine Hälfte der Folge ausgelagert wird.

void shuffle(Complex[] a, int n, int lo)
{
    int i, m=n/2;
    Complex[] b=new Complex[m];

    for (i=0; i<m; i++)
        b[i]=a[lo+i];
    for (i=0; i<m; i++)
        a[lo+i+i+1]=a[lo+i+m];
    for (i=0; i<m; i++)
        a[lo+i+i]=b[i];
}

Die inverse Fouriertransformation lässt sich mit demselben Verfahren durchführen. Aufgrund der Definition der inversen Fouriermatrix F^-1 wird jedoch statt mit der primitiven n-ten Einheitswurzel w mit der inversen n-ten Einheitswurzel w^-1 gearbeitet. In komplexe Zahlen ist dies die konjugiert komplexe n-te Einheitswurzel w (vgl. o.a. Eigenschaften (5) und (6)). Ferner werden die Elemente der invers zu transformierenden Folge zunächst durch n geteilt.

Die diskrete Fouriertransformation lässt sich interpretieren als Transformation eines Polynoms vom Grad n-1 von der Koeffizientendarstellung in die Stützstellendarstellung. Als Stützstellen werden die n-ten Einheitswurzeln w⁰, w¹, ..., w^n-1 verwendet. Im Körper komplexe Zahlen der komplexen Zahlen sind die Werte

w^k = cos(k·2π/n) + i sin(k·2π/n) = e^ik2π/n mit k Element {0, ..., n-1}

die n-ten Einheitswurzeln.

Durch Ausnutzung der Symmetrien der n-ten Einheitswurzeln lässt sich die Fouriertransformation beschleunigen; das Verfahren ist die schnelle Fouriertransformation FFT.

Mit Hilfe der FFT lässt sich die Zeitkomplexität der Polynommultiplikation von Θ(n²) auf Θ(n log(n)) verbessern.


[AHU 74]	A.V. Aho, J.E. Hopcroft, J.D. Ullman: The Design and Analysis of Computer Algorithms. Addison-Wesley (1974)
[CT 65]	J.M. Cooley, J.W. Tukey: An Algorithm for the Machine Calculation of Complex Fourier Series. Math. Comp. 19, 297-301 (1965)
[Sed 88]	R. Sedgewick: Algorithms. 2. Auflage, Addison-Wesley (1988)

[Lan 06]	H.W. Lang: Algorithmen in Java. 2. Auflage, Oldenbourg (2006) Den FFT-Algorithmus finden Sie auch in meinem Buch über Algorithmen. Weitere Themen des Buches: Sortieren, Textsuche, Parsen und Übersetzen, Graphenalgorithmen, Arithmetik, Codierung, Kryptografie, parallele Sortieralgorithmen. [Weitere Informationen] [Bestellen]