实验2FFT算法实现
实验二FFT实现信号频谱分析
0
2
4
6
4
2
0
-2
-4
-6
-4
-20246四、试验环节
4. 试验内容2旳程序运营成果如下图所示：
60
30
40
20
20
10
0
0
-10 -5
0
5
10
-40 -20
0
20 40
30
80
60 20
40 10
20
0
-40 -20
0
20 40
0
-40 -20
0
20 40
四、试验环节
|X(k)| x(n)
5. 试验内容 3旳程序运营成果如下图所示：
fft 计算迅速离散傅立叶变换
fftshift
ifft
调整fft函数旳输出顺序，将零频 位置移到频谱旳中心
计算离散傅立叶反变换
fft函数：调用方式如下
y=fft(x)：计算信号x旳迅速傅立叶变换y。当x旳长度为 2旳幂时，用基2算法，不然采用较慢旳分裂基算法。
y=fft(x,n)：计算n点FFT。当length(x)>n时，截断x，不 然补零。
【例2-11】产生一种正弦信号频率为60Hz，并用fft函数 计算并绘出其幅度谱。
fftshift函数：调用方式如下 y=fftshift(x)：假如x为向量，fftshift(x)直接将x旳左右两 部分互换；假如x为矩阵（多通道信号），将x旳左上、右 下和右上、左下四个部分两两互换。 【例2-12】产生一种正弦信号频率为60Hz,采样率为1000Hz, 用fftshift将其零频位置搬到频谱中心。
以上就是按时间抽取旳迅速傅立叶变换
实验二的应用FFT对信号进行频谱分析
实验二的应用FFT对信号进行频谱分析引言：频谱分析是通过将连续信号转换为离散信号，根据信号在频域上的强度分布来分析信号的频谱特性。
其中，FFT（Fast Fourier Transform，快速傅里叶变换）是一种常见的频谱分析算法，可以高效地计算离散信号的傅里叶变换。
实验目的：本实验旨在使用FFT算法来对一个信号进行频谱分析，从而了解FFT 的原理和应用。
实验器材：-计算机-MATLAB软件实验步骤：1.准备信号数据：首先，需要准备一个信号数据用于进行频谱分析。
可以通过MATLAB 自带的函数生成一个简单的信号数据，例如生成一个正弦信号：```Fs=1000;%采样频率T=1/Fs;%采样时间间隔L=1000;%信号长度t=(0:L-1)*T;%时间向量S = 0.7*sin(2*pi*50*t) + sin(2*pi*120*t); % 生成信号，包含50Hz和120Hz的正弦波成分```其中，Fs为采样频率，T为采样时间间隔，L为信号长度，t为时间向量，S为生成的信号数据。
2.进行FFT计算：利用MATLAB提供的fft函数，对准备好的信号数据进行FFT计算，得到信号的频谱：```Y = fft(S); % 对信号数据进行FFT计算P2 = abs(Y/L); % 取FFT结果的模值，并归一化P1=P2(1:L/2+1);%取模值前一半P1(2:end-1) = 2*P1(2:end-1); % 对非直流分量进行倍频处理f=Fs*(0:(L/2))/L;%计算对应的频率```其中，Y为FFT计算的结果，P2为对应结果的模值，并进行归一化处理，P1为P2的前一半，f为对应的频率。
3.绘制频谱图：使用MATLAB的plot函数，将频率和对应的功率谱绘制成频谱图：```plot(f,P1)title('Single-Sided Amplitude Spectrum of S(t)')xlabel('f (Hz)')ylabel('，P1(f)，')```实验结果与分析：上述实验步骤通过MATLAB实现了对一个信号的频谱分析并绘制成频谱图。
实验二用FFT做谱分析实验报告
实验二用FFT做谱分析实验报告一、引言谱分析是信号处理中一个重要的技术手段，通过分析信号的频谱特性可以得到信号的频率、幅度等信息。
傅里叶变换是一种常用的谱分析方法，通过将信号变换到频域进行分析，可以得到信号的频谱信息。
FFT（快速傅里叶变换）是一种高效的计算傅里叶变换的算法，可以大幅减少计算复杂度。
本实验旨在通过使用FFT算法实现对信号的谱分析，并进一步了解信号的频谱特性。
二、实验目的1.理解傅里叶变换的原理和谱分析的方法；2.学习使用FFT算法对信号进行谱分析；3.通过实验掌握信号的频谱特性的分析方法。
三、实验原理傅里叶变换是将信号从时域转换到频域的一种数学变换方法，可以将一个非周期性信号分解为一系列正弦和余弦函数的叠加。
FFT是一种计算傅里叶变换的快速算法，能够在较短的时间内计算出信号的频谱。
在进行FFT谱分析时，首先需要对信号进行采样，然后利用FFT算法将采样后的信号转换到频域得到信号的频谱。
频谱可以用幅度谱和相位谱表示，其中幅度谱表示信号在不同频率下的幅度，相位谱表示信号在不同频率下的相位。
四、实验装置和材料1.计算机；2.信号发生器；3.数字示波器。
五、实验步骤1.连接信号发生器和示波器，通过信号发生器产生一个周期为1s的正弦信号，并将信号输入到示波器中进行显示；2.利用示波器对信号进行采样，得到采样信号；3.利用FFT算法对采样信号进行频谱分析，得到信号的频谱图。
六、实验结果[插入频谱图]从频谱图中可以清晰地看到信号在不同频率下的幅度和相位信息。
其中，频率为2Hz的分量的幅度最大，频率为5Hz的分量的幅度次之。
七、实验分析通过对信号的频谱分析，我们可以得到信号的频率分量和其对应的幅度和相位信息。
通过分析频谱图，我们可以得到信号中各个频率分量的相对强度。
在本实验中，我们可以看到频率为2Hz的分量的幅度最大，频率为5Hz的分量的幅度次之。
这说明信号中存在2Hz和5Hz的周期性成分，且2Hz的成分更为明显。
时间抽取的基2快速傅里叶变换FFT分析与算法实现
离散时间信号的基2快速傅里叶变换FFT （时间抽取）蝶形算法实现一、一维连续信号的傅里叶变换连续函数f(x)满足Dirichlet （狄利克雷）条件下，存在积分变换：正变换：2()()()()j ux F u f x e dx R u jI u π+∞--∞==+⎰ 反变换：2()()j ux f x F u e du π+∞-∞=⎰其中()()cos(2)R u f t ut dt π+∞-∞=⎰，()()sin(2)I u f t ut dt π+∞-∞=-⎰定义幅值、相位和能量如下：幅度：1222()()()F u R u I u ⎡⎤⎡⎤=+⎣⎦⎣⎦ 相位：()arctan(()/())u I u R u ϕ= 能量：22()()(E u R u I u =+）二、一维离散信号的傅里叶变换将连续信号对自变量进行抽样得到离散信号（理想冲击抽样脉冲），利用连续信号的傅里叶变换公式得到离散时间傅里叶变换DTFT ；再利用周期性进行频域抽样，得离散傅里叶变换DFT （详情参考任何一本《数字信号处理》教材）。
DFT 变换如下：正变换：12/0()(),0,1,2,1N j ux Nx F u f x eu N π--===-∑。
反变换：12/01()(),0,1,2,1N j ux Nu f x F u ex N Nπ-===-∑。
DFT 是信号分析与处理中的一种重要变换，因为计算机等数字设备只能存储和处理离散数据（时域、频域）。
因直接计算DFT 的计算量大（与变换区间长度N 的平方成正比，当N 较大时，计算量太大），所以在快速傅里叶变换(简称FFT)出现以前，直接用DFT 算法进行谱分析和信号的实时处理是不切实际的。
直到1965年发现了DFT 的一种快速算法（快速傅里叶变换，即FFT ）以后，情况才发生了根本的变化。
FFT 有时间抽取和频率抽取两种，下面介绍时间抽取FFT 。
三、时间抽取的基2快速傅里叶变换FFT令2j NN W eπ-=，则2jkm km NNWeπ-=称为旋转因子，把DFT 正变换改写为：1[][],0,1,1N km N k F m f k W m N -===-∑将函数记作x ，变换后为X ，则为：10[][],0,1,1N kmN k X m x k W m N -===-∑时间抽取的FFT 算法利用了旋转因子的三个性质：周期性、对称性和可约性。
实验二用DFT及FFT进行谱分析
实验二用DFT及FFT进行谱分析实验二将使用DFT（离散傅里叶变换）和FFT（快速傅里叶变换）进行谱分析。
在谱分析中，我们将探索如何将时域信号转换为频域信号，并观察信号的频谱特征。
首先，我们需要了解DFT和FFT的基本概念。
DFT是一种将时域信号分解为频域信号的数学方法。
它将一个离散时间序列的N个样本转换为具有N个频率点的频率谱。
DFT在信号处理和谱分析中被广泛应用，但它的计算复杂度为O(N^2)。
为了解决DFT的计算复杂度问题，Cooley和Tukey提出了FFT算法，它是一种使用分治策略的快速计算DFT的方法。
FFT算法的计算复杂度为O(NlogN)，使得谱分析在实际应用中更加可行。
在实验中，我们将使用Python编程语言和NumPy库来实现DFT和FFT，并进行信号的谱分析。
首先，我们需要生成一个具有不同频率成分的合成信号。
我们可以使用NumPy的arange函数生成一组时间点，然后使用sin函数生成不同频率的正弦波信号。
接下来，我们将实现DFT函数。
DFT将时域信号作为输入，并返回频域信号。
DFT的公式可以表示为：X(k) = Σ(x(n) * exp(-i*2πkn/N))其中，X(k)是频域信号的第k个频率点，x(n)是时域信号的第n个样本，N是信号的长度。
我们将使用循环计算DFT，但这种方法的计算复杂度为O(N^2)。
因此，我们将在实验过程中进行一些优化。
接下来，我们将实现FFT函数。
FFT函数将时域信号作为输入，并返回频域信号。
可以使用Cooley-Tukey的分治算法来快速计算FFT。
FFT的基本思想是将一个长度为N的信号分解为两个长度为N/2的子信号，然后逐步地将子信号分解为更小的子信号。
最后，将所有子信号重新组合以得到频域信号。
实验中，我们将使用递归的方式实现FFT算法。
首先，我们将信号分解为两个子信号，然后对每个子信号进行FFT计算。
最后，将两个子信号的FFT结果重新组合以得到频域信号。
FFT算法实现
C#实现FFT算法专业：通信与信息系统学号：姓名：2014年11月一、设计原理FFT 算法就是不断把长序列的DFT 分解成几个短序列的DFT ，利用上述周期性和对称性来减少DFT 的运算次数. FFT 算法基本上分为两大类：(1)时域抽取法FFT(Decimation In Time FFT,简称DIT —FFT). (2)频域抽取法FFT(Decimation In Frequency FFT,简称DIF —FFT). 时域抽取基2FFT 算法原理设输入序列长为N=2M(M 为正整数)，将该序列按时间顺序的奇偶分解为越来越短的子序列，称为基2按时间抽取的FFT 算法。
若序列长度不满足条件N=2M ，可以加零补长使其达到 N=2M 。
算法步骤：第一步分组和变量置换：12,,1,0)()12(12)()2(2)(1,,1,0)()(2110-==+⇒+=⇒==⇒=⇒=-==∑-=Nr r x r x r n n r x r x r n n n x N k W n x k X N n nk N 其中令奇数令偶数按奇偶分成两组：将，第二步代入DFT 定义式：12,,1,0)()()()()12()2()]([)]([)(212/0212/021)12(12/012/0221-=⋅+⋅=⋅++⋅=+=∑∑∑∑-=-=+-=-=N r W Wr x Wr x Wr x Wr x r x DFT r x DFT k X kNrk N N r N r rk Nk r NN r N r rk N其中第三步求子序列的DFT蝶形图如下：CABA ＋ BCA － BC频域抽取基2FFT 算法原理设输入序列长度为N=2M(M 为正整数)，将该序列的频域输出序列X(k)按其频域顺序的奇偶分解为越来越短的子序列，称基2按频率抽取FFT 算法。
若序列长度不满足条件N=2M ，可以加零补长使其达到 N=2M 。
算法步骤：第一步分组：12,,1,0)2(12,,1,0)()()(1,,0)()(10-=+-=⇔-==∑-=Nn N n x Nn n x n x k X N k W n x k X N n nk N 后半部分：前半部分：按前后分为两部分：奇偶分组，第二步代入DFT 定义式nkNk NN N n N n kNn N Nn nk NkNn N N n nk NN n nk NW W N n x n x W N n x Wn x W N n x Wn x Wn x k X ])2()([)2()(])2()([)()(21212)2(120)2(12010++=++=++==∑∑∑∑∑-=-=+-=+-=-=第三步变量置换16,4,2,0)]2()([)(12,1,0''2112/02-=⋅++=∴-===∑-=N k W N n x n x k X N k k k Wnk N N n k N N其中：，，令蝶形运算符号：二、设计方案using System;using System.Collections.Generic; using System.Linq; using System.Text; namespace MyFFT {class Program {static void Main(string[] args) ///程序入口{const int N=128; ///定义序列最大长度，可根据需要修改Complex[] W = new Complex[N]; ///定义旋转因子Complex[] signal = new Complex[N]; ///定义序列FFT fft = new FFT(); ///定义fft 类fft.SIZE=Convert.ToInt32(Console.ReadLine());for (int i = 0; i < fft.SIZE; i++) {W[i] = new Complex(); ///初始化旋转因子signal[i] = new Complex(); ///初始化序列}fft.ShowInputMessage(signal); ///人机接口fft.initW(W); ///计算旋转因子 Console.Write("正变换FFT 请按0，逆变换IFFT 请按1：");fft.Model = Convert.ToInt32(Console.ReadLine()); if (fft.Model == 0)fft.DoFFT(W, signal);else if (fft.Model == 1)fft.DoIFFT(W,signal);Console.WriteLine("输出例如下：");fft.ShowResult(signal); ///输出变换后的序列结果Console.ReadLine();}}}namespace Fengfeng.He{///定义复数类Complex///成员：Re Imclass Complex{public double Re; ///复数实部public double Im; ///复数虚部public Complex() ///构造函数{}public Complex(double Re, double Im) ///带参的构造函数重载{this.Re = Re;this.Im = Im;}///运算符重载+号public static Complex operator +(Complex c1, Complex c2){return new Complex(c1.Re + c2.Re, c1.Im + c2.Im);}///运算符重载-号public static Complex operator -(Complex c1, Complex c2){return new Complex(c1.Re - c2.Re, c1.Im - c2.Im);}///运算符重载*号public static Complex operator *(Complex c1, Complex c2){return new Complex(c1.Re * c2.Re - c1.Im * c2.Im, c1.Re * c2.Im + c1.Im * c2.Re);}///运算符重载/号public static Complex operator /(Complex c1, Complex c2){// return new Complex(((c1.Re*c2.Re+c1.Im*c2.Im)/(c2.Re*c2.Re+c2.Im*c2.Im)),(c1. Im*c2.Re-c1.Re*c2.Im)/c2.Re*c2.Re+c2.Im*c2.Im);Complex c = new Complex();c.Re = (c1.Re * c2.Re + c1.Im * c2.Im) / (c2.Re * c2.Re + c2.Im * c2.Im);c.Im = (c1.Im * c2.Re - c1.Re * c2.Im) / (c2.Re * c2.Re + c2.Im * c2.Im);return c;}};/// <summary>/// 定义FFT类/// 成员函数：/// 1、ReverseOder()倒序变换/// 2、DoFFT()FFT变换/// 3、DoIFFT()IFFT变换/// 4、W()旋转因子计算/// </summary>class FFT{public int SIZE = 0; ///序列长度public int Model; ///变换模式private double PI = Math.Atan(1) * 4;// private double PI = 3.14159;public void ReverseOder(Complex[] x){Complex temp;int i = 0, j = 0, k = 0;double t;for (i = 0; i < SIZE; i++){k = i; j = 0;t = Math.Log(SIZE, 2);while ((t--) > 0){j = j << 1;j |= (k & 1);k = k >> 1;}if (j > i){temp = x[i];x[i] = x[j];x[j] = temp;}}}/// <summary>/// 倒序函数/// W：输入旋转因子序列数组名/// </summary>public void initW(Complex[] W){int i;#if DEBUGfor (i = 0; i < SIZE; i++){Console.WriteLine("DEBUG->W{0}+j{1}", W[i].Re, W[i].Im);}#endiffor (i = 0; i < SIZE; i++){W[i].Re = Math.Cos(2 * PI * i / SIZE); //旋转因子实部展开W[i].Im = -1 * Math.Sin(2 * PI * i / SIZE); //旋转因子虚部展开}#if DEBUGfor (i = 0; i < SIZE; i++){Console.WriteLine("DEBUG->WW{0}+j{1}",W[i].Re, W[i].Im);}#endif}public void DoFFT(Complex[] W, Complex[] x){int i = 0, j = 0, k = 0, l = 0;Complex up, down, product;ReverseOder(x); ///对输入序列进行倒序#if DEBUGfor (i = 0; i < SIZE; i++){Console.WriteLine("DEBUG->{0}+j{1}", x[i].Re, x[i].Im);}#endiffor (i = 0; i < Math.Log(SIZE, 2); i++) ///外循环，变换级数循环{l = 1 << i;for (j = 0; j < SIZE; j += 2 * l) ///中间循环，旋转因子循环{for (k = 0; k < l; k++) ///内循环，序列循环{product = x[j + k + l] * W[SIZE * k / 2 / l];up = x[j + k] + product;down = x[j + k] - product;x[j + k] = up;x[j + k + l] = down;}}}#if DEBUGfor (i = 0; i < SIZE; i++){Console.WriteLine("DEBUG->fft:{0}+j{1}", x[i].Re, x[i].Im);}#endif}public void DoIFFT(Complex[] W, Complex[] x){int i = 0, j = 0, k = 0, l = SIZE;Complex up, down;for (i = 0; i < Convert.ToInt32((Math.Log(SIZE, 2))); i++) ///外循环，级数循环{l /= 2;for (j = 0; j < SIZE; j += 2 * l) ///中间循环，旋转因子循环{for (k = 0; k < l; k++) ///内循环，序例循环{up = x[j + k] + x[j + k + l];up.Re /= 2; up.Im /= 2;down = x[j + k] - x[j + k + l];down.Re /= 2; down.Im /= 2;down = down / W[SIZE * k / 2 / l];x[j + k] = up;x[j + k + l] = down;}}}ReverseOder(x); ///倒序还原}public void ShowResult(Complex []signal){for (int i = 0; i < SIZE; i++){Console.Write("{0}", signal[i].Re);if (signal[i].Im > 0.00001) ///定义精度Console.WriteLine(" + {0}i", signal[i].Im);else if (Math.Abs(signal[i].Im) > 0.00001)Console.WriteLine(" {0}i", signal[i].Im);elseConsole.WriteLine();}}public void ShowInputMessage(Complex []signal){///人机接口，命令提示行for (int i = 0; i < SIZE; i++){三、仿真结果。
FFT算法分析实验实验报告
FFT算法分析实验实验报告一、实验目的快速傅里叶变换（Fast Fourier Transform，FFT）是数字信号处理中一种非常重要的算法。
本次实验的目的在于深入理解 FFT 算法的基本原理、性能特点，并通过实际编程实现和实验数据分析，掌握 FFT 算法在频谱分析中的应用。
二、实验原理FFT 算法是离散傅里叶变换（Discrete Fourier Transform，DFT）的快速计算方法。
DFT 的定义为：对于长度为 N 的序列 x(n)，其 DFT 为X(k) ＝∑n=0 到 N-1 x(n) e^（j 2π k n ／ N) ，其中 j 为虚数单位。
FFT 算法基于分治法的思想，将 N 点 DFT 分解为多个较小规模的DFT，从而大大减少了计算量。
常见的 FFT 算法有基 2 算法、基 4 算法等。
三、实验环境本次实验使用的编程语言为 Python，主要依赖 numpy 库来实现 FFT 计算和相关的数据处理。
四、实验步骤1、生成测试信号首先，生成一个包含不同频率成分的正弦波叠加信号，例如100Hz、200Hz 和 300Hz 的正弦波。
设定采样频率为 1000Hz，采样时间为 1 秒，以获取足够的采样点进行分析。
2、进行 FFT 计算使用 numpy 库中的 fft 函数对生成的测试信号进行 FFT 变换。
3、频谱分析计算 FFT 结果的幅度谱和相位谱。
通过幅度谱确定信号中各个频率成分的强度。
4、误差分析与理论上的频率成分进行对比，计算误差。
五、实验结果与分析1、幅度谱分析观察到在 100Hz、200Hz 和 300Hz 附近出现明显的峰值，对应于生成信号中的频率成分。
峰值的大小反映了相应频率成分的强度。
2、相位谱分析相位谱显示了各个频率成分的相位信息。
3、误差分析计算得到的频率与理论值相比，存在一定的误差，但在可接受范围内。
误差主要来源于采样过程中的量化误差以及 FFT 算法本身的近似处理。
信号实验报告( 离散傅里叶变换及其快速算法及IIR数字滤波器的设计)
信号实验一离散傅里叶变换及其快速算法一、实验目的1、掌握计算序列的离散傅里叶变换（FFT）的方法；2、掌握实现时间抽取快速傅里叶变换（FFT）编程方法；3、加深对DFT与序列的傅里叶变换和Z变换之间的关系的理解；4、复习复数序列的运算方法。
二、程序设计框图1.码位倒置程序框图2.蝶形图运算程序框图三、实验程序实验程序的源代码如下:#include"math.h"#include"stdio.h"/*------------------------------------------------------------------------------------------子函数部分------------------------------------------------------------------------------------------*/ void swap(float *a,float *b)//交换变量子函数{float T;T=*a;*a=*b;*b=T;}void fft (float A [],float B [],unsigned M)//数组A为序列的实部, 数组B为序列的虚部{unsigned long N,I,J,K,L,LE,LE1,P,Q,R;float Wr,Wi,W1r,W1i,WTr,WTi,theta,Tr,Ti;N=1<<M;J=0;for(I=0;I<N-1;I++){if(J>I){swap(&A [I],&A [J]);swap(&B [I],&B [J]);}K=N>>1;while(K>=2&&J>=K){J-=K;K>>=1;}J+=K;}for(L=1;L<=M;L++){LE=1<<L;LE1=LE/2;Wr=1.0;Wi=0.0;theta=(-1)*3.1415926536/LE1;W1r=cos (theta);W1i=sin (theta);for(R=0;R<LE1;R++){for(P=R;P<N-1;P+=LE){Q=P+LE1;//基本蝶形图的复数运算Tr=Wr*A[Q]-Wi*B[Q];Ti=Wr*B[Q]+Wi*A[Q];A[Q]=A[P]-Tr;B[Q]=B[P]-Ti;A[P]+=Tr;B[P]+=Ti;}WTr=Wr;WTi=Wi;Wr=WTr*W1r-WTi*W1i;Wi=WTr*W1i+WTi*W1r;}}return;}/*------------------------------------------------------------------------------------------主函数部分------------------------------------------------------------------------------------------*/ void main(){float A[20],B[20];char t1,t2,file_name[20];int M,N,i,iiff;FILE *fp;/*************************************数据读取部分************************************/ printf("请输入文件名:");//输入数据文件名scanf("%s",file_name);printf("FFT变换还是IFFT变换?(FFT:1,IFFT:-1):");//输入变换方式, 1为FFT, -1为IFFTscanf("%d",&iiff);while(iiff!=1&&iiff!=-1)//检错: 检验上一步的输入是否有错, 有错则重新输入{printf("输入错误, 请重新输入！ ");printf("FFT or IFFT?(FFT:1,IFFT:-1):");scanf("%d",&iiff);}fp=fopen(file_name,"r");//打开文件并读入数据fscanf(fp,"%d",&M);N=pow(2,M);//计算序列总数for(i=0;i<N;i++)//读取文件中的数据{fscanf(fp,"%f%c%c%f",&A[i],&t1,&t2,&B[i]);if(iiff==-1)//根据FFT或IFFT修正BB[i]=B[i]*-1;if(t2!='j')//检错: 检验读取格式是否有错{printf("输入格式错误\n");break;}if(t1=='+')//判断虚部的正负号B[i]=B[i];else if(t1=='-')B[i]=-B[i];}/****************************************变换部分****************************************/ fft(A,B,M);//FFT变换/**************************************数据输出部分**************************************/ fp=fopen("fft_result.txt","w"); //输出结果if(iiff==-1)fprintf(fp,"IFFT变换的输出结果是: \n");elsefprintf(fp,"FFT变换的输出结果是: \n");for(i=0;i<N;i++){if(iiff==-1) //根据FFT或IFFT修正B{B[i]=B[i]*-1/N;A[i]=A[i]/N;}if(B[i]>=0)//修正虚部的输出格式fprintf(fp,"%f+j%f\n",A[i],B[i]);else if(B[i]<0)fprintf(fp,"%f-j%f\n",A[i],-B[i]);else if(B[i]==0)fprintf(fp,"%f\n",A[i]);}fclose(fp);}四、程序运行结果检验（1） 1.对序列进行FFT变换输入文件fft_input.txt:21+j02+j0-1+j04+j0控制台输入:请输入文件名: fft_input.txtFFT变换还是IFFT变换?(FFT:1,IFFT:-1): 1输出文件fft_result.txt:FFT变换的输出结果是:6.00000+j0.000002.00000+j2.00000-6.00000+j0.000002.00000+j-2.00000运行结果分析:程序运行输出结果与计算结果相同, 表示傅里叶正变换（FFT）成功。
- 1、下载文档前请自行甄别文档内容的完整性，平台不提供额外的编辑、内容补充、找答案等附加服务。
- 2、"仅部分预览"的文档,不可在线预览部分如存在完整性等问题,可反馈申请退款(可完整预览的文档不适用该条件!)。
- 3、如文档侵犯您的权益，请联系客服反馈,我们会尽快为您处理(人工客服工作时间：9:00-18:30)。
实验2 FFT 算法实现
2.1 实验目的
1、 加深对快速傅里叶变换的理解。
2、 掌握FFT 算法及其程序的编写。
3、 掌握算法性能评测的方法。
2.2 实验原理
一个连续信号)(t x a 的频谱可以用它的傅立叶变换表示为
dt e t x j X t j a a Ω-+∞
∞-⎰=
Ω)()( （2-1）
如果对该信号进行理想采样，可以得到采样序列
)()(nT x n x a = （2-2）
同样可以对该序列进行z 变换，其中T 为采样周期
∑+∞∞--=n z n x z X )()(
（2-3）
当ωj e
z =的时候，我们就得到了序列的傅立叶变换 ∑+∞∞-=n j j e n x e X ωω)()(
（2-4）
其中ω称为数字频率，它和模拟域频率的关系为
s f T /Ω=Ω=ω （2-5）
式中的s f 是采样频率。
上式说明数字频率是模拟频率对采样率s f 的归一化。
同模拟域的情况相似，数字频率代表了序列值变化的速率，而序列的傅立叶变换称为序列的频谱。
序列的傅立叶变换和对应的采样信号频谱具有下式的对应关系。
∑+∞∞--=)2(1)(T m j X T e X a j πωω （2-6）
即序列的频谱是采样信号频谱的周期延拓。
从式（2-6）可以看出，只要分析采样序列的频谱，就可以得到相应的连续信号的频谱。
注意：这里的信号必须是带限信号，采样也必须满
足Nyquist 定理。
在各种信号序列中，有限长序列在数字信号处理中占有很重要的地位。
无限长的序列也往往可以用有限长序列来逼近。
对于有限长的序列我们可以使用离散傅立叶变换（DFT ），这一变换可以很好地反应序列的频域特性，并且容易利用快速算法在计算机上实现当序列的长度是N 时，我们定义离散傅立叶变换为：
∑-===10)()]([)(N n kn N
W n x n x DFT k X （2-7） 其中N j N e W π
2-=，它的反变换定义为：
∑-=-==10)(1)]([)(N k kn N W
k X N k X IDFT n x （2-8）
根据式（2-3）和（2-7）令k N W z -=，则有
∑-====-10)]([)(|)(N n nk N W z n x DFT W n x z X k N
（2-9） 可以得到k N j k
N e W z z X k X π2|)()(===-,k N W -是z 平面单位圆上幅角为k N
πω2=的点，就是将单位圆进行N 等分以后第k 个点。
所以，X(k)是z 变换在单位圆上的等距采样，或者说是序列傅立叶变换的等距采样。
时域采样在满足Nyquist 定理时，就不会发生频谱混淆；同样地，在频率域进行采样的时候，只要采样间隔足够小，也不会发生时域序列的混淆。
DFT 是对序列傅立叶变换的等距采样，因此可以用于序列的频谱分析。
在运用DFT 进行频谱分析的时候可能有三种误差，分析如下：
（1）混淆现象
从式（2-6）中可以看出，序列的频谱是采样信号频谱的周期延拓，周期是2π/T ，因此当采样速率不满足Nyquist 定理，即采样频率T f s /1=小于两倍的信号（这里指的是实信号）频率时，经过采样就会发生频谱混淆。
这导致采样后的信号序列频谱不能真实地反映原信号的频谱。
所以，在利用DFT 分析连续信号频谱的时候，必须注意这一问题。
避免混淆现象的唯一方法是保证采样的速率足够高，使频谱交叠的现象不出现。
这就告诉我们，在确定信号的采样频率之前，需要对频谱的性质有所了解。
在一般的情况下，为了保证高于折叠频率的分量不会出现，在采样之前，先用低通模拟滤波器对信号进行滤波。
（2）泄漏现象
实际中的信号序列往往很长，甚至是无限长序列。
为了方便，我们往往用截短的序列来近似它们。
这样可以使用较短的DFT 来对信号进行频谱分析。
这种截短等价于给原信号序列乘以一个矩形窗函数。
而矩形窗函数的频谱不是有限带宽的，从而它和原信号的频谱进行卷积以后会扩展原信号的频谱。
值得一提的是，泄漏是不能和混淆完全分离开的，因为泄露导致频谱的扩展，从而造成混淆。
为了减小泄漏的影响，可以选择适当的窗函数使频谱的扩散减到最小。
（3）栅栏效应
因为DFT 是对单位圆上z 变换的均匀采样，所以它不可能将频谱视为一个连续函数。
这样就产生了栅栏效应，从某种角度来看，用DFT 来观看频谱就好像通过一个栅栏来观看一幅景象，只能在离散点上看到真实的频谱。
这样的话就会有一些频谱的峰点或谷点被“栅栏”挡住，不能被我们观察到。
减小栅栏效应的一个方法是在源序列的末端补一些零值，从而变动DFT 的点数。
这种方法的实质是认为地改变了对真实频谱采样的点数和位置，相当于搬动了“栅栏”的位置，从而使得原来被挡住的一些频谱的峰点或谷点显露出来。
注意，这时候每根谱线多对应的频率和原来的已经不相同了。
从上面的分析过程可以看出，DFT 可以用于信号的频谱分析，但必须注意可能产生的误差，在应用过程中要尽可能减小和消除这些误差的影响。
快速傅立叶变换FFT 并不是与DFT 不相同的另一种变换，而是为了减少DFT 运算次数的一种快速算法。
它是对变换式（2-7）进行一次次的分解，使其成为若干小点数DFT 的组合，从而减小运算量。
常用的FFT 是以2为基数，其长度M
N 2 。
它的运算效率高，程序比较简单，使用也十分地方便。
当需要进行变换的序列的长度不是2的整数次方的时候，为了使用以2为基的FFT ，可以用末尾补零的方法，使其长度延长至2的整数次方。
IFFT 一般可以通过FFT 程序来完成，比较式（2-7）和（2-8），只要对X(k)取共轭，进行FFT 运算，然后再取共轭，并乘以因子1/N ，就可以完成IFFT 。
2.3 实验内容
1、 编制自己的FFT 算法。
2、 选取实验1中的典型信号序列验证算法的有效性。
3、 对所编制FFT 算法进行性能评估。
2.4 实验报告要求
1、 总结自己实现FFT 算法时候采用了哪些方法减小了运算量。
2、 给出自己的FFT 算法与实验1中自己的DFT 算法的性能比较结果。
3、 给出自己的FFT 算法与Matlab 中FFT 算法的性能比较结果。
4、 总结实验中根据实验现象得到的其他个人结论。