Forum: Digitale Signalverarbeitung / DSP / Machine Learning Einschwingzeit von IIR Filtern


von Manfred M. (bittbeisser)


Angehängte Dateien:

Lesenswert?

Ich habe da wieder mal eine Anfängerfrage. Und ja, die Überschrift ist 
etwas zu kurz gegriffen.

Was ich wissen möchte ist, nach welcher Zeit der Einschwingvorgang unter 
einem vorgegebenen Wert abgeklungen ist. In meinem Fall ist dieser Wert 
1/10000000 (-140dB).

Ich benötige das für ein Testprogramm, welches das tatsächliche 
Verhalten eines Filters testen soll. Ich will also sehen, wie sich ein 
Filter verhält, wenn die berechneten Koeffizienten nicht exakt 
eingehalten werden. Ich versuche also die Werte zu messen.

Dabei bekomme ich aber immer wieder Fehler, durch das 
Einschwingverhalten.
Als Beispiel: Ich habe hier ein Tschebyscheff-Filter, Tiefpass mit 
Grenzfrequenz 0.0125 (100Hz bei Fa=8000kHz) und 0.1% Ripple. Die Messung 
sieht schlimm aus. Aber bei 0% Ripple (Butterworth) fällt der 
Frequenzgang perfekt ins Nirwana. Siehe Bild. Der linke Teil ist zum 
Vergleich  mit Ripple 0 gemessen.

Bisher habe ich festgestellt, das dieses Problem verstärkt auftritt, 
wenn die Grenzfrequenz sich den Bandenden nähert, also gegen f=0 oder 
f=fa/2 geht oder wenn der zugelassene Ripplefaktor größer wird.

Gibt es da eine einfache Methode, wie man diese Verhalten abschätzen 
kann?

von Christoph db1uq K. (christoph_kessler)


Lesenswert?

Einschwingzeit und auch die Filter-Durchlaufzeit sind vor allem von der 
Filtercharakteristik abhängig. Die gutmütigen Typen Bessel, Gauss oder 
noch Butterworth sind schneller beruhigt als Tschebyscheff oder 
Cauer/Elliptic
so habe ich das jedenfalls bisher verstanden.

Etwas anderes sind die Störungen durch Rundungsfehler, die ja hier 
anscheinend gemessen werden sollen. Das langanhaltende Gezappel weit 
unten kommt vermutlich eher daher. Aber bitte, Experten vortreten, ich 
weiss dazu nichts genaues.

von lrep (Gast)


Lesenswert?

Manfred M. schrieb:
> Ich versuche also die Werte zu messen.

Das wird nicht ganz leicht werden.
140dB bedeutet immerhin ein Spannungsverhältnis 1:10 Millionen.
Wen du also ein Meßgerät hast, das 10µV in der letzten Stelle auflöst, 
müsstest du mit 100V anfangen.

von Manfred M. (bittbeisser)


Lesenswert?

Vielleicht hätte ich besser schreiben sollen, das ich die Messung 
simuliere. Also 32 Bit oder double Variablen ist das noch gut 
darstellbar.

Und das "gezappel" sind definitiv keine Rundungsfehler, Ich kriege das 
weg, wenn ich die Wartezeit bis zum Start der simulierten Messung 
verlängere. Die Störungen kommen meiner Meinung nach daher, das die 
Überschwinger gemessen werden, die zu dieser Zeit offenbar größer sind 
als das eigentliche Nutzsignal sein sollte.

Ich starte für jedes Pixel in x-Richtung eine neue Messung mit der 
zugehörigen Frequenz. Damit der Bildaufbau nicht ewig dauert, möchte ich 
mit der Messung nicht länger als unbedingt nötig warten.

Ich hatte gehofft, dass es dafür eine Faustregel gibt, wo dann am Ende 
ein Multiplikator für die Wartezeit rauskommt.

Bei FIR Filtern ist das ja einfach. Da warte ich einfach so viele 
Taktzyklen, wie Taps vorhanden sind.

von Markus B. (russenbaer)


Lesenswert?

Hallo,

Ich nehme an das Dein Filter ein digitales ist und Dein Programm keine 
Simulation eines LC-Filters.

