数字信号处理上机作业.docx
《数字信号处理上机作业.docx》由会员分享，可在线阅读，更多相关《数字信号处理上机作业.docx（27页珍藏版）》请在冰豆网上搜索。
数字信号处理上机作业
数字信号处理上机作业
学院:
电子工程学院
班级:
021215
组员:
实验一:
信号、系统及系统响应
1、实验目的
（1）熟悉连续信号经理想采样前后的频谱变化关系,加深对时域采样定理的理解。
（2）熟悉时域离散系统的时域特性。
（3）利用卷积方法观察分析系统的时域特性。
（4）掌握序列傅里叶变换的计算机实现方法,利用序列的傅里叶变换对连续信号、离散信号及系统响应进行频域分析。
2、实验原理与方法
（1）时域采样。
（2）LTI系统的输入输出关系。
3、实验内容及步骤
（1）认真复习采样理论、离散信号与系统、线性卷积、序列的傅里叶变换及性质等有关内容,阅读本实验原理与方法。
（2）编制实验用主程序及相应子程序。
①信号产生子程序,用于产生实验中要用到的下列信号序列:
a、xa（t）=A*e^-at*sin（Ω0t）u（t）
b、单位脉冲序列:
xb（n）=δ（n）
c、矩形序列:
xc（n）=RN（n）,N=10
②系统单位脉冲响应序列产生子程序。
本实验要用到两种FIR系统。
a、ha（n）=R10（n）;
b、hb（n）=δ（n）+2、5δ（n-1）+2、5δ（n-2）+δ（n-3）
③有限长序列线性卷积子程序
用于完成两个给定长度的序列的卷积。
可以直接调用MATLAB语言中的卷积函数conv。
conv用于两个有限长度序列的卷积,它假定两个序列都从n=0开始。
调用格式如下:
y=conv（x,h）
4、实验结果分析
①分析采样序列的特性。
a、取采样频率fs=1kHz,,即T=1ms。
b、改变采样频率,fs=300Hz,观察|X（e^jω）|的变化,并做记录（打印曲线）;进一步降低采样频率,fs=200Hz,观察频谱混叠就是否明显存在,说明原因,并记录（打印）这时的|X（e^jω）|曲线。
程序代码如下:
closeall;clearall;clc;
A=50;
a=50*sqrt
（2）*pi;
m=50*sqrt
（2）*pi;
fs1=1000;
fs2=300;
fs3=200;
T1=1/fs1;
T2=1/fs2;
T3=1/fs3;
N=100;
n=[0:
N-1];
x1=A*exp（-a*n*T1）、*sin（m*n*T1）;
x2=A*exp（-a*n*T2）、*sin（m*n*T2）;
x3=A*exp（-a*n*T3）、*sin（m*n*T3）;
w=linspace（-pi,pi,10000）;
X1=x1*exp（-j*n'*w）;
X2=x2*exp（-j*n'*w）;
X3=x3*exp（-j*n'*w）;
figure
（1）
subplot（1,3,1）
plot（w/pi,abs（X1））;
xlabel（'\omega/π'）;ylabel（'|H（e^j^\omega）|'）
title（'采样频率为1000Hz时的频谱图'）;
subplot（1,3,2）
plot（w/pi,abs（X2））;
xlabel（'\omega/π'）;ylabel（'|H（e^j^\omega）|'）
title（'采样频率为300Hz时的频谱图'）;
subplot（1,3,3）
plot（w/pi,abs（X3））;
xlabel（'\omega/π'）;ylabel（'|H（e^j^\omega）|'）
title（'采样频率为200Hz时的频谱图'）;
②时域离散信号、系统与系统响应分析。
a、观察信号xb（n）与系统hb（n）的时域与频域特性;利用线性卷积求信号xb（n）通过系统hb（n）的响应y（n）,比较所求响应y（n）与hb（n）的时域及频域特性,注意它们之间有无差别,绘图说明,并用所学理论解释所得结果。
b、观察系统ha（n）对信号xc（n）的响应特性。
程序代码如下:
closeall;clearall;clc;
xbn=[1];
xcn=ones（1,10）;
han=ones（1,10）;
hbn=[1,2、5,2、5,1];
yn=conv（xbn,hbn）;
n1=0:
length（xbn）-1;
n2=0:
length（hbn）-1;
subplot（2,1,1）;stem（n1,xbn,'、'）
xlabel（'n'）;ylabel（'xb（n）'）
title（'xb（n）的时域特性曲线'）
subplot（2,1,2）;stem（n2,hbn,'、'）
xlabel（'n'）;ylabel（'hb（n）'）
title（'hb（n）的时域特性曲线'）
figure
（2）
subplot（2,1,1）;
n1=[0:
length（xbn）-1];
w=linspace（-pi,pi,10000）;
Xbn=xbn*exp（-j*n1'*w）;
subplot（2,1,1）;
plot（w/pi,abs（Xbn））;
xlabel（'\omega/π'）;ylabel（'幅度'）
title（'DTFT[xb[n]的频谱'）;
n2=[0:
length（hbn）-1];
w=linspace（-pi,pi,10000）;
Hb=hbn*exp（-j*n2'*w）;
subplot（2,1,2）;
plot（w/pi,abs（Hb））;
xlabel（'\omega/π'）;ylabel（'幅度'）
title（'DTFT[hb（n）]的频谱'）;
figure（3）
n1=0:
length（yn）-1;
n2=0:
length（hbn）-1;
subplot（2,1,1）;stem（n1,yn,'、'）
xlabel（'n'）;ylabel（'y（n）'）
title（'y（n）的时域特性曲线'）
subplot（2,1,2）;stem（n2,hbn,'、'）
xlabel（'n'）;ylabel（'hb（n）'）
title（'hb（n）的时域特性曲线'）
figure（4）
subplot（2,1,1）;
n1=[0:
length（yn）-1];
w=linspace（-pi,pi,10000）;
Y=yn*exp（-j*n1'*w）;
subplot（2,1,1）;
plot（w/pi,abs（Y））;
xlabel（'\omega/π'）;ylabel（'幅度'）
title（'DTFT[y（n）]的频谱'）;
n2=[0:
length（hbn）-1];
w=linspace（-pi,pi,10000）;
Hb=hbn*exp（-j*n2'*w）;
subplot（2,1,2）;
plot（w/pi,abs（Hb））;
xlabel（'\omega/π'）;ylabel（'幅度'）
title（'DTFT[hb（n）]的频谱'）;
zn=conv（xcn,han）;
%以下为%系统ha（n）对信号xc（n）的频率响应特性的分析
figure（5）
n3=[0:
length（zn）-1];
w=linspace（-pi,pi,10000）;
Z=zn*exp（-j*n3'*w）;
subplot（2,1,1）;
plot（w/pi,abs（Z））;
xlabel（'\omega/π'）;ylabel（'幅度'）
title（'DTFT[zn]的幅度谱'）;
subplot（2,1,2）;
plot（w/pi,angle（Z））;
xlabel（'\omega/π'）;ylabel（'相位'）
title（'DTFT[zn]的相位谱'）;
图一:
xb（n）与hb（n）时域曲线的对比
图二:
xb（n）与hb（n）频谱函数的对比
图三:
y（n）与hb（n）时域曲线的对比
图四:
y（n）与hb（n）频谱函数的对比
图五:
系统ha（n）对信号xc（n）的频率响应特性的分析
包含z（n）幅度谱与相位谱
③卷积定理的验证。
（如下程序中Y1为直接线性卷积所得结果做DTFT,Y2为各自的DTFT相乘后所得的结果）
程序代码如下:
closeall;clearall;clc;
x1=[1,3,4,0,2];
x2=[3,2,0,0,1,5];
n1=0:
length（x1）-1;
n2=0:
length（x2）-1;
y1=conv（x1,x2）;
m1=0:
length（y1）-1;
w=linspace（-pi,pi,10000）;
X1=x1*exp（-j*n1'*w）;
X2=x2*exp（-j*n2'*w）;
Y1=y1*exp（-j*m1'*w）;
Y2=X1、*X2;
subplot（2,1,1）;
plot（w/pi,abs（Y1））;
xlabel（'\omega/π'）;ylabel（'幅度'）
title（'DTFT[y1（n）]的特性曲线'）;
subplot（2,1,2）;
plot（w/pi,abs（Y2））;
xlabel（'\omega/π'）;ylabel（'幅度'）
title（'DTFT[y2（n）]的特性曲线'）;
5、思考题
（1）在分析理想采样序列特性的实验中,采样频率不同时,相应理想采样序列的傅里叶变换频谱的数字频率度量就是否都相同？
它们所对应的模拟频率就是否相同？
为什么？
答:
数字频率度量不相同,但她们所对应的模拟频率相同。
由w=Ω*Ts得,采样间隔变化时模拟频率对应的数字频率会有相应的变化,故其度量会有所变化。
（2）在卷积定理验证的实验中,如果选用不同的频域采样点数M值,例如,选M=10与M=20,分别做序列的傅里叶变换,求得的结果有无差异？
答:
有差别。
DFT相当于序列频谱的等间隔采样,当取点少时,DFT包含序列频谱的信息会少,与序列频谱的误差会增大。
取点多时,DFT会反映更多的频谱信息,误差会小,减轻了栅栏效应。
二者因为误差而有差别。
实验二:
用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）按实验内容要求,上机实验,并写出实验报告。
3、实验内容
（1）对2中所给出的信号逐个进行谱分析。
程序代码如下:
closeall;clearall;clc;
N=32;
n=0:
999;
x1n=ones（1,4）;
x2n=[1,2,3,4,4,3,2,1];
x3n=[4,3,2,1,1,2,3,4];
x4n=cos（0、25*pi*n）;
x5n=sin（0、125*pi*n）;
x6n=cos（8*pi*n）+cos（16*pi*n）+cos（20*pi*n）;
X1=fft（x1n,N）;
X2=fft（x2n,N）;
X3=fft（x3n,N）;
X4=fft（x4n,N）;
X5=fft（x5n,N）;
X6=fft（x6n,N）;
k=0:
N-1;
figure
（1）
stem（k,abs（X1））;
xlabel（'k'）;ylabel（'|X1（k）|'）;
title（'DFT（x1n）的幅频特性图'）;
figure
（2）
stem（k,abs（X2））;
xlabel（'k'）;ylabel（'|X2（k）|'）;
title（'DFT（x2n）的幅频特性图'）;
figure（3）
stem（k,abs（X3））;
xlabel（'k'）;ylabel（'|X3（k）|'）;
title（'DFT（x3n）的幅频特性图'）;
figure（4）
stem（k,abs（X4））;
xlabel（'k'）;ylabel（'|X4（k）|'）;
title（'DFT（x4n）的幅频特性图'）;
figure（5）
stem（k,abs（X5））;
xlabel（'k'）;ylabel（'|X5（k）|'）;
title（'DFT（x5n）的幅频特性图'）;
figure（6）
stem（k,abs（X6））;
xlabel（'k'）;ylabel（'|X6（k）|'）;
title（'DFT（x6n）的幅频特性图'）;
实验结果如下:
x1n到x6n的32点DFT幅频特性图
（2）令x（n）=x4（n）+x5（n）,用FFT计算8点与16点离散傅里叶变换。
closeall;clearall;clc;
N=8;
M=16;
n=0:
10000;
x4n=cos（0、25*pi*n）;
x5n=sin（0、125*pi*n）;
xn=x4n+x5n;
X1=fft（xn,N）;
X2=fft（xn,M）;
k1=0:
N-1;
k2=0:
M-1;
figure
（1）
stem（k1,abs（X1））;
xlabel（'k'）;ylabel（'|X1（k）|'）;
title（'8点DFT[x（n）]的幅频特性图'）;
figure
（2）
stem（k2,abs（X2））;
xlabel（'k'）;ylabel（'|X2（k）|'）;
title（'16点DFT[x（n）]的幅频特性图'）;
（3）令x（n）=x4（n）+j*x5（n）,用FFT计算8点与16点离散傅里叶变换。
closeall;clearall;clc;
N=8;
M=16;
n=0:
10000;
x4n=cos（0、25*pi*n）;
x5n=sin（0、125*pi*n）;
xn=x4n+j*x5n;
X1=fft（xn,N）;
X2=fft（xn,M）;
k1=0:
N-1;
k2=0:
M-1;
figure
（1）
stem（k1,abs（X1））;
xlabel（'k'）;ylabel（'|X1（k）|'）;
title（'8点DFT[x（n）]的幅频特性图'）;
figure
（2）
stem（k2,abs（X2））;
xlabel（'k'）;ylabel（'|X2（k）|'）;
title（'16点DFT[x（n）]的幅频特性图'）;
4、思考题
（1）N=8时,序列的幅频特性相同,因为序列x2（n）就是1,2,3,4,4,3,2,1正好8个点。
x3（n）就是4,3,2,1,1,2,3,4八个点,序列x2（（n－4））8＝x3（n）,因此根据DFT循环移位,得幅频特性相同。
而N=16时,通过对x2（n）与x3（n）末端环移位特性,并且补零,两个序列不能通过移位相等,因此幅频特性不同。
（2）如果周期信号的周期预先不知道,如何用FFT进行谱分析?
可以预先设定一个较小的N,进行N点FFT,再取2N,进行2N点FFT,如果二个FFT相等,则周期可能就是N,直接取此结果作为FFT。
如果不等,取4N进行FFT,并与2N点FFT比较,依次进行,取2N直到FFT近似相等为止,求的近似的FFT。
实验三:
用窗函数法设计FIR数字滤波器
1、实验目的
（1）掌握用窗函数法设计FIR数字滤波器的原理与方法。
（2）熟悉线性相位FIR数字滤波器特性。
（3）了解各种窗函数对滤波特性的影响。
2、实验内容及步骤
（1）复习用窗函数法设计FIR数字滤波器一节内容,阅读本实验原理,掌握设计步骤。
（2）编写程序。
①编写能产生矩型窗、升余弦窗、改进升余弦窗与二阶升余弦窗的窗函数子程序。
可调用子函数如下程序给出
functionwn=wn（n,num）
switchnum
case1
wn=（boxcar（n））'
case2
wn=（hanning（n））'
case3
wn=（hamming（n））'
case4
wn=（blackman（n））'
end
end
下面用程序给出的就是四种窗的单位脉冲响应波形与幅频特性曲线
以便于研究不同窗的特点
closeall;clearall;clc;
N=64;beita=2、5;
w1=boxcar（N）;
w2=hanning（N）;
w3=hamming（N）;
w4=blackman（N）;
wvtool（w1）;
wvtool（w2）;
wvtool（w3）;
wvtool（w4）;
②编写主程序。
其中幅度特性要求用dB表示。
在此程序中设N=10;wc=pi;且分别实验四种窗函数所得的滤波结果。
closeall;clearall;clc;
n=10;
b=0:
9
a=（n-1）/2;
wc=pi;
hdn=sin（wc*（b-a））、/（pi*（b-a））;
w1=wn（n,1）;
w2=wn（n,2）;
w3=wn（n,3）;
w4=wn（n,4）;
h1=hdn、*w1;
h2=hdn、*w2;
h3=hdn、*w3;
h4=hdn、*w4;
H1=fft（h1,300）;
H2=fft（h2,300）;
H3=fft（h3,300）;
H4=fft（h4,300）;
F1=abs（H1）;
K1=20*log10（F1）;
F2=abs（H2）;
K2=20*log10（F2）;
F3=abs（H3）;
K3=20*log10（F3）;
F4=abs（H4）;
K4=20*log10（F4）;
subplot（221）;
plot（abs（K1））;
title（'矩形窗'）;
subplot（222）;
plot（abs（K2））;
title（'升余弦窗'）;
subplot（223）;
plot（abs（K3））;
title（'改进升余弦窗'）;
subplot（224）;
plot（abs（K4））;
title（'二阶升余弦窗'）;
加四种窗结果分别就是:
3、思考题
（1） 给定通带截止频率与阻带截止频率以及阻带最小衰减,用窗函数法设计线性相位低通滤波器的设计步骤:
①根据所需设计的数字滤波器类型,确定线性相位数字滤波器类型,本题为I型;
②选择合适的窗函数;
③确定理想低通数字滤波器的频率响应函数;
④计算理想低通数字滤波器的单位脉冲响应;
⑤加窗得到设计结果;
（2）用窗函数法设计带通滤波器, 且给定上、 下边带截止频率为ω1与ω2,理想带通的单位脉冲响应hd（n）的求解过程:
升级会员 