filter berechnen

wenn es vorlagen dafür gäbe würde ich die nutzen. aber du hast recht, dass meine frage scheiße ist, ein wenig mehr sollte man da wohl mitliefern.

ich benötige die mathematik für chebychev-, bessel- und eine elliptische topologie, die jeweils einen bereits bestehenden biquad mit 5 coeffizienten befeuern können.

und das in einer form, die man auch lesen und verstehen kann, wenn man nicht mathematikprofessor ist.

das ist unüblich, das so zu machen, aber ich habe gründe.

butterworth hab ich schon lange am start, und dass man keine ladder filter mit biquad bauen kann muss ich wohl akzeptieren.

bei den anderen bekomme ich die "konversion" nicht geregelt (man findet immer nur code für FIR)
 
ich hatte mal meinen butterworth in eine allgemein lesbare form gebracht, das sieht dann so aus:

A0 ->
(1.0/(1.0 + sqrt(2.)*(1.0/tanh(pi*(F/SR))) + (1.0/tanh(pi*(F/SR))) * (1.0/tanh(pi*(F/SR)))))

A1 ->
(1.0/(1.0 + sqrt(2.)*(1.0/tanh(pi*(F/SR))) + (1.0/tanh(pi*(F/SR))) * (1.0/tanh(pi*(F/SR))))) * 2.

A2 ->
(1.0/(1.0 + sqrt(2.)*(1.0/tanh(pi*(F/SR))) + (1.0/tanh(pi*(F/SR))) * (1.0/tanh(pi*(F/SR)))))

B1 ->
2.0*(1.0/(1.0 + sqrt(2.)*(1.0/tanh(pi*(F/SR))) + (1.0/tanh(pi*(F/SR))) * (1.0/tanh(pi*(F/SR))))) * (1.0 - (1.0/tanh(pi*(F/SR)))*(1.0/tanh(pi*(F/SR))))

B2 ->
(1.0/(1.0 + sqrt(2.)*(1.0/tanh(pi*(F/SR))) + (1.0/tanh(pi*(F/SR))) * (1.0/tanh(pi*(F/SR))))) * (1.0 - sqrt(2.0) * (1.0/tanh(pi*(F/SR))) + (1.0/tanh(pi*(F/SR)))*(1.0/tanh(pi*(F/SR))))

mag aber außer mir niemand. :sad:
 
Keine Ahnung ob das 'ne Option ist: Erst vorstellen wie der Filter klingen soll und dann überlegen wie er zu programmieren ist? So als unkonventionellen Ansatz.
 
ich habe keinen plan von traditioneller höherer mathematik - ich kann die nur abstrakt anwenden - und vor allem hab ich wenig disziplin.

ich bin z.b., nachdem ich (mithilfe von beispielen) vor 15 jahren meinen butterworth gebaut habe, in die matlab examples eingestiegen und dachte ich kann mir durch eine synopse die anderen filter ganz einfach da rausziehen. (wie schwer kann es sein, wenn ich beispiel A als zwei versioenn haben, und beispiel B in version, beispiel B in version 2 zu basteln? - aber ich bin damit gescheitert.)

und offenbar ist eben doch mehr zu tun als nur eine z-transform. (?)

beim cheby kommt noch hinzu, dass man da wohl ein ganzes set von berechungen braucht, weil der ja oft erst in höherer order zur praktischen anwendung kommt. (vom grundsatz her würde ich schon gerne das ripple amount auch in echtzeit modulieren können - ja, über die probleme mit modulation beim biquad bin ich mir bewusst)

das smith buch hab ich übrigens als pdf. das beispiel ist ja in der tat sehr gut versteckt! da hätte man nicht mal gesucht, wenn man gesucht hätte. :)
 
Zuletzt bearbeitet:
Die Frage ist halt, ob du in deinem vorhandenen Framework einfach die Koeffizienten A0, A1, A2, B1, B2 eingeben kannst und dann wird das Filter korrekt berechnet. Davon war ich jetzt ausgegangen. In dem Fall wäre es wirklich nur Algebra, die Koeffizienten von der einen Form in die andere umzurechnen.

Alles Weitere ist natürlich mehr Arbeit, insbesondere Filterordnungen > 2. Mit IIR gibt es ja dann auch schnell Stabilitätsprobleme.

