diff options
Diffstat (limited to 'NaiveFFT')
| -rw-r--r-- | NaiveFFT/main.c | 186 |
1 files changed, 185 insertions, 1 deletions
diff --git a/NaiveFFT/main.c b/NaiveFFT/main.c index c43c48b..fec893b 100644 --- a/NaiveFFT/main.c +++ b/NaiveFFT/main.c @@ -6,9 +6,193 @@ // #include <stdio.h> +#include <stdlib.h> +#include <unistd.h> +#include <string.h> +#include <math.h> +#include <complex.h> + +//data array +#define DATA_SIZE (8) +double data_i[DATA_SIZE] = {1.0f,0.0f,0.0f,0.0f,0.0f,0.0f,0.0f,0.0f}; +double data_q[DATA_SIZE] = {0.0f,0.0f,0.0f,0.0f,0.0f,0.0f,0.0f,0.0f}; +//calc the fft +void fft_if(double *x_i, double *x_q, int n, int inv) { + + int i=0,j=0,k=0,m=0, irem=0, sign; + double tq,ti; + + printf("FFT\n"); + k=n/2; + for (j=1,i=1;i<n; i++) { + if (i<j) { + + //swap values + ti = x_i[j]; + tq = x_q[j]; + + x_i[j] = x_i[i]; + x_q[j] = x_q[i]; + + x_i[i] = ti; + x_q[i] = tq; + printf("i=%d j=%d\n",i,j); + //find middle point + k = n/2; + while (k<j) { + j = j-k; + k=k/2; + } + } else { + printf("i=%d j=%d\n",i,j); + k = n/2; + while (k<j) { + j = j - k; + k = k/2; + } + } + j = j+k; + } + + for (i=0;i<DATA_SIZE;i++) { + printf("%.2f ",x_i[i]); + } + printf("\n"); + for (i=0;i<DATA_SIZE;i++) { + printf("%.2f ",x_q[i]); + } + printf("\n"); + + //calculate number of stages: + //m=log2(n) + m=0; + irem = n; + while(irem>1) { + irem = irem/2; + m = m+1; + } + + + //FFT or IFFT + if (inv==1) { + sign = 1; + } else { + sign = -1; + } + + // transform + for (i=0; i<m; i++) { + int le=pow(2,i); + int le1=le/2; + double ui=1.0; + double uq=0.0; + double wi = cos(M_PI/le1); + double wq = sign*sin(M_PI/le1); + + for (j=1; j<le1;j++) { + i=j; + while (i<n) { + int ip = i+le1; + + ti = x_i[ip]*ui - x_q[ip]*uq; + tq = x_q[ip]*ui + x_i[ip]*uq; + + x_i[ip] = x_i[i] - ti; + x_q[ip] = x_q[i] - tq; + + x_i[i] = x_i[i] + ti; + x_q[i] = x_q[i] + tq; + + i = i+le; + } + double temp = ui*wi - uq*wq; + uq = uq*wi + ui*wq; + ui = temp; + } + } + + //if inverse + if (inv==1) { + for (i=0;i<n;i++) { + x_i[i] = x_i[i]/n; + x_q[i] = x_q[i]/n; + } + } + +} + +//dft works fine +void dft(double *x_i, double *x_q, int n, int inv) { + double Wn,Wk; + //static array + //double Xi[DATA_SIZE],Xq[DATA_SIZE]; + //dynamic array + double *Xi, *Xq; + double c,s; + int i,j; + + Xi = malloc(n*sizeof(double)); + Xq = malloc(n*sizeof(double)); + + Wn = 2*M_PI/n; + + if (inv==1) { + Wn=-Wn; + } + + for (i=0;i<n;i++) { + Xi[i] = 0.0f; + Xq[i] = 0.0f; + + Wk = i*Wn; + for (j=0;j<n;j++) { + c = cos(j*Wk); + s = sin(j*Wk); + + //i - real, q - imaginary + Xi[i] = Xi[i] + x_i[j]*c + x_q[j]*s; + Xq[i] = Xq[i] - x_i[j]*s + x_q[j]*c; + } + + if (inv==1) { + Xi[i] = Xi[i]/n; + Xq[i] = Xq[i]/n; + } + } + + for (i=0;i<n;i++) { + x_i[i] = Xi[i]; + x_q[i] = Xq[i]; + } + +} int main(int argc, const char * argv[]) { + int i; + // insert code here... - printf("Hello, World!\n"); + printf("Hello FFT World!\n"); + + for (i=0;i<DATA_SIZE;i++) { + printf("%.2f ",data_i[i]); + } + printf("\n"); + for (i=0;i<DATA_SIZE;i++) { + printf("%.2f ",data_q[i]); + } + printf("\n"); + fft_if(data_i, data_q, DATA_SIZE, 0); + + //dft(data_i, data_q, DATA_SIZE,0); + + for (i=0;i<DATA_SIZE;i++) { + printf("%.2f ",data_i[i]); + } + printf("\n"); + for (i=0;i<DATA_SIZE;i++) { + printf("%.2f ",data_q[i]); + } + printf("\n"); + return 0; } |