数字信号处理实验报告.docx
《数字信号处理实验报告.docx》由会员分享，可在线阅读，更多相关《数字信号处理实验报告.docx（26页珍藏版）》请在冰豆网上搜索。
数字信号处理实验报告
《数字信号处理》
实验报告
课程名称：
《数字信号处理》
学院：
信息科学与工程学院
专业班级：
通信1502班
*******
学号：
**********
指导教师：
***
2017年5月28日
实验一离散时间信号和系统响应
一.实验目的
1.熟悉连续信号经理想采样前后的频谱变化关系，加深对时域采样定理的理解
2.掌握时域离散系统的时域特性
3.利用卷积方法观察分析系统的时域特性
4.掌握序列傅里叶变换的计算机实现方法，利用序列的傅里叶变换对离散信号及系统响应进行频域分析
二、实验原理
1.采样是连续信号数字化处理的第一个关键环节。
对采样过程的研究不仅可以了解采样前后信号时域和频域特性的变化以及信号信息不丢失的条件，而且可以加深对离散傅里叶变换、Z变换和序列傅里叶变换之间关系式的理解。
对连续信号
以T为采样间隔进行时域等间隔理想采样，形成采样信号：
式中
为周期冲激脉冲，
为
的理想采样。
的傅里叶变换为
：
上式表明将连续信号
采样后其频谱将变为周期的，周期为Ωs=2π/T。
也即采样信号的频谱
是原连续信号xa（t）的频谱Xa（jΩ）在频率轴上以Ωs为周期，周期延拓而成的。
因此，若对连续信号
进行采样，要保证采样频率fs≥2fm，fm为信号的最高频率，才可能由采样信号无失真地恢复出原模拟信号
计算机实现时，利用计算机计算上式并不方便，因此我们利用采样序列的傅里叶变换来实现，即
而
为采样序列的傅里叶变换
2.时域中，描述系统特性的方法是差分方程和单位脉冲响应，频域中可用系统函数描述系统特性。
已知输入信号，可以由差分方程、单位脉冲响应或系统函数求出系统对于该输入信号的响应。
本实验仅在时域求解，对于差分方程可用Matlab中的工具箱函数filter（）函数求解
一个时域离散线性时不变系统的输出与输入间的关系为：
可用Matlab中的工具箱函数conv（）函数求解
三、实验内容及步骤
1.时域采样定理的验证
给定模拟信号：
式中
。
其幅频特性如图所示：
选择三种采样频率Fs=1kHz,300Hz,200Hz,生成采样序列
分别用序列
表示。
编写程序计算三个序列的幅频特性曲线
，并绘图显示。
观察
在折叠频率附近与连续信号频谱有无明显差别，分析频谱混叠现象。
实验程序如下
%时域采样定理的验证
%Fs=1KHz
Tp=64/1000;%Tp=64ms
Fs=1000;T=1/Fs;
M=Tp*Fs;n=0:
M-1;
A=444.128;alph=pi*50*2^0.5;omega=pi*50*2^0.5;
xnt=A*exp（-alph*n*T）.*sin（omega*n*T）;
Xk=T*fft（xnt,M）;%MµãFFT
yn='xa（nT）';subplot（3,2,1）;
stem（xnt）;%» ͼ
boxon;title（'（a）Fs=1000Hz'）;
k=0:
M-1;fk=k/Tp;
subplot（3,2,2）;plot（fk,abs（Xk））;title（'（a）T*FT[xa（nT）],Fs=1000Hz'）;
xlabel（'f（Hz）'）;ylabel（'幅度'）;
axis（[0,Fs,0,1.2*max（abs（Xk））]）
%Fs=300Hz
Tp=64/1000;
Fs=300;T=1/Fs;
M=Tp*Fs;n=0:
M-1;
A=444.128;alph=pi*50*2^0.5;omega=pi*50*2^0.5;
xnt=A*exp（-alph*n*T）.*sin（omega*n*T）;
Xk=T*fft（xnt,M）;
yn='xa（nT）';subplot（3,2,1）;
stem（xnt）;
boxon;title（'（a）Fs=300Hz'）;
k=0:
M-1;fk=k/Tp;
subplot（3,2,2）;plot（fk,abs（Xk）,'r'）;title（'（a）T*FT[xa（nT）],Fs=300Hz’）;
xlabel（'f（Hz）'）;ylabel（‘幅度’）;
axis（[0,Fs,0,1.2*max（abs（Xk））]）
%Fs=200Hz
Tp=64/1000;%64ms
Fs=300;T=1/Fs;
M=Tp*Fs;n=0:
M-1;
A=444.128;alph=pi*50*2^0.5;omega=pi*50*2^0.5;
xnt=A*exp（-alph*n*T）.*sin（omega*n*T）;
Xk=T*fft（xnt,M）;
yn='xa（nT）';subplot（3,2,1）;
stem（xnt,'.'）;
boxon;title（'（a）Fs=200Hz'）;
k=0:
M-1;fk=k/Tp;
subplot（3,2,2）;plot（fk,abs（Xk））;title（'（a）T*FT[xa（nT）],Fs=200Hz'）;
xlabel（'f（Hz）'）;ylabel（'幅度'）;
axis（[0,Fs,0,1.2*max（abs（Xk））]）;
2.给定一个低通滤波器的差分方程为：
输入序列
（1）分别求出
和
的系统响应，并画出其波形
（2）求出系统的单位脉冲响应，画出其波形
A=[1,-0.9];
B=[0.05,0.05];
x1n=[ones（1,8）,zeros（1,50）]
x2n=ones（1,200）;
hn=impz（B,A,50）;
subplot（3,1,1）;stem（hn）;
title（‘
（1）系统单位脉冲响应h（n）’）;
y1n=filter（B,A,x1n）;
subplot（3,1,2）;stem（y1n）;
title（'
（2）系统对R（8）的响应y1（n）'）;
y2n=filter（B,A,x2n）;
subplot（3,1,3）;stem（y2n）;
title（‘（3）系统对u（n）的响应）y2（n）’）;
3.给定系统的单位脉冲响应为
用线性卷积法求
分别对系统
和
的输出响应，并画出波形
x1n=ones（1,8）;
h1n=[ones（1,10）zeros（1,20）];
h2n=[1,2.5,2.5,1,zeros（1,10）];
y11n=conv（h1n,x1n）;
y22n=conv（h2n,x1n）;
subplot（2,2,1）;stem（h1n,’.b’）;
title（‘（4）系统单位脉冲响应h1（n）’）;
subplot（2,2,2）;stem（y11n,'.b'）;
title（‘（5）h1（n）与R8（n）的卷积y11（n）’）;
subplot（2,2,3）;stem（h2n,'.b'）;
title（‘（6）系统单位脉冲响应h2（n）’）;
subplot（2,2,4）;stem（y22n,'.b'）;
title（‘（7）h2（n）与R8（n）的卷积y22（n）’）;
四、实验思考
1.在分析理想采样序列特性的实验中，采样频率不同时，相应理想采样序列的傅里叶变换频谱的数字频率度量是否都相同？
它们所对应的模拟频率是否相同？
为什么？
答：
当采样频率不同时，数字度量不同，但是模拟频率相同。
因为数字频率w是模拟角频率Ω用采样频率Fs归一化频率。
数字频率和模拟角频率之间的关系是w=ΩT,模拟信号的模拟角频率Ω不变，当采样频率不同时，T不同，所以数字频率Ω不同。
因此，采样频率不同时，相应理想采样序列的傅里叶变换频谱的数字频率度量不相同，但是它们所对应的模拟频率相同。
2.如果输入信号为无线长序列，系统的单位脉冲响应是有限长序列，可否用线性卷积法求系统的响应？
如何求？
答：
（1）对输入信号序列分段；
（2）求单位脉冲响应与各段的卷积；
（3）将各段卷积结果相加。
3.如果信号经过低通滤波器，把信号的高频分量滤掉，时域信号会有何变化？
用前面第二个实验结果进行分析说明
答：
把信号经过低通滤波器，把信号的高频成分滤掉，时域信号的剧烈将变得平滑。
五、实验心得及体会
通过本次实验我重新温习了MATLAB这个软件的使用方法，运行环境。
通过这款软件使我们的学习更加便利。
实验二用FFT对信号作频谱分析
一、实验目的
1.进一步加深DFT算法原理和基本性质的理解
2.掌握用FFT对连续信号和时域离散信号进行频谱分析的方法
3.了解用FFT进行频谱分析时可能出现的分析误差及其原因，以便在实际中正确应用FFT
二、实验原理
用FFT对信号作频谱分析是学习数字信号处理的重要内容，经常需要进行谱分析的信号是模拟信号和时域离散信号。
对信号进行谱分析的重要问题是频谱分辨率F和分析误差。
频谱分辨率直接和FFT的变换区间N有关，FFT能够实现的频率分辨率是2/N，因此要求2/NF。
可以根据此式选择FFT的变换区间N。
误差主要来自于用FFT作频谱分析时，得到的是离散谱，而信号（周期信号除外）是连续谱，只有当N较大时，离散谱的包络才能逼近于连续谱，因此N要适当选择大一些。
周期信号的频谱是离散谱，只有用整数倍周期的长度作FFT，得到的离散谱才能代表周期信号的频谱。
如果不知道信号周期，可以尽量选择信号的观察时间长一些。
对模拟信号进行谱分析时，首先要按照采样定理将其变为时域离散信号。
如果是模拟周期信号，也应该选取整数倍周期的长度，经过采样后形成周期序列，按照周期序列的谱分析进行。
三、实验步骤及内容
1.对以下给出的各序列进行谱分析：
选择FFT的变换区间N为8和16两种情况进行频谱分析。
分别打印其幅频特性曲线，并进行对比、分析、讨论。
x1n=[ones（1,4）];%产生R4（n）序列向量
X1k8=fft（x1n,8）;%计算x1n的8点DFT
X1k16=fft（x1n,16）;%计算x1n的16点DFT
N=8;
f=2/N*（0:
N-1）;
figure
（1）;
subplot（1,2,1）;stem（f,abs（X1k8）,'.'）;%绘制8点DFT的幅频特性图
title（'（1a）8点DFT[x_1（n）]'）;xlabel（'ω/π'）;ylabel（'幅度'）;
N=16;
f=2/N*（0:
N-1）;
subplot（1,2,2）;stem（f,abs（X1k16）,'.'）;%绘制8点DFT的幅频特性图
title（'（1a）16点DFT[x_1（n）]'）;xlabel（'ω/π'）;ylabel（'幅度'）;
%x2n和x3n
M=8;xa=1:
（M/2）;xb=（M/2）:
-1:
1;
x2n=[xa,xb];%产生长度为8的三角波序列x2（n）
x3n=[xb,xa];
X2k8=fft（x2n,8）;
X2k16=fft（x2n,16）;
X3k8=fft（x3n,8）;
X3k16=fft（x3n,16）;
figure
（2）;
N=8;
f=2/N*（0:
N-1）;
subplot（2,2,1）;stem（f,abs（X2k8）,'.'）;%绘制8点DFT的幅频特性图
title（'（2a）8点DFT[x_2（n）]'）;xlabel（'ω/π'）;ylabel（'幅度'）;
subplot（2,2,3）;stem（f,abs（X3k8）,'.'）;%绘制8点DFT的幅频特性图
title（'（3a）8点DFT[x_3（n）]'）;xlabel（'ω/π'）;ylabel（'幅度'）;
N=16;
f=2/N*（0:
N-1）;
subplot（2,2,2）;stem（f,abs（X2k16）,'.'）;%绘制8点DFT的幅频特性图
title（'（2a）16点DFT[x_2（n）]'）;xlabel（'ω/π'）;ylabel（'幅度'）;
subplot（2,2,4）;stem（f,abs（X3k16）,'.'）;%绘制8点DFT的幅频特性图
title（'（3a）16点DFT[x_3（n）]'）;xlabel（'ω/π'）;ylabel（'幅度'）;
2.对以下各周期序列进行频谱分析
选FFT的变换区间N为8和16两种情况分别对以上序列进行频谱分析。
分别打印其幅频特性曲线，并进行对比、分析、讨论。
n=0:
8;
xn4=cos（（pi.*n）/4）;
subplot（2,3,1）;stem（n,xn4,'.'）;
X8k4=fft（xn4,8）;
n21=0:
length（X8k4）-1;
subplot（2,3,2）;stem（n21,X8k4,'.'）;
X16k4=fft（xn4,16）;
n22=0:
length（X16k4）-1;
subplot（2,3,3）;stem（n22,X16k4,'.'）;%end
n=0:
16;
xn5=cos（（pi.*n）/4）+cos（（pi.*n）/8）;
subplot（2,3,4）;stem（n,xn5,'.'）;
X8k5=fft（xn5,8）;
n21=0:
length（X8k5）-1;
subplot（2,3,5）;stem（n21,X8k5,'.'）;
X16k5=fft（xn5,16）;
n22=0:
length（X16k5）-1;
subplot（2,3,6）;stem（n22,X16k5,'.'）;
3.对模拟周期信号进行频谱分析
选择样频率Fs=64Hz，对变换区间N=16，32，64三种情况进行谱分析。
分别打印其幅频特性曲线，并进行对比、分析、讨论。
程序如下：
Fs=64;T=1/Fs;
N=16;n=0:
N-1;%FFT的变换区间N=16
x6nT=cos（8*pi*n*T）+cos（16*pi*n*T）+cos（20*pi*n*T）;%对x6（t）16点采样
X6k16=fft（x6nT）;%计算x6nT的16点DFT
X6k16=fftshift（X6k16）;%将零频率移到频谱中心
Tp=N*T;F=1/Tp;%频率分辨率F
k=-N/2:
N/2-1;fk=k*F;%产生16点DFT对应的采样点频率（以零频率为中心）
subplot（3,1,1）;stem（fk,abs（X6k16）,'.'）;
boxon%绘制8点DFT的幅频特性图
title（'（6a）16点|DFT[x_6（nT）]|'）;xlabel（'f（Hz）'）;ylabel（'幅度'）;
axis（[-N*F/2-1,N*F/2-1,0,1.2*max（abs（X6k16））]）
N=32;n=0:
N-1;%FFT的变换区间N=16
x6nT=cos（8*pi*n*T）+cos（16*pi*n*T）+cos（20*pi*n*T）;%对x6（t）32点采样
X6k32=fft（x6nT）;%计算x6nT的32点DFT
X6k32=fftshift（X6k32）;%将零频率移到频谱中心
Tp=N*T;F=1/Tp;%频率分辨率F
k=-N/2:
N/2-1;fk=k*F;%产生16点DFT对应的采样点频率（以零频率为中心）
subplot（3,1,2）;stem（fk,abs（X6k32）,'.'）;
boxon%绘制8点DFT的幅频特性图
title（'（6b）32点|DFT[x_6（nT）]|'）;xlabel（'f（Hz）'）;ylabel（'幅度'）;
axis（[-N*F/2-1,N*F/2-1,0,1.2*max（abs（X6k32））]）
N=64;n=0:
N-1;%FFT的变换区间N=16
x6nT=cos（8*pi*n*T）+cos（16*pi*n*T）+cos（20*pi*n*T）;%对x6（t）64点采样
X6k64=fft（x6nT）;%计算x6nT的64点DFT
X6k64=fftshift（X6k64）;%将零频率移到频谱中心
Tp=N*T;F=1/Tp;%频率分辨率F
k=-N/2:
N/2-1;fk=k*F;%产生16点DFT对应的采样点频率（以零频率为中心）
subplot（3,1,3）;stem（fk,abs（X6k64）,'.'）;
boxon%绘制8点DFT的幅频特性图
title（'（6a）64点|DFT[x_6（nT）]|'）;xlabel（'f（Hz）'）;ylabel（'幅度'）;
axis（[-N*F/2-1,N*F/2-1,0,1.2*max（abs（X6k64））]）
四、实验思考
1.在N=8时，
和
的幅频特性会相同吗？
为什么？
N=16时呢？
答：
在N=8时，和的幅频特性相同，而N=16时不相同。
因为=,所以和的8点DFT的模相等。
但当N=16时，和不满足循环移位关系，所以两者幅频特性不相同。
2.对于周期序列，如果周期不知道，如何用FFT进行谱分析？
答：
周期信号的周期预先不知道时,可先截取M点进行DFT,再将截取长度扩大1倍截取,比较结果,如果二者的差别满足分析误差要求,则可以近似表示该信号的频谱,如果不满足误差要求就继续将截取长度加倍,重复比较,直到结果满足要求。
五、实验总结及心得体会
通过实验进一步加深DFT算法原理和基本性质的理解，掌握用FFT对连续信号和时域离散信号进行频谱分析的方法，了解用FFT进行频谱分析时可能出现的分析误差及其原因，以便在实际中正确应用FFT。
实验三用双线性变换法设计IIR数字滤波器
一、实验目的
1.熟悉用双线性变换法设计IIR数字滤波器的原理和方法
2.掌握IIR数字滤波器的Matlab实现方法
3.通过观察对实际心电图信号的滤波作用，获得数字滤波的感性认识
二、实验原理
设计IIR数字滤波器一般采用间接设计法——脉冲响应不变法和双线性变换法，应用最广泛的是双线性变换法。
脉冲响应不变法的基本思想是：
使数字滤波器的单位脉冲响应h（n）近似于模拟滤波器的单位脉冲响应ha（t），即使
其S平面和Z平面的映射关系为：
双线性变换法的基本思想是：
使描述数字滤波器的差分方程近似描述模拟滤波器的微分方程
S平面和Z平面的映射关系为：
双线性变换法中的频率变换是一种非线性变换，这种非线性引起的幅频特性畸变可通过预变形矫正法而得到校正。
设计IIR数字滤波器的一般步骤:
（1）确定所需类型数字滤波器的技术指标：
通带截止频率ωp、通带衰减αp、阻带截止频率ωs、阻带衰减αs。
（2）将所需类型数字滤波器的技术指标转换成相应类型模拟滤波器的技术指标。
（3）设计该类型模拟滤波器
（4）通过复频率变换将模拟滤波器转换成所需类型的数字滤波器。
三、实验内容
1.分别用脉冲响应不变法和双线性变换法设计一个巴特沃斯低通IIR数字滤波器，设计指标参数为：
在通带内频率低于0.2时，最大衰减小于1dB，在阻带内[0.3，]频率区间上，最小衰减大于15dB。
观察并画出所设计数字滤波器的幅频特性曲线和相频特性曲线，记录带宽和衰减量，检查是否满足要求。
比较这两种方法的优缺点。
Matlab程序为：
%脉冲响应法
T=1;wp=0.2*pi;ws=0.3*pi;rp=1;as=15;%输入低通滤波器要求
[n,wpo]=buttord（wp,ws,rp,as,'s'）;%计算阶数
[B,A]=butter（n,wpo,'s'）;%计算表达式分子分母的系数矩阵
[B1,A1]=impinvar（B,A）;
[Hk,w]=freqz（B1,A1）;
subplot（2,1,1）;
plot（w/pi,20*log10（abs（Hk）））;%画出滤波器损耗函数曲线
gridon;
title（'
（1）脉冲响应不变法衰减曲线'）;
xlabel（'频率（w/pi）'）;
ylabel（'幅度（dB）'）;
%双线性法
wp=2*tan（0.2*pi/2）;
ws=2*tan（0.3*pi/2）;
rp=1;as=15;
[n,wpo]=buttord（wp,ws,rp,as,'s'）;%计算阶数
[B,A]=butter（n,wpo,'s'）;%计算表达式分子分母的系数矩阵
[B1,A1]=bilinear（B,A,1）;
[Hk,w]=freqz（B1,A1）;
subplot（2,1,2）;
plot（w/pi,20*log10（abs（Hk）））;%画出滤波器损耗函数曲线
title（'
（1）双线性法衰减曲线'）;
gridon;
xlabel（'频率（w/pi）'）;ylabel（'幅度（dB）'）;
优缺点比较：
（1）脉冲响应不变法会产生频谱混叠，但具有很好的线性特性，其单位脉冲响应完全模仿模拟滤波器的单位冲激响应波形，时域逼近性好。
适合于带通、低通滤波器的设计。
（2）双线性法很好地消除了频谱混叠，但是其数字频率与模拟频率之间不具有线性关系。
2.用双线性变换法设计一个切比雪夫高通IIR数字滤波器，设计指标参数为：
在通带内频率高于0.3KHz,最大衰减小于1dB，在阻带内频率低于0.2KHz，最小衰减大于20dB，T=1ms。
画出所设计数字滤波器的幅频特性曲线和相频特性曲线，观察其通带损耗和阻带衰减是否满足要求。
Matlab程序如下：
%切比雪夫高通滤波器的设计
fp=2*pi*300*0.001;
fs=2*pi*200*0.001;
wp=2000*tan（fp/2）;
ws=2000*tan（fs/2）;%进行频率变换
rp=1;as=20;%输入高通滤波器要求
[n,wpo]=buttord（wp,ws,rp,as,'s'）;%计算阶数
[B,A]=butter（n,wpo,'high','s'）;%计算表达式分子分母的系数矩阵
[B1,A1]=impinvar（B,A）;
[Hk,w]=freqz（B,A）;
plot（w/pi,20*log10（abs（Hk）））;%画滤波器损耗函数曲线
gridon;
title（'（3）切比雪夫高通IIR数字滤波器'）;
xlabel（'频率（HZ）'）;ylabel（'幅度（dB）'）;
3.人体心电图信号在测量过程中往往受到工业高频干扰，所以必须经过低通滤波处理后，才能作为判断心脏功能的有用信息。
下面给出一实际心电图信号采样序列样本x（n），其中存在高频干扰。
用1所计的滤波器对心电图信号采样序列x（n）进行仿真滤波处理，画出处理前后的信号波形。
xn=[-4,-2,0,-4,-6,-4,-2,-4,-6,-6,-4,-4,-6,-6,-2,6,12,8,0,-16,-38,-60,84,-90,-66,-32,-4,-2,-4,8,12,12,10,6,6,6,4,0,0,0,0,0,-2,-4,0,0,0,-2,-2,0,0,-2,-2,-2,-2,0];
subplot（2
升级会员 