Präzise OpenCL-FFT

Gast #2459772
Ich arbeite gerade an einer blockweisen Korrelationssuche im 
Frequenzbereich.
Das Verfahren funktioniert wunderbar mit einer "klassischen" CPU-FFT, 
der Cuda-FFT - aber versagt jämmerlich bei der OpenCL-FFT, wenn die 
FFT-Größe 1024 überschreitet (ich verwende die in JavaCL RC1 angebotene 
FFT).

Die Ursache liegt in der Ungenauigkeit der FFT-Berechnung. Genauer 
gesagt liefert
 iFFT(FFT(x)*conj(FFT(y)))
mit x und y als reelle Vektoren keinen reellen Ergebnisvektor.

Deshalb die Frage: gibt es irgendwo einen brauchbaren FFT-Kernel für 
Float-Felder, der möglichst auch die reel->komplex Transformation 
beherrscht - und trotzdem in der FFT-Länge variabel ist (Zweierpotenzen 
vorausgesetzt).
Gast
Das hängt mit der internen Genauigkeit der FFT-Implementierung zusammen.

Und da muss man leider länger ausholen: CUDA unterstützt den Datentyp 
double nativ (ab Fermi-GPU). Bei OpenCL wird bis heute der Typ double 
nur optional unterstützt, die Pragma ist
1
#pragma OPENCL EXTENSION cl_khr_fp64 : enable

Wie OpenCL damit umgeht, bleibt aber den Anbeiter überlassen. Sprich: 
Selbst wenn der JavaCL-Wrapper eine double-FFT anbietet, muss OpenCL 
keineswegs mit double arbeiten und konvertiert vor der FFT zu float.

Richard schrieb:
> Deshalb die Frage: gibt es irgendwo einen brauchbaren FFT-Kernel für
> Float-Felder, der möglichst auch die reel->komplex Transformation
> beherrscht - und trotzdem in der FFT-Länge variabel ist (Zweierpotenzen
> vorausgesetzt).

Ja, viele. Der Fachausdruck heißt "Split-Radix-FFT", die Länge des 
Vektors ist dann frei wählbar und ist noch nicht mal an 2er-Potenzen 
gebunden. Etwa in FFTW3, Parallel-FFTW, CuFFT, MathNet Numerics uvam. zu 
finden. Leider weiß ich nicht, ob es hierfür Java-Wrapper gibt.
Admin Persönliche Seite
Richard schrieb:
> Die Ursache liegt in der Ungenauigkeit der FFT-Berechnung. Genauer
> gesagt liefert
>  iFFT(FFT(x)*conj(FFT(y)))
> mit x und y als reelle Vektoren keinen reellen Ergebnisvektor.

Darauf darf man sich auch niemals verlassen, weil es bei beschränkter 
Rechengenauigkeit nicht garantiert werden kann.

Ist dein Problem wirklich die Genauigkeit, oder dass dein Verfahren (was 
auch immer es ist) nicht damit klarkommt wenn ein Wert != 0 im 
Imaginärteil auftaucht?
Persönliche Seite
Richard schrieb:

> Die Ursache liegt in der Ungenauigkeit der FFT-Berechnung. Genauer
> gesagt liefert
>  iFFT(FFT(x)*conj(FFT(y)))
> mit x und y als reelle Vektoren keinen reellen Ergebnisvektor.

Das ist aber nichts Ungewöhnliches, bei Matlab hat mich das Gleiche auch 
immer wieder geärgert. Bei näherer Betrachtung sind die imaginären 
Anteile dann irgendwo im Bereich von 10-15 Größenordnungen kleiner als 
die realen Anteile, so dass man sie getrost ignorieren kann.

Wenn du weißt, dass das Ergebnis reell sein muss, warum wirfst du dann 
den Imaginärteil nicht einfach weg (Korrektheit des Algorithmus 
natürlich vorausgesetzt)?

Max

Antwort schreiben

Eine eigene Frage? Dieser Thread ist über 14 Jahre alt. Stell sie in einem eigenen Thread – dort lesen mehr Leute mit, und sie wird öfter und schneller beantwortet.

Neuen Thread starten

Eine Lösung, Ergänzung oder Rückmeldung zu diesem Thread? Dann antworte hier.

Formatierung
  • [c]C-Code[/c]
  • [code]Code in anderen Sprachen, ASCII-Zeichnungen[/code]
  • [math]Formel in LaTeX-Syntax[/math]
  • [[Titel]] – Link zu einem Artikel
  • Verweis auf einen anderen Beitrag: Rechtsklick auf dessen Beitrags-ID (#8078628), „Adresse kopieren“, in den Text einfügen

Mehr zur Formatierung

Längeren Sourcecode als Anhang, nicht in den Text. JPG nur für Fotos und Scans; Zeichnungen und Screenshots als PNG oder GIF. Siehe Bildformate.

Schon registriert? Anmelden

Wird nicht angezeigt. Wir schicken dir einen Link; sobald du damit bestätigst, erscheint dein Beitrag.

Mit dem Abschicken legst du einen kostenlosen Account an und erkennst die Nutzungsbedingungen und die Datenschutzerklärung an.