Grüsse,
synthos
 
ja, mein input is frequenz, samplingrate, und ggf Q usw, und mein output sind coeffizienten. davon bin ich dann in der lage mir zu bauen was ich brauch, nämlich die coeffizientenberechung auf 3 verschiedene arten (data, audio rate, audiorate mit vorberechneten buffern.)
leider sind die unterschiede im layout zwischen verschiedenen sprachen deutlich größer als die unterschiede zwischen verschiedenen filter topologien.

meine hoffung ist, dass jemand, der verschiedene filter wie auch immer schon selbst erstellt bzw berechnet hat, zu meinen obigen pseudo code einfach sagen kann "da und da sind der chauchy und der elliptische anders als der butter"... so extrem unterschiedlich können die ja nicht sein.
 
Zuletzt bearbeitet:
Ich hab mal das Listing 20-5 aus der PDF gezogen:

Code: 
1000 'THIS SUBROUTINE IS CALLED FROM TABLE 20-4,  LINE 340
1010 '
1020 '  Variables entering subroutine:  PI, FC, LH, PR, NP, P%
1030 '  Variables exiting subroutine:  A0, A1, A2, B1, B2
1040 '  Variables used internally:     RP, IP, ES, VX, KX, T, W, M, D, K,
1050 ' X0, X1, X2, Y1, Y2
1060 '
1070 ' 'Calculate the pole location on the unit circle
1080 RP = -COS(PI/(NP*2) + (P%-1) * PI/NP)
1090 IP =   SIN(PI/(NP*2) + (P%-1) * PI/NP)
1100 '
1110 ' 'Warp from a circle to an ellipse
1120 IF PR = 0 THEN GOTO 1210
1130 ES = SQR( (100 / (100-PR))^2 -1 )
1140 VX = (1/NP) * LOG( (1/ES) + SQR( (1/ES^2) + 1) )
1150 KX = (1/NP) * LOG( (1/ES) + SQR( (1/ES^2) - 1) )
1160 KX = (EXP(KX) + EXP(-KX))/2
1170 RP = RP * ( (EXP(VX) - EXP(-VX) ) /2 ) / KX
1180 IP  = IP * ( (EXP(VX) + EXP(-VX) ) /2 ) / KX
1190 '
1200 '            's-domain to z-domain conversion
1210 T  = 2 * TAN(1/2)
1220 W  = 2*PI*FC
1230 M  = RP^2 + IP^2
1240 D = 4 - 4*RP*T + M*T^2
1250 X0 = T^2/D
1260 X1 = 2*T^2/D
1270 X2 = T^2/D
1280 Y1 = (8 - 2*M*T^2)/D
1290 Y2 = (-4 - 4*RP*T - M*T^2)/D
1300 '
1310 ' 'LP TO LP, or LP TO HP transform
1320 IF LH = 1 THEN K = -COS(W/2 + 1/2) / COS(W/2 - 1/2)
1330 IF LH = 0 THEN K =  SIN(1/2 - W/2) / SIN(1/2 + W/2)
1340 D = 1 + Y1*K - Y2*K^2
1350 A0 = (X0 - X1*K + X2*K^2)/D
1360 A1 = (-2*X0*K + X1 + X1*K^2 - 2*X2*K)/D
1370 A2 = (X0*K^2 - X1*K + X2)/D
1380 B1 = (2*K + Y1 + Y1*K^2 - 2*Y2*K)/D
1390 B2 = (-(K^2) - Y1*K + Y2)/D
1400 IF LH = 1 THEN A1 = -A1
1410 IF LH = 1 THEN B1 = -B1
1420 '
1430 RETURN


Dieses ganze Listing ist eigentlich nur eine Funktion, die aber mehrfach und rekursiv vom Hauptprogramm (Listing 20-4) aufgerufen wird. Da Basic aber eine rein prozedurale Programmiersprache ist und mit GOTO-Sprungmarken arbeitet, ist man also einfach auf "Abstand" gegangen.
Aufgerufen wird diese Funktion in Zeile 340.
Beendet wird diese Funktion mit RETURN in Zeile 1430 und das Programm fährt mit der Schleife fort.

Es geht los mit der Zeile 1020:
Das ist eigentlich nur ein Kommentar, der aber sagt, welche Variablen an diese Funktion "übergeben" werden. Benutzereingaben hab ich gelb markiert. Im Einzelnen sind dies:

PI = einfach nur die Zahl "Pi".
FC = Cutoff-Frequenz (gültig zwischen 0 und 0,5)
LH = Lowpass- oder Highpass-Filter? (nur 0 und 1 gültig)
PR = Prozentsatz an "Ripple" (gültig 0 bis 29)
NP = Anzahl der Pole (2 bis 20, ACHTUNG: nur gerade Zahlen gültig!) [Bemerkung: im Listing 20-5 hat sich ein Tippfehler eingeschlichen und es stand dort HP statt NP, was keinen Sinn ergibt]
P% = Dies ist der konkrete Schleifendurchlauf, der von 1 bis zur Hälfte der Pol-Anzahl (deswegen muss diese Eingabe auch zwingend gerade sein!) durchrattert. Sooft wird diese ganze Subroutine von Zeile 1000 bis 1430 durchgerattert!

Dann sind zunächst die Zeilen 1070 und 1080 wichtig. Da wird eine komplexe Zahl (besteht aus RP = Realteil und IR = Imaginärteil) im Einheitskreis berechnet. ACHTUNG: Die SIN- und COS-Funktionen werden mit Bogenmaß und nicht mit Grad gefüttert!

Die Rechengrundlage hinter diesen beiden Zeilen dazu:


https://www.youtube.com/watch?v=Vhz1HnU8QSI


In Zeile 1120 wird nun geprüft, ob der Benutzer einen Ripple-Prozentsatz von genau 0 angegeben hat. Ist das der Fall, dann wird der Einheitkreis NICHT zu einer Ellipse umgeformt und das Programm springt direkt zur Zeile 1210.
Wenn aber RP > 0 ist, dann WERDEN die Zeilen 1130 bis 1180 ausgeführt und darin die Ellipse gebildet und der Realteil und der Imaginärteil neu kalkuliert.


https://www.youtube.com/watch?v=FdOQK34q2Jk



In den Zeilen 1210 bis 1290 findet die Z-Transformation statt.


https://www.youtube.com/watch?v=j0KsNoHfU8I



https://www.youtube.com/watch?v=Ii3WqRV2b1k


https://www.youtube.com/watch?v=JGviDHZQ0iU


Dann wird in 1320 bzw. 1330 die Filtercharakteristik K auf Lowpass oder Highpass festgelegt. Anschließend werden die 5 Koeffizienten A0, A1, A2, B1, und B2 berechnet. Wenn ein Hochpassfilter ausgewählt wurde, dann werden in den Zeilen 1400 und 1410 die Werte für A1 und B1 invertiert.
 
hm, die grundstruktur in diesem unteren teil erkenne ich durchaus wieder. :)

ich habs mir auch grad mal ohne zeilennummern kopiert und alles durchgelesen.

als nächstes benenne ich die variablen 110-style um, damit ich das besser lesen kann.

spontan kapieren tue ich:

- die radians konversion am eingang (ja ich weiß, das heißt anders)
- die variablen und die verknüpfungen
- die z-transform
- die tiefpass zu tiefpass transform

denn das ist erst mal genau so wie bei anderen filtertypen.

die ellipse - und die trennung des codes in HP und LP sollten auch ein nobrainer sein.


der teufel steckt dann in den details:

z.b. die aktion mit den ripples.

zum einen, ist es sicherlich sinnvoll einen subcircuit draus zu machen, weil es ja bei 0 0 ergibt - da tendiere ich zum auschalten-wollen-wenns-0-ist.

zum anderen bin ich wegen der range überfragt: warum geht das eigentlich genau bis 22%? (int? oder sollte das nicht vielmehr (1 minus 0.707...) wissenschon sein?)


unklar (im sinne von wissenslücke) ist mir auch die sache mit den poles.

ich will (32bit!) und muss (der biquad ist nicht von mir) hier kaskadieren.

warum ist das in dem basic code nicht zu finden? da muss doch der Q angepasst werden usw? oder handelt man es einfach nur so, wie man es bei anderen filtern auch machen würde.
 
Zuletzt bearbeitet:
so, dann machen wir das jetzt mal halbwegs lesbar.


unit circle


RP = -COS(PI/(poles*2) + (ripple-1) * PI/poles)

