Diskrete Fourier-Transformation (DFT)

Current version of the page has been reviewed and is approved ().


Die DFT ist durch das Folgende Transformationspaar definiert

Hintransformation:

X(n) = \sum^{N-1}_{k=0} x(k)e^{-j2\pi n k / N}
                 (5.12)

und Rücktransformation:

x(k) = \frac{1}{N} \sum^{N-1}_{n=0} X(n)e^{j2\pi n k/ N}
                 (5.13)

Entwicklung aus der DTFT (Fourier-Transformation für Abtastsignale):

Um aus der nicht an einem Computer zu berechnenden DTFT (DTFT = discrete-time Fourier transform, diese gilt für unendlich ausgedehnte Zeitsignale) eine berechenbare Transformation zu generieren, sind zwei Schritte notwendig:

  1. Beschränkung des Zeitvektors auf N Werte
  2. Abtastung des Spektrums mit N Werten (N ist hier nur die übliche Lösung, auch ein anderer Wert ist möglich)

Der erste Teil wird durch eine Multiplikation mit einer so genannten Fenster- oder Window Funktion umgesetzt. Der Name ergibt sich aus der Eigenschaft, dass nur noch ein Ausschnitt gezeigt wird, in Anlehnung an ein Glasfenster, das uns nur einen Ausschnitt der Wirklichkeit zeigt (siehe Abbildung). Die Länge
N
ist im Folgenden die Blockgröße, da nur noch ein Block mit Daten verarbeitet wird.


Fig. 5.3 Veranschaulichung der Wirkung einer Rechteck-Fensterfunktion der Länge N = 16.

Zur Lösung des zweiten Problems bei der DTFT (kontinuierliches Spektrum) nutzt man die Abtastung des Spektrums  Durch die periodische Wiederholung des Spektrums ist es aber ausreichend, nur den Bereich von
0 \leq \omega \leq 2 \pi
genauer zu betrachten. Um ein diskretes Spektrum zu erhalten, unterteilen wir das Spektrum in  gleichförmige Abschnitte (es kann eine beliebige Anzahl verwendet werden, eine Zweierpotenz führt aber zu besonders effizienten Lösungen).

Zusammenhang zwischen DTFT und DFT

Um von der exakten Darstellung des Spektrums mittels DTFT auf die computerlösbare DFT zu kommen, wurden zwei Veränderungen vorgenommen. 

Zunächst ist es interessant, die vorgenommene Diskretisierung im Frequenzbereich zu untersuchen. Die Konsequenz der Abtastung ist bereits aus der normalen Abtastung im Zeitbereich bekannt. Diese Abtastung führt zu einer periodischen Wiederholung des Spektrums. Die Abtastung im Frequenzbereich hingegen führt zu einer periodischen Wiederholung im Zeitbereich. Wenn man also ein Signal analysiert (Fourier-Transformation) und zurücktransformiert (inverse Fourier-Transformation), ergibt sich eine periodische Wiederholung. Auch hier gibt es eine andere Interpretationsmöglichkeit. Bei der Analyse von periodischen Signalen mittels Fourier-Reihen ergeben sich diskrete Spektren. Ein diskretes Spektrum führt im Umkehrschluss also zu einer periodischen Zeitfunktion, wobei bei der DFT die Periode genau
N
ist.

Fig. 5.4 Veranschaulichung der erzwungenen Signalperiodizität durch die DFT.

Der Einfluss der Fenster-Funktion lässt sich zunächst nur durch die DTFT beschreiben. Es ist bekannt, dass eine Faltung im Zeitbereich zu einer Multiplikation im Frequenzbereich führt. Das Umgekehrte gilt aber auch. Eine Multiplikation im Zeitbereich führt im Frequenzbereich zu einer Faltung der Spektren, wobei die Faltung hier kontinuierlich als Integral zu definieren ist, da wir ja ein kontinuierliches periodisches Spektrum haben.

H(e^{j\Omega}) |_{0 \leq n \leq N } = \int^{\pi}_{ \theta = -\pi} H(e^{j \theta}) H^W (e^{j(\Omega - \theta)}) \text{d} \theta

wobei
H^W
die Übertragungsfunktion der Fensterfunktion ist.

