[Ayuda] Con una FFT en C

[LorElectronica]

Hola,
Me pueden ayudar con el siguiente código, no comprendo la función de int iK en esta FFT ya que los valores cuando me dan el resultado no son los imaginarios y reales que deberian ser entonces no se si estoy ejecutando mal el codigo o como es que deberia leerlo. Ademas debo hacerle la conjugada a los numeros ingresados en signal, esto para hacer la conjugada de la FFT (abs(fft)). Muchas Gracias de Ante Mano.

Código [Seleccionar] Expandir

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

#define M_PI 3.14159265358979323846264338327

void calcImag(double radians_signal[], double signal[], int k);
void calcReal(double radians_signal[], double signal[], int k);
double sum(double c[]);


double signal[] = {2.28025, 1.32888, 0.39326, -0.49619, -1.31121, -2.02672, -2.62174,
   -3.08015, -3.39124, -3.55077, -3.55763, -3.42069, -3.15151,
   -2.76733, -2.28963, -1.74326, -1.15541, -0.55456, 0.03068,
   0.57271, 1.04606, 1.42835, 1.7122, 1.85105, 1.86948, 1.75376,
   1.50688, 1.13742, 0.65924, 0.09094, -0.54489, -1.22254,
   -1.91419, -2.59102, -3.22433, -3.78672, -4.25312, -4.60188,
   -4.81556, -4.8817, -4.79336, -4.54939, -4.15456, -3.61943,
   -2.95996, -2.19702, -1.35554, -0.46368, 0.44824, 1.34888,
   2.20699, 2.99264, 3.6783, 4.23987, 4.65761, 4.91685, 5.00855,
   4.92967, 4.68323, 4.27822, 3.72929, 3.05615, 2.28286, 1.43692,
   0.54824, -0.35197, -1.23239, -2.06272, -2.81485, -3.46385,
   -3.9889, -4.37404, -4.60875, -4.68828, -4.61378, -4.39223,
   -4.03611, -3.56286, -2.99417, -2.3551, -1.67308, -0.97682,
   -0.29515, 0.3441, 0.91522, 1.3956, 1.76664, 2.01448, 2.13052,
   2.1118, 1.96104, 1.68661, 1.30214, 0.826, 0.28055, -0.3087,
   -0.91416, -1.50721, -2.05935, -2.54337, -2.93442, -3.21099,
   -3.35583, -3.35666, -3.20667, -2.90487, -2.45619, -1.87128,
   -1.16626, -0.36208, 0.51622, 1.44039, 2.38004, 3.30375,
   4.18023, 4.97949, 5.67398, 6.23957, 6.65652, 6.91015, 6.99145,
   6.89738, 6.631, 6.20138, 5.62319, 4.91623, 4.10463, 3.216};

// size of table
int iN = 128;
// harmonics
int iK = 4;

double imc[4][128];
double rec[4][128];

/*
*
*/
int main()
{
   double amplitude[4];
   double phase[4];
   double degreesSignal[128];
   double imcSum[4];
   double radiansSignal[128];
   double recSum[4];

   int i = 0;

   for (i = 0; i < iN; i++)
   {
       degreesSignal[i] = (360.0 / (double) iN) * (double) i;

       //printf("Degree Signal: %f \n", degreesSignal[i]);

       radiansSignal[i] = degreesSignal[i] * (M_PI / 180.0);

       //printf("Radians Signal: %f \n", radiansSignal[i]);
   }

   // get the real coefficients
   for (i = 0; i < iK; i++)
   {
       calcReal(radiansSignal, signal, i + 1);
   }

   // sum the real coefficients
   for (i = 0; i < iK; i++)
   {
       recSum[i] = sum(rec[i+1]);

       //printf("Real Sum: %f \n", recSum[i]);
   }

   // get the imag coefficients
   for (i = 0; i < iK; i++)
   {
       calcImag(radiansSignal, signal, i + 1);
   }

   // sum the imag coefficients
   for (i = 0; i < iK; i++)
   {
       imcSum[i] = sum(imc[i + 1]);

       //printf("Imaginary Sum: %f \n", recSum[i]);
   }

   // calculate amplitude (Fourier coefficients)
   // a[i] = 2(abs(RECSum[i]*abs(IMCSum[i])
   for (i = 0; i < iK; i++)
   {
       amplitude[i] = 2.0 * (sqrt(pow(recSum[i], 2.0) + pow(imcSum[i], 2.0)));

       //printf("Amplitude: %f \n", amplitude[i]);
   }

   // calculate phase
   // atan(abs(b/a))
   for (i = 0; i < iK; i++)
   {
       double temp = (sqrt(pow(imcSum[i], 2)
               / pow(recSum[i], 2)));
       phase[i] = atan(temp);
       
       //printf("Phase: %f \n", phase[i]);
   }

   for (i = 0; i < iK; i++)
   {
       printf("REC %d: %f\n", (i + 1), recSum[i]);
       printf("IMC  %d: %f\n", (i + 1), imcSum[i]);
       printf("Amplitude(K=%d): %f\n", (i + 1), amplitude[i]);
       printf("Phase(K=%d): %f\n", (i + 1), phase[i]);
       printf("Degrees(K=%d): %f\n", (i + 1), phase[i]*(180.0/M_PI));
       printf("\n");
   }

   return (EXIT_SUCCESS);
}

void calcImag(double radians_signal[], double signal[], int k)
{
   imc[k][0] = signal[0] * sin(-radians_signal[0]);

   for (int i = 1; i < iN; i++)
   {
       // for when the modular value is 0
       // prevents array out of bounds
       // IMC[array.length - 0]...
       // use the first sin value instead
       if ((i * k) % iN == 0)
       {
           imc[k][i] = signal[i] * sin(radians_signal[0]);
       }
       else
       {
           imc[k][i] = (signal[i] * sin(radians_signal[iN - (i * k) % iN]));
       }

       //printf("Imaginary Signal: %f \n", imc[k][i]);
   }
}

void calcReal(double radians_signal[], double signal[], int k)
{
   rec[k][0] = signal[0] * cos(radians_signal[0]);

   for (int i = 1; i < iN; i++)
   {
       // for when the modular value is 0
       // prevents array out of bounds
       // REC[array.length - 0]
       // use the first sin value instead
       if (((i * k) % iN) == 0)
       {
           rec[k][i] = signal[i] * cos(radians_signal[0]);
       }
       else
       {
           rec[k][i] = +(signal[i] * cos(radians_signal[iN - (i * k) % iN]));
       }

       //printf("REC: %d (n=%d) %f \n", k, i, rec[k][i]);
   }
}

double sum(double c[])
{
   double dCSum = 0;

   for (int i = 0; i < iN; i++)
   {
       dCSum += c[i];
   }

   dCSum = dCSum / iN;

   return dCSum;
}