IP = SIN(PI/(poles*2) + (ripple-1) * PI/poles)


ellipse transform


IF ripples = 0 THEN GOTO z-transform


ES = SQR( (100 / (100-PR))*(100 / (100-PR)) -1 )

VX = (1./NP) * LOG( (1/ES) + SQR( (1/ES*ES) + 1.) )

KX = (1./NP) * LOG( (1/ES) + SQR( (1/ES*ES) - 1.) )

KX = (EXP(KX) + EXP(-KX)) / 2.

RP = RP * ( (EXP(VX) - EXP(-VX) ) /2. ) / KX

IP = IP * ( (EXP(VX) + EXP(-VX) ) /2. ) / KX


z-domain transform


T = 2 * TAN(1/2)

W = 2*PI*F

M = RP*RP + IP*IP

D = 4 - 4*RP*T + M*T*T

X0 = T*T / D

X1 = 2*T*T / D

X2 = T*T / D

Y1 = (8 - 2*M*T*T) / D

Y2 = (-4 - 4*RP*T - M*T*T) / D


lp to lp transform


highpass:
K = -COS(W/2 + 1/2) / COS(W/2 - 1/2)

lowpass:
K = SIN(1/2 - W/2) / SIN(1/2 + W/2)


D = 1 + Y1*K - Y2*K^K


A0 = (X0 - X1*K + X2*K*K) / D

A1 = (-2*X0*K + X1 + X1*K*K - 2*X2*K) / D

A2 = (X0*K*K - X1*K + X2) / D

B1 = (2*K + Y1 + Y1*K*K - 2*Y2*K) / D

B2 = (-(K*K) - Y1*K + Y2) / D


highpass:
A1 = -A1
B1 = -B1


das hier checke ich noch nicht.

