Inhaltsverzeichnis

1  Aliasing oder Verschiebung der zeitdiskreten Signale im Frequenzbereich  1

1.1  Mehrdeutigkeit der zeitdiskreten Signale  1

1.2  MATLAB Simulink-Experimente zum  

  Aliasing-Effekt  4

2  Aliasing-Experiment mit dem Digitalen-Oszilloskop TDS 220  8

2.1  Durchführung des Experiments  8

3  Audio-Aliasing Experiment  11

  A Spektrum der zeitdiskreten Signalen  14

  A 1 Die Dirac oder Delta-Funktion  14

  A 2 Die Abtastung als Multiplikation mit einer periodischen Delta-Folge  17

  Index  20

  Index der Abbildungen  22

Kapitel 1

Aliasing oder Verschiebung der

zeitdiskreten Signale im

Frequenzbereich

Die Abtastwerte eines gleichmäßig abgetasteten Signals kann man nicht direkt mit einem kontinuierlichen Signal assoziieren [1], [2], weil unendlich viele Signale die gleichen Abtastwerte ergeben können [3]. Nur wenn die Bandbreite des kontinuierlichen Signals bis f max = f s /2 begrenzt ist, gibt es eine eindeutige Zuordnung zwischen dem kontinuierlichen Signal und den Abtastwerten. Mit f s wurde die Abtastfrequenz bezeichnet.

Mit anderen Worten ausgedrückt, gibt es bei diskreten, gleichmässig abgetasteten Signalen eine Mehrdeutigkeit, die im nächsten Abschnitt untersucht wird.

1.1 Mehrdeutigkeit der zeitdiskreten Signale

Ein sinusförmiges Signal der Frequenz f 0 und Amplitude A x der Form

x(t) = A x sin(2πf 0 t)

wird mit einer Abtastfrequenz f s = 1/T s abgetastet, was zu einer Zeitdiskretisierung führt:

t = kT s mit k = −∞, . . . , −2, −1, 0, 1, 2, . . . , ∞ Die zeitdiskreten Abtastwerte werden dann:

x[kT s ] = A x sin(2πf 0 kT s ) = A x sin(2πf 0 kT s + m2π)

