FFT的算法

FFT是一种DFT的高效算法，称为快速傅立叶变换（fast Fourier transform），它根据离散傅氏变换的奇、偶、虚、实等特性，对离散傅立叶变换的算法进行改进获得的。FFT算法可分为按时间抽取算法和按频率抽取算法，先简要介绍FFT的基本原理。从DFT运算开始，说明FFT的基本原理。

DFT的运算为：

式中

由这种方法计算DFT对于 的每个K值，需要进行4N次实数相乘和（4N-2）次相加，对于N个k值，共需4N*4N次实数相乘和（4N-2)(4N-2）次实数相加。改进DFT算法，减小它的运算量，利用DFT中 的周期性和对称性，使整个DFT的计算变成一系列迭代运算，可大幅度提高运算过程和运算量，这就是FFT的基本思想。

FFT对傅氏变换的理论并没有新的发现，但是对于在计算机系统或者说数字系统中应用离散傅立叶变换，可以说是进了一大步。

设x(n)为N项的复数序列，由DFT变换，任一X（m）的计算都需要N次复数乘法和N-1次复数加法，而一次复数乘法等于四次实数乘法和两次实数加法，一次复数加法等于两次实数加法，即使把一次复数乘法和一次复数加法定义成一次“运算”（四次实数乘法和四次实数加法），那么求出N项复数序列的X（m）,即N点DFT变换大约就需要N^2次运算。当N=1024点甚至更多的时候，需要N2=1048576次运算，在FFT中，利用WN的周期性和对称性，把一个N项序列（设N=2k,k为正整数），分为两个N/2项的子序列，每个N/2点DFT变换需要（N/2）2次运算，再用N次运算把两个N/2点的DFT变换组合成一个N点的DFT变换。这样变换以后，总的运算次数就变成N+2*（N/2)^2=N+（N^2）/2。继续上面的例子，N=1024时，总的运算次数就变成了525312次，节省了大约50%的运算量。而如果我们将这种“一分为二”的思想不断进行下去，直到分成两两一组的DFT运算单元，那么N点的DFT变换就只需要Nlog2N次的运算，N在1024点时，运算量仅有10240次，是先前的直接算法的1%，点数越多，运算量的节约就越大，这就是FFT的优越性。

fft的概念最好看书，这个细细分析一下还是能够理解的。你要是没有相关的书的话，这个网站上的也可以看看：http://blog.csdn.net/sshcx/archive/2007/06/14/1651616.aspx

c程序我以前写过一个，注释得还是挺清楚的，你看看吧

(我这里写的是1024个点，你简单修改为64个点就行。程序思路比较简单清晰，建议你自己看懂吸收了比较好，这才是学习)

/**************************************************************/

#include<stdio.h>

#include<math.h>

#ifndef pi

#define pi 3.14159265

#endif

int N=1024//采集的数据点数，注意为2的n次方

int n=10 //2^n=N

double input[1024] //原始输入信号

double IO[1024] //经过重新排序的信号（实数）

int code[1024] //0～（2^n－1）的倒置码自然数组

struct complex c_IO[1024] //经过重新排序的信号（复数）

struct complex out[1024] //输出变换的结果（复数）

struct complex temp[1024] //保存中间变换结果（复数）

void fft(struct complex *c_IO,int m) //fft函数,输入c_IO[]为经过了码位转换的信号数据；m＝n-1表示需要分解的层数

void trans_code(int n ) //产生0～（2^n－1）的倒置码自然数组

struct complex complex_mult(struct complex a, struct complex b) //结构体复数相“乘”函数，输入为两个结构体类型的复数

struct complex complex_add(struct complex a, struct complex b)//结构体复数相“加”函数，输入为两个结构体类型的复数

struct complex complex_remove(struct complex a, struct complex b) //结构体复数相“减”函数，输入为两个结构体类型的复数

struct complex W_n_k (int n,int k) //产生权值的函数W(n,k)=exp(-j*2*k*pi/N)

struct complex doule2complex (double input) //把实数转换为结构体复数的函数

void main()

{

FILE *fp

///////////////////////////////// 获取输入信号，这里是手动（程序）输入，也可以简单修改为 “读入文本文件”

int k=0

for(k=0k<Nk++)

{

input[k]=sin(pi*k/3)+cos(pi*k/20)

}

////////////////////////////////// 把输入数据进行码位倒置

for(k=0k<Nk++)

{

code[k]=0

}

trans_code(n)

for(k=0k<Nk++)

{

IO[k]=input[ code[k] ]

}

///////////////////////////////// 把输入数据转化为复数，为下面的FFT做准备

for(k=0k<Nk++)

{

c_IO[k]=doule2complex(IO[k])

}

///////////////////////////////// 进行FFT变换，结果保存在全局变量 out[N]中

fft(c_IO, n-1)

for(k=0k<Nk++)

{

printf("%f +%f j\n",out[k].x,out[k].y)

}

////////////////////////////////////////////////把结果输入到文本文件中去，其中实部输入到dx.txt文件中，虚部输入到dy.txt中

fp=fopen("dx.txt","w")

for(k=0k<Nk++)

{

fprintf(fp,"%f\n ",out[k].x)

}

fclose(fp)

fp=fopen("dy.txt","w")

for(k=0k<Nk++)

{

fprintf(fp,"%f\n ",out[k].y)

}

fclose(fp)

///////////////////////////////////////////////////

} //main函数截至

