Forum: Digitale Signalverarbeitung / DSP / Machine Learning Filtereigenschaften Bestimmen


von Kollege von oben (Gast)


Lesenswert?

Hallo Leute,

ich benötige die Filter eigenschaften eines FIR Filters.
ich habe die Koeffizienten und die Abtastrate.

Wie kann aus denn Koeffizienten die Grenzfrequenz bestimmen?

Habe natürlich auch ein Beispiel leider nur in Pascal:

filter wird in einem 2mSec Interrupt aufgerufen
1
     pmax_filter : fdaten = (n: 1;
2
                     a:(0.2  ,      0,     0,   0,   0);
3
                     b:(    1,   -0.8,     0,   0,   0);
4
                     x:(    0,      0,     0,   0,   0);
5
                     y:(    0,      0,     0,   0,   0));
6
7
function filter(var f : fdaten; inp: single) : single;
8
 var i : integer;
9
 begin
10
  with f do
11
  begin
12
   for i:=n downto 1 do
13
   begin
14
    x[i]:=x[i-1];
15
    y[i]:=y[i-1];
16
   end;
17
   x[0]:=inp;
18
   y[0]:=a[0] * x[0];
19
   for i:=1 to n do
20
    y[0]:=y[0] + a[i] * x[i] - b[i] * y[i];
21
   filter:=y[0];
22
  end; { with f }
23
 end;

von Dergute W. (derguteweka)


Lesenswert?

Moin,

Aus den Koeffizienten und der Verschaltung des Filters oder der Software 
musst du die Uebertragungsfunktion H(z) bestimmen.
Grenzfrequenz bei Filtern ist so ne Sache. Das ist nicht ganz eindeutig. 
Oft ist's die Frequenz, bei der der Amplitudengang von H(z) gegenueber 
dem Amplitudengang bei der Frequenz 0, pi/2 oder sonstirgendeiner 
Referenzfrequenz auf irgendeinen Bruchteil (oft -3dB, aber manchmal auch 
-6dB) zurueckgegangen ist.
Such' dir eine dir genehme Definition aus und loese dann H(z) danach 
auf.
Kann aber unangenehme, nicht analytisch loesbare Gleichungen ergeben.

Alternativ kannst du auch durch deine Software einen Diracimpuls laufen 
lassen und das Ergebnis fouriertransformieren. Das ergibt zumindest in 
der Theorie auch den Frequenzgang.

Gruss
WK

von Reinhard M. (Gast)


Lesenswert?

Aus der Differenzengleichung die Übertragungsfunktion H(z) bestimmen.

Die Grenzfrequenz liegt bei -3dB ( 1/sqrt(2) )

Tip: http://www.dspguide.com/pdfbook.htm

Kapitel33 zeigt wie es geht.

von Kollege von oben (Gast)


Lesenswert?

Danke euch schon mal.

