离散傅立叶变换DFT实验报告

1、实验一离散时间系统的时域分析一、实验目的.运用 MATLAB仿真一些简单的离散时间系统，并研究它们的时域特性。 . 运用 MATLAB中的卷积运算计算系统的输出序列，加深对离散系统的差分方程、冲激响应和卷积分析方法的理解。二、实验原理离散时间系统其输入、输出关系可用以下差分方程描述：NMdk ynkpk x nkk0k0当输入信号为冲激信号时，系统的输出记为系统单位冲激响应 nh n ，则系统响应为如下的卷积计算式：ynxn hnxmhn mm当 hn 是有限长度的（n：0 ， M）时，称系统为FIR 系统；反之，称系统为IIR 系统。在 MATLAB中，可以用函数y=Filter(p,d,x

2、)求解差分方程，也可以用函数y=Conv(x,h)计算卷积。例 1clf;n=0:40;a=1;b=2;x1= *n;x2=sin(2*pi*n);x=a*x1+b*x2;num=1, ,3;den=2 -3 ;ic=0 0;%设置零初始条件y1=filter(num,den,x1,ic); %计算输入为x1(n) 时的输出y1(n)y2=filter(num,den,x2,ic); %计算输入为x2(n) 时的输出y2(n)y=filter(num,den,x,ic); %计算输入为x (n)时的输出y(n)yt= a*y1+b*y2;%画出输出信号subplot(2,1,1)stem(n,

3、y);ylabel(振幅 ) ；title(加权输入a*x1+b*x2 的输出 ) ；subplot(2,1,2)stem(n,yt);ylabel(振幅 ) ；title(加权输出a*y1+b*y2 ) ；( 一 ) 、线性和非线性系统对线性离散时间系统，若 y1 (n) 和 y2 (n) 分别是输入序列x1 (n) 和 x2 (n) 的响应， 则输入x(n)ax1 (n)bx2 (n) 的输出响应为y( n)ay1 (n)by2 (n) ，即符合叠加性， 其中对任意常量 a 和 b 以及任意输入x1 ( n) 和 x2 (n) 都成立，否则为非线性系统。( 二 ) 、时不变系统和时变系统对

4、离散时不变系统，若y1 ( n) 是 x1 ( n) 的响应，则输入x(n)=x1(n-n0)的输出响应为y(n)=y1(n-n0),式中 n0 是任意整数。该输入输出关系，对任意输入序列及其相应的输出成立，若对至少一个输入序列及其相应的输出序列不成立，则系统称之为时变的。( 三 ) 、线性卷积假设待卷积的两个序列为有限长序列，卷积运算符在MATLAB中可命令 conv 实现。例如，可以把系统的冲激响应与给定的有限长输入序列进行卷积，得到有限长冲激响应系统的输出序列。下面的MATLAB程序实现了该方法。例 2clf;h=3 2 1 -210-4 0 3;%冲激x=1 -2 3 -4321;%输

5、入序列y=conv(h,x);n=0:14;stem(n,y);xlabel(时间序号n );ylabel(振幅 ) ；title(用卷积得到的输出) ； grid;三、实验内容与步骤.假定一因果系统为y(n)(n-1)+(n-2)=(n)+(n-1)+(n-2)用 MATLAB程序仿真该系统，输入三个不同的输入序列：x1 (n)cos(20.1n) ， x2 (n)cos(20.4n) ， x2x1 ( n)3x2 (n)计算并并显示相应的输出y1 ( n) ，y2 ( n) 和 y( n) 。.用 MATLAB程序仿真步骤1 给出的系统，对两个不同的输入序列x(n) 和 x(n-10)，计

6、算并显示相应的输出序列y3(n) 和 y4(n) 。3用 MATLAB程序仿真计算下列两个有限长序列的卷积和并显示图形。x1 (n)(n)3 ( n1)2 (n2)x2 (n)u(n)u(n3)四、实验仪器设备计算机， MATLAB软件五、实验要求给出理论计算结果和程序计算结果并讨论。六、实验结果实验 1：clf;n=0:40;a=2;b=-3;x1= cos(2*pi*n);x2=cos(2*pi*n);x=a*x1+b*x2;den=1, ,;num = ;%分子系数ic=0 0;%设置零初始条件y1=filter(num,den,x1,ic); %计算输入为x1(n) 时的输出y1(n)

