mikrocontroller.net

Forum: Digitale Signalverarbeitung / DSP FFT: Numerische Genauigkeit?


Autor: Theodor (Gast)
Datum:

Bewertung
0 lesenswert
nicht lesenswert
Hi, wenn ich eine FFT mit 2^x Stützstelen nach dem Radix-2-Alogrithmus 
auf Basis von 16-Bit-Zahlen berechnen lasse, wie kann ich dann die 
numerische
Genauigkeit berechnen? Anders: wieviele Bits im Ergebnis kann ich 
wegschmeißen?

Autor: Joe Joe (j_955)
Datum:

Bewertung
0 lesenswert
nicht lesenswert
Hi,

genau darüber zerbreche ich mir auch den Kopf. Die Frage ist, wie du die 
FFT implementierst(Sprache,Plattform). Ich bin dabei, das ganze in C# zu 
machen, und dies ist wesentlich bequemer, weil dort Werkzeuge vorhanden 
sind, welche z.B. die Rundungsfehler dokumentieren.(zumindest im VS)

Grüße

Autor: Zwölf Mal Acht (hacky)
Datum:

Bewertung
0 lesenswert
nicht lesenswert
>Die Frage ist, wie du die FFT implementierst (Sprache,Plattform).

Eher nicht. Zuerst muss man die Theorie begriffen haben. Dann kann man 
voraussagen welche Fehler wie propagieren.
Alternativ kann man duch Tests herausfinden wie ein 2bit Rauschen das 
Spektrum veraendert.

Autor: Detlef _a (detlef_a)
Datum:

Bewertung
0 lesenswert
nicht lesenswert
Für die FFT muß Du die (16 Bit-) Eingangswerte ja mit cos/sin 
Koeffizienten multiplizieren. Die sind in Integer-Darstellung (aber auch 
in float-Darstellung) ja notwendigerweise quantisiert, üblich sind 
Coeffizientenbreiten, die die Hardware verarbeiten kann, also z.B. 16Bit 
oder 18Bit. Das Ergebnis wird dann wieder runterskaliert, also die 
unteren Bits werden weggeschmissen.

Eine 2^x - sample FFT hat x 'stages', die jeweils die Multiplikation mit 
sin/cos machen und eine Addition. Falls die Multiplikation wieder auf 
die ursprüngliche Bitbreite skaliert kann trotzdem pro Stage durch die 
Addition ein Bit 'hinzukommen'. Dann besteht die Gefahr des Überlaufs. 
Also entweder mit höheren Bitbreiten rechnen oder nach jeder stage 
abhängig von den Ergebnissen skalieren. Dieser Bitgewinn reflektiert, 
dass sich zB für einen Sinus die Energie in einer Spektrallinie 
konzentriert.

Also: die Koeffizienten geeignet skalieren, pro stage ein Bit 'Gewinn' 
einrechnen, bißchen Reserve lassen, ggf. dynamisch skalieren. Bei 16Bit 
Eingangsdaten beträgt Dein SNR sowieso maximal ca. 96dB, da machts dann 
keinen Sinn mit 32 oder sogar 40 Bit Zahlen zu rechnen.

Cheers
Detlef

Autor: Andreas Schwarz (andreas) (Admin) Benutzerseite Flattr this
Datum:

Bewertung
0 lesenswert
nicht lesenswert
Jo-hann Ha-nnes schrieb:
> Ich bin dabei, das ganze in C# zu machen

Da kannst du einfach float oder double verwenden und brauchst dir keine 
Gedanken mehr über die Genauigkeit machen.

Ein Trick den man bei Festkommazahlen anwenden kann ist, die Daten vor 
der FFT durch Schieben so zu skalieren dass der Wertebereich voll 
ausgenutzt wird. Auf diese Weise kann man auch bei Signalen mit großem 
Dynamikbereich oft mit relativ kurzen Wortlängen auskommen, und man 
erhält einen Fehler der näherungsweise proportional zur 
Nutzsignalleistung ist.

Autor: Joe Joe (j_955)
Datum:

