working dft

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;
 }