7、y2=filter(num,den,x2,ic); %计算输入为x2(n) 时的输出y2(n)yn=filter(num,den,x,ic); %计算输入为x (n)时的输出y(n)% 画出输出信号subplot(2,2,1)stem(n,y1);ylabel('振幅 ');title('y1输出 ');subplot(2,2,2)stem(n,y2);ylabel('振幅 ');title('y2输出 ');subplot(2,2,3)stem(n,yn);ylabel('振幅 ');title('yn输

8、出 ');实验 2：clf;n=0:40;n1=0:50;a=2;b=-3;x1= cos(2*pi*n);x2=cos(2*pi*n);x3=a*x1+b*x2;x4=zeros(1,10), x3;den=1, ,;num= ;ic=0 0;%设置零初始条件y3=filter(num,den,x3,ic);y4=filter(num,den,x4,ic);%计算输入为x (n)时的输出y(n)%画出输出信号subplot(2,1,1)stem(n,y3);ylabel('振幅 ');title('yn输出 ');subplot(2,1,2)stem(

9、n1,y4);ylabel('振幅 ');title('y1输出 ');实验 3：clf;x=132;%冲激u=111; %输入序列y=conv(u,x);n=0:4;stem(n,y);xlabel('时间序号 n');ylabel('振幅 ');title('用卷积得到的输出 ');grid;实验二（ 1） 离散时间信号的DTFT一、实验目的.运用 MATLAB理解 Z 变换及其绘制H(z) 的零极点图 。.运用 MATLAB计算逆 Z 变换。二、实验原理（一）、 MATLAB在 ZT 中的应用 。线性时不变离

10、散时间系统的冲激响应h(n) 的 z 变换是其系统函数 H(z) ， 在 MATLAB中可以利用性质求解Z 变换，例如可以利用线性卷积求的Z 变换。若 H(z) 的收敛域包含单位圆，即系统为稳定系统，即系统在单位圆上z e j处计算的是系统的频率响应。（二）、逆 Z 变换Z 变换对于分析和表示离散线性时不变系统具有重要作用。但是在MATLAB中不能直接计算 Z 变换，但是对于一些序列可以进行逆Z 变换。已知序列的 Z 变换及其收敛域，求序列称为逆Z 变换。 序列的 Z 变换及共逆 Z 变换表示如下：X ( z)x( n) z n , Rxz Rxnx(n)1X ( z) zn 1dz,c( R

11、x, Rx )2jc通常，直接计算逆Z 变换的方法有三种：围线积分法、长除法和部分分式展开法。在实际中， 直接计算围线积分比较困难，往往不直接计算围线积分。 由于序列的 Z 变换常为有理函数，因此采用部分分式展开法比较切合实际，它是将留数定律和常用序列的Z 变换相结合的一种方法。设 x(n) 的 Z 变换 X(z) 是有理函数，分母多项式是N 阶，分子多项式是M阶，将 X(z) 展成一些简单的常用的部分分式之和，通过常用序列的Z 变换求得各部分的逆变换， 再相加即得到原序列 x(n) 。在 MATLAB中提供了函数residuez来实现上述过程，调用格式如下：R ， P， K= residue

12、z（ B，A）其中 B、A 分别是有理函数分子多项式的系数和分母多项式的系数，输出 R 是留数列向量， P是极点列向量。如果分子多项式的阶数大于分母多项式的阶数，则K 返回为常数项的系数。三、实验内容与步骤选做一个实验：1、 . 运行下面程序并显示它，验证离散时间傅立叶变换DTFT的时移性。已知两个线性时不变的因果系统，系统函数分别为H 1(z) 1 z N ， H 2( z)1 z N1aN z N分别令 N=8， a=，计算并图示这两个系统的零、极点图及幅频特性。程序：2、运行下面程序并显示它，验证离散时间傅立叶变换DTFT的频移性。四、实验仪器设备计算机， MATLAB软件五、实验注意事