Dann stellt sich mir die Frage warum Du für Deinen Frequenzgang, den Du 
ja über Koeffizientenvariationen rechnen willst (d.h. verschiedene 
Koeffizientensets), nicht einfach einen Dirac-Puls in Dein Filter 
schickst?

Wenn Du Deine Impulsantwort hinreichend lange rechnest und danach eine 
FFT aller Ausgangssamples machst wirst Du Dein Ergebnis genau genug 
bekommen.

Ein "Messaufbau" mit einem stepped sweep ist im digitalen Fall nicht 
unbedingt notwendig.

Der (digitale) Diracpuls
h[n] = 1 für n = 0 und
h[n]=0 für n!=0

lg
Markus

von Manfred M. (bittbeisser)


Lesenswert?

Ja, hatte ich wohl vergessen zu erwähnen. Das ist mein erstes Experiment 
mit digitalen IIR Filtern. Die Koeffizienten habe ich nach den hier 
schon mehrfach zitiertem DSP Guide (Chapter 20, Table 20-4) 
http://www.dspguide.com/ch20/4.htm berechnet.

> Dann stellt sich mir die Frage warum Du für Deinen Frequenzgang, den Du
> ja über Koeffizientenvariationen rechnen willst (d.h. verschiedene
> Koeffizientensets), nicht einfach einen Dirac-Puls in Dein Filter
> schickst?

Das habe ich jetzt leider nicht wirklich verstanden. Ich sagte ja 
bereits, das ich in diesem Bereich Anfänger bin. Ich bin nur ein 
Funkamateur und baue an einen Morsedecoder. Da ich mit den bisherigen 
Ergebnissen nicht ganz zufrieden bin, versuche ich etwas tiefer in die 
Signal Filterung einzusteigen.

> Ein "Messaufbau" mit einem stepped sweep ist im digitalen Fall nicht
> unbedingt notwendig.

Diesen Weg habe ich aber bewusst gewählt, da ich damit auch testen kann, 
was im Falle fehlerhafter oder unpräziser Koeffizienten passiert. Ich 
weiss, das dies sowas wie ein "brute force" Angriff ist, aber der 
funktioniert im Prinzip auch, wenn ich die Übertragungsfunktion nicht 
wirklich kenne.

Das mit dem Dirac-Puls habe ich jetzt auch nicht verstanden. Ich dachte 
das sei nur ein mathematisches Hilfsmittel, das ich so praktisch nicht 
einsetzen kann. Daher habe ich bisher in meinem Testprogramm nur die 
Sprungantwort verwirklicht. Die ist zwar auch interessant, hilft mir 
jetzt aber auch nicht weiter.

Aber wenn es alternative Testmöglichkeiten gibt, würde ich diese gerne 
testen, falls ich deren Methode verstehen und in ein Programm übersetzen 
kann.

von C Programmierer (Gast)


Lesenswert?

Manfred M. schrieb:
> Das mit dem Dirac-Puls habe ich jetzt auch nicht verstanden. Ich dachte
> das sei nur ein mathematisches Hilfsmittel, das ich so praktisch nicht
> einsetzen kann.

Die Fourier-Tranformierte der Impulsantwort ist der Frequenzgang des 
Systems.

M.E. musst du hier allerdings den analogen Dirac Impuls mit der Höhe 
1/tA = fA nehmen.

von Fpgakuechle K. (Gast)


Lesenswert?

Manfred M. schrieb:
> Vielleicht hätte ich besser schreiben sollen, das ich die Messung
> simuliere. Also 32 Bit oder double Variablen ist das noch gut
> darstellbar.

Ist die "Simulation" noch realistisch?

lrep schrieb:
> 140dB bedeutet immerhin ein Spannungsverhältnis 1:10 Millionen.
> Wen du also ein Meßgerät hast, das 10µV in der letzten Stelle auflöst,
> müsstest du mit 100V anfangen.


Oder anders ausgedrückt, das durch die "Ripple" erzeugte Rauschen ist so 
gering das es wohl kaum die Schaltung beeinflußt. 105 db Dämpfung!

