Forum: Digitale Signalverarbeitung / DSP / Machine Learning FFT für Anfänger


von Michael (Gast)


Lesenswert?

Hallo

Ich möchte auch mal eine FFT programmieren. Hierbei kommt es mir erst 
mal nicht auf die Geschwindigkeit an. Ich habe versucht soviel wie 
möglich zu lesen, jedoch beschreiben die Meisten das ganze immer so 
kompliziert, das es keiner verstehen kann, daher versuche ich es über 
diesen Weg.

Ich möchte für den Anfang ein FFT mit 8 Punkten erstellen. (diese kann 
man vielleicht noch leicht verstehen)

Als µC benutze ich einen ATMEGA32 von Atmel. Mit dem integrierten µC 
digitalisiere ich einen Rechteckimpulse, der durch einen 
Funktionsgenerator eingespeist wird. Die 8 digitalisierten Werte 
speicher ich im internen SRAM ab. Nachdem ich nun meine 8 Werte habe, 
kann ich daraus die FFT berechnen.

Meine digitalisierten Werte lauten wie folgt:

x[0]= 0;  x[1]= 0;  x[2] = 50;  x[3] = 50; x[4] =50; x[5] = 50; x[6] = 
0;  x[7] = 0;

auf einigen Seiten habe ich dann gelesen, das zuerst folgene Summen 
gebildet werden

Summe 1 = x[0] + x[4]
Summe 2 = x[1] + x[5]
Summe 3 = x[2] + x[6]
Summe 4 = x[3] + x[7]

Anschließend Subtrahiert man Werte und multipliziert das Ergebnis mit 
der n-te Einheits­wurzel w^n

Differenz1 1 = (x[0] + x[4])*w^1
Differenz1 2 = (x[1] + x[5])*w^2
Differenz1 3 = (x[2] + x[6])*w^3
Differenz1 4 = (x[3] + x[7])*w^4


Jedoch komme ich mit der Einheits­wurzel nicht klar. Zudem weiß ich 
nicht wie es weitergeht.

Kann mir vielleicht jemand anhand dieses kleinen Beispiels vorrechnen 
wie die FFT funktioniert.

Grüsse

Michael

von Andreas S. (andreas) (Admin) Benutzerseite


Lesenswert?

Hast du verstanden wie die DFT (diskrete Fouriertransformation) an sich 
funktioniert?

von Michael (Gast)


Lesenswert?

Nein habe ich nicht. Wenn es eine gute Erklärung gibt dann wäre ich über 
einen Link dankbar.

Ich will doch nur die Umsetzung. Dies kann doch nicht so schwer sein. 
Sind doch nur ein paar Additionen und Multiplikationen

von pumpkin (Gast)


Lesenswert?

Es ist im Grunde genommen auch nicht schwer.

Was bitte ist die Einheitswurzel? Ich denke damit ist der Twiddle-Faktor 
gemeint. Dieser Twiddle-Faktor ist lediglich ein Zeiger in der komplexen 
Ebene der auf dem Einheitskreis wandert. Stell' dir den Einheitskreis 
vor wie einen Kuchen, der nun aufgeteilt werden soll. Und zwar in N 
gleichgroße Teile. N steht für die Ordnung deiner FFT, in diesem Fall 
ist N=8. Nun schreiben die Techniker aber ungerne immer wieder

1
  e^(-j*2*PI * k / N) = e^(-j*2*PI / N)^k

sondern kürzen sich das Ganze ab mit

1
  (W_N)^k  ("W-hoch-k-Index_N")


Im Fall mit N=8 bekommt man also die Twiddle-Faktoren

1
  W^0 =  1     + j 0
2
  W^1 =  0.707 + j 0.707
3
  W^2 =  0     + j 1
4
  W^3 = -0.707 + j 0.707
5
  W^4 = -1     + j 0
6
  W^5 = -0.707 - j 0.707
7
  W^6 =  0     - j 1
8
  W^7 =  0.707 - j 0.707
9
  W^8 = W^0
10
  W^9 = W^1
11
   .
12
   .
13
   .