Für das Rechteckfenster ergibt sich das Spektrum, indem wir zunächst die z-Transformation der Fensterfolge
w(k)
ansetzen. Es gilt

H^W (z) = \sum ^\infty _{k= -\infty} w(k) z^{-k} = \sum^{N-1}_{k=0} z^{-k}
.

Durch einsetzen von
z = e^{j \Omega}
ergibt sich das Spektrum

H^W (z) |_{z = e^{j\Omega}} = \sum ^{N-1}_{k=0} e^{-j \Omega k}
.

Mit der Formel der endlichen geometrischen Reihe erhält man das Spektrum der Rechteckfunktion

H^W(e^{j\Omega}) = \frac{1 - e^{-j\Omega N}}{1 - e^{-j\Omega}}


                   
= \frac{e^{-j\Omega N/2} \left( \underbrace{e^{j\Omega N/2} - e^{-j\Omega N/2}}_{2j\sin(\Omega N/2)} \right)} {e^{-j\Omega/2} \left( \underbrace{e^{j\Omega/2} - e^{-j\Omega/2}}_{2j\sin(\Omega/2)} \right)}


                   
= e^{-j(N-1)\Omega/2} \frac{\sin(N\Omega/2)}{\sin(\Omega/2)}


Der vordere Teil dieser Funktion entspricht einer linearen Phasenverschiebung. Der zweite Teil mit den Sinustermen stellt eine sogenannte Dirichlet-Funktion dar. Der Betrag dieser Funktion ist in Abbildung 5.5 gezeigt.
Fig. 5.5


Beispiel: DTFT und DFT bei einer Cosinus-Schwingung

Um die Unterscheide zwischen DTFT und DFT bei der Spektrumsberechnung eines Cosinus aufzuzeigen, muss zunächst das Spektrum einer abgetasteten Cosinusschwingung berechnet werden:
H(e^{j \Omega}) = \sum^\infty_{-\infty} \text{cos}(k \Omega_0)e^{-j\Omega k}

= \sum^\infty _{-\infty} \frac{1}{2} (e ^{j k \Omega_0} + e^{-j k \Omega_0}) e ^{-j \Omega k}

= \sum ^ \infty _ {- \infty} \frac{1}{2} ( e^{jk( \Omega_0 - \Omega)} + e ^{-jk(\Omega_0 + \Omega)} )

= \sum ^{\infty} _{- \infty} \frac{1}{2} e^ {-jk(\Omega - \Omega_0)} + \sum^{\infty} _{- \infty}\frac{1}{2} e^ {-jk(\Omega + \Omega_0)}

= \frac{1}{2} \delta^D ( \Omega - \Omega_0 ) + \frac{1}{2} \delta^D ( \Omega + \Omega_0)
            ,

mit
\Omega = 2 \pi f / f_s
. Das heißt, das Spektrum des abgetasteten Cosinus hat nur zwei definierte Frequenzwerte, an den Frequenzen
\Omega = \pm \Omega_0
.

Um den Einfluss der Annäherung durch die DFT zu veranschaulichen, begrenzen wir den betrachteten Signalausschnitt auf
N
Datenwerte. Dies führt wie in Gleichung (5.14) gezeigt zu einer Faltung mit dem Spektrum des Rechtecks. Durch die Siebeigenschaft der
\delta
-Funktion ergibt sich

H ( e^{j \Omega}) = ( \frac{1}{2} \delta ( \Omega - \Omega_0) + \frac{1}{2} \delta (\Omega + \Omega_0)) * H^W (e(j \Omega))

               
= \frac{1}{2} H ^W ( e^ {j(\Omega - \Omega_0)}) + \frac{1}{2} H^W ( e^{j(\Omega + \Omega_0)}).

Das Spektrum des Rechtecks wird also um
\pm \Omega_0
verschoben. Abbildung 5.6 zeigt das resultierende Spektrum für
N = 16
bei einer Grundfrequenz von
f_0 = 155
Hz und einer Abtastrate von
f_s = 1
kHz.
Fig. 5.6 Zeitsignale und Spektren eines abgetasteten Cosinus-Signals

