www.mikrocontroller.net

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


Autor: Mr Bean (Gast)
Datum:

Bewertung
0 lesenswert
nicht 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
Autor: Mark Brandis (markbrandis)
Datum:

Bewertung
0 lesenswert
nicht 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/physikalische...

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.

Autor: Mark Brandis (markbrandis)
Datum:

Bewertung
0 lesenswert
nicht lesenswert
Außerdem: Ordnung d. FIR-Filters = Anzahl d. Filterkoeffizienten minus 1

Autor: Mr Bean (Gast)
Datum:

Bewertung
0 lesenswert
nicht 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ß!

Autor: Mark Brandis (markbrandis)
Datum:

Bewertung
0 lesenswert
nicht 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?

Autor: Mr Bean (Gast)
Datum:
Angehängte Dateien:

Bewertung
0 lesenswert
nicht 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ß!

Autor: Mark Brandis (markbrandis)
Datum:

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

Autor: Mr Bean (Gast)
Datum:

Bewertung
0 lesenswert
nicht 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ß!

Autor: Mark Brandis (markbrandis)
Datum:

Bewertung
0 lesenswert
nicht 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.

Autor: fileter (Gast)
Datum:

Bewertung
0 lesenswert
nicht 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.

Autor: Mr Bean (Gast)
Datum:

Bewertung
0 lesenswert
nicht 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:
  long mittelwert_filter(long newval)
{
   static short n = 1;          //Zählvariable anlegen
   static unsigned long avgsum;      //Variable zur Summenspeicherung anlegen
   unsigned long avg;          //Variable zur Mittelwertspeicherung anlegen
   
   if (n<32)               //Mit dieser if Schleife wird die Fenstrgröße festgelegt
                       //hier sollte eine 2er potenz für die Größe gewählt werden.
   {
      avgsum += newval;          //Wert n-mal aufaddieren
      avg = avgsum/n;          //Summer durch n teilen
      n++;                //n um eins erhöhen
   }
   else 
   {
      avgsum -= avgsum/32;
      avgsum += newval;
      avg = avgsum/32;
   }
   return avg;

}

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

Autor: Stefan (Gast)
Datum:

Bewertung
0 lesenswert
nicht 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!

Autor: Mark Brandis (markbrandis)
Datum:

Bewertung
0 lesenswert
nicht 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.

Autor: Detlef _a (detlef_a)
Datum:

Bewertung
0 lesenswert
nicht 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

Autor: Mark Brandis (markbrandis)
Datum:

Bewertung
0 lesenswert
nicht 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 :-)

Autor: fred maier (Gast)
Datum:

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

Autor: Detlef _a (detlef_a)
Datum:

Bewertung
0 lesenswert
nicht 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

Autor: Mr Bean (Gast)
Datum:

Bewertung
0 lesenswert
nicht 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

Autor: Mr Bean (Gast)
Datum:

Bewertung
0 lesenswert
nicht lesenswert
sorry, sollte

heißen.

Gruß Bean

Autor: gast (Gast)
Datum:

Bewertung
0 lesenswert
nicht 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.

Autor: Mr Bean (Gast)
Datum:

Bewertung
0 lesenswert
nicht 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

Autor: Stefan (Gast)
Datum:
Angehängte Dateien:

Bewertung
0 lesenswert
nicht 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.

Autor: Mr Bean (Gast)
Datum:

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

Grüße

Dominik

Autor: gast (Gast)
Datum:

Bewertung
0 lesenswert
nicht 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.

Autor: Mr Bean (Gast)
Datum:

Bewertung
0 lesenswert
nicht 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

Autor: gast nun wach (Gast)
Datum:

Bewertung
0 lesenswert
nicht 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.

Autor: Mr Bean (Gast)
Datum:

Bewertung
0 lesenswert
nicht 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

Autor: gast nun wach (Gast)
Datum:

Bewertung
0 lesenswert
nicht 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.

Autor: Mr Bean (Gast)
Datum:

Bewertung
0 lesenswert
nicht 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

Autor: gast (Gast)
Datum:

Bewertung
0 lesenswert
nicht 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.

Antwort schreiben

Die Angabe einer E-Mail-Adresse ist freiwillig. Wenn Sie automatisch per E-Mail über Antworten auf Ihren Beitrag informiert werden möchten, melden Sie sich bitte an.

Wichtige Regeln - erst lesen, dann posten!

  • Groß- und Kleinschreibung verwenden
  • Längeren Sourcecode nicht im Text einfügen, sondern als Dateianhang

Formatierung (mehr Informationen...)

  • [c]C-Code[/c]
  • [avrasm]AVR-Assembler-Code[/avrasm]
  • [code]Code in anderen Sprachen, ASCII-Zeichnungen[/code]
  • [math]Formel in LaTeX-Syntax[/math]
  • [[Titel]] - Link zu Artikel
  • Verweis auf anderen Beitrag einfügen: Rechtsklick auf Beitragstitel,
    "Adresse kopieren", und in den Text einfügen




Bild automatisch verkleinern, falls nötig
Bitte das JPG-Format nur für Fotos und Scans verwenden!
Zeichnungen und Screenshots im PNG- oder
GIF-Format hochladen. Siehe Bildformate.
Hinweis: der ursprüngliche Beitrag ist mehr als 6 Monate alt.
Bitte hier nur auf die ursprüngliche Frage antworten,
für neue Fragen einen neuen Beitrag erstellen.

Mit dem Abschicken bestätigst du, die Nutzungsbedingungen anzuerkennen.