Ich habe diskrete Werte (Zeit und dazugehöriger Wert) die aus einem Einheitssprung entstanden sind. Wie bekomme ich daraus nun die Impulsantwort/Sprungantwort oder die Übertragungsfunktion? Geht da was mit Laplace? Oder mit Hilfe von Matlab? Kenne mich mit Regelungstechnik leider nicht aus.
:
Verschoben durch Moderator
Wenn du die Werte schon hast: -Excel -Gnuplot -Matlab such dir was raus.
Kommen jetzt alle Leute hier ausm Kindergarten?! -Googeln. Soll man das für dich machen? Vorallem du sagst nicht einmal was du überhaupt für ein Programm zur Vergügung hast. Viel Spaß.
Habe zugriff auf Matlab und MS Excel. Gegoogelt habe ich bereits aber nichts treffendes gefunden deshalb frage ich hier.
https://www.mikrocontroller.net/articles/Digitalfilter_mit_ATmega#Frequenzgang_ma.C3.9Fgeschneidert speziell der den Artikel zu FDLS: http://www2.units.it/ramponi/teaching/DSP/materiale/Ch5%28x%293e.pdf#page=29 "the FDLS algorithm produces a transfer function that approximates an arbitrary frequency response" also aus Messwerten wird eine Übertragungsfunktion (ein IIR-Filter wählbarer Stufenzahl) angenähert. Nein, das war die falsche Methode da müssen ja unterschiedliche Frequenzen gemessen werden also im Frequenzbereich, dann ist es vielleicht der "prony"-Algorithmus
:
Bearbeitet durch User
derneue schrieb: > Ich habe diskrete Werte (Zeit und dazugehöriger Wert) die aus einem > Einheitssprung entstanden sind. Wie bekomme ich daraus nun die > Impulsantwort/Sprungantwort oder die Übertragungsfunktion? Sry, ich habe dich falsch verstanden. Du hättes mal die Folge aufmalen sollen. Letztendlich musst du nur den Quotienten bilden: w(z)=1z^0 + 1z^-1 + 2^-2 (Eingangsprung) x(z)=4z^0 + 3z^-1 + 2^-2 (z.B. gegebene Folge) G(z)=x/w jetzt nur noch mit den größten Exponent erweitern. Hier mit z^2/z^2. Im Simulink braucht man dies nichteinmal zu machen. Einfach einen zwei Blöcke mit Eingangsprung und Transfer Function erstellen + scope. In Matlab ist dieser Link nützlich: https://de.mathworks.com/help/symbolic/compute-z-transforms-and-inverse-z-transforms.html Da muss man erstmal die Übtragung in den Zeitbereich transformieren. mfg
aSma>> schrieb: > w(z)=1z^0 + 1z^-1 + 2^-2 (Eingangsprung) > x(z)=4z^0 + 3z^-1 + 2^-2 (z.B. gegebene Folge) > G(z)=x/w jetzt nur noch mit den größten Exponent erweitern. > Hier mit z^2/z^2. > > Im Simulink braucht man dies nichteinmal zu machen. Einfach einen zwei > Blöcke mit Eingangsprung und Transfer Function erstellen + scope. Ja abe wie kommst du auf die Gleichungen, die sehen bei mir bestimmt anders aus? Ich weiss leider immernoch nicht wie ich vorzugehen habe. Simulink ist auch adrauf bei Matlab bei mir wie könnte ich nun vorangehen. Ich kann meine Tabelle mit den Werten in eine cell variable oder struct variable laden. Und dann?
derneue schrieb: > Geht da was mit Laplace? Nennt sich im diskreten Z-Transformation. Ist aber schon fast 30Jahre her bei mir, ich müsste genauso nachlesen wie du. Also Suche nach Z-Transformation
derneue schrieb: > Ja abe wie kommst du auf die Gleichungen, die sehen bei mir bestimmt > anders aus? Mit Verlaub, ja. Stell deine Folge hier rein, dann sehen wir weiter. > Ich weiss leider immernoch nicht wie ich vorzugehen habe. Habe ich dir schon zuvor an einen Beispiel gezeigt.
Links die spalte ist die Zeit rechts die dazugehörigen Werte.
Ich nehme alles wieder zurück. Du bist doch unfähig!!! Speicher die deine datei in den current folder als: werte.dat dann tippst du das rein: >> load werte.dat >> plot(werte)
Und was soll das? Ich weiss wie die Kurve aussieht, ich weiss auch wie ich sie plotten kann. Ich brauche doch die Übertragungsfunktion aus diesen Daten.
Hilf dir das weiter ab S.8 http://staff.ltam.lu/feljc/school/asser/2_Regelstrecken.pdf Man kann eine Pt2 Strecke als Pt1-Strecke mit einer Totzeit nachbilder oder durch einen allpass glied aus Pt1 strecken. Oder willst du das digital bearbeiten können?
Ich glaube, dass es ihm darum geht, die Pol- und Nullstellen des Systems anhand einer aufgenommenen Folge zu ermitteln. Das Stichwort heißt hier: Systemidentifikation. Ich habe mir die Messungen nicht angesehen, wenn es allerdings eine zeitliche Messung ist (Sprungantwort etc.), bietet sich z.B. die Methode der kleinsten Fehlerquadrate (Least Square) an. Gruß,
Angehängte Dateien:
-
pt2.png
12 KB
Hier das ist eine PT1 Strecke. Die geringe Verzugszeit ist eigentlich vernachlässigbar.
1 | clear all; |
2 | clc; |
3 | |
4 | load data.dat |
5 | |
6 | s = tf('s'); |
7 | K = 100; |
8 | Ts = 0.22; |
9 | t=0:0.01:1.5; |
10 | |
11 | Gs = K / (1 + Ts*s); % Transferfunktion PT1-Glied |
12 | y = step(Gs,t); |
13 | |
14 | figure
|
15 | plot(data(:,1), data(:,2),'r','LineWidth',2); |
16 | hold on; |
17 | plot(t,y,'b','LineWidth',2); |
18 | grid on; |
Du erhälst mit der Division der Fouriertransformierten des Ausgangssignals durch die Fouriertransformierte des Eingangssignals die Fouriertransformierte der Impulsantwort. Das hat bei mir, als ich das damals gebraucht habe, ganz gut funktioniert. Eine weitere Möglichkeit wäre die Verwendung eines adaptiven Filters mit LMS, NLMS oder RLS Algorithmus. Das alles geht aber nur, wenn dein System linear ist!
Wenn man diskrete Werte hat und eine Funktion davon moechte, muss man ein rationales Polynom fitten. Das Einfachste, aber Aufwendigste ist N Gleichungen mit N unbekannten. Es kann sein, dass das System einen zu hohen Grad aufweist und die Koeffizienten nicht stabil sind. Das System ist dann schlecht konditioniert. Dann macht man ein einfacheres Modell, das nun aber ueberbestimmt ist. N Gleichungen mit M Unbekannten, mit N>M. Dann muss man einen Least Square Fit machen. Ich wuerde das zweite verwenden, auch wenn es mathematisch aufwendiger ist. Noch einfacher ist wenn man schon eine bessere Vorstellung hat, dann kann man die Funktion graphisch fitten. Eine Programmieruebung.
Ich staune ueber den mathematischen Level hier. :-) Kuckt mal hier: https://de.wikipedia.org/wiki/Polynominterpolation Habt ihr das nicht in der Schule gelernt? Für die Faulen unter euch wurde es unter dem namen polyfit in Octave eingebaut. Olaf
zuerst ist das Signal zu differenzieren, da es von einem Sprung und njchf von ejnem Dirac kommt. Danach Jnverse fft. Die ANZAHL der Punkte DER inv. fft soll viel größer sein als die Messwerte, (4096 zu 100) da sonst die Messwerte zyklisch ausgewertet werden.
Hallöchen, seid ihr teilweise kirre? Sorry Leute, aber man differenziert keine Sprungantwort, um daraus ne Impulsantwort zu bekommen, viel zu rauschempfindlich (Differentiation ist ein super Hochpass). Ebenso soll hier kein Polynom gefittet werden, das ist ein dynamisches System, nix statisches. Das Zauberwort heißt Systemidentifikation, wurde oben (zwischen den absurden Antworten) schon genannt. Dazu passt Kapitel 8 "Methode der kleinsten Quadrate für dynamische Prozesse" aus dem Isermann, "Identifikation dynamischer Systeme 1". Davon gibts noch diverse Modifikationen (rekursiv, gewichtet, totale ls usw usw). Du möchtest also eine Übertragungsfunktion eines zeitdiskreten Modells haben?
mit dem Vektor der Eingangssignale u(z), dem Vektor der gemessenen Ausgangssignale y(z), der Ordnung m und N der Anzahl der Samples (im Buch kommt noch eine Totzeit dazu, lasse ich der Einfachheit halber jetzt weg) Sind einige Voraussetzungen erfüllt: - Ordnung m und Totzeit d bekannt oder zumindest grob bekannt (immer schön niedrig anfangen...) - Eingangssignal muss Messbar sein (ist erfüllt) - Störsignal ist mittelwertfrei und stationär (z.B. weißes Rauschen durch die Messungen). Kp ob das erfüllt wird bei dir, geht aber häufig auch ohne ;-) Dann kannst du deine Parameter ai und bi schätzen, wenn du
rechnest. Mit
,
, du stapelst also einfach alle deine Messungen übereinander. Dann brauchst du noch die Regressormatrix
, stapelst also auch deine Eingangs und Messgrößen in eine Matrix. Dann berechnest du wie oben angegeben die Pseudoinverse und hast deine Parameter ai und bi. Weitere Stichworte sind Parameteridentifikation, output error model, least squares, armax, system identification. Das ist alles in der System Identification Toolbox in Matlab eingebaut (und noch viel viel viel mehr). Schöne Grüße, Jan
:
Bearbeitet durch User
Servus, @ Jan K. (jan_k) ja ich habe gerade dein Ansatz mal soeben ausprobiert, damit kriege ich folgende Fkt raus: G(z)= 4.6z/(z-0.95). Diese passt mit der Übertragungsfkt. von werte.txt perfekt überein. Danke.
Hallöchen Jan, anzumerken bei deinem Ansatz wäre noch, dass dies nur für nicht-rekursive Systeme funktioniert. Der Grund ist der, dass im Vektor y das gemessene Ausgangssignal ist, die Pseudoinverse minimiert also den Fehler, wenn an deinem System genau dieses Ausgangssignal raus kommt. Leider ist das Ausgangssignal aber bei einem rekursiven System eine Funktion vom Eingangssignal und von sich selbst (zeitverschoben). Ein IIR-System kannst du so nicht identifizieren. Jan K. schrieb: > Hallöchen, seid ihr teilweise kirre? > > > Das Zauberwort heißt Systemidentifikation
Hallo, auch von mir ein Danke an Jan! Trotzdem verstehe ich die Matrix nicht. Das Element unten links -y[m+N-1] ergibt für z.B. m=2 -y[N+1]. Das höchste Element in y hat aber den Index N.
Angehängte Dateien:
-
zplane.png
5,8 KB -
step.png
18 KB -
fit.png
7,8 KB
Uhoh, da sind leider 1, 2 Fehler drin sorry. u(m) in der vorletzten Spalte ist nicht richtig, da muss natürlich u(m-1) hin. Auch y ist falsch wie du richtig erkannt hast, y muss bei y(m) starten, nicht bei y(0). Die Gesamtanzahl deiner Samples ist N+m. Essenziell für dein Problem ist übrigens, dass die Werte vorher mit einer konstanten Abtastrate resamplet werden (dein diff(t) ist nicht konstant!), sonst funktioniert das ganze nicht. Ich bekomme übrigens
für eine first order approximation heraus. Scheint nicht zu reichen, siehe Bild #3 ;-) @Tobias Doch, man kann auch rekursive Systeme mit obigem Ansatz identifizieren:
Das ist ein reines IIR System 3. Position der Pole und Nullstellen kann man im ersten Bild sehen.
1 | [y, tt] = stepz(1,[1 -.95 .001 .1],100); % daten generieren |
2 | u = ones(size(y)); % eingangssignal |
3 | m = 1; % approximation 1. Ordnung |
4 | psi = [cell2mat(arrayfun(@(i)-y(m-i+1:end-i),1:m,'unif',false)) cell2mat(arrayfun(@(i)u(m-i+1:end-i),1:m,'unif',false))]; % regressormatrix |
5 | p = psi\y(m+1:end); % parameter |
6 | |
7 | figure, |
8 | stem(y) |
9 | hold all |
10 | stepz([0 p(floor(length(p)/2)+1:end)'],[1 p(1:floor(length(p)/2))']) |
11 | legend('original','approximated')
|
First order approximation resultiert also in einem Modell
Beachte, dass im Nenner der erste Koeffizient auf 1 normiert ist (wird nicht mit identifiziert) und dieser im Zähler == 0 ist. Du kannst auch die Modellordnung höher wählen, nur dann ist das Gleichungssystem nicht mehr gut konditioniert. Funktionieren tut es trotzdem, probier es mal mit m=2 oder m=3. @asma wie bist du auf die positiven Exponenten gekommen? Wenn ich das umrechne bleibt nur ein konstanter gain im Zähler übrig, das kann diese Methode eigentlich nicht berechnen.
:
Bearbeitet durch User
@Jan man muss vielleicht präzisieren: man kann ein IIR-System so identifizieren und bekommt dann eine Übertragungsfunktion, welche den Gleichungsfehler (so hiess das glaube ich) minimiert. In der Regel möchte man aber den Ausgangsfehler minimiert haben, und das geht so nicht! Grund: dein Vektor y enthält die gemessenen Werte des Ausgangssignals. Wenn das ein IIR-System sein soll, dann sind das aber auch gleichzeitig Funktionen der Koeffizienten a. Da du die a aber nicht kennst, .... du siehst das Problem.
Jan K. schrieb: > @asma > wie bist du auf die positiven Exponenten gekommen? Wenn ich das umrechne > bleibt nur ein konstanter gain im Zähler übrig, das kann diese Methode > eigentlich nicht berechnen. Servus, ich schreibe nie etwas im Forum ab. Ist einfach zu Fehleranfällig. Ich habe einen Ansatz genommen bei dem man den Exponenten sowohl im Zähler als auch im Nenner vorgeben kann. m=0 (Zähler), n=2 (Nenner). Dies hat auch "Oh Doch (jetztnicht)" vorgeschlagen. Mit Verlaub, habe ich meiner Formel G(z)= 4.6z/(z-0.95) die Abtastzeit vergessen anzugeben: Ta=10ms.
Ich bin mir nicht sicher, ob "Oh Doch (jetztnicht)" das meinte, oder ob er einfach nur ein statisches Polynom fitten wollte. Das Modell ist zwar polynomial, aber eben im Kontext einer Übertragungsfunktion. Woher hast du deine Quelle? :-)
Jan K. schrieb: > y muss bei > y(m) starten, nicht bei y(0). Hallo Jan, danke nochmal für die Antwort. In dem Fall passt aber das Element oben links nicht mehr. Wenn y bei m Anfängt, was ist dann y(m-1)?
Da ist die Namensgebung sehr blöd von mir und beruht auf Faulheit :-/ Wenn y alle deine m+N gemessenen Werte sind, so sei Y ein Vektor mit allen Messungen, bei y(m-1) beginnend. Y ist dann der Vektor auf der rechten Seite des LGS oben. Versuche das Matlab Skript oben nachzuvollziehen. Bei p = psi\y(m+1:end); wird nicht der gesamte Vektor y, sondern eben ab dem m+1-ten Element genommen. Wichtig: Die Schreibweise in der Matrix oben geht von der Indizierung von 0 aus, y(0) ist also dein erstes gemessene Ausgangselement zum Zeitpunkt k=0. In Matlab muss alles um +1 verschoben werden, da das erste Element am Index 1 ist.
Jan K. schrieb: > Da ist die Namensgebung sehr blöd von mir und beruht auf Faulheit :-/ Deswegen macht man einfach ein Foto. Dieses rumgetippe. Anscheinend hasst du nichts bessres zu tun. Such dir mal ne Freundin. Jan K. schrieb: > Woher hast du deine Quelle? :-) Quid pro quo (kirre ausgesprochen). Gibt erstmal eine fehlerfreie Version von dir ab. Dann will ich deine mathematik Künste mal testen: >>p = psi\y. Wie berechnet man dies (ohne Matlab)?
Wat bist du denn für einer? Danke, dass du dich um mein Privatwohl kümmerst. Nein, man macht natürlich kein Foto, schonmal was von Copyright gehört? Außerdem ist es in dem Buch noch etwas komplizierter dargestellt. Die Matlab Version funktioniert und sollte mit der Matrix oben übereinstimmen. Wie das LGS gelöst wird steht oben schon. Matlab bildet beim "\" bzw "mldivide" automatisch die Pseudoinverse [bzw berechnet die implizit über eine QR Zerlegung], aber man kann das natürlich auch zu Fuß ausrechnen. Wie gesagt, wie die Pseudoinverse definiert ist steht in meinem ersten Post. Wie man Matrizen transponiert und invertiert wissen hoffentlich alle, wenn nicht ist der Gauß-Jordan Algo eine erste Anlaufstelle. Bin mir bei dir nicht sicher, ob du trollst oder einfach irgendwo was gegoogelt hast.
:
Bearbeitet durch User
Servus, hier das gewünschte Skript ohne Kommentar: https://www.tu-ilmenau.de/fileadmin/media/systemanalyse/Lehre/SI_PA2/Skript-Systemidentifikation-V19.pdf
Hallo Jan, vielen Dank für die Hilfe! Eine Frage hätte ich noch: Wird in dem Buch / in den Büchern auch die Systemidentifikation für Mehrgrößensysteme beschrieben?
Thread beobachten
Seitenaufteilung abschaltenBitte melde dich an um einen Beitrag zu schreiben. Anmeldung ist kostenlos und dauert nur eine Minute.
Bestehender Account
Schon ein Account bei Google/GoogleMail? Keine Anmeldung erforderlich!
Mit Google-Account einloggen
Mit Google-Account einloggen
Noch kein Account? Hier anmelden.