按时间抽取的基2FFT算法分析

1、第四章快速傅里叶变换有限长序列可以通过离散傅里叶变换(DFT)将其频域也离散化成有限长序列但其计算量太大，很难实时地处理问题，因此引出了快速傅里叶变换(FFT). 1965年，Cooley和Tukey提出了计算离散傅里叶变换( DFT)的快 速算法，将DFT的运算量减少了几个数量级。从此，对快速傅里叶变换(FFT)算法的研究便不断深入，数字信号处理这门新兴学科也随FFT的出现和发展而迅速发展。根据对序列分解与选取方法的不同而产生了FFT的多种算法，基本算法是基2 DIT和基2 DIF。 FFT在离散傅里叶反变换、线性卷积 和线性相关等方面也有重要应用。快速傅里叶变换(FFT)是计算离散傅里叶变

2、换(DFT )的快速算法。DFT的定义式为NX(k)八 x(n)WRN(k)n=0在所有复指数值 W,n的值全部已算好的情况下，要计算一个X(k)需要N次复数乘法和N 1次复数加法。算出全部 N点X(k)共需N 2次复数乘法 和N(N -1)次复数加法。即计算量是与 N2成正比的。FFT的基本思想：将大点数的DFT分解为若干个小点数 DFT的组合，从而减少运算量。Wn因子具有以下两个特性，可使DFT运算量尽量分解为小点数的DFT运算：(1) 周期性：WNkhWNT =WNn N)k(2) 对称性：WNiN/2)=-W,利用这两个性质，可以使 DFT运算中有些项合并，以减少乘法次数。例子: 求当

3、N = 4时，X(2)的值3X(2)二嘉 x(n)W42n =x(O)W： x(1)W42 x(2)W44 x(3)w44n zO= x(0) x(2)W40 x(1) x(3)W：(周期性)= x(0)x(2)-x(1)x(3)W40(对称性)通过合并，使乘法次数由4次减少到1次，运算量减少。FFT的算法形式有很多种，但基本上可以分为两大类：按时间抽取(DIT )和按频率抽取(DIF )。4.1按时间抽取(DIT )的FTT为了将大点数的 DFT分解为小点数的 DFT运算，要求序列的长度 N为 复合数，最常用的是 N =2M的情况(M为正整数)。该情况下的变换称为 基2FFT。下面讨论基2情

4、况的算法。先将序列x(n)按奇偶项分解为两组：x(2r)=xi(r)、x(2r +1) =x2(r)r =0,1/ ,- -12将DFT运算也相应分为两组N-1X(k)二 DFTx(n) =、x(n)Wn =0N 4N A=Z x(n)WNkn+ E x(n)W,nn =0n n为偶数n为奇数N /2 4N/2 4 xJwN x(2r1)WN(2r 1)kr =0r =0N /2 AN /2 J八 xOWjk WN、X2(r)Wf2rkr Z0r z0N /2 JN/2八 xJWN W, v X2(r)wNk/2(因为 W：rk 二W,/2)r z0r =0k卡你)WNX2(k)其中 X,k)

5、、X2(k)分别是 x,n)、X2(n)的 N/2 点的 DFTN /2 JN /2d rkrkNXi(k) = x Xi(r)WN/2 二 x(2r)WN/2,0 乞 k1r =Qr z02N /2 dN / 2 J rk_rkNX2(k) =、X2(JWn/2 二 x(2r 1)Wn/2,0 乞 k1r =0r=02至此，一个N点DFT被分解为两个 N/2点的DFT。上面是否将全部 N点的X(k)求解出来了？N分析：Xi(k)和 X2(k)只有 N/2 个点(k=0,1,_，2-1)，则由kX(k) =Xi(k)，WNX2(k)只能求出X(k)的前N/2个点的DFT，要求出全部N点的X(k

6、)，需要找出X1(k) X2(k)和X(k N /2)的关系，其中0,1节-1。由式子X(k) = X1(k)叽的可得X(k N /2X1(k N/2) W, N/2X2(k N/2)化简得kNX(k N /2) = = XMk)-WzX2(k) , k = 0,1，12这样N点DFT可全部由下式确定出来：k =0,1,罟 T2严)X(k) =X,k) +W,X2(k)kX(k +N /2) Xk) WNX2(k)上式可用一个专用的碟形符号来表示，这个符号对应一次复乘和两次复加运算。图蝶形运算符号a wNba-wNb通过这样的分解以后，每一个2N 2 N 2N/2点的DFT只需要（一）2工24

