2017-08-13 3 views
1

J'essaye d'implémenter la définition "Naive" de la Transformée de Fourier Discrète en n'utilisant rien d'autre que C++ et ses libs standard, pour ma propre compréhension du matériel mathématique, et je continue à obtenir une sortie incorrecte, malgré que mon code semble être une traduction directe de la définition mathématique de la DFT.Calcul de la DFT selon la définition (C++ avec <complex>)

Voici le code (modifié de façon à nettoyer la sortie mieux que l'original):

void DFT(std::complex<double>* outputs, int N, std::complex<double>* inputs) 
{ 
    for (int k = 0; k < N; ++k) 
    { 
    outputs[k] = std::complex<double>(0.0, 0.0); 
    for (int t = 0; t < N; ++t) 
    { 
     // tk2π/N 
     double angle = t * k * TWOPI/double(N); 
     // e^(-jtk2π/N) 
     std::complex<double> dft_sino = std::polar(1.0, -angle); 

     if(fabs(real(dft_sino)) < std::numeric_limits<float>::epsilon()) 
     dft_sino = std::complex<double>(0.0, imag(dft_sino)); 

     if(fabs(imag(dft_sino)) < std::numeric_limits<float>::epsilon()) 
     dft_sino = std::complex<double>(real(dft_sino), 0.0); 

     std::cout << "DFT sinosoid (" << k << "), Sample " << t << ": " << dft_sino << std::endl; 
     outputs[k] += inputs[t] * dft_sino; 
    } 
    } 
} 

L'ensemble d'entrée est:

1+1i 
-1+1i 
1-1i 
-1-1i 
-3-1i 
4-2i 
0-3i 
2-6i 

Le (vérifié) sortie Je pense est la suivante:

3-12i 
10.24+5.243i 
3+2i 
9.071+10.07i 
-5+4i 
1.757-3.243i 
-9+6i 
-5.071-4.071i 

Et la sortie du programme dans sa forme actuelle est:

3+2i 
3.293+6.364i 
-2-1i 
6.121+6.95i 
-5-4.91e-15i 
4.707-6.364i 
-4+3i 
1.879-2.95i 

Y at-il quelque chose d'évident que je manque ici?

Édition: ci-dessous est la sortie des échantillons pour les échantillons sinusoïdaux DFT correspondants qui sont multipliés par l'entrée pendant la sommation.

DFT sinusoid (0): 
Sample 0: 1+0i 
Sample 1: 1+0i 
Sample 2: 1+0i 
Sample 3: 1+0i 
Sample 4: 1+0i 
Sample 5: 1+0i 
Sample 6: 1+0i 
Sample 7: 1+0i 

DFT sinusoid (1): 
Sample 0: 1+0i 
Sample 1: 0.7071-0.7071i 
Sample 2: 0-1i 
Sample 3: -0.7071-0.7071i 
Sample 4: -1+0i 
Sample 5: -0.7071+0.7071i 
Sample 6: 0+1i 
Sample 7: 0.7071+0.7071i 

DFT sinusoid (2): 
Sample 0: 1+0i 
Sample 1: 0-1i 
Sample 2: -1+0i 
Sample 3: 0+1i 
Sample 4: 1+0i 
Sample 5: 0-1i 
Sample 6: -1+0i 
Sample 7: 0+1i 

DFT sinusoid (3): 
Sample 0: 1+0i 
Sample 1: -0.7071-0.7071i 
Sample 2: 0+1i 
Sample 3: 0.7071-0.7071i 
Sample 4: -1+0i 
Sample 5: 0.7071+0.7071i 
Sample 6: 0-1i 
Sample 7: -0.7071+0.7071i 

DFT sinusoid (4): 
Sample 0: 1+0i 
Sample 1: -1+0i 
Sample 2: 1+0i 
Sample 3: -1+0i 
Sample 4: 1+0i 
Sample 5: -1+0i 
Sample 6: 1+0i 
Sample 7: -1+0i 

DFT sinusoid (5): 
Sample 0: 1+0i 
Sample 1: -0.7071+0.7071i 
Sample 2: 0-1i 
Sample 3: 0.7071+0.7071i 
Sample 4: -1+0i 
Sample 5: 0.7071-0.7071i 
Sample 6: 0+1i 
Sample 7: -0.7071-0.7071i 

DFT sinusoid (6): 
Sample 0: 1+0i 
Sample 1: 0+1i 
Sample 2: -1+0i 
Sample 3: 0-1i 
Sample 4: 1+0i 
Sample 5: 0+1i 
Sample 6: -1+0i 
Sample 7: 0-1i 

DFT sinusoid (7): 
Sample 0: 1+0i 
Sample 1: 0.7071+0.7071i 
Sample 2: 0+1i 
Sample 3: -0.7071+0.7071i 
Sample 4: -1+0i 
Sample 5: -0.7071-0.7071i 
Sample 6: 0-1i 
Sample 7: 0.7071-0.7071i 
+0

Comme même votre terme de sortie à fréquence zéro est faux, cela devrait être facile à déboguer. Le terme de fréquence zéro est littéralement juste la somme des termes d'entrée, donc vous devriez être en mesure de franchir ligne par ligne pour voir où le comportement va mal. –

+0

J'ai essayé de le traverser plusieurs fois. Malheureusement, je ne suis pas assez expérimenté en mathématiques pour savoir ce qui serait considéré comme un comportement inhabituel ou incorrect. Que devrais-je rechercher? –

+0

Vous pouvez décomposer l'expression mathématique complexe en une série d'instructions distinctes, en stockant les résultats intermédiaires dans des variables locales et en les vérifiant une par une. – Frank

Répondre

0

Les bonnes nouvelles sont qu'il n'y a rien de mal à votre mise en œuvre de DFT, et je suis en mesure de confirmer que, avec l'entrée appropriée, la fonction produit la sortie que vous attendez. Cela signifie cependant qu'il y a un problème avec les entrées qu'il reçoit. La mise en œuvre de la transformation inverse de façon similaire suggère que l'entrée suivante a été effectivement alimenté à la fonction DFT:

1+1i 
-1+1i 
1+0i 
-1+0i 
-3+0i 
4+0i 
0+0i 
2+0i 

De cela, il semble que vous manquez la partie imaginaire de la plupart des entrées (à l'exception de inputs[0] et inputs[1]). Comme vous n'avez pas indiqué comment les données de inputs sont générées et transmises à DFT, il est impossible de dire exactement pourquoi. Cependant, voici un certain nombre de causes communes qui peuvent conduire à une telle question, que vous pouvez utiliser comme points de départ alors qu'il enquêtait sur la raison pour laquelle la inputs à la fonction est pas ce que vous pensiez qu'il était:

  • Accidentellement (ou pas) écraser les parties imaginaires du tampon d'entrée avec des zéros. Cela pourrait à son tour être due à:
    • offrant le même tableau pour les inputs et outputs arguments à DFT (cette mise en œuvre, elles doivent être des tableaux distincts);
    • Utilisation de tableaux d'entrée et de sortie plus petits que la taille de transformation N et écriture ultérieure des limites lors de l'écriture de valeurs dans ces tableaux;
  • En effet, oublier de définir les parties imaginaires;
  • Obtenir de manière incorrecte des parties réelles/imaginaires en les lisant à partir d'un fichier ou d'autres parties du programme (par exemple en supposant qu'un fichier contient 8 parties réelles suivies de 8 parties imaginaires lorsque le fichier contient réellement 10 valeurs).