Hi Ich hoffe, dass ich mit meinem Anliegen in diesem Unterforum den richtigen Platz gefunden habe. Es geht um folgendes. An einen Cortex M3 (80MHz) wird ein Mikrofon angebunden, das mit 48kHz abtastet, 24 Bit pro Samlple. Nun möchte ich in C den Goertzel Algorithmus implementieren. Hier stoße ich an einige Probleme was die berechnung angeht. Die Verwendung der Fließkommadarstellung entspricht meinem natürlichen Zahlengefühl von Rechnungen auf dem Papier, hat aber, der verbreiteten Meinung nach, deutliche Geschwindigkeitsnachteile auf einem Festkommaprozessor ohne FPU. Gibt es hier eine Faustformel, mit der man abschätzen kann, ab wann die Nutzung von Fließkommazahlen auf dem M3 "zu langsam" wird? Ein Cortex M4 kann nicht verwendet werden, um diesem Vorschlag schoneinmal zu entgegnen. Bei der Verwendung mit Festkommazahlen blicke ich noch nicht so ganz durch. Im Goertzel Algorithmus werden Samples, die ich in Variablen vom Typ long speichere (32 Bit), mit trigonometrischen Funktionen Multipliziert. Interpretiere ich die 24 Bit breiten Samples als Integerzahlen und die Ergebnisse z.B. von cos(x) als Festkommazahl des Typs 1Q31, dann erhalte ich im Endeffekt ein Ergebnis der Breite 64 Bit. Bei meiner dürftigen Literaturrecherche kam dann heraus, dass das Ergebnis einer solchen Rechnung um eine Bitstelle nach Links geschoben werden solle und dann müsse man die oberen 32 Bit als Ergebnis nehmen. Setze ich eine solche 32x32 Bit Multiplikation in meinem Programmcode um, dann wird die Multiplikation als MUL Anweisung in Assembler realisiert. Das Ergebnis dieser Operation enthält dann aber leider nur die unteren 32 Bit. Der SMULL Befehl, der ein 64 Ergebnis auf zwei 32 Bit Register verteilt, wäre scheinbar die optimalere Lösung. Leider gibt es für den Compiler (TMS430, Code Composer Studio) keine entsprechende Intrinsic-Funktion oder ähnliches für den M3 (für den M4 schon...). Würde ich dann am Schreiben einer eigenen Inline Assembler Funktion nicht vorbei kommen? Die CMSIS DSP Bibliothek bietet zwar auch Funktionen zur Multiplikation an, diese sind aber für Festkomma Datentypen gedacht (q31, q15 und q7) und sehen keine ganzzahligen Werte vor. Wie gehe ich desweiteren mit Brüchen um, die keinen ganzzahligen Anteil besitzen? Ich habe mal beispielhaft mit der Berechnung des Wertes 2787/48000 herumgespielt. Das Ergebnis soll als Q31 interpretiert werden. Wenn ich die 2787 vorher mit 2^31 multipliziere kommt es natürlich zu einem Überlauf. Wenn ich ((2787<<15)/48000)<<16 rechne und diesen Wert als Argument für die Cosinus-Funktion der CMSIS DSP Lib nehme, dann bekomme ich eine Abweichung von knapp 6,6% im Vergleich zu dem Wert, den mir z.B. mein Taschenrechner für cos(2787/48000) ausgibt. Mit gewissen Ungenauigkeiten bei endlicher Darstellungsbreite rechnen zu müssen ist mir klar. Aber 6,6% kommen mir doch etwas arg groß vor, zumal sich dieser Fehler bei nachfolgenden Multiplikationen ja weiter fortpflanzt. Ich bin für jeden etwas ausführlicher erklärten Gedanken, Hilfsanstoß, etc. sehr dankbar. Schöne Grüße
Moin, Der Verständnistrick ist dabei wohl, dass man die Audiowerte als 1.23-Format auffasst, also von -1.0 bis 0.999... normiert betrachtet (2^23-1 ist die 0.9999). Der Goertzel ist ein klassischer IIR-Filter, bei dem Du ein multiply-accumulate machst. Dabei treten dann beim Addieren natürlich höhere Werte also 0.999... auf. Ein typischer DSP-Akkumulator federt das ab, indem er ein paar Bits mehr vor dem Komma spendiert, z.B. 16.32-Format (48 Bit). Am Schluss musst Du das ganze - sobald Du es in ein Register speicherst - normieren. Die Operationen kannste ja im Prinzip in C ganz "stupide" mit long long Datentypen und bitshift-Operationen nachbilden, nur wird ev. nicht der optimale Code generiert. Bei den verschiedenen ARM cores kenne ich mich jetzt nicht generell mit den Bitformaten aus. Würde mich da also auch auf die Libraries verlassen, mit 1.31 solltest Du auf jeden Fall ohne relevante Rundungsprobleme den Goertzel erschlagen können. Wichtig ist, dass die Library einen Overflow signalisieren kann. Oder Du kannst dein Berechnungs-Intervall für den Goertzel begrenzen (siehe auch "Sliding Goertzel"). Die 6.6% Abweichung kann ich so nich nachvollziehen. Aus der Fix-point-Division sollte sie auf jeden Fall nicht kommen. Also: N = 2^m, dann sollte [ ((N * 2787) / 48000)) >> m ] (als Integer-Division) nahe an derselben Berechnung in Float liegen. m wäre hier 31. Grüsse, - Strubi
Hi Sorry, dass ich mich jetzt erst melde. Deine Antwort hat mit sehr gut geholfen und mit einem (long long) Cast und dem passenden Schieben des Ergebnisses komme ich zumindest in Sachen Multiplikation zu dem gewünschten Ergebnis. Zwei Fragen sind nur noch offen: 1) Wie gehe ich mit Überläufen bei der Addition bzw. Subtraktion um? Also wie sieht dann die normierung sinnvollerweise aus? Könntest du mir ein kleines Beispiel zeigen? 2) Wenn ich die Samlpewerte als 1.23 Format interpretiere und in einer 32 Bit Variabel speichere, dann würde ich die Abtastwerte ja in den oberen 24 Bit Speichern und die unteren 8 Bit sind halt null, korrekt? Vielen Dank und schöne Grüße
Hi André, > > Sorry, dass ich mich jetzt erst melde. Deine Antwort hat mit sehr gut > geholfen und mit einem (long long) Cast und dem passenden Schieben des > Ergebnisses komme ich zumindest in Sachen Multiplikation zu dem > gewünschten Ergebnis. Zwei Fragen sind nur noch offen: > > 1) Wie gehe ich mit Überläufen bei der Addition bzw. Subtraktion um? > Also wie sieht dann die normierung sinnvollerweise aus? Könntest du mir > ein kleines Beispiel zeigen? > In C gibt's keinen plattformunabhängigen Overflow-Check oder ein Carry-Bit (siehe Wiki), geht also nur so (Beispiel 1.15 Format-Emulation):
1 | long l; |
2 | ... |
3 | if (l > 32767 || l < -32768) overflow = 1; |
Unter Umst. ist es schneller, obigen Overflow ohne if-Abfrage (mit entsprechendem Bit-Gefuddel) zu erledigen. > 2) Wenn ich die Samlpewerte als 1.23 Format interpretiere und in einer > 32 Bit Variabel speichere, dann würde ich die Abtastwerte ja in den > oberen 24 Bit Speichern und die unteren 8 Bit sind halt null, korrekt? > Das hängt davon ab, wo du den Kompromiss zwischen Genauigkeit und 'Headroom' ansetzt -- und ob's Du ein einfaches Register oder einen Akkumulator "emulierst". Für die Multiplikation von Werten mit Betrag < 1 reicht das so, aber wenn Du akkumulierst, wird's ja grösser. Wenn Du mehr MSB headroom spendierst, sagen wir 4 (=4.28 Format), kannst Du min 16 Werte ohne Overflow addieren. Dann hast Du noch 4 Bit für die hinteren Nachkommastellen, die schlussendlich abgeschnitten werden. Hast Du von letzeren noch weniger (zugunsten Headroom), treten die berüchtigen Rundungsfehler in den Vordergrund. Einfach dabei nicht vergessen, dass eine x.y-Format (Fixpoint-) Multiplikation eine klassische Integer-MUL mit anschliessendem Shift ist. Und am Schluss noch die Normierung (also Shift) vom Akku in ein Ergebnis-Register. Die üblichen DSPs haben dafür ne Menge Rundungs- und Shift-Modifier schon eingebaut, die Du in C teils recht mühsam implementieren musst. Das Shiften ist Pipifax, aber die Abfrage einzelner Bitsequenzen innerhalb einer dicken Schleife kann ganz schön Zeit verbraten (if-Abfragen sind fast immer schlecht für die Pipeline) Es scheiden sich ja oft die Geister, aber ich bin immer noch der Meinung, dass sich in einem solchen Fall für die "inner loops" auszahlt, die Geschichte in Assembler umzusetzen. Dann muss man sich halt auf die Architektur bisserl festlegen. Dafür hat man die Arithmetik im Griff und keine garstigen Artefakte wie bei nicht ganz konform arbeitenden FPUs (davon liefen mir schon n paar übern Weg). Inzwischen gehen übrigens einige (mich eingeschlossen) auch den Weg, gewisse überschaubare DSP-Ops gleich per VHDL in Hardware zu giessen, da man sich da seine x.y-Fixkommaformate und Pipeline beliebig zusammenstellen kann. Grüsse & viel Erfolg, - Strubi
Hallo, hatte mal meinen C-Code für Goertzel hier Beitrag "1750 Hz Ton Auswertung mit TMS320VC5505" gepostet. Der geht. Cheers Detlef
Erstmal danke für eure Antworten Detlef _a schrieb: > hatte mal meinen C-Code für Goertzel hier > > Beitrag "1750 Hz Ton Auswertung mit TMS320VC5505" > > gepostet. Der geht. Den Code habe ich im Grunde schon, wobei es natürlich nicht schadet den Vergleich zu haben. Martin S. schrieb: >> 2) Wenn ich die Samlpewerte als 1.23 Format interpretiere und in einer >> 32 Bit Variabel speichere, dann würde ich die Abtastwerte ja in den >> oberen 24 Bit Speichern und die unteren 8 Bit sind halt null, korrekt? >> > > Das hängt davon ab, wo du den Kompromiss zwischen Genauigkeit und > 'Headroom' ansetzt -- und ob's Du ein einfaches Register oder einen > Akkumulator "emulierst". Für die Multiplikation von Werten mit Betrag < > 1 reicht das so, aber wenn Du akkumulierst, wird's ja grösser. > Wenn Du mehr MSB headroom spendierst, sagen wir 4 (=4.28 Format), kannst > Du min 16 Werte ohne Overflow addieren. Dann hast Du noch 4 Bit für die > hinteren Nachkommastellen, die schlussendlich abgeschnitten werden. Hast > Du von letzeren noch weniger (zugunsten Headroom), treten die > berüchtigen Rundungsfehler in den Vordergrund. So, jetzt ist die Frage, wenn ich mit dem 1.31 Format arbeite und die 24 Bit der Samples an die untere Grenze lege, so dass das LSB des Eingangssignal auch das LSB im 1.31 Format ist, dann habe ich doch keine Verluste der Bits und damit keine Verschlechterung des SNR aber trotzdem nach oben "dicke" Platz Platz für Summationen. Da der Goerztel sowohl positive als auch negative Summanden enthalten kann, bewege ich mich statistisch gesehen ja immer im Rahmen einer Standardabweichung nach oben und nach unten. Nach unten sind es im schlimmsten Fall 24 bit und nach oben bestenfalls 7 Bit, bis es zu einem Überlauf kommt. Oder irre ich mich da? So hätte ich ja genug Spielraum und die Multiplikation verhält sich wie sonst auch. Martin S. schrieb: > Die üblichen DSPs haben dafür ne Menge Rundungs- und Shift-Modifier > schon eingebaut, die Du in C teils recht mühsam implementieren musst. > Das Shiften ist Pipifax, aber die Abfrage einzelner Bitsequenzen > innerhalb einer dicken Schleife kann ganz schön Zeit verbraten > (if-Abfragen sind fast immer schlecht für die Pipeline) > Es scheiden sich ja oft die Geister, aber ich bin immer noch der > Meinung, dass sich in einem solchen Fall für die "inner loops" auszahlt, > die Geschichte in Assembler umzusetzen. Dann muss man sich halt auf die > Architektur bisserl festlegen. Dafür hat man die Arithmetik im Griff und > keine garstigen Artefakte wie bei nicht ganz konform arbeitenden FPUs > (davon liefen mir schon n paar übern Weg). > > Inzwischen gehen übrigens einige (mich eingeschlossen) auch den Weg, > gewisse überschaubare DSP-Ops gleich per VHDL in Hardware zu giessen, da > man sich da seine x.y-Fixkommaformate und Pipeline beliebig > zusammenstellen kann. Assembler möchte ich eher vermeiden, wobei ich mir natürlich immer wieder angucke, was der Compiler für Code erzeugt und schau, ob sich über den C Code etwas optimieren lässt. Und klar, für eine optimale Implementierung sollte man vielleicht auch mit anderen Hardwareplattformen liebäugeln. Die Rahmenbedingungen sind aber, wie sie sind und damit lässt sich nicht an ihnen rütteln. Ich schätze aber, dass die Leistung des Cortex bei der Abtastrate von 48kHz seine Arbeit locker erledigen können sollte. Martin S. schrieb: > Grüsse & viel Erfolg, Danke, ich hoffe, ich werde ihn haben ;-)
Hi André, > So, jetzt ist die Frage, wenn ich mit dem 1.31 Format arbeite und die 24 > Bit der Samples an die untere Grenze lege, so dass das LSB des > Eingangssignal auch das LSB im 1.31 Format ist, dann habe ich doch keine > Verluste der Bits und damit keine Verschlechterung des SNR aber trotzdem > nach oben "dicke" Platz Platz für Summationen. Da der Goerztel sowohl > positive als auch negative Summanden enthalten kann, bewege ich mich > statistisch gesehen ja immer im Rahmen einer Standardabweichung nach > oben und nach unten. Nach unten sind es im schlimmsten Fall 24 bit und > nach oben bestenfalls 7 Bit, bis es zu einem Überlauf kommt. Oder irre > ich mich da? So hätte ich ja genug Spielraum und die Multiplikation > verhält sich wie sonst auch. Beim Goertzel bin ich mit Dir nicht ganz einig: Wenn Du die gesuchte Frequenz triffst, "schwingt" das Ding beliebig an, heisst, Du kriegst Sequenzen a la -0.1, 0.2, -0.3, 0.4, ... Wenn Du die LSBs von 32 und 24 bit so wie oben alignierst, wäre es das 8.24-Format. Heisst, vor dem (Fix-)Komma geht Dein Wertebereich von -128.00 bis +127.99.. Aber es geht Dir eben Genauigkeit verloren, da Du bei der 32-Bit-Fixkomma-MUL die unteren Bits durch den unmittelbaren "Shift-Rechts" (y >> x) abschneidest. Gruss, - Strubi
Strubi schrieb: > Beim Goertzel bin ich mit Dir nicht ganz einig: Wenn Du die gesuchte > Frequenz triffst, "schwingt" das Ding beliebig an, heisst, Du kriegst > Sequenzen a la -0.1, 0.2, -0.3, 0.4, ... Wie meinst du das genau? Strubi schrieb: > Aber es geht Dir eben Genauigkeit verloren, da Du bei der > 32-Bit-Fixkomma-MUL die unteren Bits durch den unmittelbaren > "Shift-Rechts" (y >> x) abschneidest. Da hast du allerdings recht. Zumal ich davon ausgehe, dass sich die Messwerte eher im unteren Bereich bewegen und weniger bei der Vollaussteuerung. Das würde enorme Fehler verursachen. Die Frage, die bleibt, ist, wie ich nun aber in der Praxis auf dem uC mit Überläufen bei Summen umgehe? Bei meiner Recherche bin ich darüber gestolpert, dass im 2er Komplement bei Zwischensummen ruhig zu Überläufen kommen kann, das Endergebnis aber trotzdem korrekt ist. Mir ist das aber noch nicht so ganz verständlich.
Und wo wir eh schon beim Goertzel sind...wie schaut das eigentlich aus? Das zu messende Signal ist eine Sinusschwingung mit Rauschen. Für eine normale DFT sollte man ja am besten die Eingangssequenz mit einer Fensterfunktion gewichten, da es ja praktisch kaum möglich ist Vielfache der Periodendauer abzutasten. Zumal der Goertzel-Bandpass ja eh nach N Samples aufhört, egal bei welcher Periode des Eingangssignals. Wird man da trotzdem Fenstern ode wie geht man da vor? In der Literatur habe ich zu der Kombination Goertzel + Fenstern nichts gefunden.
Thread beobachten
Seitenaufteilung abschaltenBitte 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
Mit Google-Account einloggen
Noch kein Account? Hier anmelden.