时间抽选基2FFT及IFFT算法C语言实现.rar_WINDOWS 编程_fft_fft c语言_信号与系统_信号处理

/*时间抽选基2FFT及IFFT算法C语言实现*/ /*Author :Junyi Sun*/ /*Copyright 2004-2005*/ /*Mail:ccnusjy@yahoo.com.cn*/ #include <stdio.h> #include <math.h> #include <stdlib.h> #define N 1000 /*定义复数类型*/ typedef struct{ double real; double img; }complex; complex x[N], *W; /*输入序列,变换核*/ int size_x=0; /*输入序列的大小，在本程序中仅限2的次幂*/ double PI; /*圆周率*/ int main(){ int i,method; void fft(); /*快速傅里叶变换*/ void ifft(); void initW(); /*初始化变换核*/ void change(); /*变址*/ void add(complex a,complex b,complex *c); /*复数加法*/ void mul(complex a,complex b,complex *c); /*复数乘法*/ void sub(complex a,complex b,complex *c); /*复数减法*/ void divi(complex a,complex b,complex *c);/*复数除法*/ void output(); /*输出结果*/ system("cls"); PI=atan(1)*4; printf("Please input the size of x:\n"); scanf("%d",&size_x); printf("Please input the data in x[N]:\n"); for(i=0;i<size_x;i++) scanf("%lf%lf",&x.real,&x.img); initW(); printf("Use FFT(0) or IFFT(1)?\n"); scanf("%d",&method); if(method==0) fft(); else ifft(); output(); return 0; } /*快速傅里叶变换*/ void fft(){ int i=0,j=0,k=0,l=0; complex up,down; change(); for(i=0;i< (int)( log(size_x)/log(2) );i++){ /*一级蝶形运算*/ l=( 1<<i ); for(j=0;j<size_x;j+= (1<<l) ){ /*一组蝶形运算*/ for(k=0;k<l;k++){ /*一个蝶形运算*/ mul(x[j+k+l],W[size_x*k/2/l],&up); add(x[j+k],up,&up); mul(x[j+k+l],W[size_x*k/2/l],&down); sub(x[j+k],down,&down); x[j+k]=up; x[j+k+l]=down; } } } } /*快速傅里叶逆变换*/ void ifft(){ int i=0,j=0,k=0,l=size_x; complex up,down; for(i=0;i< (int)( log(size_x)/log(2) );i++){ /*一级蝶形运算*/ l/=2; for(j=0;j<size_x;j+= (1<<l) ){ /*一组蝶形运算*/ for(k=0;k<l;k++){ /*一个蝶形运算*/ add(x[j+k],x[j+k+l],&up); up.real/=2;up.img/=2; sub(x[j+k],x[j+k+l],&down); down.real/=2;down.img/=2; divi(down,W[size_x*k/2/l],&down); x[j+k]=up; x[j+k+l]=down; } } } change(); } /*初始化变换核*/ void initW(){ int i; W=(complex *)malloc(sizeof(complex) * size_x); for(i=0;i<size_x;i++){ W.real=cos(2*PI/size_x*i); W.img=-1*sin(2*PI/size_x*i); } } /*变址计算，将x(n)码位倒置*/ void change(){ complex temp; int i=0,j=0,k=0,t; for(i=0;i<size_x;i++){ k=i;j=0; t=(unsigned) (log(size_x)/log(2)); while(t--){