Transformationen

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

Außerdem hat die Fourier­transformation Anwendungen in der Numerik; z.B. lässt sich eine Polynom­multiplikation schneller durchführen, wenn die Polynome in Stützstellen­darstellung statt in Koeffizientendarstellung vorliegen.

Das Verfahren der schnellen Fourier­transformation (engl.: Fast Fourier Transform – FFT) hat eine Zeit­komplexität von Θ(n log(n)). Dadurch ist die Polynom­multiplikation sogar einschließ­lich Trans­formation in Stützstellen­darstellung und Rück­trans­formation noch schneller als die direkte Multi­plikation in Koeffizientendarstellung.

Grundlagen

Allgemein gilt in ℂ:

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

ist n-te Einheits­wurzel; ist k teilerfremd zu n, so ist w primitiv.

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

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

Eigen­schaften von n-ten Einheits­wurzeln

1. Es gibt in ℂ genau n verschiedene n-te Einheits­wurzeln, diese sind darstellbar als die Potenzen einer primitiven n-ten Einheits­wurzel w:

   w⁰, w¹, w², ..., w^n-1

    

2. Jede ganzzahlige Potenz w^k einer n-ten Einheits­wurzel w ist wieder n-te Einheits­wurzel, denn

   (w^k)ⁿ  =  w^k·n   =  (wⁿ)^k  =  1^k  =  1

   Dies gilt auch für negative k.

    

3. Ist n gerade, so gilt für jede primitive n-te Einheits­wurzel w:

   w^n/2  =  -1,

   denn (w^n/2)² = wⁿ = 1, d.h. w^n/2 ist 2-te Einheits­wurzel, also 1 oder -1. Da aber  w^n/2 ≠ 1 ist, da w primitiv ist, gilt w^n/2 = -1.

    

4. Das Quadrat w² einer primitiven n-ten Einheits­wurzel (n gerade) ist primitive n/2-te Einheits­wurzel, 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.

    

5. Ist w primitive n-te Einheits­wurzel, so ist w^-1 ebenfalls primitive n-te Einheits­wurzel, 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.

    

6. In ℂ gilt für die zu einer n-ten Einheits­wurzel 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)2 + sin(k·2π/n)2
   		=		1.

Diskrete Fouriertransformation

Die diskrete Fourier­transformation ist definiert als eine Vektor-Matrix-Multi­plikation.

Beispiel:  Sei n = 4. Dann ist i primitive n-te Einheits­wurzel. Die zugehörige Fourier­matrix 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 Fourier­transformation des Vektors a = [1 1 1 0] ergibt

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

Inverse Fourier­transformation

Beispiel:  Sei n = 4. Die inverse Fourier­matrix ist

F -1  =  1/4 ·

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

Die inverse Fourier­transformation des Vektors y = [3 i 1 -i] ergibt

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

Inter­pretation als Polynom­auswertung

Offenbar ist die Fourier­matrix F die Vandermonde-Matrix des Vektors w⁰, ..., w^n-1. Die Matrix-Vektor-Multi­plikation y = a·F lässt sich somit als Polynom­auswertung an den Stellen w⁰, ..., w^n-1 auf­fassen, wobei der Vektor a die Koeffizienten des Polynoms enthält (Koeffizient von x⁰ zuerst). Das Ergebnis y₀, ..., y_n-1 ist die Stützstellen­darstellung des Polynoms.

Mit der inversen Fourier­transformation wird das Polynom von der Stützstellen­darstellung wieder in die Koeffizientendarstellung zurück­verwandelt.

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

Schnelle Fouriertransformation (FFT)

Die Idee des Verfahrens der schnellen Fourier­transformation ist, die einzelnen Berechnungen der Matrix-Vektor-Multi­plikation y = a·F in einer speziellen Reihenfolge auszuführen, sodass jeweils auf schon berechnete Zwischen­ergebnisse zurück­gegriffen werden kann. Dabei werden die o.a. Eigen­schaften (3) und (4) ausgenutzt, nämlich dass w^n/2 = -1 ist und dass das Quadrat w² der primitiven n-ten Einheits­wurzel w primitive n/2-te Einheits­wurzel ist. Das Verfahren setzt voraus, dass n eine Zweierpotenz ist.

Komponenten von y mit geradem Index

Zunächst werden die Komponenten von y mit geradem Index berechnet, indem der Vektor a mit den ent­sprechenden Spalten der Fourier­matrix multipliziert wird. Es gilt für alle k ∈ {0, ..., n/2-1}:

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

Die Summe wird in zwei Hälften aufgespalten:

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

Nun ist aber

w^(i+n/2)·2k  =  w^i·2k + nk  =  w^i·2k·w^nk  =  w^i·2k

