FIR aus Impulsantwort generieren, um Frequenzgang zu begradigen

[Joe F. (easylife)]

Ich habe eine gemessene Impulsantwort (blau) eines Lautsprechers und 
möchte jetzt einen FIR Filter (grün) generieren, um den Frequenzgang zu 
begradigen.
Dazu würde ich die gemessene Impulsantwort durch eine FFT jagen, und 
müsste dann das Spektrum auf irgendeine Weise "horizontal spiegeln", um 
dann damit das Eingangssignal falten zu können.
Was muss man nach der FFT der gemessenen Impulsantwort mit den realen 
und imaginären Teilen tun, um diese Spiegelung zu erreichen?

Ich möchte das ganze in C umsetzen, also kein Mathlab o.Ä...

[Dergute W. (derguteweka)]

Moin,

Joe F. schrieb:
> Was muss man nach der FFT der gemessenen Impulsantwort mit den realen
> und imaginären Teilen tun, um diese Spiegelung zu erreichen?

Wenn deine Lautsprecherimpulsantwort h1(t) ist; dann ist die 
Fouriertransformierte z.B. H1(s). Also brauchst du noch eine Funktion 
z.B. H2(s), so das gilt:

H1(s)*H2(s)=exp(-sT) //(das * ist eine Multiplikaton, keine Faltung)

Da steht deshalb exp(-gedoens) und nicht 1, weil durch jedes Filter das 
Signal ja verzoegert wird. Zumindest, wenn die Filter kausal sein 
sollen, was immer eine gute Idee ist, wenn man's in echt nachbauen will.
Also hast du dann idealerweise ein Filter H1, den Lautsprecher und ein 
Filter H2, das Korrekturfilter, und beide zusammen machen nur eine 
Verzoegerung des Signals, keine weiteren Amplituden oder 
Phasenveraenderungen.
Kannst ja mal versuchen das T=0 zu setzen, damit vereinfacht sich das zu
H1(s)*H2(s)=1;
H2(s)=1/H1(s)
vielleicht haste Glueck und wenn du dein H2 wieder 
zuruecktransformierst, kriegste ein paar Koeffizienten, die nicht ganz 
daneben liegen und das Filter hat dann die fft-Laenge an Verzoegerung. 
Bin mir aber nicht sicher, ob das immer so funktioniert.


Wuerdest du das nicht in C, sondern in Matlab oder Octave machen wollen, 
wuerde ich da die Funktionen invfreqz() und deconv() empfehlen.

Gruss
WK

[Joe F. (easylife)]

> H2(s)=1/H1(s)

genau da hakt's bei mir.
Heisst das H2(s).re = 1/H1(s).re und H2(s).im = 1/H1(s).im?
Mal abgesehen vom Fall dass H1(s) Null sein kann, wenn H1(s) bei mir im 
Bereich von 0 bis 1 liegt, wird H2(s) ja beliebig gross...

[Manfred M. (bittbeisser)]

So wie ich das verstanden habe, muss man dazu eine inverse FFT (oder 
inverse DFT wenn N keine Zweierpotenz) mit dem abgetasteten 
Wunsch-Frequenzgang durchführen (also mit der Spiegelung). Details dazu 
hier: http://www.dspguide.com/ch17/1.htm

 > Ich möchte das ganze in C umsetzen, also kein Mathlab o.Ä...
Bezieht sich das jetzt auf die Ermittlung der Koeffizienten oder nur auf 
den Code für das Filter?

Der Code für die IFFT bzw. IDFT ist im verlinkten Buch zwar nur in einem 
rudimentären BASIC Dialekt vorhanden, aber das sollte sich problemlos in 
C übersetzen lassen.

[Dr. Sommer (Gast)]

Versuch's mal mit dem LMS Algorithmus und dem Stichwort "Channel 
Equalization". Gut dafür ist das Buch "Fundamentals of adaptive 
Filtering". Wenn du das so machst kann sich das System automatisch an 
verschiedene Lautsprecher anpassen. Ob du das in C oder Matlab 
implementierst ist letztlich egal. Es werden im Endeffekt nur 
Vektor-Multiplikationen und Grundrechenarten gebraucht

[Joe F. (easylife)]

Manfred M. schrieb:
> Bezieht sich das jetzt auf die Ermittlung der Koeffizienten oder nur auf
> den Code für das Filter?

Der Code ist kein Problem. Mir geht es darum zu verstehen, wie ich aus 
der gemessenen Impulsantwort h1 -> FFT -> H1 den passenden 
(kompensierenden) Filter H2 konstruiere.

Und ich denke, ich habe bei 1/H(s) einfach nur zu kurz gedacht.
Habe gerade das hier gefunden:

http://www.mathe-online.at/materialien/Andreas.Pester/files/ComNum/images/Eqn4.gif

und werde morgen dementsprechend mal mit

a = H1(s).re * H1(s).re + H1(s).im * H1(s).im

H2(s).re = H1(s).re / a
und
H2(s).im = - H1(s).im / a

probieren.

[Dergute W. (derguteweka)]

Moin,

