Forum: Digitale Signalverarbeitung / DSP / Machine Learning Berechnung/Hintergrund von Digitalem Tiefpass (gleitender Mittelwert)


von Mr Bean (Gast)


Lesenswert?

Hallo

Ich habe mit meinem µController eine gleitende Mittelwertberechnung 
realisiert. Wenn ich mir die Werte die errechnet werden in Excel 
importiere, ist das Tiefpassverhalten gut zu sehen.
Wie kann ich diesen Tiefpass nun berechnen? Also Grenzfrequenz, 
Eckfrequenz, Ordnung des Tiefpasses...
Generell würde mich der theoretische Hintergrund interessieren.
Sicher könnt ihr mir da weiter helfen.

Grüße!

Bean

: Verschoben durch Admin
von Mark B. (markbrandis)


Lesenswert?

Die Grenzfrequenz berechnet man nicht, die gibt man vor, wenn man ein 
Filter bastelt. Also genauer: Wenn man die Koeffizienten des Filters 
berechnet. Und dafür gibt's z.B. MATLAB :-)

Für die Theorie siehe z.B. hier: 
http://wwwex.physik.uni-ulm.de/lehre/physikalischeelektronik/phys_elektr/node65.html

Ach so, jetzt versteht ich die Frage erst so richtig glaube ich: Die 
Koeffizienten eines FIR-Filters erhält man ja, wenn man sich die 
Impulsantwort anschaut. Von dieser kommt man zum Frequenzgang durch eine 
inverse diskrete Fouriertransformation.

von Mark B. (markbrandis)


Lesenswert?

Außerdem: Ordnung d. FIR-Filters = Anzahl d. Filterkoeffizienten minus 1

von Mr Bean (Gast)


Lesenswert?

Hallo ja genau,

Also ich habe mir jetzt im Excel eine Sprungantwort des Filters zeichnen 
lassen. Aber ich habe diese Sprungantwort eben nur als "Bild". Wie 
bekomme ich jetzt die Werte raus mit denen ich weiterrechnen kann?

Gruß!

von Mark B. (markbrandis)


Lesenswert?

Ausgehend von der Sprungantwort müsstest Du die Ableitung bilden, um die 
Impulsantwort zu erhalten. Aber wozu die Mühe, wenn eine einzelne 1 ins 
Filer reingeschickt doch schon die Koeffizienten ergibt :-)
Die Zahlenwerte an sich hast Du doch, oder was steht in den Zellen 
Deiner Excel-Tabelle?

von Mr Bean (Gast)


Angehängte Dateien:

Lesenswert?

hmm ja, in den Zellen stehen Zahlenwerte. Ich messe einen Strom. Dieser 
steigt dann Sprunghaft an. So "simuliere" ich einen Sprung. (Siehe Bild) 
Die Werte des Mittelwerts nähern sich dann an den eigentlichen Strom an.
Ich brauch doch aber irgendwie die Zeitkonstante oder täusche ich mich 
da?
Wie kann ich mit den Werten die ich in Excel habe weiter rechnen?

Gruß!

von Mark B. (markbrandis)


Lesenswert?

Was willst Du jetzt berechnen, ein digitales FIR-Filter oder ein 
analoges RC-Glied? Das ist mir im Moment nicht so ganz klar.

von Mr Bean (Gast)


Lesenswert?

Einen digitalen Filter.
Die rote linie in dem Diagramm sind die Mittelwerte die in meinem 
Controller berechnet wurden. Ist also denke ich digital ;-)

Gruß!

von Mark B. (markbrandis)


Lesenswert?

Am einfachsten würdest Du die Grenzfrequenz erhalten, wenn Du einen 
Funktionsgenerator benutzt und einen Sinus mit variabler Frequenz duch 
Deinen µC schickst. Dann siehst Du mit einem Oszilloksop am Ausgang, bei 
welcher Frequenz sich Dein Eingangssignal mehr oder weniger schnell 
verabschiedet ;-)

