两点DFT简化
假设输入为x[0],x[1];输出为X[0],X[1]. 伪代码如下:
// ------------------------------------------------------------------
#define N 2
#define PI 3.1415926
//tw为旋转因子
// ------------------------------------------------------------------
int i, j
for(i=0, X[i]=0.0; i for(j=0; j X[i] += x[j] * ( cos(2*PI*i*j/N) - sin(2*PI*i*j/N) ); 注意到 j=0j=1 i=0cos 1 sin 0 tw 1 cos 1 sin 0 tw 1 i=1cos 1 Sin 0 tw 1 cos -1 sin 0 tw -1 X[0]= x[0]*(1-0) + x[1]*(1-0) = x[0] + 1*x[1]; X[1]= x[0]*(1-0) + x[1]*(-1-0) = x[0] - 1*x[1]; 这就是单个2点蝶形算法. FFT实现流程图分析(N=8,以8点信号为例) FFT implementation of an 8-point DFT as two 4-point DFTs and four 2-point DFTs 8点FFT流程图(Layer表示层, gr表示当前层的颗粒)下面以LayerI为例. LayerI部分,具有4个颗粒,每个颗粒2个输入 (注意2个输入的来源,由时域信号提供) 将输入x[k]分为两部分x_r[k], x_i[k].具有实部和虚部,时域信号本没有虚部的,因此可以让x_i[k]为0.因为LayerII, LayerIII的输入是复数,为了编码统一而强行分的.当然编码时可以判断当前层是否为1来决定是否分. 旋转因子tw = cos(2*PI*k/N)-j*sin(2*PI*k/N);也可以分为实部和虚部,令其为tw_r, tw_i; 则tw = tw_r - j*tw_i; X[k] = (x_r[k] + j*x_i[k]) + (tw_r–j*tw_i) * (x_r[k+N/2]+j*x_i[k+N/2]) 则 X_R[k] = x_r[k] + tw_r*x_r[k+N/2] + tw_i*x_i[k+N/2]; X_I[k] = x_i[k] - tw_i*x_r[k+N/2] + tw_r*x_i[k+N/2]; (注意4个输入的来源,由LayerI提供) (注意8个输入的来源,由LayerII提供) LayerI, LayerII, LayerIII从左往右,蝶形信号运算流非常明显! 假令输入为x[k], x[k+N/2],输出为X[k], X[k+N/2]. x[k]分解为x_r[k], x_i[k]部分则该蝶形运算为 X[k] = (x_r[k]-j*x_i[k]) + (x_r[k+N/2]-j*x_i[k+N/2])*(cos(2*PI*k/N)-j*sin(2*PI*k/N));再令cos(2*PI*k/N)为tw1, sin(2*PI*k/N)为tw2则 X[k] = (x_r[k]-j*x_i[k]) + (x_r[k+N/2]-j*x_i[k+N/2])*(tw1-j*tw2); X_R[k] = x_r[k] + x_r[k+N/2]*tw1 - x_i[k+N/2]*tw2; X_I[K] = x_i[k] x_r[k] = x_r[k] + x_r[k+b]*tw1 + x_i[k+b]*tw2; x_i[k] = x_i[k] - x_r[k+b]*tw2 + x_i[k+b]*tw1; 譬如8点输入x[8] 1. 先分割成2部分: x[0], x[2], x[4], x[6]和x[1], x[3], x[5], x[7] 2. 信号x[0], x[2], x[4], x[6]再分割成x[0], x[4]和x[2], x[6] 信号x[1], x[3], x[5], x[7]再分割成x[1], x[5]和x[3], x[7] 如上图: 在LayerI的时候,我们是对2点进行DFT.(一共4次DFT ) 输入为x[0]&x[4]; x[2]&x[6]; x[1]&x[5]; x[3]&x[7] 输出为y[0],y[1]; Y[2],y[3]; Y[4],y[5]; Y[6],y[7]; 流程: I.希望将输入直接转换为x[0], x[4], x[2], x[6], x[1], x[5], x[3], x[7]的顺序II.对转换顺序后的信号进行4次DFT 步骤I代码实现 /** *反转算法.这个算法效率比较低!先用起来在说,之后需要进行优化. */ static void bitrev(void) { int p=1, q, i; int bit_rev[ N ]; float xx_r[ N ]; bit_rev[ 0 ] = 0; while( p < N ) { for(q=0; q { bit_rev[ q ] = bit_rev[ q ] * 2; bit_rev[ q + p ] = bit_rev[ q ] + 1; } p *= 2; for(i=0; i for(i=0; i } // ------------------------此刻序列x重排完毕------------------------ 步骤II代码实现 int j; float TR; //临时变量 float tw1; //旋转因子 /*两点DFT */ for(k=0; k { //两点DFT简化告诉我们tw1=1 TR = x_r[k]; // TR就是A, x_r[k+b]就是B. x_r[k] = TR + tw1*x_r[k+b]; x_r[k+b] = TR - tw1*x_r[k+b]; } 在LayerII的时候,我们希望得到z,就需要对y进行DFT. y[0],y[2]; y[1],y[3]; y[4],y[6]; y[5],y[7]; z[0], z[1]; z[2],z[3]; z[4],z[5]; z[6],z[7]; 在LayerIII的时候,我们希望得到v,就需要对z进行DFT. z[0],z[4]; z[1],z[5]; z[2],z[6]; z[3],z[7]; v[0],v[1]; v[2],v[3]; v[4],v[5]; v[6],v[7]; 准备 令输入为x[s], x[s+N/2],输出为y[s], y[s+N/2] 这个N绝对不是上面的8,这个N是当前颗粒的输入样本总量 对于LayerI而言N是2;对于LayerII而言N是4;对于LayerIII而言N是8 复数乘法:(a+j*b) * (c+j*d) 实部= a*c – bd; 虚部= ad + bc; 旋转因子: 实现(C描述) #include #include #include //#include "complex.h" // -------------------------------------------------------------------------- #define N 8 //64 #define M 3 //6 //2^m=N #define PI 3.1415926 // -------------------------------------------------------------------------- float twiddle[N/2] = {1.0, 0.707, 0.0, -0.707}; float x_r[N] = {1, 1, 1, 1, 0, 0, 0, 0}; float x_i[N]; //N=8 /* float twiddle[N/2] = {1, 0.9951, 0.9808, 0.9570, 0.9239, 0.8820, 0.8317, 0.7733, 0.7075, 0.6349, 0.5561, 0.4721, 0.3835, 0.2912, 0.1961, 0.0991, 0.0000,-0.0991,-0.1961,-0.2912,-0.3835,-0.4721,-0.5561,-0.6349, -0.7075,-0.7733, 0.8317,-0.8820,-0.9239,-0.9570,-0.9808,-0.9951}; //N=64 float x_r[N]={1,1,1,1,1,1,1,1, 1,1,1,1,1,1,1,1, 1,1,1,1,1,1,1,1, 1,1,1,1,1,1,1,1, 0,0,0,0,0,0,0,0, 0,0,0,0,0,0,0,0, 0,0,0,0,0,0,0,0, 0,0,0,0,0,0,0,0,}; float x_i[N]; */ FILE *fp; // ----------------------------------- func -----------------------------------/** *初始化输出虚部 */ static void fft_init(void) { int i; for(i=0; i } /** *反转算法.将时域信号重新排序. *这个算法有改进的空间 */ static void bitrev(void) { int p=1, q, i; int bit_rev[ N ]; // float xx_r[ N ]; // bit_rev[ 0 ] = 0; while( p < N ) { for(q=0; q { bit_rev[ q ] = bit_rev[ q ] * 2; bit_rev[ q + p ] = bit_rev[ q ] + 1; } p *= 2; } for(i=0; i for(i=0; i } /* ------------ add by sshc625 ------------ */ static void bitrev2(void) { return; } /* */ void display(void) { printf("\n\n"); int i; for(i=0; i printf("%f\t%f\n", x_r[i], x_i[i]); } /** * */ void fft1(void) { fp = fopen("log1.txt","a+"); int L, i, b, j, p, k, tx1, tx2; float TR, TI, temp; //临时变量 float tw1, tw2; /*深M.对层进行循环. L为当前层,总层数为M. */ for(L=1; L<=M; L++) { fprintf(fp,"----------Layer=%d----------\n", L); /* b的意义非常重大,b表示当前层的颗粒具有的输入样本点数*/ b = 1; i = L - 1; while(i > 0) { b *= 2; i--; } // --------------是否外层对颗粒循环,内层对样本点循环逻辑性更强一些呢! --------------/* * outter对参与DFT的样本点进行循环 * L=1,循环了1次(4个颗粒,每个颗粒2个样本点) * L=2,循环了2次(2个颗粒,每个颗粒4个样本点) * L=3,循环了4次(1个颗粒,每个颗粒8个样本点) */ for(j=0; j { /*求旋转因子tw1 */ p = 1; i = M - L; // M是为总层数, L为当前层. while(i > 0) { p = p*2; i--; } p = p * j; tx1 = p % N; tx2 = tx1 + 3*N/4; tx2 = tx2 % N; // tw1是cos部分,实部; tw2是sin部分,虚数部分. tw1 = ( tx1>=N/2)? -twiddle[tx1-N/2] : twiddle[ tx1 ]; tw2 = ( tx2>=N/2)? -twiddle[tx2-(N/2)] : twiddle[tx2]; /* * inner对颗粒进行循环 * L=1,循环了4次(4个颗粒,每个颗粒2个输入) * L=2,循环了2次(2个颗粒,每个颗粒4个输入) * L=3,循环了1次(1个颗粒,每个颗粒8个输入) */ for(k=j; k { TR = x_r[k]; // TR就是A, x_r[k+b]就是B. TI = x_i[k]; temp = x_r[k+b]; /* *如果复习一下(a+j*b)(c+j*d)两个复数相乘后的实部虚部分别是什么 *就能理解为什么会如下运算了,只有在L=1时候输入才是实数,之后层的 *输入都是复数,为了让所有的层的输入都是复数,我们只好让L=1时候的 *输入虚部为0 * x_i[k+b]*tw2是两个虚数相乘 */ fprintf(fp,"tw1=%f, tw2=%f\n", tw1, tw2); x_r[k] = TR + x_r[k+b]*tw1 + x_i[k+b]*tw2; x_i[k] = TI - x_r[k+b]*tw2 + x_i[k+b]*tw1; x_r[k+b] = TR - x_r[k+b]*tw1 - x_i[k+b]*tw2; x_i[k+b] = TI + temp*tw2 - x_i[k+b]*tw1; fprintf(fp,"k=%d, x_r[k]=%f, x_i[k]=%f\n", k, x_r[k], x_i[k]); fprintf(fp,"k=%d, x_r[k]=%f, x_i[k]=%f\n", k+b, x_r[k+b], x_i[k+b]); }// }// }// } /** * ------------ add by sshc625 ------------ *该实现的流程为 * for( Layer ) * for( Granule ) * for( Sample ) * * * * */ void fft2(void) { fp = fopen("log2.txt","a+"); int cur_layer, gr_num, i, k, p; float tmp_real, tmp_imag, temp; //临时变量,记录实部 float tw1, tw2;//旋转因子,tw1为旋转因子的实部cos部分, tw2为旋转因子的虚部sin部分. int step; //步进 int sample_num; //颗粒的样本总数(各层不同,因为各层颗粒的输入不同) /*对层循环*/ for(cur_layer=1; cur_layer<=M; cur_layer++) { /*求当前层拥有多少个颗粒(gr_num) */ gr_num = 1; i = M - cur_layer; while(i > 0) { i--; gr_num *= 2; } /*每个颗粒的输入样本数N' */ sample_num = (int)pow(2, cur_layer); /*步进.步进是N'/2 */ step = sample_num/2; /* */ k = 0; /*对颗粒进行循环*/ for(i=0; i { /* *对样本点进行循环,注意上限和步进 */ for(p=0; p { //旋转因子,需要优化... tw1 = cos(2*PI*p/pow(2, cur_layer)); tw2 = -sin(2*PI*p/pow(2, cur_layer)); tmp_real = x_r[k+p]; tmp_imag = x_i[k+p]; temp = x_r[k+p+step]; /*(tw1+jtw2)(x_r[k]+jx_i[k]) * * real : tw1*x_r[k] - tw2*x_i[k] * imag : tw1*x_i[k] + tw2*x_r[k] *我想不抽象出一个 * typedef struct { * double real; //实部 * double imag; //虚部 * } complex;以及针对complex的操作 *来简化复数运算是否是因为效率上的考虑! */ /*蝶形算法*/ x_r[k+p] = tmp_real + ( tw1*x_r[k+p+step] - tw2*x_i[k+p+step] ); x_i[k+p] = tmp_imag + ( tw2*x_r[k+p+step] + tw1*x_i[k+p+step] ); /* X[k] = A(k)+WB(k) * X[k+N/2] = A(k)-WB(k)的性质可以优化这里*/ //旋转因子,需要优化... tw1 = cos(2*PI*(p+step)/pow(2, cur_layer)); tw2 = -sin(2*PI*(p+step)/pow(2, cur_layer)); x_r[k+p+step] = tmp_real + ( tw1*temp - tw2*x_i[k+p+step] ); x_i[k+p+step] = tmp_imag + ( tw2*temp + tw1*x_i[k+p+step] ); printf("k=%d, x_r[k]=%f, x_i[k]=%f\n", k+p, x_r[k+p], x_i[k+p]); printf("k=%d, x_r[k]=%f, x_i[k]=%f\n", k+p+step, x_r[k+p+step], x_i[k+p+step]); } /*开跳!:) */ k += 2*step; } } } /* *后记: *究竟是颗粒在外层循环还是样本输入在外层,好象也差不多,复杂度完全一样. *但以我资质愚钝花费了不少时间才弄明白这数十行代码. *从中我发现一个于我非常有帮助的教训,很久以前我写过一部分算法,其中绝大多数都是递归. *将数据量减少,减少再减少,用归纳的方式来找出数据量加大代码的规律 *比如FFT * 1.先写死LayerI的代码;然后再把LayerI的输出作为LayerII的输入,又写死代码; ...... * 大约3层就可以统计出规律来.这和递归也是一样,先写死一两层,自然就出来了! * 2.有的功能可以写伪代码,不急于求出结果,降低复杂性,把逻辑结果定出来后再添加. * 比如旋转因子就可以写死,就写1.0.流程出来后再写旋转因子. */ void dft(void) { int i, n, k, tx1, tx2; float tw1,tw2; float xx_r[N],xx_i[N]; /* * clear any data in Real and Imaginary result arrays prior to DFT */ for(k=0; k<=N-1; k++) xx_r[k] = xx_i[k] = x_i[k] = 0.0; // caculate the DFT for(k=0; k<=(N-1); k++) { for(n=0; n<=(N-1); n++) { tx1 = (n*k); tx2 = tx1+(3*N)/4; tx1 = tx1%(N); tx2 = tx2%(N); if(tx1 >= (N/2)) tw1 = -twiddle[tx1-(N/2)]; else tw1 = twiddle[tx1]; if(tx2 >= (N/2)) tw2 = -twiddle[tx2-(N/2)]; else tw2 = twiddle[tx2]; xx_r[k] = xx_r[k]+x_r[n]*tw1; xx_i[k] = xx_i[k]+x_r[n]*tw2; } xx_i[k] = -xx_i[k]; } // display for(i=0; i printf("%f\t%f\n", xx_r[i], xx_i[i]); } // ---------------------------------------------------------------------------int main(void) { fft_init( ); bitrev( ); // bitrev2( ); //fft1( ); fft2( ); display( ); system("pause"); // dft(); return1; } 创作编号： GB8878185555334563BT9125XW 创作者： 凤呜大王* 《测试信号分析及处理》课程作业 快速傅里叶变换 一、程序设计思路 快速傅里叶变换的目的是减少运算量，其用到的方法是分级进行运算。全部计算分解为M 级，其中N M 2log =；在输入序列()i x 中是按码位倒序排列的，输出序列()k X 是按顺序排列；每级包含 2N 个蝶形单元，第i 级有i N 2 个群，每个群有12-i 个蝶形单元； 每个蝶形单元都包含乘r N W 和r N W -系数 的运算，每个蝶形单元数据的间隔为12-i ，i 为第i 级； 同一级中各个群的系数W 分布规律完全相同。 将输入序列()i x 按码位倒序排列时，用到的是倒序算法——雷德算法。 自然序排列的二进制数，其下面一个数总比上面的数大1，而倒序二进制数的下面一个数是上面一个数在最高位加1并由高位向低位仅为而得到的。 若已知某数的倒序数是J ，求下一个倒序数，应先判断J 的最高位是 否为0，与2 N k =进行比较即可得到结果。如果J k >，说明最高位为0， 应把其变成1，即2 N J +，这样就得到倒序数了。如果J k ≤，即J 的最高 位为1，将最高位化为0，即2N J -，再判断次高位；与4N k =进行比较， 若为0，将其变位1，即4 N J +，即得到倒序数，如果次高位为1，将其化 为0，再判断下一位……即从高位到低位依次判断其是否为1，为1将其变位0，若这一位为0，将其变位1，即可得到倒序数。若倒序数小于顺序数，进行换位，否则不变，防治重复交换，变回原数。 注：因为0的倒序数为0，所以可从1开始进行求解。 二、程序设计框图 （1）倒序算法——雷德算法流程图 (2)FFT算法流程 FFT相关原理及使用注意事项 在信号分析与处理中，频谱分析是重要的工具。FFT（Fast Fourier Transform，快速傅立叶变换）可以将时域信号转换至频域，以获得信号的频率结构、幅度、相位等信息。该算法在理工科课程中都有介绍，众多的仪器或软件亦集成此功能。FFT实用且高效，相关原理与使用注意事项也值得好好学习。 一、何为FFT 对于模拟信号的频谱分析，首先得使用ADC（模拟数字转换器）进行采样，转换为有限序列，其非零值长度为N，经DFT（离散傅立叶变换）即可转化为频域。DFT变换式为： 在上式中，N点序列的DFT需要进行N2次复数乘法和次复数加法，运算量大。 FFT是DFT的快速算法，利用DFT运算中的对称性与周期性，将长序列DFT分解为短序列DFT 之和。最终运算量明显减少，使得FFT应用更加广泛。 FFT基于一个基本理论：任何连续的波形，都可以分解为不同频率的正弦波形的叠加。FFT将采样得到的原始信号，转化此信号所包含的正弦波信号的频率、幅度、相位，为信号分析提供一个创新视觉。 例如在日常生活中有使用到的AM（Amplitude Modulation，幅度调制）广播，其原理是将人的声音（频率约20Hz至20kHz，称为调制波）调制到500kHz～1500kHz正弦波上（称为载波）中，载波的幅度随调制波的幅度变化。声音经这样调制后，可以传播得更远。在AM的时域波形（波形电压随时间的变化曲线），载波与调制波特征不易体现，而在FFT后的幅频曲线中则一目了然。如下图为1000kHz载波、10kHz调制波的AM调制信号，时域信号经FFT后其频率能量出现在990kHz、1.01MHz频率处，符合理论计算。 图 1 调制波10kHz、载波1000kHz的AM时域与频域曲线 二、FFT相关知识 现实生活中的模拟信号，大多都是连续复杂的，其频谱分量十分丰富。正如在数学中常量π，其真实值是个无理数。当用3.14来替代π时，计算值与真实值就会有偏差。在使用FFT这个工具时，受限于采样时的频率Fs、采样点长度N、ADC的分辨率n bit等因素的制约，所得到的信息会有所缺失与混淆。 1.奈奎斯特区与波形混叠 FFT分析结果中，存在一个那奈奎斯特区的概念，其宽度为采样率的一半Fs/2，信号频 实验七 快速傅立叶变换（FFT ）实验 一 实验目的 1． 熟悉CCS 集成开发环境； 2． 了解FFT 的算法原理和基本性质； 3． 熟悉DSP 中cmd 文件的作用及对它的修改； 4． 学习用FFT 对连续信号和时域信号进行频谱分析的方法； 5． 利用DSPLIB 中现有的库函数； 6． 了解DSP 处理FFT 算法的特殊寻址方式； 7． 熟悉对FFT 的调试方法。 二 实验内容 本实验要求使用FFT 变换对一个时域信号进行频谱分析，同时进行IFFT 。这里用到时域信号可以是来源于信号发生器输入到CODEC 输入端，也可以是通过其他工具计算获取的数据表。本实验使用Matlab 语言实现对FFT 算法的仿真，然后将结果和DSP 分析的结果进行比较，其中原始数据也直接来自Matlab 。 三 实验原理 一个N 点序列][k x 的DFT ][m X ,以及IDFT 分别定义为： 1,,1,0,][][1 0-==∑-=N m W k x m X km N N k 1,,1,0,][1 ][1 -== --=∑ N k W m X N k x km N N m 如果利用上式直接计算DFT,对于每一个固定的m,需要计算N 次复数乘法，N-1次加法,对于N 个不同的m,共需计算N 的2次方复数乘法,N*(N-1)次复数加法.显然,随着N 的增加,运算量将急剧增加, 快速傅里叶算法有效提高计算速度（本例使用基2 FFT 快速算法）,利用FFT 算法只需(N/2)logN 次运算。 四 知识要点 . 1、 CMD 文件的功能及编写 2、 一种特殊的寻址方式：间接寻址 间接寻址是按照存放在某个辅助寄存器的16位地址寻址的。C54x 的8个辅助寄存器 （AR0—AR7）都可以用来寻址64K 字数据存储空间中的任何一个存储单元。 3、 TMS320C54x DSPLIB 中关于FFT 变换的一些函数的调用（SPRA480B.pdf ） 利用DSPLIB 库时，在主程序中要包含头文件：54xdsp.lib 4、 FFT 在CCS 集成开发环境下的相关头文件 #include FFT算法 FFT是离散傅立叶变换的快速算法，可以将一个信号变换到频域。有些信号在时域上是很难看出什么特征的，但是如果变换到频域之后，就很容易看出特征了。这就是很多信号分析采用FFT变换的原因。另外，FFT可以将一个信号的频谱提取出来，这在频谱分析方面也是经常用的。虽然很多人都知道FFT是什么，可以用来做什么，怎么去做，但是却不知道FFT之后的结果是什意思、如何决定要使用多少点来做FFT。现在圈圈就根据实际经验来说说FFT结果的具体物理意义。一个模拟信号，经过ADC采样之后，就变成了数字信号。采样定理告诉我们，采样频率要大于信号频率的两倍，这些我就不在此罗嗦了。 采样得到的数字信号，就可以做FFT变换了。N个采样点，经过FFT之后，就可以得到N个点的FFT结果。为了方便进行FFT运算，通常N取2的整数次方。 假设采样频率为Fs，信号频率F，采样点数为N。那么FFT之后结果就是一个为N点的复数。每一个点就对应着一个频率点。这个点的模值，就是该频率值下的幅度特性。具体跟原始信号的幅度有什么关系呢假设原始信号的峰值为A，那么FFT的结果的每个点（除了第一个点直流分量之外）的模值就是A的N/2倍。而第一个点就是直流分量，它的模值就是直流分量的N倍。而每个点的相位呢，就是在该频率下的信号的相位。第一个点表示直流分量（即0Hz），而最后一个点N的再下一个点（实际上这个点是不存在的，这里是假设的第N+1个点，也可以看做是将第一个点分做两半分，另一半移到最后）则表示采样频率Fs，这中间被N-1个点平均分成N等份，每个点的频率依次增加。例如某点n所表示的频率为：Fn=(n-1)*Fs/N。由上面的公式可以看出，Fn所能分辨到频率为为Fs/N，如果采样频率Fs为1024Hz，采样点数为1024点，则可以分辨到1Hz。1024Hz的采样率采样1024点，刚好是1秒，也就是说，采样1秒时间的信号并做FFT，则结果可以分析到1Hz，如果采样2秒时间的信号并做FFT，则结果可以分析到。如果要提高频率分辨力，则必须增加采样点数，也即采样时间。频率分辨率和采样时间是倒数关系。 假设FFT之后某点n用复数a+bi表示，那么这个复数的模就是An=根号a*a+b*b，相位就是Pn=atan2(b,a)。根据以上的结果，就可以计算出n点（n≠1，且n<=N/2）对应的信号的表达式为：An/(N/2)*cos(2*pi*Fn*t+Pn)，即2*An/N*cos(2*pi*Fn*t+Pn)。对于n=1点的信号，是直流分量，幅度即为A1/N。由于FFT结果的对称性，通常我们只使用前半部分的结果，即小于采样频率一半的结果。好了，说了半天，看着公式也晕，下面圈圈以一个实际的信号来做说明。 假设我们有一个信号，它含有2V的直流分量，频率为50Hz、相位为-30度、幅度为3V的交流信号，以及一个频率为75Hz、相位为90度、幅度为的交流信号。用数学表达式就是如下： S=2+3*cos(2*pi*50*t-pi*30/180)+*cos(2*pi*75*t+pi*90/180) 式中cos参数为弧度，所以-30度和90度要分别换算成弧度。 我们以256Hz的采样率对这个信号进行采样，总共采样256点。按照我们上面的分析，Fn=(n-1)*Fs/N，我们可以知道，每 《测试信号分析及处理》课程作业 快速傅里叶变换 一、程序设计思路 快速傅里叶变换的目的是减少运算量，其用到的方法是分级进行运算。全部计算分解为M 级，其中N M 2log =；在输入序列()i x 中是按码位倒序排列的，输出序列()k X 是按顺序排列；每级包含2N 个蝶形单元，第i 级有i N 2 个群，每个群有12-i 个蝶形单元； 每个蝶形单元都包含乘r N W 和r N W -系数的运算，每个蝶形 单元数据的间隔为12-i ，i 为第i 级； 同一级中各个群的系数W 分布规律完全相同。 将输入序列()i x 按码位倒序排列时，用到的是倒序算法——雷德算法。 自然序排列的二进制数，其下面一个数总比上面的数大1，而倒序二进制数的下面一个数是上面一个数在最高位加1并由高位向低位仅为而得到的。 若已知某数的倒序数是J ，求下一个倒序数，应先判断J 的最高位是否为0，与2 N k =进行比较即可得到结果。如果J k >，说明最高位为0，应把其变成1，即2 N J +，这样就得到倒序数了。如果J k ≤，即J 的最高位为1，将最高位化为0，即2N J -，再判断次高位；与4N k =进行比较，若为0，将其变位1，即4 N J +，即得到倒序数，如果次高位为1，将其化为0，再判断下一位……即从高位到低位依次判断其是否为1，为1将其变位0，若这一位为0，将其变位1，即可得到倒序数。若倒序数小于顺序数，进行换位，否则不变，防治重复交换，变回原数。 注：因为0的倒序数为0，所以可从1开始进行求解。 二、程序设计框图 （1）倒序算法——雷德算法流程图 (2)FFT算法流程 FFT原理与实现 在数字信号处理中常常需要用到离散傅立叶变换(DFT)，以获取信号的频域特征。尽管传统的DFT算法能够获取信号频域特征，但是算法计算量大，耗时长，不利于计算机实时对信号进行处理。因此至DFT被发现以来，在很长的一段时间内都不能被应用到实际的工程项目中，直到一种快速的离散傅立叶计算方法——FFT，被发现，离散是傅立叶变换才在实际的工程中得到广泛应用。需要强调的是，FFT并不是一种新的频域特征获取方式，而是DFT的一种快速实现算法。本文就FFT的原理以及具体实现过程进行详尽讲解。 DFT计算公式 其中x(n)表示输入的离散数字信号序列，WN为旋转因子，X(k)一组N点组成的频率成分的相对幅度。一般情况下，假设x(n)来自于低通采样，采样频率为fs，那么X(k)表示了从-fs/2率开始，频率间隔为fs/N,到fs/2-fs/N截至的N个频率点的相对幅度。因为DFT计算得到的一组离散频率幅度只实际上是在频率轴上从成周期变化的，即X(k+N)=X(k)。因此任意取N个点均可以表示DFT的计算效果，负频率成分比较抽象，难于理解，根据X(k)的周期特性，于是我们又可以认为X(k)表示了从零频率开始，频率间隔为fs/N,到fs-fs/N截至的N个频率点的相对幅度。 N点DFT的计算量 根据(1)式给出的DFT计算公式，我们可以知道每计算一个频率点X(k)均需要进行N次复数乘法和N-1次复数加法，计算N各点的X(k)共需要N^2次复数乘法和N*(N-1)次复数加法。当x(n)为实数的情况下，计算N点的DFT需要2*N^2次实数乘法，2*N*(N-1)次实数加法。 旋转因子WN的特性 1. W的对称性 N W的周期性 2. N W的可约性 3. N 快速傅里叶变换的原理及其应用 摘要 快速傅氏变换（FFT），是离散傅氏变换的快速算法，它是根据离散傅氏变换的奇、偶、虚、实等特性，对离散傅立叶变换的算法进行改进获得的。它对傅氏变换的理论并没有新的发现，但是对于在计算机系统或者说数字系统中应用离散傅立叶变换，可以说是进了一大步。傅里叶变换的理论与方法在“数理方程”、“线性系统分析”、“信号处理、仿真”等很多学科领域都有着广泛应用,由于计算机只能处理有限长度的离散的序列,所以真正在计算机上运算的是一种离散傅里叶变换. 虽然傅里叶运算在各方面计算中有着重要的作用，但是它的计算过于复杂，大量的计算对于系统的运算负担过于庞大，使得一些对于耗电量少，运算速度慢的系统对其敬而远之，然而，快速傅里叶变换的产生，使得傅里叶变换大为简化，在不牺牲耗电量的条件下提高了系统的运算速度，增强了系统的综合能力，提高了运算速度，因此快速傅里叶变换在生产和生活中都有着非常重要的作用，对于学习掌握都有着非常大的意义。 关键词快速傅氏变换；快速算法；简化；广泛应用 Abstract Fast Fourier Transform (FFT), is a discrete fast Fourier transform algorithm, which is based on the Discrete Fourier Transform of odd and even, false, false, and other characteristics of the Discrete Fourier Transform algorithms improvements obtained. Its Fourier transform theory has not found a new, but in the computer system or the application of digital systems Discrete Fourier Transform can be said to be a big step into. Fourier transform theory and methods in the "mathematical equation" and "linear systems analysis" and "signal processing, simulation," and many other areas have a wide range of applications, as the computer can only handle a limited length of the sequence of discrete, so true On the computer's operation is a discrete Fourier transform. Fourier Although all aspects of computing in the calculation has an important role, but its calculation was too complicated, a lot of computing system for calculating the burden is too large for some Less power consumption, the slow speed of operation of its system at arm's length, however, have the fast Fourier transform, Fourier transform greatly simplifying the making, not in power at the expense of the conditions to increase the speed of computing systems, and enhance the system The comprehensive ability to improve the speed of operation, the Fast Fourier Transform in the production and life have a very important role in learning to master all have great significance. Key words Fast Fourier Transform; fast algorithm; simplified; widely used 2010-10-07 21:10:09| 分类：数字信号处理 | 标签：fft dft |字号订阅 在数字信号处理中常常需要用到离散傅立叶变换(DFT)，以获取信号的频域特征。尽管传统的DFT算法能够获取信号频域特征，但是算法计算量大，耗时长，不利于计算机实时对信号进行处理。因此至DFT被发现以来，在很长的一段时间内都不能被应用到实际的工程项目中，直到一种快速的离散傅立叶计算方法——FFT，被发现，离散是傅立叶变换才在实际的工程中得到广泛应用。需要强调的是，FFT并不是一种新的频域特征获取方式，而是DFT的一种快速实现算法。本文就FFT的原理以及具体实现过程进行详尽讲解。 DFT计算公式 本文不加推导地直接给出DFT的计算公式： 其中x(n)表示输入的离散数字信号序列，WN为旋转因子，X(k)一组N点组成的频率成分的相对幅度。一般情况下，假设x(n)来自于低通采样，采样频率为fs，那么X(k)表示了从-fs/2率开始，频率间隔为fs/N,到fs/2-fs/N截至的N个频率点的相对幅度。因为DFT 计算得到的一组离散频率幅度只实际上是在频率轴上从成周期变化的，即X(k+N)=X(k)。因此任意取N个点均可以表示DFT的计算效果，负频率成分比较抽象，难于理解，根据X(k)的周期特性，于是我们又可以认为X(k)表示了从零频率开始，频率间隔为fs/N,到fs-fs/N截至的N个频率点的相对幅度。 N点DFT的计算量 根据(1)式给出的DFT计算公式，我们可以知道每计算一个频率点X(k)均需要进行N次复数乘法和N-1次复数加法，计算N各点的X(k)共需要N^2次复数乘法和N*(N-1)次复数加法。当x(n)为实数的情况下，计算N点的DFT需要2*N^2次实数乘法，2*N*(N-1)次实数加法。 旋转因子WN的特性 1.WN的对称性 2.WN的周期性 3.WN的可约性 快速傅里叶变换(FFT)的原理及公式 原理及公式 非周期性连续时间信号x(t)的傅里叶变换可以表示为 式中计算出来的是信号x(t)的连续频谱。但是，在实际的控制系统中能够得到的是连续信号x(t)的离散采样值x(nT)。因此需要利用离散信号x(nT)来计算信号x(t)的频谱。 有限长离散信号x(n)，n=0，1，…，N-1的DFT定义为： 可以看出，DFT需要计算大约N2次乘法和N2次加法。当N较大时，这个计算量是很大的。利用WN的对称性和周期性，将N点DFT分解为两个N／2点 的DFT，这样两个N／2点DFT总的计算量只是原来的一半，即(N／2)2+(N／2)2=N2／2，这样可以继续分解下去，将N／2再分解为N／4点DFT等。对于N=2m点的DFT都可以分解为2点的DFT，这样其计算量可以减少为(N／2)log2N 次乘法和Nlog2N次加法。图1为FFT与DFT-所需运算量与计算点数的关系曲线。由图可以明显看出FFT算法的优越性。 将x(n)分解为偶数与奇数的两个序列之和，即 x1(n)和x2(n)的长度都是N／2，x1(n)是偶数序列，x2(n)是奇数序列，则 其中X1(k)和X2(k)分别为x1(n)和x2(n)的N／2点DFT。由于X1(k)和X2(k)均以N／2为周期，且WN k+N/2=-WN k，所以X(k)又可表示为： 上式的运算可以用图2表示，根据其形状称之为蝶形运算。依此类推，经过m-1次分解，最后将N点DFT分解为N／2个两点DFT。图3为8点FFT的分解流程。 FFT算法的原理是通过许多小的更加容易进行的变换去实现大规模的变换，降低了运算要求，提高了与运算速度。FFT不是DFT的近似运算，它们完全是等效的。 关于FFT精度的说明： 因为这个变换采用了浮点运算，因此需要足够的精度，以使在出现舍入误差时，结果中的每个组成部分的准确整数值仍是可辨认的。为了FFT的舍入误差，应该允许增加几倍log2(log2N)位的二进制。以256为基数、长度为N字节的数 哈!经过连续几个晚上的奋战,终于弄懂了FFT推导过程及实现! Happy? 基2 FFT总的思想是将输入信号对半分割,再对半分割,再再对半分割(以下省略10000个再再...?)直至分割到2点. 两点DFT简化 假设输入为x[0],x[1];输出为X[0],X[1]. 伪代码如下: // ------------------------------------------------------------------ #define N 2 #define PI 3.1415926 // ------------------------------------------------------------------ int i, j for(i=0, X[i]=0.0; i 串行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变换就只需要Nlog(2)(N)次的运算，N在1024点时，运算量仅有10240次，是先前的直接算法的1%，点数越多，运算量的节约就越大，这就是FFT的优越性。 关键字：FFT 蝶式计算傅里叶变换 目录 一．题目及要求 (1) 1.1题目 (1) 二．设计算法、算法原理 (1) 2.1算法原理与设计 (1) 2.2设计步骤 (2) 三．算法描述、设计流程 (4) 3.1算法描述 (4) 3.2流程图 (6) 四．源程序代码及运行结果 (8) 4.1源程序代码 (8) 4.2运行结果 (13) 五．算法分析、优缺点 (15) 5.1算法分析 (15) 5.2优缺点 (16) 六．总结 (17) 七．参考文献 (18) FFT原理与实现 2010-10-07 21:10:09| 分类：数字信号处理 | 标签：fft dft |举报|字号订阅 在数字信号处理中常常需要用到离散傅立叶变换(DFT)，以获取信号的频域特征。尽管传统的DFT算法能够获取信号频域特征，但是算法计算量大，耗时长，不利于计算机实时对信号进行处理。因此至DFT被发现以来，在很长的一段时间内都不能被应用到实际的工程项目中，直到一种快速的离散傅立叶计算方法——FFT，被发现，离散傅立叶变换才在实际的工程中得到广泛应用。需要强调的是，FFT并不是一种新的频域特征获取方式，而是DFT的一种快速实现算法。本文就FFT的原理以及具体实现过程进行详尽讲解。 DFT计算公式 本文不加推导地直接给出DFT的计算公式： 其中x(n)表示输入的离散数字信号序列，WN为旋转因子，X(k)为输入序列x(n)对应的N个离散频率点的相对幅度。一般情况下，假设x(n)来自于低通采样，采样频率为fs，那么X(k)表示了从-fs/2率开始，频率间隔为fs/N,到fs/2-fs/N截至的N个频率点的相对幅度。因为DFT计算得到的一组离散频率幅度值实际上是在频率轴上从成周期变化的，即X(k+N)=X(k)。因此任意取连续的N个点均可以表示DFT的计算效果，负频率成分比较抽象，难于理解，根据X(k)的周期特性，于是我们又可以认为X(k)表示了从零频率开始，频率间隔为fs/N,到fs-fs/N 截至的N个频率点的相对幅度。 N点DFT的计算量 根据(1)式给出的DFT计算公式，我们可以知道每计算一个频率点X(k)均需要进行N次复数乘法和N-1次复数加法，计算N各点的X(k)共需要N^2次复数乘法和N*(N-1)次复数加法。当x(n)为实数的情况下，计算N点的DFT需要2*N^2次实数乘法，2*N*(N-1)次实数加法。 旋转因子WN的特性 1.WN的对称性 2.WN的周期性 3.WN的可约性 根据以上这些性质，我们可以得到式(5)的一系列有用结果 基-2 FFT算法推导 假设采样序列点数为N=2^L，L为整数，如果不满足这个条件可以人为地添加若干个0以使采样序列点数满足这一要求。首先我们将序列x(n)按照奇偶分为两组如下： 于是根据DFT计算公式(1)有： 基4FFT原理及MATLAB实现 一．时域抽取法基4FFT基本原理： 有限长序列错误！未找到引用源。的N点DFT为： 错误！未找到引用源。 k=0,1,2…,N-1 设序列错误！未找到引用源。长度为N，且满足错误！未找到引用源。,M 为自然数，可把错误！未找到引用源。按下列方法分解为4个长为N/4的子序列：错误！未找到引用源。 r=0,1,…,N/4-1 错误！未找到引用源。 r=0,1,…,N/4-1 错误！未找到引用源。 r=0,1,…,N/4-1 错误！未找到引用源。 r=0,1,…,N/4-1 则错误！未找到引用源。的DFT为X(k)=错误！未找到引用源。+错误！未找到引用源。 =错误！未找到引用源。 =错误！未找到引用源。 因为：错误！未找到引用源。 所以：原式可表示为： X(k)=错误！未找到引用源。 =错误！未找到引用源。 k=0,1,…,N/1 其中：错误！未找到引用源。分别为错误！未找到引用源。的N/4点DFT，即：错误！未找到引用源。=错误！未找到引用源。 错误！未找到引用源。= 错误！未找到引用源。 错误！未找到引用源。= 错误！未找到引用源。 错误！未找到引用源。= 错误！未找到引用源。 由于错误！未找到引用源。均以N/4为周期，且有： 错误！未找到引用源。，错误！未找到引用源。，错误！未找到引用源。，错误！未找到引用源。 所以： X(k)=错误！未找到引用源。 k=0,1,…,N/4-1 X(k+ N/4)=错误！未找到引用源。 k=0,1,…,N/4-1 X(k+ 2N/4)=错误！未找到引用源。 k=0,1,…,N/4-1 X(k+ 3N/4)=错误！未找到引用源。 k=0,1,…,N/4-1 二．运算规律及编程思想： 1.按照上述分解法，再对错误！未找到引用源。进行反复分解，直到每个序列的长度等于4为止，这个4点的DFT的数学表达式为： 错误！未找到引用源。错误！未找到引用源。错误！未找到引用源。错误！未找到引用源。 2.旋转因子错误！未找到引用源。与运算级数的关系（参考《数字信号处理（第三版）》西安电子科技大学出版社第115页）如下： 其中，L表示运算级数（L=1，2，…,M）(M=错误！未找到引用源。)(J=0,1,2,…,错误！未找到引用源。) 错误！未找到引用源。 ,(错误！未找到引用源。0,1,2…,错误！未找到引用源。) 3.序列的倒序： 与基2FFT的倒序相似（参考《数字信号处理（第三版）》西安电子科技大学出版社第116页） 由于错误！未找到引用源。，因此顺序数可用M位4进制数（错误！未找到引用源。）表示。M次时域抽取如下： 第一次按最低位错误！未找到引用源。的0,1,2,3将x(n)分解为4组，第二次又按次低位错误！未找到引用源。的0,1,2,3值分别对上面所得的4组分组；以此类推，第M次按错误！未找到引用源。位分解，最后得到4进制倒序数。最终可以得到这样的规律：只要将顺序数（错误！未找到引用源。）的4进制数倒置，就能得到对应的4进制倒序数（错误！未找到引用源。）。 4.运算流程图： FFT原理 FFT是离散傅立叶变换的快速算法，可以将一个信号变换 到频域。有些信号在时域上是很难看出什么特征的，但是如 果变换到频域之后，就很容易看出特征了。这就是很多信号 分析采用FFT变换的原因。另外，FFT可以将一个信号的频谱 提取出来，这在频谱分析方面也是经常用的。 虽然很多人都知道FFT是什么，可以用来做什么，怎么去 做，但是却不知道FFT之后的结果是什意思、如何决定要使用 多少点来做FFT。 现在圈圈就根据实际经验来说说FFT结果的具体物理意义。 一个模拟信号，经过ADC采样之后，就变成了数字信号。采样 定理告诉我们，采样频率要大于信号频率的两倍，这些我就 不在此罗嗦了。 采样得到的数字信号，就可以做FFT变换了。N个采样点， 经过FFT之后，就可以得到N个点的FFT结果。为了方便进行FFT 运算，通常N取2的整数次方。 假设采样频率为Fs，信号频率F，采样点数为N。那么FFT 之后结果就是一个为N点的复数。每一个点就对应着一个频率 点。这个点的模值，就是该频率值下的幅度特性。具体跟原始 信号的幅度有什么关系呢？假设原始信号的峰值为A，那么FFT 的结果的每个点（除了第一个点直流分量之外）的模值就是A 的N/2倍。而第一个点就是直流分量，它的模值就是直流分量 的N倍。而每个点的相位呢，就是在该频率下的信号的相位。 第一个点表示直流分量（即0Hz），而最后一个点N的再下一个 点（实际上这个点是不存在的，这里是假设的第N+1个点，也 可以看做是将第一个点分做两半分，另一半移到最后）则表示 采样频率Fs，这中间被N-1个点平均分成N等份，每个点的频率 依次增加。例如某点n所表示的频率为：Fn=(n-1)*Fs/N。 由上面的公式可以看出，Fn所能分辨到频率为为Fs/N，如果 采样频率Fs为1024Hz，采样点数为1024点，则可以分辨到1Hz。1024Hz的采样率采样1024点，刚好是1秒，也就是说，采样1秒时间的信号并做FFT，则结果可以分析到1Hz，如果采样2秒时 间的信号并做FFT，则结果可以分析到0.5Hz。如果要提高频率 分辨力，则必须增加采样点数，也即采样时间。频率分辨率和 采样时间是倒数关系。 假设FFT之后某点n用复数a+bi表示，那么这个复数的模就是 An=根号a*a+b*b，相位就是Pn=atan2(b,a)。根据以上的结果， 就可以计算出n点（n≠1，且n<=N/2）对应的信号的表达式为： An/(N/2)*cos(2*pi*Fn*t+Pn)，即2*An/N*cos(2*pi*Fn*t+Pn)。 对于n=1点的信号，是直流分量，幅度即为A1/N。 快速傅里叶变换的原理及其应用 ————————————————————————————————作者：————————————————————————————————日期： 快速傅里叶变换的原理及其应用 摘要: 快速傅氏变换（FFT ），是离散傅氏变换的快速算法，它是根据离散傅氏变换的奇、偶、虚、实等特性，对离散傅立叶变换的算法进行改进获得的。它对傅氏变换的理论并没有新的发现，但是对于在计算机系统或者说数字系统中应用离散傅立叶变换，可以说是进了一大步。 傅里叶变换的理论与方法在“数理方程”、“线性系统分析”、“信号处理、仿真”等很多学科领域都有着广泛应用,由于计算机只能处理有限长度的离散的序列,所以真正在计算机上运算的是一种离散傅里叶变换. 虽然傅里叶运算在各方面计算中有着重要的作用，但是它的计算过于复杂，大量的计算对于系统的运算负担过于庞大，使得一些对于耗电量少，运算速度慢的系统对其敬而远之，然而，快速傅里叶变换的产生，使得傅里叶变换大为简化，在不牺牲耗电量的条件下提高了系统的运算速度，增强了系统的综合能力，提高了运算速度，因此快速傅里叶变换在生产和生活中都有着非常重要的作用，对于学习掌握都有着非常大的意义。 关键词：快速傅氏变换；图像处理；matlab 前言： 傅里叶变换在信号处理中具有十分重要的作用，但是基于离散时间的傅里叶变换具有很大的时间复杂度，根据傅里叶变换理论，对一个有限长度且长度为N 的离散信号，做傅里叶变换的时间复杂度为)(2N O ，当N 很大时τ，其实现的时间是相当惊人的（比如当N 为410时，其完成时间为τ810 (τ为计算机的时钟周期） ），故其实现难度是相当大的，同时也严重制约了DFT 在信号分析中的应用，故需要提出一种快速的且有效的算法来实现。 正是鉴于DFT 极其复杂的时间复杂度，1965年..JWCooley 和..JWTukey 巧妙地利用 NW 因子的周期性和对称性，提出了一个DFT 的快速算法，即快速傅里叶变换（FFT ），从 而使得DFT 在信号处理中才得到真正的广泛应用。 傅立叶变化的原理； （1）原理 《数字信号处理》课程论文 离散傅里叶变换(DFT) 算法简介及其应用举例 姓名：安昱 学号：12011243986 专业：通信工程 班级：2011级（1）班 指导老师：武永峰 学院：物理电气信息学院 完成日期：2013年11月11 日 离散傅里叶变换(DFT)算法简介及其应用举例 （安昱12011243986 2011级1班） [摘要]：离散傅立叶变换（DFT）实现了信号首次在频域表示的离散化，使得频域也能够用 计算机进行处理。并且这种DFT变换可以有多种实用的快速算法。使信号处理在时、频域的处理和转换均可离散化和快速化。最后就该项目做了总结。 [关键词]DFT算法Matlab语言频域采样DFT应用 一、DFT算法原理 快速傅氏变换（FFT）是离散傅氏变换的快速算法，它是根据离散傅氏变换的奇、偶、虚、实等特性，对离散傅立叶变换的算法进行改进获得的。它对傅氏变换的理论并没有新的发现，但是对于在计算机系统或者说数字系统中应用离散傅立叶变换，可以说是进了一大步。 设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变换。这样变换以后，总的运算次数就变成N+2（N/2）2=N+N2/2。继续上面的例子，N=1024时，总的运算次数就变成了525312次，节省了大约50%的运算量。而如果我们将这种“一分为二”的思想不断进行下去，直到分成两两一组的DFT运算单元，那么N点的DFT变换就只需要Nlog2N次的运算，N在1024点时，运算量仅有10240次，是先前的直接算法的1%，点数越多，运算量的节约就越大，这就是FFT的优越性。 二、FFT的算法原理 FFT算法的输出X(K)为自然顺序，但为了适应原位计算，其输入序列不是按x(n)的自然顺序排序，这种经过M-1次奇偶抽选后的排序为序列的倒序。因此，在运算之前应先对序列x(n)进行倒序。倒序的规律就是把顺序数的二进制位倒置，即可得到倒序值。倒序数是在M位二进制数最高位加一，逢2向右进位。M 位二进制数最高位的权值为N/2，且从左到右二进制位的权值依次为你N/4， 快速傅里叶变换的原理及其应用 摘要: 快速傅氏变换（FFT），是离散傅氏变换的快速算法，它是根据离散傅氏变换的奇、偶、虚、实等特性，对离散傅立叶变换的算法进行改进获得的。它对傅氏变换的理论并没有新的发现，但是对于在计算机系统或者说数字系统中应用离散傅立叶变换，可以说是进了一大步。傅里叶变换的理论与方法在“数理方程”、“线性系统分析”、“信号处理、仿真”等很多学科领域都有着广泛应用,由于计算机只能处理有限长度的离散的序列,所以真正在计算机上运算的是一种离散傅里叶变换. 虽然傅里叶运算在各方面计算中有着 重要的作用，但是它的计算过于复杂，大量的计算对于系统的运算负担过于庞大，使得一些对于耗电量少，运算速度慢的系统对其敬而远之，然而，快速傅里叶变换的产生，使得傅里叶变换大为简化，在不牺牲耗电量的条件下提高了系统的运算速度，增强了系统的综合能力，提高了运算速度，因此快速傅里叶变换在生产和生活中都有着非常重要的作用，对于学习掌握都有着非常大的意义。 关键词：快速傅氏变换；图像处理；matlab 前言： 傅里叶变换在信号处理中具有十分重要的作用，但是基于离散时间的傅里叶变换具有很大的时间复杂度，根据傅里叶变换理论，对一个有限长度且 长度为N的离散信号，做傅里叶变换的时间复杂度为) O，当N很大时τ， (2 N 其实现的时间是相当惊人的（比如当N为4 10 10时，其完成时间为τ8 (τ为计算机的时钟周期）），故其实现难度是相当大的，同时也严重制约了DFT在信号分析中的应用，故需要提出一种快速的且有效的算法来实 现。正是鉴于DFT极其复杂的时间复杂度，1965年..JWCooley 和..JWTukey巧妙地利用 NW因子的周期性和对称性，提出了一个DFT的快速算法，即快速傅里叶变换（FFT），从而使得DFT在信号处理中才得到真正的广泛应用。 傅立叶变化的原理； （1）原理 两点DFT简化 假设输入为x[0],x[1];输出为X[0],X[1]. 伪代码如下: // ------------------------------------------------------------------ #define N 2 #define PI 3.1415926 //tw为旋转因子 // ------------------------------------------------------------------ int i, j for(i=0, X[i]=0.0; i 8点FFT流程图(Layer表示层, gr表示当前层的颗粒)下面以LayerI为例. LayerI部分,具有4个颗粒,每个颗粒2个输入 (注意2个输入的来源,由时域信号提供) 将输入x[k]分为两部分x_r[k], x_i[k].具有实部和虚部,时域信号本没有虚部的,因此可以让x_i[k]为0.因为LayerII, LayerIII的输入是复数,为了编码统一而强行分的.当然编码时可以判断当前层是否为1来决定是否分. 旋转因子tw = cos(2*PI*k/N)-j*sin(2*PI*k/N);也可以分为实部和虚部,令其为tw_r, tw_i; 则tw = tw_r - j*tw_i; X[k] = (x_r[k] + j*x_i[k]) + (tw_r–j*tw_i) * (x_r[k+N/2]+j*x_i[k+N/2]) 则 X_R[k] = x_r[k] + tw_r*x_r[k+N/2] + tw_i*x_i[k+N/2]; X_I[k] = x_i[k] - tw_i*x_r[k+N/2] + tw_r*x_i[k+N/2];快速傅里叶变换(FFT)原理及源程序
FFT相关原理及使用注意事项
实验八 快速傅立叶变换(FFT)实验
fft_原理详解
快速傅里叶变换(FFT)原理及源程序
FFT原理与实现
快速傅里叶变换原理及其应用(快速入门)
FFT算法原理
快速傅里叶变换(FFT)的原理及公式
FFT原理及实现
串行FFT递归算法(蝶式递归计算原理)求傅里叶变换
快速傅里叶变换FFT原理与实现
基4fft原理及matlab实现
fft_原理详解
快速傅里叶变换的原理及其应用
DFT算法原理、FFT的算法原理
快速傅里叶变换的原理及其应用
FFT原理及C++实现
相关主题
文本预览