Bewertung
0 lesenswert
nicht lesenswert
Detlef _a schrieb:
> Eine 2^x - sample FFT hat x 'stages', die jeweils die Multiplikation mit
> sin/cos machen und eine Addition.

Ich glaube pro Stage sind das N/2 kompl. MUL und eine kompl. ADD, wobei 
N die Anzahl meiner Punkte ist, die sich ergeben aus 2^x

Detlef _a schrieb:
> Bei 16Bit
> Eingangsdaten beträgt Dein SNR sowieso maximal ca. 96dB

Bedeutet das automatisch, dass ich mir jetzt das SNR ausrechnen kann und 
darüber eine Aussage habe, wie sich die Quantisierung auf mein 
Ergebnis(-> sprich Spektrum) auswirkt? Im Prinzip sind das ja dann nur 
die von Null versch. Werte, die im Spektrum auftauchen, richtig?

Meine mich zu erinnern, dass mir die Quantisierung beim AD/DA Wandler 
mal über den Weg gelaufen ist, kann ich dies dann soweit einfach 
übertragen und hier anwenden? Im Prinzip habe ich ja ein Signal und ein 
Quantisierungsfehlersignal...Hmm.

Autor: Joe Joe (j_955)
Datum:

Bewertung
0 lesenswert
nicht lesenswert
Jo-hann Ha-nnes schrieb:
> Im Prinzip sind das ja dann nur
> die von Null versch. Werte, die im Spektrum auftauchen, richtig?


Geimeint sind natürlich die Werte ausser den Spektralanteilen selbst....

Autor: Detlef _a (detlef_a)
Datum:

Bewertung
0 lesenswert
nicht lesenswert
Zum Quantisierungsrauschen des AD-Wandlers:
http://de.wikipedia.org/wiki/Quantisierungsrauschen

Zum möglichen SNR bei integer-Zahlen: n-Bit erlauben Dir eine 
Quantisierung in 2^n Stufen, dh. Du machst einen 'Fehler' von 1/(2^n). 
Der 1 entspricht 0dB: 20*log10(1)=0. Der Fehler in dB beträgt 
20*log10(1/(2^n))= -n*20*log10(2) ~ -n*6.02 , also pi mal Daumen 6dB pro 
Bit.

Oder anders: Wenn Du einen Wert von 0 bis 1 mit einem Bit quantisierts, 
liegt Du im Mittel 0.5 daneben: 20*log10(0.5)=6.02dB SNR

Cheers
Detlef

Autor: Joe Joe (j_955)
Datum:

Bewertung
0 lesenswert
nicht lesenswert
Wie verhält es sich dann mit anderen Datentypen, glaube das mit dem 
Integer ist mir jetzt klar geworden, aber so kann ich das ja dann nun 
nicht auf Fließkommadatentypen übertragen, nicht?

Autor: Theodor (Gast)
Datum:

Bewertung
0 lesenswert
nicht lesenswert
Da hat es ja doch nochgefunkt ;-)

Um mit Float oder Double alles durch zu rechnen reicht der Speicher 
meines MCs leider nicht aus.

Andreas Schwarz schrieb:
-----------------------------------------------------------------------
Da kannst du einfach float oder double verwenden und brauchst dir keine
Gedanken mehr über die Genauigkeit machen.

Ein Trick den man bei Festkommazahlen anwenden kann ist, die Daten vor
der FFT durch Schieben so zu skalieren dass der Wertebereich voll
ausgenutzt wird. Auf diese Weise kann man auch bei Signalen mit großem
Dynamikbereich oft mit relativ kurzen Wortlängen auskommen, und man
erhält einen Fehler der näherungsweise proportional zur
Nutzsignalleistung ist.
-----------------------------------------------------------------------

Kann ich das auch erreichen, indem ich vor den 
Multiplikationen/Additionen meine S16 auf double caste, die Rechnungen 
ausführe und dann wieder in S16 caste um den Wert abzulegen?
Gschwindigkeit spielt bei mir keine große Rolle.

Theodor

Autor: Detlef _a (detlef_a)
Datum:

Bewertung
0 lesenswert
nicht lesenswert
Du weißt, dass float bzw. double beim AVR 4Byte lang sind?

Cheers
Detlef

Autor: Theodor (Gast)
Datum:

Bewertung
0 lesenswert
nicht lesenswert
>Du weißt, dass float bzw. double beim AVR 4Byte lang sind?
Bei mir hat double 8 Byte und das beschränkt sich nicht auf den AVR.
Wie kommst du darauf, dass ich einen solchen verwende?

Autor: Andreas Schwarz (andreas) (Admin) Benutzerseite Flattr this
Datum:

Bewertung
0 lesenswert
nicht lesenswert
Theodor schrieb:
> Kann ich das auch erreichen, indem ich vor den
> Multiplikationen/Additionen meine S16 auf double caste, die Rechnungen
> ausführe und dann wieder in S16 caste um den Wert abzulegen?

Klar, kannst du, wenn Laufzeit und Speicherverbrauch nicht stört.

Autor: Andreas Schwarz (andreas) (Admin) Benutzerseite Flattr this
Datum:

Bewertung
0 lesenswert
nicht lesenswert
Jo-hann Ha-nnes schrieb:
> Wie verhält es sich dann mit anderen Datentypen, glaube das mit dem
> Integer ist mir jetzt klar geworden, aber so kann ich das ja dann nun
> nicht auf Fließkommadatentypen übertragen, nicht?

Fehler bei Floating-Point-Rechnungen lassen sich praktisch nicht 
berechnen, nur abschätzen. In der Regel geht man davon aus dass die 
Fehler bei double (64 Bit) so klein sind dass sie keine Rolle spielen. 
Um herauszufinden ob float (32 Bit) ausreichend ist führt man die 
gleiche Rechnung mit float und double durch und vergleicht das Ergebnis. 
Die Differenz ist das Quantisierungsrauschen das durch die Ungenauigkeit 
von float erzeugt wird, damit kann man das SNR für das Signal berechnen 
und entscheiden ob es ausreichend ist.

Autor: Theodor (Gast)
Datum:

Bewertung
0 lesenswert
nicht lesenswert
>Klar, kannst du, wenn Laufzeit und Speicherverbrauch nicht stört.
In diesem Fall stört es nicht, ich habe ein Datenvolumen
von 16KB an 16-Bit-Werten im Speicher liegen, was dem halben RAM 
entspricht.
Ein paar hundert Byte zusätzlich für die Rechnungen mit double stören 
nicht,
ich kann nur nicht alles gleichzeitig mit double berechnen.
Danke, das hat mir sehr geholfen!

Theodor

Autor: Joe Joe (j_955)
Datum:

Bewertung
0 lesenswert
nicht lesenswert
@Autor:  Andreas Schwarz

thx!! Aber bei Integer ist es dann normales QUantisierungsrauschen, 
welches nur von meiner Bitbreite abhängt, wenn ich das richtig 
verstanden habe....Gut zu wissen!

Autor: Helmut S. (helmuts)
Datum:

Bewertung
0 lesenswert
nicht lesenswert
> Ein paar hundert Byte zusätzlich für die Rechnungen mit double stören
> nicht,...

Wie soll das gehen?
In deiner FFT muss dann alles in double sein.
Du hast ein Riesenproblem mit deinem kleinen Speicher wenn du eine 
16k-FFT in long oder gar in double machen willst. Das mit int16 zu 
machen kann man vergessen.

Autor: Theodor (Gast)
Datum:

Bewertung
0 lesenswert
nicht lesenswert
>Wie soll das gehen?
Indem nur der Butterfly mit double arbeitet.

>In deiner FFT muss dann alles in double sein.
Nein, warum den? Der Butterfly arbeitete vorher mit S16,
jetzt mit double, zweimal gecastet eben.

>Du hast ein Riesenproblem mit deinem kleinen Speicher wenn du eine 16k-FFT
>in long oder gar in double machen willst.
Es ist eine 8K-FFT mit 16-Bit-Werten.
Ich habe kein Problem damit, der Butterfly ist nun auf double geändert
und ich erhalte das richtige Spektrum zu meinen Testsignalen, wie
vorher auch, nur wird es jetzt wohl etwas genauer sein.

