Fourier Reihe einer harmonischen Schwingung

Übung 5
Analyse von Zeitreihen in der Umweltphysik und Geophysik
1
Fourier Reihe einer harmonischen Schwingung
In dieser Übung geht es darum, ein periodisches Signal mittels Fouriertransformation in harmonische Schwingungen unterschiedlicher Frequenzen zu zerlegen und die Amplitude und Phase
der einzelnen Komponenten zu bestimmen.
15
10
5
0
−5
−10
0
0.5
1
1.5
Zeit (s)
2
2.5
3
Abbildung 1: Periodisches Signal bestehend aus harmonischen Schwingungen.
Das abgebildete Signal besteht aus einem Mittelwert, einer Grundschwingung mit der Periode
T = 1 s und mehren Oberschwingungen. Die Datei harmonics.txt enthält den Zeitvektor t
in der ersten Spalte und das Signal x in der zweiten Spalte. Laden Sie zunächst das Skript
uebung5.m, welches den Zeitvektor t und das Signal x initialisiert:
>> uebung5
>> whos
Name
Size
Bytes
Class
dt
1x1
8
double
nt
1x1
8
double
t
1000x1
8000
double
x
1000x1
8000
double
Attributes
>> plot(t, x)
Der Zeitvektor t hat ein Abtastinverval ∆t von 0.001 s. Wie Sie im Zusammenhang mit der
diskreten Fourier Transformation noch erfahren werden, hängt die höchste Frequenz, die man
auflösen kann, direkt von ∆t ab. Da die Abtastrate 1/∆t in diesem Fall viel höher ist als die
kleinste Frequenz von 10 Hz in der harmonischen Schwingung, betrachten wir das Signal für
diese Übung als quasi kontinuierlich.
Übung 5
Analyse von Zeitreihen in der Umweltphysik und Geophysik
2
Bestimmung der Fourier Koeffizienten
Die Fourier Reihe eines Signals x(t) ist definiert als
∞
x(t) =
a0 X
+
an cos(nωt) + bn sin(nωt)
2
(1)
n=1
mit ω = 2π/T und
Z
2 T
x(t) cos(nωt)dt
an =
T 0
Z
2 T
bn =
x(t) sin(nωt)dt
T 0
(2)
Da x(t) in diesem Fall in der diskreten Form xn mit N Werten vorliegt, werden die Integrale in
2 zu Summen:
an =
bn =
N
2X
xn cos(nωtn )∆t
T
2
T
1
N
X
(3)
xn sin(nωtn )∆t
1
Berechnen Sie die Fourier Koeffizienten an und bn für n = 0, 1, 2, 3 . . . , 10 aus dem periodischen
Signal in x. Ein mögliches Vorgehen ist (Befehle für omega und b(n) ergänzen):
nfreq=10
a=zeros(nfreq, 1);
b=zeros(nfreq, 1);
T=max(t);
% Grundperiode
omega = ?
for n=1:nfreq
a(n) = 2./T * sum(x .* cos(omega*n*t)) * dt;
b(n) = ?
end
Die Amplitude cn und Phase φn für jede Frequenz n lassen sich aus an und bn berechnen:
p
bn
−1
2
2
cn = an + bn und φn = tan
an
Ergänzen Sie die obenstehende Schlaufe, um jeweils auch die Amplitude c(n) und die phase
ph(n) zu berechnen:
c=zeros(nfreq, 1)
ph=zeros(nfreq, 1)
...
c(n) = sqrt( ? );
ph(n) = atan2(b(n), a(n));
Übung 5
Analyse von Zeitreihen in der Umweltphysik und Geophysik
3
Aus den Werten in den Vektoren c und ph können Sie die Amplitude- und Phase der einzelnen
Oberschwingungen im Signal ablesen, oder Sie können die Amplituden mit stem(c) plotten.
Zur Illustration der Fourier Transformation können Sie die einzelnen Schwingungen darstellen,
aus denen das Signal besteht. Plotten Sie dazu in einer neuen Figur jeweils die einzelnen Komponenten der harmonischen Reihe in der rechten Spalte und die Summe der Komponenten in
der linken Spalte:
figure(2)
s2=zeros(nt,1);
offset=5;
for n=1:nfreq
subplot(121)
nsig = c(n)*cos(omega*n*t - ph(n));
plot(t, nsig + offset*n)
hold on
subplot(122)
s2=s2+nsig;
plot(t, s2 + offset*n)
hold on
end
Exponentielle Form
In exponentieller Schreibweise lautet die Fourier Transformation:
x(t) =
∞
X
Xn einωt
(4)
x(t)e−inωt
(5)
n=−∞
mit
1
Xn =
T
Z
T
0
Berechnen Sie die Werte Xn für n = 0, ±1, ±2, . . . , ±10 aus dem Signal in Abbildung 1. Ein
mögliches Vorgehen ist (benutzen Sie die Funktion exp zum ergänzen bei Xn):
Xn=zeros(nfreq*2+1, 1);
for n=-nfreq:nfreq
f2(n+nfreq+1)=n;
Xn(n+nfreq+1)= ?
end
Was stellen Sie beim Betrachten von Xn fest (Vergleich von negativen und positiven Frequenzen)?
Aus dem komplexen Vektor Xn können Sie die Amplitude mit abs und die Phase mit phase
bestimmen:
Übung 5
Analyse von Zeitreihen in der Umweltphysik und Geophysik
4
c2=abs(Xn)
ph2=phase(Xn)
stem(f2, c2)
Vergleichen Sie die Amplituden in c2 mit denen in c, und die Phasen in ph mit den Phasen in
ph2 (da diese im Bogenmass angegeben sind, teilen Sie am besten durch π).

Download Report