Aber ich bin nicht so das mathe Ass.... :-(..

Könnt ihr vielleicht noch weiter helfen

von Dergute W. (derguteweka)


Lesenswert?

Moin,

Kollege von oben schrieb:
> Aber ich bin nicht so das mathe Ass.... :-(..

Und ich hab keinen Pascalcompiler da... Unn nu?

Geb' mal einen Dirac drauf und poste die Impulsantwort. Oder zumindest 
den Anfang.

Gruss
WK

von Kollege von oben (Gast)


Angehängte Dateien:

Lesenswert?

Hallo WK,

ich habs gerade mal auf die schnelle in C# umgebaut.

Könntest du mir weiterhelfen

von Dergute W. (derguteweka)


Lesenswert?

Moin,

Du kannsts auch noch auf Java, Python, Postscript oder Assembler umbauen 
- es wird nicht wirklich weiterhelfen.
Die Impulsantwort ist das Mittel zur Filteranalyse, nicht die 
Programmiersprache.

Kollege von oben schrieb:
> Könntest du mir weiterhelfen

Ungern, so wie's aussieht.

Gruss
WK

von Reinhard M. (Gast)


Lesenswert?

Angenommen deine Differenzengleichung im Programm ist

   y(n) = a0 * x(n) + b1 * y(n-1) + b2 * y(n-2)


mit z -Transformation: X(z) = Summe [x(n) * z^(-n)]

y(n) * z^0  = a0 * x(n) * z^0  + b1 * y(n) * z^(-1)  + b2 * y(n) * 
z^(-2)

// z^0 = 1

y(n)  = a0 * x(n)  + b1 * y(n) * z^(-1)  + b2 * y(n) * z^(-2)

Y(z) =  a0 * X(z)+   b1 * Y(z) * z^(-1)  + b2 * Y(z) * z^(-2)

. . .

Y(z)             a0
---- = ---------------------------- = H(z)
X(z)    1 - b1*z^(-1) - b2* z^(-2)


Bei Grenzfreq.  H(z) = 1/sqrt(2)

hilft das ein bisschen ?

von Reinhard M. (Gast)


Lesenswert?

Nachtrag:
ist ein IIR Filter, kein FIR, weil rekursiv
-- von vorangegangenen Ausgangswerten abhängig --

Sollte ein schmaler Bandpass sein,
mit den Eckfrequenzen 70Hz und 85Hz.

von Burkhard K. (buks)


Lesenswert?

Reinhard M. schrieb:
> Nachtrag:
> ist ein IIR Filter, kein FIR, weil rekursiv
> -- von vorangegangenen Ausgangswerten abhängig --

Wirklich?

Wenn nur a0 <> Null mit x(n) multipliziert wird gibt es kein Feedback 
und damit auch keinen "infinite response". (Normalerweise würde man a0 
dann eher mit b0 bezeichnen).

In der Schleife
1
   for i:=1 to n do
2
    y[0]:=y[0] + a[i] * x[i] - b[i] * y[i];
3
   filter:=y[0];
mit (nur a0 <> Null):
1
     pmax_filter : fdaten = (n: 1;
2
                   a:(0.2  ,      0,     0,   0,   0);
3
     ...

ist das Glied (a[i] * x[i]) = 0 und könnte somit auch weggelassen 
werden.

von W.S. (Gast)



Lesenswert?

Kollege von oben schrieb:
> ich habe die Koeffizienten und die Abtastrate.

Ach so.

Es wäre durchaus hilfreich, eben diese Koeffizienten mal zu posten.

Aus deinem seltsamen Stück Pascal kann ich nicht erkennen, was da 
eigentlich wie gefiltert werden soll.

Also, poste mal deinen Filterkernel in normaler Form, dann kann ich dir 
das mal simulieren, aber so wird das nix. Ich häng dir mal ne Textdatei 
und das Wobbelergebnis dran, wo du ein SSB-Filter sehen kannst. Kommt 
auf diesem PC bissel mies, da das Programm noch mit D6 gemacht wurde und 
unter Win7 nicht gut skaliert. Aber wie Filtertaps so aussehen, kannst 
du draus entnehmen - und wenn du den Zinnober in Excel o.ä. importierst, 
kannst du dir auch den Filterkernel mal grafisch angucken.

W.S.

von Reinhard M. (Gast)


Lesenswert?

Burkhard K. schrieb:
> Wirklich?
>
> Wenn nur a0 <> Null mit x(n) multipliziert wird gibt es kein Feedback
> und damit auch keinen "infinite response". (Normalerweise würde man a0
> dann eher mit b0 bezeichnen).

Ja, wirklich

y[0]:=a[0] * x[0];  // aktueller Ausgang = a0 * aktuellem Eingang
   for i:=1 to n do
    y[0]:=y[0] + a[i] * x[i] - b[i] * y[i];
//          ^          ^          ^
            |          |          |______ b[i] * alte y- werte
                       |                     REKURSIVTEIL (IIR) !
                       |
                       |_________ NICHT REKURSIVER (FIR) TEIL
                                  kann in diesem Fall weggelassen werden
                                  weil alle a[i] mit i>0  gleich null 
sind
   filter:=y[0];

von Burkhard K. (buks)


Lesenswert?

Du hast recht, hier gibt es (ausschliesslich) einen Feedback-Anteil - 
ich hab mich von den vertauschten Koeffizientenbezeichnern a vs. b 
verwirren lassen.

Worauf ich ursprünglich hinweisen wollte: eine rekursive Struktur allein 
ist noch kein hinreichendes Merkmal für einen IIR-Filter, es gibt 
durchaus FIR-Filter mit rekursiver Struktur, z.B. Moving Averager, CIC- 
oder Frequency-Sampling-Filter.

von Reinhard M. (Gast)


Lesenswert?

das C#-Teil mal mit Dirac erweitert

1
using System;
2
using System.Collections.Generic;
3
//using System.Linq;
4
using System.Text;
5
//using System.Threading.Tasks;
6
7
namespace FilterCsharp
8
{
9
    class Program
10
    {
11
        public struct filter
12
        {
13
            public int n;
14
            public float[] a;
15
            public float[] b;
16
            public float[] x;
17
            public float[] y;
18
        };
19
20
        static void Main(string[] args)
21
        {
22
            filter divFilter;
23
24
            divFilter.a = new float[4];
25
            divFilter.b = new float[4];
26
            divFilter.x = new float[4];
27
            divFilter.y = new float[4];
28
29
            divFilter.n = 1;
30
            divFilter.a[0] = 0.2F;
31
            divFilter.a[1] = 0;
32
            divFilter.a[2] = 0;
33
            divFilter.a[3] = 0;
34
35
            divFilter.b[0] = 1;
36
            divFilter.b[1] = -0.8F;
37
            divFilter.b[2] = 0F;
38
            divFilter.b[3] = 0;
39
40
            for ( int  samples = 0; samples < 16; samples++)
41
            {
42
                
43
                float returnValue;
44
                
45
                if(samples == 0)
46
                  returnValue = filterFunction(divFilter, 1);
47
                else
48
                  returnValue = filterFunction(divFilter, 0);
49
50
                Console.WriteLine("Filter Value {0}", returnValue);
51
            }
52
53
            Console.Write("Press any key to continue . . . ");
54
            Console.ReadKey(true);
55
        }
56
57
        public static float filterFunction(filter fil, float newSample)
58
        {
59
60
            int i = fil.n;
61
62
            for (i = fil.n; i >= 1; i--)
63
            {
64
                fil.x[i] = fil.x[i - 1];
65
                fil.y[i] = fil.y[i - 1];
66
            }
67
68
            fil.x[0] = newSample;
69
            fil.y[0] = fil.a[0] * fil.x[0];
70
71
            for (i = 1; i <= fil.n; i++)
72
            {
73
                fil.y[0] = fil.y[0] + fil.a[i] * fil.x[i] - fil.b[i] * fil.y[i];
74
            }
75
76
            return fil.y[0];
77
        }
78
    }
79
}


Ausgabe:
Filter Value 0,2
Filter Value 0,16
Filter Value 0,128
Filter Value 0,1024
Filter Value 0,08192001
Filter Value 0,06553601
Filter Value 0,05242881
Filter Value 0,04194305
Filter Value 0,03355444
Filter Value 0,02684355
Filter Value 0,02147484
Filter Value 0,01717987
Filter Value 0,0137439
Filter Value 0,01099512
Filter Value 0,008796096
Filter Value 0,007036877

von Reinhard M. (Gast)


Lesenswert?

nach fft
schaut doch nach TP aus

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.