基于MATLAB的频谱分析和信号去噪.docx
《基于MATLAB的频谱分析和信号去噪.docx》由会员分享，可在线阅读，更多相关《基于MATLAB的频谱分析和信号去噪.docx（22页珍藏版）》请在冰豆网上搜索。
基于MATLAB的频谱分析和信号去噪
基于MATLAB的频谱分析和信号去噪
摘要
滤波器设计在数字信号处理中占有极其重要的地位，FIR数字滤波器和IIR滤波器是滤波器设计的重要组成部分。
利用MATLAB信号处理工具箱可以快速有效地设计各种数字滤波器。
课题基于MATLAB有噪音语音信号处理的设计与实现，综合运用数字信号处理的理论知识对加噪声语音信号进行时域、频域分析和滤波。
通过理论推导得出相应结论，再利用MATLAB作为编程工具进行计算机实现。
在设计实现的过程中，使用窗函数法来设计FIR数字滤波器，用巴特沃斯、切比雪夫和双线性变法设计IIR数字滤波器，并利用MATLAB作为辅助工具完成设计中的计算与图形的绘制。
通过对对所设计滤波器的仿真和频率特性分析，可知利用MATLAB信号处理工具箱可以有效快捷地设计FIR和IIR数字滤波器，过程简单方便，结果的各项性能指标均达到指定要求。
关键词MATLAB/滤波器/频域分析/滤波
BasedontheMATLABspectrumanalysisandsignaldenoising
ABSTRACT
Filterdesignindigitalsignalprocessinginoccupiesanextremelyimportantstatus,FIRdigitalfiltersandfilterIIRfilteristheimportantcomponent.MATLABsignalprocessingtoolboxcaneffectivelydesignofdigitalfilter.SubjecthavenoisebasedonMATLABspeechsignalprocessing,thedesignandimplementationofintegrateduseofdigitalsignalprocessingtheoryknowledgeofspeechsignalwithnoiseaddingtimedomainandfrequencydomainanalysisandfiltering.Throughtheoreticalderivationcorrespondingconclusionasaprogrammingtool,usingMATLABforcomputerimplementation.Inthedesignandimplementationofprocess,USESwindowfunctionmethodtodesignFIRdigitalfilters,withbartleavenworth,chebyshevandbilinearanshiIIRdigitalfilter,anddesignMATLABastheauxiliarytoolcompletedesignofcomputingandgraphicsdrawing.Throughthedesignofthefilterselvesorsimulationandfrequencycharacteristicswereanalyzed,usingMATLABsignalprocessingtoolboxthatcaneffectivelyquicklydesignFIRandIIRdigitalfilters,processissimpleandconvenient,resultsofvariousperformanceindicatorsareasspecifiedrequirements.
KEYWORDSMATLAB,filter,frequencydomainanalysis,filtering
1绪论
数字滤波器,是数字信号处理中及其重要的一部分。
随着信息时代和数字技术的发展，受到人们越来越多的重视。
数字滤波器可以通过数值运算实现滤波，所以数字滤波器处理精度高、稳定、体积小、重量轻、灵活不存在阻抗匹配问题，可以实现模拟滤波器无法实现的特殊功能。
数字滤波器种类很多，根据其实现的网络结构或者其冲激响应函数的时域特性，可分为两种，即有限冲激响应（FIR，FiniteImpulseResponse）滤波器和无限冲激响应（IIR，InfiniteImpulseResponse）滤波器。
FIR滤波器结构上主要是非递归结构，没有输出到输入的反馈，系统函数H（z）在处收敛，极点全部在z=0处（因果系统），因而只能用较高的阶数达到高的选择性。
FIR数字滤波器的幅频特性精度较之于IIR数字滤波器低，但是线性相位，就是不同频率分量的信号经过fir滤波器后他们的时间差不变，这是很好的性质。
FIR数字滤波器是有限的单位响应也有利于对数字信号的处理，便于编程，用于计算的时延也小，这对实时的信号处理很重要。
FIR滤波器因具有系统稳定，易实现相位控制，允许设计多通带（或多阻带）滤波器等优点收到人们的青睐。
IIR滤波器采用递归型结构，即结构上带有反馈环路。
IIR滤波器运算结构通常由延时、乘以系数和相加等基本运算组成，可以组合成直接型、正准型、级联型、并联型四种结构形式，都具有反馈回路。
同时，IIR数字滤波器在设计上可以借助成熟的模拟滤波器的成果，如巴特沃斯、契比雪夫和椭圆滤波器等，有现成的设计数据或图表可查，在设计一个IIR数字滤波器时，我们根据指标先写出模拟滤波器的公式，然后通过一定的变换，将模拟滤波器的公式转换成数字滤波器的公式。
本文研究的主要内容中是利用一段语音信号对其进行加噪处理，然后利用FIR和IRR两种滤波器对加有噪声的语音信号进行过滤和频谱分析，画出滤波之后的频谱图与时域波形，并对信号滤波前后进行分析比较，分析信号的变化，来进行比较IIR与FIR两种滤波器对语音信号的处理的特点。
2基于滤波器的语音信号的处理
选择一个语音信号作为分析的对象，或录制一段语音信号，对其进行频谱分析；利用MATLAB中的随机函数产生噪声加入到语音信号中，模仿语音信号被污染，并对其进行频谱分析；设计FIR和IIR数字滤波器，并对被噪声污染的语音信号进行滤波，分析滤波后信号的时域和频域特征，回放语音信号；最后，设计一个信号处理系统界面。
在所设计的系统界面上可以选择滤波器的类型。
2.1语音信号的采集
利用PC机上的声卡和WINDOWS操作系统可以进行数字信号的采集。
将话筒输入计算机的语音输入插口上,启动录音机。
按下录音按钮,接着对话筒说话“语音信号处理”,说完后停止录音,屏幕左侧将显示所录声音的长度。
点击放音按钮,可以实现所录音的重现。
以文件名“speech”保存入g:
\MATLAB\work中。
可以看到,文件存储器的后缀默认为.wav,这是WINDOWS操作系统规定的声音文件存的标准。
2.2设计FIR和IIR数字滤波器
IIR滤波器和FIR滤波器的设计方法完全不同。
IIR滤波器设计方法有间接法和直接法，间接法是借助于模拟滤波器的设计方法进行的。
其设计步骤是：
先设计过渡模拟滤波器得到系统函数H（s），然后将H（s）按某种方法转换成数字滤波器的系统函数H（z）。
FIR滤波器比鞥采用间接法，常用的方法有窗函数法、频率采样发和切比雪夫等波纹逼近法。
对于线性相位滤波器，经常采用FIR滤波器。
对于数字高通、带通滤波器的设计，通用方法为双线性变换法。
可以借助于模拟滤波器的频率转换设计一个所需类型的过渡模拟滤波器，再经过双线性变换将其转换策划那个所需的数字滤波器。
具体设计步骤如下：
（1）确定所需类型数字滤波器的技术指标。
（2）将所需类型数字滤波器的边界频率转换成相应的模拟滤波器的边界频率，转换公式为Ω=2/Ttan（0.5ω）
（3）将相应类型的模拟滤波器技术指标转换成模拟低通滤波器技术指标。
（4）设计模拟低通滤波器。
（5）通过频率变换将模拟低通转换成相应类型的过渡模拟滤波器。
（6）采用双线性变换法将相应类型的过渡模拟滤波器转换成所需类型的数字滤波器。
我们知道，脉冲响应不变法的主要缺点是会产生频谱混叠现象，使数字滤波器的频响偏离模拟滤波器的频响特性。
为了克服之一缺点，可以采用双线性变换法。
下面我们总结一下利用模拟滤波器设计IIR数字低通滤波器的步骤：
（1）确定数字低通滤波器的技术指标：
通带边界频率、通带最大衰减，阻带截止频率、阻带最小衰减。
（2）将数字低通滤波器的技术指标转换成相应的模拟低通滤波器的技术指标。
（3）按照模拟低通滤波器的技术指标设计及过渡模拟低通滤波器。
（4）用双线性变换法，模拟滤波器系统函数转换成数字低通滤波器系统函数。
如前所述，IIR滤波器和FIR滤波器的设计方法有很大的区别。
下面我们着重介绍用窗函数法设计FIR滤波器的步骤。
如下：
（1）根据对阻带衰减及过渡带的指标要求，选择串窗数类型（矩形窗、三角窗、汉宁窗、哈明窗、凯塞窗等），并估计窗口长度N。
先按照阻带衰减选择窗函数类型。
原则是在保证阻带衰减满足要求的情况下，尽量选择主瓣的窗函数。
（2）构造希望逼近的频率响应函数。
（3）计算h（n）.。
（4）加窗得到设计结果。
接下来，我们根据语音信号的特点给出有关滤波器的技术指标：
①低通滤波器的性能指标：
fp=1000Hz，fc=1200Hz，As=50db,Ap=1dB
②高通滤波器的性能指标：
fp=3500Hz，fc=4000Hz，As=50dB，Ap=1dB；
在Matlab中，可以利用函数fir1设计FIR滤波器，利用函数butter,cheby1和ellip设计IIR滤波器，利用Matlab中的函数freqz画出各步步器的频率响应。
hn=fir1（M，wc，window），可以指定窗函数向量window。
如果缺省window参数，则fir1默认为哈明窗。
其中可选的窗函数有RectangularBarlrttHammingHannBlackman窗，其相应的都有实现函数。
MATLAB信号处理工具箱函数buttpbuttorbutter是巴特沃斯滤波器设计函数，其有5种调用格式，本课程设计中用到的是[N,wc]=butter（N,wc,Rp,As,’s’）,该格式用于计算巴特沃斯模拟滤波器的阶数N和3dB截止频率wc。
MATLAB信号处理工具箱函数cheblap,cheblord和cheeby1是切比雪夫I型滤波器设计函数。
我们用到的是cheeby1函数，其调用格式如下：
[B,A]=cheby1（N,Rp,wpo,’ftypr’）
[B,A]=cheby1（N,Rp,wpo,’ftypr’,’s’）
函数butter,cheby1和ellip设计IIR滤波器时都是默认的双线性变换法，所以在设计滤波器时只需要代入相应的实现函数即可。
下面我们将给出FIR和IIR数字滤波器的主要程序。
=========================IIR低通滤波器=======================
IIR低通滤波器程序见附录3
图1IIR低通滤波器
======================IIR高通==========================
IIR高通滤波器程序见附录4
图2IIR高通滤波器
**************************FIR低通滤波器*****************************
FIR低通滤波器程序见附录5
图3FIR低通滤波器
***************************%FIR高通滤波器*******************
FIR高通滤波器程序见附录6
图4FIR高通滤波器
2.3用滤波器对加噪语音信号进行滤波
用自己设计的各滤波器分别对加噪的语音信号进行滤波，在Matlab中，FIR滤波器利用函数fftfilt对信号进行滤波，IIR滤波器利用函数filter对信号进行滤波。
函数fftfilt用的是重叠相加法实现线性卷积的计算。
调用格式为：
y=fftfilter（h,x,M）。
其中，h是系统单位冲击响应向量；x是输入序列向量；y是系统的输出序列向量；M是有用户选择的输入序列的分段长度，缺省时，默认的输入向量的重长度M=512。
函数filter的调用格式：
yn=filter（B，A.xn），它是按照直线型结构实现对xn的滤波。
其中xn是输入信号向量，yn输出信号向量。
3仿真及其结果分析
3.1语音信号的时频分析
利用MATLAB中的“wavread”命令来读入（采集）语音信号，将它赋值给某一向量。
再对其进行采样，记住采样频率和采样点数。
下面介绍Wavread函数几种调用格式。
（1）y=wavread（file）
功能说明：
读取file所规定的wav文件，返回采样值放在向量y中。
（2）[y,fs,nbits]=wavread（file）
功能说明：
采样值放在向量y中，fs表示采样频率（hz），nbits表示采样位数。
（3）y=wavread（file，N）
功能说明：
读取钱N点的采样值放在向量y中。
（4）y=wavread（file，[N1,N2]）
功能说明：
读取从N1到N2点的采样值放在向量y中。
接下来，对语音信号OriSound.wav进行采样。
其程序如下：
>>[y,fs,nbits]=wavered（‘OriSound’）;%把语音信号加载入Matlab仿真软件平台中
然后，画出语音信号的时域波形，再对语音信号进行频谱分析。
MATLAB提供了快速傅里叶变换算法FFT计算DFT的函数fft,其调用格式如下：
Xk=fft（xn,N）
参数xn为被变换的时域序列向量，N是DFT变换区间长度，当N大于xn的长度时，fft函数自动在xn后面补零。
，当N小于xn的长度时，fft函数计算xn的前N个元素，忽略其后面的元素。
在本次设计中，我们利用fft对语音信号进行快速傅里叶变换，就可以得到信号的频谱特性。
其程序见附录1。
程序结果如下图：
图5原始信号波形及频谱
3.2加噪后的语音信号及其频谱分析
本文中，利用MATLAB中的随机函数（rand或randn）产生噪声加入到语音信号中，模仿语音信号被污染，并对其频谱分析。
Randn函数有两种基本调用格式：
Randn（n）和Randn（m,n）,前者产生n×n服从标准高斯分布的随机数矩阵，后者产生m×n的随机数矩阵。
在这里，我们选用Randn（m,n）函数。
语音信号添加噪声及其频谱分析的主要程序见附录2.
程序结果如下图：
图6加噪后的波形及频谱分析
3.3通过滤波器处理的加噪语音信号的分析
*************************IIR低通滤波器************************
IIR低通滤波器主要程序见附录7
图7IIR低通滤波
%****************************IIR高通滤波器*****************************
IIR高通滤波器主要程序见附录8
图8IIR高通滤波
分析：
通过IIR滤波器滤波前后的频谱比较相近，但仔细观察仍可以看到途中结果滤掉了相应的成分。
通过声音回放能明显的听到前后两者声音的不同，IIR高通滤波后只能听到少量的杂音，IIR低通能听到原来音乐，但是声音低沉发闷，但是很接近原来的声音。
**************************FIR低通滤波器****************************
FIR低通滤波器主要程序见附录9
图9FIR低通滤波
****************************FIR高通滤波器**************************
FIR高通滤波器主要程序见附录10
图10FIR高通滤波
分析：
通过FIR高通滤波器滤波前后的频谱比较较大，FIR低通滤波器则是不太明显，结果滤掉了相应的成分。
通过声音回放能明显的听到前后两者声音的不同，FIR高通滤波后只能听到的是杂音，原来音乐的声音很模糊低，FIR低通滤波后能听到原来音乐，但是其中还掺杂少量杂音且低沉发闷。
3.4IIR与FIR滤波器处理结果的的比较分析
IIR数字滤波器采用递归型结构，即结构上带有反馈环路。
IIR滤波器运算结构通常由延时、乘以系数和相加等基本运算组成，可以组合成直接型、正准型、级联型、并联型四种结构形式，都具有反馈回路。
由于运算中的舍入处理，使误差不断累积，有时会产生微弱的寄生振荡。
（1）IIR数字滤波器的相位特性不好控制，对相位要求较高时，需加相位校准网络。
FIR滤波器则要求较低。
（2）与FIR滤波器的设计不同，IIR滤波器设计时的阶数不是由设计者指定，而是根据设计者输入的各个滤波器参数（截止频率、通带滤纹、阻带衰减等），由软件设计出满足这些参数的最低滤波器阶数。
在MATLAB下设计不同类型IIR滤波器均有与之对应的函数用于阶数的选择。
（3）IIR单位响应为无限脉冲序列FIR单位响应为有限的
（4）IIR幅频特性精度很高，不是线性相位的，可以应用于对相位信息不敏感的音频信号上；
（5）FIR幅频特性精度较之于iir低，但是线性相位，就是不同频率分量的信号经过FIR滤波器后他们的时间差不变。
这是很好的性质。
（6）IIR滤波器有噪声反馈，而且噪声较大，FIR滤波器噪声较小。
（7）IIR滤波器运算误差大，有可能出现极限环振荡，FIR相比之下运算误差较小，不会出现极限环振荡。
IIR幅频特性精度很高，不是线性相位的，可以应用于对相位信息不敏感的音频信号上；
FIR幅频特性精度较之于iir低，但是线性相位，就是不同频率分量的信号经过FIR滤波器后他们的时间差不变。
这是很好的性质。
第四章结论
语音信号处理是语音学与数字信号处理技术相结合的交叉学科，课题在这里不讨论语音学，而是将语音当做一种特殊的信号，即一种“复杂向量”来看待。
也就是说，课题更多的还是体现了数字信号处理技术。
从课题的中心来看，课题“基于MATLAB的频谱分析和信号去噪”是希望将数字信号处理技术应用于某一实际领域，这里就是指对语音及加噪处理。
作为存储于计算机中的语音信号，其本身就是离散化了的向量，我们只需将这些离散的量提取出来，就可以对其进行处理了。
这一过程的实现，用到了处理数字信号的强有力工具MATLAB。
通过MATLAB里几个命令函数的调用，很轻易的在实际语音与数字信号的理论之间搭了一座桥。
课题的特色在于它将语音信号看作一个向量，于是就把语音数字化了。
那么，就可以完全利用数字信号处理的知识来解决语音及加噪处理问题。
我们可以像给一般信号做频谱分析一样，来对语音信号做频谱分析，也可以较容易的用数字滤波器来对语音进行滤波处理。
通过比较加噪前后，语音的频谱和语音回放，能明显的感觉到加入噪声后回放的声音与原始的语音信号有很大的不同,前者随较尖锐的干扰啸叫声。
从含噪语音信号的频谱图中可以看出含噪声的语音信号频谱,在整个频域范围内分是布均匀。
其实，这正是干扰所造成的。
通过滤波前后的对比，低通滤波后效果最好，高通滤波后的效果最差。
由此可见，语音信号主要分布在低频段，而噪声主要分布在高频段。
致谢
在此论文撰写过程中，要特别感谢我的导师老师的指导与督促，同时感谢她的谅解与包容。
没有老师的帮助也就没有今天的这篇论文。
求学历程是艰苦的，但又是快乐的。
感谢我的辅导员段老师，谢谢他在这四年中为我们全班所做的一切，他不求回报，无私奉献的精神很让我感动，再次向他表示由衷的感谢。
在这四年的学期中结识的各位生活和学习上的挚友让我得到了人生最大的一笔财富。
在此，也对他们表示衷心感谢。
本文参考了大量的文献资料，在此，向各学术界的前辈们致敬！
参考文献
[1]WolfgangR.Podilchuk,C.Delp,EPerpetualwatermarksfordigitalimagesandvideo1999（07）
[2]MallatS.Zhong,SCharacterizationofsignalsformmultiscaleedges1992（07）
[3]WangH-JM.Su,P-C.Kuo,C-C.JWavelet-Baseddigitalimagewatermarking1998（12）
[4]高西全，丁玉美．数字信号处理．第3版．北京：
西安电子科技大学出版社，2008
[5]刘泉，阙大顺．数字信号处理原理与实现．北京：
电子工业出版社，2005
[6]张磊，毕靖，郭莲英．MATLAB实用教程．北京：
人民邮电出版社，2008
[7]张威．MATLAB基础与编程入门．西安：
西安电子科技大学出版社，2006
[5]周利清，苏菲．数字信号处理基础．北京：
北京邮电大学出版社，2005
附录
附录1：
clf;
Ft=8000;
Fp=1000;
Fs=1200;
wp=2*pi*Fp/Ft;
ws=2*pi*Fs/Ft;
fp=2*Ft*tan（wp/2）;
fs=2*Fs*tan（wp/2）;
[n11,wn11]=buttord（wp,ws,1,50,'s'）;%求低通滤波器的阶数和截止频率
[b11,a11]=butter（n11,wn11,'s'）;%求S域的频率响应的参数
[num11,den11]=bilinear（b11,a11,0.5）;%双线性变换实现S域到Z域的变换
[h,w]=freqz（num11,den11）;%根据参数求出频率响应
plot（w*8000*0.5/pi,abs（h））;
legend（'用butter设计'）;
grid;
附录2：
clf;
Ft=8000;
Fp=4000;
Fs=3500;
wp1=tan（pi*Fp/Ft）;%高通到低通滤波器参数转换
ws1=tan（pi*Fs/Ft）;
wp=1;
ws=wp1*wp/ws1;
[n13,wn13]=cheb1ord（wp,ws,1,50,'s'）;%求模拟的低通滤波器阶数和截止频率
[b13,a13]=cheby1（n13,1,wn13,'s'）;%求S域的频率响应的参数
[num,den]=lp2hp（b13,a13,wn13）;%将S域低通参数转为高通的
[num13,den13]=bilinear（num,den,0.5）;%利用双线性变换实现S域到Z域转换
[h,w]=freqz（num13,den13）;
plot（w*21000*0.5/pi,abs（h））;
title（'IIR高通滤波器'）;
legend（用cheby1设计'）;
附录3：
clf;
Ft=8000;
Fp=1000;
Fs=1200;
wp=2*Fp/Ft;
ws=2*Fs/Ft;
rp=1;
rs=50;
p=1-10.^（-rp/20）;%通带阻带波纹
s=10.^（-rs/20）;
fpts=[wpws];
mag=[10];
dev=[ps];
[n21,wn21,beta,ftype]=kaiserord（fpts
升级会员 