NFFT = 2^nextpow2(length(y)); % Next power of 2 from length of y
fft_y = fft(y,NFFT)/length(y); %Fouriertransformation Ausgang
fft_u = fft(u,NFFT)/length(u); %Fouriertransformation Eingang
die Funktion wurde an bekannten Systemen getestet und hat auch für ein VZ1-Glied funtioniert wenn ich ChirpSignal drauf gegeben habe. für ein Integrator aber schon nicht mehr.
Beim Eingang von Sprung funktioniert es aber auch beim VZ1 nicht mehr.
zweiter Ansatz:
gegebene Sprungantwort fouriertransformieren, im Bildbereichableiten...Ergebnis müsste ja eigentlich die Übertragungsfunktion sein!?
Code:
fs = 1/(t(2)-t(1)); %Abtastfrequenz
NFFT = 2^nextpow2(length(y)); % Next power of 2 from length of y
fft_y = fft(y,NFFT)/length(y); % Fourertransformieren
diese Code bringt keine guten Ergebnisse und ich weiß nicht warum.
Bin für jede hilfe dankbar. egal ob Korrektion an meinem Code, Erklärungen warum mein Vorgehen nicht funktionieren kann oder komplett neue Lösungsforschläge (auch für Simulink)
ich bin mir nicht sicher aber deine ableitung im bildbereich sieht merkwürdig aus. die multiplikation mit (jw) bedeutet eine differentiation im zeitbereich. du möchtest aber den bildbereich ableiten. hier müsste man wahrscheinlich mit der funktion diff arbeiten, und dann durch den frequenzabstand teilen.
Gast
Beiträge: ---
Anmeldedatum: ---
Wohnort: ---
Version: ---
Verfasst am: 19.01.2012, 13:07
Titel:
hey, danke für deine Antwort...Sorry, habe ich schlecht formuliert. Ich möchte im Bildbereich die Operation durchführen, die im Zeitbereich eine Ableitung nach t darstellt, weil die Sprungantwort im Zeitbreich abgeleitet meine Stoßantwort ist. und die wiederum ist im Bildbereich meine Übertragungsfunktion. überführe ich die Sprungantwort aber vorher in den Bildbereich müsste ich sie doch dort mit (jw) Multiplizieren um auf meine Übertragungsfunktion zu kommen.
und warum führst du die ableitung nicht im zeitbereich durch und transformierst dann?
Gast
Beiträge: ---
Anmeldedatum: ---
Wohnort: ---
Version: ---
Verfasst am: 19.01.2012, 13:34
Titel:
hm... , danke für den Tipp, ich werde es mal ausprobieren.
Gast
Beiträge: ---
Anmeldedatum: ---
Wohnort: ---
Version: ---
Verfasst am: 20.01.2012, 11:46
Titel:
ich hab es dann jetzt mal ausprobiert die Sprungantwort vorher abzuleiten.
Das scheint auch ganz gut zu funktionieren, nur stimmt die Amplitude irgendwie nicht. Die amplitudenverstärkungskurve ist nach unten verschoben.
hat jemand eine Idee woran das liegen könnte? muss ich irgendwo eine Amplitudenkorrektur einbauen?
auch wenn sich nicht viel geändert hat füge ich den Code nochmal an:
Code:
fs = 1/(t(2)-t(1)); %Abtastfrequenz
NFFT = 2^(nextpow2(length(g))); % Next power of 2 from length of y
fft_g = fft(g,NFFT)/length(g);
% g ist jetzt die zeitliche Ableitung meiner Sprungantwort
das könnte an der weise liegen wie du ableitest....
wie löst du denn die ableitung? beziehst du auch die Abtastzeit mit ein
Gast
Beiträge: ---
Anmeldedatum: ---
Wohnort: ---
Version: ---
Verfasst am: 23.01.2012, 11:20
Titel:
danke, wiedermal für den Vorschlag. ich hab es auf drei verschiedene Arten abgeleitet:
1) mit simulink-block: du/dt
2) dy3 = diff(y(:,2))./diff(y(:,1));
3) mit einem selbst geschriebenen Abeleitschätzer
(Abtastzeit wird überall berücksichtigt)
leider hat es nichts gebracht
die diagramme waren alle ähnlich:
phase ok, nur die Amplitudeverstärkung ist um ca. 20db nach unten verschoben!
des wegen dachte ich auch, dass ich an der Amplitude was korrigieren muss.
zum beispiel verstehe ich nämlich auch im Matlab-help-beispiel für die fft nicht warum sie hier mit 2 multipliziert wird!
Code:
Fs = 1000; % Sampling frequency
T = 1/Fs; % Sample time
L = 1000; % Length of signal
t = (0:L-1)*T; % Time vector % Sum of a 50 Hz sinusoid and a 120 Hz sinusoid
x = 0.7*sin(2*pi*50*t) + sin(2*pi*120*t);
y = x + 2*randn(size(t)); % Sinusoids plus noise plot(Fs*t(1:50),y(1:50)) title('Signal Corrupted with Zero-MeanRandom Noise') xlabel('time (milliseconds)')
NFFT = 2^nextpow2(L); % Next power of 2 from length of y
Y = fft(y,NFFT)/L;
f = Fs/2*linspace(0,1,NFFT/2+1);
Hab nochmal weiter gesucht und bin auf die chirp z-transform (czt) gestoßen.
kennt sich jemand damit aus und kann mir sagen ob ich damit ggf auf bessere Ergebnisse komme ?
Das ist doch annähernd richtig was hier raus kommt (letzter Code mit dem verrauschten Testsignal). Lediglich der Gleichsignalanteil sowie bei der Nyquistfreq. wird nicht mit L/2 skaliert sondern mit L. Aber das Spektrum ist doch korrekt...die Abweichung der Amplituden 0.6 statt 0.7 bei 50Hz und 0.93 statt 1 bei 120 Hz (ich beziehe mich auf die lineare Darstellung - nicht dB skaliert) kommen durch den Leakage-effekt. Die Auflösung des Spektrums df ist kein vielfaches der beiden Sinusfreq.
df = FS/NFFT = 1000/1024
Somit liegen die Frequenzlinien auch nicht exakt bei 50 und 120 Hz und die Amplitude ist etwas kleiner
Edit: Wegen dem multiplizieren mit 2...die fft Funktion berechnet ein zweiseitiges Spektrum. Die Energie, die in den Spektrallinien dargestellt wird, wird geteilt...für den positiven und negativen Freq.-bereich. Du willst ja nur den Positiven darstellen, weshalb du nun die Amplitudenwerte mit 2 mult. musst, damit die richtige Energie rauskommt. Der Gleichanteil und bei der Nyquistfreq. sind in dem Ergebnisvektor Y=fft(...) aber nur ein mal vorhanden, weshalb diese nicht mit 2 multipliziert werden dürfen.
Zum Thema FFT und richtige Darstellung empfehle ich dir das umfassende FFT Bsp. hier aus der Skriptecke.
Wenn dein System einen relativ einfachen Aufbau hat (nicht sehr hohe Ordnung) könntest du auch mit einem Parameterschätzverfahren die Übertragungsfuntkion ermitteln. Hier gibt es z.B. ARX oder ARMAX Modelle, womit Übertragunsglieder einfacher Art (z.B. PTn) relativ gut geschätzt werden können. Vorraussetzung ist jedoch, dass du sowohl das Eingangs- als auch das Ausgangssignal kennst/aufgezeichnet hast.
Im Zeitbreich wird ja das Eingangssignal mit dem Systemverhalten gefaltet und ergibt das Ausgangssignal. Im Frequenzbereich wird aus der Faltung eine Multiplikation. Ich bin mir momentan aber nicht sicher, ob die Rechenvorschrift in Matlab so richtig umgesetzt wird.
Schau dir dazu mal die Formel unter Punkt 9.2 an...
genau nach dieser Formel den Real- und Imaginärteil deines Systems berechnet. Ich würde dir daher empfehlen, die Ergebnisvektoren der FFT mal in RT und ImT zu trennen und nach der Formel 9.2 berechnen zu lassen. Da müsste ja das gleiche herauskommen...
Gast
Beiträge: ---
Anmeldedatum: ---
Wohnort: ---
Version: ---
Verfasst am: 08.02.2012, 11:33
Titel:
Hey Danke für deine Antworten.... konnte leider nicht eher schreiben.
habe das:
Code:
fft_g = fft_y./fft_u;
% äquivalent zu
H(f) = Y(F)./X(F)% siehe Sktipt wobei H die Frequenzantowrt des Systems ist, im Zeitbereich = Impulsantwort
auffällig hier bei ist, dass
length(g) im diskreten ja quasi meinem Zeitfenster im Kontinuierlichen entspricht.
nur wirklich verstehen tue ich es nicht.
hat jemand eine Ahnung warum ich hier diese Korrektur durchführen muss?
Du kannst Beiträge in dieses Forum schreiben. Du kannst auf Beiträge in diesem Forum antworten. Du kannst deine Beiträge in diesem Forum nicht bearbeiten. Du kannst deine Beiträge in diesem Forum nicht löschen. Du kannst an Umfragen in diesem Forum nicht mitmachen. Du kannst Dateien in diesem Forum posten Du kannst Dateien in diesem Forum herunterladen
MATLAB, Simulink, Stateflow, Handle Graphics, Real-Time Workshop, SimBiology, SimHydraulics, SimEvents, and xPC TargetBox are registered trademarks and The MathWorks, the L-shaped membrane logo, and Embedded MATLAB are trademarks of The MathWorks, Inc.
, danke für den Tipp, ich werde es mal ausprobieren.
RSS