13、项课前预先阅读并理解实验程序；六、实验结果clearnum1=1 0 0 0 0 0 0 0 -1;%分子系数高阶到低阶den1=1 0 0 0 0 0 0 0 0;subplot(2,2,1)zplane(num1, den1)grid;title('H1零极点分布图') ;H,w=freqz(num1,den1,200,'whole'); %中 B 和 A 分别为离散系统的系统函数分子、分母多项式的系数向量，HF=abs(H);%返回量 H则包含了离散系统频响在0pi 范围内 N个频率等分点的值（其中N 为正整数）subplot(2,2,2);%w则包含了范

14、围内N 个频率等分点。plot(w,HF)title('H1幅频响应特性曲线');a=;A=a8;num2=1 0 0 0 0 0 0 0 -1;%分子系数高阶到低阶den2=1 0 0 0 0 0 0 0 A;subplot(2,2,3)zplane(num2, den2);grid;title('H2零极点分布图') ;H,w=freqz(num2,den2,200,'whole'); %中 B 和 A 分别为离散系统的系统函数分子、分母多项式的系数向量，HF=abs(H);%返回量 H则包含了离散系统频响在0pi 范围内 N个频率等分点的值

15、（其中N 为正整数）subplot(2,2,4);%w则包含了范围内N 个频率等分点。plot(w,HF)title('H2幅频响应特性曲线');实验二（ 2） 离散傅立叶变换DFT一、实验目的 . 运用 MATLAB计算有限长序列的 DFT和 IDFT。. 运用 MATLAB验证离散傅立叶变换的性质。3 . 运用 MATLAB计算有限长序列的圆周卷积。二、实验原理（一）、离散傅立叶变换DFT的定义一个有限长度的序列x(n) （ 0 n<N-1 ） ,它的 DFT X(k)可以通过在 轴（ 02 ）上对 X (e j) 均匀采样得到X (k ) X (ej)x(n)e j

