für die Auswertung von Schallmessungen bin ich auf der Suche nach Auswertefunktionen und habe ein paar Fragen zur Mittelwertbildung.
1)
In der Akustik werden sehr viele Größen als Pegel definiert. Gibt es hierzu in Matlab entsprechende Toolboxen oder muss das Rad wirklich neu erfunden werden?
2)
Für die Darstellung eines Schalldruckpegels müssen digitalisierte Rohdaten eines Schalldrucksignals zeit- und frequenzbewertet werden. Gibt es hierzu fertige Funktionen in Matlab?
3)
3a)
Wie lässt sich in Matlab für die Schalldruckpegelbildung eine Zeitbewertung - also eine gleitende Mittelwertbildung (Effektivwert) mittels "Quadratischen Mittelwerts" - für große Datenmengen effektiv realisieren? (http://de.wikipedia.org/wiki/Zeitbewertung)
3b)
Welche Art von Gleiten wird hierbei benötigt bzw. angewendet (vorauseilend, mittig, oder nacheilend)?
die benutzen dort die standardfunktionen conv und filter (wobei filter schneller arbeitet)
quadratischer mittelwert ist glaube ich hypot...
die frequenzbewertungen musst du auch selber berechnen. da ich mit matlab nur terzen verarbeite und ggf in oktaven umwandle, konnte ich aus entsprechender literatur die pegelanpassungen übernehmen.
vielen Dank für die Hinweise und den Link. Ich habe versucht das dort Gesagte zu testen und umzusetzen. Hierbei haben sich folgende Fragen und Erkenntnisse ergeben:
0) Für die Mittelwertbildung stehen verschiedene Methoden bzw. Funktionen zur Verfügung. Diese werden in verschiedenen Stellen des Forums vorgeschlagen:
Leider wird in oben angegebenen Diskussionen selten detailliert, welche Art von Mittelwert verwendet werden soll. Daher mein Versuch, für die weiteren Fragen erst einmal eine Diskussionsbasis zu schaffen:
Definition der verschiedenen Mittelwerte: http://de.wikipedia.org/wiki/Mittelwert. Die vielfältigen weiterführenden Links und Begriffe wie Effektivwert, RMS etc. zeigen, dass diese mathematischen Grunddefinitionen in einzelnen Disziplinen unter anderen Namen auftauchen (versteckt werden).
Nach dem Versuch der Zusammenführung nun zu meinen Fragen:
1) Ich habe es geschafft, mittels for-Schleife und
einen (zentriert)-gleitendes quadratisches Mittel zu definieren. Hierfür habe ich genaugenommen die for-Schleife in drei Schleifen unterteilt, damit das Indexing am Anfang und Ende auch funktioniert.
Code:
function Vektor_GMW = fkt_GleitenderMittelwert(Vektor,N_vor,N_nach)
for k = 1 : N_vor
Vektor_GMW(k) = fkt_QuadratischerMittelwert...
(Vektor(1:k+N_nach));
end
for k = N_vor + 1 : length(Vektor) - N_nach
Vektor_GMW(k) = fkt_QuadratischerMittelwert...
(Vektor(k-N_vor : k+N_nach));
end
for k = length(Vektor) - N_nach + 1 : length(Vektor)
Vektor_GMW(k) = fkt_QuadratischerMittelwert...
(Vektor(k-N_vor : length(Vektor)));
end
function QuadratischerMittelwert = fkt_QuadratischerMittelwert(Vektor)
QuadratischerMittelwert = norm(Vektor)/sqrt(length(Vektor));
Es bleibt die Frage, ob mein Umgang mit den Anfangs- und Endwerten richtig ist. Oder anders gesagt, wie ist das ganze mathematisch zu definieren? Meiner Ansicht nach sind die nichtvorhanden Werte NICHT durch Null zu ersetzen, sondern sie fallen weg und die Anzahl N (Zählerterm unter der Wurzel) wird auch auf die vorhandene Anzahl an Elementen reduziert. Kann mir hier jemand zustimmen?
2) Allerdings ist die selbstgebaute Funktion für die Anwenung in der Akustik mit einer sehr hohen Anzahl an Werten (da hohe Samplingrate, z.B. f_s=50000 1/s) und dazu noch einer hohen Anzahl an Mittelungswerten (Mittelungszeitraum für die Bewertung 'Fast' ist t_int=1/8s also z.B. 6250 Werte!) enorm langsam. Für einen schnellen Ansatz wird in diesem Forum die Funktion
ausgewählt werden können. Für ein alpha = 0 (, und nur für dieses!) entspricht das gefilterte Ergebnis dem der Funktion filter, bzw. conv für Polynome (Fenster) von 1/N*[1,1,...].
Ist es richtig, dass für das quadratische Mittel alpha = 2 ist? (Definition in der Hilfe: "2 = square weighted moving average").
4) Was BEWIRKT das Argument "lag" in der Funktion movavg? Es soll sich ja nach meinem Verständnis der Hilfe um den nachlauffenden Anteil des Fensters handeln. Führt man die Mittelwertbildung für ein beliebiges Signal durch, so führt die Veränderung des Arguments "lag" allerdings zu keinerlei Veränderung, wie man an folgendem Beispiel überprüfen kann:
5) Weiterhin führen alle oben genannten Funktionen, die sich (alle) der Faltung bedienen (?), zu anderen ersten Werten als meine selbstgebastelte langsame Funktion oder eine händische Berechnung des quadratischen Mittels. Dies habe ich am Beispiel eines Kamm-Vektors festgestellt:
Konkret führt eine ungerade Anzahl an Elementen zu einer Verschiebung des Ergebnisses der Funktion movavg. Ebenso weichen vorallem die Anfangswerte voneinander ab. Woran liegt das, bzw. wie lässt sich dies Abweichung der Anfangswerte und die Verschiebung für movavg, conf und filter in den Griff bekommen?
6) Nein, ich habe keine weiteren Fragen mehr, aber
... Fragen über Fragen... und sicherlich hat es der Ein- oder Andere nicht bis hierhin geschafft...
Trotzdem freue ich mich über alle Antworten, die sicherlich nicht nur mir weiterhelfen werden.
da ich leider keine Antworten erhalten habe, inzwischen aber ein wenig schlauer geworden bin, gebe ich an dieser Stelle mal selber die Antwort. Gerne könnt ihr mich korrigieren, falls ich irgendwo noch falsch denke.
Leider wird die Mittelwertbildung je nach Anwendungsgebiet (Regelungstechnik, Messtechnik, Signaltechnik, Stochastik, Finanzmathematik etc.) völlig unterschiedlich bezeichnet, daher musste ich auf der Suche nach der Antwort verschiedene Wege gehen.
Folgende von mir gestellte Frage soll nun nachfolgend beantwortet werden:
1) Wie wird ein gleitendes quadratische Mittel gebildet
2) Was hat das mit der function movavg zu tun (Eingabeargumente alpha)
lässt sich entnehmen, dass diese Funktion dem Bereich der Finanzmathematik zugeordnet ist. Daher sind Literaturstellen auch vorangig in diesem Bereich zu suchen. Stichworte sind "MA (Moving Average)" sowie ähnliche Abkürzungen, "Gleitendes Mittel" etc. Ein verständlicher Überblick und das "Warum + Wofür" habe ich sehr gut in dem Buch "Chartanalyse für Dummies" gefunden (googlebooks).
Dieser und anderen Quellen lässt sich entnehmen, das movavg für meinen Fall nicht zielführend ist, da es um gewichtete (nicht um sonst steht da der Zusatz weigthed) Mittelwerte geht. Die Werte im gleitenden Fenster sollen in meinem Fall aber nicht mit einem Vorfaktor (z.b. 2 = quadratisch) gewichtet werden, sondern quadratisch eingehen. Somit müssen die Werte einfach vorher einzeln quadriert werden...
Die Lösung lautet:
Code:
function Signal = GleitenderQuadratischerMittelwert...
(Signal, Fensterlaenge)
window = rectwin(Fensterlaenge);
power = Signal.^2;
rms = conv(power,window,'same'); % Zentriertes Fenster mit Fehler zu Beginn und Ende
wie angekündigt sind leider bei mir ein paar Fragen übrig geblieben. Konkret stellt sich für mich die Frage warum in movavg ein Faktor "2" verwendet wird?
Die oben schon teilweise aufgezeigten unterschiedlichen Ansätze habe ich weiterverfolgt, um einige DIN-Auswerteverfahren zu realisieren. Allerdings bin ich leider über besagten Faktor gestolpert.
Zum Hintergrund:
Der gleitende quadratische Mittelwert (beidseitig gleitend) ist bei hohen Abtastraten und großer Fensterbreite (z.B. Fs = 50kHz und tau,slow = 1s ) schon ganz schön rechenintensiv. Daher wird in der Fachliteratur (siehe: GoogleBooks, Stichwort "Gleitender Effektivwert" z.B. in "Gerhard Müller, Taschenbuch der Technischen Akustik -> Beurteilungsgrößen") als Alternative der einseitig gleitende Mittelwert exponentielller Bauart, kurz EMA, vorgeschlagen, da dieser in eine rekursive Form überführt werden kann. Letzterer weicht laut verschiedener Quellen nur geringfügig vom beidseitigigen gleitenden Mittelwert ab (siehe: "Dr. A. Ziegler, Digitale Signalverarbeitung bei Erschütterungsmessungen
Dr. A. Ziegler; ZIEGLER CONSULTANTS, Zürich
1", http://www.z-c.ch/Publikationen/Digitale%20Signalverabeitung.pdf).
Problemannäherung:
Mit diesem Wissen habe ich mir den Quellcode zur passenden matlab-Funktion angeschaut, da ich mir nicht sicher war, wie Matlab arbeitet und wie die Fensterbreite (zeitbewertung 'fast' oder 'slow') für diesen Fall zu definieren ist. Siehe da, mit edit('movavg') findet man dann besagte rekursive Formel mit vorheriger Berechnung der 'Dämpfungskonstante':
Da ich zuerst nicht wusste wie die schlauen Programmierer von Matlab da nun wieder drauf gekommen sind, bin ich der in Wiki zitierten Quelle gefolgt und habe mir die relevanten Passagen in "Jens-Peter Kreiß, Georg Neuhaus: Einführung in die Zeitreihenanalyse. Springer, 2006. ISBN 978-3-540-25628-1" angetan. Wenn man dann bemerkt, dass das dortige alpha nicht exakt das Gleiche ist, hat man es schon fast geschafft. Mit "alpha* = 1 - alpha" kann man die beiden rekursiven Formen ineinader überführen, wenn da nicht der Faktor 2 im Quellcode stehen würde!
Konkret sähe ein Quellcode der Kenrfunktion wie folgt aus, wobei mir nicht klar ist, ob X=1 ("Wiki bzw. Kreiß") oder X=2 ("movavg"). (Anmerkung: Für Auswertungen von akustischen Signalen vorher quadrieren und später wurzelziehen...):
Code:
function Signal_Neu = movavg2(Signal, Zeitbereich, Messtastfrequenz) % Konstante bestimmen:
fensterlaenge = Zeitbereich * Messtastfrequenz;
alpha = X/(fensterlaenge + 1);
a = 1 - alpha; % Umrechnung für Schreibweise nach "Kreiß"
% Neues Signal vorbelegen:
Signal_Neu = [zeros(length(Signal),1)];
% Initialisierung des ersten Werts:
Signal_Neu(1) = Signal(1);
% Itterative Berechnung aller weiteren Werte: for s = 2:length(Signal) % Schreibweise siehe Quelle (anders als in Quelltext von movavg):
Signal_Neu(s) = (1 - a) * Signal(s) + a * Signal_Neu(s-1);
end
Um die Verwirrung noch rund zu machen folgt man dann einfach noch den Hinweisen aus "DIN EN 61672-1:2003 Seite 9, Abb. 1" und schreibt diese exponentielle Glättung als diskretes Filter des PT1-Glied, oder übernimmt dafür einfach die Lösung aus http://de.wikipedia.org/wiki/PT1-Glied im letzen Abschnitt (Zeitdiskretes PT1-Glied). Wie nicht weiter verwunderlich ist dies wieder eben oben zitierte rekursive Formel, aber diesmal wieder mit T* = 1/...., also X=EINS.
Die Frage lautet: EINS oder ZWEI?
Vielen Dank fürs Mitverfolgen und Lösen des Problems bzw. meines Gedankenfehlers,
Germanus
Einstellungen und Berechtigungen
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.
RSS