Joe F. schrieb:
> Mal abgesehen vom Fall dass H1(s) Null sein kann, wenn H1(s) bei mir im
> Bereich von 0 bis 1 liegt, wird H2(s) ja beliebig gross...

Jepp. Kann passieren. Dann haste bisschen verkackt. Is ja auch irgendwie 
nachvollziehbar: Wenn du ein Filter hast, was irgendeine Frequenz 
komplett sprerrt, ist die ja danach nicht mehr im Signal und du weisst 
ums verrecken nicht mehr ob und mit welchem Pegel/Phasenlage sie vorher 
drinnen war...

Joe F. schrieb:
> und werde morgen dementsprechend mal mit
> ...
> probieren.

Jepp. Halt komplexe Rechnung. In C dann eben schulbuchmaessig mit 
konjungiert komplexer Erweiterung etc. bla.
Ich schaetz' aber mal, dass das auch nicht immer hinhauen wird, weil ja 
beide Filter zusammen eben nicht mehr 1 als Gesamtfrequenzgang haben 
koennen wegen der Kausalitaet, sondern irgendeine Verzoegerung haben.
Und je groesser diese Verzoegerung sein darf, desto "besser" wird sich 
das Korrekturfilter berechnen lassen.

Gruss
WK

[Oh (Gast)]

Musik rein und nach ner halben Stund Musik aus den Lautsprechern wieder 
raus? :(

[Dergute W. (derguteweka)]

Moin,

Oh schrieb:
> Musik rein und nach ner halben Stund Musik aus den Lautsprechern
> wieder
> raus? :(

Ja. Aber dafuer mit 1A Qualitaet aus'm 1000W P.M.P.O 
Plastik-Bruellwuerfel aus dem 1 EUR Shop. :)

Gruss
WK

[Joe F. (easylife)]

Dergute W. schrieb:
> Ich schaetz' aber mal, dass das auch nicht immer hinhauen wird, weil ja
> beide Filter zusammen eben nicht mehr 1 als Gesamtfrequenzgang haben
> koennen wegen der Kausalitaet, sondern irgendeine Verzoegerung haben.

Klar. Mir reichen nach unten ca. 50-60 Hz.

Dergute W. schrieb:
> Ja. Aber dafuer mit 1A Qualitaet aus'm 1000W P.M.P.O
> Plastik-Bruellwuerfel aus dem 1 EUR Shop. :)

Ich will es nicht 100% gerade machen. Der Effekt die größeren Dellen 
etwas auszugleichen und den Bass (80 Hz) um 6dB anzuheben ist aber nicht 
zu unterschätzen. Erste Experimente mit IIR Filtern haben schon recht 
überzeugende Hör-Ergebnisse gebracht. Die "Box" hat eher 5W und ca. 1.5L 
Volumen...

[Tobias P. (hubertus)]

Du könntest auch quasi offline einen LMS, NLMS oder RLS Filter drüber 
laufen lassen. Das geht auch mit einem kleinen C-Programm. Danach die 
Koeffizienten vom Filter abgreifen und in ein 'festes' Filter eingeben.

Oder eine interessante Alternative wär gleich ein adaptiver Lautsprecher 
:-)

[Joe F. (easylife)]

Vermelde Erfolg.
Man muss die komplexen Zahlen halt korrekt dividieren (ach)...
schwarz: Messung, rot: gemittelter Freqenzverlauf, blau: errechneter 
Filter
Anhören konnte ich es noch nicht, aber so falsch sieht's nicht aus.

[Dogbert (Gast)]

Glückwunsch zur ersten Dekonvolution.

Praktische Implemenationen gibt es open source:

http://drc-fir.sourceforge.net/

[Edi M. (Gast)]

Was war da jetzt das Problem? Im logarithmischen ist das doch nur eine 
Spiegelung des Frequenzganges an der 0dB-Linie? Etwas Skalierung 
abgerechnet ...

[Joe F. (easylife)]

Edi M. schrieb:
> Was war da jetzt das Problem?

Das Problem war die Spiegelung im Frequenzbereich nach der FFT korrekt 
zu berechnen. Das Spektrum liegt dann ja als Array von komplexen Werten 
vor, wenn man dort nur die Amplitude ändern würde, wird's nicht 
richtig... ;-)
Allerdings habe ich mich inzwischen gegen FIR und für eine Kombination 
aus 8 wohl balancierten IIRs entschieden. Gründe waren bessere 
Beherrschbarkeit der Korrektur, geringere Welligkeit der Filter, 
geringere Latenz des Systems.

[Edi M. (Gast)]

Joe F. schrieb:
> Das Spektrum liegt dann ja als Array von komplexen Werten
> vor, wenn man dort nur die Amplitude ändern würde, wird's nicht
> richtig... ;-)

hm, das ist wohl in der Tat ein Punkt, den Ich übersehen zu haben 
scheine. Dann müsste mn komplexe invertieren? Vom Betrag müsste es aber 
genau so aussehen, denke Ich und vom Winkel einfach negiert? Stimmt das 
?

[Joe F. (easylife)]

Edi. M. schrieb:
> Stimmt das ?

Nein, dann würdest du ja nur die Phase drehen.

x * 1/x

ist da schon richtig, um überall auf '1' zu kommen.