快速傅里叶变化算法及应用

快速傅里叶变化算法及应用

课程数字信号处理实验

系（院）物理与电子科学系

专业

班级

学生姓名

学号

一、 离散傅里叶变换的性质

离散傅里叶变换（DFT ），是连续傅里叶变换在时域和频域上都离散的形式，将时域信号的采样变换为在离散时间傅里叶变换（DTFT ）频域的采样。在形式上，变换两端（时域和频域上）的序列是有限长的，而实际上这两组序列都应当被认为是离散周期信号的主值序列。即使对有限长的离散信号作DFT ，也应当将其看作经过周期延拓成为周期信号再作变换。在实际应用中通常采用快速傅里叶变换以高效计算DFT 。

离散傅里叶变换具有下列性质： 1) 线性

{}{}{}][DFT ][DFT ][][DFT 2121k x b k x a k bx k ax +=+

需将较短序列补零后，再按长序列的点数做DFT 2) 循环位移(Circular shift of a sequence)

][])[(][k R n k x k y N N += 3) 对称性 (symmetry)

周期共轭对称(Periodic conjugate symmetry)定义为

][*][])[(*][k N x k R k x k x N N -=-=

周期共轭反对称(Periodic conjugate antisymmetry)定义为

][*][])[(*][k N x k R k x k x N N --=--=

4) 循环卷积定理

][])[(])[(][*][210121k R n k x n x k x k x N N N n N ⎭

⎬⎫

⎩⎨⎧-=∑-=

二、快速傅里叶变换的算法原理

1. 将长序列DFT 分解为短序列的DFT

2. 利用旋转因子的周期性、对称性、可约性。将时域序列逐次分解为一组子序列，利用旋转因子的特性，由子序列的DFT 来实现整个序列的DFT 。

其中：快速傅里叶变换分为两种，分为基2时间抽取算法和基2频率抽取算法 1) 基2时间抽取(Decimation in time)FFT 算法

⎩⎨

⎧+=]

12[]

2[][r x r x k x 其中:r=0,1,2…

12

-N

2) 基2频率抽取(Decimation in frequency)FFT 算法

⎩⎨

⎧+=]

12[]

2[][m X m X m X

三、FFT 的特点和规律

1) FFT 算法特点：(r

N 2=)

（1）共需r 次迭代；

（2）第)1(r L L ≤≤次迭代对偶结点的偶距为L

L r L L N K 2/2==--，因此一组结点覆盖

的序号个数是

12)(2-=

-L L L N K K 。

（3）第)1(r L L ≤≤次迭代结点的组数为[]1

2)(2/-=-L L L K K N 。

（4）

L

P

N W 可以预先计算好，而且L P 的变化范围是

12~

0-N 。

2) FFT 算法流程：(r

N 2=)

（1）初始化：

1

0),()(0-≤≤←N n n x n x ；

（2）第)1(r L L ≤≤次迭代：

（a ）下标控制变量初始化0=L K ； （b ）“结点对”的个数初始化0=num ; （c ）

DO N

num WHILE L )2(<

四、蝶形运算的MATALAB 的程序设计

下面以16点为例：

disp('请输入一个16点序列'); for ii=1:16 %自由输入序列 x(ii) = input(['x(',num2str(ii),')=']); end

%整体运用原位计算

m=nextpow2(x);N=2^m; % 求x 的长度对应的2的最低幂次m if length(x)

x=[x,zeros(1,N-length(x))]; % 若x的长度不是2的幂，补零到2的整数幂

end

nxd=bin2dec(fliplr(dec2bin([1:N]-1,m)))+1; % 求1:2^m数列序号的倒序

y=x(nxd); % 将x倒序排列作为y的初始值

for mm=1:m % 将DFT作m次基2分解,从左到右，对每次分解作DFT运算,共做m 级蝶形运算，每一级都有2^(mm-1)个蝶形结

Nz=2^mm;u=1; % 旋转因子u初始化为WN^0=1

WN=exp(-i*2*pi/Nz); % 本次分解的基本DFT因子WN=exp(-i*2*pi/Nz)

for j=1:Nz/2 % 本次跨越间隔内的各次蝶形运算,在进行第mm级运算时需要2^(mm-1)个蝶形

for k=j:Nz:N % 本次蝶形运算的跨越间隔为Nz=2^mm

kp=k+Nz/2; % 蝶形运算的两个因子对应单元下标的关系

t=y(kp)*u; % 蝶形运算的乘积项

y(kp)=y(k)-t; % 蝶形运算

y(k)=y(k)+t; % 蝶形运算

end

u=u*WN; % 修改旋转因子,多乘一个基本DFT因子WN

end

end

y

y1=fft(x)

五、快速傅里叶变换MATALAB程序设计

function samples(f,Fs,M)

N=2^M; % fft点数=取样总点数

Ts=1/Fs; % 取样时间间隔

T=N*T s; % 取样总时间=取样总点数*取样时间间隔

n=0:N-1;

t=n*Ts;

Xn=cos(2*f*pi*t);