Ah ja und Deinen Code könntest Du auch mal posten, zumindest den für die 
Berechnung relevanten Teil.

von fileter (Gast)


Lesenswert?

Wie bereits gesagt:

Filterordnung = Koeffizienten   -1 = Fensterlänge -1

Ein Mittelwertfilter mit gerader Anzahl Koeffizienten blockt alle 
Signale mit 2, 4, 6, ..., halbe Anzahl Koffizienten x Sampling Frequenz.

Am einfachsten kann man die Grenz und Eckfrequenzen ermitteln, indem man 
die Amplituden über Frequenz plottet.

von Mr Bean (Gast)


Lesenswert?

Hallo

Vielen Dank für die Antworten!

Hier mal noch der Code den ich zur Mittelwertbildung verwende. Stammt 
auch aus einem Beitrag hier im Forum:
1
  long mittelwert_filter(long newval)
2
{
3
   static short n = 1;          //Zählvariable anlegen
4
   static unsigned long avgsum;      //Variable zur Summenspeicherung anlegen
5
   unsigned long avg;          //Variable zur Mittelwertspeicherung anlegen
6
   
7
   if (n<32)               //Mit dieser if Schleife wird die Fenstrgröße festgelegt
8
                       //hier sollte eine 2er potenz für die Größe gewählt werden.
9
   {
10
      avgsum += newval;          //Wert n-mal aufaddieren
11
      avg = avgsum/n;          //Summer durch n teilen
12
      n++;                //n um eins erhöhen
13
   }
14
   else 
15
   {
16
      avgsum -= avgsum/32;
17
      avgsum += newval;
18
      avg = avgsum/32;
19
   }
20
   return avg;
21
22
}

Zum Messen der Grenzfrequenz mit dem Oszi noch eine Frage. Ich habe ja 
in dem Sinne keinen Ausgang. HIerfür bräuchte ich ja wieder einen DAC. 
Oder wie hast Du das gemeint @Marc Brandis ;-)

Gruß

Bean

von Stefan (Gast)


Lesenswert?

Dein TP per Mittelwertbildung sollte durch ein PT1-Glied gut angenähert 
sein.
In der Sprungantwort in Deinem Diagramm sucht Du 63% vom Maximalwert 
(Strom). Genau an dieser Stelle findest Du Deine Zeitkonstante!

von Mark B. (markbrandis)


Lesenswert?

Mr Bean schrieb:
> Zum Messen der Grenzfrequenz mit dem Oszi noch eine Frage. Ich habe ja
> in dem Sinne keinen Ausgang. HIerfür bräuchte ich ja wieder einen DAC.
> Oder wie hast Du das gemeint @Marc Brandis ;-)

Ja, so in der Art. Sinus-Signal mit variabler Frequenz --> ADC --> Rein 
ins das "Filter" --> Digitale Ausgangswerte aus Deinem Algorithmus --> 
DAC --> Eingang Oszilloskop.

von Detlef _. (detlef_a)


Lesenswert?

Das 'gleitende Mittelwertfilter', das Du beschrieben hast, ist ein IIR 
Filter mit der Differenzengleichung

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

Das Ding hat die z-Übertragungsfkt.

H(z)= 1/(z+1/32)

Ersetzen des z durch exp(j*w) liefert Dir den Frequenzgang, nen Tiefpaß 
1.Ordnung mit einer maximalen Sperrdämpfung bei der halben 
Abtastfrequenz von ca. 0.5dB.

Ich hoffe, dass ich mich nicht verrechnet habe.

Cheers
Detlef

von Mark B. (markbrandis)


Lesenswert?

Detlef _a schrieb:
> Das 'gleitende Mittelwertfilter', das Du beschrieben hast, ist ein IIR
> Filter mit der Differenzengleichung
>
> y(k)=-1/32*y(k-1)+x(k)

Nehmen wir an, diese Differenzengleichung sei korrekt. Nun die 
z-Transformation:

Mit z durchmultiplizieren:

Die Terme mit Y(z) auf eine Seite bringen:

