数字信号处理实验Word格式文档下载.docx
《数字信号处理实验Word格式文档下载.docx》由会员分享，可在线阅读，更多相关《数字信号处理实验Word格式文档下载.docx（25页珍藏版）》请在冰点文库上搜索。
五、思考题
（1）在分析理想采样序列特性的实验中，采样频率不同时，相应理想采样序列的傅里叶变换频谱的数字频率度量是否都相同?
它们所对应的模拟频率是否相同?
为什么?
（2）在卷积定理验证的实验中，如果选用不同的频域采样点数M值，例如，选M=10和M=20，分别做序列的傅里叶变换，求得的结果有无差异？
为什么？
六、实验程序及结果
本实验使用自定义函数的方法产生信号：
门函数：
function[y,n]=gate（np,ns,nf）
ifns>
np|ns>
nf|np>
nf
error（'
Ê
ä
È
ë
Î
»
Ö
Ã
²
ý
Â
ú
×
ã
ns<
=np<
=nf'
）;
else
n=[ns:
nf];
y=（n>
=0&
n<
np）;
end
冲击函数：
function[y,n]=impseq（np,ns,nf）
y=[（n-np）==0];
时域函数：
function[y,ts]=sig（t0,tp,t1）;
a=50.*sqrt
（2）.*pi;
w0=50.*sqrt
（2）.*pi;
ts=t0:
tp:
t1;
y=A.*exp（-a.*ts）.*sin（w0.*ts）.*（ts>
=0）;
;
hb=impseq（0,0,10）+2.5.*impseq（1,0,10）+2.5.*impseq（2,0,10）+impseq（3,0,10）;
[ha,N]=gate（10,0,10）;
第一局部主程序：
[y1,t1]=sig（0,0.001,1）;
%t1=0:
0.001:
1;
%y1=0.5.*sin（2.*pi.*15.*t1）;
figure;
plot（t1,y1）;
axis（[00.15-10150]）;
title（'
f=1000Ê
±
Ó
ò
Ð
Å
º
'
f=-500:
500;
yy1=fftshift（fft（y1））;
plot（f,abs（yy1））;
f=1000Æ
µ
Æ
[y2,t2]=sig（0,1/300,1）;
plot（t2,y2）;
f=300Ê
yy2=fftshift（fft（y2））;
plot（-150:
150,abs（yy2））;
f=300Æ
[y3,t3]=sig（0,0.02,1）;
plot（t3,y3）;
f=200Ê
axis（[00.15-85]）;
yy3=fftshift（fft（y3））;
plot（-25:
25,abs（yy3））;
f=200Æ
第二局部主程序：
xb=impseq（0,0,100）;
stem（0:
length（xb）-1,xb）;
axis（[-1,3,0,2]）;
y1=fft（xb）;
plot（0:
length（y1）-1,y1）;
hb=impseq（0,0,100）+2.5.*impseq（1,0,100）+2.5.*impseq（2,0,100）+impseq（3,0,100）;
length（hb）-1,hb）;
axis（[-1,5,0,3]）
figure
y2=fftshift（fft（hb））;
plot（-（length（y2）-1）/2:
（length（y2）-1）/2,abs（y2））;
y3=conv（double（xb）,double（hb））;
length（y3）-1,y3）;
y4=fftshift（fft（y3））;
plot（-（length（y4）-1）/2:
（length（y4）-1）/2,abs（y4））;
xc=gate（10,0,100）;
ha=xc;
y5=fftshift（fft（conv（double（xc）,double（ha））））;
plot（-（length（y5）-1）/2:
（length（y5）-1）/2,abs（y5））;
实验结果：
第一局部：
采样频率为200Hz时时域恢复200Hz时频谱
采样频率为300Hz时时域恢复300Hz时频谱
采样频率为1000Hz是时域恢复1000Hz时频谱
第二局部：
Xb频域分析xb时域分析
hb频域分析hb时域分析
信号通过xb频域信号通过xc频域
卷积定理的验证：
%%%%%卷积定理的验证
y1=[0,1,2,3,4,5,6,7,8];
y2=[8,7,6,5,4,3,2,1,0];
y3=conv（y1,y2）;
f1=fft（y1,20）;
f2=fft（y2,20）;
f3=fft（y3,20）;
f=f1.*f2;
subplot（2,1,1）;
length（f）-1,f,'
.'
subplot（2,1,2）;
length（f3）-1,f3,'
时域卷积后求频谱和频域相乘
可见，f=conv（y1,y2）的频谱和y1,y2的频谱相乘后结果相同。
即满足卷积定理
七、答复以下问题。
（1）在分析理想采样序列特性的实验中，采样频率不同时，相应理想采样序列的傅里叶变换频谱的数字频率度量是否都相同？
它们所对应的模拟频率是否相同？
答：
数字频率度量不相同，但他们所对应的模拟频率相同。
由w=Ω*Ts公式得，采样间隔变化时模拟频率对应的数字频率会有相应的变化，故其度量会有所变化。
而且采样频率的大小直接关系到能否将能否将原始信号恢复出来。
（2）在卷积定理验证的实验中，如果选用不同的频域采样点数M值，例如，选M=10和M=20，分别做序列的傅里叶变换，求得的结果有无差异？
有差异，所到的结果点数不同。
八、总结。
本实验主要是后续实验的根底。
涉及内容也比拟浅显。
单位冲击序列与hb卷积后得到的频谱与hb原频谱相同。
原因很简单，是因为单位冲击序列卷积任何函数仍然是原函数。
实验二：
用FFT作谱分析
一、实验目的
（1）进一步加深DFT算法原理和根本性质的理解（因为FFT只是DFT的一种快速算法，所以FFT的运算结果必然满足DFT的根本性质）。
（2）熟悉FFT算法原理和FFT子程序的应用。
（3）学习用FFT对连续信号和时域离散信号进行谱分析的方法，了解可能出现的分析误差及其原因，以便在实际中正确应用FFT。
二、实验步骤
（1）复习DFT的定义、性质和用DFT作谱分析的有关内容。
（2）复习FFT算法原理与编程思想，并对照DIT-FFT运算流图和程序框图，读懂本实验提供的FFT子程序。
（3）编制信号产生子程序，产生以下典型信号供谱分析用：
（4）编写主程序。
以下图给出了主程序框图，供参考。
本实验提供FFT子程序和通用绘图子程序。
（5）按实验内容要求，上机实验，并写出实验报告。
三、实验内容
（1）对2中所给出的信号逐个进行谱分析。
（2）令x（n）=x4（n）+x5（n），用FFT计算8点和16点离散傅里叶变换，
X（k）=DFT［x（n）］
（3）令x（n）=x4（n）+jx5（n），重复
（2）。
四、思考题
（1）在N=8时，x2（n）和x3（n）的幅频特性会相同吗?
N=16呢?
（2）如果周期信号的周期预先不知道，如何用FFT进行谱分析?
五、实验程序及结果
x1=[1,1,1,1,0,0,0,0,0,0];
x2=[1,2,3,4,4,3,2,1,0,0];
x3=[4,3,2,1,1,2,3,4,0,0];
n=0:
1000;
x4=cos（n）;
x5=sin（n）;
t=0:
x6=cos（8*pi*t）+cos（16*pi*t）+cos（20*pi*t）;
%figure;
%title（'
ð
¸
ö
·
）
subplot（3,2,1）;
length（x1）-1,x1）;
subplot（3,2,2）;
length（x2）-1,x2）;
subplot（3,2,3）;
length（x3）-1,x3）;
subplot（3,2,4）;
20,x4（1:
21））;
subplot（3,2,5）;
20,x5（1:
subplot（3,2,6）;
20,x6（1:
200,x6（1:
201））;
y1=abs（fft（x1,8））;
y2=abs（fft（x2,8））;
y3=abs（fft（x3,8））;
y4=abs（fft（x4,8））;
y5=abs（fft（x5,8））;
y6=abs（fft（x6,8））;
length（y2）-1,y2）;
length（y5）-1,y5）;
length（y6）-1,y6）;
%%%figure;
%%%stem（0:
15,y6（1:
16））;
%%%%%%%%%%%%%%%
y1=abs（fft（x1,16））;
y2=abs（fft（x2,16））;
y3=abs（fft（x3,16））;
y4=abs（fft（x4,16））;
y5=abs（fft（x5,16））;
y6=abs（fft（x6,16））;
length（y6）-1,y4）;
length（y6）-1,y5）;
%figure;
%stem（0:
x7=x4+j*x5;
y71=abs（fft（x7,8））;
y72=abs（fft（x7,16））;
subplot（1,2,1）;
length（y71）-1,y71）;
subplot（1,2,2）;
length（y72）-1,y72）;
程序结果如下所示：
每个信号时域分析
因为第六个信号不易在上图中展示〔包含信息不到一个周期〕所以单独画出图形如下：
x6的时域波形
逐个进行8点离散傅里叶变换
逐个进行16点离散傅里叶变换
x（n）8点和16点离散傅里叶变换比拟
六、答复以下问题
（1）在N=8时，x2（n）和x3（n）的幅频特性会相同吗?
N=8时幅频特性一样，N=16时幅频特性不一样。
Fn所能分辨到频率为为Fs/N，如果采样频率Fs为1024Hz，采样点数为1024点，那么可以分辨到1Hz。
1024Hz的采样率采样1024点，刚好是1秒，也就是说，采样1秒时间的信号并做FFT，那么结果可以分析到1Hz，如果采样2秒时间的信号并做FFT，那么结果可以分析到。
如果要提高频率分辨力，那么必须增加采样点数，也即采样时间。
频率分辨率和采样时间是倒数关系。
总结：
采样点N的不同，分辨率不同，在N=8，俩函数不同的局部没有分辨出来，因此特性相同，而N=16时已经足以分辨出不同的频率，因此频谱不同。
（2）如果周期信号的周期预先不知道，如何用FFT进行谱分析?
设一个定长的m值,先取2m,看2m/m的误差是否大，如大的话再取4m,看4m/2m的误差是否大，如不大，4m（4倍的m值）那么可近似原来点的谱分析。
七、总结
本实验主要掌握fft函数的用法，包括采样点的意义。
这个意义已经在七、答复以下问题中描述，不再赘述。
实验三：
用窗函数法设计FIR数字滤波器
（1）掌握用窗函数法设计FIR数字滤波器的原理和方法。
（2）熟悉线性相位FIR数字滤波器特性。
（3）了解各种窗函数对滤波特性的影响。
二、实验内容及步骤
（1）复习用窗函数法设计FIR数字滤波器一节内容，阅读本实验原理，掌握设计步骤。
（2）编写程序。
①编写能产生矩型窗、升余弦窗、改良升余弦窗和二阶升余弦窗的窗函数子程序。
②编写主程序。
其中幅度特性要求用dB表示。
窗函数法设计滤波器主程序框图
三、思考题
（1）如果给定通带截止频率和阻带截止频率以及阻带最小衰减，如何用窗函数法设计线性相位低通滤波器?
写出设计步骤。
（2）如果要求用窗函数法设计带通滤波器，且给定上、下边带截止频率为ω1和ω2，试求理想带通的单位脉冲响应hd（n）。
四、实验程序及结果
myfreqz（b,a）函数，是自编的改良freqz（）的函数。
function[db,mag,pha,grd,w]=myfreqz（b,a）;
%advancedfreqz
%[db,mag,pha,grd,w]=myfreqz（b,a）;
[H,w]=freqz（b,a,1000,'
whole'
H=（H（1:
1:
501））'
w=（w（1:
mag=abs（H）;
db=20*log10（（mag+eps）/max（mag））;
pha=angle（H）;
grd=grpdelay（b,a,w）;
以下是主程序局部
wp=0.2*pi;
ws=0.3*pi;
deltaw=ws-wp;
N0=ceil（6.6*pi/deltaw）;
N=N0+mod（N0+1,2）;
wdham=（hamming（N））'
wc=（ws+wp）/2;
%hd=ideallp（wc,N）;
tao=（N-1）/2;
n=[0:
（N-1）];
m=n-tao+eps;
hd=sin（wc*m）./（pi*m）;
h=hd.*wdham;
[db,mag,pha,grd,w]=myfreqz（h,[1]）;
dw=2*pi/1000;
Rp=-（min（db（1:
wp/dw+1）））;
As=-round（max（db（ws/dw+1:
501）））;
subplot（2,2,1）;
length（hd）-1,hd,'
subplot（2,2,2）;
length（wdham）-1,wdham,'
subplot（2,2,3）;
length（h）-1,h,'
subplot（2,2,4）;
plot（1/length（db）:
1/length（db）:
1,db）;
%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%
wn1=boxcar（N）;
N-1,wn1,'
boxcar'
wn2=bartlett（N）;
N-1,wn2,'
bartlett'
wn3=hanning（N）;
N-1,wn3,'
hanning'
%wn4=hamming（N）;
wn4=blackman（N）;
N-1,wn4,'
blackman'
%wn=（N,beta）;
得到的结果如以下图：
其中，图一是理想脉冲响应，第一行图二是汉明窗；
图三是实际脉冲响应，图四是幅度响应。
一下是几个典型窗函数的图像。
，汉明窗上例已画出，不再画。
每个窗的名字已在每个子图上方标出。
典型窗函数
五、答复以下问题
（1）给定通带截止频率和阻带截止频率以及阻带最小衰减，用窗函数法设计线性相位低通滤波器的设计步骤：
技术指标Wp=0.2*pi,Ws=0.4*pi,Ap=0.25dB,As=50dB
方法一
选择海明窗
clearall;
Wp=0.2*pi;
Ws=0.4*pi;
tr_wide=Ws-Wp;
%过渡带宽度
N=ceil（6.6*pi/tr_wide）+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:
Wp/delta_w+1）））%实际通带纹波
As=-round（max（db（Ws/delta_w+1:
501）））%实际阻带纹波
subplot（221）
stem（n,hd）
理想单位脉冲响应hd（n）'
subplot（222）
stem（n,w_ham）
海明窗'
subplot（223）
stem（n,h）
实际单位脉冲响应hd（n）'
subplot（224）
plot（wi/pi,db）
幅度响应（dB）'
axis（[0,1,-100,10]）
方法二
Window=blackman（16）;
b=fir1（15,0.3*pi,'
low'
Window）;
freqz（b,128）
（2）用窗函数法设计带通滤波器，且给定上、下边带截止频率为ω1和ω2，理想带通的单位脉冲响应hd（n）的求解过程：
由hd（n）=1/2π∫w1woe-jwaejnwdw+1/2π∫w2w3e-jwaejnwdw
〔其中w0=-w0-wc,w1=-w0+wc,w2=w0-wc,w3=w0+wc〕
计算整理后可得：
hd〔n〕=2/（（n-a）*π）*sin[（n-a）wc]*cos[（n-a）w0]
=2wc/π*sa[（n-a）wc]*cos[（n-a）w0]
升级会员 