Matlab语音信号处理
目录
1绪 论 ............................................................................................ 2 2设计分析 ...................................................................................... 3
2.1设计思路 ............................................................................. 3
2.2
................................................................................. 4 3详细设计内容与调试分析 ........................................................... 5
3.1选取语音 ............................................................................. 5
3.2原始信号 ............................................................................. 6
3.3加噪信号 ............................................................................. 7
3.4双线性法设计滤波器(巴特沃思) .................................. 8
3.5窗函数法设计滤波器(汉明窗) ........................................ 10
3.6界面的设计 ....................................................................... 12 小结 ............................................................................................... 13 谢辞 ............................................................................................... 14 参考文献 ....................................................................................... 15 附录A主程序................................................................................ 16 附录B界面程序 ............................................................................ 20
1绪 论
数字信号处理是一门理论和技术发展十分迅速、应用非常广泛的交叉性、前沿性的学科,它的理论性和实践性都很强, 即它的理论基础涉及到众多的学科,但其成果又为这些学科起着重要的促进作用;它的应用范围甚广,如通信、航天航空、雷达、地震勘探,以及生物医学等领域。因此它已成为理论与实践并重,在高新科技领域中占有重要地位的新兴学科。
数字信号处理的特点是:
(1)要求的数学知识多,包括高等代数、数值分析、概率统计、随机过程等。
(2)要求掌握的基础知识强,网络理论、信号与系统是本课程的理论基础。
(3)与其他学科密切相关,即与通信理论、计算机、微电子技术不可分,又是人工智能、模式识别、神经网络等新兴学科的理论基础之一。
MATLAB是数字信号处理技术实现的重要手段。随着信号处理应用领域的不断扩大,一般的信号处理
难以完全满足实际工程的需要,从而促使人们在信号处理理论和方法方面不断创新与探索,另一方面,MATLAB的信号处理工具箱也随着信号处理理论与方法的发展而丰富。即MATBAB信号处理工具箱随着信号处理理论与方法的发展而发展,同时又为信号处理理论与方法的实现及工具应用提供了有效的帮助。
Matlab语言是一种广泛应用于工程计算及数值分析 领域的新型高级语言,Matlab功能强大、简单易学、编程效率高,深受广大科技工作者的欢迎。特别是Matlab还具有信号分析工具箱,不需具备很强的编程能力,就可以很方便地进行信号分析、处理和设计。因此,本次的数字信号处理选择用Matlab进行课程设计。
2设计分析
2.1设计思路
该课程设计的主要目的是根据具体语音情况利用双线性变换法和窗函数法设计出不同的滤波器将语音中的噪声滤除掉,并进行对比分析,最后为了使操作更加清晰和简洁,在此可以设计一个界面以便操作。因此先应利用电脑中Windows下的录音机,录制一段清晰的语音,时间不要太长,然后在Matlab软件平台下,利用函数wavread对语音信号进行采样,记住采样频率和采样点数。画出语音信号的时域波形;然后对语音号进行快速傅里叶变换,得到信号的频谱特性,然后在
添加余弦信号的噪声程序,即将噪声加到语音信号中,原信号程序中
绘制出加噪后的时域图和频谱图,再根据加噪后信号的频谱特性确定各性能指标,选择合适的滤波器,用双线性变换法和窗函数法设计滤波器,接着用自己设计的各滤波器分别对采集的信号进行滤波,在Matlab中,FIR滤波器利用函数fftfilt对信号进行滤波,IIR滤波器利用函数filter对信号进行滤波,将滤波前后的信号通过图像和语音进行对比分析,最后设计一个界面以便更容易对程序的操作。
具体思路为:录音
采样
加噪声
设计滤波器
滤波
设计界面
2.2流程图
录制语音
采样 时域波形 频谱图
、
添加噪声 频谱图 时域波形
确定性能指标
低通巴特沃思 设计滤波器 海明窗
滤波器
设计界面 滤波 滤波
频谱图 时域波形 时域波形 频谱图
频谱图
注:在设计过程中应将所得出的的语音和图像(时域图和频谱图)进行对比。
3详细设计内容与调试分析 3.1选取语音
利用Windows下的录音机,录制一段清晰的语音,时间不应过长,在此录制时间为4.22秒,窗口图形如下:
音频属性:
3.2原始信号
在Matlab软件平台下,利用函数wavread对语音信号进行采样, x1=wavread('f:\zmzf1.wav');在此语音信号采样频率为8kHz,利用语音回放语句sound(x1,22050);实现播放。利用plot(x1)绘制出原始信号时域波形。
利用语句freqz(x1)绘制出原始语音信号频率响应图,
注:上方为幅频图,下方为相位图
通过y1=fft(x1,1024)语句对原始语音信号做1024点FFT变换,以及通过plot(f,abs(y1(1:512)))和plot(abs(y1(1:512)))语句分别绘出原始语音信号的频谱图和原始语音信号FFT频谱图,并将其进行对比。
3.3加噪信号
利用同维数信号相加法将噪音信号d加到原始信号x1中,形成
x2=x1+d,再利用plot(t,x2)语句绘制出加噪后的信号时域图,并利用
sound(x2,8000);回放加噪后的语音。
所要添加的高频余弦信号为d=[0.3*cos(2*pi*3200*t)]';由于所要添加的信号是一维信号,因此在选取语音信号时要使其为单声道信号,因为只有维数相同的信号才能相加。除此之外,还应将所加噪声信号的点数调整到与原始信号相同,所以在此有t=0:1/8000:(size(x1)-1)/8000; 再利用x2=x1+d;将噪声添加到原始信号中形成两种信号的混合,从对比频谱图中可以发现加了噪音信号以后有一定的区别。
3.4双线性法设计滤波器(巴特沃思)
巴特沃思滤波器设计分析:
巴特沃思滤波器属于IIR滤波器,使用时一般要用到的语句有[N,Wn]=buttord(Wp,Ws,Rp,Rs)和N,Wn]=buttord(Wp,Ws,Rp,Rs,’s’) 其中Wp,Ws分别为通带和阻带截止频率,为归一化频率,Rp,Rs分别为通带和阻带波纹。
根据加入噪声后的频谱特性确定各参数,通过观察频谱图发现在此取通带截止频率wp=0.5π,阻带截止频率ws=0.6π; 通带波纹Rp=1dB; 阻带波纹Rs=15dB; 由于这所取的都是模拟指标,
在此可以先将模拟指标转换成数字指标,即使wp1=2/Ts*tan(wp/2);
ws1=2/Ts*tan(ws/2);再利用语句[N,Wn]=buttord(wp1,ws1,Rp,Rs,'s')确定滤波器的最小阶数N,通过[Z,P,K]=buttap(N) 创建butterworth模拟滤波器,再用双线性变换法实现模拟滤波器到数字滤波器的转换,即[bz,az]=bilinear(b,a,fs),它由模拟滤波器H(S)得到数字滤波器H(z),其中b,a,fs分别为H(s)的分子、分母多项式的系数向量和抽样频率,bz,az分别为H(z) 的分子、分母多项式的系数向量。最后通过[H,W]=freqz(bz,az) 绘制出低通巴特沃思滤波器频率响应曲线。
滤波前和滤波后对比分析:
利用plot( )语句分别绘制出滤波前和滤波后的时域图和频谱图,并进行对比。在此所添加的是频率为3200Hz的余弦信号,通过图像发现通过巴特沃思低通滤波器滤波前后图像都发生了一定的变化,可以将噪声基本上滤除掉,但在滤除的同时也会将原始信号的一部分滤除掉,利用f1=filter(bz,az,x2);sound(f1,22050);回放出滤波后的语音,并进行对比发现语音基本上不会受到太大的影响。其中bz,az分别为H(z) 的分子、分母多项式的系数向量,而x2为前面提到的加噪信号,
3.5窗函数法设计滤波器(汉明窗)
窗函数法设计分析:
海明窗的函数表达式为w(n)={0.54-0.46cos[2nπ/(N-1)] }RN(n),主瓣宽度为4×2π/N=8π/N,过度带宽Δw=3.3×2π/N,最后利用语句b=fir1(N,wn/pi,hamming(N+1)); 实现函数的选择与设计,上式中wn/pi是频率向量,其值在0~1之间,hamming(N+1)是与wn/pi相对应的所希望的幅频响应,不指定窗函数类型时,则自动选择海明窗。
该设计中首先通过阻带衰减来确定窗形状,由过度带宽确定阶次N,通过观察加噪后的信号频谱图发现可选通带截止频率 wp=0.5π; 阻带截止频率 ws=0.6π; 则Δw=wp-ws=0.1π, 所以利用公式得N=6.6π/0.1π=66,又wn=(0.2+0.3) π/2,所以在此选择海明窗设计滤波器。
滤波前和滤波后对比分析:
利用plot( )语句分别绘制出滤波前和滤波后的时域图和频谱图,并进行对比。通过图像发现通过海明窗滤波器滤波前后图像同样发生了一定的变化,可以将噪声基本上滤除掉,但在滤除的同时也会将原始信号的一部分滤除掉,其滤波效果和巴特沃思低通滤波器基本上差不多。利用f2=filter(bz,az,x2);sound(f2,22050);对经过汉明窗滤波器滤波后的语音,其中bz,az分别为H(z) 的分子、分母多项式的系数向量,而x2为前面提到的加噪信号。
3.6界面的设计
界面设计分析:先将四个主程序保存在m文件中,并分别对它们进行命名,再通过MATLAB中的File->new->GUI->OK,便出现设计界面的界面图,在此我将语音与图像分开设计按钮,通过双击按钮输入m文件名和按钮名,便可完成界面的设计。在此过程中发现如果m文件不能顺利通过运行则无法完成界面的设计,在输入汉字按钮名后无法在上面显示按钮名,但在运行时却能出现汉字,如图所示。界面的设计使操作更加方便和清晰。
小结
在设计之前,在图书馆参考了许多相关的资料,同时也从网络上查看和下载了一些有关的程序。在设计前又回顾了MATLAB的操作和滤波器的有关知识,从而有了基本的思路。
但着手设计时,又出现了许多未预料到的问题,首先便是语音信号的选取,这不仅仅是通过电脑中Windows下的录音机录制一段信号便可,好的语音能使设计更加清晰简单,在该录制中要求录制的语音必须为8位单声道的,且语音的清晰度越高越好,录制时间不应过长,光这一项便花费了很长一段时间。录音完成后利用函数wavread对语音信号进行采样,在采样时采样频率为8Hz,但在做程序时发现可以使用其它的采样频率,但采样后声音与原声明显不同,在设计中选择不同点数进行FFT变换时所得出的结果不同,我在设计中尝试了很多种情况。在绘制原始信号的时域波形和频谱图遇到的困难不是很多,但在添加噪声信号时却出现了很多问题,例如开始时始终无法将两个信号加在一块,后来在老师的解说下才发现录音的格式不符合,即所录的音必须为单声道的,这样才能和一维的高频余弦信号相加,解决完这点后又发现加噪声的程序会出现警告,这点会影响到M文件的使用,即到了后面设计界面时无法进行,后经过自己的检查和同学的指导最后还是解决了这一问题。在加完噪声后开始确定各参数设计滤波器时,这时需要前面学到的知识再次应用,根据具体的参数选择设计合适的滤波器对加噪后的信号进行滤波,在该设计中我选择的是巴特沃思低通滤波器和汉明窗滤波器,在设计这两种滤波器时要求对它们的原理都要有一定的了解,这同时给了我们复习课程的机会,通过学到的和课本最终还是按要求得出了合适的滤波器,通过对比滤波前后的频谱特性图发现滤波后的信号不能完全达到原来信号那么完整和清晰,这也许是误差,也许是滤波器的选择还是存在一定的问题,但我认为这应该是误差吧~在设计界面的时要将文件分别保存到M文件中,发现汉字不能显示在最外层,设计完界面后利用界面来调取程序感觉简洁了好多,使人很容易就可以了解设计的思路。
通过这次滤波器的设计,我发现了以往学习中的许多不足,也让我掌握了以往许多掌握的不太牢的知识,感觉学到了很多东西,相信在以后的学习中会有很大的帮助~
谢辞
能够顺利的完成本次数字信号处理的课程设计,首先要感谢余慧娟老师的指导,以及通信教研室对同学们上机的支持;感谢学校教育部给我们这次课程设计的机会;感谢学校图书馆能给我们提供丰富的资料。同时还要感谢多位同学给于我的帮助,因为他们不光在课程设计中指导我,还有当发现了有用的资料时,总能共同分享,大家在一块讨论问题,解决问题,让我再次感觉到团队的力量是无比强大的。最后感谢余慧娟老师在本次课程设计中对我们的严格要求,使我不仅回顾到了过去学到的知识,同时也学到了很多从未接触过的以及从课本上无法学到的东西。
参考文献
1、 程佩青 . 数字信号处理教程(第二版). 北京:清华大学出版社 ,2001 2、 刘敏 魏玲(Matlab通信仿真与应用,M,(北京:国防工业出版社, 2001 3、 陈亚勇. MATLAB信号处理详解. 北京:人民邮电出版社, 2001 4、 黄文梅 熊桂林 杨勇 . 信号分析与处理——MATLAB语言及应用 . 西
安:西安交通大学出版社, 2002
5、 刘波 文忠 曾涯 . MATLAB信号处理 . 北京:电子工业出版社,20006 6、 李勇 除震. MATLAB辅助现代工程数字信号处理. 西安:西安电子科技大
学出版社,2002
附录A主程序
程序1:用MATLAB对原始语音信号进行分析,画出它的时域波形和频谱
fs=8000; %语音信号采样频率为8000 x1=wavread('f:\zmzf1.wav'); %读取语音信号的数据,赋给变量x1 sound(x1,8000); %播放语音信号
y1=fft(x1,1024); %对信号做1024点FFT变换 f=fs*(0:511)/1024;
figure(1)
plot(x1) %绘制原始语音信号的时域图形 title('原始语音信号时域图');
xlabel('time n');ylabel('fuzhi n'); figure(2)
freqz(x1) %绘制原始语音信号的频率响应图 title('原始语音信号频率响应图');
figure(3)
subplot(2,1,1);
plot(f,abs(y1(1:512))); %绘制原始语音信号的频谱图 title('原始语音信号频谱')
xlabel('Hz');ylabel('fuzhi'); subplot(2,1,2);
plot(abs(y1(1:512))) %做原始语音信号的FFT频谱图 title('原始语音信号FFT频谱');
程序2:给原始的语音信号加上一个高频余弦噪声,频率为3200Hz。绘制出
加噪后的语音信号时域和频谱图,与原始信号对比,可以很明显的看出区
别。
fs=8000; %语音信号采样频率为8000 x1=wavread('f:\zmzf1.wav'); %读取语音信号的数据,赋给变量x1 f=fs*(0:511)/1024;
t=0:1/8000:(size(x1)-1)/8000; %将所加噪声信号的点数调整到与原始信号的点数相同
d=[0.3*cos(2*pi*3200*t)]'; %噪声为3200Hz的余弦信号 x2=x1+d; %将两信号相加
sound(x2,8000); %播放加噪声后的语音信号
y2=fft(x2,1024); %对加噪信号做1024点FFT变换 figure(1)
plot(t,x2) %绘制加噪后语音信号的时域波形 title('加噪后的信号时域波形');
xlabel('time n');ylabel('fuzhi n'); figure(2)
subplot(2,1,1);
plot(f,abs(y1(1:512))); %绘制原始语音信号频谱
title('原始语音信号频谱');
xlabel('Hz');ylabel('fuzhi'); subplot(2,1,2);
plot(f,abs(y2(1:512))); %绘制加噪后的信号频谱
title('加噪后的信号频谱');
xlabel('Hz');ylabel('fuzhi');
程序3:双线性变换法设计Butterworth滤波器
fs=8000; %语音信号采样频率为8000Hz x1=wavread('f:\zmzf1.wav'); %读取语音信号的数据,赋给变量x1 t=0:1/8000:(size(x1)-1)/8000; %将所加噪声信号的点数调整到与原始信号的点数相同
d=[0.3*cos(2*pi*3200*t)]'; %噪声为3200Hz的余弦信号 x2=x1+d; %将两信号相加
wp=0.5*pi; %通带截止频率
ws=0.6*pi; %阻带截止频率
Rp=1; %通带波纹
Rs=15; %阻带波纹
fs=8000;
Ts=1/fs; %采样周期
wp1=2/Ts*tan(wp/2); %将模拟指标转换成数字指标 ws1=2/Ts*tan(ws/2);
[N,Wn]=buttord(wp1,ws1,Rp,Rs,'s'); %选择滤波器的最小阶数
[Z,P,K]=buttap(N); %创建butterworth模拟滤波器 [Bap,Aap]=zp2tf(Z,P,K); [b,a]=lp2lp(Bap,Aap,Wn); [bz,az]=bilinear(b,a,fs); %用双线性变换法实现模拟滤波器到数字滤波器的转换 [H,W]=freqz(bz,az);
figure(1)
plot(W*fs/(2*pi),abs(H)) %巴特沃思低通滤波器绘制频率响应曲线 grid
xlabel('频率,Hz')
ylabel('频率响应幅度')
title('Butterworth')
f1=filter(bz,az,x2);
figure(2)
subplot(2,1,1)
plot(t,x2) %绘制出滤波前的时域图 title('滤波前的时域波形');
subplot(2,1,2)
plot(t,f1); %绘制出滤波后的时域图 title('滤波后的时域波形');
sound(f1,8000); %播放滤波后的信号 F0=fft(f1,1024);
f=fs*(0:511)/1024;
figure(3)
y2=fft(x2,1024);
subplot(2,1,1);
plot(f,abs(y2(1:512))); %绘制出滤波前的频谱图 title('滤波前的频谱')
xlabel('Hz');ylabel('fuzhi'); subplot(2,1,2)
F1=plot(f,abs(F0(1:512))); %绘制出滤波后的频谱图 title('滤波后的频谱')
xlabel('Hz');ylabel('fuzhi');
程序4:窗函数法设计滤波器:
fs=8000; %语音信号采样频率为8000Hz x1=wavread('f:\zmzf1.wav'); %读取语音信号的数据,赋给变量x1 t=0:1/8000:(size(x1)-1)/8000; %将所加噪声信号的点数调整到与原始信号的点数相同
d=[0.3*cos(2*pi*3200*t)]'; %噪声为3200Hz的余弦信号 x2=x1+d; %将两信号相加
wp=0.5*pi; %通带截止频率
ws=0.6*pi; %阻带截止频率
wdelta=ws-wp; N=ceil(6.6*pi/wdelta); %取整
wn=(0.2+0.3)*pi/2;
b=fir1(N,wn/pi,hamming(N+1)); %选择窗函数,并归一化截止频率 figure(1)
freqz(b,1,512)
f2=filter(bz,az,x2)
figure(2)
subplot(2,1,1)
plot(t,x2)
title('滤波前的时域波形'); %绘制滤波前的时域波形 subplot(2,1,2)
plot(t,f2); %绘制滤波后的时域波形 title('滤波后的时域波形');
sound(f2,8000); %播放滤波后的语音信号 F0=fft(f2,1024);
f=fs*(0:511)/1024;
figure(3)
y2=fft(x2,1024);
subplot(2,1,1);
plot(f,abs(y2(1:512))); %绘制滤波前信号的频谱 title('滤波前的频谱')
xlabel('Hz');ylabel('fuzhi'); subplot(2,1,2)
F2=plot(f,abs(F0(1:512))); %绘制滤波后信号的频谱 title('滤波后的频谱')
xlabel('Hz');ylabel('fuzhi');
附录B界面程序
function varargout = zmz(varargin)
% ZMZ M-file for zmz.fig
% ZMZ, by itself, creates a new ZMZ or raises the existing % singleton*.
%
% H = ZMZ returns the handle to a new ZMZ or the handle to % the existing singleton*.
%
% ZMZ('CALLBACK',hObject,eventData,handles,...) calls the local % function named CALLBACK in ZMZ.M with the given input arguments. %
% ZMZ('Property','Value',...) creates a new ZMZ or raises the % existing singleton*. Starting from the left, property value pairs are % applied to the GUI before zmz_OpeningFunction gets called. An % unrecognized property name or invalid value makes property application % stop. All inputs are passed to zmz_OpeningFcn via varargin. %
% *See GUI Options on GUIDE's Tools menu. Choose "GUI allows only one % instance to run (singleton)".
%
% See also: GUIDE, GUIDATA, GUIHANDLES
% Edit the above text to modify the response to help zmz
% Last Modified by GUIDE v2.5 25-Mar-2002 12:19:15
% Begin initialization code - DO NOT EDIT
gui_Singleton = 1;
gui_State = struct('gui_Name', mfilename, ...
'gui_Singleton', gui_Singleton, ...
'gui_OpeningFcn', @zmz_OpeningFcn, ...
'gui_OutputFcn', @zmz_OutputFcn, ...
'gui_LayoutFcn', [] , ...
'gui_Callback', []);
if nargin & isstr(varargin{1})
gui_State.gui_Callback = str2func(varargin{1}); end
if nargout
[varargout{1:nargout}] = gui_mainfcn(gui_State, varargin{:}); else
gui_mainfcn(gui_State, varargin{:});
end
% End initialization code - DO NOT EDIT
% --- Executes just before zmz is made visible.
function zmz_OpeningFcn(hObject, eventdata, handles, varargin) % This function has no output args, see OutputFcn.
% hObject handle to figure
% eventdata reserved - to be defined in a future version of MATLAB % handles structure with handles and user data (see GUIDATA) % varargin command line arguments to zmz (see VARARGIN)
% Choose default command line output for zmz
handles.output = hObject;
% Update handles structure
guidata(hObject, handles);
% UIWAIT makes zmz wait for user response (see UIRESUME) % uiwait(handles.figure1);
% --- Outputs from this function are returned to the command line. function varargout = zmz_OutputFcn(hObject, eventdata, handles) % varargout cell array for returning output args (see VARARGOUT); % hObject handle to figure
% eventdata reserved - to be defined in a future version of MATLAB % handles structure with handles and user data (see GUIDATA)
% Get default command line output from handles structure varargout{1} = handles.output;