La taille du vecteur de sortie FFTW est incorrecte après l'exécution du plan FFT.

Question

J'ai des problèmes pour exécuter un plan FFT avec la bibliothèque fftw. Je déclare fftIn et fftOut vecteurs avec la taille 2048, et en définissant un fftwf_plan_dft_1d avec eux.

Mais lorsque j'exécute le plan, le vecteur fftOut se met tout à coup à mal. Dans mon code, je peux voir en déboguant que sa taille va de 2048 à 0. Quand j'ai écrit le petit exemple exécutable ci-dessous, je reçois juste de mauvaises tailles comme 17293858104588369909.

Bien sûr, quand j'essaie d'accéder à la premier élément du vecteur, un SIGSEGV arrive.

Mon code:

#include <complex> 
#include <iostream> 
#include <fftw3.h> 
#include <vector> 
using namespace std; 

typedef std::vector<std::complex<float>> fft_vector; 

int main() { 
    const unsigned int VECTOR_SIZE = 2048; 
    fft_vector* fftIn = new fft_vector(VECTOR_SIZE); 
    fft_vector* fftOut = new fft_vector(VECTOR_SIZE); 

    fftwf_plan fftPlan = fftwf_plan_dft_1d(VECTOR_SIZE, reinterpret_cast<fftwf_complex*>(fftIn), reinterpret_cast<fftwf_complex*>(fftOut), FFTW_FORWARD, FFTW_ESTIMATE); 

    fftwf_execute(fftPlan); 

    std::cout << fftOut->size() << std::endl; 
    std::cout << fftOut->at(0).real() << fftOut->at(0).imag() << std::endl; 

    return 0; 
} 

Bien sûr, je sais fftIn vecteur est vide dans cet exemple, mais la sortie est cassée quand il est pas de toute façon. Dans ce cas, le SIGSEGV se passe dans le second cout comme décrit précédemment. Mon code complet a des threads (mais FFT se passe tous dans le même thread, donc les conditions de course ne devraient pas s'appliquer), et c'était une des raisons pour essayer d'isoler le code dans ce petit exemple, juste au cas où, mais il semble avoir quelque chose de mal de toute façon.

Une idée?

Answer 1

Le principal problème est que vous passez un vecteur, ce qui ne fonctionne pas; vous devez passer le contenu du vecteur: & ((* fftIn) [0]) et & ((* fftOut) [0]) ou quelque chose de comparable. En l'état, vous dites à fftw de taper sur les métadonnées de l'objet vectoriel (y compris la longueur, ce qui explique pourquoi c'était parfois 0 et parfois du charabia). Il écrivait littéralement au début de la structure vectorielle parce que c'est ce à quoi pointe le pointeur. En outre, il utilisait les métadonnées de fftIn comme partie de l'entrée de la fft, ce qui n'est pas non plus ce que vous voulez. au lieu

Vous pouvez envisager d'utiliser fftw_complex et fftw_malloc, qui veillera à ce que vos données sont stockées aligné la façon dont a besoin fftw: http://www.fftw.org/doc/SIMD-alignment-and-fftw_005fmalloc.html Vous pouvez mettre fftOut dans un vecteur plus tard si vous avez vraiment besoin. Cela vous permettra d'obtenir les avantages de SIMD, si disponible, et d'éviter tout comportement spécifique au compilateur ou à la plate-forme avec des types complexes et l'allocation de mémoire. L'utilisation du type et de l'allocateur de fftw signifie que votre code fonctionnera toujours de manière optimale et comme prévu, quelle que soit la plateforme sur laquelle vous l'exécutez et quel compilateur vous utilisez.

Voici un exemple de votre transformation de la documentation de FFTW (http://www.fftw.org/doc/Complex-One_002dDimensional-DFTs.html):

#include <fftw3.h> 
... 
{ 
    fftw_complex *in, *out; 
    fftw_plan p; 
    ... 
    in = (fftw_complex*) fftw_malloc(sizeof(fftw_complex) * N); 
    out = (fftw_complex*) fftw_malloc(sizeof(fftw_complex) * N); 
    p = fftw_plan_dft_1d(N, in, out, FFTW_FORWARD, FFTW_ESTIMATE); 
    ... 
    fftw_execute(p); /* repeat as needed */ 
    ... 
    fftw_destroy_plan(p); 
    fftw_free(in); fftw_free(out); 
} 

Essayez de faire correspondre cet exemple et votre code devrait faire ce que vous attendez.

EDIT: Si vous devez utiliser C++ complexe et vecteur, cela devrait fonctionner:

const unsigned int VECTOR_SIZE = 2048; 
fft_vector* fftIn = new fft_vector(VECTOR_SIZE); 
fft_vector* fftOut = new fft_vector(VECTOR_SIZE); 

fftwf_plan fftPlan = fftwf_plan_dft_1d(VECTOR_SIZE, 
    reinterpret_cast<fftwf_complex*>(&(*fftIn)[0]), 
    reinterpret_cast<fftwf_complex*>(&(*fftOut)[0]), 
    FFTW_FORWARD, FFTW_ESTIMATE); 

fftwf_execute(fftPlan); 

std::cout << fftOut->size() << std::endl; 
std::cout << fftOut->at(0).real() << fftOut->at(0).imag() << std::endl; 

Notez que vous devez déréférencer le pointeur de vecteur pour obtenir le [] opérateur de travailler; Obtenir l'index 0 vous donne le premier élément des données vectorielles réelles et donc l'adresse de cet élément est l'adresse (pointeur) que vous devez donner à fftw_plan.

Answer 2

vous ne pouvez pas reinterpret_cast un vecteur vers un pointeur sur des éléments complexes: vous accédez aux variables internes vectorielles (comme la taille), pas nécessairement aux données vectorielles brutes, ce qui explique probablement pourquoi l'objet est jeté en sortie, et tout les mauvaises choses qui arrivent ensuite.

D'autre part, la coulée de l'adresse des données vectorielles fonctionne.

Il suffit de faire (pas besoin de new les vecteurs aussi, utilisez la déclaration standard):

fft_vector fftIn(VECTOR_SIZE); 
fft_vector fftOut(VECTOR_SIZE); 
// load your data in fftIn 
fftwf_plan_dft_1d(VECTOR_SIZE, static_cast<const fftwf_complex*>(&fftIn[0]), static_cast<fftwf_complex*>(&fftOut[0]), FFTW_FORWARD, FFTW_ESTIMATE);