关于DFT与FFT运算速度的比较Word下载.docx
《关于DFT与FFT运算速度的比较Word下载.docx》由会员分享，可在线阅读，更多相关《关于DFT与FFT运算速度的比较Word下载.docx（17页珍藏版）》请在冰豆网上搜索。
个抽样值
，这N个点的宽度为N的DFT为：
。
与之对应的反变换，即IDFT的定义：
是连续频率函数
，这N个点的宽度为N的IDFT为：
其中
称为N点DFT的变换核函数，
称为N点IDFT的变换核函数。
它们互为共轭。
如果引入
，正逆变换的核函数分别可以表示为
和
DFT可以表示为：
，IDFT可以表示为：
用FFT实现DFT：
先设序列点数为N=2M，M为整数。
如果不满足这个条件，可以人为地加上若干零值点，使之达到这一要求。
这种N为2的整数幂的FFT称基2FFT。
设输入序列长度为N=2M（M为正整数），将该序列按时间顺序的奇偶分解为越来越短的子序列，称为按时间抽取（DIT）的FFT算法。
下面给出8点基
2FFT的运算流图：
图1-1
下面对两种变化的运算量进行比较。
设x（n）为N项的复数序列，由DFT变换，任一X（m）的计算都需要N次复数乘法和N-1次复数加法，而一次复数乘法等于四次实数乘法和两次实数加法，一次复数加法等于两次实数加法，即使把一次复数乘法和一次复数加法定义成一次“运算”（四次实数乘法和四次实数加法），那么求出N项复数序列的X（m）,即N点DFT变换大约就需要N2次运算。
当N=1024点甚至更多的时候，需要N2=1048576次运算，在FFT中，利用WN的周期性和对称性，把一个N项序列（设N=2k,k为正整数），分为两个N/2项的子序列，每个N/2点DFT变换需要（N/2）2次运算，再用N次运算把两个N/2点的DFT变换组合成一个N点的DFT变换。
这样变换以后，总的运算次数就变成2*（N/2）2=N2/2。
继续上面的例子，N=1024时，总的运算次数就变成了525312次，节省了大约50%的运算量。
而如果我们将这种“一分为二”的思想不断进行下去，直到分成两两一组的DFT运算单元，那么N点的DFT变换就只需要Nlog2N次的运算，N在1024点时，运算量仅有10240次，是先前的直接算法的1%，点数越多，运算量的节约就越大，这就是FFT的优越性。
2.DFT与FFT等价性验证
如上面论述的，FFT只是DFT的快速实现算法，它仍是以DFT为理论基础的，因此他们两者完全等价。
这里将在VC6.0下编制DFT与FFT的程序，用其对两组离散序列进行变化，并输出变换结果。
矩形脉冲信号，序列长64，其中前16点为1，其余为0。
其DFT与FFT前5点输出如下图：
图2-1
正弦脉冲信号，序列长61，数字角频率为PI/2。
图2-2
比较可以发现用DFT和FFT对同一序列变化得到的结果是完全一样的。
这说明DFT与FFT是完全等价的，而区别仅在于其具体算法。
3.DFT与FFT运算速度比较
如上面介绍的，N点DFT需进行N2次运算，而N（N=2M）点FFT只需要Nlog2N次运算。
因此DFT的运算时间与FFT的运算时间比应为：
N/log2N。
这里将对上述结论在数量级层次进行验证。
为了严谨，这里先验证DFT与FFT运算时间仅与N有关，而与具体信号类型无关。
这里分别对128点、256点、512点矩形脉冲信号和正弦脉冲信号进行DFT与FFT变换，运算时间如图：
图3-1
图3-2
图3-3
可以看出，对不同类型信号进行DFT与FFT变换，只要其点数相同，其运算时间也是大致相同的。
考虑到C程序计时并不如汇编等机器语言准确，且100次重复计算会对每次运算时间之差进行累加，上面3次试验反映出的时间差是可以忽略的。
由此可见，DFT与FFT运算时间仅与运算点数N有关，而与具体信号类型无关。
当信号点数较多时，对不同类型信号进行DFT与FFT变化将花费大量时间，而上面已经验证DFT与FFT运算时间仅与N有关，而与具体信号类型无关。
故下面将仅对同一类型不同点数信号的DFT与FFT的运算时间的比较，而这也具有一般代表性。
下面将仅对不同点数正弦信号做DFT与FFT变化，并进行运算时间比较。
由于C程序中计时受机器滴答时间限制，故对点数较小信号进行较多次数DFT与FFT变换。
而当点数较多时，DFT运算时间将按N2规律增长，故对点数较多信号进行较少次数DFT与FFT变换。
各次DFT与FFT运算时间如下诸图：
图3-4
图3-5
图3-6
图3-7
图3-8
图3-9
图3-10
图3-11
图3-12
图3-13
图3-14
图3-15
图3-16
图3-17
下面对以上数据进行整理分析。
（注：
e+n表示*10+n）
N
变换次数
TDFT/s
TFFT/s
mean（TDFT）/s
mean（TFFT）/s
TDFT/TFFT
理论比值
2
10000
0.032
0.031
3.2e-6
3.1e-6
1.03
4
0.156
0.063
1.56e-5
6.3e-6
2.48
8
1000
0.062
0.016
6.2e-5
1.6e-5
3.875
2.67
16
0.25
0.047
2.5e-4
4.7e-5
5.32
32
1.016
0.109
1.016e-3
1.09e-4
9.32
6.4
64
100
0.406
4.06e-3
1.6e-4
25.38
10.67
128
1.703
1.703e-2
6.3e-4
27.03
18.29
256
6.609
0.141
6.609e-2
1.41e-3
46.87
512
26.218
0.313
2.6218e-1
3.13e-3
83.76
56.89
1024
10
10.531
1.0531
6.3e-3
167.16
102.4
2048
42.172
4.2172
1.56e-2
270.33
186.18
4096
167.047
0.344
1.67047e+1
3.44e-2
485.60
341.33
8192
1
67.141
0.078
6.7141e+1
7.8e-2
860.78
630.15
16384
267.485
2.67485e+2
1.56e-1
1714.65
1170.29
表3-1
在Matlab中对以上数据绘图如下：
图3-18
理论分析中已经指出，N点DFT运算次数将按N2规律增长，而N点FFT以Nlog2N规律增长，两者运算时间比值为N/log2N。
当N较小时，两者运算时间尚可比拟。
但随着N的增大，FFT相对DFT减少的运算量相对值就越大。
在N=1024时，运算次数已相差两个数量级，而当N=16384时，运算次数更是相差达三个数量级。
故信号序列点数越多，FFT相对DFT越有优势。
这一点在上述试验中得到了充分证实。
当N<
256时，一次DFT运算时间还较小，而当N>
1024后一次DFT运算时间就已经超过1秒。
事实上，在实际工程应用时不可能只进行一次DFT变化，因此即便在短序列分析中FFT的优势也是明显的。
上述实验中，DFT与FFT运算时间比与理论值有一定偏差。
其原因有以下两点：
一，程序与理论推导不同，程序难免在基于理论的同时存在一定的冗余，而两者冗余度的不同将直接影响其运算时间比；
二，本次试验的实验环境是VC6.0，其在计时方面本身就劣于汇编等机器语言。
尽管如此，在数量级层面试验反映的结果还是和理论分析一致的。
此外，由于N=2M，M为整数，DFT的运算时间以N2规律增长，即本实验DFT运算时间应以4为倍数增长，这一点在实验数据中有很好的体现。
4.总结
由上面的分析我们可以看出，快速傅里叶变化是离散傅里叶变化的一种高速算法，它的理论基础仍然是离散傅里叶变换，而并无任何创新。
但是由于快速傅里叶变化充分利用了离散傅里叶变换的奇、偶、虚、实等特性，并对离散傅立叶变换的算法进行改进，使得离散傅里叶变化运算次数由N2减少为Nlog2N次，使得离散傅里叶变换应用于实际变成现实。
参考文献
[1]AlanV.OppenheimAlanS.WillskywithS.HamidNawab.Signals&
Systems.Prentice-Hall,1997;
[2]AlanV.OppenheimRonaldW.ShaferwithJhonR.Buck.Discrete-TimeSignalProcessing.Prentice-Hall,1999;
[3]赵淑清.郑薇.随机信号分析.哈尔滨工业大学出版社.1999;
附录
程序
#include<
time.h>
stdlib.h>
stdio.h>
#defineN256
#defineM_PI3.14159
#defineEXP2.7182818288
voidFFT（float*xr,float*xi,float*yr,float*yi,intn,intinv）;
voidDFT（intn,float*x,intm,float*yr,float*yi）;
voidSignalRect（intn,intm,float*yr,float*yi）;
voidSignalCos（intn,intf,float*yr,float*yi）;
voidtra（float*x,float*y）;
main（）
{
/*以不同程序实现不同功能，具体见下文*/
}
voidSignalRect（intn,intm,float*yr,float*yi）
inti;
for（i=0;
i<
n;
i++）
yr[i]=yi[0]=0;
m;
yr[i]=1;
voidSignalCos（intn,intf,float*yr,float*yi）
yr[i]=cos（2*M_PI*f*i/n）;
voidDFT（intn,float*x,intm,float*yr,float*yi）
intj,k;
for（k=0;
k<
k++）
{
for（j=0;
j<
j++）
yr[k]=yr[k]+x[j]*cos（2*M_PI*k*j/m）;
yi[k]=yi[k]+x[j]*sin（2*M_PI*k*j/m）;
}
return;
voidFFT（float*xr,float*xi,float*yr,float*yi,intn,intinv）
inti,j,a,b,k,m;
intep,arg,mt,s0,s1;
floatsign,pr,pi,ph;
float*c,*s;
c=（float*）calloc（n,sizeof（float））;
if（c==NULL）exit
（1）;
s=（float*）calloc（n,sizeof（float））;
if（s==NULL）exit
（1）;
j=0;
if（inv==0）
sign=1.0;
yr[i]=xr[i]/n;
yi[i]=xi[i]/n;
elsesign=-1.0;
n-1;
if（i<
j）
tra（&
yr[i],&
yr[j]）;
yi[i],&
yi[j]）;
k=n/2;
while（k<
=j）
j=j-k;
k=k/2;
j=j+k;
ep=0;
i=n;
while（i!
=1）
ep=ep+1;
i=i/2;
ph=2*M_PI/n;
s[i]=sign*sin（ph*i）;
c[i]=cos（ph*i）;
a=2;
b=1;
for（mt=1;
mt<
=ep;
mt++）
s0=n/a;
s1=0;
b;
i=k;
while（i<
n）
arg=i+b;
if（k==0）
pr=yr[arg];
pi=yi[arg];
else
pr=yr[arg]*c[s1]-yi[arg]*s[s1];
pi=yr[arg]*s[s1]+yi[arg]*c[s1];
yr[arg]=yr[i]-pr;
yi[arg]=yi[i]-pi;
yr[i]=yr[i]+pr;
yi[i]=yi[i]+pi;
i=i+a;
s1=s1+s0;
a=2*a;
b=b*2;
free（c）;
free（s）;
voidtra（float*x,float*y）
floatt;
t=（*x）;
（*x）=（*y）;
（*y）=t;
验证DFT与FFT等价性主函数：
floatyr[N]={0.0},yi[N]={0.0},xr[N]={0.0},xi[N]={0.0};
inti;
SignalCos（N,16,xr,xi）;
/*产生64点cos（pi*i/2）离散信号*/
/*SignalRect（64,16,xr,xi）;
产生矩形脉冲*/
DFT（N,xr,N,yr,yi）;
printf（"
theoutput（byDFT,thefirst5）is"
）;
for（i=0;
5;
%f+%fj"
yr[i],yi[i]）;
\n"
FFT（xr,xi,yr,yi,N,0）;
theoutput（byFFT,thefirst5）is"
getch（）;
验证DFT与FFT运算时间仅与N有关主程序：
inti=100;
clock_tstart,end;
SignalCos（N,3,xr,xi）;
/*SignalCos（N,3,xr,xi）orSignalRect（64,16,xr,xi）;
*/
start=1*clock（）;
while（i--）
end=1*clock（）;
ThetimeforSignalCosbyFFT（N=%d,100times）was:
%fseconds\n"
N,（double）（end-start）/CLK_TCK）;
i=100;
ThetimeforSignalCosbyDFT（N=%d,100times）was:
SignalRect（64,16,xr,xi）;
ThetimeforSignalRectbyFFT（N=%d,100times）was:
ThetimeforSignalRectbyDFT（N=%d,100times）was:
对不同点数正弦信号DFT与FFT运算时间比较主程序（需对N，i赋值进行相应调整）：
inti=10000;
ThetimeforSignalCosbyFFT（N=%d,10000times）was:
i=10000;
ThetimeforSignalCosbyDFT（N=%d,10000times）was:
升级会员 