Das sind die einzigsten komplexen Faktoren die man hier brauchen wird. 
Zeichne dir das am besten mal auf, dann wird du schnell begreifen wie es 
sich verhält. Zeichne den Einheitskreis auch einmal mit N=16 oder N=4. 
Achte darauf, dass der Zeiger bei 1 + j0 anfängt und mathematisch 
positiv dreht, d.h. gegen den Uhrzeigersinn.

Wenn du das begriffen hast, dann schaue dir ein Butterfly-Plot an um zu 
sehen, wie man eine FFT beispielsweise implementiert.

von Tom (Gast)


Lesenswert?

Bisschen Grundlagen in komplexer Mathematik wären auch hilfreich.

Ich würde erstmal eine DFT in einer Hochsprache programmieren, um das 
Prinzip verstehen zu lernen. Dann eine FFT, was ja eigentlich nur eine 
Folge von 2-Punkt-DFT ist. Auch in einer Hochsprache, denn da hast du 
mehr Möglichkeite, z.B. float-Variablen. Und dann das Ganze in 
Assembler.

von Michael (Gast)


Lesenswert?

Ich habe schon mal mit der Labview-Demo gearbeitet. Dort möchte ich das 
ganze zuerst ausprobieren, bevor ich es in den µC implementiere. Man 
kann da schöne Signale erzeugen und einfach ein Rauschen draufgeben und 
sich anschließend die FFT anschauen. Der Vorteil ist man kann die 
berechnete FFT von Labview mit der eigenen vergleichen. Sollten ja beide 
identisch sein.

von Tommi H. (drmota)


Lesenswert?

Michael bist du vielleicht der Newbie von 
Beitrag "Probleme mit Aufgabe zu digitale Signalverarbeitung"

von Michael (Gast)


Lesenswert?

Nein bin ich nicht. Ich denke mal Michael heißen wohl sehr viele in 
diesen Forum

von Michael (Gast)


Lesenswert?

So ich habe mal ein wenig im Internet nach dem Twiddle-Faktor gesucht. 
Jedoch ist mir dann nicht ganz klar wie Du auf die folgenden Ergebnisse 
gekommen bist:

  W^0 =  1     + j 0
  W^1 =  0.707 + j 0.707
  W^2 =  0     + j 1
  W^3 = -0.707 + j 0.707
  W^4 = -1     + j 0
  W^5 = -0.707 - j 0.707
  W^6 =  0     - j 1
  W^7 =  0.707 - j 0.707
  W^8 = W^0
  W^9 = W^1

e^(-j*2*pi*n*k/N) = cos(2*pi*n*k/N) -j sin(2*pi*n*k/N).

W^0 = cos (2*pi*0*0/8) -j sin(2*pi*0*0/8) => cos(0) - j sin(0) => 1
W^1 = cos (2*pi*1*1/8) -j sind(2*pi*1*1/8) => cos(pi/4) -j sin(pi/4) => 
0,0137

Was mache ich denn falsch denn Du kommst ja bei W^1 auf 0.707 + j 0.707

von Alex (Gast)


Lesenswert?

Hi Newbie,

das für W^0 = 1 + j0 = 1 ist siehst du hoffentlich ein.

Dass cos(pi/4) = 0.707 und sin(pi/4) = 0.707 ist kriegt ein 
Fünftklässler raus, dem man einen Taschenrechner in die Hand drückt.

Wie du es geschafft hast, Real- und Imaginärteil zu 0.0137 zu vermehren 
wird wohl dein Geheimnis bleiben (du weißt was eine komplexe Zahl ist?).

Dass pumpkin die Vorzeichen etwas verdreht hat steht auf einem anderen 
Blatt, er hat sich gegen anstatt mit dem Urzeigersinn auf dem 
Einheitskreis in der komplexen Ebene bewegt.

W^1 = .707 - j.707
...

von Michael (Gast)


Lesenswert?

Ich habe den taschenrechner von Windos verwendet. Vielleicht hat dieser 
sich aber verrechnet, da ich die CPU stark übertaktet habe.

pi/4 = 0,78539816339744830961566084581988

cos(0,785398) ergibt 0,99990
sin(0,785398) ergibt 0,01370

habe leider gerade keinen anderen Taschenrechner zur Hand

von HildeK (Gast)


Lesenswert?

>cos(pi/4) -j sin(pi/4) => 0,0137
Im Taschenrechner RAD als Basis eingeben! Dann kommt das raus, was Alex 
schrieb.
Bei meinem Win-Rechner kommt 0,70710678118654752440084436210485 heraus!