Gerade die Periodizität der stationären, sinusförmigen (oder cosinusförmigen) Komponenten, in dieser Gleichung mit dem Term m2π, m ∈ Z dargestellt, zusammen mit der gleichmässigen Abtastung führt zu einer Mehrdeutigkeit. Die obige Gleichung kann auch in folgender Form m x[kT s ] = A x sin(2π(f 0 +

1.1 Mehrdeutigkeit der zeitdiskreten Signale 2

geschrieben werden, wobei für jeden Wert k ein entsprechender Wert für m gewählt werden kann, so dass das Verhältnis n = m/k immer eine ganze Zahl wird. Die innere Klammer der letzten Form in Gl. (1.3) f n = f 0 + nf s für n = 0, ±1, ±2, ±3, . . .

stellt die Frequenzen f n aller Komponenten dar, die zu den gleichen Abtastwerten wie die Komponente der Frequenz f 0 führen.

Mit einem Zahlenbeispiel soll dieses Ergebnis näher erläutert werden. Angenommen eine kontinuierliche sinusförmige Komponente der Frequenz f 0 = 125 kHz wird mit einer Abtastfrequenz von f s = 100 kHz abgetastet. Die Tabelle 1.1 zeigt die Frequenzen einiger Komponenten, die zu gleichen Abtastwerte führen, mit der Voraussetzung, dass alle die gleiche Amplitude A x besitzen.

Die negativen Frequenzen sollen hier nicht stören. Eine negative Frequenz bedeutet nur, dass das Argument der entsprechenden sinusförmigen Komponente negativ ist. Mit einer Nullphase von ±π kann das Argument so geschrieben werden, dass die Frequenz positiv wird. Als Beispiel, für die Komponenten der Frequenz −175000 Hz würde das zu A x sin(2π(−175000)kT s ) = A x sin(2π(175000)kT s ± π) führen.

Abb. 1.1: Amplitudenspektrum der Komponenten mit gleichen Abtastwerten

Abb. 1.1 stellt in Form eines Amplitudenspektrums die Komponenten der gleichen Abtastwerten aus der Tabelle 1.1 dar und Abb. 1.2 zeigt die drei Komponenten für n = −1, 0, 1 und deren gleiche Abtastwerte. Die tiefste Frequenz entspricht der verschobenen Komponente der Frequenz 25 kHz, die höchste Frequenz gehört der Komponente

1.1 Mehrdeutigkeit der zeitdiskreten Signale 3

Abb. 1.2: Drei Komponenten (n=-1,0,1) und deren gleiche Abtastwerte

der Frequenz 225 kHz und die verbleibende entspricht der ursprünglichen Komponente der Frequenz von 125 kHz.

Der Frequenzbereich von f = 0 bis f = f s /2 definiert den so genannten Nyquist-Bereich. Aus gegebenen Abtastwerten kann immer eine kontinuierliche Komponente rekonstruiert werden. Wenn die Frequenz des ursprünglichen, kontinuierlichen Signals im Nyquist-Bereich liegt, dann ergibt die Rekonstruktion aus dessen Abtastwerten das ursprüngliche Signal.

Das Signal der höchsten Frequenz, das mit diskreten Abtastwerten darstellbar ist, besteht aus einem positiven Abtastwert gefolgt von einem negativen Abtastwert. Die Periode dieses Signals ist T min = 2T s und somit ist die höchste Frequenz eines kontinuierlichen Signals, das eindeutig mit Abtastwerten darstellbar ist, durch

f max ≤ f s /2

gegeben (Abb. 1.3). Diese Bedingung stellt das Abtasttheorem dar. Es wurde angenommen, dass die Abtastwerte nicht an den Nulldurchgängen des kontinuierlichen Signals erfasst werden. Damit das nie geschehen kann, muss die Bedingung immer mit f max < f s /2 erfüllt sein.

Abb. 1.3: Signal der höchsten Frequenz, das mit Abtastwerten darstellbar ist

Für alle anderen kontinuierlichen Komponenten mit Frequenzen außerhalb des Nyquist-Bereichs, welche die gleichen Abtastwerte besitzen, stellt die Rekonstruktion ein kontinuierliches Signal dar, das der Aliased- oder der verschobenen Komponente entspricht.

1.2 MATLAB/Simulink-Experimente zum Aliasing-Effekt

In dem gezeigten Zahlenbeispiel, wird die aus den Abtastwerten rekonstruierte, kontinuierliche Komponente der Frequenz 25 kHz nicht mehr die ursprüngliche Komponente der Frequenz 125 kHz darstellen, sondern sie bildet die verschobene Komponente.

Auch im alltäglichen Leben kommen solche Sachverhalte vor. Bei einem Formel-1-Rennen sieht man plötzlich, dass die Räder stehen bleiben oder sich in umgekehrter Richtung drehen, obwohl sich der Wagen richtig bewegt. Abb. 1.4 zeigt wie das zustande kommt.

Die Bildsequenzen stellen Abtastbilder bei t, t + T s , t + 2T s , t + 3T s , t + 4T s dar, und in einer Sequenz wird das Rad mit den Speichern an einer Stelle und in der nächsten Sequenz an einer anderen Stelle erfasst.

Die Stellen können langsam vorrücken, obwohl sich die Räder vielmals gedreht haben. Es entsteht der Eindruck, dass die Räder zu langsam drehen (Abb. 1.4a). Wenn die Speicher in derselben Lage in jeder Sequenz erfasst werden, hat man den Eindruck die Räder drehen sich nicht (Abb. 1.4b). Kleine Abweichungen nach hinten führen zum Eindruck, dass die Räder umgekehrt drehen (Abb. 1.4c).

Es gibt auch andere ähnliche Beispiele, alle verbunden mit Abtastereignissen. So z.B. funktioniert das Stroboskop, das man für die Einstellung der Vorzündung eines Fahrzeugs einsetzt. Das pulsierende Licht des Stroboskops wird mit der Zündung synchronisiert, so dass die Stelle am Schwungrad, bei der die Zündung stattfindet, als stehend erscheint.

1.2 MATLAB/Simulink-Experimente zum

Aliasing-Effekt

Mit MATLAB- bzw. Simulink-Simulationen kann der Aliasing-Effekt einfach reproduziert werden [4]. Abb. 1.5 zeigt das Modell abtast_12.mdl das über das Programm abtast12.m initialisiert und Aufgerufen wird. Die Blöcke Sine Wave aus der Sources- Bibliothekakzeptieren keine negativen Frequenzen und um das zu umgehen, werden die sinusförmigen Signale mit Hilfe von Fcn-Blöcken (Function-Blöcken), die als Eingang die Zeit haben, erzeugt.

1.2 MATLAB/Simulink-Experimente zum Aliasing-Effekt

Der obere Generator erzeugt ein Signal der Frequenz f 0 und der untere erzeugt eine Komponente der Frequenz f n = f 0 + nf s (mit beliebigen n-Werten), die gleiche Abtastwerte besitzt. Mit einem analogen Tiefpassfilter wird aus den Abtastwerten des Signals der Frequenz f 0 ein kontinuierliches Signal rekonstruiert. Ein ideales Tiefpassfilter mit Bandbreite bis f s /2 müsste daraus das ursprüngliche kontinuierliche Signal rekonstruieren, wenn dessen Frequenz im Nyquist-Bereich liegt, oder die Aliased-Komponente, wenn die Frequenz außerhalb liegt.

Abb. 1.5: Modell zur Simulation des Aliasing-Effekts (abtast12.m, abtast_12.mdl)

Abb. 1.6 zeigt oben das Signal der Frequenz f 0 = 15000 Hz und das Signal der Frequenz f n = −185000, das einem Wert n = −2 bei f s = 100000 Hz entspricht. Zusätzlich werden die Abtastwerte dieser Signale gezeigt. Sie erscheinen wegen der Zero-Order Hold-Blöcke (Halteglied Nullter-Ordnung), die zur Abtastung dienen, als Treppen. Wie man sieht, haben beide Signale dieselben Abtastwerte. Im Programm abtast12.m werden mit

% ------- Signal (Grundsignal)

f0 = 15000 phi_sig = 0;

% ------- Abtastfrequenz fs = 100000; % Abtastfrequenz Ts = 1/fs; % Abtastperiode % ------- Signal mit gleichen Abtastwerten n = -2 fn = f0+n*fs

1.2 MATLAB/Simulink-Experimente zum Aliasing-Effekt

Abb. 1.6: Signale mit gleichen Abtastwerten (abtast12.m, abtast_12.mdl) (abtast12.m, abtast_12.mdl)

phin = 0;

die Parameter der zwei Signalgeneratoren initialisiert. Danach wird die Zeit für die Simulation so bestimmt, dass drei Perioden des Signals mit der kleinsten Frequenz erfasst werden und es wird die Simulation aufgerufen: fd = min(f0,abs(fn));

delta_t = 1/(10*abs(fn));

n_perioden = 3;

sim(’abtast_12’,[0,n_perioden/fd]); In der Senke To Workspace werden die Signale als Structure with Time mit Name y eingefangen. Die Darstellung aus Abb. 1.6 ist mit folgender Programmsequenz entstanden:

nsig = length(y.signals.values); nx = fix((n_perioden/fd)/delta_t); nd = nsig-nx:nsig;

figure(1); clf;

subplot(211), plot(y.time(nd), y.signals.values(nd,1:4)); title([’Komponenten gleicher Abtastwerte und ... Abtastwerte (n = ’,...

num2str(n),’, f0 = ’,num2str(f0),’ Hz)’]); xlabel(’Zeit in s’); grid;

La = axis;

axis([min(y.time(nd)), max(y.time(nd)), -1.2, 1.2]);

1.2 MATLAB/Simulink-Experimente zum Aliasing-Effekt

subplot(212), plot(y.time(nd), y.signals.values(nd,5:6)); title([’Abtastwerte und rekonstruierte Komponente ... (Aliased-Komponente)’]);

xlabel(’Zeit in s’); grid;

La = axis;

axis([min(y.time(nd)), max(y.time(nd)), -1.2, 1.2]);

Abb. 1.7: Signale mit gleichen Abtastwerten (abtast12.m, abtast_12.mdl)

Im unteren Teil der Abb. 1.6 ist das kontinuierliche, über das analoge (nicht ideales) Tiefpassfilter rekonstruierte Signal, zusammen mit den Abtastwerten gezeigt. Wenn die Abtastwerte ein Signal mit einer Frequenz darstellen, die näher zur Grenze des Nyquist-Bereichs (f s /2) liegt, dann sind die Fehler der Rekonstruktion wegen des nicht idealen Tiefpassfilters mit Bandbreite f s /2 viel größer.

Abb. 1.7 zeigt die Signale für f 0 = 13000 Hz, f s = 100000 Hz und n = -1, die zu einer Aliased-Komponente von f n = −30000 Hz ergibt, die viel näher an f s /2 = 50000 Hz liegt. Das kontinuierliche rekonstruierte Signal im unteren Teil unterscheidet sich viel mehr vom Aliased-Signal aus dem oberen Teil der Abbildung.

Mit dem Block Transport Delay versucht man die Verspätung, die das Tiefpassrekonstruktionsfilter erbringt zu kompensieren, so dass das diskrete Signal und das rekonstruierte, kontinuierliche Signal überlappt dargestellt werden können. Die Variable delta_t wird im Fenster der Parametrierung der Simulation benutzt, um die maximale Schrittweite der Simulation zu begrenzen. Für die Darstellung aus Abb. 1.7 muss diese Schrittweite noch kleiner gemacht werden, (delta_t = 1/(10*abs(fn));

um schöne sinusförmige Signale zu erhalten. Die Anzahl der Perioden für die Darstellung wurde ebenfalls geändert auf n_perioden=5.

Kapitel 2

Aliasing-Experiment mit dem

Digitalen-Oszilloskop TDS 220

Das Tektronix TDS 220 Oszilloskop ist ein so genanntes digitales Speicher-Oszilloskop, das die Verläufe in Echtzeit aufnimmt. Wenn Einkanal-Betrieb verwendet wird, werden 2500 Abtastwerte des Eingangssignals aufgenommen, und die Hälfte davon, wenn beide Kanäle eingesetzt werden. Nachdem diese Werte im Speicher hinterlegt sind, werden sie bearbeitet und am Display angezeigt.

Abb. 2.1 zeigt diesen Sachverhalt für den Einkanal-Betrieb. Die neue Erfassung muss mit dem Signal synchronisiert sein, um ein stabiles Bild der Anzeige zu erhalten.

2.1 Durchführung des Experiments

Das Oszilloskop hat eine Bandbreite von 100 MHz und eine Abtastfrequenz bis 1 GS/s ¹ oder 1 GHz. Die höchste Abtastfrequenz von 1 GHz ist zehn mal größer als die Bandbreite von 100 MHz und somit gibt es 10 Abtastwerte in einer Periode für das Signal der höchsten Frequenz.

Mit Hilfe der Abtastfrequenz wird die Zeitablenkung bei diesen Oszilloskopen eingestellt. Der Bereich für diesen Parameter liegt zwischen 5 ns/Div bis 5 s/Div. Bei 10

2.1 Durchführung des Experiments 9

Divisionen in x-Richtung bedeutet das 50 ns bis 50 s für den ganzen Bildschirm. In dieser Zeit werden die schon erwähnten 2500 Abtastwerte erfasst und somit lässt sich die Abtastfrequenz berechnen:

f s = 2500/(10 ∗ a b ), Hz

Durch a b wurde die Ablenkungszeit in s/Div bezeichnet. Als Beispiel, bei a b = 5 ms/-Div erhält man eine Abtastfrequenz von f s = 50 kHz.

Für die Durchführung des Experiments benötigt man einen Signalgenerator mit sinusförmigen Ausgang, der eine gute Stabilität der Frequenz besitzt und gleichzeitig eine sehr feine Einstellung der Frequenz erlaubt.

Das Experiment wird mit dem HAMEG-Signalgenerator Typ HM 8130 durchgeführt. Man wählt z.B. eine Frequenz für das Signal des Generators von 1 MHz und die Möglichkeit diese in Schritte von 25 kHz zu ändern.

Das Oszilloskop wird zuerst so eingestellt, das man 5 oder 10 Perioden des Signals von 1 MHz am Bildschirm anzeigt. Die ganze Ablenkungszeit muss somit (5 ÷ 10)/10 ⁶ Sekunden sein. Bei 10 Divisionen bedeutet das eine Ablenkungszeit von 0, 5 ÷ 1µs/Div. Es wird danach die Ablenkungszeit vergrößert bis man den Wert von 5 ms/Div erreicht, der einer Abtastfrequenz von 50 kHz entspricht. Da die Ablenkungszeit und dadurch implizit auch die Abtastfrequenz nicht fein einstellbar ist, muss man jetzt die Frequenz des Generators so ändern, bis auf dem Bildschirm Aliased-Komponenten erscheinen, deren Frequenz vom Oszilloskop ermittelt und angezeigt wird. Die Frequenz dieser Komponente ist laut Gl. 1.4 durch f n = f 0 + nf s für n = 0, ±1, ±2, ±3, . . .

gegeben.

Für f s = 50 kHz (Zeitablenkung 5 ms/Div) kann man viele Aliased-Komponenten erzeugen. Mit n = −20 erhält man eine Scheinkomponente der Frequenz null. Dafür muss man die Frequenz des Generators um die Frequenz von 1 MHz sehr fein ändern können, weil der Ablenkungsoszillator des Oszilloskops nicht genau ist.

Abb. 2.2: Scheinfrequenzen für f s = 50 kHz und f 0 in der Umgebung von 1 MHz

Danach ergeben Frequenzänderungen in Schritte von 50 kHz beim Generator immer wieder diese Scheinkomponente der Frequenz null, weil diese Änderungen bei konstanter Abtastfrequenz von 50 kHz einer Änderung von n mit eins entspricht.

2.1 Durchführung des Experiments 10

Wenn man die Änderungsschritte jetzt auf 25 kHz reduziert, dann erhält man Scheinfrequenzen von null und von 25 kHz. Bei einer Frequenz des Generators von 975 kHz und n = −20 ist die Scheinfrequenz gleich 975 − 20 ∗ 50 = −25 kHz bzw. bei einer Frequenz des Generators von 1025 kHz ist die Scheinfrequenz 1025 − 20 ∗ 50 = 25 kHz. Wenn man jetzt die Frequenz des Generators mit mehr als der halben Abtastfrequenz ändert, z.B. nimmt man f 0 = 1030 kHz, dann ändert sich auch n auf n = −21 und die Scheinfrequenz wird −20 kHz. Dieser Sachverhalt ist in Abb. 2.2 dargestellt. Man kann sich viele andere Frequenzen des Generators f 0 und andere Abtastfrequenzen f s vorstellen, die zu Aliasing führen und mit einer ähnlichen Darstellung zu erklären sind.

Wenn man mit Simulink den beschriebenen Fall simuliert, stellt man fest, dass die Aliased-Komponenten relativ wenig Abtastwerte enthalten und somit ist nicht so leicht zu erkennen, dass sie von sinusförmigen Signalen hervorgehen. Das Oszilloskop stellt sie aber sinusförmig dar und zeigt dadurch, das die Abtastwerte der Komponenten die Aliased sind, interpoliert werden.

Kapitel 3

Audio-Aliasing-Experiment

Mit Hilfe eines Chirp-Signals, das seine Frequenz linear verändert, kann das Aliasing der abgetasteten Signale im Zeitbereich anschaulich gezeigt werden. Die Frequenzen der Signale sind so gewählt worden, dass die resultierenden Signale hörbar sind, wenn der PC, mit dem das Experiment durchgeführt wird, eine Sound-Karte enthält. Die Verschiebung ist dann leicht zu erkennen.

Abb. 3.1 zeigt das Simulink-Modell des Experiments. Das Chirp-Signal wird mit Hilfe des Blocks Fcn erzeugt, der als Eingang die Zeit t hat, die mit u im Block bezeichnet ist. Die Frequenz des Signals f (t) muss linear steigen, ausgehend von einer minimalen Frequenz f min = 0 Hz:

f (t) = µ f t + f ^min Die Geschwindigkeit der Änderung ist durch den Faktor µ f = 50 Hz/s gegeben. In 6 Sekunden erreicht die Frequenz 300 Hz.

Um ein Signal mit dieser Frequenzänderung zu erzeugen, wird die entsprechende Phase als Zeitveränderliche ermittelt:

Das Signal erhält man dann durch:

x(t) = ˆ

Dieser Ausdruck ist in dem Block Fcn für µ f = 50 Hz/s, f min = 0 Hz und ϕ(0) = π/3 eingetragen:

cos(pi*50*u^2+pi/3)

Mit einer Abtastfrequenz von f s = 100 Hz wird nur der erste Teil des Chirp-Signals bis zur Frequenz f = 50 Hz korrekt mit den Abtastwerten dargestellt. Dieser Teil entspricht der Zeit von null bis zu einer Sekunde. Danach für f > 50 Hz, wie z.B. f = 70 Hz, entsteht Aliasing und die Scheinfrequenz ist f n = 70 + n100 = −30 Hz bei n = −1. Das Minusvorzeichen stellt nur eine Phasenverschiebung von 180° relativ zum Teil in dem die Frequenz f = 30 Hz ist (f < f s /2) dar.

Wenn die Frequenz des Chirp-Signals sich dem Wert f s = 100 Hz bei 2 Sekunden nähert, dann nähert sich die Scheinfrequenz dem Wert null. Für Frequenzen oberhalb von f s , wie z.B. für f = 120 Hz, erhält man eine Scheinfrequenz von f n = 120+nf s = 20

Abb. 3.1: Simulink-Modell des Audio-Aliasing-Experiments (audio_alias1.mdl, audio_alias2.mdl)

Hz bei n = −1, erneut mit einer relativen Phasenlage gleich dem Zustand für f = 20 Hz.

Eine Frequenz des Chirp-Signals, die den Wert 150 Hz überschreitet, wie z.B. f = 170 Hz, führt zu einer Scheinfrequenz von f n = 170 + nf s = −30 Hz mit n = −2. Bei f = 200 Hz erhält man erneut eine Scheinfrequenz von null. Danach wiederholt sich die Scheinfrequenz für jedes Intervall nf s < f < (n + 1)f s .

Abb. 3.2 zeigt oben das Chirp-Signal mit einer Frequenzänderung von 0 Hz bis 300 Hz in 6 Sekunden, in der Mitte ist das abgetastete Signal gezeigt und unten das aus den Abtastwerten interpolierte Signal. Das letztere kann als das rekonstruierte, kontinuierliche Signal gelten. Der Interpolationsfaktor ist 10, was auch bedeutet, dass zwischen den Abtastwerten mit f s = 100 Hz noch 9 Werte interpoliert werden, die einer Abtastfrequenz von f s = 1000 Hz entspricht.

Die relativ kleinen Werte des abgetasteten Signals zu den Zeitpunkten 1, 3 und 5 Sekunden (Abb. 3.2 in der Mitte) ergeben sich durch die Phasenänderung von 180° des Signals von 50 Hz bei diesen Zeitmomenten.

Das digitale Tiefpassfilter für die Interpolation muss eine relative Bandbreite von f p /f s = 1/10 oder in der MATLAB-Konvention eine relative Frequenz von 2/10 = 0, 2 besitzen.

Das ganze Modell arbeitet digital mit einer Schrittweite der Simulation von (1/f s )/10 = 1/1000 Sekunden. In der Senke To Workspace werden die Signale in dem Feld y eingefangen und können danach beliebig bearbeitet werden. Mit

wavwrite(y(:,1), 1000, ’audio1.wav’);

wavwrite(y(:,2), 1000, ’audio2.wav’);

wavwrite(y(:,3), 1000, ’audio3.wav’);

Abb. 3.2: Chirp-Signal, abgetastetes und rekonstruiertes Signal

(audio_alias1.mdl, audio_alias2.mdl)

werden z.B. Wav ¹ -Dateien gebildet, die man danach abhören kann. Die gewählten Frequenzen sind relativ niedrig, um die Signale anschauend darstellen zu können. Das Modell kann einfach für höhere Frequenzen umgewandelt werden. Mit µ f = 250 Hz/s erhält man eine Änderung der Frequenz von 0 Hz bis 1500 Hz in 6 Sekunden. Für eine Abtastfrequenz von f s = 500 Hz erhält man jede 2 Sekunden eine Scheinfrequenz gleich null.

Die Interpolation mit dem gleichen Faktor von 10 führt zu einer Abtastfrequenz für die Interpolation von 500 × 10 = 5000 Hz, bei der das FIR-Tiefpassfilter (Block Discrete Filter) arbeiten muss. Die relative Bandbreite dieses Filters bleibt die gleiche. Die Wav-Dateien müssen dann durch wavwrite(y(:,1), 5000, ’audio1.wav’);

wavwrite(y(:,2), 5000, ’audio2.wav’);

wavwrite(y(:,3), 5000, ’audio3.wav’);

gebildet werden.

Anhang A

Spektrum der zeitdiskreten

Signalen

Es wird das Spektrum der zeitdiskreten Signalen über die Fourier-Transformation ermittelt. Der partikuläre Fall der periodisch, stationären, sinusförmigen Komponenten wurde im ersten Kapitel mit einfachen mathematischen Mitteln abgeleitet. Die Fourier-Transformation [5], [6] erweitert die Betrachtungen im Frequenzbereich auch für zeitdiskrete aperiodische Signale. Wenn man die Werkzeuge, die man bei kontinuierlichen Signalen einsetzt, auch für die zeitdiskreten Signale anwenden möchte, dann müssen die Abtastwerte als Dirac-Funktionen [5] angenommen werden. Nur so kann man alle Integralfunktionen (wie z.B. die Fourier-Transformation oder das Faltungsintegral) auch für die zeitdiskreten Signale benutzen.

A.1 Die Dirac oder Delta-Funktion

Ingenieurmäßig wird die Dirac- oder Delta-Funktion als Grenzfunktion eines Rechteckimpulses der Form (Abb. A.1a)

∆(t) =

eingeführt. Er ist ein Impuls der Dauer τ und Höhe 1/τ und somit besitzt er eine Fläche gleich eins. Die Grenzfunktion δ(t) = lim

bildet die so genannte Dirac- oder Delta-Funktion δ(t). Die Delta-Funktion hat also für alle t = 0 den Wert 0 und trotzdem gilt:

∞

Im Sinne der klassischen Funktionen und des Riemann-Integralls ist das nicht möglich. Mathematisch korrekt wird die Delta-Funktion als Distribution [1] eingeführt. Dieser Aspekt wird hier nicht mehr weiter verfolgt.

A.1 Die Dirac oder Delta-Funktion 15

Abb. A.1: Rechteckimpuls und dessen Grenzfunktion als Delta-Funktion

Zwei andere Definitionen der Delta-Funktionen sind sehr nützlich [1], [7]:

Mit einfachen MATLAB-Programmen können diese Integrale annähernd nachgebildet werden, um zu zeigen, dass sie zu Delta-Funktionen führen. Als Beispiel wird die eine Definition von δ(t) mit folgender Annäherung des Integrals programmiert:

Das Programm (delta1.m) begint mit der Wahl eines Bereichs für die Zeit und für die Frequenz:

% ------ Zeitbereich

tmax = 5; tmin = -tmax;

nt = 100; dt = (tmax-tmin)/nt;

t = tmin:dt:tmax;

% ------ Frequenzbereich

nf = 100;

f_max = 1/(2*dt); df = f_max/nf;

f = -f_max:df:f_max;

n_t = length(t); n_f = length(f);

Danach wird eine Matrix delta initialisiert, die in den Zeilen die Cosinusfunktionen für verschiedenen Frequenzen enthalten und in einer for-Schleife werden diese Zeilen gebildet:

A.1 Die Dirac oder Delta-Funktion 16

% ------ Delta-Funktion als Zeitfunktion delta = zeros(n_f, n_t); % Initialisierung for k = 1:n_f delta(k,:) = cos(2*pi*f(k)*t); end;

Die Annäherung des Integrals über die Frequenz ergibt sich durch die Summe entlang der Spalten mal ∆f :

delta_g = sum(delta)*df;

% Überprüfen der Fläche der Delta-Funktion flaeche = sum(delta_g)*dt Die Darstellung einiger Cosinusfunktionen und der angenäherten Delta-Funktion wird mit folgendem Programmabschnitt erhalten (Abb. A.2). figure(1); clf;

subplot(121), plot(t, delta(nf-5:nf+5,:)’); title(’Cosinuskomponenten’); xlabel(’Zeit’); grid;

subplot(122), plot(t, delta_g); title(’Annaeherung der Delta Funktion ...

xlabel(’Zeit’);

A.2 Die Abtastung als Multiplikation mit einer periodischen Delta-Folge 17

Für die graphische Darstellung der Delta-Funktion als Limes des Rechteckimpulses oder der gezeigten Integralen, wird ein vertikaler Pfeil der Länge eins benutzt (Abb. A.1b). Dementsprechend wird die Funktion aδ(t) durch einen Pfeil der Länge a dargestellt.

Die Delta-Funktion besitzt die Extraktionseigenschaft

∞

mit deren Hilfe viele Beziehungen aus der kontinuierlichen Welt auch für die zeitdiskreten Signale verwendet werden können. Hier ist ϕ(t) eine beliebige bei t 0 stetige Funktion.

So z.B. ergibt die Faltung eines Signals ϕ(t) mit der Delta-Funktion δ(t − t 0 ) eine Verschiebung des Signals an Stelle der Delta-Funktion (Abb. A.3):

A.2 Die Abtastung als Multiplikation mit einer periodi-

schen Delta-Folge

Man kann jetzt die Abtastung durch eine Multiplikation einer periodischen Folge von Delta-Funktionen p(t) mit dem kontinuierlichen Signal x(t) als Modell annehmen, Abb. A.4. Das Spektrum des kontinuierlichen Signals wird mit X(jω) bezeichnet und die periodische Folge von Delta-Funktionen wird durch die komplexe Fourier-Reihe ausgedrückt:

A.2 Die Abtastung als Multiplikation mit einer periodischen Delta-Folge 18

Wobei T s die Abtastperiode ist. Die Koeffizienten c n dieser Fourier-Reihe werden, wie üblich, durch folgendes Integral berechnet:

Somit ist die periodische Folge von Delta-Funktionen durch

gegeben.

Das Spektrum der periodischen Folge von Delta-Funktionen als Fourier-Transformierte dieser Folge wird jetzt:

Das Endergebnis basiert auf der Definition der Delta-Funktion laut Gl. (A.4) und zeigt, dass das Spektrum der in T s periodischen Folge von Delta-Funktionen auch eine periodische Folge von Delta-Funktionen der Periode ω s ist (Abb. A.4e).

A.2 Die Abtastung als Multiplikation mit einer periodischen Delta-Folge 19

Dieses Ergebnis war zu erwarten, weil die periodische Folge der Delta-Funktionen, als Grenzfolge von schmalen Rechteckimpulsen, ein Linienspektrum besitzen muss. Die Leistung dieser Impulse ist an den Vielfachen der Grundfrequenz f s = 1/T s konzentriert und führt im Grenzfall zu den gezeigten Delta-Funktionen als Spektrum im Frequenzbereich.

Das Produkt x(t)p(t) im Zeitbereich führt zu einer Faltung X(jω) ∗ P (jω) im Frequenzbereich und ergibt (laut Gl. (A.7)) das Spektrum des abgetasteten Signals X s (jω)

mit der Darstellung in Abb. A.4f.

Es entsteht in diesem Fall Aliasing, wenn die Bedingung f s /2 ≥ f max nicht erfüllt wird. Die Spektren X(j(ω − nω s ) schneiden sich und das aus den Abtastwerten rekonstruierte, kontinuierliche Signal ist stark verzerrt.

Abb. A.5: Verschiebung eines bandpassbegrenzten Signals im Basisband (Nyquist-Bereich)

In der Kommunikationstechnik werden oft bandpassbegrenzte Signale verwendet, für die auch eine Verletzung des Abtasttheorems in Frage kommt. Abb. A.5a zeigt das Spektrum eines realen bandbegrenzten Signals mit einer Bandbreite B um die Mittenfrequenz f 0 .

Wenn die Abtastfrequenz f s kurz unterhalb der Frequenz f 0 −B/2 gewählt wird und der Bereich 0 bis f s −B/2 genügend groß ist, entsteht das Aliased-Spektrum im Nyquist-Bereich (grau hervorgehoben), das in der Sprache der Kommunikationstechnik das im Basisband verschobene Signal bildet.

Wenn die Spektren dieser verschiedenen Komponenten sich nicht schneiden, dann kann man aus dem Basisband-Signal ohne Verluste das original Bandpasssignal rekonstruieren. Vielmals betrachtet man diesen Fall als Abtastung ohne Aliasing und wenn sich die Spektren schneiden, sagt man, dass die Abtastung zu Aliasing geführt hat. Mit dieser Unterabtastung kann man hochfrequente Signale, als zeitdiskrete Signale im Basisband verschoben, leichter bearbeiten, um sie danach eventuell wieder im Bandpassbereich zu verschieben.

Literaturverzeichnis

[1] A.V. Oppenheim, R.W. Schäfer. Zeitdiscrete Signalverarbeitung. Oldenburg-Verlag, 1992.

[2] S.D. Stearns. Digitale Verarbeitung analoger Signale. Oldenburg-Verlag, 1998. [3] Richard G. Lyons. Understanding Digital Signal Processing. Addison-Wesley Publishing Company, 1997.

[4] J. Hoffmann. MATLAB und Simulink in Signalverarbeitung und Kommunikationstechnik. Addison-Wesley, 1999.

[5] A.V. Oppenheim, A.S. Willsky. Signale und Systeme. VCH-Verlag, 1992. [6] Edward W. Kame, Bonnie S. Heck. Fundamentals of Signals and Systems, Using MAT-LAB. Prentice-Hall, 1997.

[7] E. Oran Brigham. FFT Schnelle Fourier-Transformation. Oldenburg-Verlag, 1982.

Index

  Abtasttheorem  3

  Aliased 3  9

  Aliasing 1 4 8 10  19

  Amplitudenspektrum  2

  Digital-Oszilloskop  8

  Mehrdeutigkeit  1

  Negative Frequenzen  2

  Nyquist-Bereich  3

  Scheinfrequenzen 9  22

  Simulink-Modell  4

Index der Abbildungen

  Amplitudenspektrum der Komponenten  

  mit gleichen Abtastwerten  2

  Annäherung der Definition  

  δ(t) cos(2πf t)df  16

  Chirp-Signal abgetastetes und rekonstru  

  iertes Signal  13

  Die Faltung mit der Delta-Funktion  17

  Drei Komponenten (n 1 0 1) und deren  

  gleiche Abtastwerte  3

  Erfassung und Anzeigen der Abtastwer  

  te  8

  Modell zur Simulation des Aliasing  

  Effekts  5

  Rechteckimpuls und dessen Grenzfunk  

  tion als Delta-Funktion  15

  Scheinfrequenzen für f s 50 kHz und  

  f 0 in der Umgebung von 1 MHz  

  9  

  Signal der höchsten Frequenz das mit  

  Abtastwerten darstellbar ist  3

  Signale mit gleichen Abtastwerten 6  7

  Simulink-Modell des Audio-Aliasing  

  Experiments  12

  Spektrum eines zeitdiskreten Signals  18

  Stellen der Speicher eines Rades in Bild  

  sequenzen  4

  Verschiebung eines bandpassbegrenzten  

  Signals im Basisband (Nyquist  

  Bereich)  19

Ermittlung der Frequenzgänge von digitalen

Filtern mit Hilfe der FFT in MATLAB

Inhaltsverzeichnis

1  Ermittlung der Frequenzgänge von FIR- und IIR-Filtern  1

1.1  Frequenzgang von FIR-Filtern  1

1.2  Frequenzgang von IIR-Filtern  4

2  Frequenzgänge von Filterbänken  7

2.1  DCT-Filterbänke  7

  Literatur  8

Kapitel 1

Ermittlung der Frequenzgänge

von FIR- und IIR-Filtern

Es wird gezeigt, wie man mit Hilfe der FFT die Frequenzgänge von digitalen FIR- und IIR-Filtern effizient ermitteln kann. Diese Methode ist auch in der MATLAB-Funktion freqz implementiert.

1.1 Frequenzgang von FIR-Filtern

Ein FIR-Filter realisiert folgende zeitdiskrete Eingang-Ausgangsbeziehung: y[kT s ] = b 0 x[kT s ] + b 1 x[(k − 1)T s ] + · · · + b m x[(k − m)T s ]

Wobei durch T s die Abtastperiode und durch f s = 1/T s die Abtastfrequenz der diskreten Sequenzen bezeichnet wird.

Die dazugehörige Übertragungsfunktion ist:

Den Frequenzgang erhält man, wenn z durch z = e ^jωTs in H(z) ersetzt wird [1], [2]:

Es ist bekannt [3], dass der Frequenzgang eine periodische Funktion in ω der Periode ω s = 2π/T s ist. Somit muss man den Frequenzgang nur für den Bereich 0 ≤ ω < 2π/T ^s ermitteln. Dafür wird dieser Bereich in N Intervalle unterteilt und der Frequenzgang wird nur für diese Werte ermittelt:

Der Frequenzgang H(e ^jωTs ) für ω = ω i wird:

1.1 Frequenzgang von FIR-Filtern 2

Wenn die Koeffizienten des FIR-Filters b 0 , b 1 , b 2 , . . . , b m mit Nullwerten bis zur Länge N ≥ m erweitert werden, dann geht die Funktion H i in die DFT ¹ oder FFT ² [1], [2] der erweiterten Koeffizienten über:

wobei

die erweiterten Koeffizienten des FIR-Filters sind. Die Erweiterung der Koeffizienten führt dazu, dass der Frequenzgang für N Werte der Frequenz im Bereich 0 ≤ ω < ω ^s berechnet wird.

Im folgenden MATLAB-Programm [4] (freqz_fft.m) wird gezeigt, wie man die beschriebenen Sachverhalte für die Ermittlung des Frequenzgangs eines FIR-Tiefpassfilters über die DFT oder FFT anwendet.

%--------------------------------------------------------------------------------------% Programm freqz_fft.m in dem der Frequenzgang eines % FIR-Filters über die FFT ermittelt und dargestellt wird

² Fast-Fourier-Transform

1.1 Frequenzgang von FIR-Filtern 3

% ------- Entwicklung eines FIR-Tiefpassfilters % mit der Funktion fir1 nord = 64; % Ordnung des Filters m = nord + 1; % Anzahl der Koeffizienten fr = 0.3; % Relative Bandbreite (f/fs)

b = fir1(nord, 2*fr);

% ------- Frequenzgang für N Punkte % entlang des Einheitskreises N = 1024;

be = [b, zeros(1,N-m)]; % Erweiterung mit Nullwerten

Hi = fft(be); % DFT (oder FFT)

% ------- Darstellung des Frequenzgangs figure(1); clf;

subplot(221), plot((0:N-1)/N, 20*log10(abs(Hi))); title(’Amplitudengang’);

xlabel(’Relative Frequenz f/fs’); grid;

subplot(223), plot((0:N-1)/N, unwrap(angle(Hi))); title(’Phasengang’);

xlabel(’Relative Frequenz f/fs’); grid

% ------- Darstellung des Frequenzgangs mit absoluten % Frequenzen fs = 1000; % Abtastfrequenz subplot(222), plot((0:N-1)*fs/N, 20*log10(abs(Hi))); title(’Amplitudengang’);

xlabel(’Absolute Frequenz f/fs’); grid;

subplot(224), plot((0:N-1)*fs/N, (angle(Hi))); title(’Phasengang’);

xlabel(’Absolute Frequenz f/fs’); grid Das Filter der Ordnung 64 besitzt 65 Koeffizienten, die im Zeilenvektor b hinterlegt sind. Es wurde mit der MATLAB-Funktion fir1 entwickelt. Die Koeffizienten werden mit Nullwerten bis zu N = 1024 erweitert und im Vektor be gespeichert. Somit werden N Punkte des Frequenzgangs entlang des Einheitskreises berechnet. Die Abszisse der Darstellungen in relativen Frequenzen wird durch (0:N-1)/N und für absolute Frequenzen mit (0:N-1)*fs/N erhalten. Abb. 1.1 zeigt den Frequenzgang für relative und absolute Frequenzen. Es wurde eine Abtastfrequenz f s = 1000 Hz angenommen.

Mit der Funktion unwrap in subplot(223), plot((0:N-1)/N, unwrap(angle(Hi)));

1.2 Frequenzgang von IIR-Filtern 4

wird der Phasengang ohne die typischen Sprünge von 2π dargestellt. Diese Sprünge werden eingesetzt, um die Auflösung der Darstellung der Phase zu verbessern. In der Darstellung rechts unten sieht man diese Sprünge, weil hier die Phase ohne die Funktion unwrap dargestellt wird:

subplot(224), plot((0:N-1)/N, (angle(Hi)));

1.2 Frequenzgang von IIR-Filtern

Ein IIR-Filter realisiert folgende zeitdiskrete Eingang-Ausgangsbeziehung: y[kT s ] = b 0 x[kT s ] + b 1 x[(k − 1)T s ] + · · · + b m x[(k − m)T s ] (1.8)

Die entsprechende Übertragungsfunktion ist jetzt:

Den Frequenzgang erhält man auch hier, wenn z durch z = e ^jωTs in H(z) ersetzt wird:

H(e ^jωTs ) =

mit a 0 = 1

Statt der Summe aus Gl. (1.3), mit der man den Frequenzgang eines FIR-Filters berechnen kann, ist jetzt noch eine ähnliche Summe für den Nenner in Gl. (1.10) hinzugekommen. Sie kann in der selben Art für N Punkte entlang des Einheitskreises berechnet werden.

Zusammenfassend stellt man fest, dass für IIR-Filter zwei DFTs (oder FFTs) berechnet werden müssen und zwar einmal für den Zähler und einmal für den Nenner des komplexen Frequenzgangs H(e ^jωTs ). Dazu werden die Koeffizienten des Zählers und des Nenners mit Nullwerten bis zu einer Länge N erweitert und danach die DFT dieser erweiterten Koeffizienten ermittelt. Die elementweise Teilung dieser DFTs ergibt den komplexen Frequenzgang des IIR-Filters für N Punkte entlang des Einheitskreises. Wenn die ursprüngliche Anzahl der Koeffizienten im Zähler und Nenner verschieden sind, müssen sie bis zur gleichen Länge N erweitert werden. Im Programm freqz_iir.m ist ein Beispiel gezeigt:

% Programm freqz_iir.m in dem der Frequenzgang eines % IIR-Filters über die FFT ermittelt und dargestellt wird

% ------- Entwicklung eines IIR-Tiefpassfilters % mit der Funktion ellip nord = 8; % Ordnung des Filters m = nord + 1; % Anzahl der Koeffizienten in Zähler und Nenner

1.2 Frequenzgang von IIR-Filtern 5

fr = 0.3; % Relative Bandbreite (f/fs) Rp = 0.1; % Welligkeit im Durchlassbereich Rs = 60; % Dämpfung im Sperrbereich

[b,a] = ellip(nord, Rp, Rs, 2*fr);

% ------- Frequenzgang für N Punkte % entlang des Einheitskreises N = 1024;

be = [b, zeros(1,N-m)]; % Erweiterung mit Nullwerten ae = [a, zeros(1,N-m)];

Hz = fft(be); % DFT (oder FFT) Hn = fft(ae);

Hi = Hz./Hn;

% ------- Darstellung des Frequenzgangs figure(1); clf;

subplot(221), plot((0:N-1)/N, 20*log10(abs(Hi))); title(’Amplitudengang’); xlabel(’Relative Frequenz f/fs’); grid;

subplot(223), plot((0:N-1)/N, unwrap(angle(Hi)));

1.2 Frequenzgang von IIR-Filtern 6

title(’Phasengang’);

xlabel(’Relative Frequenz f/fs’); grid

% ------- Darstellung des Frequenzgangs mit absoluten % Frequenzen fs = 1000;

subplot(222), plot((0:N-1)*fs/N, 20*log10(abs(Hi))); title(’Amplitudengang’);

xlabel(’Absolute Frequenz f/fs’); grid;

subplot(224), plot((0:N-1)*fs/N, (angle(Hi))); title(’Phasengang’);

xlabel(’Absolute Frequenz f/fs’); grid Abb. 1.2 zeigt den Frequenzgang des IIR-Tiefpassfilters mit und ohne unwrap für den Phasengang.

Die Funktion fft erzwingt automatisch eine Erweiterung mit Nullwerten bis zur Länge N wenn sie, wie folgt, aufgerufen wird: Hz = fft(b,N); % statt Hz = fft(be) Hn = fft(a,N); % statt Hn = fft(ae) Es wurde angenommen, dass die Längen der Vektoren b und a kleiner als N sind.

Kapitel 2

Frequenzgänge von Filterbänken

Wenn der Frequenzgang einer Filterbank bestehend aus mehreren Bandpassfiltern zu ermitteln ist, dann kann dies mit Hilfe der FFT in MATLAB sehr leicht realisiert werden. Angenommen in der Matrix hbank stellen die Zeilen die Impulsantworten der FIR-Filter der Filterbank dar. Dann wird mit N = 1024;

Hbank = fft(hbank,N,2);

eine komplexe Matrix Hbank ermittelt, in der die komplexen Frequenzgänge der Filter in den Zeilen enthalten sind. Die Koeffizienten der Filter werden mit Nullwerten bis N erweitert und mit 2 signalisiert man, dass die DFT (oder FFT) entlang der Zeilen zu berechnen ist. Ohne diesen Wert wird defaultmäßig die Berechnung entlang der Spalten durchgeführt.

Danach können die Frequenzgänge der Filterbank direkt dargestellt werden: subplot(211), plot((0:N-1)/N, abs(Hbank.’)); ...

subplot(212), plot((0:N-1)/N, angle(Hbank.’)); ...

Mit Hbank.’ wird wegen des Punktes die Transponierung der komplexen Matrix, ohne dass die Elemente der Matrix durch ihre konjugiert komplexen Werte zu ersetzen, realisiert. Die Transponierung ist notwendig, weil die Funktion plot stellt defaultmäßig die Spalten einer Matrix dar.

2.1 DCT-Filterbänke

Als Beispiel wird im Programm freqz_filterbank.m eine DCT ¹ -Filterbank [4] mit 8 Filtern entwickelt:

% Programm freqz_filterbank.m in dem der Frequenzgang einer % FIR-Filterbank über die FFT ermittelt und dargestellt wird

% ------- Entwicklung einer FIR-DCT-Filterbank % über die Matrix der DCT

2.1 DCT-Filterbänke 8

nfilter = 8;

hbank = fliplr(dctmtx(nfilter));

% ------- Frequenzgang für N Punkte % entlang des Einheitskreises N = 1024;

Hbank = fft(hbank,N,2);

% ------- Darstellung des Frequenzgangs figure(1); clf;

subplot(211), plot((0:N-1)/N, abs(Hbank.’)); title(’Amplitudengang’);

xlabel(’Relative Frequenz f/fs’); grid;

subplot(212), plot((0:N-1)/N, unwrap(angle(Hbank.’))); title(’Phasengang’);

xlabel(’Relative Frequenz f/fs’); grid Die Funktion dctmtx erzeugt die Matrix der DCT-Transformation. Um daraus die Filter der Bank zu bilden, werden die Zeilen mit fliplr gedreht. Weil sie bei der Faltung in der Filterungsfunktion nochmals gedreht werden, stellen sie die DCT über eine Filterung dar.

Abb. 2.1a zeigt die Amplitudengänge der DCT-Filterbank mit 8 Filtern und Abb. 2.1b stellt nur die Amplitudengänge einer DCT-Filterbank mit 32 Filtern in dB dar, ausgehend von einer 32 × 32 DCT-Matrix.

2.1 DCT-Filterbänke 9

Literaturverzeichnis

[1] A.V. Oppenheim, R.W. Schäfer: Zeitdiskrete Signalverarbeitung. Oldenburg-Verlag, 1992.

[2] S.D. Sterns: Digitale Verarbeitung analoger Signale. Oldenburg-Verlag, 1998. [3] A.V. Oppenheim, A.S. Willsky: Signale und Systeme. VCH-Verlag, 1992. [4] J. Hoffmann: MATLAB und Simulink in Signalverarbeitung und Kommunikationstechnik. Addison-Wesley, 1999.