Schnelle Faltung zur Filterung von kontinuierlichen Daten

Und das schaffst Du nicht unter iOS v47 in Java auf einem 5GHz i7 ?

Will fragen: in welcher Umgebung ist diese Aufgabenstellung eine 
Herausforderung?

Wie gross ist den die Ordnung (n?) Deines FIR Filters?

Die schnelle Faltung mittels FFT lohnt sich normalerweise erst wenn n 
ziemlich gross ist. Und dann musst Du ziemlich genau rechnen. Wenn Du 
nen FPGA mit genügend DSP-Blöcken (Multiply-Add MAC) hast, kannst Du ein 
FIR-Filter durch die Parallelisierung sehr schnell berechnen.

Microcontroller-Bastler schrieb:
>> Die schnelle Faltung mittels FFT lohnt sich normalerweise erst wenn n
>> ziemlich gross ist.
> Wären 128 schon "ziemlich gross"?

Wenn ich gerade nicht total Müll gerechnet habe lohnt sich das schon 
relativ schnell. Für eine FFT braucht man laut Wikipedia log2(N)*N 
Multiplikationen und Additionen für einen N Werte großen Vektor.

Für die schnelle Faltung braucht man eine FFT, dann N Multiplikationen 
und die inverse FFT. Macht grob 2*log2(N)*N+N Rechenschritte. Die 
Faltung direkt sind etwa 2*N^2 Rechenschritte. Siehe Plot:

https://www.wolframalpha.com/input/?i=plot+2*log2(N)*N%2BN,+2*N%5E2+from+N%3D1+to+10

Gut die Multiplikationen bei der Variante über die FFT sind komplex aber 
alles in allem sollte spätestens ab 8 oder 16 Vektorlänge die schnelle 
Faltung wirklich schneller sein.

Microcontroller-Bastler schrieb:
> Marcus H. schrieb:
>> Und das schaffst Du nicht unter iOS v47 in Java auf einem 5GHz i7 ?
>
>> Will fragen: in welcher Umgebung ist diese Aufgabenstellung
>> eine Herausforderung?
>
> Es ist aus der Sicht ein Problem, als dass Ich nicht verstanden habe,
> wie die schnelle Faltung genau arbeitet und dahingehend nicht in der
> Lage bin, abzuschätzen, ob es sich lohnt, das anzusehen und anzuwenden.
>
> Danke aber für Deinen Versuch, zu helfen.

Es war tatsächlich ein Versuch, rauszubekommen, was Du willst.
Es hat schlicht Kontext gefehlt. Mit der aktuellen Info verstehe ich die 
Frage so:

- Ihr seid gerade dabei eine passende MCU für die Aufgabe zu suchen
- Ihr braucht jemanden, der Euch die Anzahl MUL/FMUL, OPS/FLOPS sowie 
Speicherbedarf(!) anhand der Filterlänge ausrechnet

Ich bin mir immer noch nicht so sicher, ob die Filterlänge dynamisch 
sein soll und welches Zahlenformat (Ganzzahl, Festkomma, Fließkomma) 
verwendet werden soll.

Mein Vorschlag, da Ihr anscheinend flexibel seid: keinen 8-Bitter.
Ein STM32F0 könnte eng werden, ein STM32F4/F7 bietet Platz zum 
experimentieren.

Moin,

Das wird hier zwar des oefteren geschrieben:


Sebastian V. schrieb:
> Wenn ich gerade nicht total Müll gerechnet habe lohnt sich das schon
> relativ schnell. Für eine FFT braucht man laut Wikipedia log2(N)*N
> Multiplikationen und Additionen für einen N Werte großen Vektor.
>
> Für die schnelle Faltung braucht man eine FFT, dann N Multiplikationen
> und die inverse FFT. Macht grob 2*log2(N)*N+N Rechenschritte. Die
> Faltung direkt sind etwa 2*N^2 Rechenschritte. Siehe Plot:
>
> 
https://www.wolframalpha.com/input/?i=plot+2*log2(N)*N%2BN,+2*N%5E2+from+N%3D1+to+10
>
> Gut die Multiplikationen bei der Variante über die FFT sind komplex aber
> alles in allem sollte spätestens ab 8 oder 16 Vektorlänge die schnelle
> Faltung wirklich schneller sein.

aber ich hab' noch nie ein konkretes Beispiel gesehen, wo das 
funktioniert und tatsaechlich Zyklen spart.

Wie schon angedeutet: Dann sind erstmal alle Multiplikationen im 
Frequenzbereich ploetzlich komplex; d.h. mehr als den 4fachen 
Rechenaufwand einer reellen Multiplikation.
Dann gibts natuerlich auch so Gimmicks wie: Wenn du mit z.b. 64k grossen 
Bloecken arbeitest fuer die FFT; dann ist bei einem normalen FIR z.B. 
5ter Ordnung im Ausgangswert aus Sample # 65537 noch die gesamte 
Vorgeschichte der Samples 65536,65535,65534, 65533,65532 enthalten. Bei 
der Version mit den 64k grossen FFT Bloecken aber nicht, denn da ging ja 
die neue FFT bei Sample 65536 los...

