《FFT在功率谱密度计算中地应用》由会员分享，可在线阅读，更多相关《FFT在功率谱密度计算中地应用（22页珍藏版）》请在金锄头文库上搜索。
1、wordFFT在功率谱密度计算中的应用一、FFT算法理论依据和编程思想FFT算法的根本思想：考察DFT与IDFT的运算发现，利用以下两个特性可减少运算量：)系数是一个周期函数，它的周期性和对称性可利用来改良运算，提高计算效率。如:因此利用这些周期性和对称性，DFT运算中有些项可合并；利用WNnk的周期性和对称性，把长度为N点的大点数的DFT运算分解为假如干个小点数的DFT。因为DFT的计算量正比于N2，N小计算量也就小。 FFT算法正是基于这样的根本思想开展起来的。它有多种形式，下面是按时间抽取的FFTN点DFT运算的分解先从一个特殊情况开始，假定N是2的整数次方，N=2M，M：正整数 1.将
2、N点的DFT分解为两个N/2点的DFT：首先将序列xn分解为两组，一组为偶数项，一组为奇数项r=0,1，N/2-1将DFT运算也相应分为两组：其中X1k和X2k分别是x1r和x2r的N/2点DFT。可见，一个N点的DFT可以分解为两个N/2点的DFT，这两个N/2点的DFT再按照上面1式合成为一个N点DFT，注意到，X1k，X2k有N/2个点，即k=0，1，N/2-1，由1式得到Xk只有N/2点，而实际上Xk有N个点，即k=0，1，N-1，要用X1k，X2k表示全部XK值，还必须应用系数w的周期性和对称性。2.X(k)的(N/2)N-1点表示：由X(k)= X1k+wkN X2k, k=0,1
3、,2,N/2-1 得：，2a因为,且同样。考虑到WNk对称性：。故2b2a式表示了Xk前半局部k=0N/2-1时的组成方式，2b式如此表示了后半局部k=N/2N-1时的组成方式。这两式所表示的运算过程可用一个称作蝶形的信号流图来表示。3.蝶形信号流图：如图1a所示，图中左面两支为输入，中间以一个小圆圈表示加、减运算，右上支为相加输出，右下支为相减输出，如果在某一支路上信号需要进展乘法运算，如此在该支路上标以箭头，并将相乘的系数标在简头边，这样2a，2b所表示的运算，可用图1b所表示的“蝶形结来表示。采用这种表示法，可将以上以讨论的分解过程用计算流图来表示。图2.6所示为N=23=8的例子。通过
4、这样分解以后，每一个N/2点DFT只需要 图2.6 N点DFT分解为2个N/2点DFTN=8(N/2)2= N2/4次复数乘法运算，两个N/2点的DFT需要2(N/2)2= N2/2 次复乘，再加上将两个N/2点DFT合成为N点DFT时，在蝶形结前的N/2次复乘，共需要(N/2)2+N/2 N2/2次复乘，由此可见，经过这样的分解处理，运算量差不多节省了一倍。4.将N/2点的DFT分解为两个N/4点的DFT：既然这样分解是有效的，由于N=2M，N/2仍然是偶数，因此可对两个N/2点的DFT再分别作进一步分解，例如对x1r和x2r可以再按其偶数局部与奇数局部分解为两个N/4点的DFT，既然这样分
5、解是有效的，由于N=2M，N/2仍然是偶数，因此可对两个N/2点的DFT再分别作进一步分解，例如对x1r和x2r可以再按其偶数局部与奇数局部分解为两个N/4点的DFT， l=0,1,，N/4-1而同样X2(k)也可这样分解，并且将系数统一为 ，这样一个8点DFT就可分解为四个2点的DFT，如图2.7所示。 图2.7 N点DFT分解为4个N/4点的DFT(N=8)5.2个点DFT的表示：最后剩下的是2点DFT，它可以用一个蝶形结表示，例如，x0，x4所组成的2点DFT就可表示式：这样，一个8点的完整的按时间抽取运算的流图如图2.8所示。 由于这样的方法每一步分解都是按输入时间序列是属于偶数还是奇
6、数来抽取的，所以称为“按时间抽取法或“时间抽取法。6.时间抽取法FFT运算特点：1蝶形运算对任何一2的整数幂N=2M，总是可以通过M次分解最后完全成为2点的DFT运算。这样的M次分解，就构成从x(n)到X(k)的M级运算过程。从上面的流图可看到，每一级运算都由N/2个蝶形运算构成。因此每一级运算都需要 (N/2次复乘和N次复加(每个结作加、减各一次)，这样，经过时间抽取后M级运算总共需要的运算：复乘复加 NM=Nlog2N实际运算量与这个数字稍有出入，因为W这几个系数实际上都不用乘法运算，因此在上面N=8的例子中，实际上只有两个系数WN1与WN3是需要乘法运算的。用时间抽取法所需的计算量，不论
7、是复乘还是复加都与Nlog2N成正比，而直接运算时如此与N2成正比。例N=2048，N2=4194304，(N/2)log2N=11264，N2/(N/2)log2N=392.4倍。FFT显然要比直接法快得多。(2原位计算当数据输入到存储器中以后，每一级运算的结果仍然储存在同一组存储器中，直到最后输出，中间无需其它存储器，这叫原位计算。例如，N=8的FFT运算，输入x(0)，x(4)，x(2)，x(6)，x(7)可分别存入A(1)，A(2)，A(8)这9个存储单元中，在第一级运算中，首先是存储单元A(1)，A(2)中x(0)，x(4)进入蝶形运算，x(0)，x(4)输入运算器后，其数值不再需要
8、保存，因此蝶形运算的结果可仍然送回存储单元A(1)，A(2)中保存，然后A(3)，A(4)中x(2)，x(6)再进入蝶形运算，其结果再送回A(3)，A(4)，一直到算完A(7)，A(8)，如此完成了第一级运算过程。第二级运算仍可采用这种原位的方式，但是进入蝶形结的组合关系不同，首先进入蝶形结的是A(1)、A(3)存储单元中的数据，运算结果仍可送回A(1)、A(3)保存，然后进入蝶形结的是A(2)、A(4)，依此类推，每一级运算均可在原位进展，这种原位运算结构可节省存储单元，降低设备本钱，还可节省找地址的时间。3序数重排对按时间抽取FFT的原位运算结构，当运算完毕时，这种结构存储单元A(1)、A
9、(2)，A(8)中正好顺序存放着X(0)，X(1)，X(2)，X(7)，因此可直接按顺序输出，但这种原位运算的输入x(n)却不能按这种自然顺序存入存储单元中，而是按X(0)，X(4)，X(2)，X(6)，X(7)的顺序存入存储单元，这种顺序看起来相当杂乱，然而它也是有规律的。当用二进制表示这个顺序时，它正好是“码位倒置的顺序。例如，原来的自然顺序应是x(1)的地方，现在放着x(4)，用二进制码表示这一规律时，如此是在x(0 0 1)处放着x(1 0 0)，x(0 1 1)处放着x(1 1 0)。即将自然顺序的二进制码位倒置过来，第一位码变成最末位码，这样倒置以后的顺序正是输入所需要的顺序。下表
10、列出N=8时按码位倒置规律所得的顺序，其结果与按时间抽取FFT流图中的输入顺序是一致的。表一 码位倒置顺序自然顺序二进码表示码位倒置码位倒置顺序0000000010011004201001023011110641000011510110156110011371111117在实际运算中，一般直接将输入数据xn按码位倒置的顺序排好输入很不方便，总是先按自然顺序输入存储单元，然后再通过变址运算将自然顺序的存储转换成码位倒置顺序的存储，然后进展FFT的原位计算。目前有许多通用DSP芯片支持这种码位倒置的寻址功能。4蝶形类型随迭代次数成倍增加观察8点FFT的三次迭代运算第一级迭代，只有一种类型的蝶形运算
11、系数W08第二级迭代，有二种类型的蝶形运算系数W08、W28，参加运算的两个数据点间隔为2。第三级迭代，有四类蝶形运算系数W08、W18、W28、W38，参加运算的两个数据点间隔为4。 所以，每次迭代的蝶形类型比上一次蝶代增加一倍，数据点间隔也增大一倍。根据相关定理与维纳-辛钦关系式可得随机信号序列x(n)的功率谱密度 (1)其估计值 (2) 如果观察到序列x(n)的N个值，即x(0), x(1),x(N-1)，就可以通过FFT直接求得X(k),再按式(2)求得Sx(k),其计算过程如图2.9所示。x(n)FFT平方平方加法除法Sx(k)X(k)2/NX(k)2X2R(k)XI2(k)XR(k
12、)图2.9 周期图法计算功率谱密度流程XI(k)二、程序设计本设计在主程序中分别调用了输入get_in()、倒序re_order()和蝶式运算butterfly(),功率谱密度计算power(),绘图on_draw五个子函数。（1） 运算主程序框图(2)、整序程序流程图LHN/2jLHN1N-2kLHij?Tx(i)x(i)x(j)x(j)TTy(i)y(i)y(j)y(j)TNNY返回i=1, N1jj-kk=k/2jj+kJk?图7 整序程序流程图3蝶形运算框图i1j0r0m12i-1m22m1m3N/ m2m4j*m2k1m4+rk2k1+m1x(k2)x(k1)-uy(k2)y(k1)
13、-vx(k1)x(k1)+uy(k1)y(k1)+vrr+1rm1?jj+1jm3?ii+1jL?返回YYYNN图8 蝶形运算框图三、程序源代码主程序/ FFT.cpp : Defines the class behaviors for the application.#include stdafx.h#include FFT.h#include FFTDlg.h/自己添加的.#includemath.h#ifdef _DEBUG#define new DEBUG_NEW#undef THIS_FILEstatic char THIS_FILE = _FILE_;#endif/ CFftDrawingAppBEGIN_MESSAGE_MAP(CFftDrawingApp, CWinApp)/AFX_MSG_MAP(CFftDrawingApp)/ NOTE - the ClassWizard will add and remove mapping macros here./ DO NO
