实验四报告.docx
《实验四报告.docx》由会员分享，可在线阅读，更多相关《实验四报告.docx（21页珍藏版）》请在冰豆网上搜索。
实验四报告
实验四离散时间系统的特性和离散傅里叶变换
任务一IIR系统的特性
某线性系统用差分方程表示为
（4.1.1）
1、求出系统函数，编程调用函数zplane画出系统函数的零极图；
对差分方程两边进行Z变换得
（4.1.2）
（4.1.3）
2、调用函数freqz，画出此系统的频率响应的幅度和相位。
编程原理、思路和公式
画零极点图用函数zplane画，调用函数freqz，出此系统的频率响应的幅度和相位。
函数的调用参考课本例题。
程序脚本，并注释
clc;
clear
b=[1,1];%系统函数分子系数向量
a=[1,-0.9,0.81];%系统函数分母系数向量
rp=roots（a）;%求极点
rz=roots（b）;%求零点
figure
（1）%图一
zplane（b,a）;%画出零极点图
[H,w]=freqz（b,a,1024,'whole'）;%计算频率响应
magX=abs（H）;angX=angle（H）;%计算幅度和相位
figure
（2）%图二
subplot（2,1,1）;
plot（w,magX）,gridon;%画幅度谱
xlabel（'\omega'）;%横轴
ylabel（'|X（e^j^\omega）|'）;%纵轴为幅度谱
title（'幅度谱'）;%给图片取名
subplot（2,1,2）;
plot（w,angX）,gridon;%画相位谱
xlabel（'\omega'）;%横轴
ylabel（'argX（e^j^\omega）'）;%纵轴
title（'相位谱'）;%给图片取名
仿真结果、图形
图4.1.1系统函数的零极图
图4.1.2系统的频率响应的幅度和相位
3、能否用编写的DTFT子函数无误差地计算此系统的频率响应特性？
不能，该系统的DTFT是连续的，但是，如果调用编写DTFT子涵是，DTFT是离散后的，所以无法无误差地计算此系统的频率响应特性
遇到的问题、解决方法及收获
1、学会了根据差分方程画零极点图，求系统频率响应
2、学会了无误差地系统的频率响应特性的方法
任务二FIR系统的特性
某系统的单位冲激响应为
（升余弦系统）（4.2.1）
1、
，画出系统的单位冲激响应图形；
编程原理、思路和公式
1、定义一个矩形序列
2、定义函数
，然后点乘
程序脚本，并注释
clc;
clear
n1=-40;n2=40;%定义观察起止点
n=n1:
n2;%观察范围n;
N=32;%举行矩阵的长度
L=length（n）;%定义n的长度
R=zeros（1,L）;%赋初值
R（-n1+1:
-n1+N）=1;%定义举行矩阵R
h=（0.5*（1-cos（2*pi*n/（N-1））））.*R;%定义单位冲击函数
stem（n,h,'.'）;%画图
xlabel（'n'）;%横轴
ylabel（'h（n）'）;%纵轴
title（'单位冲击响应h（n）=0.5（1-cos2n\pi/N-1）R_N（n）'）;%给图片取名
仿真结果、图形
图4.2.1系统的单位冲激响应
，N=32图形
2、利用编写的DTFT子函数计算此系统的频率响应特性，画出幅频特性图形、增益（dB）图形、相频特性，
。
编程原理、思路和公式
调用实验二写好的DTFT子程序画图
程序脚本，并注释
clc;
clear
n1=-40;n2=40;%定义观察起止点
n=n1:
n2;%观察范围n;
N=32;%举行矩阵的长度
M=256;
L=length（n）;%定义n的长度
R=zeros（1,L）;%赋初值
R（-n1+1:
-n1+N）=1;%定义举行矩阵R
x=（0.5*（1-cos（2*pi*n/（N-1））））.*R;%定义单位冲击函数
[X,w]=dtft2（x,n,M）;%调用函数求频率响应
magX=abs（X）;angX=angle（X）;%计算幅度和相位
G=20*log10（magX/max（magX））;%定义G为增益
subplot（3,1,1）;
plot（w,magX）,gridon;%画幅度谱
xlabel（'\omega'）;%横轴
ylabel（'|X（e^j^\omega）|'）;%纵轴为幅度谱
title（'幅频特性图形'）;%给图片取名
subplot（3,1,2）;
plot（w,angX）,gridon;%画相位谱
xlabel（'\omega'）;%横轴
ylabel（'argX（e^j^\omega）'）;%纵轴
title（'相频特性图形'）;%给图片取名
subplot（3,1,3）;
plot（w,G）,gridon;%画增益
xlabel（'\omega'）;%横轴为频率w
ylabel（'G（dB）'）;%纵轴为增益，单位为dB
title（'增益图形'）;%给图片取名增益
仿真结果、图形
图4.2.2系统的幅频特性图形、增益（dB）图形、相频特性
图4.2.3系统的衰减（dB）图形
3、判断系统的滤波特性，估计系统的第一零点带宽和第一旁瓣衰减。
由图4.2.2得是低通滤波器，
由图4.2.2得第一零点带宽为
由图4.2.3得第一旁瓣31dB衰减为
遇到的问题、解决方法及收获
1、了解了升余弦系统的滤波特性
2、了解升余弦系统第一旁瓣的衰减
任务三DFT分析信号的频谱
分析计算序列
的频谱。
1、理论计算信号的频谱（DTFT）和周期，用计算机画出其幅度谱的图形，标注坐标轴的含义和单位。
表达式
（4.3.1）
它是
的周期函数，周期为
，为了便于观察，取一个周期即可，
（4.3.2）
程序脚本，并注释
clc;
clear%清除全部
w=0:
2*pi/1024:
2*pi;%w的取值范围及间隔
L=length（w）;%定义L为w的长度
y=zeros（1,L）;%赋初值
a=pi/4/（2*pi/1024）+1;%第pi/4/（2*pi/1024）+1个点
b=7*pi/4/（2*pi/1024）+1;%第7*pi/4/（2*pi/1024）+1个点
c=pi/8/（2*pi/1024）+1;%第pi/8/（2*pi/1024）+1个点
d=14*pi/8/（2*pi/1024）+1+pi/8/（2*pi/1024）;%第15*pi/8/（2*pi/1024）+1个点
y（[abcd]）=1;%第abcd个点为一
X=abs（y）;%Xb为y的模值，即x（n）的幅度谱
stem（w,X,'.'）,gridon;%画图，为离散的
xlabel（'\omega（rad）'）;%横轴为频率w,单位为弧度
ylabel（'|X（e^j^w）|'）;%纵轴为幅度
title（'x=cos（n*\pi/4）+cos（n*\pi/8）的幅度谱'）;%给图片取名x=cos（n*\pi/4）+cos（n*\pi/8）的幅度谱
仿真结果、图形
图4.3.1
幅度谱的图形
2、将
截取有限长度
为
，编写子函数计算
的N点的DFT，画出幅度谱，标注坐标轴。
编程原理、思路和公式
将
截取有限长度
为
，则
其DFS为
（4.3.1）
截取有限长序列
的DFS则为DFT公式
（4.3.2）
将式4.3.2写成子程序需要两重循环，内层循环是对n求和的过程，循环N次，外层循环为k与X（k）一一对应形成数组。
过程与实验二写DTFT的过程一样
程序脚本，并注释
function[X,k]=dft（x,n,N）%定义函数
for（q=1:
N）%外层循环，对w的循环，循环M+1次
k（q）=q-1;%循环计算不同的w值
sum=0;%求和定初值
for（m=1:
N）%内层循环对n的求和过程，循环的次数为n的长度
sum=sum+x（m）*exp（-i*2*pi*k（q）*n（m）/N）;%求和
X（q）=sum;%把得到的和的值赋给，进行对下一个w的求和过程
end%内循环结束
end%外循环结束
1）分别取N=16，64，128，仔细观察幅度谱的图形，并与信号的理想频谱相比较，得到N的取值与信号的理想频谱的有关结论。
程序脚本，并注释
clc;
clear%清除全部
N=16;%定义N，
n1=0;%x（n）的的起点
n2=N-1;%x（n）的终点
n=n1:
1:
n2;%n的取值范围
x=cos（n*pi/4）+cos（n*pi/8）;%定义序列x
[X,k]=dft（x,n,N）%;%调用函数dft
Xb=abs（X）;%X为X的模值，即x（n）的幅度谱
subplot（3,1,1）;
stem（2*pi/N*k,Xb,'.'）,gridon;%画图，为离散的
axis（[0,7,0,10]）;
xlabel（'\omega'）;%横轴
ylabel（'dft[x（n）]'）;%纵轴为幅度
title（'N=16,x=cos（n*\pi/4）+cos（n*\pi/8）的幅度谱'）;%给图片取名幅度谱
N=64;%定义N，
n1=0;%x（n）的的起点
n2=N-1;%x（n）的终点
n=n1:
1:
n2;%n的取值范围
x=cos（n*pi/4）+cos（n*pi/8）;%定义序列x
[X,k]=dft（x,n,N）%;%调用函数dft
Xb=abs（X）;%X为X的模值，即x（n）的幅度谱
subplot（3,1,2）;
stem（2*pi/N*k,Xb,'.'）,gridon;%画图，为离散的
xlabel（'\omega'）;%横轴
ylabel（'dft[x（n）]'）;%纵轴为幅度
title（'N=64,x=cos（n*\pi/4）+cos（n*\pi/8）的幅度谱'）;%给图片取名幅度谱
N=128;%定义N，
n1=0;%x（n）的的起点
n2=N-1;%x（n）的终点
n=n1:
1:
n2;%n的取值范围
x=cos（n*pi/4）+cos（n*pi/8）;%定义序列x
[X,k]=dft（x,n,N）%;%调用函数dft
Xb=abs（X）;%X为X的模值，即x（n）的幅度谱
subplot（3,1,3）;
stem（2*pi/N*k,Xb,'.'）,gridon;%画图，为离散的
xlabel（'\omega'）;%横轴
ylabel（'dft[x（n）]'）;%纵轴为幅度
title（'N=128,x=cos（n*\pi/4）+cos（n*\pi/8）的幅度谱'）;%给图片取名幅度谱
仿真结果、图形
图4.3.2
的DFT幅度谱，从上到下N=16,64,128
图4.3.3
的DFT幅度谱，从上到下N=16,64,128
2）分别取N=4，25，75，100，仔细观察幅度谱的图形，与信号的理想频谱比较，结合1）中N的取值的结果，得到有关对正弦信号利用DFT分析频谱时，有关截取信号的长度（或DFT的点数）的结论。
（用stem画图）
clc;
clear%清除全部
N=4;%定义N，
n1=0;%x（n）的的起点
n2=N-1;%x（n）的终点
n=n1:
1:
n2;%n的取值范围
x=cos（n*pi/4）+cos（n*pi/8）;%定义序列x
[X,k]=dft（x,n,N）%;%调用函数dft
Xb=abs（X）;%X为X的模值，即x（n）的幅度谱
subplot（4,1,1）;
stem（2*pi/N*k,Xb,'.'）,gridon;%画图，为离散的
axis（[0,7,0,5]）;
xlabel（'\omega'）;%横轴
ylabel（'dft[x（n）]'）;%纵轴为幅度
title（'N=4,x=cos（n*\pi/4）+cos（n*\pi/8）的幅度谱'）;%给图片取名幅度谱
N=25;%定义N，
n1=0;%x（n）的的起点
n2=N-1;%x（n）的终点
n=n1:
1:
n2;%n的取值范围
x=cos（n*pi/4）+cos（n*pi/8）;%定义序列x
[X,k]=dft（x,n,N）%;%调用函数dft
Xb=abs（X）;%X为X的模值，即x（n）的幅度谱
subplot（4,1,2）;
stem（2*pi/N*k,Xb,'.'）,gridon;%画图，为离散的
xlabel（'\omega'）;%横轴
ylabel（'dft[x（n）]'）;%纵轴为幅度
title（'N=25,x=cos（n*\pi/4）+cos（n*\pi/8）的幅度谱'）;%给图片取名幅度谱
N=75;%定义N，
n1=0;%x（n）的的起点
n2=N-1;%x（n）的终点
n=n1:
1:
n2;%n的取值范围
x=cos（n*pi/4）+cos（n*pi/8）;%定义序列x
[X,k]=dft（x,n,N）%;%调用函数dft
Xb=abs（X）;%X为X的模值，即x（n）的幅度谱
subplot（4,1,3）;
stem（2*pi/N*k,Xb,'.'）,gridon;%画图，为离散的
xlabel（'\omega'）;%横轴
ylabel（'dft[x（n）]'）;%纵轴为幅度
title（'N=75,x=cos（n*\pi/4）+cos（n*\pi/8）的幅度谱'）;%给图片取名幅度谱
N=100;%定义N，
n1=0;%x（n）的的起点
n2=N-1;%x（n）的终点
n=n1:
1:
n2;%n的取值范围
x=cos（n*pi/4）+cos（n*pi/8）;%定义序列x
[X,k]=dft（x,n,N）%;%调用函数dft
Xb=abs（X）;%X为X的模值，即x（n）的幅度谱
subplot（4,1,4）;
stem（2*pi/N*k,Xb,'.'）,gridon;%画图，为离散的
xlabel（'\omega'）;%横轴
ylabel（'dft[x（n）]'）;%纵轴为幅度
title（'N=100,x=cos（n*\pi/4）+cos（n*\pi/8）的幅度谱'）;%给图片取名幅度谱
仿真结果、图形
图4.3.4
的4点DFT幅度谱，从上到下N为别为4,25,75,100
图4.3.5
的4点DFT幅度谱，从上到下N为别为4,25,75,100
结果分析和结论
观察幅度谱的图形4.3.3和4.3.5，与信号的理想频谱4.3.1比较，对正弦信号利用DFT分析频谱时
1、截取信号的长度（或DFT的点数）为正弦信号的周期或者周期的整数倍是，得到频谱更接近理想频谱。
如：
图4.3.3和4.3.5比较，4.3.3更接近图4.3.1，因为图4.3.3中的N的取值为16,64,128，都是x（n）的周期16的整数倍；而图4.3.5中的图形和4.3.1比较，相差很远。
2、截取信号的长度（或DFT的点数）多了比少了更接近信号的理想频谱。
图4.3.3和图4.3.5中最后一幅图比第一幅图比较更接近理想频谱。
遇到的问题、解决方法及收获
1、画4.3.3和4.3.5时，开始用k做横轴，在和理想频谱比较时，无法将各自的点一一对应，得到的的图如4.3.2所示，后来想到将横轴改为
，这样得到的图就可以和理想频谱对应比较。
2、频谱分析过程中，不是很理解N的取值造成的影响，后来查阅资料得到
3、通过此次实验，更深刻理解DFT分析频谱的方法，也直观的感受到N取值的不同对DFT分析的影响。
3、将
截取有限长度
为信号
，计算其10点和100点的DFT，分别画出幅度谱，理解栅栏效应和频率分辨率的问题（参考课本的描述）。
clc;
clear%清除全部
N=10;%定义N，
n1=0;%x（n）的的起点
n2=N-1;%x（n）的终点
n=n1:
1:
n2;%n的取值范围
R=zeros（1,N）;
R（-n1+1:
-n1+10）=1;
x=（cos（n*pi/4）+cos（n*pi/8））.*R;%定义序列x
[X,k]=dft（x,n,N）%;%调用函数dft
Xb=abs（X）;%X为X的模值，即x（n）的幅度谱
subplot（2,1,1）;
stem（k,Xb,'.'）,gridon;%画图，为离散的
xlabel（'k'）;%横轴为k
ylabel（'dft[x（n）]'）;%纵轴为幅度
title（'N=10,x=cos（n*\pi/4）+cos（n*\pi/8）的幅度谱'）;%给图片取名幅度谱
N=100;%定义N，
n1=0;%x（n）的的起点
n2=N-1;%x（n）的终点
n=n1:
1:
n2;%n的取值范围
R=zeros（1,N）;
R（-n1+1:
-n1+10）=1;
x=（cos（n*pi/4）+cos（n*pi/8））.*R;%定义序列x
[X,k]=dft（x,n,N）%;%调用函数dft
Xb=abs（X）;%X为X的模值，即x（n）的幅度谱
subplot（2,1,2）;
stem（k,Xb,'.'）,gridon;%画图，为离散的
xlabel（'k'）;%横轴为k
ylabel（'dft[x（n）]'）;%纵轴为幅度
title（'N=100,x=cos（n*\pi/4）+cos（n*\pi/8）的幅度谱'）;%给图片取名幅度谱
图4.3.6
的10点和100点DFT
图4.3.7
的10点和100点DFT
结果分析和结论
为了减小栅栏效应，可以有两个办法：
1在数据长度T不变的情况下，增加DFT变化的点数，如图4.3.3和4.3.5；2在N点有效数据的尾部增加零点，使整个数据长度为M点这就相当于DFT的变换点数为M，就如图4.3.7所示，第二个图比第一个图栅栏效应小。
提高频率分辨率的为一办法，是增加数据的有效长度T。
遇到的问题、解决方法及收获
通过此次实验，进一步理解数据长度、DFT变换点数、尾部补零、栅栏效应、频率分辨率等基本概念。
同时，对DFT频谱分析的方法有了更深一步的了解。
升级会员 