Wie groß ist die Genauigkeit der Wandlers - 10 bit? - 12 bit? Welchen 
"Ripples" würde das entsprechen? IMHO  ist das "Rauschen" durch die 
Filterkurve fern vom "Grundrauschen" und damit komplett unbedeudent.

MfG,

von Markus B. (russenbaer)


Lesenswert?

Hallo,

Abgesehen davon das Fpga Kuechle et.al. recht hat mit der Frage nach der 
Sinnhaftigkeit bei Dämpfungen >100dB den Ripple anzuschauen...


Die Fouriertransformierte der Impulsantwort ist die 
Übertragungsfunktion.
Du suchst die Übertragungsfunktion für verschiedene Koeffizientensets 
Deines Fílters.

Wenn Du nun im Zeitbereich ein Signal in Dein System schickts das wie 
folgt ausschaut:

1, 0, 0, 0, 0, 0 ....

und das aufgezeichnete Ausgangssignal Fouriertransormierst erhältst Du 
die gesuchte Übertragungsfunktion (mit je nach Aufzeichnungslänge 
größeren oder kleinerem Restfehler)


Noch genauer (begrenzt durch die Rechengenauigkeit Deines PCs und 
wahrscheinlich auch schneller) kannst Du die Übertragungsfunktion 
berechnen indem Du die Differenzengleichung Deines Filters z 
transformierst und auf die Polynombruchform der Übertragungsfunktion 
umformst. Wenn Du das ganze dann am Einheitskreis auswertest ( z = 
exp(j2Pi f/fs) hast Du die Übertragungsfunktion die Du suchst.

Hier ein Beispiel für ein Biquad-Filter, sei:

h[n] gesuchte Impulsantwort, H(z) Übertragungsfunktion des Systems
y[n] Ausgang deines Filters zum Zeitpunkt n, Y(z) z-Transformierte des 
Ausgangssignals
x[n] Eingang deines Filters zum Zeitpunkt n, X(z) z-Transformierte des 
Eingangssignals
a1, a2 Rückkopplungskoeffizienten
b0, b1, b2 Vorwärtskoeffizienten


y[n] = x[n] b0 + x[n-1] b1 + x[n-2] b2 - y[n-1] a1 - y[n-2] a2

Y(z) = X(z)b0+ X(z) z^-1 b1 + X(z) z^-2 b2 - Y(z) z^-1 a1 - Y(z) z^-2 a2 
= X(z) (b0+z^-1 b1 + z^-2 b2) - Y(z) (z^-1 a1 + z^-2 a2)

H(z) = Y(z)/X(z) = (b0+z^-1 b1 + z^-2 b2) / (1 + z^-1 a1 + z^-2 a2)


setze nun z=exp(j*2*pi*f/fs)
Dann kannst Du sowohl 20*log10(abs(H(z))) als auch arg(H(z)) berechnen 
-> gewünschte Übertragungsfunktion (für verschiedene f)

: Bearbeitet durch User
von derguteweka (Gast)


Lesenswert?

Moin,

Manfred M. schrieb:
> Ich hatte gehofft, dass es dafür eine Faustregel gibt, wo dann am Ende
> ein Multiplikator für die Wartezeit rauskommt.

Ah, ich glaub' ich ahne auf was du raus willst. Nein, zumindest mir ist 
keine Regel bekannt. Die Dinger heissen ja IIR, wo ein I fuer Infinite 
(=so lange, dass nur Chuck Norris bis dahin zaehlen kann) steht. Leider 
kann dieses Infinite auch fuer das Einschwingen der Filterantwort 
gelten.
Insbesondere dann, wenn man auf IIR dann ausweicht, weil einem die 
Filterordnung bei FIR zu hoch wird - sprich: Wenn die Polstellen dem 
Einheitskreis nahe auf die Pelle ruecken...

Gruss
WK

von Manfred M. (bittbeisser)


Lesenswert?

Also ich befinde mich noch weitgehend in der Planungsphase und will 
erstmal sehen was machbar ist. Testen tue ich mit WAV Dateien, die 
teilweise mit Recordern aus dem Musiker Bereich kommen (24Bit 
Auflösung).

Das der Berich bis -140dB mehr in dem Bereich Kosmetik fällt ist mir 
auch klar. Den habe ich aber auch soweit ausgedehnt, weil ich sonst bei 
FIR Tiefpässen bei Anwendung einiger Fensterfunktionen im Sperrbereich 
nur noch vereinzelte Kamelhöcker gesehen hatte.

> Ah, ich glaub' ich ahne auf was du raus willst.

Endlich versteht mich mal einer ;-)
Ne, im Ernst, ich wollte gerade erklären, das ich die benötigte Zeit 
wohl vorher anhand der Sprungantwort messen muss (beim Tiefpass geht das 
ja gut).
Dabei ist mir aufgefallen, das man schon ganz am Anfang der 
Sprungantwort gravierende Unterschiede erkennen kann.