Um nun das DFT-Spektrum zu berechnen, ist zusätzlich die Abtastung im Frequenzbereich notwendig. Dies hat implizit zur Folge, dass die Cosinusfolge periodisch wiederholt wird. Das abgetastete Spektrum ist auf der unteren linken Seite in Abbildung 5.6 zu sehen und die resultierende Folge auf der rechten Seite.

Das resultierende Spektrum hat nicht das erwartete Maximum bei
f_0
, da 155 Hz nicht im Abtastraster einer 16 Punkte FFT bei einer Abtastrate von 1 kHz liegt. Statt dessen ist die Leistung des Cosinus über mehrere Spektrallinien verteilt. Dieser Effekt wird auch als Leakage bezeichnet.

Symmetrien der DFT

Eine der wichtigsten Symmetrien für die diskrete sowie für die zeitdiskrete Fourier-Transformation (DFT und DTFT) ergibt sich für reelle Signale. Aus den Eigenschaften der z-Transformation ist deutlich, dass sich nur dann reelle Koeffizienten bei Signalen und Systemen ergeben, wenn die Pole und Nullstellen konjugiert komplex auftreten. Daraus folgt, dass die Spektren reeller Signale konjugiert komplex sind. Es ergibt sich also eine Symmetrie der Spektren an der Null-Hertz-Linie (Gleichstrom). Der Realteil ist dabei achsensymmetrisch (gerade) und der Imaginärteil punktsymmetrisch (ungerade). Diese Symmetrien ergeben sich auch, wenn Betrag und Phase betrachtet werden. Der Betrag reeller Funktionen ist immer achsensymmetrisch und die Phase punktsymmetrisch (siehe Abbildung 5.7).

Fig. 5.7 Beispiel eines Frequenz- und Phasenverlaufs eines reellwertigen Systems.

Die Symmetrien lassen sich bei der Berechnung der DFT und IDFT (inverse DFT) ausnutzen, da immer nur das halbe Spektrum berechnet werden muss, während die andere Hälfte durch Spiegelung erzeugt werden kann. Häufig kommt es vor, dass man im Spektralbereich etwas berechnet und an der dazugehörigen Zeitfolge interessiert ist (zum Beispiel bei einen Filterentwurf). Es reicht aus, für eine
N
-Punkte DFT,
N / 2 + 1
Spektralwerte zu kennen. Man benötigt einen Wert mehr, da die DFT für die höchste Frequenz bei
f_s / 2
nur einen reellwertigen Koeffizienten besitzt. 

Wie lässt sich die Symmetriebedingung zeigen?

Ausgehend von

X(n) = \sum^{N-1}_{k=0} x (k) e ^{-2\pi k n / N}


ergibt sich für die an der y-Achse gespiegelte Folge durch Variablensubstitution

X(-n) = \sum^{N-1}_{k=0} x(k) e^{j2 \pi k n / N}
           .

Konjugiert man dieses Signal ergibt sich

X^{*} (-n) = \sum^{N-1}_{k=0} x^{\ast}(k)e^{-j2 \pi kn / N}
             ,

wobei der negative Exponent zu beachten ist. Für ein reellwertiges Signal gilt
x^{\ast}(k) = x(k)
, woraus folgt, dass
X^{\ast}(-n) = X(n)
x

 Spektren reeller gerader Folgen

Die Symmetriebedingung reeller gerader Folgen kann aus den vorherigen Überlegungen geschlossen werden. Da gerade Funktionen achsensymmetrisch sind, der Imaginärteil aber punktsymmetrisch, können wir schließen, dass reelle, gerade Funktionen ein reelles, gerades Spektrum haben.

Rechenregeln

Linearität

Die DFT ist eine lineare Transformation. Es gilt also das

a_1 x_1(k) + a_2 x_2(k); \quad \circ\! \text{--} \!\bullet \quad a_1 X_1(n) + a_2 X_2(n).
          .      (5.15)

Theorem von Parseval

Das Theorem von Parseval sagt aus, dass man die Energie eines Signals im Zeit, oder im Frequenzbereich berechnen kann, bzw. dass die Leistung eines Signals im Zeit- und Frequenzbereich gleich ist. Für die DTFT gilt