////////////////////////////////////////////

////////////////////////////////////////////

// 以下是子函数的实现部分 //

////////////////////////////////////////////

////////////////////////////////////////////

void fft(struct complex c_IO[],int m)

{

int group,i,p,q,k

if(m==-1)

{

for(k=0k<Nk++)

{

out[k]=c_IO[k]

}

}

else

{

fft(c_IO, m-1) //递归调用

group=(int)( N/pow(2,m+1) )//每一级的分组数

for(i=0i<groupi++)

{

for(p=(int)( i*pow(2,m+1) ) p<=(int)( i*pow(2,m+1)+pow(2,m)-1 ) p++ )

{

q=p+(int)pow(2,m)

temp[p]=complex_add(out[p],complex_mult(W_n_k( (int)pow(2,m+1),(int)(p-i*pow(2,m+1)) ),out[q] ) )

temp[q]=complex_remove(out[p],complex_mult(W_n_k( (int)pow(2,m+1),(int)( p-i*pow(2,m+1) ) ),out[q]) )

}

}

for(k=0k<Nk++)

{

out[k]=temp[k]

}

}

}

void trans_code(int n ) //change code order from 1 to 2^n

{

int bit[20]={0}

int i=0

int j=0

int k=0

int mytemp=0

int num

num=(int)pow(2,n)

for(j=0j<numj++)

{

mytemp=j

////////////////////////////////////////////////////////////

for(i=0i<ni++) //十进制数二进制化

{

bit[i]=mytemp%2

mytemp=mytemp/2

}

////////////////////////////////////////////////////////////

for(k=0k<(n/2)k++) //码位倒置

{

mytemp=bit[k]

bit[k]=bit[n-1-k]

bit[n-1-k]=mytemp

}

///////////////////////////////////////////////////////////

for(i=0i<ni++) //恢复为十进制数

{

code[j]=code[j]+(int)pow(2,i)*bit[i]

}

}

}

struct complex complex_mult(struct complex a, struct complex b)

{

struct complex c

c.x=a.x*b.x-a.y*b.y

c.y=a.x*b.y+a.y*b.x

return c

}

struct complex complex_add(struct complex a, struct complex b)

{

struct complex c

c.x=a.x+b.x

c.y=b.y+a.y

return c

}

struct complex complex_remove(struct complex a, struct complex b)

{

struct complex c

c.x=a.x-b.x

c.y=a.y-b.y

return c

}

struct complex W_n_k (int n,int k)

{

struct complex c

c.x=cos(2*pi/n*k)

c.y=-sin(2*pi/n*k)

return c

}

struct complex doule2complex (double input)

{

int i=0

struct complex c

c.x=input

c.y=0

return c

}

FFT（快速傅里叶变换）是一种实现DFT（离散傅里叶变换）的快速算法，是利用复数形式的离散傅里叶变换来计算实数形式的离散傅里叶变换，matlab中的fft()函数是实现该算法的实现。

MATLAB它将数值分析、矩阵计算、科学数据可视化以及非线性动态系统的建模和仿真等诸多强大功能集成在一个易于使用的视窗环境中，为科学研究、工程设计以及必须进行有效数值计算的众多科学领域提供了一种全面的解决方案，并在很大程度上摆脱了传统非交互式程序设计语言（如C、Fortran）的编辑模式，代表了当今国际科学计算软件的先进水平。

快速傅里叶变换, 即利用计算机计算离散傅里叶变换（DFT)的高效、快速计算方法的统称，简称FFT。快速傅里叶变换是1965年由J.W.库利和T.W.图基提出的。采用这种算法能使计算机计算离散傅里叶变换所需要的乘法次数大为减少，特别是被变换的抽样点数N越多，FFT算法计算量的节省就越显著。

扩展资料：

matlab优势特点：

1、高效的数值计算及符号计算功能，能使用户从繁杂的数学运算分析中解脱出来；

2、具有完备的图形处理功能，实现计算结果和编程的可视化；

3、友好的用户界面及接近数学表达式的自然化语言，使学者易于学习和掌握；

4、功能丰富的应用工具箱(如信号处理工具箱、通信工具箱等) ，为用户提供了大量方便实用的处理工具。

参考资料来源：

百度百科-快速傅里叶变换

百度百科-MATLAB

---