数字信处理实验.docx
《数字信处理实验.docx》由会员分享，可在线阅读，更多相关《数字信处理实验.docx（12页珍藏版）》请在冰豆网上搜索。
数字信处理实验
实验一：
系统及响应时域采样及频域采样
1.实验目的
（1）掌握用卷积求系统响应及卷积定理的验证；
（2）掌握连续信号经理想采样前后的频谱变化关系，加深对时域采样定理的理解。
（3）掌握频域采样引起时域周期化概念，加深对频域采样定理的理解。
（4）掌握序列傅里叶变换的计算机实现方法，利用序列的傅里叶变换对连续信号、离散信号及系统响应进行频域分析。
3.实验内容及步骤
（1）认真复习卷积定理、时域采样和频域采样理论。
（2）编制实验用主程序及相应子程序。
①系统单位脉冲响应序列产生子程序。
有限长序列线性卷积子程序，
用于完成两个给定长度的序列的卷积。
可以直接调用MATLAB语言中的卷积函数conv。
conv用于两个有限长度序列的卷积，它假定两个序列
都从n=0开始。
调用格式如下：
y=conv（x,h）
②卷积定理的验证。
（3）时域采样定理的验证：
信号产生子程序，用于产生实验中要用到的下列信号序列：
xa（t）=Ae-atsin（Ω0t）u（t）
进行采样，可得到采样序列
xa（n）=xa（nT）=Ae-anTsin（Ω0nT）u（n）,0≤n<50
其中A为幅度因子，a为衰减因子，Ω0是模拟角频率，T为采样间隔。
这些参数都要在实验过程中由键盘输入，产生不同的xa（t）和xa（n）。
>>%1时域采样序列分析
A=400;a=200;w=200;
n=0:
50-1;fs=1000;
xa=A*exp（（-a）*n/fs）.*sin（w*n/fs）;
k=-200:
200;w=（pi/100）*k;
Xk=fft（xa,length（k））;magX=abs（Xk）;angX=angle（Xk）;
subplot（2,1,1）;
stem（n,xa,'.'）;xlabel（'n'）;ylabel（'xa（n）'）;
title（'信号的类型'）;
subplot（2,1,2）;plot（w/pi,magX）;xlabel（'w/pi'）;
ylabel（'|Yjw|'）;title（'Y（|jw|）'）;
（4）频域采样定理的验证：
>>%1时域采样序列分析fs=1000
A=400;a=200;w=200;;
ts=64*10^（-3）;fs=1000;T=1/fs;
n=0:
ts/T-1;xn=A*exp（（-a）*n/fs）.*sin（w*n/fs）;
Xk=fft（xn）;
subplot（3,2,1）;stem（n,xn）;xlabel（'n,fs=1000Hz'）;ylabel（'xn'）;title（'xn'）;
subplot（3,2,2）;plot（n,abs（Xk））;xlabel（'k,fs=1000Hz'）;title（'|X（k）|'）;
>>%频域采样定理验证
M=26;N=32;n=0:
M;n1=0:
13;x1=n1+1;
n2=14:
26;x2=27-n2;
x=[x1,x2];Xk=fft（x,512）;
X32k=fft（x,32）;
k=0:
511;w=（pi/512）*k;
subplot（321）;stem（n,x）;xlabel（'n'）;
ylabel（'xn'）;axis（[0,31,0,15]）;
subplot（322）;plot（w,abs（Xk））;xlabel（'k'）;
ylabel（'|X（k）|'）;axis（[0,1,0,200]）
X16k=X32k（1:
2:
N）;
x32n=ifft（X32k）;x16n=ifft（X16k,16）;
k1=0:
31;k2=0:
15;
subplot（323）;stem（k1,abs（X32k））;xlabel（'k'）;
ylabel（'X32k'）;axis（[0,31,0,200]）;
subplot（325）;stem（k2,abs（X16k））;xlabel（'k'）;
ylabel（'|X（k）|'）;axis（[0,15,0,200]）
n=0:
31;
subplot（324）;stem（n,abs（x32n））;xlabel（'n'）;
ylabel（'|x（n）|'）;axis（[0,31,0,15]）
n1=0:
15;
subplot（326）;stem（n1,abs（x16n））;xlabel（'n'）;
ylabel（'|x（n）|'）;axis（[0,31,0,15]）
实验二：
用FFT作谱分析
1.实验目的
（1）进一步加深DFT算法原理和基本性质的理解（因为FFT只是DFT的一种快速算法，所以FFT的运算结果必然满足DFT的基本性质）。
（2）熟悉FFT算法原理和FFT子程序的应用。
（3）学习用FFT对连续信号和时域离散信号进行谱分析的方法，了解可能出现的分析误差及其原因，以便在实际中正确应用FFT。
2.实验步骤
（1）复习DFT的定义、性质和用DFT作谱分析的有关内容。
（2）复习FFT算法原理与编程思想，并对照DIT-FFT运算流图和程序框图，读懂本实验提供的FFT子程序。
（3）编制信号产生子程序，产生以下典型信号供谱分析用：
（4）编写程序。
（5）按实验内容要求，上机实验，并写出实验报告。
>>%
x1=[11110000];
x2=[12344321];
x3=[43211234];
x4=cos*pi*n）;
N=8;n=0:
7;k=0:
7;
X1k=fft（x1,N）;
subplot（2,2,1）;stem（n,x1,'.'）;
xlabel（'n'）;ylabel（'|x1（n）|'）;
subplot（2,2,2）;stem（k,abs（X1k）,'.'）;
xlabel（'k'）;ylabel（'|X1（k）|'）;
X2k=fft（x2,N）;
subplot（2,2,3）;stem（n,x2,'.'）;
xlabel（'n'）;ylabel（'|x2（n）|'）;
subplot（2,2,4）;stem（k,abs（X2k）,'.'）;
xlabel（'k'）;ylabel（'|X2（k）|'）;
>>x1=[11110000];
x2=[12344321];
x3=[43211234];
x4=cos*pi*n）;
N=8;n=0:
7;k=0:
7;
figure
（2）
X3k=fft（x3,N）;
subplot（2,2,1）;stem（n,x3,'.'）;
xlabel（'n'）;ylabel（'|x3（n）|'）;
subplot（2,2,2）;stem（k,abs（X3k）,'.'）;
xlabel（'k'）;ylabel（'|X3（k）|'）;
X2k=fft（x4,N）;
subplot（2,2,3）;stem（n,x4,'.'）;
xlabel（'n'）;ylabel（'|x4（n）|'）;
subplot（2,2,4）;stem（k,abs（X2k）,'.'）;
xlabel（'k'）;ylabel（'|X4（k）|'）;
>>%ex3
（2）main
fs=64;N=16;
n=0:
N-1;k=n;
x5=cos（n*pi/4）+cos（n*pi/8）;
x6=cos（8*pi*n/fs）+cos（16*pi*n/fs）+cos（20*pi*n/fs）;
X5k=fft（x5,N）;
X6k=fft（x6,N）;
figure（3）
subplot（2,2,1）;stem（n,x5,'.'）;
xlabel（'n'）;ylabel（'|x5（n）|'）;
subplot（2,2,2）;stem（abs（X5k）,'.'）;
xlabel（'k'）;ylabel（'|X5（k）|'）;
subplot（2,2,3）;stem（n,x6,'.'）;
xlabel（'n'）;ylabel（'|x6（n）|'）;
subplot（2,2,4）;stem（abs（X6k）,'.'）;
xlabel（'k'）;ylabel（'|X6（k）|'）;
总结
通过这次实验我学到了：
1、MATLAB中程序的调试：
M文件中，按下F5设置断点，然后F10运行就可以调试自己需要的程序了。
2、时域采样定理，采样频率s必须大于等于模拟信号最高频率的两倍以上，才能使采样信号的频谱不产生频谱混叠；频域采样原理，频域采样点数N必须大于等于时域离散信号的长度M（即N≥M），才能使时域不产生混叠。
在数字信号处理的应用中，只要涉及时域或者频域采样，都必须服从这两个采样理论。
并且可知“时域采样频谱周期延拓，频域采样时域信号周期延拓”。
3、周期信号的频谱是离散谱，只有用整数倍周期的长度作FFT，得到的离散谱才能代表周期信号的频谱。
如果不知道信号周期，可以尽量选择信号的观察时间长一些。
对模拟信号进行谱分析时，首先要按照采样定理将其变成时域离散信号。
如果是模拟周期信号，也应该选取整数倍周期的长度，经过采样后形成周期序列，按照周期序列的谱分析进行。
误差主要来自于用FFT作频谱分析时，得到的是离散谱，而信号（周期信号除外）是连续谱，只有当N较大时离散谱的包络才能逼近于连续谱，因此N要适当选择大一些。
实验三：
用双线性变换法设计
IIR数字滤波器
1.实验目的
（1）熟悉用双线性变换法设计IIR数字滤波器的原理与方法。
（2）掌握数字滤波器的计算机仿真方法。
（3）通过观察对实际心电图信号的滤波作用，获得数字滤波的感性知识。
2.实验内容
（1）用双线性变换法设计一个巴特沃斯低通IIR数字滤波器。
设计指标参数为：
在通带截止频率为π，最大衰减为1dB；阻带截止频率为π，最小衰减为15dB。
（2）以π为采样间隔，打印出数字滤波器在频率区间［0，π/2］上的幅频响应特性曲线。
（3）用所设计的滤波器对实际心电图信号采样序列（在本实验后面给出）进行仿真滤波处理，并分别打印出滤波前后的心电图信号波形图，观察总结滤波作用与效果。
（4）设计一个工作于采样频率5MHz的椭圆数字带通滤波器,要求通带边界频率为560kHz和780kHz，通带最大衰减为1dB，阻带边界频率为375kHz和1MHz，阻带最小衰减为50dB，调用MATLAB工具箱ellipord和ellip设计，并显示数字滤波器的系统函数H（z）的系数，绘出幅频特性和相频特性。
（5）设计一个工作于采样频率2500kHz的椭圆高通数字滤波器,要求通带边界频率为325kHz，通带最大衰减为1dB，阻带边界频率为225kHz，阻带最小衰减为40dB，调用MATLAB工具箱ellipord和ellip设计，并显示数字滤波器的系统函数H（z）的系数，绘出幅频特性和相频特性。
3.实验步骤
（1）复习有关巴特沃斯模拟滤波器设计和用双线性变换法设计IIR数字滤波器的内容
（2）编写滤波器仿真程序,计算H（z）对心电图信号采样序列x（n）的响应序列y（n）
（3）在通用计算机上运行仿真滤波程序，并调用通用绘图子程序，完成实验内容
（2）到（5）。
4.思考题
用双线性变换法设计数字滤波器过程中，变换公式
中T的取值，对设计结果有无影响?
为什么?
5.实验报告要求
（1）简述实验目的及原理。
（2）由所打印的|H（ejω）|特性曲线及设计过程简述双线性变换法的特点。
（3）对比滤波前后的心电图信号波形，说明数字滤波器的滤波过程与滤波作用。
（4）简要回答思考题。
>>wp=*pi;
ws=*pi;
Rp=1;
As=15;T=2;
wp1=（2/T）*tan（wp/2）;%PrewarpPrototypePassbandfreq
ws1=（2/T）*tan（ws/2）;
[N,wc]=buttord（wp1,ws1,Rp,As,'s'）
[B,A]=butter（N,wc,'s'）
[b,a]=bilinear（B,A,T）;
x=[-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];
k=1;
closeall;
figure
（1）
subplot（2,2,1）
n=0:
55;
stem（n,x,'.'）;
axis（[056-10050]）;
holdon;
n=0:
60;
m=zeros（61）;
plot（n,m）;
xlabel（'n'）;
ylabel（'x（n）'）;
title（'心电图信号采样序列x（n）'）;
y=filter（b,a,x）;
subplot（2,2,3）
n=0:
55;
stem（n,y,'.'）;
axis（[056-155]）;
holdon;
n=0:
60;
m=zeros（61）;
plot（n,m）;
xlabel（'n'）;
ylabel（'y（n）'）;
title（'三级滤波后的心电图信号'）;
axis（[056-10050]）;
[H,w]=freqz（b,a,100）;
mag=abs（H）;
db=20*log10（（mag+eps）/max（mag））;
subplot（2,2,2）
plot（w/pi,db）;
axis（[0,,-50,10]）;
title（'滤波器的幅频响应'）;
N=
6
wc=
B=
000000
A=
>>
实验四FIR数字滤波器的设计
1.实验目的
（1）掌握用窗函数法设计FIR数字滤波器的原理和方法。
（2）熟悉线性相位FIR数字滤波器特性。
（3）了解各种窗函数对滤波特性的影响。
2.实验原理与方法
如果所希望的滤波器的理想频率响应函数为
Hd（ejω），则其对应的单位脉冲响应为
用窗函数w（n）将hd（n）截断，并进行加权处理，得到:
h（n）就作为实际设计的FIR数字滤波器的单位脉
冲响应序列，其频率响应函数H（ejω）为
如果要求线性相位特性，则h（n）还必须满足：
根据上式中的正、负号和长度N的奇偶性又将线性相位FIR滤波器分成四类。
要根据所设计的滤波特性正确选择其中一类。
例如，要设计线性相位低通特性，可选择h（n）=h（N-1-n）一类，而不能选h（n）=-h（N-1-n）一类。
3.实验内容及步骤
（1）复习用窗函数法设计FIR数字滤波器一节内容，阅读本实验原理，掌握设计步骤。
（2）编写程序。
①要求不调用fir1函数，直接按照窗函数设计法编程，用矩型窗、hanning窗、hamming窗和blackman窗设计FIR低通滤波器，要求编写求理想hd（n）的子程序和主程序，技术指标：
通带截止频率wp=*pirad;阻带截止频率ws=*pirad，分别画出各h（n），幅频特性、相频特性及衰减特性。
用窗函数法设计滤波器主程序框图
②调用fir1函数设计上述低通FIR滤波器，分别画出各h（n），幅频特性、相频特性及衰减特性。
③调用remezord和remez设计FIR高通滤波器，要求：
采样频率为16kHz，通带截止频率为，通带衰减为1dB，过渡带小于，阻带衰减为75dB，分别画出各h（n），幅频特性、相频特性及衰减特性。
④用频率采样法设计FIR带通滤波器，N＝33，理想幅度为：
>>%ex52hanningwindow
>>clearall
>>wp=*pi;ws=*pi;As=40;
Bt=ws-wp;
N0=ceil*pi/Bt）;%向上取整
N=N0+mod（N0+1,2）;%确保N为奇数
wc=（wp+ws）/2/pi;
n=0:
N-1;
hn=fir1（N-1,wc,hanning（N））;
fh0=fft（hn,1024）;
fh=20*log10（abs（fh0））;
wk=2*[0:
1023]/1024;
subplot（221）;
stem（n,hn）
title（'h（n）'）;
subplot（222）;
plot（wk,fh）;grid;
title（'衰减特性'）;xlabel（'w/pi'）
axis（[0,1,-150,0]）
pha=angle（fh0）;
subplot（223）
plot（wk,abs（fh0））
title（'幅频特性'）;xlabel（'w/pi'）
axis（[0,1,0,]）
subplot（224）;
plot（wk,pha）
axis（[0,1,-4,4]）
title（'相频特性'）;xlabel（'w'）
axis（[0,1,-4,4]）
思考题
实验三
用双线性变换法设计数字滤波器过程中，变换公式
中T的取值，对设计结果有无影响?
为什么?
答：
无影响。
依靠双线性变换是建立起来s平面和z平面的单值映射关系，因此可以有效避免频谱混叠现象，无论T取何值都是单值映射关系，对设计结果不会有影。
实验四
（1）如果给定通带截止频率和阻带截止频率以及阻带最小衰减，如何用窗函数法设计线性相位低通滤波器?
写出设计步骤。
技术指标Wp=*pi,Ws=*pi,Ap=,As=50dB
选择海明窗
clearall;
Wp=*pi;
Ws=*pi;
tr_wide=Ws-Wp;%过渡带宽度
N=ceil*pi/tr_wide）+1;%滤波器长度
n=0:
1:
N-1;
Wc=（Wp+Ws）/2;%理想低通滤波器的截止频率
hd=ideal_lp1（Wc,N）;%理想滤波器的单位冲击响应
w_ham=（hamming（N））';%海明窗
h=hd.*w_ham;%实际海明窗的响应
[db,mag,pha,w]=freqz_m2（h,[1]）;%计算实际滤波器的幅度响应
delta_w=2*pi/1000;
Ap=-（min（db（1:
1:
Wp/delta_w+1）））%实际通带纹波
As=-round（max（db（Ws/delta_w+1:
1:
501）））%实际阻带纹波
subplot（221）
stem（n,hd）
title（'理想单位脉冲响应hd（n）'）
subplot（222）
stem（n,w_ham）
title（'海明窗'）
subplot（223）
stem（n,h）
title（'实际单位脉冲响应hd（n）'）
subplot（224）
plot（wi/pi,db）
title（'幅度响应（dB）'）
axis（[0,1,-100,10]）
（2）如果要求用窗函数法设计带通滤波器，且给定上、下边带截止频率为ω1和ω2，试求理想带通的单位脉冲响应hd（n）。
解：
理想线性相位带通滤波器的频率响应为
求单位冲激响应hd（n）,即
升级会员 