RP = RP * ( (EXP( ...

IP = IP * ( (EXP( ...

wie soll man eine variable anhand sich selbst errechnen? oder ist damit sample t-1 gemeint oder was?

und wozu berechnen wir eine variable KX wenn die später nirgends benutzt wird?

vorher brauch ich gar nicht weitermachen.
 
Zuletzt bearbeitet:
einseinsnull schrieb:
das hier checke ich noch nicht.

RP = RP * ( (EXP( ...

IP = IP * ( (EXP( ...

wie soll man eine variable anhand sich selbst errechnen? oder ist damit sample t-1 gemeint oder was?
Zum Vergrößern anklicken....

Der Gleichzeichen-Operator ist in BASIC eine sogenannte Zuweisung. Das ist NICHT mit einem arithmetischen Gleichheitszeichen zu verwechseln!

In BASIC wird zunächst der Teil RECHTS vom Gleichzeichen-Operator berechnet und das Ergebnis dann der Variable LINKS vom Gleichzeichen-Operator zugewiesen.

Deswegen wäre eine solche Reihenfolge auch gültig:

A = 5
B = 2
C = A + B

Der Inhalt der Variable "C" ist jetzt "7".

Weiter:

C = C - (A - B)

Der Inhalt der Variable "C" ist jetzt "4"! Klammern werden von innen nach außen hin zuerst ausgewertet und danach dann der Rest rechts vom Gleichzeichen-Operator.
 
2bit schrieb:
In BASIC wird zunächst der Teil RECHTS vom Gleichzeichen-Operator berechnet und das Ergebnis dann der Variable LINKS vom Gleichzeichen-Operator zugewiesen.
Zum Vergrößern anklicken....

ja, aber wie denn, wenn man nicht weiss, was RP ist? ich steh da echt auf dem schlauch.

da steht ja nur:

1. wir führen eine interne variable B ein

2. B = B + 5

???

dass die reihenfolge nicht von oben nach unten sein muss hab ich an anderen stellen schon gemerkt. das ist okay. :) die anderen stellen bekomme ich alle zusammen. da sehe ich auf einen blick, wie ich das in meine welt übertragen kann.
 
Zuletzt bearbeitet:
jetzt häng ich wieder.

M = RP*RP + IP*IP

welches real und imaginary ist das jetzt? das obere (zeile 1080)? oder das erneut umgerechnete (zeile 1170)?
 
einseinsnull schrieb:
welches real und imaginary ist das jetzt? das obere (zeile 1080)? oder das erneut umgerechnete (zeile 1170)?
Zum Vergrößern anklicken....

Beides. :D

Zeile 1170 wird dann ausgeführt, wenn die Bedingung in Zeile 1120 NICHT erfüllt ist, also wenn PR (= Prozentsatz an "Ripple") ungleich 0 ist. Ansonsten wird der ganze Block einfach per "GOTO 1210" übersprungen.
 
beides! gute antwort.

jetzt hat der kleine 110 auch verstanden, warum der buchmensch da die gleiche variable 2 mal benutzt.

es war also doch nicht nur um mich zu verwirren.
 
einseinsnull schrieb:
hehe, erst macht er W = 2*pi*F und später macht er dann überall W/2., was für ein noob. :lollo:
Zum Vergrößern anklicken....

Das ist für die bessere Lesbarkeit, ist ja schließlich Beispiel-Code.
2*Pi*f ist die Kreisfrequenz. Den "noobs" hilft es manchmal beim Verstehen des Codes, wenn dieser Wert wiedererkennbar auftaucht.
Das W ist auch nicht willkürlich als Variablenname gewählt: in der Literatur wird die Kreisfrequenz in der Regel mit einem kleinen Omega abgekürzt.
 
Literaturtip: Tieze-Schenk (Halbleiter-Schaltungstechnik), da sind die ganzen Formeln drin und kurz aber gut erklärt, und auch einige Tabellen mit Koeffizienten für die verschiedenen Typen. Zu IIR und FIR sind auch die Basics drin sowie zur Z-Trafo.
Das kann man sich gut unters Kopfkissen legen, wenn man Nackenprobleme hat wegen zu viel Bildschirmarbeit, ganzschön dick.
Ansonsten gibt es für die Berechnung der Koeffizienten Software, das macht heute kaum einer von Hand, zumal es da einige Dinge zu optimieren gibt beim Design, die man mit Bleistift und Papier nicht mal schnell berechnet. Matlab ist dafür recht gängig. Matlab hat eine eigene Toolbox für das Filterdesign und andere Toolboxen, um mit dem optimierten Filterdesign automatisch Code für verschiedene Zielplattformen zu generieren.

http://www.winfilter.20m.com/ ist auch ganz gut und recht leicht für den Einstieg, wenn man einigermaßen weiss was man macht und die Koeffizienten braucht.
 
Gibt es 'ne Möglichkeit das Prinzip das hinter den Filter Algorithmen steckt zu erklären, was mit den Daten passiert, ohne dass man sich mit Formel belasten muss oder wird da wirklich versucht irgendwelche Schaltungen bzw. Bauteile in Formeln abzubilden?
 
Zuletzt bearbeitet:
Duch einen Tiefpassfilter sollen schnelle Änderungen gedämpft werden. Ein Tiefpassfilter erhält man z.B. durch die Bildung des Mittelwertes. Also z.B. der Mittelwert aus dem vorigen Wert und dem neuesten Wert (nennt man wenn man das kontinuierlich macht moving average).
Dieser Mittelwertfilter ist allerdings nicht besonders steil, er sieht aus wie ein Teil von einem Sinus im Frequenzbereich bzw. eine sinc-Funktion.
Wenn man es besser machen will muss man schnelle Änderungen stärker unterdrücken, langsame Änderungen hingegen weniger. Dazu muss man geschickt werte aus der Vergangenheit und den aktuellen Wert miteinander verwurschteln (mit Koeffizienten multiplizieren und addieren). Divisionen versucht man typischerweise zu vermeiden, weil es vergleichsweise rechenintensiv bzw. aufwändig ist, ganz kann man es aber meist nicht vermeiden.

Einen Hochpassfilter erreicht man im Prinzip, indem man das Originalsignal nimmt, und das Tiefpass-Signal (das von dem Originalsignal abgeleitet wurde) davon phasenrichtig abzieht.

Tieze-Schenk gibt es in Auszügen als google-buch
books.google.de

Halbleiter-Schaltungstechnik

Nicht ohne Grund ist der "Tietze/Schenk" seit langem das Standard-Lehrbuch und -Nachschlagewerk der Halbleiter-Schaltungstechnik. Das Werk wurde seit dem Erscheinen der ersten Auflage im Jahre 1969 ständig überarbeitet und bietet daher nicht nur Einführung und Überblick über die heutige...
books.google.de
 
Zuletzt bearbeitet:
haesslich schrieb:
Das ist für die bessere Lesbarkeit, ist ja schließlich Beispiel-Code.
Zum Vergrößern anklicken....

is schon klar - bei dem ganzen gewusel kommts eh nicht mehr auf "optimierung" an. das ist ja echt drei mal mehr als bei gescheiten filtern.^^

Summa schrieb:
Gibt es 'ne Möglichkeit das Prinzip das hinter den Filter Algorithmen steckt zu erklären, was mit den Daten passiert, ohne dass man sich mit Formel belasten muss oder wird da wirklich versucht irgendwelche Schaltungen bzw. Bauteile in Formeln abzubilden?
Zum Vergrößern anklicken....

nicht wirklich. nur ein mathematiker würde behaupten, dass eine elektronische schaltung eigentlich das gleiche macht wie eine reihe von assembler befehlen oder ein max patch.

ansonsten sprechen programmierer und elektroniker schon sehr unterschiedliche sprachen. vor allem ist digital das spektrum begrenzt.

wie corc richtig sagt macht man sowas heutzutage normalerweise auch nicht selbst, denn das ist alles z.b. in matlab vorhanden. und übrigens auch in max, supercollider, reaktor, in diversen c++ libraries usw, da müsste man theoretisch nur noch benutzen, was schon da ist.

es gibt aber immer auch gründe, sowas auch mal selbst zu bauen. (vertrauen ist gut - kontrolle ist besser! - hab genug schwachsinn im source code von max und cubase entdeckt die letzten 20 jahre)

mein problem damit ist hier, dass ich bislang keine guten beispiele für einen cheby als IIR gefunden hatte, die ich in meine welt übersetzen kann. dieses übersetzen ist das, was die arbeit macht.
 
Zuletzt bearbeitet:
C0r€ schrieb:
Duch einen Tiefpassfilter sollen schnelle Änderungen gedämpft werden. Ein Tiefpassfilter erhält man z.B. durch die Bildung des Mittelwertes. Also z.B. der Mittelwert aus dem vorigen Wert und dem neuesten Wert (nennt man wenn man das kontinuierlich macht moving average).
Dieser Mittelwertfilter ist allerdings nicht besonders steil, er sieht aus wie ein Teil von einem Sinus im Frequenzbereich bzw. eine sinc-Funktion.
Wenn man es besser machen will muss man schnelle Änderungen stärker unterdrücken, langsame Änderungen hingegen weniger. Dazu muss man geschickt werte aus der Vergangenheit und den aktuellen Wert miteinander verwurschteln (mit Koeffizienten multiplizieren und addieren). Divisionen versucht man typischerweise zu vermeiden, weil es vergleichsweise rechenintensiv bzw. aufwändig ist, ganz kann man es aber meist nicht vermeiden.
Zum Vergrößern anklicken....
In wie fern spielt die Filterfrequenz dabei eine Rolle?

C0r€ schrieb:
Einen Hochpassfilter erreicht man im Prinzip, indem man das Originalsignal nimmt, und das Tiefpass-Signal (das von dem Originalsignal abgeleitet wurde) davon phasenrichtig abzieht.
Zum Vergrößern anklicken....
Den Trick kenne ich von div. Samplern, das Sample invertieren und überlagern ;-)
 
Du musst dich einloggen oder registrieren, um hier zu antworten.

Similar threads

Tommi
Pure Data Experiment: Berechnung Intensität Obertöne
Antworten
12
Aufrufe
360
Tommi
Tommi
Koch
Arduino als Solina, welches Filter ist das?
Antworten
1
Aufrufe
229
fanwander
fanwander
Edobot
ROLAND D 50 Filter Mapping
Antworten
13
Aufrufe
341
Edobot
Edobot
goorooj
Midi & Daw: Sub37 nimmt erst auf wenn ich am Filter drehe...
Antworten
12
Aufrufe
322
goorooj
goorooj
Nick Name
Für die Bearbeitung von komplexen Signalen mit Filtern - Gibt es da Unterschiede zu beachten!?
Antworten
0
Aufrufe
160
Nick Name
Nick Name


Neueste Beiträge

News

Zurück
Oben