2fft.txt

matlab实现基2fft按时间抽取

源代码如下：

#include<math.h>
#include<stdio.h>
//#include<stdlib.h>
#include<windows.h>
#define N 1024
void FFT2(float x[],float y[],int n,int sign)
{ 
int i,j,k,l,m,n1,n2;
     float c,c1,e,s,s1,t,tr,ti;
    for(j=1,i=1;i<16;i++)
{
   m=i;
   j=j<<1;
   if(j==n)break;
}

//位序颠倒，二进制位换位
n1=n-1;
for(j=0,i=0;i<n1;i++)
{
   if(i<j)
   {
    tr=x[j];
    ti=y[j];
    x[j]=x[i];
    y[j]=y[i];
    x[i]=tr;
    y[i]=ti;
   }
   k=n>>1;
   while(k<(j+1))
   {
    j=j-k;
    k=k>>1;
   }
   j=j+k; 
}


   //FFT核心，不断将数据分为两部分，依次求解，降低了复杂度。
n1=1;
l=1;
while(l<=m)    //外循环，log2N
{//初始赋值。
     n1=n1<<1;
     n2=n1>>1;
     e=3.14159265359/n2;
     c=1;
     s=0;
     c1=cos(e);
     s1=-sign*sin(e);   
     j=0;
     while(j<n2)//for(j=0;j<n2;j++)//每个蝶形的两个输入数据相距B= 2^(L-1)个点
     {    for(i=j;i<n;i+=n1)//同一旋转因子对应着间隔为2^L点的2^(M-L)个蝶形
      {    k=i+n2;
       tr=c*x[k]-s*y[k];
       ti=c*y[k]+s*x[k];
       x[k]=x[i]-tr;
        y[k]=y[i]-ti;
       x[i]=x[i]+tr;
       y[i]=y[i]+ti;
              
      }
      t=c;
      c=c*c1-s*s1;   //三角递归公式
      s=t*s1+s*c1;
      j++;
     }
     l++;   
}

}
void main()
{   
    struct   high_resolution_timer   
    {   
   high_resolution_timer()   {::QueryPerformanceFrequency(&liFreq);}   
   void   start(){::QueryPerformanceCounter(&liStart);}   
   void   end(){::QueryPerformanceCounter(&liEnd);}   
   double get_duration(){return(liEnd.LowPart-liStart.LowPart)/(double)liFreq.QuadPart;}   
     private:   
    LARGE_INTEGER   liStart, liEnd, liFreq;   
    };  
   high_resolution_timer   hrt;   

  
//clock_t start, finish;   //计时参数
     //double duration; 
int n;
int i=0;
float x[N],y[N];
printf("-----------------时间抽取FFT2的计算机实现及验证-------------\n");
n=1024;    //用1024个数据进行测试，对程序运行时间进行测试
for(i=0;i<n;i++)
{
   x[i]=i+1;
   y[i]=0;

}
printf("--------------需要变换的数如下：-------------\n");
       for(i =0;i <n; i++)
    {
    printf("%f",x[i]);
    printf("+");
    printf("%f",y[i]);
    printf("i\n");  
    }
   //计时开始
hrt.start(); 
FFT2(x,y,n,1);//调用FFT2函数，利用原位运算，即输出结果还存放在输入数组中，可以节省存储单元，降低设备成本。
hrt.end(); 
     printf("————————傅里叶变换为——————-\n");
         for(i =0;i <n; i++)
    {
    printf("%f",x[i]);
    printf("+");
    printf("%f",y[i]);
    printf("i\n");  
    }
     FFT2 (x, y, n,-1);//利用原函数直接求逆函数，在三角变换中变换符号
     printf("————————逆傅里叶变换为——————-\n");
         for(i =0;i <n; i++)
   {
    printf("%f",x[i]/n);
    printf("+");
    printf("%f",y[i]/n);
    printf("i\n");  
   }
//输出计时时间        
     printf("%lf    s",    hrt.get_duration());   
     system("pause");  
}