数字信号处理课程设计报告25770.docx
《数字信号处理课程设计报告25770.docx》由会员分享,可在线阅读,更多相关《数字信号处理课程设计报告25770.docx(30页珍藏版)》请在冰豆网上搜索。
数字信号处理课程设计报告25770
课程设计报告书
数字信号处理课程设计报告书
基于Matlab对语音信号进行频谱分析及滤波
学院:
专业:
设计者:
学号:
设计周数:
完成日期:
[作者]钟伟雄
[机构]广东工业大学在校本科生
摘要:
数字信号处理是通信工程专业的一门相当重要的学科,对日后就业和科研有重大意义,通过MATLAB,我们可以清晰地理解数字信号处理中难以理解的一面,对理论的知识加以深化。
关键字:
MATLAB数字信号处理GUI频谱相位滤波器
一、课程设计题目
应用Matlab对语音信号进行频谱分析及滤波
二、课程设计目的
数字信号处理是一门以算法为核心,理论和实践性较强的学科。
是电子信息工程、通信工程专业、电子信息科学与技术专业的一门重要的专业技术基础课。
数字信号处理课程是在学习完数字信号处理的相关理论后,进行的综合性训练课程,其目的是:
1.使学生进一步巩固数字信号处理的基本概念、理论、分析方法和实现方法;
2.增强学生应用Matlab语言编写数字信号处理的应用程序及分析、解决实际问题的能力;
三、课程设计内容描述和要求
为了巩固所学的数字信号处理理论知识,使学生对信号的采集、处理、传输、显示和存储等有一个系统的掌握和理解,安排了以下的课程设计的内容:
录制一段个人自己的语音信号,并对录制的信号进行采样;画出采样后语音信号的时域波形和频谱图;给定滤波器的性能指标,采用窗函数法和双线性变换设计滤波器,并画出滤波器的频率响应;然后用自己设计的滤波器对采集的信号进行滤波,画出滤波后信号的时域波形和频谱,并对滤波前后的信号进行对比,分析信号的变化;回放语音信号;最后,设计一个信号处理系统界面。
下面对各步骤加以具体说明。
2.1语音信号的采集
利用Windows下的录音机,录制一段自己的话音,时间在1s内。
然后基于Matlab软件平台下,利用函数[x,fs,bits]=wavread(filename)对语音信号进行采样,记住采样频率和采样点数。
通过wavread函数的使用,理解采样频率、采样位数等概念。
2.2语音信号的频谱分析
首先播放读入的语音信号,画出语音信号的时域波形;然后调用FFT函数对语音号进行快速傅里叶变换,得到信号的频谱特性,调用plot,从而加深对频谱特性的理解。
2.3设计数字滤波器和画出其频率响应
给出各滤波器的性能指标:
(1)低通滤波器性能指标fb=1000Hz,fc=1200Hz,As=100dB,Ap=1dB。
(2)高通滤波器性能指标fc=4800Hz,fb=5000HzAs=100dB,Ap=1dB。
(3)带通滤波器性能指标fb1=1200Hz,fb2=3000Hz,fc1=1000Hz,fc2=3200Hz,As=100dB,Ap=1dB。
要求用窗函数法和双线性变换法设计上面要求的3种滤波器。
在Matlab中,可以利用函数fir1设计FIR滤波器,可以利用函数butte,cheby1和ellip设计IIR滤波器;利用Matlab中的函数freqz画出各滤波器的频率响应。
2.4用滤波器对信号进行滤波
要求用自己设计的各滤波器分别对采集的信号进行滤波,在Matlab中,FIR滤波器利用函数fftfilt对信号进行滤波,IIR滤波器利用函数filter对信号进行滤波。
2.5比较滤波前后语音信号的波形及频谱
要求在一个窗口同时画出滤波前后的波形及频谱。
2.6回放语音信号
在Matlab中,调用函数sound对声音进行回放。
其调用格式:
sound(x,fs,bits);可以感觉滤波前后的声音有变化。
2.7设计系统界面
为了使编制的程序操作方便,设计处理系统的用户界面。
在所设计的系统界面上可以选择滤波器的类型,滤波方式,输入滤波器的参数,显示滤波器的频率响应,选择音频信号等。
这里选用matlab的GUI工具进行界面设计。
四、课程设计进度安排
序号
设计内容
所用时间
1
熟悉Matlab程序设计方法,了解数字信号处理工具箱使用
1天
2
分析题目,设计程序框图,编写程序代码
1天
3
上机调试程序,修改并完善设计,并完成设计报告
1天
合计
3天
五、课程设计基本流程
5.1熟悉MATLAB基本操作,基本语法,以及GUI界面设计中各控件的调用方式。
MATLAB是一款功能超强的数学软件,应用于各个行业。
而其基本操作却很大众化,操作起来很人性化。
其中所附带的GUI面向对象的友好编程方式更便于学术交流和设计开发。
5.2熟悉各种滤波器设计方案以及和收集与滤波器相关的MATLAB函数及其调用格式。
一.IIR数字滤波器
(一)IIR数字滤波器的传递函数及特点
设IIR滤波器的输入序列为x(n),则IIR滤波器的输入序列x(n)与输出序列y(n)之间的关系可以用下面的方程式表示:
其中,
和
是滤波器的系数,其中
中至少有一个非零。
与之相对应的差分方程为:
由传递函数可以发现无限常单位冲激响应滤波器有如下特点:
a.单位冲激响应h(n)是无限长的。
b.系统传递函数H(z)在有限z平面上有极点存在。
c.结构上存在着输出到输入的反馈,也就是结构上是递归型的。
(二)IIR数字滤波器的设计与实现
IIR数字滤波器的设计有多种方法,如频率变换法、数字域直接设计以及计算辅助设计等。
下面只介绍频率变换设计法。
首先考虑由模拟低通滤波器到数字低通滤波器的转换,其基本的设计过程如下:
A.将数字滤波器的技术指标转换为模拟滤波器的技术指标;
B.设计模拟滤波器G(S);
C.将G(S)转换成数字滤波器H(Z);
在低通滤波器的设计基础上,可以得到数字高通、带通、带阻滤波器的设计流程如下:
A.给定数字滤波器的设计要求(高通、带阻、带通);
B.转换为模拟(高通、带阻、带通)滤波器的技术指标;
C.转换为模拟低通滤波器的指标;
D.设计得到满足第三步要求的低通滤波器传递函数;
E.通过频率转换得到模拟(高通、带阻、带通)滤波器;
F.变换为数字(高通、带阻、带通)滤波器。
(三)在matlab中设计IIR滤波器的方法及其它们所用到的函数如表所示。
表1.基于matlab中设计IIR滤波器的方法列表
方法
描述
函数
模拟原型法
采用经典低通滤波器作为
连续域上的设计模型,通
过频率变换得到IIR数字
滤波器,最后进行离散化处理
完整设计函数:
Beself,butter,cheby1,cheby2,ellip
滤波器的阶估计函数:
Buttord,cheb1ord,cheb2ord,ellipord
低通模拟滤波器原型函数:
beselap,buttap,cheb1ap,cheb2ap,ellipap
频域变换函数:
Lp2bp,lp2bs,lp2hp,lp2lp
其他函数:
Bilinear,impinvar
直接设计方法
直接在离散时域上估计线
性的幅度响应
yulewalk
通用butterworth
设计方法
使用butterworth设计低通
数字滤波器
Maxflat
参数建模方法
寻找接近于所需要设计的
滤波器的通用原型
时域上的建模函数:
Lpc,prony,stmcb
频域上的建模函数:
Invfreqs,invfreqz
表2频率转换函数列表
频率转换
转换函数
低通到低通
[numt,dent]=lp2lp(num.den,w0)
[At,Bt,Ct,Dt]=lp2lp(A,B,C,D,w0)
低通到高通
[numt,dent]=lp2hp(num.den,w0)
[At,Bt,Ct,Dt]=lp2hp(A,B,C,D,w0)
低通到带通
[numt,dent]=lp2bp(num.den,w0)
[At,Bt,Ct,Dt]=lp2bp(A,B,C,D,w0)
低通到带阻
[numt,dent]=lp2bs(num.den,w0)
[At,Bt,Ct,Dt]=lp2bs(A,B,C,D,w0)
(四)标准数字滤波器设计函数
Matlab提供了一组标准的数字滤波器设计函数,大大简化了滤波器的设计过程。
butter
功能:
Butterworth模拟/数字滤波器设计
格式:
[b,a]=butter(n,wn,'ftype',’s’)
[b,a]=butter(n,wn,'ftype')
说明:
选项中加入‘S’用于设计各种模拟Butterworth滤波器;不加设计各种数字Butterworth滤波器
Ftype为缺省,设计低通滤波器
Ftype=hign,设计高通滤波器
Ftype=stop,设计带阻滤波器
【实例1】设计一个5阶Butterworth数字高通滤波器,阻带截止频率为250Hz。
设采样频率为1000Hz。
源代码如下:
[b,a]=butter(5,250/500,'high')
[z,p,k]=butter(5,250/500,'high')
freqz(b,a,512,1000)
程序运行后,产生结果如下所示。
b=
0.0528-0.26390.5279-0.52790.2639-0.0528
a=
1.0000-0.00000.6334-0.00000.0557-0.0000
z=
11111
p=
0.0000+0.7265i0.0000-0.7265i0.0000+0.3249i
0.0000-0.3249i0.0000
k=0.0528
图15阶Butterworth数字高通滤波器
Cheby1、Cheby2
功能:
chebyshevI、chebyshevII型模拟/数字滤波器设计
格式:
[b,a]=cheby1(n,Rp,wn,'ftype',)
[b,a]=cheby2(n,Rs,wn,'ftype')
【实例2】设计一个7阶chebyshevII型数字低通滤波器,截止频率为3000Hz,
Rs=30dB。
设采样频率为1000Hz。
解:
源程序如下:
[b,a]=cheby2(7,30,300/500');
[z,p,k]=butter(5,250/500,'high');
freqz(b,a,512,1000)
程序运行后,产生如图6-2所示的波形。
图27阶chebyshevII型数字低通滤波器
(五)冲激响应不变法
一般来说,在要求时域冲激响应能模仿模拟滤波器的场合,一般使用冲激响应不变法。
冲激响应不变法一个重要特点是频率坐标的变换是线性的,因此如果模拟滤波器的频响带限于折叠频率的话,则通过变换后滤波器的频响应可不失真的反映原响应与频率的关系。
【实例3】设计一个中心频率为500HZ,带宽为600Hz的数字带通滤波器,采样频率为1000Hz。
解:
源代码如下:
[z,p,k]=buttap(3);
[b,a]=zp2tf(z,p,k);
[bt,at]=lp2bp(b,a,500*2*pi,600*2*pi);
[bz,az]=impinvar(bt,at,1000);%将模拟滤波器变换成数字滤波器
freqz(bz,az,512,'whole',1000)
程序运行后,产生如图6-3所示的波形。
图3用冲激响应不变法设计数字低通滤波器
(六)双线性变换法
与冲激响应不变法比较,双线性变换的主要优点是靠频率的非线性关系得到S平面与Z平面的单值一一对应关系,整个值对应于单位圆一周。
所以从模拟传递函数可直接通过代数置换得到数字滤波器的传递函数。
【实例4】设计一个截止频率为200Hz的数字低通滤波器,采用频率为1000Hz。
解:
源代码如下:
[z,p,k]=buttap(3);
[b,a]=zp2tf(z,p,k);
[bt,at]=lp2lp(b,a,200*2*pi);
[bz,az]=bilinear(bt,at,1000);
freqz(bz,az,512,1000)
结果如图6-4所示:
bz=
0.07530.22590.22590.0753
az=
1.0000-0.82660.5154-0.0865
图4用双线性变换法设计数字低通滤波器
【实例5】基于Butterworth模拟滤波器原型,使用双线性转换设计数字滤波器,其中参数指标为:
通带截止频率:
通带波动值:
阻带截止频率:
阻带波动值:
解:
首先确定滤波器的阶数N,同时根据
确定
=0.5。
接着使用bilinear进行双线性转换,最后绘制在频域上的各种图像,其源代码如下:
%数字滤波器指标
wp=0.2*pi;
ws=0.3*pi;
Rp=1;
As=15;
%将数字滤波器指标反转变化为模拟滤波器的参数
T=1;
fs=1/T;
omegap=(2/T)*tan(wp/T);
omegas=(2/T)*tan(ws/T);
ep=sqrt(10^(Rp/10)-1);
Ripple=sqrt(1/(1+ep*ep));
Attn=1/(10^(As/20));
%butterworth原型模拟滤波器的设计
[cs,ds]=afd_butt(omegap,omegas,Rp,As);
%双线性变换
[b,a]=bilinear(cs,ds,T,fs);
%频域图像的绘制
freqz(b,a);
程序运行后,产生4阶的butterworth数字滤波器,频率响应如图5所示的波形。
图54阶Butterworth数字滤波器
二.FIR低通数字滤波器
(一)设计原理
如果系统的冲激响应
为已知,则系统的输入/输出关系为:
对于低通滤波器,只要设计出低通滤波器的冲激响应函数,就可以由上式得到系统的输出了。
假设所希望的数字滤波器的频率响应为
,它是频域的周期函数,周期为2
,那么它与
相对应的傅立叶系数为
以
为冲激响应的数字滤波器将具有频域响
。
但是将
作为滤波器脉冲响应有两个问题:
(1)它是无限长的,与FIP滤波器脉冲响应有限长这一前提不一致
(2)它是非因果的,
对此,要采取以下的措施,
(1)将
截短
(2)将其往右平移,
由此得到
的实际频域响应
,与理想频域响应
相近,但不完全一致。
理论证明上述现象是对
进行简单截短处理的必然结果,一般称为吉布斯现象,为尽可能的减少吉布斯现象,应对
进行加窗截取,即以
作为FIR滤波器的系数。
常用的窗函数有矩形窗、海明窗和布莱克曼窗等。
(二)用窗函数法设计FIR滤波器
Matlab设计FIR滤波器有多种方法和对应的函数,见表7-1。
表7-1matlab设计FIR滤波器的方法和函数
方法
描述
函数
窗方法
使用窗函数和逆傅立叶变换实现
fir1,fir2,kaiserord等
多带方法
包含子带频率域
firls,remez等
最小二乘法
使用最小二乘法将整个频率域上的错误几率压缩到最小
fircls,fircls1等
任意响应法
使用任意响应,包括非线性相位以及复滤波器
cremez等
余弦法
使用三角函数的低通响应
firrcos等
窗函数方法不仅在数字滤波器的设计中占有重要的地位,同时可以用于功率谱的估计,从根本上讲,使用窗函数的目的就是消除由无限序列的截短而引起的Gibbs现象所带来的影响。
窗函数设计线性相位FIR滤波器步骤如下:
(1)确定数字滤波器的性能要求,临界频率
,滤波器单位脉冲响应长度N
(2)根据性能要求,合理选择单位脉冲响应h(n)的奇偶对称性,从而确定理想频率响应
的幅频特性和相频特性
(3)求理想单位脉冲响应
,在实际计算中,可对
采样,并对其求IDFT的
,用
代替
(4)选择适当的窗函数w(n),根据
求所需设计的FIR滤波器单位脉冲响应
(5)求
,分析其幅频特性,若不满足要求,可适当改变窗函数形式或长度N,重复上述设计过程,以得到满意的结果
【实例1】设计一个34阶的高通滤波器,截止频率为0.48,使用具有30dB波纹的chebyshev窗。
解:
源程序如下:
b=fir1(34,0.48,'high',chebwin(35,30));
freqz(b,1,512)
其响应波形如图6所示。
图6带通FIR滤波器
【实例2】设计具有下面指标的低通FIR滤波器
由于其最小阻带衰减为50dB,因此可以选择hamming窗来实现这个滤波器,因为它具有较小的过渡带。
解:
MATLAB源程序为
%数字滤波器指标
wp=0.2*pi;
ws=0.3*pi;
tr_width=ws-wp;
M=ceil(6.6*pi/tr_width)+1;
n=[0:
1:
M-1];
wc=(ws+wp)/2;
hd=ideal_lp(wc,M);
%生成hamming窗
w_ham=(hamming(M))';
%频域图像的绘制
h=hd.*w_ham;
freqz(h,[1])
figure
(2);
subplot(2,2,1),stem(n,hd);title('理想脉冲响应')
axis([0M-1-0.30.3]);xlabel('n');ylabel('hd(n)')
xa=0.*n;
holdon
plot(n,xa,'k');
holdoff
subplot(2,2,2),stem(n,w_ham);title('hamming窗')
axis([0M-1-0.31.2]);xlabel('n');ylabel('w(n)')
subplot(2,2,3),stem(n,h);title('实际脉冲响应')
axis([0M-1-0.30.3]);xlabel('n');ylabel('h(n)')
holdon
plot(n,xa,'k');
holdoff
其响应波形如图7所示。
图7hamming窗函数设计FIR滤波器
【实例3】设带通滤波器的指标为
选择Blackman窗来实现这个滤波器。
解:
MATLAB源程序为
%数字滤波器指标
ws1=0.2*pi;wp1=0.35*pi;
ws2=0.65*pi;wp2=0.8*pi;
As=60;
tr_width=min((wp1-ws1),(wp2-ws2));
M=ceil(11*pi/tr_width)+1;
n=[0:
1:
M-1];
wc1=(ws1+wp1)/2;
wc2=(ws2+wp2)/2;
hd=ideal_lp(wc2,M)-ideal_lp(wc1,M);
%生成blackman窗
w_bla=(blackman(M))';
h=hd.*w_bla
%频域图像的绘制
freqz(h,[1])
figure
(2);
subplot(2,2,1),stem(n,hd);title('idaelimpulseresponse')
axis([0M-1-0.40.5]);xlabel('n');ylabel('hd(n)')
xa=0.*n;
holdon
plot(n,xa,'k');
holdoff
subplot(2,2,2),stem(n,w_bla);title('blackmanwindow')
axis([0M-101.1]);xlabel('n');ylabel('w(n)')
subplot(2,2,3),stem(n,h);title('actualimpulseresponse')
axis([0M-1-0.40.5]);xlabel('n');ylabel('h(n)')
holdon
plot(n,xa,'k');
holdoff
其响应波形如图8所示。
图8blackman窗函数设计FIR滤波器
5.3基于MATLAB的GUI界面设计。
5.4基于MATLAB的GUI滤波器设计方案实现。
5.5运行调试和界面美化以及方案优化改进。
六、基于MATLAB的滤波器结果设计分析
6.1滤波器实现展示:
1.FIR低通滤波器频率相位特性响应图及滤波图:
2.FIR高通滤波器频率相位特性响应图及滤波图:
3.FIR带通滤波器频率相位特性响应图及滤波图:
4.IIR低通滤波器频率相位特性响应图及滤波图:
5.IIR高通滤波器频率相位特性响应图及滤波图:
6.IIR带通滤波器频率相位特性响应图及滤波图:
6.2结果分析:
(详见设计流程)。
本次数字信号处理课程设计,GUI只是一个工具,一个托体,通过GUI面向对象的人性化简化了信号处理的抽象性;便于理解。
对MATLAB函数的调用设计滤波器以及对读入音频处理滤波,得到其时域,频域图。
七、设计经验总结
Malta无疑是一款功能超强的数学软件,基于它可实现多领域的科学研究和开发。
在本次课程设计中,收获最多的是如何看待问题以及如何解决问题。
一个很现实的问题就是搞科研开发,遇到问题是必不可免,关键是如何看待这样的问题,持怎样的心态,持怎样的目的。
再就是问题是遇到了,要如何解决以及用怎样最优的方案去处理,都得下心思去思考和摸索。
再一个很现实的问题是,要学会学习,学会自学。
毕竟社会日新月异,科技蓬勃发展,新领域新知识不断涌现,怎么样才能在社会保持有效的立足,怎样才能翻弄潮流,面对社会变迁做到运筹帷幄之中呢?
在我看来,关键要保持良好心态,保持一种积极永不妥协永不放弃的信念,不断进取不断开拓。
在学好基本课内的专业知识之外,应尽量多补充课外的和本专业相协调的专业知识。
强势未来就业利剑。
参考文献:
【1】MATLAB与科学计算
【2】精通基于MATLAB的GUI设计
【3】MATLAB在电子类学科的应用
参考网址:
【1】
【2】
【3】
附录:
方案主体代码:
globaly1;
globalfs;
globalfe;
globalbits;
[fname,dirpath]=uigetfile('*.wav');
myfile=[dirpathfname]
[y1,fe,bits]=wavread(myfile);
[H,f]=freqz(b,1,512,fs);%
axes(handles.axes4);
plot(f,20*log10(abs(H)));
xlabel('频率/Hz');
ylabel('幅值');
gridon;
title('滤波器频率响应曲线');
fb=str2num(get(handles.edit1,'string'));
fc=str2num(get(handles.edit2,'string'));
FB=str2