>Das mit int16 zu machen kann man vergessen.
Ich weiß nicht, was du für eine Implementierung verwendest, bei mir
klappt das super.

Autor: Theodor (Gast)
Datum:

Bewertung
0 lesenswert
nicht lesenswert
P.S.: Wenn ich ein Signal mit 86db Dynamikbereich (16-Bit) reingebe, 
kann mein Spektrum auch höchstens diesem Bereich aufweisen.
Insofern macht die Verwendung von double für alle Ergebniswerte nicht 
viel Sinn. Da spätestens bei der 86-db-Schwelle der Rauschpegel liegt.

Autor: Detlef _a (detlef_a)
Datum:

Bewertung
0 lesenswert
nicht lesenswert
Mein Beitrag vom 24.8.2010 ist fehlerhaft und ich muß ihn vollinhaltlich 
zurücknehmen. Die 6dB stimmen natürlich, aber die begründung war falsch. 
Richtig lautet sie so:


Ein Sinus 2*A Spitze/Spitze wird mit n Bit quantisiert, dann beträgt die 
Quantisierungsstufe q=2A/2^N . Der Quntisierungsfehler sei im Intervall 
-q/2 bis q/2 gleichverteilt, also ist die Wahrscheinlichkeitsdichte im 
Intervall konstant = 1/q .

Jetzt: Die Leistung eines stochastischen Signals der 
Wahscheinlichkeitsdichte p(e) beträgt :

sq=Integral -unendlich bis + unendlich von e^2*p(e)de, also bei uns:

integral -q/2 bis q/2 von e^2*1/q de =  (1/q) *(1/3)* q^3 von -q/2 bis 
q/2 =q^2/12

Der Effektivwert des Sinus beträgt A/sqrt(2), die mittlere Leistung also 
A^2/2

Mittlere Leistung im Verhältnis zu Rauschleistung also:

(A^2/2)/(q^2/12)=(A^2/2)/(4A^2/2^2N)=3/2*2^2N.

10*log(3/2*2^2N)= 1.76dB+N*6.02dB.

Schlüssel ist die Leistung eines stochastischen Signals der 
Amplitudendichteverteilung 1/q .

http://de.wikipedia.org/wiki/Quantisierungsrauschen
http://en.wikipedia.org/wiki/Quantization_error
https://ccrma.stanford.edu/~jos/mdft/Round_Off_Err...

So ist's richtig
math rulez!

Cheers
Detlef

Antwort schreiben

Die Angabe einer E-Mail-Adresse ist freiwillig. Wenn Sie automatisch per E-Mail über Antworten auf Ihren Beitrag informiert werden möchten, melden Sie sich bitte an.

Wichtige Regeln - erst lesen, dann posten!

  • Groß- und Kleinschreibung verwenden
  • Längeren Sourcecode nicht im Text einfügen, sondern als Dateianhang

Formatierung (mehr Informationen...)

  • [c]C-Code[/c]
  • [avrasm]AVR-Assembler-Code[/avrasm]
  • [code]Code in anderen Sprachen, ASCII-Zeichnungen[/code]
  • [math]Formel in LaTeX-Syntax[/math]
  • [[Titel]] - Link zu Artikel
  • Verweis auf anderen Beitrag einfügen: Rechtsklick auf Beitragstitel,
    "Adresse kopieren", und in den Text einfügen




Bild automatisch verkleinern, falls nötig
Bitte das JPG-Format nur für Fotos und Scans verwenden!
Zeichnungen und Screenshots im PNG- oder
GIF-Format hochladen. Siehe Bildformate.
Hinweis: der ursprüngliche Beitrag ist mehr als 6 Monate alt.
Bitte hier nur auf die ursprüngliche Frage antworten,
für neue Fragen einen neuen Beitrag erstellen.

Mit dem Abschicken bestätigst du, die Nutzungsbedingungen anzuerkennen.