数字信号处理报告Word文档格式.docx
《数字信号处理报告Word文档格式.docx》由会员分享，可在线阅读，更多相关《数字信号处理报告Word文档格式.docx（26页珍藏版）》请在冰点文库上搜索。
%放入数据位
[x,fs,Nbits]=wavread（-g:
\Matlab\work\chaotian.wav）;
%把语音信号进行加载入Matlab仿真软件平台中
sound（x,fs）;
%对加载的语音信号进行回放
stem（x）;
title（‘语音信号的时域波形’）;
1.2语音信号的时域分析
因录音采样频率无法准确。
故用余弦序列表示信号的迟延，反转，上采样下采样。
迟延：
题目要求时延d为0.4s.取点数k为0.4*Fs。
产生一组新序列将原序列依次推后k点。
【四】程序清单与运行结果
程序清单：
y=load（'
E:
\erf1s1t0'
）;
fs=8192;
n=length（y）;
m=n+round（0.4*fs）;
q=（m-1）/fs;
l=round（0.4*fs）;
z=zeros（1,m）;
x=z;
fori=1:
nx（i）=y（i）;
end%原信号
nz（l+i）=y（i）;
end%信号的延拖
o=x+0.4*z;
%回声信号
sound（o）%声音
t=0:
1/fs:
q;
subplot（3,1,1）;
plot（t,x）
title（'
原信号'
subplot（3,1,2）;
plot（t,z）
时延信号'
subplot（3,1,3）;
plot（t,o）
回声信号'
运行结果：
翻转：
K为序列长度，将序列第i个点值依次与第k-i点值进行调换产生新序列，即为翻转序列。
采样：
下采样，也就是对信号的抽取。
其实，上采样和下采样都是对数字信号进行重采，重采的采样率与原来获得该数字信号的采样率比较，大于原信号的称为上采样，小于的则称为下采样。
上采样的实质也就是内插或插值。
%对上述语音信号进行时间反转x（-k）、上采样x（2k）、下采样x（k/2）操作
data=load（'
x=data;
%原声信号
d=0.4;
%d为时延
Fs=8192;
%Fs为采样频率
k=0.4*8192;
a=0.5;
%a为时延信号的衰减幅度
z=linspace（0,0,k）;
%linspace为线性间隔的向量
y=[zx'
];
x=[x'
z];
f=x+a*y;
f1=fliplr（x）;
%翻转信号
f2=downsample（x,2）;
%下采样信号
f3=upsample（x,2）;
%上采样信号
subplot（4,1,1）;
plot（x）;
原声信号'
）
subplot（4,1,2）;
plot（f1）;
翻转信号'
subplot（4,1,3）;
plot（f2）;
下采样'
subplot（4,1,4）;
plot（f3）;
上采样'
sound（f1）
设计二正余弦信号的谱分析
【一】设计目的
用DFT实现对正余弦信号的谱分析；
观察DFT长度和窗函数长度对频谱的影响；
对DFT进行谱分析中的误差现象获得感性认识。
【二】设计内容
1.对一个频率为10Hz，采样频率为64Hz的32点余弦序列进行谱分析，画出其频谱图；
若将频率改为11Hz，其他参数不变，重新画出该序列的频谱图，观察频率泄漏现象，分析原因。
2.考察DFT的长度对双频率信号频谱分析的影响。
设待分析的信号为
（1.2）
令两个长度为16的正余弦序列的数字频率为
及
。
取N为四个不同值16，32，64，128。
画出四个DFT幅频图，分析DFT长度对频谱分辨率的影响。
3.在上题中若把两个正弦波的频率取得较近，令
，
，试问怎样选择FFT参数才能在频谱分析中分辨出这两个分量？
【三】设计过程，调试，分析。
连续信号的频谱分析在工程上有着广泛的应用，计算机的第一步就是把连续信号离散化，一是采样，二是截断。
由此会产生频率混叠误差、栅栏效应和截断误差。
（一）对一个频率为10Hz，和11Hz，采样频率为64Hz的32点余弦序列用Matlab计算它的DFT程序为
N=32;
Fs=64;
f=10;
n=0:
N-1;
x=cos（2*pi*10*n/64）
y=cos（2*pi*11*n/64）
subplot（2,1,1）
X=abs（fft（x,N））
stem（n,X,'
.'
xlabel（'
n'
10HZ32点变化幅频曲线'
subplot（2,1,2）
Y=abs（fft（y,N））
stem（n,Y,'
11HZ32点变化幅频曲线'
通过图可看出：
频率为10Hz的余弦曲线DFT只有两个点不等于零，位于k=5和k=27处，k=5对应于频率10Hz，k=27对应于频率54Hz（也就是-10Hz）。
这样DFT确实正确的分辨了余弦信号的频率。
但是这样的理想结果是恰巧得到的，此时我们借去了五个完整的余弦周期（f*N/Fs=5）.
将频率改为11Hz，采样频率和窗长度依然为32点，计算图像可看出：
频谱图上k=5和k=27处都有较大的峰值，而其它的点上幅度不再为零。
这两个峰值对应的频率为10Hz和12Hz，所以，信号的峰值位于两者之间，本来是单一的11Hz频率的能量会分不到许多DFT频率上，这种现象叫频率泄露，来源于截断效应。
（二）对双信号的频谱分析程序为：
16点，32点，64点，128点的频谱分析程序：
N=16;
15;
f1=0.22;
f2=0.34;
x=0.5*sin（2*pi*f1*n）+sin（2*pi*f2*n）;
subplot（4,2,1）,stem（n,x,'
%xlable（'
ylable（'
x1（n）'
余弦序列'
X=abs（fft（x,N））;
subplot（4,2,2）
k=0:
stem（k,X,'
k'
ylable（'
X（k）'
string=[num2str（N）,'
点FFT幅频曲线'
title（string）;
subplot（4,2,3）,stem（n,x,'
subplot（4,2,4）
N=64;
subplot（4,2,5）,stem（n,x,'
subplot（4,2,6）
string=[num2str（64）,'
N=128;
subplot（4,2,7）,stem（n,x,'
subplot（4,2,8）
string=[num2str（128）,'
运行结果:
分析：
DFT样本值就是其DTFT在相应位置的采样。
在图一中很难看出两个峰值，因此要提高它的分辨率，故把N增大，逐渐可以看出它有两个峰值，将k换算成数字频率f=w/2*pi=k/N.这样可确定峰值的位置大体在f=0.21和0.35之附近，与信号的给定频率有一定的误差，这也是截断和泄露带来的问题，在这图上还可以看到一些较小的峰，这是很难判断是输入信号固有的，还是由泄露引起的。
说明了增加DFT长度N减小了相邻样本间的频率间距，提高频谱的视在分辨率，因而可以提高样本未知的测定精度。
（三）把两个正弦波的频率取得较近，另fr=[0.22,0.25]应怎样选择FFT参数
要能分清两个频率，分辨率至少应达到f=0.03.因为此处的数字频率是对采样频率Fs进行归一化后的，几fr最大为1，因此总的样本数至少要达到1/0.03=33
加窗以后可以使频谱函数更加光滑，便于分辨峰值位置和准确的数值，为了提高实际分辨率，应该尽量增加信号的长度n及DF长度N,当受到条件限制不能提高n，则单独提高N可以提高视在分辨率。
当
时，16点，32点，64点，128点的频谱运行结果：
（四）谱分析中的误差现象
①混叠误差。
产生混叠误差的原因是：
信号的离散化是通过抽样实现的，而抽样频率再高总是有限的。
除带限信号外，如果信号的最高频率Ω趋于无穷，则实际器件无法满足抽样定理。
而抽样过程如果不满足抽样定理，就会产生频谱的混叠，即混叠误差。
要减少或避免混叠误差，应该提高抽样频率，以设法满足抽样定理，或者采用抗混叠滤波这样的信号预处理措施。
②栅栏效应。
对于非周期信号来说，理论上应具有连续的频谱，但数字谱分析是用DFT来近似的，是用频谱的抽样值逼近连续频谱值，只能观察到有限个频谱值，每一个间隔中的频谱都观察不到。
如同通过“栅栏”观察景物一样，一部分被“栅栏”所阻挡，看不见，把这种现象称为栅栏效应。
连续时间信号只要采用数字谱分析的方法，就必定产生栅栏效应，栅栏效应只能减小而无法避免。
能够感受的频谱最小间隔值，称为频谱分辨率，而频率分辨率与信号截断长度成反比，且频率分辨率越小，插值效果越好。
因此即使连续信号是有限时间的，也应选择一个大的截断长度来改善栅栏效应。
③截断误差。
截断误差是由于对信号进行截断，把无限长的信号限定为有限长，即令有限区间以外的函数值均为零值的近似处理而产生的，这种处理相当于用一个矩形（窗）信号乘待分析的连续时间信号。
必然会引起吉布斯效应（波动），也会把窗函数的频谱引入信号频谱,造成混叠，所以需要考虑其误差的问题。
设计三数字滤波器的设计及实现
1.熟悉IIR数字滤波器和FIR数字滤波器的设计原理和方法；
2.学会调用MATLAB信号处理工具箱中的滤波器设计函数设计各种IIR和FIR数字滤波器，学会根据滤波要求确定滤波器指标参数；
3.掌握用IIR和FIR数字滤波器的MATLAB实现方法，并能绘制滤波器的幅频特性、相频特性；
4.通过观察滤波器的输入、输出信号的时域波形及其频谱，建立数字滤波的概念。
1.调用信号产生函数mstg产生由三路抑制载波调幅信号相加构成的复合信号st，观察st的时域波形和幅频特性曲线；
2.通过观察st的幅频特性曲线，分别确定可以分离st中三路抑制载波单频调幅信号的三个滤波器（低通滤波器、带通滤波器、高通滤波器）的通带截止频率和阻带截止频率，假定要求滤波器的通带最大衰减为0.1dB，阻带最小衰减为60dB；
3.编程调用MATLAB滤波器设计函数分别设计这三个数字滤波器，并绘图显示其幅频特性曲线；
4.用所设计的三个滤波器分别对复合信号st进行滤波，分离出st中的三路不同载波频率的调幅信号，并绘图显示滤波后信号的时域波形和频谱，观察分离效果。
（一）IIR数字滤波设计方法
数字滤波是数字信号处理的一种重要算法，广泛用于对信号的过滤、检测与参数的估计等信号处理中。
数字滤波器按其单位脉冲响应的长度可分为有限脉冲响应（FIR）滤波器和无限脉冲响应（IIR）滤波器两类。
数字滤波器按其通频带分为低通滤波器，高通滤波器，带通滤波器和带阻滤波器。
在本次实验中旨在用椭圆滤波器实现IIR低通，带通，高通滤波。
IIR滤波器基本概念
一个N阶递归型数字滤波器的差分方程
H（Z）=
MEbrz-r
r=0
N
1+K=Eakz-k
从以上的系统函数可知,设计IIR滤波器的任务就是通过计算寻求一个因果、物理上可实现的系统函数H（Z）,使其频率响应H（ejw）满足所希望得到的频域指标,即符合给定的通带截止频率、阻带截止、通带衰减和阻带衰减。
数字滤波器的设计过程
数字滤波器的一般设计过程可分为以下四个步骤:
1）按照实际的需要,确定滤波器的性能要求或指标;
2）用一个因果稳定的离散线性时不变系统,去逼近这一性能指标;
3）用有限精度的运算实现所设计的系统;
4）通过模拟,验证所设计的系统是否符合性能要求
（二）
1：
调用信号产生函数mstg产生由三路抑制载波调幅信号相加构成的复合信号st,画出st的时域波形和幅频特性，三路信号在时域混叠无法在时域分离，但在频域分离，所以可以通过滤波的方法在频域分离。
抑制载波单频调幅信号的数学表达式为
这三路调幅信号的载波频率分别为250Hz、500Hz和1000Hz
2：
通过观察st的幅频特性，分别可以确定三路信号的低通，带通，高通滤波器。
的通带截止频率和阻带截止频率，要求滤波器的通带最大衰减为0.1db,阻带最小衰减为60db
分离混合信号st的滤波器指标参数选取如下：
载波为250Hz的调幅信号，用低通滤波器分离。
指标为：
通带截止频率fp=300HZ,阻带截止频率fs=400Hz,通带最大衰减ap=0.1db,阻带最小衰减as=60db
对载波频率为500Hz的调幅信号，用带通滤波器分离，指标为：
通带截止频率fl=400HZ,fu=500Hz,阻带截止频率fl=600Hz,fu=700Hz,通带最大衰减ap=0.1db,阻带最小衰减as=60db
对载波频率为1000Hz的调幅信号，用高通滤波器分离，其指标为：
通带截止频率fp=900HZ,阻带截止频率fs=800Hz,通带最大衰减ap=0.1db,阻带最小衰减，60db
3：
用椭圆滤波器实现
椭圆滤波器特点：
幅值响应在通带和阻带内都是等波纹的。
对于给定的阶数和给定的波纹要求，椭圆滤波器能获得较其它滤波器更窄的过渡带宽，就这点而言，椭圆滤波器是最优的。
通带和阻带内波纹固定时，阶数越高，过渡带越窄
阶数固定，通带和阻带波纹越小，过渡带越宽
%1.调用信号产生函数mstg产生由三路抑制载波调幅信号相加构成的复合信号st，观察%st的时域
%波形和幅频特性曲线；
functionst=mstg
%产生信号序列st，并显示st的时域波形和频谱
%st=mstg返回三路调幅信号相加形成的混合信号，长度N=800
N=800;
%信号长度N=800
Fs=10000;
%采样频率Fs=10kHz,Tp为采样时间
T=1/Fs;
Tp=N*T;
T:
（N-1）*T;
f=k/Tp;
fc1=Fs/10;
%第一路调幅信号载波频率fc1=1000Hz
fm1=fc1/10;
%第一路调幅信号的调制信号频率fm1=100Hz
fc2=Fs/20;
%第二路调幅信号载波频率fc2=500Hz
fm2=fc2/10;
%第二路调幅信号的调制信号频率fm2=50Hz
fc3=Fs/40;
%第三路调幅信号载波频率fc3=250Hz
fm3=fc3/10;
%第三路调幅信号的调制信号频率fm3=25Hz
xt1=cos（2*pi*fm1*t）.*cos（2*pi*fc1*t）;
%产生第一路调幅信号
xt2=cos（2*pi*fm2*t）.*cos（2*pi*fc2*t）;
%产生第二路调幅信号
xt3=cos（2*pi*fm3*t）.*cos（2*pi*fc3*t）;
%产生第三路调幅信号
st=xt1+xt2+xt3;
%三路信号相加，得到复合信号
fxt=fft（st,N）;
%计算信号st的频谱
%以下为绘图命令
figure
（1）
plot（t,st）;
grid;
t/s'
ylabel（'
s（t）'
axis（[0,Tp,min（st）,max（st）]）;
（a）s（t）的波形'
stem（f,abs（fxt）/max（abs（fxt））,'
（b）s（t）的频谱'
axis（[0,Fs/8,0,1.2]）;
f/Hz'
幅度'
%低通滤波器
fp=320;
fs=350;
wp=2*fp/Fs;
ws=2*fs/Fs;
rp=0.1;
rs=60;
%DF指标（低通滤波器的通、阻带边界频）
[N,wp]=ellipord（wp,ws,rp,rs）;
%调用ellipord计算椭圆DF阶数N和通带截止频率wp
[B,A]=ellip（N,rp,rs,wp）;
%调用ellip计算椭圆带通DF系统函数系数向量B和A
y1t=filter（B,A,st）;
%滤波器软件实现
y1=fft（y1t）;
%绘图部分
figure
（2）;
subplot（3,3,1）;
[H,W]=freqz（B,A,800）;
plot（W*Fs/2/pi,abs（H））;
Hz'
H（w）'
低通滤波器'
axis（[0,2000,0,1.2]）;
subplot（3,3,2）;
plot（t,y1t）;
t'
y（t）'
分离出的250Hz的波形'
subplot（3,3,3）;
stem（f,abs（y1）/max（abs（y1））,'
|H|'
频谱'
axis（[0,1200,0,1]）;
%带通滤波器
fpl=420;
fpu=580;
fsl=300;
fsu=700;
wp=[2*fpl/Fs,2*fpu/Fs];
ws=[2*fsl/Fs,2*fsu/Fs];
y2t=filter（B,A,st）;
y2=fft（y2t）;
%绘图部分
subplot（3,3,4）;
带通滤波器'
subplot（3,3,5）;
plot（t,y2t）;
分离出的500Hz的波形'
subplot（3,3,6）;
stem（f,abs（y2）/max（abs（y2））,'
%高通滤波器
fp=870;
fs=780;
[B,A]=ellip（N,rp,rs,wp,'
high'
y3t=filter（B,A,st）;
y3=fft（y3t）;
subplot（3,3,7）;
高通滤波器'
升级会员 