Y(z) ausklammern:

Durch X(z) und durch (z+1/32) dividieren:

Dies ist die Übertragungsfunktion. Es gibt eine Nullstelle (z=0) und 
eine Polstelle (z=-1/32). Dies erscheint mir für ein IIR-Filter ehrlich 
gesagt auch plausibler als nur eine Polstelle :-)

von fred maier (Gast)


Lesenswert?

Dein Filter berechnet einen "modifizierten gleitenden Mittelwert" nicht 
den "gleitenden Mittelwert" im eigentlichen Sinne. (FIR Filter mit 
Mittelwert über die letzten N Werte)

von Detlef _. (detlef_a)


Lesenswert?

ja, meiner Überragungsfkt. fehlt das z.

Das macht aber für den Betragsfrequenzgang keinen Unterschied, weil die 
Entfernung von jedem Punkt des Einheitskreises zum Nullpunkt immer eins 
ist (das ist ja gerade das Schöne am Einheitskreis).

Die Eins steht dann im Zähler der Übertragungsfkt. wenn der 'alte wert 
des Zustandsspeichers' (also im Beispiel avgsum vorm update) 
Ausgangswert wäre.

Cheers
Detlef

von Mr Bean (Gast)


Lesenswert?

Hallo!

Vielen Dank für eure Antworten. Hat mir egentlich wirklich sehr 
geholfen.
Ich hab jetzt aber nochmal für mich die Formel aufgestellt und komme 
irgendwie auf was anderes.

Das ist was ich heraus bekomme. Wo liegt da mein Rechenfehler? Ich hab 
halt noch den Therm mit

drinnen.

Bin von folgender Formel ausgegangen:

Oder bin ich einfach von der falschen Anfangsformel ausgegangen?

Gruß!
Bean

von Mr Bean (Gast)


Lesenswert?

sorry, sollte

heißen.

Gruß Bean

von gast (Gast)


Lesenswert?

y ist der modifizierte, gleitende Mittelwert und u ist der Eingang.


Die Differenzengleichung ist dann

y_k+1=31/32 y_k+1/32 u_k+1

(oder evtl  y_k+1=31/32 y_k+1/32 u_k  was aber eher ungewöhnlich wäre.)

Im ersten Fall hat man  die Übertragungsfunktion

zY =31/32  Y +1/32 zU

U/Y= (1/32 )z/ (z-31/32)

Also eine Nullstelle bei z=0 und ein Pol bei z=31/32<1.


Somit ist das Filter stabil und zudem schwingt es nicht. Es hat eine 
statinäre Verstärkung von 1.
In der "Filterdarstellung.

U/Y= 1/ (1-31 z^1)

Es gibt also keine Verzögerung.

von Mr Bean (Gast)


Lesenswert?

Guten Morgen!

@Gast: Danke für die Antwort. Allerdings habe ich Schwierigkeiten 
nachzuvollziehen wie Du auf die Differentialgleichung kommst. Kannst Du 
das noch etwas näher erklären?
Hast Du die aus meinem C Code abgeleitet?
Außerdem wären ein paar Klammern glaube ich hilfreich.

Gruß

Bean

von Stefan (Gast)


Angehängte Dateien:

Lesenswert?

Wie ich schon am 06.07.2009 schrieb, kann man die Mittelwertbildung 
näherungsweise durch ein PT1-Glied beschreiben.
Sicher ist das Aufstellen und Berechnen der Übertragungsfunktion die 
100%-ige Variante... und wenn's gefordert ist, führt kein Weg daran 
vorbei!
Wenn's aber nur darum geht, das Verhalten des Mittelwertfilters in guter 
Näherung zu beschreiben, dann würde ich auch die wesentlich einfachere 
Näherungslösung bevorzugen. Exemplarisch mal an einem Mittelwertfilter 
mit N=4 dargelegt.

von Mr Bean (Gast)


Lesenswert?

Ja das stimmt, klar kann man das ganze annähern. Ich wollt aber trotzdem 
gerne die Übertragungsfunktion. :-)