7、次复数乘法，两个 N/2点的DFT需要2N 2 N22（）2次复乘，再加上将两个N/2点2 2DFT合并成为 N点DFT时有N / 2次与 W因子相乘，一共需要N2次复乘。可见，通过这样的分解，运算量节省了近一半。2 2因为N =2M，N/2仍然是偶数，因此可以对两个 N/2点的DFT再分别作进一步的分解，将两个 N/2点的DFT分解成两个N/4点的DFT。例如对xdr），可以在按其偶数部分及奇数部分进行分解:X(2I)=X3(I)N “丿I =0,1,，一 -1Xi (2I +1)*1)4则的运算可相应分为两组：N /4 4N /4 4Xi(k)Xi(2I)wN/k2 工 Xi(2IIWN%

8、1沐0I -0N/4 二N/44八 X3(I)WNk/4 Wj/2 X4(I)WNk/4I z0I =0kN”3(k) W,/2X4(k) k =0,1,1将系数统一为以N为周期，即W,/2 =WN2k，可得k =0,1,， -14Xi（k） =X3（k） +Wn X4（k）Xi（k N/4）=X3（k） -WfN!kX4（k）同样，对X2（k）也可进行类似的分解。一直分解下去，最后是2点的3DFT, 2点DFT的运算也可用碟形符号来表示。这样，对于一个N = 2 =8的DFT运算，其按时间抽取的分解过程及完整流图如下图所示。J * * *1-忆伍 n p317 x(x(x(xlx(x(xXX

9、 X X X X X X XJ u u n 0 2 4 6 fl fl flX X X Xu n n u 13 5 7 fl fl flX X X XDFT占E00S11 )P1I3)Ml151161PI x(x(xrx(x(x N1+N2-1 ,并且N =2M (M为整数)，即x(n)二；x(n), n = 0,1,，N11Q n = N1, N1 +1,N 1一、h(n), n= 0,1,，N21 h(n) = *0, n = N 2, N2+1，N1(2) 用FFT计算x(n)与h(n)的离散傅里叶变换x(n)= (FFT)二 X(k)(N 点)h(n)= (FFT )二 H (k)(N

10、 点)(3) 计算 Y(k)=X(k)H(k)(4) 用IFFT计算Y(k)的离散傅里叶反变换得：y(n)=IFFTY(k) (N点)4.2.2 利用FFT求相关一快速相关互相关及自相关的运算已广泛的应用于信号分析与统计分析，应用于连续时间系统也用于离散事件系统。用FFT计算相关函数称为快速相关，它与快速卷积完全类似，不同的是一个应用离散相关定理，另一个应用离散卷积定理。同样都要注意到离散傅里叶变换固有的周期性，也同样用补零的方法来绕过这个障碍。设两个离散时间信号 x(n)与y(n)为已知，离散互相关函数记作Rxy (n)，定义为QORxy (n)二 x(m)y(n m)m =O如果x(n)与

11、y(n)的序列长度分别为 N1和N2，则用FFT求相关的计算步骤 如下：(1 )对序列x(n)与y( n)补零至长为N ,使N N1+N2-1 ,并且 N =2M (M为整数)，即(、|x(n), n= 0,1,，N1-1x(n) = *0, n = N1,N1+1,，N 1/、y(n), n=0,1，N2-1y(n)=丿0, n = N2,N2+1，N-1(2) 用FFT计算x(n)与y(n)的离散傅里叶变换x(n)= (FFT)二 X(k)(N 点)y(n)= (FFT)u Y(k)(N 点)(3) 将X(k)的虚部lmX(k)改变符号，求得其共轭 X*(k)(4) 计算 Rxy(k)=x

12、*(k)Y(k)(5) 用IFFT求得相关序列Rxy(n)Rxy(n) = IFFT Rxy(k)(N 点)如果x(n) = y(n)，则求得的是自相关序列Rxx(n)4.2.3 Chirp - Z 变换采用FFT算法可以很快的计算出全部DFT值，也即Z变换在单位圆上的全部等间隔采样值。但是，很多场合下，并非整个单位圆上的频谱都是很有意义的，例如对于窄带信号过程，往往只需要对信号所在的一段频带进 行分析，这是，希望采样能密集在这段频带内，而对频带以外的部分，则可 以完全不管。另外，有时也希望采样能不局限于单位圆上，例如，语音信号 处理中，往往需要知道其Z变换的极点所在频率，如果极点位置离单位圆较远，则其单位圆上的频谱就很平滑，如图所示，这是很难从中识别出极 点所在的频率。要是采样不是沿单位圆而是沿一条接近这些极点的弧线进 行，则所得的结果将会在极点所在频率上出现明显的尖峰，从而可以较准确 的测定极点频率。螺线采样就是一种适用于这种需要的变换，并且可以采用FFT来快速计算，这种变换也称为 Chirp-Z变换，它是沿Z平面上的一段螺线作等分 角的采样，这些采样点可表达为kzk 二 AW ,k =0,1, ,M -1