\sum^{\infty}_{k=-\infty} x^2(k) = \frac{1}{2\pi} \int^{\pi}_{-\pi} |X(e(j\Omega))|^2 d\Omega
         ,     (5.16)

bzw. für die DFT

\sum^{N-1}_{k=0} x^2(k) = \frac{1}{N} \sum^{N-1}_{n=0} |X(n)|^2
         .     (5.17)


Effiziente Implementierung (Fast Fourier Transformation)

Die DFT lässt sich durch Ausnutzung unterschiedlicher Symmetrien sehr effizient berechnen. Um dies zu verdeutlichen, soll das sogenannten Decimation in Time-Verfahren zur drastischen Reduktion der benötigten Rechenleistung genauer gezeigt werden. Andere Verfahren können in der sehr umfangreichen Literatur zur Entwicklung der Fast Fourier Transform (FFT) gefunden werden.

Um eine effiziente Realisierung zu finden, legen wir die Länge der FFT als Zweierpotenz fest. Besonders häufig in der Audio- und Sprachsignalverarbeitung genutzte FFT-Längen sind 256, 512, 1024 und 2048. Durch diese Forderung ist es möglich, die Folge in zwei Teilfolgen zu zerlegen, wobei wir immer abwechselnd die Elemente der Folge den jeweiligen neuen Teilfolgen zuordnen. Es ergeben sich die Folgen
x(2k)
und
x(2k+1)
(siehe Abbildung 5.9).

Fig. 5.9 Aufteilung einer Sequenz in zwei Teilsequenzen zur Erklärung der FFT (Decimation in Time)

Die Konsequenzen für die DFT lassen sich nun wie folgt berechnen:

\sum_{k=0}^{N-1} x(k)e^{-j 2\pi nk/N} = \sum_{k=0}^{N/2 - 1} \underbrace{x(2k)}_{u(k)} e^{-j 2\pi n 2k/N} + \sum_{k=0}^{N/2 - 1} \underbrace{x(2k+1)}_{v(k)} e^{-j 2\pi n (2k+1)/N}


= \underbrace{\sum_{k=0}^{N/2 - 1} u(k) e^{-j 2\pi nk/(N/2)}}_{U(n) \quad \textbf{N/2 DFT}} + e^{-j 2\pi n /N} \underbrace{\sum_{k=0}^{N/2 - 1} v(k) e^{-j 2\pi nk/(N/2)}}_{V(n) \quad \textbf{N/2 DFT}}


X(n) = U(n) + e^{-j 2\pi n /N} V(n)


Die N-Punkte DFT lässt sich also in zwei N/2-Punkte DFT zerlegen. Hierbei tritt jetzt das Problem auf, dass die Folgen
U(n)
und
V(n)
nur bis
N/ 2
definiert sind, da ja auch nur eine N/2 DFT durchgeführt wurde. Die Lösung für dieses Problem ist durch die Periodizität der DFT aber sehr einfach zu umgehen, da sich das Spektrum immer wiederholt.

Aber warum stellt es einen Vorteil dar, wenn man die DFT so zerlegen kann? Um eine diskrete Frequenz zu berechnen, sind
N
komplexe Multiplikationen nötig. Dieser Schritt muss für alle diskreten Frequenzen durchgeführt werden. Die Berechnung des vollständigen Spektrums benötigt also
N^2
   Multiplikationen. Teilen wir die Aufgabe in zwei Teilspektren benötigt man
2(N/2)^2 +N
Multiplikationen für den Drehfaktor vor
V(n)
. Im Vergleich ergeben sich für
N = 8
   einmal 64 Multiplikationen und für die aufgeteilten Spektren 40 Multiplikationen. Der Schritt der Aufteilung kann nun solange wiederholt werden, bis die Folge nicht weiter aufgeteilt werden kann
(N = 2)
. Zusätzlich können einige Multiplikationen vernachlässigt werden, da
e^{j0} = 1
   ist. Eine Reduktion auf
\frac{N}{2}(\text{ld}(\frac{N}{2}))
ist so möglich (ld = Logarithmus zur Basis 2, logarithmus dualis). Somit ergibt sich eine im Rechenaufwand stark reduzierte DFT, die als Fast Fourier Transform (FFT) bekannt ist.