Discrete Fourier Transform in C

Complex Numbers in C

• Support from ISO C99
  
• 3.0 + 4.0 * I
  
• Lib: require complex.h
  

types

• float complex
  
• double complex
  
• long double complex
  

PI Number
pi radian = 180°
tan(45°) = 1
arctan(1) = 45°
in C
printf(“%f”,atan2(1,1)); ⇒» 0.7853

 #include <stdio.h>  
 #include <stdlib.h>  
 #include <math.h>  
 #include <complex.h>  
 float pinum;  
 typedef double complex comp;  
 int main(void)  
 {  
  FILE *sin_tab;  
  sin_tab = fopen("sin_tab.txt","w+") ;  
  FILE *dft_tab;  
  dft_tab = fopen("dft_tab.txt","w+") ;  
  FILE *mag_tab;  
  mag_tab = fopen("mag_tab.txt","w+");  
  pinum = atan2(1.0,1.0) * 4 ;  
  int Fs = 1000; // Sampling Rate(Hz)  
  double time_interval = (double)1/Fs; //Sampling Period.. i.e. 1ms  
  int data_len = 1000; // Input data length  
  comp sin_arr[data_len],fft_out[data_len];  
  comp buff = 0 + 0 * I;  
  int f1 = 111;  
  int f2 = 300; // sine frequencies  
 // Generate input test signals  
  for(int L=0;L<data_len;L++)  
      {  
           sin_arr[L] = (10000 * sin(2 * pinum * f1 * L * time_interval ) + (2000 * sin(2 * pinum * f2 * L * time_interval )));  
           fprintf(sin_tab,"%g\n",creal(sin_arr[L]));  
      }  
 // Get DFT    
  for(int k=0;k<data_len;k+=1)  
      {  
           for(int n=0;n<data_len;n+=1)  
           {  
                buff += sin_arr[n] * cexp(-2 * I * pinum * n * k / data_len) ;  
           }  
           fft_out[k]= buff ;  
           buff = 0 + 0 * I ;  
          mag[k] = sqrt(pow(creal(fft_out[k]),2.0) + pow(cimag(fft_out[k]),2.0) ) ;  
           fprintf(dft_tab,"%d----%g + j %g ve mag : %lf\r\n",k,creal(fft_out[k]),cimag(fft_out[k]),mag[k]);            
      }  
 // Write Magnitudes to a file   
  for(int k=0;k<(data_len/2);k+=1)  
      {  
           fprintf(mag_tab,"%lf\n",mag[k]);  
      }   
  return 0;  
  }   

Sonuç: Kodda decimation in time yapılmamıştır. Yani giriş sinyali cooley-tukey alg. olduğu gibi radix-2 olarak tek ve çift olarak ayrılmamıştır. Bu yüzden DFT kodu diyebiliriz. N boyutlu giriş için N*N tane kompleks çarpma yapıyor. Bu da hızı düşürmektedir. Kod 3 fazdan oluşuyor. 1-Giriş dizisi belirleniyor(Burada sinüs).2-DFT alınıyor; Giriş dizisiyle aynı boyutta.3-Sonuçlar txt dosyasına kaydediliyor.

-----------------------
A Good document: http://www.robots.ox.ac.uk/~sjrob/Teaching/SP/l7.pdf