Grüße

Dominik

von gast (Gast)


Lesenswert?

Hi Mr Bean,

also ich habe mir den C Code nochmal angeschaut.

Du hast

      avgsum -= avgsum/32;
      avgsum += newval;
      avg = avgsum/32;

Was abgesehen vom Rundungsverhalten, Rechenaufwand (?) und der 
Initialisierung sich wie

      avgsum -= avgsum/32; (1)
      avgsum += newval/32;  (2)
      avg = avgsum;     (3)
verhält. Du willst nur die Ü-Funktion, also ist dies egal.
Nimm

y - avgsum
u - newval

Also ist

y_(k+1)=y_(k)-y_(k)/32+u_(k+1)/32


(die ungewöhnliche Variante wäre   Vertauschung von (2), (3))



Dann die z-Trafo zur Ü-Funktionbestimmung (ohne Berücksichtigung von 
Anfangswerten)

zY(z)=31/32 Y(z)+1/32 U(z)

Umstellen, fertig.

Matlab liefert ein monotones Tiefpassverhalte mit -36 db bei der halben 
Abtastfrequnz. Grenzfrequenz -3db ca bei 1/200 der Abtastfrequenz.

von Mr Bean (Gast)


Lesenswert?

Also sehe ich das richtig dass sich die Differentialgleichung zu:

und die Übertragungsfunktion zu :

??

Wenn ich das dann durch Y ausklammern etc. umstelle, komme ich auf:

Fehlt mir da nicht noch ein ein z im Zähler?

Gruß und Danke

Bean

von gast nun wach (Gast)


Lesenswert?

Ah ich habe da ein Tipp-Fehler gemacht

zY(z)=31/32 Y(z)+1/32 zU(z) <-- z vor U(z)

Es ist keine Differentialgleichung, sondern eine Differenzengleichung.

von Mr Bean (Gast)


Lesenswert?

Hallo

Ah ja, jetzt passt es.
Aber der Pol bei 31/32 macht mich doch etwas stutzig, ich meine mich 
dunkel daran erinnern zu können dass bei einem stabilen System kein Pol 
in der rechten s-Halbebene vorkommen darf.
Ist aber auch schon ne ganze Weile her... :-)

Gruß

Bean

von gast nun wach (Gast)


Lesenswert?

Es ist ein zeitdiskretes System (Differenzengleichung) und kein 
zeitdiskretes (DGL)!

Für Stabilität müssen alle Pole in der z-Ebene (!) im Einheitskreis sein 
- 31/32 ist das.

Da der Pol auf der rechten Seite auf der Realteil-Achse ist, schwingt 
das System nicht.

von Mr Bean (Gast)


Lesenswert?

Hallo

Ich hab mir die Sache noch einmal angeschaut.  Ich habe leider noch 
immer Probleme mit der Formel ansich die sich aus dem Code ergibt. Muss 
ich die Summe die sich in der ersten if Schleife ergibt nicht auch bei 
der Transformation berücksichtigen?

Ich hab ja eigentlich sowas:

Wie kann ich eine Summe transformieren?

Sorry dass ich den Thread nochmal vorkram und nochmal frag, aber ich 
stolper immer wieder über diese Formel und bin mir nicht sicher ob die 
Transformation so stimmt.

Gruß
Bean

von gast (Gast)


Lesenswert?

Hi Bean,


Der Einfachheitshalber würde ich nur den Teil ab n=32 anschauen (mit 
z-Trafo). Nach einiger Zeit ist der Anfangsteil vernachlässigbar. 
(Filter ist asymptotisch stabil)

Auch kann man mit der z-Trafo nur Übertragungsfunktionen von 
zeitinvarianten Systemen anschauen.

Fazit: die ersten 31 Schritte ändern das Übertragungsverhalten nur am 
Anfang des "Filterns".

Summe transformieren: wie oben erwähnt musst du dazu die Filterformeln 
umschreiben.

Bitte 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
Noch kein Account? Hier anmelden.