Also ein TP mit fg=1000 Hz und 0.5% Ripple erreicht den Einheitswert 
erstmalig bereits nach 1.4 ms. Der gleiche TP mit fg=150 Hz aber erst 
nach 9.4 ms.

Wenn es da keine Faustregel gibt, hilft vielleicht eine Tabelle.

von Manfred M. (bittbeisser)


Lesenswert?

Ich habe mal angefangen einige Messreihen aufzunehmen. Die Darstellung 
ist allerdings schwierig, da sie eigentlich 4-dimensional sein müsste 
und deutlich aufwändiger ist als ich anfangs vermutete.

Aber es scheint so, als wenn die Einschwingzeit deutlich stärker vom 
Frequenzverhältnis abhängt als von allen anderen Faktoren. Mein erster 
Verdacht war, das die Einschwingzeit zu niedrigen Frequenzen nach 
e-Funktion zunimmt, aber eine erste Nachrechnung tendiert eher zu einer 
hyperbolischen Abhängigkeit. Die Abhängigkeit von der Welligkeit 
erscheint eher gering, da die beiden ersten Kurven, die ich aufgenommen 
habe, mit geringem Abstand parallel verlaufen.

Die Abhängigkeit von der Anzahl der Filterpole habe ich noch nicht 
grafisch untersucht, aber sie scheint auch deutlich stärker zu sein, als 
die vom zugelassenen Ripple-Faktor.

von C Programmierer (Gast)


Lesenswert?

Du brauchst die Einschwingzeit nicht zu berechnen oder zu schätzen, 
schau dir den Post von Markus nochmal an. Genauer und für den Rechner um 
ein großes Vielfaches schneller wirst du den Frequenzgang nicht 
berechnen können.

Markus B. schrieb:
> y[n] = x[n] b0 + x[n-1] b1 + x[n-2] b2 - y[n-1] a1 - y[n-2] a2
>
> Y(z) = X(z)b0+ X(z) z^-1 b1 + X(z) z^-2 b2 - Y(z) z^-1 a1 - Y(z) z^-2 a2
> = X(z) (b0+z^-1 b1 + z^-2 b2) - Y(z) (z^-1 a1 + z^-2 a2)
>
> H(z) = Y(z)/X(z) = (b0+z^-1 b1 + z^-2 b2) / (1 + z^-1 a1 + z^-2 a2)
>
> setze nun z=exp(j*2*pi*f/fs)

von derguteweka (Gast)


Lesenswert?

Moin,

Bloss mal so ins Blaue reinfabuliert, ohne 's wirklich zu wissen: 
Koennt's sein, dass diese Einschwingzeit irgendwas mit der 
Gruppenlaufzeit des jeweiligen Filters zu tun hat? Also so in der Art: 
Je groesser die Gruppenlaufzeit(differenz im Durchlassbereich) desto 
laenger die Einschwingzeit?

Gruss
WK

von Gört (Gast)


Lesenswert?

derguteweka schrieb:
> Koennt's sein, dass diese Einschwingzeit irgendwas mit der
> Gruppenlaufzeit des jeweiligen Filters zu tun hat

Zusammenhängen tut es schon, das selbe ist es nicht. Siehe
Beitrag "Gruppenlaufzeit vs. Einschwingvorgang (Verständnisproblem)"

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.