Mein Eindruck ist ganz stark, dass das mit der schnellen Falterei 
hoechstens fuer so blockorientierte Dinger, wie z.B. Bilder 
funktioniert, aber nicht vernuenftig fuer kontinuierliche Signale wie 
z.B. Audio oder irgendwelches HF-Zeugs, was keinen definierten Anfang 
und Ende hat.

Aber ich lass' mich gerne mit konkreten Beispielen eines Besseren 
belehren.

Ansonsten klingt ein Prozessor mit 16Mhz Taktfrequenz doch voellig im 
Bereich des Machbaren, auch ohne irgendwelche 
"FFT-hin-und-zurueck-und-trotzdem-zeit-spar-tricks"

Gruss
WK

Dergute W. schrieb:
> aber ich hab' noch nie ein konkretes Beispiel gesehen, wo das
> funktioniert und tatsaechlich Zyklen spart.

Ich hab das mal für eine Auswertung genutzt. Dort war es aber eine 
3-dimensionale Faltung mit Kantenlänge N=256. Wenn man die Faltung naiv 
macht sind das N^6 = 2,81E14 Multiplikationen. Das hätte Ewigkeiten 
gedauert. Über die FFT nur ~1 Sekunde auf einer Grafikkarte.

Edit: Und es war auch nichts kontinuierliches. Einfach einmal den 256^3 
Block falten um die räumliche Korrelation rauszukriegen.

Ich glaub bei der schnellen Faltung braucht man eine längere FFT.
Mal angenommen wir falten zwei reelle Sequenzen mit N Punkten.
Dann ist das Ergebnis 2N Punkte lang. Evtl. wirds wegen der Symmetrie er 
FFT sogar noch schlimmer.

Hat jemand nen Literaturhinweis wo jemand schon mal ne schnelle Faltung 
auf nem FPGA mit ner Implementation mit Hilfe von DSP-MAC Blöcken auf 
demselben FPGA verglichen hat?

Martin O. schrieb:
> Hat jemand nen Literaturhinweis wo jemand schon mal ne schnelle Faltung
> auf nem FPGA mit ner Implementation mit Hilfe von DSP-MAC Blöcken auf
> demselben FPGA verglichen hat?

Wie ist das gemeint? Die DSP-Blöcke sind die Multiplier und die würden 
bei einer FPGA-Lösung in jedem Fall benutzt, egal ob Faltung oder 
flacher FIR. Du kannst Dich eigentlich nur entscheiden, was parallel und 
was sequenziell läuft und die FFT ist im Vergleich zu einer DFT schon 
sehr stark parallelisiert, nutzt also bereits ausgeführte 
Zwischenergebnisse mehrfach.

Ich sehe nicht, dass man bei der Faltung etwa Multiplikationen sparen 
kann, wenigstens nicht aus Implementierungssicht.

So, ich hab jetzt mal ne "schnelle Faltung" in Java implementiert, und 
nun zähle ich die MAC (Multiply Add) Operationen.

Die Impulsantwort des Filters hat Nh Samples (z.B. Nh=64).
Die lange Eingangssequenz hat M Samples (M sehr gross, z.B. M=65536).

Man faltet jetzt immer Nh Samples der Eingangsfolge mit der 
Impulsantwort.
Das sind M/Nh Blöcke, bei jedem Block hat man ne Nh * Nh Faltung mit 
einem Resultat der Länge 2*Nh. Dazu braucht man 2 FFTs der Länge Nf=2Nh.
Beim Resultat macht man "overlap and add".

Eine FFT der Länge Nf kostet 4*Nf*log2(Nf) reelle MACs.
Wir brauchen eine FFT hin, und eine Zurück.

Alles zusammen:
Rechenaufwand für M samples:
M/Nh * 2 *4*Nf*log(Nf)=M/Nh*8*2*Nh*log2(2*Nh)
Rechenaufwand pro Sample
16*log2(2*Nh)=16*7=112 MACs/Sample

Eine direkte Implementation brauch Nh=64 MACs

Für Nh=64 ist die direkte Impplementation auf einer CPU vermutlich noch 
schneller, bei Nh=128 ist gerade der Break/even.

So, ich hoffe ich habe nichts vergessen.
Hat man auf nem FPGA natürlich genügend DSP/MAC Blöcke kann man beide 
Verfahren parallelisieren, da müsste man dann Aufwand/Geschwindigkeit 
genauer analysieren.

MAC - Multiply ACummulate

Ist z.B. im Befehlssatz bei ARM CPU als INT schon implementiert, ggf. 
auch in ausreichender Auflösung in der FPU beim F4/F7.

Ich hab's im Thread immer noch nicht gefunden - in welchem Zahlenformat 
willst Du rechnen? Int/Float  Komplex  Bittiefe?

Die Auflösung ist wohl der entscheidende Punkt bei der Betrachtung, will 
man das direkt vergleichen.

