Transformée de Fourier discret donnant des résultats incorrects

CAVEAT: EM> Je suis un peu rouillé sur ce

Comme je me souviens de cela, la partie COS / Real donne le spectre de fréquence et la partie SIN / IMAG donne le spectre de puissance. Donc, ce que je pense que vous voulez imprimer, c'est le spectre de la fréquence (qui est juste outraal [i] code>) plutôt que ce que vous avez fait (c'est-à-dire ajouter nombre inférieur code> et et indimag code> est faux). Vous pouvez parcourir les deux couleurs différentes, mais je commencerais simplement. P>

aussi, je commencerais avec une entrée de données plus simple. P>

J'ai changé le théorique code> faire une seule vague de cosinus [ainsi que votre original]. Ceci est prévisible et vous devriez avoir une seule pointe, que vous le faites, donc je pense que dftriqueal code> est correct. P>

J'ai ajouté des options de ligne de commande afin que vous puissiez expérimenter. p>

J'ai aussi trouvé que diviser i / linecount code> est parfois tronqué, donc j'ai créé un frag code> macro pour illustrer / corriger ceci. Sans cela, la pointe d'une entrée de COS (2 * pi * t) code> a fini dans outraal [0] code> (partie DC?) Au lieu de Outral [ 1] code> p>

En outre, à des fins de test, vous pouvez obtenir avec beaucoup moins d'échantillons (par exemple 1000). Cela peut également aider à étendre votre intrigue horizontalement pour être mieux capable de voir les choses. P>

Quoi qu'il en soit, voici le code [Veuillez pardonner le nettoyage de style gratuit]: P>

#include <math.h>
#include <stdio.h>
#include <stdlib.h>
#include <string.h>
#include <time.h>

//#define M_PI 3.14159265358979323846
#define M_2PI (M_PI * 2.0)

int opt_A;                              // algorithm to use
int linecount;                          // number of samples
int opt_a;                              // use absolute value on output
int opt_s;                              // use squared value on output
int opt_i;                              // use integer output index
int opt_j;                              // output imaginary part
int opt_c;                              // .csv compatibility

time_t todlast;

// the first was your original and would truncate
#if 0
#define FRAG(_i)    ((_i) / linecount)
#else
#define FRAG(_i)    ((double) (_i) / linecount)
#endif

void
pgr(int k,int n,int count)
{
    time_t todnow = time(NULL);

    if ((todnow - todlast) >= 1) {
        todlast = todnow;
        printf("\r%d %d ",count,k);
        fflush(stdout);
    }
}

// the test wavefunction
double
theoretical(double t)
{
    double a;

    switch (opt_A) {
    case 0:
        a = 0.0;
        a += sin(M_PI * t);
        a += 2.0 * sin(2.0 * M_PI * t);
        a += 4.0 * sin(4.0 * M_PI * t);
        break;
    default:
        a = cos(M_2PI * t * opt_A);
        break;
    }

    return a;
}

//-------------------------------------------------------------------------
void
dftreal(double inreal[], double outreal[], double outimag[], int linecount)
{
    int n;
    int k;

    for (k = 0; k < linecount; k++) {
        double sumreal = 0;
        double sumimag = 0;
        double omega_k = (M_2PI * k) / (double) linecount;

        for (n = 0; n < linecount; n++) {
            // omega k:
#if 0
            double angle = M_2PI * n * FRAG(k);
#else
            double angle = omega_k * n;
#endif

            sumreal += inreal[n] * cos(angle);
            sumimag += inreal[n] * sin(angle);

            pgr(k,n,linecount);
        }

        outreal[k] = sumreal;
        outimag[k] = sumimag;
    }
}

//=========================================================================
int
main(int argc,char **argv)
{
    char *cp;

    --argc;
    ++argv;

    for (;  argc > 0;  --argc, ++argv) {
        cp = *argv;
        if (*cp != '-')
            break;

        switch (cp[1]) {
        case 'a':  // absolute value
            opt_a = ! opt_a;
            break;
        case 's':  // square the output value
            opt_s = ! opt_s;
            break;
        case 'c':  // .csv output
            opt_c = ! opt_c;
            break;
        case 'i':  // integer output
            opt_i = ! opt_i;
            break;
        case 'j':  // imaginary output
            opt_j = ! opt_j;
            break;
        case 'A':
            opt_A = atoi(cp + 2);
            break;
        case 'N':
            linecount = atoi(cp + 2);
            break;
        }
    }

    if (linecount <= 0)
        linecount = 44100 / 2;

    // creates all necessary arrays
    double inreal[linecount];
    double outreal[linecount];
    double outimag[linecount];
#if 0
    double p[linecount];
#endif

    FILE *fout = fopen("Output.txt", "w");

    for (int i = 0; i < linecount; ++i)
        inreal[i] = theoretical(FRAG(i));

    todlast = time(NULL);

    // actually computes the transform
    dftreal(inreal, outreal, outimag, linecount);

    for (int i = 0; i < linecount; ++i) {
#if 0
        p[i] = 2 * (outreal[i] * outreal[i] + outimag[i] * outimag[i]);
        fprintf(fout, "%f %f\n", (i / (double) linecount), p[i]);
#endif

#if 0
        fprintf(fout, "%f %f %f\n", (i / (double) linecount),
            outreal[i] * outreal[i],
            outimag[i] * outimag[i]);
#endif

#if 0
        fprintf(fout, "%f %f\n",
            (i / (double) linecount),outreal[i] * outreal[i]);
#endif

        if (opt_a) {
            outreal[i] = fabs(outreal[i]);
            outimag[i] = fabs(outimag[i]);
        }

        if (opt_s) {
            outreal[i] *= outreal[i];
            outimag[i] *= outimag[i];
        }

        if (! opt_c) {
            if (opt_i)
                fprintf(fout, "%d ",i);
            else
                fprintf(fout, "%f ",FRAG(i));
        }

        fprintf(fout, "%f",outreal[i]);
        if (opt_j)
            fprintf(fout, "%f",outimag[i]);
        fprintf(fout, "\n");
    }
    fclose(fout);

    printf("\nEnd of program\n");

    //getchar();
    return 0;
}

La bonne nouvelle est qu'il n'y a rien de mal à votre implémentation dftreal .

Le problème, si on pouvait l'appeler, est que la forme d'onde de test que vous utilisez comprend des composants de fréquence qui sont des composants de fréquence qui sont de très basses fréquences relatives à votre taux d'échantillonnage linocount . En conséquence, le DFT montre que l'énergie est concentrée dans les premiers bacs, avec des fuites spectrales dans des bacs supérieurs car les composants de fréquence de forme d'onde ne sont pas des multiples exacts de la fréquence d'échantillonnage.

Si vous augmentez le test. Les fréquences de forme d'onde en faisant la fréquence relative à la fréquence d'échantillonnage avec par exemple: xxx

Vous devez obtenir un parcelle telle que: