quadratischer mittelwert

[tobias hofer (Gast)]

Hallo

Ich möchte mit einem DSP den koninuierlichen quadratischen mittelwert 
eines eingangs signales berechnen.
d.h.

ich habe einen ad wander der mit 10kHz sample frequenz läuft. nun möchte 
ich mit diesem kontinuierlichen eingangsstrom von analog werten den 
laufenden quadratischen mittelwert berechnen. ich möchte nicht z.b. 1024 
sample einlesen und dann in einer schlaufe den rms wert berechnen.

kann mir hier jemand auf die sprünge helfen?

[Detlef _. (detlef_a)]

Hallo,

kontinuierlich das qaudratische Mittel der Wandlerwerte berechnen. Wenn 
x(k+1) der neue Wandlerwert ist, dann berechnest Du 
y(k+1)=((n-1)*y(k)+x(k+1)*x(k+1))/n .y(k+1) folgt dann dem Mittelwert 
der Quadrate von x, je größer n ist, umso langsamer. Nach Wuzelziehen 
und Umskalieren kommt hast Du dann RMS.

Cheers
Detlef

[tobias hofer (Gast)]

Hallo,

Kannst du das ein bisschen näher ausführen, damit ich das auch 
mathematisch verstehe?
Wie wird n festgelegt?

Gruss

Tobias

[Franko P. (sgssn)]

Hallo Tobias

da ist noch ein Punkt: Was hast du für einen A/D-Wandler ? Ist das ein 
externer bipolarer A/D-Wandler, der ein vorzeichenbehaftetes Ergebnis 
liefert, oder ist das ein unipolarer A/D-Wandler ?

Gerhard

[Detlef _. (detlef_a)]

Hallo,

wenn ich das richtig verstanden habe, kommen kontinuierlich Werte von 
Deinem AD Wandler und Du möchtest laufend den quadratischen Mittelwert 
berechnen. Ich schlage vor,das folgendermaßen zu machen:

Die Werte von Deinem AD-Wandler nenne ich x(k). Zunächst brauchst Du für 
den RMS das durchschnittliche Quadrat der x(k). Also berechnest Du 
x(k)*x(k); Die einfachste Möglichkeit, einen laufenden Durchschnitt der 
x(k)*x(k) zu bilden ist die Anwendung folgender Formel:

y(k+1)=((n-1)*y(k)+x(k+1)*x(k+1))/n

Die Folge der y(k) folgt dann langsam dem mittleren Quadrat der x(k). Je 
größer n ist, umso 'träger' folgt y(k) den x(k)*x(k). Das neue y(k+1) 
ist ja die Summe aus dem alten y(k) mit dem Vorfaktor(n-1)/n und dem 
neuen x(k)*x(k) mit dem Vorfaktor 1/n. In der Summe geben die 
Vorfaktoren den Wert 1. Signalverarbeitungstechnisch ist die Rechnung 
nen Tiefpaß erster Ordnung, je größer n ist umso größer ist dessen 
Zeitkonstante.

RMS bekommst Du, indem Du die Wurzel aus den y(k) ziehst.  Für eine 
Integerrechnung ist es immer günstig n als Potenz von 2 zu wählen weil 
dann die Division einfach wird. Vorsicht auch bei Integerrechnung mit 
dem Zahlenbereich, 10 Bit AD Werte haben nach der Multiplikation 20Bit, 
für den Vorfaktor auch noch genügend Raum lassen, mit 32 Bit rechnen.

Cheers
Detlef

[Klaus (Gast)]

Hi Tobias,

schau mal unter
http://de.wikipedia.org/wiki/Mittelwert
nach.

Bye
Klaus

[Detlef _. (detlef_a)]

Ähm, Gerhard hat natürlich Recht, ich bin bei meiner Einlassung davon 
ausgegangen, daß der AD-Wandler 'DC-freie' Werte liefert, also z.B -1023 
für -1V, 0 für 0V und 1023 für 1V. Falls der AD Wandler das anders 
macht, brauchst Du nicht nur das quadratische Mittel sondern auch den 
laufenden Mittelwert Deiner Eingangsfolge, den Du dann subtrahierst.

Cheers
Detlef

[tobias hofer (Gast)]

Hallo Detlef

Mein AD Wandler zahlenformat ist signed integer.

Nun habe ich das ganze vestanden.

besten dank, werde das mal programmieren.

gruss tobias

[tobias hofer (Gast)]

hallo Detlef

Noch eine letzte Frage, gerne würde ich wissen wie die mathematische 
herleitung deiner berechnung funktioniert, ebenso wie ich n zur cutoff 
frequenz des filters setze.

hast du gute lieteratur oder links?

gruss
tobias

[Detlef _. (detlef_a)]

Mit den einfachen x(k) als Eingangswerte lautet die Differenzengleichung 
für dieses IIR-Filter 1. Ordnung
y(k+1)=((n-1)y(k)+x(k))/n =((n-1)/n)*y(k) +x(k)/n