von Michael (Gast)


Lesenswert?

Oh ich hatte bei Windows nicht auf Rad gestell ^^

von Michael (Gast)


Lesenswert?

Jetzt komme ich auch auf das Ergebnis.

Ich habe noch Probleme mit der Formel e^(-j*2*pi*n*k/N)

j ist klar
pi auch
N = Anzahl der abtastwerte

nur was ist denn k und n? Welchen Bereich hat k und n?

von pumpkin (Gast)


Lesenswert?

> Dass pumpkin die Vorzeichen etwas verdreht hat steht auf einem anderen
> Blatt, er hat sich gegen anstatt mit dem Urzeigersinn auf dem
> Einheitskreis in der komplexen Ebene bewegt.
>
> W^1 = .707 - j.707
> ...

Konform zu meiner Definition, ich habe es sogar angemerkt. Ich kenne 
keine Literatur in der es anders ist. Stifte also bitte keine unnötige 
Verwirrung.


> Ich habe den taschenrechner von Windos verwendet. Vielleicht hat dieser
> sich aber verrechnet, da ich die CPU stark übertaktet habe.

GAS GAS GAS!

Ich bin raus.

von pumpkin (Gast)


Angehängte Dateien:

Lesenswert?

Na gut, einen noch.

> nur was ist denn k und n? Welchen Bereich hat k und n?

k und n sind Indizes der DFT. Wenn man die Formel (1) aus dem Anhang 
annimmt, dann steht n für den Frequenzindex und k für den Zeitindex. 
Beides sind ganze Zahlen bzw. beide sind diskret.

von pumpkin (Gast)


Lesenswert?

an Alex:

Entschuldige, jetzt weiss sehe ich erst was du meinst! Das ist mir 
wirklich unangenehm. Korrektur:


1
  W^0 =  1     + j 0
2
  W^1 =  0.707 - j 0.707
3
  W^2 =  0     - j 1
4
  W^3 = -0.707 - j 0.707
5
  W^4 = -1     + j 0
6
  W^5 = -0.707 + j 0.707
7
  W^6 =  0     + j 1
8
  W^7 =  0.707 + j 0.707
9
  W^8 = W^0
10
  W^9 = W^1
11
   .
12
   .
13
   .
14
15
Achte darauf, dass der Zeiger bei 1 + j0 anfängt und mathematisch
16
positiv dreht, d.h. mit dem Uhrzeigersinn.

von Alex (Gast)


Lesenswert?

Mathematisch positiv heißt gegen den Uhrzeigersinn, sorry ;-)

von pumpkin (Gast)


Lesenswert?

Bin ich jetzt total bescheuert? Wenn ich schon das eine ändere dann muss 
ich auch das andere ändern:

1
Mathematisch negativ, d.h. im Uhrzeigersinn.

Sonntägliches Kopfkino! Ich halt' lieber den Rand!   ;^)

von Michael (Gast)


Lesenswert?

Habe ich schon gemerkt. Ich werde das morgen einmal alles durchrechnen 
und posten. Dann habe ich sicher auch noch eine Frage.

von Michael (Gast)


Lesenswert?

Wir sollten hier ein FAQ für FFT einrichten. Mal sehen ob ich es mache.

von Jemand (Gast)


Lesenswert?

Guck mal die gnumeric-Datei da an:
Beitrag "Re: Schnelle FFT in Assembler"
Da hat es jemand ne FFT in Gnumeric nachgebaut. vieleicht hilfts.

von Tobias W. (tobis)


Lesenswert?

also das mit dem drehzeiger ist völlig klar...ganz normaler drehzeiger 
eben.

Doch wie gehts dann weiter???

Butterfly-Plot???

von Jeso (Gast)


Lesenswert?

Hi,

hab auch noch so meine Probleme bei der FFT vor allem weil ich den 
mathematischen zusammenhang net so ganz kapier fehlt halt viel Mathe.

Aber jetzt zum Thema versuch grad mal Testweise in Excel den Butterfly 
nachzumachen. Soll dann später mal auf nen AVR oder 8051 umgesetzt 
werden.