16、 2kn / N2 k / N0k N 1n可以看到 X (k) 也是频域上的有限长序列，长度为 N。序列 X (k ) 称为序列 x(n) 的 N 点 DFT。N称为 DFT变换区间长度。通常表示WNe j 2/ N可将定义式表示为X (k)x(n)W knkN1n0X(k) 的离散傅里叶逆变换(IDFT) 为x(n)1X ( k)W kn0nN1N n（二）、 DFT的性质1圆周移位定义序列x(n) 的 m单位的圆周移位y(n) 为：m)RN (n) x(n m) N RN (n)y(n) x(n( x(n m) N 即对 x(n) 以 N 为周期进行周期延拓的序列RN (n) 表示对x

17、(n) 的 m点移位，此延拓移位后再取主值序列)1 圆周卷积设 x1 (n)DFTX1( k)0 k N 1Nx2 (n)DFTX 2 ( k)0 k N 1N则 x1 ( n) x2 (n)DFTX 1( k) X2(k ) 0 k N 1N这里x1 ( n)x2 (n) 表示 x1 (n) 与x2 (n) 的 N 点循环卷积。N1x1 (n)x2 (n)x2 (m) x1 (nm) N RN (n), n0,1, N1m02 共轭对称性x(n)xep ( n) xop (n),0n N1xep (n)1 x(n)x* (Nn)n N 12,0xop (n)1 x(n)x* (Nn)2x(n

18、)DFTX (k )NDFTxep (n)N1* X (k)X (k)Re X (k )X r (k)实际应用中，利用上述对称性质可以减少DFT的运算量，提高运算效率。三、实验内容与步骤：（2，3 选做一个）.构造离散傅立叶正、反变换函数的MATLAB程序，其中dft(xn,N)为离散傅立叶正变换，idft(xn,N)为离散傅立叶反变换。2、如果 x(n)sin(n/ 8)sin( n/ 4) 是一个 N=16的有限长序列， 利用离散傅立叶变换函数求其 16 点 DFT。3、如果 x(n)cos(0.82 n)2sin(0.43n) 是一个 0n100 的有限长序列， 绘制 x(n) 及其离散

19、傅立叶变换X（ K）的幅度、相位图。四、实验仪器设备计算机， MATLAB软件五、实验注意事项课前预先阅读并理解实验程序；六、实验结果Dft: 程序function xk=dft(xn,N)%dftn=0:1:N-1;k=n;WN=exp(-i*2*pi/N); %旋转因子nk=n'*k;WNnk=WN.nk;xk=xn*WNnk;endidft：程序function xn=idft(xk,N)%idftn=0:1:N-1;k=n;WN=exp(-j*2*pi/N);nk=n'*k;%矩阵的转制 *KWNnk=WN.(-nk);xn=xk*WNnk/N;end实验程序：选做 2

20、k=16; % 序列长N=16;%dft 点数n1=0:1:15;xn1=sin(pi/8*n1/k)+sin(pi/4*n1/k); %抽样信号xk1=dft(xn1,N);subplot(1,2,1);stem(n1,xn1);xlabel('t/T');ylabel('x(n)');subplot(1,2,2);stem(n1,xk1);grid;xlabel('k');ylabel('x(k)');实验二（ 3） 快速傅立叶变换FFT 及其应用一、实验目的 . 利用 MATLAB的快速傅立叶变换来计算信号的离散傅立叶变换。

21、. 利用 MATLAB程序，理解进一步离散傅立叶变换的物理意义。3. 利用 MATLAB程序，理解快速卷积算法。二、实验原理在 MATLAB中，使用函数 fft可以很容易地计算有限长序列x(n) 的离散傅立叶变换Xk 。此函数有两种形式， fft(x)计算序列 x(n)的离散傅立叶变换X(k) ，这里 X(k) 的长度与x(n)的长度相等。 fft(x， L) 计算序列x(n)的 L 点离散傅立叶变换，其中L N。若 L>N，在计算离散傅立叶变换之前，对 x(n) 尾部的 L-N 个值进行补零。 同样，离散傅立叶变换序列X(k)的离散傅立叶逆变换x(n) 用函数 ifft计算，它也有两种

22、形式。（一）、基本序列的离散傅立叶变换计算N 点离散傅立叶变换的一种物理解释就是，Xk是 x(n) 以 N 为周期的周期延拓序列的离 散 傅 立 叶 级 数 系 数 X (k) 的 主 值 区 间 序 列 ， 即 X (k )X (k ) RN ( k) 。 例 如 序 列cos(n)RN (n) , 当 N=16时 ， cos( n)RN (n) 正 好 是 cos( n) 的 一 个 周 期 ， 所 以888cos(n)RN (n) 的 周 期 延 拓 序 列 就 是 这 种 单 一 频 率 的 正 弦 序 列 。 而 当 N=8 时 ，8cos(n)RN (n) 正好是 cos(n) 的

23、半个周期， cos(n) RN (n) 的周期延拓就不再是单一频888率的正弦序列， 而是含有丰富的谐波成分， 其离散傅立叶级数的系数与N=16时的差别很大，因此对信号进行谱分析时，一定要截取整个周期，否则得到错误的频谱。（二）、验证 N 点 DFT的物理意义假如 x(n) 非周期、有限长，则傅立叶变换存在，那么对X (e j) 在 N 个等间隔频率k =2 k/N, k=0,1, ,N-1取样 , 则可得 X(k) 。X (k )X ( )2 k / Nx( n)e j 2 kn / N0 k N 1n序列 x(n) 的 N点 DFT的物理意义是对X( ) 在 0 ， 2 上进行 N点的等间

24、隔采样。（三）、利用 FFT 计算序列的线性卷积直接计算线性卷积计算量大，并且计算机无法判断 y(n)的长度，需要计算多少的 y(n)值，若输入为无限长，就更无法计算，其运算量随长度成级数增长。由于可以利用FFT 对DFT进行有效的计算，我们希望能够利用DFT来计算线性卷积。设 x(n) 和 h(n)是长度分别为 M和 N 的有限长序列，令 L=M+N-1，定义两个长度L 的有限长序列：x(n),0nM1x'( n)MnL10,h(n),0nN1h'( n)0,NnL1通过对 x(n)和 h(n) 补充零样本值得到上面两个序列。那么：yl (n) x( n) h(n) yc (

25、n) x' ( n) h' (n)(3.4.10)上面的过程如下图所示：计算线性卷积也可以直接调用函数con 来计算，因为MATLAB中的计时比较粗糙，所以只有M和 N较大的时候，才能比较两种方法的执行时间快慢。三、实验内容与步骤（选做一个）j n . 对复正弦序列 x(n) e 8 RN (n) ，利用 MATLAB程序求当 N=16 和 N=8 时的离散傅立叶变换，并显示其图形。 . 已知 x( n)1eR4 (n) , X ( )e1j 4,绘制相应的幅频和相频曲线，并计算N=8 和jN=16 时的 DFT。四、实验仪器设备计算机， MATLAB软件五、实验注意事项课前预

26、先阅读并理解实验程序；六、实验结果k1=16; % 序列长N1=16;%dft 点数n1=0:1:15;xn1=exp(j*pi/8*n1/k1); %抽样信号xk1=dft(xn1,N1);subplot(2,2,1);stem(n1,xn1);xlabel('t/T');ylabel('x(n)');subplot(2,2,2);stem(n1,xk1);grid;xlabel('k');ylabel('x(k)');k2=8; % 序列长N2=8;%dftn2=0:1:7;点数xn2=exp(j*pi/8*n2/k2); %

27、xk2=dft(xn2,N2);subplot(2,2,3);stem(n2,xn2);xlabel('t/T');ylabel('x(n)');subplot(2,2,4);stem(n2,xk2);grid;xlabel('k');ylabel('x(k)');抽样信号实验三基于MATLAB的IIR数字滤波器设计一、实验目的.进一步熟悉IIR数字滤波器的理论知识。.熟悉与IIR数字滤波器设计有关的MATLAB函数。3 .学会通过MATLAB，利用脉冲响应不变法和双线性变换法设计IIR数字滤波器，加深对数字滤波器的常用指标和设计

28、过程的理解。二、实验原理（一）、低通滤波器的常用指标：1PG(e j )1P , forPG(ej)S , forS通带边缘频率：p ，阻带边缘频率：s，通带起伏：p ，通带峰值起伏：p20 log 10 (1p ) dB ，阻带起伏： s最小阻带衰减：S20 log10 ( s ) dB 。（二）、 IIR 数字滤波器设计目前，设计 IIR 数字滤波器的通用方法是先设计相应的低通滤波器，然后再通过双线性变换法和频率变换得到所需要的数字滤波器。模拟滤波器从功能上分有低通、高通、带通及带阻四种，从类型上分有巴特沃斯滤波器、切比雪夫滤波器、椭圆滤波器以及贝塞尔滤波器等。1、利用模拟滤波器设计IIR

29、 数字低通滤波器的步骤。(1) 确定数字低通滤波器的技术指标：通带截止频率p、通带衰减 p、阻带截止频率s、阻带衰减 s。(2) 将数字低通滤波器的技术指标转换成模拟低通滤波器的技术指标。脉冲响应不变法：T2 tan( 1)双线性变换法：T2(3) 按照模拟低通滤波器的技术指标设计模拟低通滤波器。(4) 将模拟滤波器Ha(s) ，从 s 平面转换到z 平面，得到数字低通滤波器系统函数H(z) 。2、下面给出与IIR 数字滤波器设计有关的MATLAB文件。（ 1）用来确定数字低通或模拟低通滤波器的阶次，其调用格式分别是a. N,Wn=buttord(Wp,Ws,Rp,Rs)b. N,Wn=but

30、tord(Wp,Ws,Rp,Rs, s )格式 a 对应数字滤波器，式中 Wp,Ws分别是通带和阻带的截止频率，实际上它们是归一化频率，其值在0-1 之间， 1 对应（即对的归一化）。Rp,Rs 分别是通带和阻带衰减，单位为 dB。 N是求出的相应低通滤波器的阶次，Wn是求出的3dB 频率。格式 b 对应模拟滤波器，式中各个变量的含义和格式a 相同，但Wp,Ws及 Wn是模拟角频率，单位为rad/s 。（ 2）用来设计模拟低通原型（归一化）滤波器Ha(p) ，其调用的格式为z , p, k=buttap(N)N是欲设计的低通原型（归一化） 滤波器的阶次， z, p 和 k 分别是设计出Ha(p

31、) 的极点、零点及增益。（ 3）将模拟低通原型 （归一化） 滤波器 Ha(p) 转换为实际的低通滤波器Ha(s) 。（去归一化） ，其调用格式为：B，A=lp2lp(b,a,Wn)b,a 分别是模拟低通原型滤波器Ha(p) 的分子、 分母多项式的系数向量，其中 B，A 是去归一化后Ha(s)的分子、分母多项式的系数向量, Wn 为截止频率。（ 4）实现双线性变换，即由模拟滤波器Ha(s) 得到数字滤波器H(z) 。其调用格式是：Bz,Az=bilinear(B,A,Fs)B， A 是去归一化后Ha(s)的分子、分母多项式的系数向量，Bz，Az 是 H (z)的分子、分母多项式的系数向量, Fs

32、是抽样频率。（ 4）由脉冲响应不变法将模拟滤波器Ha(s) 转换为数字滤波器H(z) 。其调用格式是：Bz,Az= impinvar(B,A,Fs)B， A 是去归一化后Ha(s)的分子、分母多项式的系数向量，Bz，Az 是 H (z)的分子、分母多项式的系数向量, Fs是抽样频率。(5)用来直接设计巴特沃斯数字滤波器（双线性变换法） ，实际上它把，及等文件都包含进去，从而使设计过程更简捷，其调用格式为：a.B， A=butter(N, Wn)b.B， A=butter(N，Wn， s)格式 a 是设计低通数字滤波器，格式b 是设计低通模拟滤波器。B，A 是 H (z)的分子、分母多项式的系数

33、向量，Wn是截止频率。三、实验内容与步骤以下选做一个.设计 MATLAB程序，采用脉冲响应不变法设计一个巴特沃斯低通数字滤波器，其通带上限临界频率为400Hz，阻带临界频率为600Hz，抽样频率是1000Hz，在通带内的最大衰减为 ,阻带内的最小衰减为60dB，并绘出幅频特性曲线。.设计 MATLAB程序，采用双线性变换法设计一个巴特沃斯低通数字滤波器，要求在通带0 ， 内衰减不大于3dB, 在阻带 , 内衰减不小于40dB，并绘出幅频特性曲线。四、实验仪器设备计算机， MATLAB软件五、实验要求根据要求独立编程设计，并根据程序运行结果写出滤波器的系统函数六、实验结果选做 1：fp=400;

34、 % 通带上限临界频率fs=600; % 阻带临界频率Rp=; %通带允许的最大衰减Rs=60; % 阻带允许的最小衰减Fs=1000; % 采样频率Wp=2*pi*fp; %通带截止平率Ws=2*pi*fs; %阻带截止平率%Nn=256;n=(0:100-1);%采样点数N,Wn=buttord(Wp,Ws,Rp,Rs,'s'); %用于计算阶数和截止平率b,a=butter(N,Wn,'s'); %计算分子向量b，分母向量aw=linspace(1,400,100)*2*pi;%起始值，终止值，元素个数H=freqs(b,a,w); %在 0,2 上进行采

35、样，采样频率点由矢量w 指定figure(1);plot(w/(2*pi),20*log10(abs(H);title('巴特沃斯模拟滤波器幅频特性');xlabel('频率 /Hz');ylabel('幅度 /db');%bz,az=impinvar(b,a,Fs); %caiyong冲击响应不变法转换为数字滤波器实验四基于 MATLAB的 FIR 数字滤波器设计一、实验目的.进一步熟悉FIR 数字滤波器的理论知识。.熟悉与 FIR 数字滤波器设计有关的MATLAB函数。3. 学会通过 MATLAB，利用窗函数法设计 FIR 数字滤波器 。二、

36、实验原理设计FIR 滤波器实际上是要在满足线性相位的条件下，实现幅度响应的逼近。而一个FIR 滤波器若是符合线性相位，则必须满足一定的条件，即：一个 FIR 滤波器若是线性相位的，则其单位冲激响应必然满足h( n)h(N1n)n=0,1, ,N-1h(n) 是关于 (N-1)/2对称（奇对称或偶对称）即,(1) h(n) 是偶对称序列1N21 n,0 n N 1h n h N(1) h(n)是奇对称（反对称）序列N12h nh N 1n设滤波器要求的理想频率响应为Hd(ejw) ,那么 FIR 滤波器的设计问题在于寻找N 1一系统函数 H ( z)h(n) z n，使其频率响应H (ejw ) H ( z) |z ejw 逼近 Hd(ejw) 。若要n 0求 FIR 滤波器具有线性相位特性，则 h(n) 必须满足上节所述的对称条件。逼近的方法有三种：窗口设计法（时域逼近） ；频率采样法 Frequency-sampling （频域逼近）；最优化设计 Optimum Equiripple（等波纹逼近）