Die z-Übertragungsfunktion davon lautet:
H(z)=(1/n)/(z-(n-1)/n)

für z=exp(j*w) bekommt man den Frequenzgang des Filters, hier mal für n 
von 2 bis 100 gezeichnet, Matlab script hängt auch hinten dran.

Literatur/links: Digitale Filtertheorie ist nen weites Feld, da gibts 
unendlich viel Zeug.

Mein Handbuch aufm Schreibtisch ist
Lacroix, A.: Digitale Filter, Oldenburg Verlag München Wien 1996
Standard ist auch:
Oppenheim; Schafer: Zeitdiskrete Signalverarbeitung,Oldenburg Verlag 
München Wien 1992

Cheers
Detlef

clear
nn=(2:5:100);
HH=[];

for(k=1:length(nn))
n=nn(k);
H=freqz(1/n,[1 -(n-1)/n],256);
HH=[HH H];
end;
semilogx((0:255)/256,20*log10(abs(HH)))
grid

title('IIR 1. Ordnung')
xlabel('f/f_N_y_q_u_i_s_t')
ylabel('Dämpfung [dB]')

return

[tobias hofer (Gast)]

Hallo detlef

Kannst du mir sagen wie du von

y(k+1)=((n-1)y(k)+x(k))/n =((n-1)/n)*y(k) +x(k)/n auf
H(z)=(1/n)/(z-(n-1)/n) kommst?

Ich kenne mich nur mit laplace aus, leider nicht mit der 
z-transformation

danke,
gruss tobias

[Detlef _. (detlef_a)]

Normalerweise schlage ich die z-Übertragungsfunktionen der 
Filterstrukturen in der angegebenen Literatur nach. Herleitung explizit: 
Für die angegebene Differenzengleichung lautet die Anwortfolge auf einen 
Einheitsimpuls (eine 1, danach nur noch Nullen)

h(k)= (1/n)*((n-1)/n)^k

Die z-Transformierte dieser Impulsantwort ist die Übertragungsfkt. des 
Filters. Die Tansformierte lautete (siehe Definition der 
z-Transformation)

H(z)= (1/n)*(Summe k=0 bis unendlich) (((n-1)/n)^k*z(^-1)))
     =(1/n)*(Summe k=0 bis unendlich) (((n-1)/(z*n))^k)

Mein Bronstein/Semendjajew Mathehandbuch sagt dazu:

H(z) = (1/n)/(1-((n-1)/(z*n))
H(z) = (z/n)/(z-(n-1)/n)

Dieses H(z) ist um den Faktor z unterschiedlich zu dem von mir oben 
angegebenen, die zugehörigen Folgen sind um einen Wert gegeneinander 
verzögert (siehe dazu auch Verschiebungsatz der z-Tranformation). Das 
spielt keine wesentliche Rolle, einmal ist der Ausgang des Filters 'vor 
dem Summierer', einmal dahinter.

Danke für die interessante Frage, Antwort war Spaß.

Cheers
Detlef

[Detlef _. (detlef_a)]

Ähm, so

H(z)= (1/n)*(Summe k=0 bis unendlich) (((n-1)/n)^k*z^(-k)))

[Rumpel (Gast)]

Sieh Dir mal die Begriffe "mittlere quadratische Kontingenz" an.

[tobias hofer (Gast)]

Hallo Delef_a

Wie kann ich in deinem Matlab beispiel vom ausgelegten filter
mit hilfe matlab die sprungantwort berechnen, und die pole und 
nullstellen berechnen.

für die step response kenne ich den befehl step(sys). nur ist mir hier 
nicht klar wie ich sys bekomme.

gruss tobias

[Detlef _. (detlef_a)]

Hi,

die Nullstellen des Filters sind die Nullstellen des Zählerpolynoms, die 
Polstellen sind die Nullstellen des Nennerpolynoms.
Für H(z) = (z/n)/(z-(n-1)/n) also denkbar schlicht Nullstelle bei z=0 
und Polstelle bei z= (n-1)/n. Die Nullstellen von Polynomen berechnest 
Du bei Matlab mit dem Befehls 'roots'.

'step' berechnet anscheinend die Sprungantwort, 'impulse' berechnet die 
Impulsantwort. Die beiden Befehle kenne ich nicht, die gehören zu 
irgendeiner toolbox oder zu Simulink. Das geht aber viel einfacher mit 
dem Befehl 'filter', der berechnet digitale filter, mal in dessen 
Beschreibung reinkucken.

Also (für n=100):
impulsantwort= filter([ 0 1/100],[1 -99/100],[1 zeros(1,1000)]);
sprungantwort= filter([ 0 1/100],[1 -99/100],[   ones(1,1000)]);

Impulsantwort springt auf 1/100 und geht dann exponentiell gegen 0, 
Sprungantwort geht exponentiell gegen 1, muß sie ja auch, als 
Mittelwertbildner.

Cheers
Detlef