数字信号处理实验报告完整版.docx
《数字信号处理实验报告完整版.docx》由会员分享，可在线阅读，更多相关《数字信号处理实验报告完整版.docx（47页珍藏版）》请在冰豆网上搜索。
数字信号处理实验报告完整版
实验1利用DFT分析信号频谱
、实验目的
1.加深对DFT原理的理解
2.应用DFT分析信号的频谱。
3.深刻理解利用DFT分析信号频谱的原理，分析实现过程中出现的现象及解决方法。
二、实验设备与环境
计算机、MATLA软件环境三、实验基础理论
1.DFT与DTFT的关系
有限长序列？
?
（?
?
）（g?
?
?
-1）的离散时间傅里叶变换？
?
（?
?
?
?
在频率区间
（0
?
<2?
?
的N个等间隔分布的点？
?
?
=2?
?
?
?
/?
?
（凑?
?
?
-1）上的N个取样
值可以由下式表示：
由上式可知，序列？
?
（?
?
）N点DFT?
?
?
?
实际上就是？
?
（?
序列的DTFT在N个等间隔频率点？
?
?
=2?
?
?
?
/?
?
（僕？
客?
?
-1）上样本?
?
?
?
。
2.利用DFT求DTFT
方法1:
由恢复出的方法如下：
由图2.1所示流程可知:
由上式可以得到:
N
X（ej）X（k）（
k1
其中Jil为内插函数
（x）
sin（N/2）
e
Nsin（/2）
方法2:
实际在MATLAB+算中，上述插值运算不见得是最好的办法。
由于
DFT是DTFT的取样值，其相邻两个频率样本点的间距为2n/N，所以如果我们增加数据的长度N,使得到的DFT谱线就更加精细，其包络就越接近DTFT的结果,这样就可以利用DFT计算DTFT如果没有更多的数据，可以通过补零来增加数据长度。
3.利用DFT分析连续信号的频谱
采用计算机分析连续时间信号的频谱，第一步就是把连续信号离散化，这里需要进行两个操作：
一是采样，二是截断。
对于连续时间非周期信号，按采样间隔T进行采样，阶段长度M那么：
M1
Xa（j）Xa（t）ejtdtTXa（nT）e
n0
对进行N点频域采样，得到
因此，可以将利用DFT分析连续非周期信号频谱的步骤归纳如下:
（1）确定时域米样间隔T，得到离散序列•廿
（2）确定截取长度M,得到M点离散序列这里・匚为窗函数。
（3）确定频域采样点数N,要求NA
（4）利用FFT计算离散序列的N点DFT得到広/k）.
（5）根据上式由-y计算采样点，二的近似值。
采用上述方法计算信号的频谱需要注意如下三个问题：
（1）频谱混叠。
如果不满足采样定理的条件，频谱会出现混叠误差。
对于频谱无限宽的信号，应考虑覆盖大部分主要频率分量的范围。
（2）栅栏效应和频谱分辨率。
使用DFT计算频谱，得到的结果只是N个频谱样本值，样本值之间的频谱是未知的，像通过一个栅栏观察频谱，称为“栅栏效应”。
频谱分辨率与记录长度成反比，要提高频谱分辨率，就要增加记录时间。
（3）频谱泄露。
对信号截断会把窗函数的频谱引入信号频谱，造成频谱泄露。
解决这个问题的主要办法是采用旁瓣小的窗函数，频谱泄露和窗函数均会引
因此，要合理选取采样间隔和截取长度，必要时还需考虑加适当的窗。
对于连续时间周期信号，我们在采用计算机进行计算时，也总是要进行截断，序列总是有限长的，仍然可以采用上述方法近似计算。
4.可能用到的MATLA函数与代码
实验中DFT运算可采用MATLA中提供的函数fft来实现。
DTFT可采用MATLA矩阵运算的方法进行计算，如下式所示：
jn1
jnN
四、实验内容
x（n）{2,1,1,1}
1、已知，完成如下要求：
（1）计算其DTFT并画出，区间的波形。
（2）计算4点DFT并把结果显示在
（1）所画的图形中。
（3）对x（n）补零，计算64点DFT并显示结果。
（4）根据实验结果，分析是否可以由DFT计算DTFT如果可以，如何实现
解：
（1）计算其DTFT并画出，区间的波形。
>>n=0:
3;
>>x=[2-111];
>>w=-pi:
0.01*pi:
pi;
>>X=x*exp（-j*n'*w）;
>>subplot（211）;
>>plot（w,abs（X））;
>>xlabel（'\Omega/\pi'）;
>>title（'Magnitude'）;
>>axistight;
>>subplot（212）;
>>plot（w,angle（X）/pi）;
>>xlabel（'\Omega/\pi'）;
>>title（'Phase'）;
>>axistight;
Magnitude
Phase
（2）计算4点DFT并把结果显示在
（1）所画的图形中
>>Xk=fft（x）;
>>subplot（211）;
>>holdon;
>>stem（n,abs（Xk）,'filled'）;
>>plot（w,abs（X））;
>>axistight;
>>xlabel（'\Omega/\pi'）;
>>title（'Magnitude'）;
>>subplot（212）;
>>holdon;
>>plot（w,angle（X）/pi）;
>>stem（n,angle（Xk）,'filled'）;
>>axistight;
>>xlabel（'\Omega/\pi'）;
>>title（'Phase'）;
运行结果如下：
Magnitude
3
2
1
0
/
Phase
/
（3）对x（n）补零，计算64点DFT并显示结果
>>x=[2-111zeros（1,60）];
>>n=0:
63;
>>Xk=fft（x）;
>>subplot（211）;
>>stem（n,abs（Xk）,'filled'）;
>>subplot（212）;
>>stem（n,angle（Xk）,'filled'）;
（4）
根据实验结果，分析是否可以由DFT计算DTFT如果可以，如何实现
的波形在一定的分辨率下已经相同
2、考察序列
x（n）=cos（0.48nn）+cos（0.52nn）
（1）0<=n<=10时，用DFT估计x（n）的频谱；将x（n）补零加长到长度为100点
序列用DFT估计x（n）的频谱，要求画出相应波形。
（2）0<=*=100时，用DFT估计x（n）的频谱。
并画出波形。
（3）根据实验结果，分析怎样提高频谱分辨率
解：
（1）0<=n<=10时，用DFT估计x（n）的频谱；将x（n）补零加长到长度为100
点序列用DFT估计x（n）的频谱，要求画出相应波形。
>>n=0:
10;
>>x=cos（0.48*pi.*n）+cos（0.52*pi.*n）;
>>Xk=fft（x）;
>>subplot（211）;
>>stem（n,Xk,'filled'）;
>>x=[cos（0.48*pi.*n）+cos（0.52*pi.*n）zeros（1,89）];
>>Xk=fft（x）;
>>subplot（212）;
>>stem（Xk,'filled'）;
10
-10[1[[11[[[
0102030405060708090100
（2）0<=*=100时，用DFT估计x（n）的频谱。
并画出波形
程序代码如下：
>>n=0:
100;
>>x=cos（0.48*pi.*n）+cos（0.52*pi.*n）;
>>Xk=fft（x）;
>>stem（Xk,'filled'）;
>>axistight;
30
20
10
0
-10
-20
-30
102030405060708090100
（3）根据实验结果，分析怎样提高频谱分辨率可以通过如下三种方式来增加分辨率。
a、增加时域内信号采样时间
b、提高采样频率
c、补零
3、已知信号x（t）=0.15sin（2nf1t）+sin（2nf2t）-0.1sin（2nfst）,其
f1=1Hz,f2=2Hz,f3=3Hz。
从x（t）的表达式可以看出，它包含三个频率的正弦波，但是，从其时域波形来看，似乎是一个正弦信号，利用DFT做频谱分析，确定
适合的参数，使得到的频谱的频率分辨率符合需要。
n=0:
10;
x=0.15*sin（2*pi.*n）+sin（4*pi.*n）-0.1*sin（6*pi.*n）;
Xk=fft（x）;
stem（abs（Xk）,'filled'）;
n=0:
100;
x=0.15*sin（0.2*pi.*in）+sin（0.4*pi.*in）-0.1*sin（0.6*pi.*n）;
Xk=fft（x）;
stem（abs（Xk）,'filled'）;
n=0:
200;
x=0.15*sin（0.1*pi.*n）+sin（0.2*pi.*n）-0.1*sin（0.3*pi.*n）;
Xk=fft（x）;
stem（abs（Xk）,'filled'）;
结果分析：
上图为x（t）信号截取过后的连续时间信号的傅里叶变换幅频特性曲
线，截取周期为100S（即为采样时间M）,根据X（j）MCk的特点，知1hz
处的幅值为X（j2*1）1000.157.5，X（j2*2）5，X（j2*3）50与图
2
像相符。
4、利用DFT近似分析连续时间信号x（t）=e-0.1u（t）的频谱（幅度值）。
分析采用不同的采样间隔和截取长度进行计算的结果，并最终确定合适的参数。
n=0:
10;
x=exp（-0.1.*n）;
Xk=fft（x）;
stem（Xk,'filled'）;
n=0:
20;
x=exp（-0.05.*n）;
Xk=fft（x）;
stem（Xk,'filled'）;
14
12
10
8
6
斗
2
0510152025
n=0:
40;
x=exp（-0.025.*n）;
Xk=fft（x）;
stem（Xk,'filled'）;
五、心得体会
通过本次实验，加深了对DFT原理的理解。
学会了应用DFT分析信号的频谱。
深刻理解到利用DFT分析信号频谱的原理，能够分析实现过程中出现的现象及解决方法。
实验二利用FFT计算线性卷积
一、实验目的
1.掌握利用FFT计算线性卷积的原理及具体实现方法。
2.加深理解重叠相加法和重叠保留法。
3.考察利用FFT计算线性卷积各种方法的适用范围。
二、实验基础理论
1.线性卷积与圆周卷积
设x（n）为L点序列，h（n）为M点序列，x（n）和h（n）的线性卷积为
OO
yi（n）=x（n）?
h（n）=刀x（m）h（n-m）
m=-s
yi（n）的长度为L+M-1。
x（n）和h（n）的圆周卷积为
N-1
y（n）=x（n）Oh（n）=刀x（m）h（n-m）NRN（n）
m=0
圆周卷积与线性卷积相等而不产生交叠的必要条件为
N>L+M-1
圆周卷积定理：
根据DFT性质，x（n）和h（n）的N点圆周卷积的DFT等于它们的
DFT的乘积：
DFT[x（n）Oh（n）]=X（k）H（k）
2.快速卷积
快速卷积发运用圆周卷积实现线性卷积，根据圆周卷积定理利用FFT算法实现圆
周卷积。
可将快速卷积运算的步骤归纳如下：
（1）必须选择N>L+M-1;为了能使用基-2算法，要求N=2Y。
采用补零的办法使得x（n）和h（n）的长度均为N。
⑵计算x（n）和h（n）的N点FFT。
FFT
x（n）—X（k）
FFT
h（n）—H（k）
（3）组成乘积
Y（k）=X（k）H（k）
⑶利用IFFT计算丫（k）的IDFT,得到线性卷积y（n）
IFFT
Y（k）-y（n）
3.分段卷积
我们考察单位取样响应为h（n）的线性系统,输入为x（n）,输出为y（n）,则y（n）=x（n）?
h（n）
如果x（n）极长时,如果要等x（n）全部集齐时再开始进行卷积,会使输出有较大延时；如果序列太长,需要大量存储单元。
为此,我们把x（n）分段,为别求出每段的卷积,合在一起得到最后的总输出。
这称为分段卷积。
分段卷积可以细分为重叠保留法和重叠相加法。
重叠保留法：
设x（n）的长度为Nx,h（n）的长度为M把序列x（n）分成多段N点序列Xi（n），每段雨前一段重写M-1个样本。
并在第一个输入段前面补M-1个零。
计算每一段与h（n）的圆周卷积,其结果中前M-1个不等与线性卷积,应当舍去,只保留后面N-M+1个正确的输出样本，把它们合起来得到总的输出。
利用FFT实现重叠保留法的步骤如下：
（1）在x（n）前面填充M-1个零，扩大以后的序列为
x?
（n）={0,0,?
0,x（n）}
⑵将x（n）分为若干段N点子段，设L=N-M+1为每一段的有效长度，则第i段的数据为：
xi（n）=?
（m）iL0,0⑶计算每一段与h（n）的N点圆周卷积，利用FFT计算圆周卷积
FFT
xi（n）-Xi（k）
FFT
h（n）-H（k）
FFT
Y（k）-Xi（k）H（k）
IFFT
Y（k）-y（n）
（4）设每一段卷积结果的前M-1个样本，连接剩下的样本得到卷积结果y（n）
重叠相加法：
设
没断信号，则：
h（n）长度为M将信号x（n）分解成长为L的子段。
以xi（n）表示
OO
x（n）=刀Xi（n）
i=0
x（n+iL）,x（n）0xi（n）={0其他
x（n）?
h（n）=刀Xi（n）?
h（n）
i=0
每一段卷积yi（n）的长度为L+M-1,所以在做求和时，相邻两段序列由M-1个样本重叠，即前一段的最后M-1个样本和下一段前M-1个样本序列重叠，这个重叠部分相加，再与不重叠的部分共同组成y（n）。
利用FFT实现重叠保留法的步骤如下：
（1）将x（n）分为若干L点子段Xi（n）。
（2）计算每一段与h（n）的卷积，根据快速卷积法利用FFT计算卷积。
（3）将各段相加，得到输出y（n）。
OO
y（n）=刀yi（n-iL）
i=0
三、实验内容
假设要计算序列x（n）u（n）u（nL）,0nL和h（n）cos（0.2n）,0nM的线
性卷积，完成以下实验内容：
设L=M根据线性卷积的表达式和快速卷积的原理，分别编程实现计算两个序列线性卷积的方法，比较当序列长度分别为8，16，32,64,256,512,1024时，两种
计算方法计算线性卷积所需时间。
>>fori=1:
7
L=input（'L:
'）;n=0:
L;
x=heaviside（n）-heaviside（n-L）;h=cos（0.2*pi.*n）;
tic
y=conv（x,h）;
tocendL:
8时间已过0.000050秒。
L:
16
时间已过0.000037秒。
L:
32
时间已过0.000056秒。
L:
64
时间已过0.000065秒。
L:
256
时间已过0.000101秒。
L:
512
时间已过0.000188秒。
L:
1024
时间已过0.000254秒。
>>fori=1:
7
L=input（'L:
'）;n=0:
L;x=heaviside（n）-heaviside（n-L）;h=cos（0.2*pi.*n）;
tic
Xk=fft（x,L+1）;
Hk=fft（h,L+1）;
Yk=Xk.*Hk;
y=ifft（Yk）;
toc
end
L:
8时间已过0.000048秒。
L:
16
时间已过0.000044秒。
L:
32
时间已过0.000042秒。
L:
64时间已过0.000055秒。
L:
256
时间已过0.000134秒。
L:
512时间已过0.000147秒。
L:
1024时间已过0.000176秒。
当L=2048且M=256时，比较计算线性卷积和快速卷积所需的时间，进一步考察当L=4096且M=256时两种算法所需时间。
>>fori=1:
2
L=input（'L:
'）;
M=input（'M:
'）;
n1=0:
L;
n2=0:
M;
x=heaviside（n1）-heaviside（n1-L）;
h=cos（0.2*pi.*n2）;
tic
Xk=fft（x,L+1）;Hk=fft（h,L+1）;
Yk=Xk.*Hk;y=ifft（Yk）;
tocendL:
2048M:
256时间已过0.001139秒。
L:
4096M:
256时间已过0.001503秒。
>>fori=1:
2
L=input（'L:
'）;M=input（'M:
'）;
n1=0:
L;n2=0:
M;
x=heaviside（n1）-heaviside（n1-L）;h=cos（0.2*pi.*n2）;
tic
y=conv（x,h）;
tocendL:
2048M:
256时间已过0.000263秒。
L:
4096
M:
256时间已过0.000368秒。
3.
>>L=input（'L:
'）;M=input（'M:
'）;n2=0:
M;k=L/M;
h=cos（0.2*pi.*n2）;
Y=0;
ticfori=1:
k
n1=（i-1）*M:
i*M-1;x=heaviside（n1）-heaviside（n1-L）;
Xk=fft（x,M）;Hk=fft（h,M）;
Yk=Xk.*Hk;y=ifft（Yk,M）;
Y=Y+y;endtocL:
2048M:
256时间已过0.025870秒。
>>L=input（'L:
'）;M=input（'M:
'）;n2=0:
M;k=L/M;
h=cos（0.2*pi.*n2）;
Y=0;ticfori=1:
k
n1=（i-1）*M:
i*M-1;x=heaviside（n1）-heaviside（n1-L）;
Xk=fft（x,M）;Hk=fft（h,M）;
Yk=Xk.*Hk;y=ifft（Yk,M）;
Y=Y+y;endtocL:
4096M:
256时间已过0.022300秒。
L=4096;n1=0:
L;x=heaviside（n1）-heaviside（n1-L）;n2=0:
256;
h=cos（0.2*pi.*n2）;tic
y1=conv（x,h）;tocticXk=fft（x）;Hk=fft（h）;Yk=Xk.*Hk;
y2=ifft（Yk）;toc
L=2048且M=256时计
编程实现利用重叠相加法计算两个序列的线性卷积，
算线性卷积的时间，与第二题结果进行比较。
clear
tic
N=512;
m=0:
256;
h=cos（0.2*pi*m）;
n=0:
2048;
x=heaviside（n）-heaviside（n-2048）;
Lenx=length（x）;
M=length（h）;
M1=M-1;
L=N-M1;
h=fft（h,N）;
K=ceil（Lenx/L）;
fori=Lenx:
K*L-1
x（i+1）=0;
end
Y=zeros（K,N）;
YY=zeros（1,（K-1）*L+N）;
fork=0:
K-1
xk=[x（k*L+1:
k*L+L）,zeros（1,M1）];
Y（k+1,:
）=（ifft（fft（xk）.*h））;
YY（k*L+1:
k*L+N）=YY（k*L+1:
k*L+N）+Y（k+1,:
）;
endtoc
时间已过0.028816秒
L=2048且M=256时计
编程实现利用重叠保留法计算两个序列的线性卷积，算线性卷积的时间，与第二题结果进行比较。
clc
clear
tic
N=512;m=0:
256;
h=cos（0.2*pi*m）;n=0:
2048;
x=heaviside（n）-heaviside（n-2048）;Lenx=length（x）;
M=length（h）;
M1=M-1;
L=N-M1;h=fft（h,N）;
K=floor（（Lenx+M1-1）/L）+1;p=（K）*L-Lenx;
x1=[zeros（1,M1）,x,zeros（1,p）];Y=zeros（K,N）;
fork=0:
K-1xk=fft（x1（k*L+1:
k*L+N））;
Y（k+1,:
）=（ifft（xk.*h））;
end
Z=reshape（Y（:
M:
N）',1,[]）;
toc
时间已过0.044424秒。
四、实验心得
加深了理
通过本次实验，掌握了利用FFT计算线性卷积的原理及具体实现方法,解重叠相加法和重叠保留法。
实验3IIR数字滤波器设计
一、实验目的
掌握利用脉冲响应不变法设计IIR数字滤波器的原理及具体方法。
加深理解数字滤波器和模拟滤波器之间的技术指标转化。
掌握脉冲响应不变法设计IIR数字滤波器的优缺点及适用范围。
二、实验内容
设采样频率为fs10kHz,设计数字低通滤波器，满足如下指标
通带截止频率：
fp1kHz,通带波动：
Rp1dB:
阻带截止频率：
fs1.5kHz，阻带衰减：
R15dB:
要求分别采用巴特沃斯、切比雪夫I型、切比雪夫U型和椭圆模拟原型滤波器并分别结合脉冲响应不变法和双线性变换法进行设计。
结合实验结果，分别讨论采用上述方法设计的数字滤波器是否都能满足给定指标要求，分析脉冲响应不变法和双线性变换法设计IIR数字滤波器的优缺点及适用范围。
1巴特沃斯
脉冲响应不变法：
fp=1000
fs=1500
f=10000
Wp=2*pi*fp
Ws=2*pi*fs
wp=2*f*tan（2*pi*fp/（2*f））
ws=2*f*tan（2*pi*fs/（2*f））
Rp=1
As=15
N仁ceil（（log10（（10A（Rp/10）-1”（10A（As/10）-1）））/（2*log10（Wp/Ws）））
N2=ceil（（log10（（10A（Rp/10）-1”（10A（As/10）-1）））/（2*log10（wp/ws）））
Omegac仁Wp/（（10A（Rp/10）-1）A（1/（2*N1）））
Omegac2=wp/（（10A（Rp/10）-1F（1/（2*N1）））[z,p,k]=buttap（N1）
b1=k*Omegac1AN1a1=poly（p*Omegac1）
[H,w]=freqs（b1,a1）subplot（221）plot（w/pi,abs（H））gridonsubplot（223）
plot（w/pi,ang
升级会员 