Twiddlefaktor hab ich so langsam auch mal kapiert und hab mir sie auch 
ausgerechnet. Gut so weit.
Jetzt hab ich mir ne Eingangsfolge überlegt. Soll en Rechteck werden.
Hab 16 Samples da sind jeweils 2 0 und 2 = 255. Nun stellt sich bei mir 
die Frage wie verrechne ich jetzt den Twiddlefaktor der ja komplex ist 
mit ner rein reelen eingangszahl z.B. die 255. Wenn das jemand 
anschaulich erklären könnte wäre auch von großem nutzen für die 
Umsetzung auf nem AVR.

Wenn ich mir den Butterfly anschaue und ich den dann bei mir die 4 
Stufen ausgerechnet hab dann hab ich am Ende direkt meine Spektren?

Danke für die Antworten.

von Günter -. (guenter)


Lesenswert?

Die Berechnung ist immer komplex. Selbst wenn dein Eingangssignal real 
ist, "wandelst" du es in eine komplex Zahl um.

Also 255 wird 255 +j0 den Wert multiplizierst du jetzt mit dem 
Twiddelfaktor.

von matthias (Gast)


Lesenswert?

Ein bisschen Googeln würde gut tun.

Fixfertig auf dem Silbertablar:
http://www.codeproject.com/KB/recipes/howtofft.aspx

von JESO (Gast)


Lesenswert?

0   0   0
1   0   0
2  255  0
3  255  0
4   0  255
5   0  255
6  255 255
7  255 255
8   0   0
9   0   0
10 255  0
11 255  0
12  0  255
13  0  255
14 255 255
15 255 255

So hier hab ich meine Eingangsfolge und danach die Bit Reversal folge 
von meinen 16 Werten.
in der 1. Stufe ist ja der Twiddle Faktor 1 und j 0.
Also addiere ich in Zeile null die Zeile  0 mit Zeile 1 und 
Multipliziere mit Twiddle Faktor Real. und in zeile 1 Subtrahiere ich 
Zeile 0 mit Zeile 1 und multipliziere mit Twiddlefaktor Imaginär. Oder 
hab ich das falsch verstanden. Und in Stufe 2 sind es dann 2 
Twiddlefaktoren und so weiter.

von Roland B. (rolandb)


Lesenswert?

Versuche das ganze vlt nicht an einem Rechteck. Da kommt dann was 
komsiches vlt für dich raus (irgendetwas si oder sinc) förmiges.. Nehme 
erstmal einen Sinus oder Cosinus.


mfg

von Jeso (Gast)


Lesenswert?

Werde es auch noch mit einem anderen Eingangsignal versuchen rechteck 
war erstmal das leichteste. Mir war erstmal nur nicht ganz klar wie ich 
den Butterfly dann rechne vorallem wenn ich nur die Samplewerte habe die 
ja keine Komplexe Zahlen sind. Mal ne dumme Frage wieso soll ich kein 
Rechteck nehmen da sollte doch ein relativ glecihmäßig verteiltes 
Spektrum geben. Werde aber mal en paar sinus werte berechnen.

von Jeso (Gast)


Lesenswert?

Irgendwie hab ich komplexe Zahlen noch nicht ganz verstanden.

Wenn ich jetzt 250 + 5j habe und diese mit Twiddelfaktor 0707 - 0707j 
multipliziere kommt raus:  176,75 - 3,53j stimmt dies.

 (250*0,707) + (5j * -0,707j)

von Günter -. (guenter)


Lesenswert?

Nein, die Berechnung ist so:
1
(250 + 5j) * (0707 - 0707j)
2
3
= 250 * 0707 - 250 * 0707j + 5j * 0707 - 5j * 0707j
4
5
j * j = -1

Das Ergebnis ist dann:
1
(250*0707 + 5*0707 +j 5*0707 - 250*0707)
Neuer Realwert ist also 250*0707 + 5*0707 und der imaginäre Wert ist 
5*0707 - 250*0707.

von geri (Gast)


Lesenswert?

hallo leute,

frage zum link: http://www.mikrocontroller.net/attachment/28612/dft.JPG

X(n) = sum(x(k) * exp(-j2*Pi*n*k/M), 0<=n<=M-1) ist ja auch schon im 
frequenzbereich periodisch?

bye

Bitte 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
Noch kein Account? Hier anmelden.