Hi,
ich habe den Spannungsverlauf an einer Asynchronmaschine mit dSpace
aufgenommen und möchte jetzt gerne die Grundschwingung / 1. Harmonische
in MATLAB berechnen, da ich für weitere Arbeiten den
Grundschwingungseffektivwert usw. berechnen will.
%Speichern der Messwerte in Vektor
t=mw(:,1);
u1=mw(:,2);
Die beiden Vektoren habe eine Länge von 1024000, wobei genau 20 Perioden
gespeichert wurden (2 s Aufnahmezeit bei 10 Hz).
Wie kann man nun die 1. Harmonische ermitteln? Ist das richtig, dass ich
zuerst einmal den Datenvektor auf 1 Periode beschneiden sollte?
% Daten auf 1 Periode zuschneiden
u1 = u1(1:length(u1)/20);
t = t(1:length(t)/20);
Ist meine Berechnung der Grundschwingung korrekt?
% Berechnung der Grundschwingung
N = length(u1);
for n=1:N
u1gs(n) = u1(n)*exp(-1i*((2*pi)/N)*n);
end
Ich würde euch bitten mir da zu helfen.
Alex schrieb:> Ist das richtig, dass ich zuerst einmal den Datenvektor auf 1 Periode> beschneiden sollte?
Das kannst du tun, must es aber nicht. Wichtig ist nur, dass du eine
ganzzahlige Anzahl von Perioden bearbeitest. Mehr Perioden bedeutet mehr
Rechenzeit, aber auch eine bessere Unterdrückung des Messrauschens.
> u1gs(n) = u1(n)*exp(-1i*((2*pi)/N)*n);
Stimmt, allerdings ist das noch nicht die Grundschwingung selbst. Von
u1gs musst du noch den Mittelwert bilden, diesen mit 2 multiplizieren,
und schon hast du die komplexe Amplitude der Grundschwingung. Diesen
kannst du dann in Betrag und Phase zerlegen. Für die Mittelwertbildung
gibt es in Matlab sicher eine fertige Funktion, so dass du das Ganze
ohne (ineffiziente) Schleife hinschreiben kannst.
Ich habe jetzt die Zeilen so geändert (noch ineffizient, aber für's
erste für mich nicht wichtig), wie du es mir vorgeschlagen hast:
N = length(u1);
for n=1:N
u1gs(n) = 2/N*(u1(n)*exp(-1i*((2*pi)/N)*n));
end
Im Anhang siehst du die Messdaten in einem Plot, die ich auswerte. Ich
würde jetzt gerne die Grundschwingung als Zeitfunktion in den Plot
einzeichnen. Ist das mit der komplexen Zahl, die ich bei u1gs
rausbekomme, überhaupt möglich?
Wenn ich Betrag und Phase mit "abs" und "angle" ausrechne, wie komme ich
dann auf die Zeitfunktion für die Grundschwingung?
Wenn die Frequenz bekannt ist kann man gleich eine Vereinfachte fourier
mit nur einer Frequenz machen. Die Grundschwingung ist dann :
realteil = 1/T integral f(t) sin ft dt
Imaginaer = 1/T integral f(t) cos ft dt
Fuer die Oberwellen N = 2..M dann :
realteil = 1/T integral f(t) sin Nft dt
Imaginaer = 1/T integral f(t) cos Nft dt
Alex schrieb:> Wenn ich Betrag und Phase mit "abs" und "angle" ausrechne, wie komme ich> dann auf die Zeitfunktion für die Grundschwingung?
u(t) = abs(a) * cos(omega * t + angle(a))
a ist dabei die zuvor berechnete komplexe Amplitude, omega die
(bekannte) Kreisfrequenz der Grundschwingung.
Ich habe jetzt die Grundschwingung mit der Formel von dir berechnet.
u1 = u1(1:length(u1)/20);
t = t(1:length(t)/20);
N = length(u1);
for n=1:N
u1gs(n) = 2/N*(u1(n)*exp(-1i*((2*pi)/N)*n));
i1gs(n) = 2/N*(i1(n)*exp(-1i*((2*pi)/N)*n));
end
for n=1:length(abs(u1gs))
u1gs_t(n) = abs(u1gs(n))*cos(2*pi*10*t(n)+angle(u1gs(n)));
i1gs_t(n) = abs(i1gs(n))*cos(2*pi*10*t(n)+angle(i1gs(n)));
end
Interessanterweise habe ich aber nur so eine geringe Amplitude im
Millibereich erhalten. Wenn ich mir aber die FFT von u1 anschaue, so
habe ich aber einen viel, viel größeren Peak bei dieser Frequenz.
Da muss irgendwas falsch sein bei der Berechnung der Grundschwingung.
Alex schrieb:> Da muss irgendwas falsch sein bei der Berechnung der Grundschwingung.
Ja, da fehlt die Mittelwertbildung.
Weil ich Matlab weder habe noch kann, habe ich die Berechnung mal in
Python programmiert, was du sicher ganz leicht in Matlab übersetzen
kannst:
1
from math import sin, cos, pi, radians, degrees
2
from cmath import exp, phase
3
from sys import stderr
4
5
f = 10 # Frequenz
6
T = 1 / f # Periodendauer
7
omega = 2 * pi * f # Kreisfrequenz
8
n = 1024000 // 20 # Anzahl Samples in 1 Periode
9
dt = T / n # Zeit zwischen zwei aufeinanderfolgenden Samples
10
11
# Funktion für die Messreihe mit 7 Harmonischen
12
# Die Grundschwingung hat den Betrag 3.0 und die Phasenverschiebung 72°
13
14
def f(t):
15
return 3.0 * cos(1 * omega * t + radians(72)) \
16
- 2.4 * cos(2 * omega * t + radians(10)) \
17
+ 1.9 * cos(3 * omega * t - radians(20)) \
18
+ 1.7 * cos(4 * omega * t + radians(30)) \
19
- 1.5 * cos(5 * omega * t - radians(40)) \
20
- 1.3 * cos(6 * omega * t + radians(50)) \
21
+ 1.1 * cos(7 * omega * t + radians(60))
22
23
# gen_t1_u1 generiert eine Beispielsmessreihe aus der obigen Funktion f
Diese Werte stimmen genau mit den Parametern des ersten Summanden in der
Funktion f überein. Im Anhang sind u1 und die berechnete Grundschwingung
über die Zeit geplottet.
