FFT的计算机实现
- 格式:doc
- 大小:1016.50 KB
- 文档页数:14
快速傅立叶变换〔FFT〕算法试验一.试验目的1.加深对DFT 算法原理和根本性质的理解;2.生疏FFT 算法原理和FFT 子程序的应用;3.学习用FFT 对连续信号和时域信号进展谱分析的方法,了解可能消灭的分析误差及其缘由,以便在实际中正确应用FFT。
二.试验设备计算机,CCS 3.1 版软件,E300 试验箱,DSP 仿真器,导线三.根本原理1.离散傅立叶变换DFT 的定义:将时域的采样变换成频域的周期性离散函数,频域的采样也可以变换成时域的周期性离散函数,这样的变换称为离散傅立叶变换,简称DFT。
2.FFT 是DFT 的一种快速算法,将DFT 的N2 步运算削减为〔N/2〕logN 步,极大2的提高了运算的速度。
3.旋转因子的变化规律。
4.蝶形运算规律。
5.基2FFT 算法。
四.试验步骤1.E300 底板的开关SW4 的第1 位置ON,其余置OFF。
其余开关不用具体设置。
2.E300 板子上的SW7 开关的第1 位置OFF,其余位置ON3.阅读本试验所供给的样例子程序;4.运行CCS 软件,对样例程序进展跟踪,分析结果;记录必要的参数。
5.填写试验报告。
6.供给样例程序试验操作说明A.试验前预备用导线连接“Signal expansion Unit”中2 号孔接口“SIN”和“A/D 单元”的2 号孔接口“AD_IN0”。
〔试验承受的是外部的AD模块〕B.试验1.正确完成计算机、DSP 仿真器和试验箱的连接后,系统上电。
2.启动CCS3.1,Project/Open 翻开“algorithm\01_fft”子名目下“fft.pjt”工程文件;双击“fft.pjt”及“Source”可查看各源程序;加载“Debug\fft.out”;3.单击“Debug\Go main”进入到主程序,在主程序“flag=0;”处设置断点;4.单击“Debug \ Run”运行程序,或按F5 运行程序;程序将运行至断点处停顿;5.用View / Graph / Time/Frequency 翻开一个图形观看窗口;设置该观看图形窗口变量及参数;承受双踪观看在启始地址分别为px 和pz,长度为128,数值类型为16 位整型,p x:存放经A/D 转换后的输入信号;p z:对该信号进展FFT 变换的结果。
fft 计算频率和幅值傅里叶变换(Fourier Transform)是一种将时域信号转换为频域信号的数学工具,可以用于分析信号的频率成分和幅值。
离散傅里叶变换(Discrete Fourier Transform,DFT)是傅里叶变换在数字信号处理中的离散形式,可以通过计算机算法进行实现。
快速傅里叶变换(Fast Fourier Transform,FFT)是一种高效计算离散傅里叶变换的算法,可以在计算机中快速地对信号进行频域分析。
要计算信号的频率和幅值,可以使用FFT算法进行以下步骤:1. 选择信号并确定采样率:首先,选择要分析的信号。
这可以是任何连续或离散的信号,例如音频、振动传感器数据等。
确定信号的采样率,即每秒采集的样本数。
2. 做零填充(Zero-padding):为了获得更准确的频率分辨率,可以在信号末尾添加零填充。
零填充是指在信号的末尾添加一些零值样本,通常将信号长度扩展到2的幂次方。
这对于使用FFT算法进行频谱分析时是有益的。
3. 应用FFT算法:使用FFT算法对零填充后的信号进行离散傅里叶变换。
FFT算法将信号从时域转换为频域,并返回频域上的复数结果。
FFT算法的输入是时域信号的离散样本,输出是频域上的复数结果,包括幅度和相位信息。
4. 计算频率和幅值:从FFT的输出中提取频率和幅值信息。
FFT 的输出是一个复数数组,其中每个元素对应于信号的一个频率分量。
频率分量的幅度可以通过计算复数模数(或绝对值)来获得,表示该频率分量的强度或振幅。
频率分量的相位可以通过计算复数的相角来获得,表示该频率分量的相位偏移。
5. 频率分辨率和频率范围:频率分辨率是指能够分辨两个频率分量之间的最小频率差。
在FFT中,频率分辨率取决于采样率和信号长度。
频率范围是指在频域上可以得到的有效频率范围,它与采样率和信号长度有关。
频率范围从零频率开始,一直延伸到采样率的一半。
通过对FFT输出进行适当的缩放和归一化,可以将幅度转换为能量表示或相对强度表示。
快速傅立叶变换(FFT )的实现一、实验目的1.了解FFT 的原理及算法;2.了解DSP 中FFT 的设计及编程方法;3.熟悉FFT 的调试方法;二、实验原理FFT 是一种高效实现离散付立叶变换的算法,把信号从时域变换到频域,在频域分析处理信息。
对于长度为N 的有限长序列x (n ),它的离散傅里叶变换为:(2/)j N nk N W e π-=,称为旋转因子,或蝶形因子。
在x (n )为复数序列的情况下,计算X (k ):对某个k 值,需要N 次复数乘法、(N -1)次复数加法;对所有N 个k 值,需要2N 次复数乘法和N (N -1)次复数加法。
对于N 相当大时(如1024)来说,直接计算它的DFT 所作的计算量是很大的,FFT 的基本思想在于: 利用2()j nk N N W e π-=的周期性即:k N k N N W W +=对称性:/2k k N N N W W +=-将原有的N 点序列分成两个较短的序列,这些序列的DFT 可以很简单的组合起来得到原序列的DFT 。
按时间抽取的FFT ——DIT FFT 信号流图如图5.1所示:图5.1 时间抽取的FFT —DIT FFT 信号流图FFT 算法主要分为以下四步。
第一步 输入数据的组合和位倒序∑=-=10)()(N n nk N W n x k X把输入序列作位倒序是为了在整个运算最后的输出中得到的序列是自然顺序。
第二步 实现N 点复数FFT第一级蝶形运算;第二级蝶形运算;第三级至log2N 级蝶形运算;FFT 运算中的旋转因子N W 是一个复数,可表示:为了实现旋转因子N W 的运算,在存储空间分别建立正弦表和余弦表,每个表对应从0度到180度,采用循环寻址来对正弦表和余弦表进行寻址。
第三步 功率谱的计算X (k )是由实部()R X k 和虚部()I X k 组成的复数:()()()R I X k X k jX k =+;计算功率谱时只需将FFT 变换好的数据,按照实部()R X k 和虚部()I X k 求它们的平方和,然后对平方和进行开平方运算。
fft函数
FFT(快速傅里叶变换)是一种实现DFT(离散傅里叶变换)的快速算法,是利用复数形式的离散傅里叶变换来计算实数形式的离散傅里叶变换,matlab中的fft()函数是实现该算法的实现。
MATLAB它将数值分析、矩阵计算、科学数据可视化以及非线性动态系统的建模和仿真等诸多强大功能集成在一个易于使用的视窗环境中,为科学研究、工程设计以及必须进行有效数值计算的众多科学领域提供了一种全面的解决方案,并在很大程度上摆脱了传统非交互式程序设计语言(如C、Fortran)的编辑模式,代表了当今国际科学计算软件的先进水平。
快速傅里叶变换, 即利用计算机计算离散傅里叶变换(DFT)的高效、快速计算方法的统称,简称FFT。
快速傅里叶变换是1965年由J.W.库利和T.W.图基提出的。
采用这种算法能使计算机计算离散傅里叶变换所需要的乘法次数大为减少,特别是被变换的抽样点数N越多,FFT算法计算量的节省就越显著。
FFT及IFFTC语言实现FFT(快速傅里叶变换)和IFFT(快速傅里叶逆变换)是傅里叶变换在计算机科学和信号处理中的高效实现方法。
它们广泛应用于图像处理、音频处理、通信系统等领域。
下面将详细介绍FFT和IFFT的C语言实现。
首先,让我们了解一下DFT(离散傅里叶变换)。
DFT将一个离散的时间域序列转换为离散的频域序列,其定义如下:其中,N表示序列的长度,x(n)是输入序列,X(k)是输出序列。
FFT是DFT的一种高效实现方法,它利用了序列的对称性质,将操作的复杂度从O(N^2)降低到O(NlogN)。
IFFT则是FFT的逆过程,可以将频域序列恢复为时间域序列。
以下是FFT的C语言实现代码:```c#include <stdio.h>#include <math.h>typedef structdouble real;double imag;result.real = a.real * b.real - a.imag * b.imag;result.imag = a.real * b.imag + a.imag * b.real;return result;result.real = a.real + b.real; result.imag = a.imag + b.imag; return result;result.real = a.real - b.real; result.imag = a.imag - b.imag; return result;if (N <= 1)return;}for (int i = 0; i < N / 2; i++) even[i] = x[2 * i];odd[i] = x[2 * i + 1];}fft(even, N / 2);fft(odd, N / 2);for (int k = 0; k < N / 2; k++)x[k] = add(even[k], t);x[k + N / 2] = subtract(even[k], t); }//完整的IFFT实现代码略,与FFT相似int mai//输入序列x(n)int N = sizeof(x) / sizeof(x[0]);//调用FFTfft(x, N);//输出频域序列X(k)for (int k = 0; k < N; k++)printf("X(%d) = %.2f + %.2fi\n", k, x[k].real, x[k].imag);}//调用IFFT//...return 0;```IFFT的实现与FFT类似,只需将其中几处计算公式略作变换即可。
FFT的算法原理应用FFT(Fast Fourier Transform)是一种高效计算离散傅里叶变换(Discrete Fourier Transform,DFT)的算法。
通过使用FFT算法,可以将DFT的计算时间从O(N^2)降低到O(NlogN),其中N是离散序列的长度。
FFT的算法原理基于Radix-2分治策略,将一个长序列分解为两个较短序列,并重复此过程,直到仅剩两个元素相乘为止。
FFT的算法主要应用于信号处理和频谱分析等领域。
其在频谱分析中的应用可以帮助我们了解信号的频率内容以及频率分量的强度。
在信号处理中,FFT可以用于将时域数据转换为频域数据,使得信号处理更加简化和高效。
下面将详细介绍FFT的算法原理和主要应用。
1.FFT算法原理:具体步骤如下:1)通过对输入序列进行重新排列,将序列按照奇偶位进行分组,分为两个长度为N/2的子序列。
2)对这两个子序列分别进行DFT计算,得到两个长度为N/2的频域序列。
3)将这两个序列分别与旋转因子进行乘积,得到两个长度为N/2的频域子序列。
4)将这两个频域子序列连接起来,得到长度为N的频域序列。
5)递归地将这个过程应用于每个子序列,直到序列长度为2,此时不需要再进行分解。
6)将分解后的频域序列进行合并和重排,得到最终的频域序列。
通过这种分治策略,FFT能够将DFT的复杂度从O(N^2)降低到O(NlogN),大大提高了计算效率。
2.FFT的应用:(1)频谱分析:FFT算法可以将时域信号转换为频域信号,分析信号的频率成分和强度。
通过FFT,可以得到信号的频谱信息,帮助我们了解信号的频率特点和分布情况。
常见的应用包括音频分析、图像处理、通信信号分析等。
(2)信号处理:FFT在信号处理中广泛应用,例如滤波、模式识别、降噪等。
通过将信号转换为频域,在频域进行处理后再进行逆变换,可以实现对信号的特定频率的增强或者抑制。
(3)图像处理:FFT在图像处理中的应用主要是基于频率域滤波。
FFT原理与实现FFT(Fast Fourier Transform)是一种快速傅里叶变换算法,用于将时域信号转换为频域信号。
它是傅里叶变换的一种优化算法,通过分解和递归计算可以大幅缩短计算时间。
FFT在信号处理、图像处理、音频处理、通信等领域得到广泛应用。
FFT的原理可以通过以下几个步骤进行解释:1.输入信号:FFT的输入是时域信号,即在时间上采样得到的离散信号。
通常,这个信号是一个在计算机中存储的数字序列。
2.傅里叶变换:FFT将时域信号转换为频域信号,使用了傅里叶变换的数学公式。
傅里叶变换将信号分解为一系列正弦和余弦函数的和,每个函数对应一种频率成分。
3.分解和递归计算:FFT将输入信号分解成多个较小的子问题,然后通过递归计算解决。
这个分解过程是基于信号的二进制反转,将数字序列重新排列为二进制反转的次序。
4.合并和求和:计算出每个子问题的频域结果后,FFT通过合并和求和操作将子问题的结果合并为最终的频域结果。
这个操作需要对结果进行两两合并,并加上对应的旋转因子,以得到更高频率成分的结果。
5.输出频域信号:最终得到的频域结果可以表示为幅度和相位,幅度表示在对应频率上的能量大小,相位表示信号在对应频率上的相位关系。
FFT的实现可以分为两种常见的算法:基2FFT和基4FFT。
基2FFT是将输入信号按照二进制表示进行分解和合并,而基4FFT是将输入信号按照四进制进行分解和合并。
两种算法时间复杂度相同,但基4FFT的计算量更少。
在实现FFT时,需要注意以下几点:1.选择适当的FFT长度:FFT的计算结果是根据输入信号的采样个数来计算的,因此需要根据输入信号的长度选择合适的FFT长度。
一般来说,FFT长度选择为2的幂次方,这样可以充分利用FFT的分解和合并特性。
2.使用合适的内存分配:FFT需要使用数组来存储输入信号和计算结果,因此需要根据输入信号的长度分配合适的内存空间。
同时,为了提高计算效率,可以使用一些优化策略,如使用位操作和循环展开。
数字信号处理实验报告 实验二 FFT 算法的MATLAB 实现(一)实验目的:理解离散傅立叶变换时信号分析与处理的一种重要变换,特别是FFT 在数字信号处理中的高效率应用。
(二)实验原理:1、有限长序列x(n)的DFT 的概念和公式:⎪⎪⎩⎪⎪⎨⎧-≤≤=-≤≤=∑∑-=--=101010)(1)(10)()(N k kn N N n kn N N n W k x N n x N k W n x k x)/2(N j N eW π-=2、FFT 算法调用格式是 X= fft(x) 或 X=fft(x,N)对前者,若x 的长度是2的整数次幂,则按该长度实现x 的快速变换,否则,实现的是慢速的非2的整数次幂的变换;对后者,N 应为2的整数次幂,若x 的长度小于N ,则补零,若超过N ,则舍弃N 以后的数据。
Ifft 的调用格式与之相同。
(三)实验内容1、题一:若x(n)=cos(n*pi/6)是一个N=12的有限序列,利用MATLAB 计算它的DFT 并画出图形。
源程序: clc; N=12; n=0:N-1; k=0:N-1;xn=cos(n*pi/6); W=exp(-j*2*pi/N); kn=n'*kXk=xn*(W.^kn) stem(n,Xk); xlabel('k'); ylabel('Xk'); grid on ;也可用FFT 算法直接得出结果,程序如下: clc; N=12; n=0:N-1;xn=cos(n*pi/6);Xk=fft(xn,N); stem(n,Xk); xlabel('k'); ylabel('Xk'); grid on ;实验结果:24681012kX k分析实验结果:用DFT 和用FFT 对序列进行运算,最后得到的结果相同。
但用快速傅立叶变换的运算速度可以快很多。
2、题二:一被噪声污染的信号,很难看出它所包含的频率分量,如一个由50Hz 和120Hz 正弦信号构成的信号,受均值随机噪声的干扰,数据采样率为1000Hz ,通过FFT 来分析其信号频率成分,用MA TLAB 实现。
定点fft matlab代码1.引言1.1 概述在文章的引言部分,我们首先要概述一下所要讨论的主题,即定点FFT (快速傅里叶变换)算法的Matlab代码实现。
定点FFT算法是一种计算机快速傅里叶变换的算法。
傅里叶变换是一种重要的信号处理工具,在很多领域中都有广泛的应用,如通信、图像处理、音频处理等。
传统的傅里叶变换算法复杂度较高,需要进行大量的复数运算,导致计算时间较长。
而快速傅里叶变换算法通过巧妙地利用对称性和周期性的特点,在计算复杂度上有很大的优势,能够快速地对信号进行频域分析。
Matlab是一种功能强大的数学软件,广泛应用于科学计算、数据分析等领域。
在Matlab中,有很多已经实现好的函数可以方便地进行FFT 计算。
然而,这些函数通常是基于浮点数运算的,即使用双精度浮点数进行计算。
在某些应用场景下,我们可能需要使用定点数进行傅里叶变换,如在一些嵌入式系统中由于硬件限制无法支持浮点数运算。
因此,我们需要对FFT算法进行定点化的实现。
本文将介绍定点FFT算法的原理和在Matlab中的实现。
在实现过程中,我们将讨论如何进行定点数的表示和运算,并给出详细的代码实现。
同时,我们还将分析定点FFT算法在不同精度下的计算性能和结果精度,并进行相关的讨论和总结。
通过本文的阅读,读者将能够了解到定点FFT算法的原理和编程实现,以及在Matlab中如何使用定点数进行傅里叶变换。
这对于需要在嵌入式系统中进行傅里叶变换的工程师和研究人员来说,将是一份有价值的参考资料。
1.2 文章结构文章将分为三个主要部分:引言、正文和结论。
在引言部分,我们将给出本文的概述,简要介绍定点FFT算法,并明确文章的目的。
首先,我们将解释FFT算法的基本原理以及其在信号处理中的应用。
接着,我们将介绍定点FFT算法的原理和特点,包括其对计算资源的要求和性能优化方面的研究。
最后,我们将明确文章的目的,即在Matlab中实现定点FFT算法,并对实验结果进行分析与讨论。
快速傅立叶变换(FFT )的计算机实现邓凯 电气07061.FFT 运算的原理DFT 运算过程中如果有限长序列的长度很长时,即N 很大时,在运算过程中所做的乘法和加法运算将很多。
即便是采用高速的计算机进行运算,也很难达到信息实时处理的要求。
由库利、图基等科学家提出的快速傅里叶变换FFT 算法大大减少了运算过程中的乘法和加法次数,适合信息处理对实时性的要求,从而得到广泛应用。
按时间抽取的FFT 算法的基本思想是将输入的有限长序列首先分成奇数序列和偶数序列,分别计算出奇、偶序列的DFT ,然后根据DFT 的周期性和对称性质,将其化简,接着将已分成的奇、偶序列再次分别划分成奇、偶序列,即前面的奇序列按其长度再划分成奇数序列和偶数序列;前面的偶序列按其长度再划分成奇数序列和偶数序列。
分别计算其DFT 后,再按上述方法进行化简。
如此反复,直至被划分成的奇、偶序列长度为1为止。
1.1基二FFT 算法原理若设X(k)=DFT[x(n)],且0≤n, k ≤N-1,N 为偶数(如果N 为奇数,则添上一个零值点使长度N 为偶数)。
把它分为奇序列和偶序列: 21i) x((i) x = )i x((i) x 122+= )120(-≤≤Ni 又 (i )] D F T [x (k ) (i )] X D F T [x (k )X 2211==)120(-≤≤Ni则有 ⎪⎩⎪⎨⎧-=++=kkNNW k X k X N k X W k X k X k X )()()2()()()(2121 )120(-≤≤N i 其中,][1k X ][2k X 分别是][1i x 和][2i x 的2N点DFT 即要求][k X 的值仅仅只需要求][1k X ][2k X 在)120(-≤≤Ni 部分的值即可。
这样,我们就可以将][k X n 一直分为奇序列和偶序列来求其值,直到奇偶序列长度为1为止。
1.2蝶形图对于FFT 运算,我们通常用蝶形图来表示比如⎪⎩⎪⎨⎧-=++=kkN N W k X k X N k X W k X k X k X )()()2()()()(2121 我们表示为[1X (2X kN W k X k X k )()()(21+=k N W k X k X Nk X )()(2(21-=+图1.21下面以N=8为例,运用FFT 算法原理计算DFT ,如图 1.22、图1.23、图1.24所示。
图1.22 按时间抽取,将一个N=8点DFT分解为两个 4 点DFTX(0)图1.23 按时间抽取,将一个N=8点DFT分解为四个 2点DFTX(0)-14图1.23 按时间抽取,将一个N=8点DFT 分解为8个 1点DFT21.3 FFT 算法特点(1)原位运算。
从图2.9的FFT 运算流图中我们可以看出,这种运算是很有规律可循的。
其中的每级(每列)运算都是由N/2个蝶形运算所构成,如式(1.31)所示。
⎩⎨⎧-=+=----kNm m m Nm m m W j X k X j X W j X k X k X )()()()()()(11110(式1.31) 其中,m 表示第m 列迭代,k, j 表示数据所在的行数,该式所对应的蝶形运算结构如图1.31所示。
图1.31x x rNm m m W j x k x k )()()(11--+=rNm m m W j x k x j )()()(11---=由图1.31的流程图我们可以看出,对于任一个蝶形中的任何两个节点k 和j ,经过运算后所得的结果作为下一列(级)k, j 两个节点的节点变量,同其他节点变量无关,这样就可以采用原位运算方法。
即某一列的N 个数据送入存储器经蝶形运算后,结果为另一级(列)数据,而这些数据结果可以以蝶形为单位仍然存储在这同一组存储器中,直到最后输出,中间无需另加存储器。
也就是说,蝶形运算的两个输出值仍放回到蝶形的两个输入所在的存储器中(即将原先的输入值覆盖掉)。
前一级蝶形运算完成后才进行下一级蝶形原位运算,因而整个运算结构由于采用这种原位运算而大大减少了存储单元的个数,有效地降低了成本。
(2)倒位序规律。
由图2.9我们还可看出,由于采用了原位计算方法,FFT 的输出X(k)是按k 值的自然顺序排列在存储单元中的,即排列序列为X(0), X(1), …, X(7),而输入的时间序列都不是按照自然顺序排列,而是按x(0), x(4), x(2), …, x(7),输入到存储单元中。
这种排列初看起来像是“混乱无序”的,但实质是有序的。
如果我们用二进制数来表示0~7个数,设为n(n2 n1 n0)2,然后再来观察n (n0 n1 n2)2,如表1.1所示。
由表我们可以看出,只要将2)(012n n n n 转换成22)(10n n n n ,即将二进制的最高位与最低位相交换、次高位与次低位相交换……所得的n 就是以自然顺序排列的,故通常n 的这种排列规律我们称之为倒位序。
这两条性质给我们编程带来了极大的方便。
2.C 语言实现FFT 算法2.1算法基本原理如上所述,FFT 算法给编程带来了极大大方便,对比于DFT 不仅运算次数大大降低,就连所需的内存空间在运算中也不会增加。
主要算法有两部分,一个是倒二进制排序,还有最重要的蝶形算法。
倒二进制算法对数组下标进行倒二进制运算,把输入数组按倒二进制排序即可,算法很简单,具体可参见代码蝶形算法核心即为如下公式,只需进行几次循环分级进行蝶形算法即可,具体可参见代码。
⎩⎨⎧-=+=----kN m m m Nm m m W j X k X j X W j X k X k X )()()()()()(111102.2程序流程图否FFT 变换核心流程其中的N 表示运算的总级数,即将FFT 时域序列分出奇偶序列,直至每个序列点数为1所需的次数。
若有1024个点,则N 为10。
3.FFT的Matlab实验3.1 改变序列的长度N对计算结果的影响x=0.5sin(2*pi*15*t)+2*sin(2*pi*40*t)。
采样频率fs=100Hz,分别绘制N=128、1024点幅频图。
在Matlab中输入如下代码:>> clf;fs=100;N=128; %采样频率和数据点数n=0:N-1;t=n/fs; %时间序列x=0.5*sin(2*pi*15*t)+2*sin(2*pi*40*t); %信号y=fft(x,N); %对信号进行快速Fourier变换mag=abs(y); %求得Fourier变换后的振幅f=n*fs/N; %频率序列subplot(2,2,1),plot(f,mag); %绘出随频率变化的振幅xlabel('频率/Hz');ylabel('振幅');title('N=128');grid on;subplot(2,2,2),plot(f(1:N/2),mag(1:N/2)); %绘出Nyquist频率之前随频率变化的振幅xlabel('频率/Hz');ylabel('振幅');title('N=128');grid on;%对信号采样数据为1024点的处理fs=100;N=1024;n=0:N-1;t=n/fs;x=0.5*sin(2*pi*15*t)+2*sin(2*pi*40*t); %信号y=fft(x,N); %对信号进行快速Fourier变换mag=abs(y); %求取Fourier变换的振幅f=n*fs/N;subplot(2,2,3),plot(f,mag); %绘出随频率变化的振幅xlabel('频率/Hz');ylabel('振幅');title('N=1024');grid on;subplot(2,2,4)plot(f(1:N/2),mag(1:N/2)); %绘出Nyquist频率之前随频率变化的振幅xlabel('频率/Hz');ylabel('振幅');title('N=1024');grid on;>>x=0.5*sin(2*pi*15*t)+2*sin(2*pi*40*t)。
采样频率fs=100Hz,当N=128、1024点幅频图如下fs=100Hz,Nyquist频率为fs/2=50Hz。
整个频谱图是以Nyquist频率为对称轴的。
并且可以明显识别出信号中含有两种频率成分:15Hz和40Hz。
由此可以知道FFT变换数据的对称性。
因此用FFT对信号做谱分析,只需考察0~Nyquist频率范围内的福频特性。
若没有给出采样频率和采样间隔,则分析通常对归一化频率0~1进行。
另外,振幅的大小与所用采样点数有关,采用128点和1024点的相同频率的振幅是有不同的表现值,但在同一幅图中,40Hz 与15Hz振动幅值之比均为4:1,与真实振幅0.5:2是一致的。
为了与真实振幅对应,需要将变换后结果乘以2除以N。
由上图可得,点数取得越多,由FFT所得的频谱越精确。
4.总结傅里叶变换在物理学、数论、组合数学、信号处理、概率论、统计学、密码学、声学、光学、海洋学、结构动力学等领域都有着广泛的应用(例如在信号处理中,傅里叶变换的典型用途是将信号分解成幅值分量和频率分量)。
而FFT则实现了傅里叶变换在计算机领域的广泛应用,学习FFT的原理和实现,对后续课程如电力系统的分析等是很有必要的。
本次课程设计,我选择了FFT在计算机上的实现这一个题目,采用C语言编程实现FFT算法,同时用Matlab对C程序进行验证和分析。
在完成课设过程中,不仅弄清楚了FFT的原理,还学到了Matlab相关的知识,C语言编程能力也得到了较大的提高。
附录程序清单1.1复数的算法#include "complex.h"CComplex::CComplex(void){Realpart = 0;ImagPart = 0;}CComplex::CComplex(double Real,double Imag){Realpart = Real;ImagPart = Imag;}CComplex::CComplex(double real){Realpart = real;ImagPart = 0;}CComplex::~CComplex(void){}//////////////////////////////////////////////////////////////////////////////////* 重载数学运算符*////////////////////////////////////////////////////////////////////////////////////+重载CComplex CComplex::operator +(CComplex OtherComplex){CComplex temp;temp.Realpart=Realpart+OtherComplex.Realpart;temp.ImagPart=ImagPart+OtherComplex.ImagPart;return temp;}//-重载CComplex CComplex::operator -(CComplex OtherComplex){CComplex temp;temp.Realpart=Realpart-OtherComplex.Realpart;temp.ImagPart=ImagPart-OtherComplex.ImagPart;return temp;}//*重载CComplex CComplex::operator*(CComplex OtherComplex){CComplex temp;temp.Realpart=Realpart*OtherComplex.Realpart-ImagPart*OtherComplex.ImagPart;temp.ImagPart=Realpart*OtherComplex.ImagPart+ImagPart*OtherComplex.Realpart;return temp;}1.2FFT运算#include "complex.h"#define FFTPOINT 8///////////////////////////////////////////////////////////////////////////////////////*输入时域取样数据*/void InPut(double* pData){double Ipt;for (int i=0 ;i<FFTPOINT; i++){cin>>Ipt;*(pData+i)=Ipt;}}/*输出变换后的频域数据,复数形式*/void Output(CComplex* pData){CComplex Out;cout<<"FFT变换后,频域数据为"<<endl;for (int i=0;i<FFTPOINT;i++){ Out=*(pData+i);if ((Out.Realpart)>=0){printf("%.4fj+%.4f\n",Out.ImagPart,Out.Realpart);}else{printf("%.4fj%.4f\n",Out.ImagPart,Out.Realpart);}}}/////////////////////////////////////////////////////////////////////////////////////// /*将输入的数据按倒二进制排序*/void FFTArry(double* pData1,CComplex* pData2){int i=0,m=2,t=0;for (N=0;m <= FFTPOINT;N++){m*=2;}for(i=0;i<FFTPOINT;i++){int temp=i;t=0;for (int j=N-1;j>=0;j--){if(temp/(int)(pow(2,j)))t+=pow(2,N-j-1);temp=temp%(int)(pow(2,j));}pData2[t].Realpart=*(pData1+i);}}////////////////////////////////////////////////////////////////////////////////////////////////////extern int N=0;void main(){void FFT(CComplex* pData);double FFTDataIn[FFTPOINT]={0.0};CComplex FFTDataOut[FFTPOINT];cout<<"请输入原数据,最大个数为:"<<FFTPOINT<<endl;InPut(FFTDataIn);//输入数据FFTArry(FFTDataIn,FFTDataOut);//将输入数据按到二进制排序FFT(FFTDataOut);//进行FFT变换Output(FFTDataOut);//输出}///////////////////////////////////////////////////////////////////////////////////////void FFT(CComplex* pData){int i,r,steplen,k,j;for(i=0;i<N;i++) /*实现N级运算*/{steplen=(int)(pow(2,i)); /*计算步长*/for(r=0;r<pow(2,i);r++){CComplex W(cos(2*PI*r/pow(2,i+1)),-sin(2*PI*r/pow(2,i+1)));/*旋转因子的值*/for(j=0;j<(FFTPOINT/pow(2,i+1));j++)/*r相同的值的组的计算*/{CComplex temp1,temp2;k=(int)( r+j*pow(2,i+1)); /*每一组成员的下标*/temp1=pData[k]+W*pData[k+steplen];temp2=pData[k]-W*pData[k+steplen];pData[k]=temp1;pData[k+steplen]=temp2;}}}}//////////////////////////////////////////////////////////////////////////////////////////////////。