按频率抽取的快速傅里叶变换
快速傅里叶变换
N X 2 ( k ) X 2 (k ) 2
第4章 快速傅里叶变换(FFT)
(N k ) 2 又由于WN
k WN WN
N 2
k WN
,所以
N N N k N 2 X (k ) X 1 (k ) WN X 2 (k ) 2 2 2
k X 1 (k ) WN X 2 (k ),
X 1 (k ) x1 (r )W x(2r )W
r 0 rk 4 r 0
3
3
rk 4
k 0,1,2,3
第4章 快速傅里叶变换(FFT)
(2) n为奇数时,分别记作:
x2 (0) x (1), x2 (1) x (3), x2 ( 2) x (5), x2 (3) x (7);
k N
1 1
k WN
-1
N X ( k ) X 1 (k ) WNk X 2 (k ) (后一半) 2
5.计算工作量分析
按奇、偶分组后的计算量:
第4章 快速傅里叶变换(FFT)
由上图可知,N点DFT的复乘为N2 ;复加N(N-1); 与分解后相比可知,计算工作点差不多减少 一半。
第4章 快速傅里叶变换(FFT)
一个X(k)的值的工作量,如X(1)
0 1 X (1) x(0)WN x(1)WN x(2)WN2 x( N 1)WNN 1
nk 通常x(n)和 W 都是复数,所以计算一个 N X(k)的值需要N次复数乘法运算,和N 1 次 复数加法运算.那么,所有的X(k)就要N2次复 数乘法运算,N(N-1)次复数加法运算.当N很 大时,运算量将是惊人的,如N=1024,则要完 成1048576 次(一百多万次)运算.这样,难 以做到实时处理.
快速傅里叶变换 FFT
因为
,所以:
上式中X1(k)和X2(k)分别为x2(r)和x2(r)的N/2点DFT, 即
由于X1(k)和X2(k)均以N/2为周期,且 X(k)又可表示为:
,以
即将一个N点的DFT分解成为两个N/2点的DFT。 上述运算可用右下图来表示,称为蝶形运算符号。
从右图可知,要完成一 个蝶形运算需要进行一 次复数相乘和两次复数 相加运算。
对于象雷达、通信、声纳等需要实时处理的信号, 因为其运算量更大,所以无法满足信号处理的实时 性要求。迫切需要有新的算法。
二、DFT运算的特点
实际上,DFT运算中包含有大量的重复运算。在WN
矩阵中,虽然其中有N2个元素,但由于WN的周期
性,其中只有N个独立的值,即
,且
这N个值也有一些对称关系。总之,WN因子具有如 下所述周期性及对称性:
N 1
X (k) x(n)WNkn n0
x(n)
1 N
N 1
X (k)WNkn
k 0
k 0,1,2, , N 1 n 0,1,2, , N 1
计算X(k)的运算量:需要N2次复数乘法,N(N-1)
次复数加法。在N较大时计算量很大。
例如:N=1024时, 需要1,048,576次复数乘法, 即 4,194,304次实数乘法
1.对称性
2.周期性 由上述特性还可得出:
利用上述对称特性,可使DFT运算中有些项可以合 并,这样,可使乘法次数减少大约一半;利用WN 矩阵的对称性及周期性,可以将长序列的DFT分解 为短序列的DFT,N越小,运算量能够减少。 例如,对于四点的DFT,直接计算需要16次复数乘 法,根据上述特性可以有以下形5年,J. W. Cooley和J. W. Tukey巧妙应用DFT中 W因子的周期性及对称性提出了最早的FFT,这是 基于时间抽取的FFT。具有里程碑式的贡献(运算量 缩短两个数量级)
《快速傅里叶变换》课件
FFT的历史背景
01
1960年代,Cooley和Tukey提 出了基于“分治”思想的FFT 算法,为快速傅里叶变换的实 用化奠定了基础。
02
随后,出现了多种FFT算法的 变种和优化,如Radix-2、 Radix-4等。
03
随着计算机技术的发展,FFT 算法在硬件实现上也得到了广 泛应用,如FPGA、GPU等。
《快速傅里叶变换》ppt课件
contents
目录
• FFT简介 • FFT基本原理 • FFT实现 • FFT的应用 • FFT的优化与改进 • FFT的挑战与未来发展
01 FFT简介
FFT的定义
快速傅里叶变换(FFT):一种高效计算离散傅里叶变换(DFT)及其逆变换的 算法。它将复杂度为$O(N^2)$的DFT计算降低到$O(Nlog N)$,大大提高了计 算效率。
详细描述
混合基数FFT算法结合了基数-2和基数-4算法的特点,利用两者在计算过程中的 互补性,减少了计算量,提高了计算效率。同时,该算法在处理大规模数据时 ,能够保持较高的精度。
分段FFT算法
总结词
分段FFT算法将输入数据分成若干段,对每一段进行快速傅里叶变换,以降低计算复杂度和提高计算效率。
详细描述
02 FFT基本原理
离散傅里叶变换(DFT)
定义
应用
DFT是时间域信号到频域的变换,通 过计算信号中各个频率成分的幅度和 相位,可以分析信号的频谱特性。
DFT在信号处理、图像处理、频谱分 析等领域有广泛应用。
计算量
DFT的计算量随着信号长度N的增加 而呈平方关系增长,因此对于长信号 ,计算量巨大。
05-第五章-快速傅里叶变换(蝶形运算)
W
r N
的确定
W
r N
的确定
以N=8为例:
m 1 时 W N r W , N j/4 W 2 j m W 2 0 ,j 0
m 2 时 W N r W , N j/2 W 2 j m W 4 j,j 0 , 1
m 3 时 W N r W N , j W 2 j m W 8 j,j 0 , 1 , 2 , 3 N2M,第L级:
8.0 1024 1 048 576 5 120 204.8
12.8 2048 4 194 304 11 264 372.4
64 4049 192
21.4
24
5.3.3 按时间抽取的FFT算法的特点
序列的逆序排列
同址运算(原位运算)
蝶形运算两节点间的距离
W
r N
的确定
25
序列的逆序排列
序列的逆序排列
以N/2点序列x1(r)为例
x1x (1 2 (l2 l)1 ) x3 x(4 l()l)
l0,1, ,N1 4
则有
N21
N4 1
N4 1
X1(k) x1(r)WNrk2 x1(2 l)W N 2l2k x1(2 l 1 )W N (2 2 l 1 )k
r0
l 0
l 0
N41
N41
x3(l)W N lk 4W N k2 x4(l)W N lk 4
FFT并不是一种与DFT不同的变换,而是 DFT的一种快速计算的算法。
3
5.2 直接计算DFT的问题及改进的途径
DFT的运算量
设复序列x(n) 长度为N点,其DFT为
N1
X(k) x(n)WNnk n0
k=0,,…,N-1
快速傅里叶变换发展史(可编辑)
快速傅里叶变换发展史快速傅立叶变换(FFT)个人日记2010-04-16 12:24:48 阅读163 评论0 字号:大中小订阅近十多年来数字信号处理技术同数字计算机、大规模集成电路等先进技术一样,有了突飞猛进的发展,日新月异,已经形成了一门具有强大生命力的技术科学。
由于它本身具有一系列的优点,所以能有效地促进各工程技术领域的技术改造和学科发展,应用领域也更加广泛、深入,越来越受到人们的重视。
在数字信号处理中,离散傅里叶变换(Discrete Fourier Transform,DFT)是常用的变换方法,它在各种数字信号处理系统中扮演着重要的角色。
傅里叶变换已有一百多年的历史了,我们知道频域分析常常比时域分析更优越,不仅简单,且易于分析复杂信号。
但用较精确的数字方法,即DFT进行谱分析,在FFT出现以前是不切实际的。
这是因为DFT计算量太大。
直到1965年出现了DFT]运算的一种快速方法以后,情况才发生了根本的变化。
快速傅里叶变换〔Fast Fourier Transfonn,FFT〕并不是与离散傅里叶变换不同的另一种变换,而是为了减少DFT计算次数的一种快速有效的算法。
当时Garwin在自己的研究中极需要一个计算傅立叶变换的快速方法,而正在写有关傅里叶变换的文章,Tukey概括地对Garwin介绍了一种方法,它实质上就是后来著名的Cooley-Tukey算法。
在Garwin的迫切要求下,1963年,IBM公司的Cooley根据Tukey的想法编写了第一个FFT算法程序。
在FFT算法中,Tukey主要利用了旋转因子的周期性和对称性。
这两个性质使DFT运算中的某些项可以合并,使DFT运算尽量分解为更少点数的DFT运算。
因为DFT的运算量与Pow(N,2)成比例,所以如果将一个大点数的DFT分解为若干个小点数的DFT 的组合,将有效地减少运算量。
Cooley在计算机上实现该算法时,为节省存储空间和减少寻址时间,采用了3维标号映射方法和在算法内部的循环结构,这些结构和技巧对后来的FFT算法研究及实现同样产生了很大影响。
《快速傅里叶变换(FFT) 第四章》
方法: 分解N为较小值:把序列分解为几个较短的 序列,分别计算其DFT值; 利用旋转因子WNk的周期性、对称性、可 约性进行合并、归类处理,以减少DFT的运 算次数。 k ( kn WN m WNN m WN ( nlN ) WNk lN ) n WN 周期性: N m m N m N m m m m 对称性:Wm WNm [W WN N WNN [WNNN m ]] WN WN 2 WN WN 可约性:W mN N W knmW kn / m W kn m kmn ,m 2 2
x ( r ) W x ( r )W x ( r ) W x ( r )W e (r W x r) xxr) W( r ) W (WW (r )W W e W (2 ) x x x(2 r 1)
W e
2 j 2 kr 2 kr N N /2
N 2
2 这样将N点DFT分解为两个N/2点的DFT
N X (k ) X 1 (k ) W X 2k(k ) k 0,1, 1 N X (k ) X 1 (k ) WN X 2 (k ) k 0,1, N 1 2 k X (kN X 1 (k ) WN X 2 (k ) k 0,1, 2 1 ) N2 k X (k N X 1 (k ) WN X 2 (k ) k 0,1, 1 N ) k X (k 2 N X 1 (k ) WN X 2 (k ) k 0,1, N 1 ) 2 k 2 X (k ) X (k ) W X (k ) k 0,1, 2 1
4.1 离散傅里叶变换的高效计算思路 DFT是信号分析与处理中的一种重要变换。但直接 计算DFT的计算量与变换区间长度N的平方成正比, 当N较大时,计算量太大,直接用DFT算法进行谱分 析和信号的实时处理是不切实际的。
快速傅里叶变换
快速傅⾥叶变换快速傅⾥叶变换快速傅⾥叶变换(FFT )是根据计算量的最⼩化原理来设计和实施离散傅⾥叶变换(DFT)计算的⽅法。
1965年,库利(T.W.Cooley )和图基(J.W.tukey )发表了著名的《计算机计算傅⾥叶级数的⼀种算法》论⽂。
从此掀起了快速傅⾥叶变换计算⽅法研究的热潮。
快速傅⾥叶变换(FFT )的出现,实现了快速、⾼效的信号分析和信号处理,为离散傅⾥叶变换(DFT)的⼴泛应⽤奠定了基础。
1.1离散傅⾥叶变换(DFT)的计算设x(n)是⼀个长度为M 的有限长序列,则定义x(n)的N 点离散傅⾥叶变换为∑-===10)()]([)(N n kn NW n x n x DFT k X 其中由于计算⼀个X(k)值需要N 次复乘法和(N-1)次复数加法,因⽽计算N 个X(k)值,共需N2次复乘法和N(N-1)次复加法。
每次复乘法包括4次实数乘法和2次实数加法,每次复加法包括2次实数加法,因此计算N 点的DFT 共需要4N2次实数乘法和(2N2+2N ·(N-1))次实数加法。
当N 很⼤时,这是⼀个⾮常⼤的计算量。
1.2减少DFT 计算量的⽅法减少DFT 的计算量的主要途径是利⽤k N W 的性质和计算表达式的组合使⽤,其本质是减少DFT 计算的点数N 以便减少DFT 的计算量。
k N W 的性质:(1)对称性: (2)周期性: (3) 可约性: (4) 特殊点:选择其中⼀个证明N N j k N j N k N j N k N e e e W 222)2(22πππ--+-+==ππj k N j e e --=2k N j e π2--=k N W -=FFT 算法是基于可以将⼀个长度为N 的序列的离散傅⾥叶变换逐次分解为较短的离散傅⾥叶变换来计算这⼀基本原理的。
这⼀原理产⽣了许多不同的算法,但它们在计算速度上均取得了⼤致相当的改善。
0,1,,1k N =-()*nk nk N N W W -=()()nk N n k n N k N N NW W W ++==nk mnk N mN W W =//nk nk m N N mW W =01N W =/21N N W =-(/2)k N k N NW W +=-在这⾥讨论两类基本的FFT 算法。
快速傅里叶变换(含详细实验过程分析)
[实验2] 快速傅里叶变换 (FFT) 实现一、实验目的1、掌握FFT 算法和卷积运算的基本原理;2、掌握用C 语言编写DSP 程序的方法;3、了解利用FFT 算法在数字信号处理中的应用。
二、实验设备 1. 一台装有CCS 软件的计算机; 2. DSP 实验箱的TMS320C5410主控板; 3. DSP 硬件仿真器。
三、实验原理 (一)快速傅里叶变换傅里叶变换是一种将信号从时域变换到频域的变换形式,是信号处理的重要分析工具。
离散傅里叶变换(DFT )是傅里叶变换在离散系统中的表示形式。
但是DFT 的计算量非常大, FFT 就是DFT 的一种快速算法, FFT 将DFT 的N 2 步运算减少至 ( N/2 )log 2N 步。
离散信号x(n)的傅里叶变换可以表示为∑=-=10][)(N N nk N W n x k X , Nj N e W /2π-=式中的W N 称为蝶形因子,利用它的对称性和周期性可以减少运算量。
一般而言,FFT 算法分为时间抽取(DIT )和频率抽取(DIF )两大类。
两者的区别是蝶形因子出现的位置不同,前者中蝶形因子出现在输入端,后者中出现在输出端。
本实验以时间抽取方法为例。
时间抽取FFT 是将N 点输入序列x(n) 按照偶数项和奇数项分解为偶序列和奇序列。
偶序列为:x(0), x(2), x(4),…, x(N-2);奇序列为:x(1), x(3), x(5),…, x(N-1)。
这样x(n) 的N 点DFT 可写成:()()∑++∑=-=+-=12/0)12(12/02122)(N n kn NN n nkNW n x Wn x k X考虑到W N 的性质,即2/)2//(22/)2(2][N N j N j N W e e W ===--ππ因此有:()()∑++∑=-=-=12/02/12/02/122)(N n nkN k NN n nkN W n x WWn x k X或者写成:()()12()kN X k X k W X k =+由于X 1(k) 与X 2(k) 的周期为N/2,并且利用W N 的对称性和周期性,即:k N N k N W W -=+2/可得:()()12(/2)kN X k N X k W X k +=-对X 1(k) 与X 2(k)继续以同样的方式分解下去,就可以使一个N 点的DFT 最终用一组2点的DFT 来计算。
程佩青《数字信号处理教程》(第4版)(课后习题详解 快速傅里叶变换(FFT))
4.2 课后习题详解4-1 如果一台通用计算机的速度为平均每次复乘40ns ,每次复加5ns ,用它来计算512点的DFT[x (n )],问直接计算需要多少时问?用FFT 运算需要多少时间?若做128点快速卷积运算,问最低抽样频率应是多少?解:①直接利用DFT 计算:复乘次数为N 2,复加次数为N (N-1)。
②利用FFT计算:复乘次数为,复加次数为N㏒2N 。
(1)直接计算复乘所需时间复加所需时间所以(2)用FFT 计算复乘所需时间复加所需时间所以4-2 N =16时,画出基-2按频率抽选法的FFT 流图采用输入自然顺序,输出倒位序),统计所需乘法次数(乘±1,乘±j 都不计在内)。
根据任一种流图确定序列x (n )=4cos (n π/2)(0≤n ≤15)的DFT 。
解:按频率抽取法的FFT 流图中的复数乘法出现在减法之后,其运算量为复数乘法:;复数加法:;由于N =16,有,,,不需要乘法。
按频率抽取,见图4-1(a )。
图4-1(a )运算量:复数乘法:由于,,,不需要乘法。
由图P4.2(a )可知,共有的个数为1+2+4+8=15有的个数为1+2+4=7所以总的乘法次数为32-15-7=10(个)复数加法:举例:对序列x (n )=4cos (n π/2)(0≤n ≤15)可表示为由于N =16,可采用P4.2(b )的流图。
设Xi (k )=(i =1,2,3,4)分别为第i 级蝶形结构的输出序列,则由P4.2(b )的流图可知由于采用的是顺序输入、逆序输出的结构,因此输出X (k )与X 4(k )为逆序关系,即,为k 二进制逆序值由此可知,x (n )的DFT 为X (4)=X 4(2)=32,X (12)=X 4(3)=12图4-1(b )4-3 用MATLAB 或C 语言编制以下几个子程序。
(1)蝶形结运算子程序;(2)求二进制倒位序子程序;(3)基-2 DIT FFT 流程图,即迭代次数计算子程序。
第5章 快速傅立叶变换
2
N 2
N 2
−1
r =0
2
第 5章 快速傅立叶变换
式中 X 1 (k ) 和 X 2 (k ) 分别是 N / 2 点DFT,得到的只是X (k ) 的 , 前一半项数的结果。 前一半项数的结果。 X 性质可知, 由DFT性质可知, 1 (k ) 、 2 (k ) 都是周期为N / 2 的周期函数 性质可知 X
x3 (0) = x1 (0) = x(0)
x3 (1) = x1 (2) = x(4)
X 3 (0)
N/4点 点 DFT X 3 (1)
X 1 (0)
0 N /2
x4 (0) = x1 (1) = x(2)
x4 (1) = x1 (3) = x(6)
N/4点 点 X 4 (1) DFT
W X 4 (0)
W
k N
-1
X(
(后一半) 后一半)
例如讨论N=8的蝶形流图表示方法 的蝶形流图表示方法 例如讨论 x(2r ) = x1 (r ) , r = 0,1, 2,3 x(2r + 1) = x2 (r )
X (k ) = X 1 (k ) + W8k X 2 (k ), k = 0,1, 2,3 X (k + 4) = X 1 (k ) − W8k X 2 (k ), k = 0,1, 2,3
lk N /4
= ∑ x6 (l )W
l =0
N , k = 0,1,..., − 1 4
2 将系数统一为: k 将系数统一为: N / 2 = WN k 则N=8点DFT就可以分解为四 W N=8点DFT就可以分解为四
个N/4=2点DFT,就可以得到如下流程图 点 ,
频率抽取FFT
W
3 N
x[0]
2点
DFT
0 WN
X[0] X[4] X[2] X[6] X[1]
x[1]
x[2]
-1
0 N
x[3]
x[4] x[5] x[6]
W
W
-1
-1
-1
2 N
2点 DFT
2点 DFT
W W
1 N
2 N 3 N
X[5]
X[3] X[7]
W
-1 -1
-1
-1
0 N
x[7]
W
W
2 N
2点 DFT
基4时间抽取FFT算法推导
且 W
k+N / 4 N
jW
k N
W
k+2N / 4 N
1W
k N
k+3N / 4 k k+4N / 4 k WN jWN WN WN
k 2k 3k X(k ) X 1 (k ) WN X 2 (k ) WN X 3 (k ) WN X 4 (k )
将长度为N的序列分解为N/4的短序列,可进一步提
高计算效率。 当长序列分解的短序列越短时,其分解的短序列越多, 相应的合成运算越复杂。 需要综合平衡分解与合成过程中的综合计算复杂度!
基4时间抽取FFT算法
将长度为N的序列x(n)分解为4个短序列
x1 (n) x(4n), x2 (n) x(4n 1), x3 (n) x(4n 2),
1、原位运算(同址运算)
X m1 (k )
X m1 ( j )
1
X m (k ) X m1 (k ) X m1 ( j )
r WNr X m ( j ) [ X m1 (k ) X m1 ( j )]WN
快速傅里叶变换
r 0,1,
, N 2 1
N / 2 1
r 0
x(2r )WN
2 rk
N / 2 1
r 0
x(2r 1)WN
(2 r 1) k
N / 2 1
r 0
rk k x(2r )WN W /2 N
N / 2 1
r 0
rk x(2r 1)WN /2
N 2 N N ( N 1) N 2 所以8点FFT求X (k )共需 36次复数乘法, 2 2 2 2 N N2 N ( 1) N 32次复数加法, 共计为68次。看出仅做一次 2 2 分解就可以节省约一半计算量。
(2)N/2(4点)-->N/4(2点)FFT
N/2点 DFT
X1(1)
X1(2)
X1(3) X2(0)
WN
0
X(0) X(1) X(2) X(3)
1 2
x2(r) x2(0)=x(1)
奇 数 序 列
x2(1)=x(3) x2(2)=x(5) x2(3)=x(7)
N/2点 DFT
X2(1)
X2(2)
WN
WN
X2(3)
3 WN
X(4) -1 X(5) -1 X(6) -1 X(7)
-1 X(4)~X(7)
如:X(0) X 1 (0) X 2 (0)W80 X(2) X 1 (2) X 2 (2)W82
X( 1) X 1 (1) X 2 (1)W81
X(3) X 1 (3) X 2 (3)W83
同学们自已写
比较直接计算N=8点DFT 与分解2个4点DFT的
DSP--FFT-深入浅出-详细讲解快速傅里叶变换
第一节 引言
一、迅速付里叶变换FFT
• 有限长序列经过离散傅里叶变换 (DFT)将其频 域离散化成有限长序列.但其计算量太大(与N 旳平方成正比), 极难 实时地处理问题 , 因 此 引 出 了 快 速 傅 里 叶 变 换(FFT) .
• 一种复数乘法涉及4个实数乘法和2个实数相 法。
(a+jb)(c+jd)=(ac-bd)+j(bc+ad)
2次实数加法
4次实数乘法
4.计算DFT需要旳实数运算量
N 1
X (k) {(Re[x(n)]Re[WNkn ] Im[x(n)]Im[WNkn ]) n0
j(Re[x(n) Im[WNkn ] Im[x(n)]Re[WNkn ])}
4
4
X (k) N X (k) N
(
N
)
4 2
4
+
(
N 4
4
)2
=
N 4
2
这么一直分下去,剩余两点旳变换。
2、将长序列DFT利用对称性和 周期性分解为短序列DFT--结论
• 迅速付里时变换(FFT)就是在此特征基础上 发展起来旳,并产生了多种FFT算法,其基 本上可提成两大类:
• 按抽取措施分: 时间抽取法(DIT);频率抽取法(DIF)
r 0
r 0
W 2
j 2 2
e N
j 2
e N/2
W
3.求出子序列旳DFT
上式得:
N / 21
N / 21
X(k)
x1(r)WNrk/ 2
x2 (r)WNrk/ 2WNk
FFT原理
2
X 1 (k )
同理:X 2 ( k
N 2
) X 2 (k )
再利用W系数的对称性:
WN
(k
N ) 2
k WN
X (k
N
2 2 可见,一个N点的DFT被分解为前后两个N/2点的DFT,
) X 1 k WN X 2 k
k
k 0,1,
N
1
这两个N/2点的DFT再合成为一个N点DFT。
0 wN
0 wN 2 wN1
w11 N
2 w1( N 1) wN (N 1) N
x (0) (N 1) 1 x (1) wN (N 1)(N 1) wN x ( N 1)
(2)原位运算结构(同址运算)
FFT运算结构很有规律,它具有强烈的对称性, 它的所有运算都可以由左下图所示的基本运算构成。 由于运算结构规律,所以硬件实现起来简单,软件 编程方便。 FFT的蝶形运算结构很经济,它可以采用原位运 算。原位运算—— 指当数据输入到存储器中以后, 每一级运算的结果仍然储存在同 一组存储器中,直到最后输出, 中间无需其它任何存储器,称为 原位运算。
2
rk wN
2
N / 21
r 0
x (2r ) WNrk
2
N / 2 1
r 0
k x(2r 1) wN WNrk
2
故 : k X
N / 2 1
r 0 2
x1 (r ) W rk WNk
N
N / 2 1
r 0 2
N 1
}
快速傅里叶变换
前半部分
k0,1, ,N21
X ( N /2 k ) X 1 ( k ) W N k X 2 ( k )
后半部分
蝶形运算流图符号
说明:
X1(k)
X1(k)W N kX2(k)
(1) 左边两路为输入 (2) 右边两路为输出
W
k N
X2(k)
X1(k)W N kX2(k)
(3) 中间以一个小圆表示加、
X(k)X1(k)W N kX2(k)
k0, ,N/21
X(kN/2)X1(k)W N kX2(k)
得到X1(k)和X2(k)需要: 复乘:(N/2)2+ (N/2)2次 = 32次 复加:N/2(N/2-1)+N/2(N/2-1) =12+12 =24次
此外,还有4个蝶形结,每个蝶形结需要1次复乘,2次复加。 一共是:复乘4次,复加8次。
W
N N
k
W
N
n
k
W
n N
N
W
N
n
k
周 期 性 W N n k W N ( N n ) k W N n ( N k )
可 约 性 W N n k W m m N n k WNnk WNnk/m /m
j 2 m nk
e mN
j2N
e N2
e j
1
特 殊 点 : W N 0 1 W N N / 2 1 W N ( k N / 2 ) W N k
偶序列 奇序 (l0 ..列 N 4 . 1 )此l处 0 ,1 ,
那么,X1(k)又可表示为
N /4 1
N /4 1
X 1 (k ) x 1 (2 l)W N 2 l/2 k x 1 (2 l 1 )W N (2 /l2 1 )k
快速傅立叶变换(FFT)算法_DSP实验
快速傅立叶变换(FFT)算法实验摘要:FFT(Fast Fourier Transformation),即为快速傅里叶变换,是离散傅里叶变换的快速算法,它是根据离散傅里叶变换的奇、偶、虚、实等特性,对离散傅立叶变换的算法进行改进获得的。
这种算法大大减少了变换中的运算量,使得其在数字信号处理中有了广泛的运用。
本实验主要要求掌握在CCS环境下用窗函数法设计FFT快速傅里叶的原理和方法;并且熟悉FFT快速傅里叶特性;以及通过本次试验了解各种窗函数对快速傅里叶特性的影响等。
引言:快速傅里叶变换FFT是离散傅里叶变换DFT的一种快速算法。
起初DFT的计算在数字信号处理中就非常有用,但由于计算量太大,即使采用计算机也很难对问题进行实时处理,所以并没有得到真正的运用。
1965年J.W.库利和T.W.图基提出快速傅里叶变换,采用这种算法能使计算机计算离散傅里叶变换所需要的乘法次数大为减少,特别是被变换的抽样点数N越多,FFT算法计算量的节省就越显著。
从此,对快速傅里叶变换(FFT)算法的研究便不断深入,数字信号处理这门新兴学科也随FFT的出现和发展而迅速发展。
根据对序列分解与选取方法的不同而产生了FFT的多种算法,基本算法是基2DIT和基2DIF。
FFT 的出现,使信号分析从时域分析向频域分析成为可能,极大地推动了信号分析在各领域的实际应用。
FFT在离散傅里叶反变换、线性卷积和线性相关等方面也有重要应用。
一、 实验原理:FFT 并不是一种新的变换,它是离散傅立叶变换(DFT )的一种快速算法。
由于我们在计算DFT 时一次复数乘法需用四次实数乘法和二次实数加法;一次复数加法则需二次实数加法。
每运算一个X (k )需要4N 次复数乘法及2N+2(N-1)=2(2N-1)次实数加法。
所以整个DFT 运算总共需要4N^2次实数乘法和N*2(2N-1)=2N(2N-1)次实数加法。
如此一来,计算时乘法次数和加法次数都是和N^2成正比的,当N 很大时,运算量是可观的,因而需要改进对DFT 的算法减少运算速度。
7-快速傅里叶变换
将N点序列分解成2个N/2点序列
x(0) x(0)
x(3) x(5)
将N/2点DFT分为两个N/4点DFT
2点DFT
x(0) + +
x(4)
+ -
8点碟2FFT
x(0) x(0)
5 3
按频率抽取(DIF)
x1[n]=x[n]
x2[n]=x[n+N/2]
n=0,1,…,N/2-1
按频率抽取(DIF)
x[2n 1]W
N / 2 1 n 0
( 2 n 1) k N
N / 2 1
nk k x1[n]W N / 2 W N
nk x 2[n]W N / 2
k X 1(k ) W N X 2(k )
上式的最后一步是因为
W [e
2 N
2 j( ) 2 N
]
e
N
4 16 16
2
N 2 - N ( N / 2) log 2 N N log N 2
12 240 4 32 8 64
256
64 256 1024
4096 65536 1048576
4032 65280 1047552
192 1024 5120
384 2048 10240
将N点序列分解成2个N/2点序列
X (k ) X 1(k N / 2) W X 2(k N / 2)
k N k X 1(k ) W N X 2(k )
这是因为:
W
(N、 k 2) N
W
( N / 2) N
W W
k Nห้องสมุดไป่ตู้
k N
将N点序列分解成2个N/2点序列
快速傅里叶变换_蝶形运算_按频率抽取基2-fft算法_MATLAB代码
function y=MyFFT_FB(x,n)%MYFFT_TB:My Fast Fourier Transform Frequency Based%按频率抽取基2-fft算法%input:% x -- 输入的一维样本% n -- 变换长度,缺省时 n=length(x) 当n小于x数据长度时,x数据被截断到第n个数据% 当n大于时,x数据在尾部补0直到 x 含n个数据%output:% y -- 1*n的向量,快速傅里叶变换结果%variable define:% N -- 一维数据x的长度% xtem -- 临时储存x数据用% m,M -- 对N进行分解 N=2^m*M,M为不能被2整除的整数% two_m -- 2^m% adr -- 变址,1*N的向量% l -- 当前蝶形运算的级数% W -- 长为 N/2的向量,记录 W(0,N),W(1,N),...W(N/2-1,N)% d -- 蝶形运算两点间距离% t -- 第l级蝶形运算含有的奇偶数组的个数% mul -- 标量,乘数% ind1,ind2 -- 标量,下标% tem -- 标量,用于临时储存%参考文献:% 81c 输入参数个数检查msg=nargchk(1,2,nargin);error(msg);%% 输入数据截断或加0N=length(x);if nargin==2if N<n % 加0xtem=x;x=zeros(1,n);x(1:N)=xtem;N=n;else % 截断xtem=x;x=xtem(1:n);N=n;endend%% 对N进行分解 N=2^m*M[m,M]=factorize(N);two_m=N/M;%% 变换if m~=0%% 如果N可以被2整除adr=address(m,M,two_m);y=x; % 蝶形运算级数 l=m 时%% 计算W向量W=exp(-2*pi*i* ( 0:N/2-1 ) /N);%% 蝶形运算d=N/2;t=1;for l=1:m% 加for ii=0:t-1ind1=ii*2*d+1;ind2=ind1+d;for r=0:d-1tem=y(ind1)+y(ind2);y(ind2)=y(ind1)-y(ind2);y(ind1)=tem;ind1=ind1+1;ind2=ind2+1;endend% 乘for r=0:d-1mul=W(r*t+1);for ii=0:t-1y(ii*2*d+d+1+r) = y(ii*2*d+d+1+r)*mul;endendd=d/2;t=t*2;end%% 直接傅里叶变换if M~=1 % N 分解含有非2因数M时,对y中每M个数据做直接傅里叶变换 for ii=1:two_my((ii-1)*M+1 : ii*M ) = DDFT( y((ii-1)*M+1 : ii*M ) );end%% 变址输出y=y(adr+1);else%% 如果N 不能被2整除y=DDFT(x);endend%% 内嵌函数 ====================================================== function y=DDFT(x)%% 直接离散傅里叶变换%input:% x -- 样本数据,N维向量%output:% y -- N维向量%参考文献:% 结构动力学,克拉夫,P82% variable define% s -- sum,用于求和N=length(x);y=zeros(size(x));for n=1:Ns=0;for m=1:Ns=s+x(m)*exp( -i*2*pi*(m-1)*(n-1)/N );endy(n)=s;endendfunction [m,M]=factorize(N)%% 对N分解m=0;while trueif mod(N,2)==0m=m+1;N=N/2;elsebreak;endendendfunction adr=address(m,M,two_m)%% 变址% b -- 2^m * m 的矩阵,用来存储二进制数据% ds -- 数,公差adr=zeros(two_m,M);b=de2bi(0:two_m-1,m);%转换为2进制注:matlab中二进制[0 1 1]=6 b=b(:,end:-1:1);% 逆序adr(:,1)=bi2de(b);%2进制转换为10进制if M~=1ds=two_m;adr=adr(:,1)*ones(1,M);adr=adr+ds*ones(size(adr,1),1)*(0:M-1);adr=reshape(adr',1,[]);endend。
快速傅里叶变换2
WNk H
(k)
后N/2点
12
3、计算量分析: 按奇、偶分组后的计算量
▪ 一个N/2点的DFT运算量G(k) 或 H:(k)
复乘次数: ( N )2 复N加2 次数:
2
4
▪ 两个N/2点的DFT运算量:
N ( N 1) 22
复乘次数: N 2 复加次数: ▪ 一个蝶形运算 2
N ( N 1) 2
DIT FFT属于原位运算
流图中输入‘乱序’ ,输出‘顺序’。
21
‘乱序’原因
n =0 1
n0=0(偶) n1=1
x(n2, n1, n00)
n1=0 n0=1 (奇)
n1=1
(n2)
0 x(000) 0 1 x (100) 4 0 x (010) 2 1 x (110) 6 0 x (001) 1 1 x (101) 5
0 x (011) 3 1 x (111) 7
22
‘整序’规 律 将输入序号按自然顺序排序后,用相应位数的 二进制码表示,在进行反序,即可实现输入端的乱序。
自然顺序n 二进制n2n1n0 反序二进制n0n1n2 倒位顺序n
0
000
000
0
1
001
100
4
2
010
010
2
3
011
110
6
4
100
001
x2(2l) x5(l) x2(2l 1) x6(l)
l
0,1,L
,
N 4
1
l
0,1,L
,
N 4
1
16
x1 (r ) 的DFT:
x1(r) 分解为
x1(2l) x3(l),
快速傅里叶变换fft变换
快速傅里叶变换FFT的C语言算法彻底研究LED音乐频谱显示的核心算法就是快速傅里叶变换,FFT的理解和编程还是比较难的,特地撰写此文分享一下研究成果。
一、彻底理解傅里叶变换快速傅里叶变换(Fast Fourier Transform)是离散傅里叶变换的一种快速算法,简称FFT,通过FFT可以将一个信号从时域变换到频域。
模拟信号经过A/D转换变为数字信号的过程称为采样。
为保证采样后信号的频谱形状不失真,采样频率必须大于信号中最高频率成分的2倍,这称之为采样定理。
假设采样频率为fs,采样点数为N,那么FFT结果就是一个N点的复数,每一个点就对应着一个频率点,某一点n(n 从1开始)表示的频率为:fn=(n-1)*fs/N。
举例说明:用1kHz的采样频率采样128点,则FFT结果的128个数据即对应的频率点分别是0,1k/128,2k/128,3k/128,…,127k/128 Hz。
这个频率点的幅值为:该点复数的模值除以N/2(n=1时是直流分量,其幅值是该点的模值除以N)。
二、傅里叶变换的C语言编程1、对于快速傅里叶变换FFT,第一个要解决的问题就是码位倒序。
假设一个N 点的输入序列,那么它的序号二进制数位数就是t=log2N.码位倒序要解决两个问题:①将t位二进制数倒序;②将倒序后的两个存储单元进行交换。
如果输入序列的自然顺序号i用二进制数表示,例如若最大序号为15,即用4位就可表示n3n2n1n0,则其倒序后j对应的二进制数就是n0n1n2n3,那么怎样才能实现倒序呢?利用C语言的移位功能!程序如下,我不多说,看不懂者智商一定在180以下!复数类型定义及其运算#define N 64 //64点#define log2N 6 //log2N=6/*复数类型*/typedef struct{float real;float img;}complex;complex xdata x[N]; //输入序列/*复数加法*/complex add(complex a,complex b){complex c;c.real=a.real+b.real;c.img=a.img+b.img;return c;}/*复数减法*/complex sub(complex a,complex b){complex c;c.real=a.real-b.real;c.img=a.img-b.img;return c;}/*复数乘法*/complex mul(complex a,complex b){complex c;c.real=a.real*b.real - a.img*b.img;c.img=a.real*b.img + a.img*b.real;return c;}/***码位倒序函数***/void Reverse(void){unsigned int i,j,k;unsigned int t;complex temp;//临时交换变量for(i=0;i<N;i++)//从第0个序号到第N-1个序号{k=i;//当前第i个序号j=0;//存储倒序后的序号,先初始化为0for(t=0;t<log2N;t++)//共移位t次,其中log2N是事先宏定义算好的{j<<=1;j|=(k&1);//j左移一位然后加上k的最低位k>>=1;//k右移一位,次低位变为最低位}if(j>i)//如果倒序后大于原序数,就将两个存储单元进行交换(判断j>i是为了防止重复交换){temp=x[i];x[i]=x[j];x[j]=temp;}}}2、第二个要解决的问题就是蝶形运算①第1级(第1列)每个蝶形的两节点“距离”为1,第2级每个蝶形的两节点“距离”为2,第3级每个蝶形的两节点“距离”为4,第4级每个蝶形的两节点“距离”为8。
- 1、下载文档前请自行甄别文档内容的完整性,平台不提供额外的编辑、内容补充、找答案等附加服务。
- 2、"仅部分预览"的文档,不可在线预览部分如存在完整性等问题,可反馈申请退款(可完整预览的文档不适用该条件!)。
- 3、如文档侵犯您的权益,请联系客服反馈,我们会尽快为您处理(人工客服工作时间:9:00-18:30)。
《数字信号处理》课程设计报告按频率抽取的DFT快速算法分析及MATLAB实现专业:通信工程班级:组次:姓名:学号:目录摘要 (1)关键字 (1)0 引言 (1)1 按频率抽取的DFT快速算法原理 (1)2 DIF-FFT的运算规律及编程思想 (2)2.1 原位计算 (2)2.2 序列的倒序 (2)2.3 旋转因子的变换规律 (2)2.4 蝶形运算规律 (4)2.5 编程思想及程序框图 (4)3 DIF-FFT算法运算量分析 (5)4 MATLAB程序实现 (5)5 结束语 (7)参考文献 (7)按频率抽取的DFT 快速算法分析及MATLAB 实现摘要:DFT 是数字信号分析与处理中的一种重要变换。
但直接计算DFT 的计算量与变换区间长度N 的平方成正比,计算量非常大。
DFT 的快速算法使运算效率提高了1~2个数量级,为数字信号处理技术应用于各种信号的实时处理创造了条件。
为了对FFT 有更加深入的了解,本文对DIF-FFT 的原理进行了分析,并给出MATLAB 程序实现的方法与步骤。
关键词:DFT;DIF-FFT;MATLAB;0 引言DFT 是数字信号分析与处理中的一种重要变换。
但直接计算DFT 的计算量与变换区间长度N 的平方成正比,计算量非常大。
DFT 的快速算法使运算效率提高了1~2个数量级,为数字信号处理技术应用于各种信号的实时处理创造了条件。
本文通过对按频率抽取的DFT 快速算法原理介绍与MATLAB 实现以期使我们对傅里叶快速算法有更全面的理解,为我们以后更复杂的快速算法学习打下基础。
1 按频率抽取的DFT 快速算法原理设序列x(n)的长度为M N 2=,将序列前后对半分开,得到两个子序列,如下:式中:k kN NW )1(2/-=将x(k)分解成偶数组与奇数组,当k 取偶数(k=2m,m=0,1,…,N/2-1)时:∑∑-=-=++=++=12/02/12/02)]2()([)]2()([)2(N n mnN N n mn N W N n x n x W N n x n x m x (1)当k 取奇数(k=2m+1,m=0,1,…,N/2-1)时,∑∑-=-=+⋅+-=+-=+12/02/12/0)12()]2()([)]2()([)12(N n nmN n N N n m n N W W N n x n x W N n x n x m x (2)nk N N n Nk N kN n NN n nkN N n N N n nkN N n nk NN n nk NW W N n x n x W N n x W n x W n x Wn x Wn x k X ∑∑∑∑∑∑-=⎪⎭⎫ ⎝⎛+-=-=-=-=-=⎥⎦⎤⎢⎣⎡⎪⎭⎫ ⎝⎛++=⎪⎭⎫ ⎝⎛++=+==1202/212012012101202)(2)()()()()(令:⎭⎬⎫++=+-=)2()()()]2()([)(12Nn x n x n x W N n x n x n x n N 其中, n=0,1,2,…,N/2-1 将)()(21n x n x 和分别代入(1)、(2)式,可得:⎪⎭⎪⎬⎫∑∑-=-===+12/02/112/02/2)()2()()12(N n mnN N n nmN W n x m X W n x m X (3)(3)式表明,X(k)按奇偶k 值分为两组,其偶数组是)(1n x 的N/2点DFT ,奇数组则是)(2n x 的N/2点DFT 。
)(1n x 、)(2n x 和x (n )之间的关系可以用图1所示的蝶形运算流图符号表示。
图2表示N=8的DIF-FFT 运算流图。
图1 DTF-FFT 蝶形运算流图符号图2 DIF-FFT 的运算流图(N=8)-1x (n )x (n +N / 2)nN W x (n )+x (n +N / 2)[x (n )-x (n +N / 2)]n N W -1-10NW 2N W x (0)x (1)x (2)x (3)-1-1x (4)x (5)x (6)x (7)0NW 1N W 2N W 3N W X (0)X (4)X (2)X (6)X (1)X (5)X (3)X (7)N W 2N W -1-1-1-1-1-1-1-10NW 0N W 0NW 0NW2 DIF-FFT 的运算规律及编程思想 2.1 原位计算M N 2=点的FFT 共进行M 级运算,每级由N/2个蝶形运算组成。
同一级中,每个蝶形的两个输入数据只对计算本蝶形有用,而且每个蝶形的输入、输出数据结点又同在一条水平线上,这就意味着计算完一个蝶形后,所得输出数据可立即存入原输入数据所占用的存储单元。
这样,经过M 级运算后,原来存放输入序列数据的N 个存储单元中便依次存放X (k )的N 个值。
原位计算可节省大量内存,从而使设备陈本降低。
2.2 序列的倒序由图2可知,DIF-FFT 算法输入序列为自然序列,而输出为倒序排列。
因此M 级运算完后,要对输出数据进行倒序才能得到自然顺序的X(k)。
图3为顺序与倒序二进制对照图。
图3 顺序与倒序二进制对照图2.3 旋转因子的变换规律N 点的DFT 快速傅里叶运算流图中,每级都有N/2个蝶形。
每个蝶形都要乘以因子PN W ,称其为旋转因子,P 为旋转因子的指数。
但各级的旋转因子和循环方式都有所不同。
为了编写计算程序,下面列出旋转因子P N W 与运算级数的关系。
用L 表示从左到右的运算级数(L=1,2,…,M ),第L 级共有12-L 个不同的旋转因子。
顺序倒序十进制数I二进制数 二进制数 十进制数J0 000 000 0 1 001 100 4 2 010 010 2 3 011 110 6 4 100 001 1 5 101 101 5 6 110 011 3 71111117对N=M 2的一般情况,第L 级的旋转因子为:J PN L W W 2= J=0,1,2,…,12-L -1因为 M L M L M LN --⋅=⨯=2222所以 LM M L J NJ N P N W W W --⋅⋅==22 J=0,1,2,12-L -1L M J P-⋅=22.4 蝶形运算规律对M N 2=点FFT ,输入倒位序,输出自然序,设第L 级运算每个蝶形的两节点距离为B行,则第L 级运算:{)()()(*)]()([)(1111B J A J A J A W B J A J A B J A L L L PNL L L ++⇐+-⇐+----2.5 编程思想及程序框图观察图2可以归纳出一些对编程有用的运算规律:第L 级中,每个蝶形的两个输入 数据相距B=12-L 个点;每级有B 个不同的因子;同一旋转因子对应着间隔为L2点的LM -2个蝶形。
频率抽取法输入为自然顺序,输出为倒序。
图4为大致流程图。
图5 为DIF-FFT 运算程序框图6倒序程序框图图4 大致流程图计算x 的长度n,不到2的数幂,补0开始读入x(n)蝶形运算序列倒序输出 结束DIF-FFT运算程序框图倒序程序框图图 5 图63 DIF-FFT算法运算量分析按频率抽取的 FFT算法是将频域信号序列X(K)分解为奇偶两部分,但算法仍是由时域信号序列开始逐级运算,把 N点分成N/2点计算FFT,可以把直接计算离散傅里叶变换所需的2N次乘法缩减到次。
4 MATLAB程序实现clc;clear all;close all;%关闭程序,清屏x=[1,2,3,4,5,6,7];x1=x;%本程序对输入序列实现DIF-FFT基2算法,点数取大于等于长度的2的幂次m=nextpow2(length(x)); %求的x长度对应的2的最低幂次mN=2^m;if length(x)<Nx=[x,zeros(1,N-length(x))]; %若的长度不是2的幂,补0到2的整数幂endfor L=m:-1:1 %将DFT做m次基2分解,从左到右,对每次分解作DFT运算D=2^L;u=1; %旋转因子u初始化WN=exp(-1i*2*pi/D); %本次分解的基本DFT因子WN=exp(-i*2*pi/D)for j=1:D/2 %本次跨越间隔内的各次蝶形运算for k=j:D:N %本次蝶形运算的跨越间隔为Dkp=k+D/2; %确定蝶形运算的对应单元下标temp=x(k); %保存x(k)的值x(k)=x(k)+x(kp); %加法运算x(kp)=(temp-x(kp))*u; %乘法运算endu=u*WN; %修改旋转因子,多乘一个基本DFT因子WNendendnxd=bin2dec(fliplr(dec2bin([1:N]-1,m)))+1; %求1:2^m数列的倒序disp('自编程序结果:');y=x(nxd)disp('系统自带函数结果:');y1=fft(x1,N)自编程序结果:y =Columns 1 through 528.0000 -9.6569 + 4.0000i -4.0000 - 4.0000i 1.6569 - 4.0000i 4.0000 Columns 6 through 81.6569 + 4.0000i -4.0000 + 4.0000i -9.6569 - 4.0000i系统自带函数结果:y1 =Columns 1 through 528.0000 -9.6569 + 4.0000i -4.0000 - 4.0000i 1.6569 - 4.0000i 4.0000Columns 6 through 81.6569 + 4.0000i -4.0000 + 4.0000i -9.6569 - 4.0000i5 结束语通过这次课程设计,能够提高我独立思考,解决学习问题的能力,并且重新温习了DIT-FFT运算,自学了DIF-FFT运算,对以前学过的知识掌握得更加牢固,同时也增加了MATLAB编程的信心。
参考文献[1] 唐向宏,岳恒立,郑雪峰.MATLAB及在电子信息类课程中的应用(第2版)[M].北京:电子工业出版社, 2009.6[2] 高西全,丁玉美.数字信号处理(第三版)[M].西安:西安电子科技大学出版社,2008。