da w^nk = 1 ist, und damit

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

Mit m = n/2 sowie v = w², v primitive m-te Einheits­wurzel, gilt:

y_k'  =   _{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 Fourier­transformation des Vektors

(a_i + a_i+m) _{i = 0, ..., m-1}

der Länge m.

Komponenten von y mit ungeradem Index

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

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

Wiederum wird die Summe in zwei Hälften aufgespalten:

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

Nun ist aber

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

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

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

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

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

y_k''  =   _{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 Fourier­transformation des Vektors

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

Rekursive Ausführung

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

Analyse

Um die Fourier­transformation 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 Multi­plikationen erforderlich sowie nochmals m Multi­plikationen, um jeweils wⁱ aus w^i-1 zu berechnen. Insgesamt ergibt dies also 2m = n Additionen und 2m = n Multi­plikationen.

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

Die Ergebnis­vektoren 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 Zeit­komplexität T(n) der FFT ergibt sich somit als

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

T(1) = 0.

Die Auf­lösung dieser Rekursion ergibt

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

T(n) ∈ Θ(n log(n)).

Mit der Zeit­komplexität von T(n) ∈ Θ(n log(n)) ist die schnelle Fourier­transformation erheblich schneller als die diskrete Fourier­transformation in Form der Vektor-Matrix-Multi­plikation mit der Zeit­komplexität von T(n) ∈ Θ(n²).

Programm

Es folgt eine Implementierung in der Programmier­sprache Python.

Die Funktion fft berechnet die Fourier­transformation eines Vektors a, beginnend beim Index lo und der Länge n. Der Parameter w steht für die primitive n-te Einheits­wurzel.

Die Funktion eignet sich für die Berechnung der Fourier­transformation sowohl im Körper ℂ der komplexen Zahlen als auch in einem endlichen Körper ℤ_p. Für Berechnungen im Körper ℤ_p wird die Python-Klasse ModInt verwendet.

In welchem Zahlentyp gerechnet wird, also complex oder ModInt, wird vom Typ des Parameters w bestimmt,. Die Zahl 1 in dem jeweiligen Zahlentyp wird durch Berechnung von w⁰ erzeugt.

Die Funktion shuffle_ verschränkt die beiden, von den rekursiven Aufrufen von fft erzeugten Hälften reiß­verschluss­artig ineinander. Die ent­sprechende Permutation für n = 8 lautet

Mit der folgenden Funktion fftmain wird ein Vektor a der Länge n trans­formiert.

Die folgende Funktion complexroot erzeugt eine komplexe primitive n-te Einheits­wurzel.

Inverse Fouriertransformation

Die inverse Fouriertransformation lässt sich mit demselben Verfahren durchführen. Aufgrund der Definition der inversen Fourier­matrix F ^-1 wird jedoch statt mit der primitiven n-ten Einheits­wurzel w mit der inversen n-ten Einheits­wurzel w^-1 gearbeitet. In ℂ ist dies die konjugiert komplexe n-te Einheits­wurzel w (vgl. o.a. Eigen­schaften (5) und (6)). Ferner werden die Elemente der invers zu trans­formierenden Folge zunächst durch n geteilt.

In der folgenden Implementation werden verschiedene Vorausberechnungen durchgeführt, um die Zahlen 0, 1 und n im jeweiligen Zahlentyp (complex oder ModInt) zu erzeugen.

Anwendungsbeispiel

Das folgende Programm trans­formiert einen Vektor der Länge 8 von ModInt-Zahlen hin und wieder zurück.

Zusammenfassung

Die diskrete Fourier­transformation lässt sich inter­pretieren als Trans­formation eines Polynoms vom Grad n-1 von der Koeffizientendarstellung in die Stützstellen­darstellung. Als Stützstellen werden die n-ten Einheits­wurzeln w⁰, w¹, ..., w^n-1 verwendet. Im Körper ℂ der komplexen Zahlen sind die Werte

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

die n-ten Einheits­wurzeln.

Durch Ausnutzung der Symmetrien der n-ten Einheits­wurzeln lässt sich die Fourier­transformation beschleunigen; das Verfahren ist die schnelle Fourier­transformation FFT.

Mit Hilfe der FFT lässt sich die Zeit­komplexität der Polynom­multiplikation von Θ(n²) auf Θ(n log(n)) verbessern.

Literatur

[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 12]   H.W. Lang: Algorithmen in Java. 3. Auflage, Oldenbourg (2012)

Den FFT-Algorithmus finden Sie auch in meinem Buch über Algorithmen.
Weitere Themen des Buches: Sortieren, Textsuche, Graphenalgorithmen, Arithmetik, Codierung, Kryptografie, parallele Sortieralgorithmen.

[Weitere Informationen]