Es gibt aber noch einen anderen Punkte und wir hatten die Thematik bei 
kontinuierliche Signalen wie Audio schon mehrfach:

Bei jeder Filterung mit einem FIR entsteht eine Fensterung, welche das 
Signal mit einem Spektrum beaufschlagt. Lässt man das Fenster weg, hat 
man ein Rechteck. Das ist bei Audio und den meisten anderen Anwendungen 
untauglich. Daher nutzt man ein zusätzliches Fenster, das dan 
seinerseits eine Welle aufprägt.

Will man das richtig machen, muss sich das Fenster sehr weit über die 
Periode der tiefsten Signalkomponente hinaus erstrecken. Also braucht 
man einen langen Filter.

Das geht aber nicht, will man nicht oder es ist nicht effizient, also 
arbeitet man mit Halbbandfiltern, Kaskaden oder parallelen überlappenden 
Strukturen. Damit entsteht das Problem der Rekonstruktion des Signals 
aus mehren Pfaden und auch da braucht es wieder ein Fenster für die 
Resynthese - insbesondere bei der IFFT!

Dieses muss aber möglichst schnell und kurz sein, weil Ich theoretisch 
mit jedem sample switchen oder smoothen muss. Also braucht es einen 
anwendungspezifischen Kompromiss zwischen den Filterlängen, den 
Fensterbreiten und den parallelen Strukturen. Und dann muss man eben 
schauen, wo es Identitäten gibt, wo man gleiche Koeffizienten hat und wo 
man was zusammenfassen kann. Meistens sind diese Lösungen dann sehr 
stark mit dem Fokus auf Rechenzeitersparnis ausgelegt, besonders bei 
Software / Audiosoftware und den einschlägigen Plugins. Ich habe das 
Thema gerade intensivst mit einem Kunden durchgearbeitet.

Beim FPGA sind die Randbedingungen etwas anders, weil man 
Zwischenergebnisse nicht beliebig lange speichern kann, um sie in der 
pipeline weiternutzen zu können, weil jedes Ergebnis Registerspeicher 
kostet und man da balancieren muss. Umgekehrt bietet eine Rechenpipiline 
die Möglichkeit, durch ihre Quasiparallelität Rechenleistung zu 
gewinnen, die sich auf die Abtastfrequenz aufmultipliziert.  Damit ist 
das Rechenverhältnis wieder 1:1 zzgl. der echten Parallelität durch 
mehrere Einheiten.

Wenn man also z.B. 20KHz Audio mit 192kHz-Rate verarbeiten will (und 
auch das passiert mit komplexen Zahlen!) benötigt man für (2xLambda/4 
für 24 Hz) =12Hz etwa 16000 TAPs. Bei rund 200Mhz wären das 128 
zeitgleiche Multiplikationen, die sich durch Spiegelung der 
Filterkoeffizienten und Benutzung eines Speichers mit der halben 
Filterlänge für das overlapping jeweils 2x teilen lassen.

Macht dann 16 Vor-Addierer, 32 Multiplier, 32 Nach-Addierer in den 
Akumulatoren der pipelines und später 8 Multiplizierer + 2 Addieren für 
die Reko, jeweils mit voller Taktrate für einen Monokanal. Hinzukommt 
etwas Verwaltung für Zähler, Addressumrechnung und enables. Das Fenster 
ist den Filterkoeffizienten integriert und arbeitet auf 12Hz. Davor 
klebt noch ein hochgenauer IIR als 3Hz Hochpass. Wenn man den hinteren 
Filter verdoppelt, kann man es sich auch leisten nur jedes zweite Sample 
zu verwenden und mit 96kHz zu arbeiten. Der Platzbedarf geht dann vorne 
quadratisch!
In der ersten Version des Filter waren es noch 48kHz, was für die 
damalige Zeit aktuell war. Für Stereo lief das dann sogar in einem 
Spartan 3.

http://96khz.org/htm/precisionaudiofilter.htm

Wenn da jemand einen Vorschlag hat, wie man das noch kleiner kriegt, bin 
Ich offen.

Mir eschlisst sich nicht, was das heissen soll:

"Bei jeder Filterung mit einem FIR entsteht eine Fensterung, welche das
Signal mit einem Spektrum beaufschlagt. Lässt man das Fenster weg, hat
man ein Rechteck. Das ist bei Audio und den meisten anderen Anwendungen
untauglich. Daher nutzt man ein zusätzliches Fenster, das dan
seinerseits eine Welle aufprägt."

Ich designe meine Filter immer so, wie ich anhand der Forderungen 
festgelegt habe. Bei hohen Anforderungen wird das FIR Filter dann shr 
lang. Mit "Fenstern und Welle" habe ich dabei erstmal nichts zu tun 
(evtl. das Filterdesignprogramm intern). Ich kontrolliere immer nur 
Frequenz- und Phasengang (Gruppenlaufzeit). Und wenn ich das Filter 
lauen lasse, filtert es genauso we im Design festgelegt.