语音信号滤波去噪使用TUKEYWIN窗设计的FIR滤波器文档格式.docx
《语音信号滤波去噪使用TUKEYWIN窗设计的FIR滤波器文档格式.docx》由会员分享，可在线阅读，更多相关《语音信号滤波去噪使用TUKEYWIN窗设计的FIR滤波器文档格式.docx（14页珍藏版）》请在冰豆网上搜索。
该软件有三大特点：
一是功能强大；
二是界面友善、语言自然；
三是开放性强。
目前，Mathworks公司已推出30多个应用工具箱。
MATLAB在线性代数、矩阵分析、数值及优化、数理统计和随机信号分析、电路与系统、系统动力学、信号和图像处理、控制理论分析和系统设计、过程控制、建模和仿真、通信系统、以及财政金融等众多领域的理论研究和工程设计中得到了广泛应用。
MATLAB在信号与系统中的应用主要包括符号运算和数值计算仿真分析。
由于信号与系统课程的许多内容都是基于公式演算，而MATLAB借助符号数学工具箱提供的符号运算功能能基本满足信号与系统课程的需求。
例如，解微分方程、傅里叶正反变换、拉普拉斯正反变换、z正反变换等。
MATLAB在信号与系统中的另一主要应用是数值计算与仿真分析，主要包括函数波形绘制、函数运算、冲激响应与阶跃响应仿真分析、信号的时域分析、信号的频谱分析、系统的S域分析、零极点图绘制等内容。
数值计算仿真分析可以帮助学生更深入理解信号与系统的理论知识，并为将来使用MATLAB进行信号处理领域的各种分析和实际应用打下基础[2]。
2基本理论
2.1FIR滤波器
FIR滤波器：
有限长单位冲激响应滤波器，是数字信号系统中最基本的元件，它可以在保证任意幅频特性的同时具有严格的线性相频特性，同时其单位抽样响应是有限长的，因而滤波器是稳定的系统。
因此，FIR滤波器在通信、图像处理、模式识别等领域都有着广泛的应用。
FIR数字滤波器设计的基本步骤如下：
（1）确定指标
在设计一个滤波器之前，必须首先根据工程实际的需要确定滤波器的技术指标。
在很多实际应用中，数字滤波器常常被用来实现选频操作。
因此，指标的形式一般在频域中给出幅度和相位响应。
幅度指标主要以两种方式给出。
第一种是绝对指标。
它提供对幅度响应函数的要求，一般应用于FIR滤波器的设计。
第二种指标是相对指标。
它以分贝值的形式给出要求。
（2）逼近
确定了技术指标后，就可以建立一个目标的数字滤波器模型。
通常用理想的数字滤波器模型。
之后，利用数字滤波器的设计方法，设计出一个实际滤波器模型来逼近给定的目标。
（3）性能分析和计算机仿真
上两步的结果是得到以差分或冲激响应描述的滤波器。
根据这个描述就可以分析其频率特性和相位特性，以验证设计结果是否满足指标要求；
或者利用计算机仿真实现设计的滤波器，再分析滤波结果来判断。
FIR滤波器的设计问题实际上是确定能满足所要求的转移序列或脉冲响应的常数的问题，设计方法主要有窗函数、频率响应法和等波纹最佳逼近法等。
2.2窗口设计法
窗口设计法是一种通过截断和计权的方法使无限长非因果序列成为有限长脉冲响应序列的设计方法。
通常在设计滤波器之前，应该先根据具体的工程应用确定滤波器的技术指标。
在大多数实际应用中，数字滤波器常常被用来实现选频操作，所以指标的形式一般为在频域中以分贝值给出的相对幅度响应和相位响应。
窗口设计法步骤如下：
（1）根据过渡带宽及阻带衰减要求，选择窗函数的类型并估计窗口长度N。
窗函数的类型可根据最小阻带衰减AS独立选择。
（2）根据待求滤波器的理想频率响应求出理想单位脉冲响应hd（n）。
（3）由性能指标确定窗函数W（n）和长度N。
（4）求得实际滤波器的单位脉冲响应h（n），h（n）即为所设计FIR滤波器系数向量b（n）。
（2.1）
常见的窗函数性能表如下图1所示：
图1常见窗函数性能表
名称
滤波器
过渡带宽
最小阻带衰减
矩形
1.8π/M
21dB
PARZENWIN
6.6π/M
56dB
巴特利特
6.1π/M
25dB
FLATTOPWIN
19.6π/M
108dB
汉宁
6.2π/M
44dB
GAUSSWIN
5.8π/M
60dB
汉明
51dB
BARTHANNWIN
3.6π/M
40dB
布莱克曼
11π/M
74dB
BLACKMANHARRIS
16.1π/M
109dB
BOHMANWIN
51.5dB
CHEBWIN
15.2π/M
113dB
NUTTALLWIN
15.4π/M
TUKEYWIN
2.4π/M
22dB
2.3TUKEYWIN窗
W=tukeywin（M）在数组w中产生M点的Tukeywin窗函数。
由上表可知：
Tukeywin窗的滤波器的过滤带宽为2.4π/M，最小阻带衰减22db。
利用wvtool可绘制出Tukeywin窗函数波形如图2：
N=128;
t=tukeywin（N）;
Wvtool（t）;
图2Tukeywin窗函数波形
3设计步骤
3.1设计流程图
根据设计的要求，首先自己录制一段语音信号，对语音信号加入噪声干扰，再利用TUKEYWIN窗设计合理的FIR滤波器。
最后用滤波器对干扰后的语音信号进行滤波去噪。
具体设计流程图如下图3所示：
NO
YES
图3设计流程图
3.2录制音乐信号
从电脑上录制一段语音信号，并命名为“speech_dft.wav”。
在Matlab平台下调用wavread函数可采集到语音信号。
运行程序后得到的波形图如图4所示：
图4原始音乐信号时域波形图
画出语音信号的时域波形；
然后对语音号进行快速傅里叶变换，得到信号的频谱特性，并将原始音乐信号的波形图与加干扰后的波形图进行比较
运行程序后我们所得的波形图如图5所示：
图5干扰前后信号时域图与频谱图的比较
3.3滤波器设计
滤波器的设计就是要找到一组能满足特定滤波要求的系数向量a和b，其中b为系统函数的分子系数，a为系统函数分母系数。
在本次的课程设计中我所采用的就是利用TUKEYWIN窗函数来设计FIR滤波器如图6所示。
图6设计的滤波器图形
3.4信号滤波处理
滤波器设计完成后，在MATLAB平台上用函数filter实现滤波，
得到的滤波前后音乐信号的时域波形图和频谱图对比图如图7所示：
图7滤波前后语音信号的时域波形图和频谱图对比图
3.5结果分析
在MATLAB中，经过sound（x,fs,bits）函数，对经过TUKEYWIN窗设计的FIR滤波器之后的音乐信号进行回放，可以听出滤波之后的信号跟原始信号一样清晰，完全滤除掉了噪声的干扰，通过下面语句来进行语音信号回放比较：
>
sound（y_fil,fs,bits）;
播放滤波后的音乐信号
sound（x,fs,bits）;
所得结果证明了用TUKEYWIN窗设计的FIR滤波器和音乐信号去噪设计是成功的。
4出现的问题及解决方法
在本次课程设计中我遇到的问题如下：
1、不知道如何将音乐信号截断成自己所需的长度。
2、对利用TUKEYWIN窗函数设计FIR滤波器的设计步骤很生疏。
3、在采用TUKEYWIN窗函数设计的FIR滤波器时得不到理想的滤波器，因而信号的恢复不是特别理想。
针对以上问题，相应的解决方案如下：
1、向同学学习如何截取音乐。
2、自己上网查阅资料，或向图书馆借阅资料，掌握利用TUKEYWIN窗函数设计FIR滤波器的方法和步骤。
3、通过不断设置参数的值，最终达到最理想的值，设计出理想的滤波器，使信号得到理想恢复。
5结束语
本次的课程设计，我的任务是利用TUKEYWIN窗函数设计FIR滤波器对音乐信号滤波去噪。
在本次课程设计之前，我对TUKEYWIN窗函数完全没有了解，因此在看到这个题目时，我是一头雾水。
但是通过自己翻阅资料和询问同学，我掌握了用TUKEYWIN窗函数设计FIR滤波器的方法步骤，了解了窗函数的基本设计流程。
经过这两周的课程设计，我学会了很多东西。
我们通信工程专业是个实践性很强的专业，而我们在校大部分的学习时间都是花在理论学习上面，实践的机会很少。
因而我对很多所学的理论知识如何跟实践联系的概念很模糊，这次的课程设计给了我这个机会，加深了我对理论联系实际的理解，增强了自己独立分析问题和解决问题的能力，开阔了自己的思维。
还有让我看到了自己的不足，自己对本专业的相关知识掌握的还很少，还有很多知识都没掌握，还让我认识到解决问题的方法、途径很多，做事要开阔自己的思维，看待问题要从多个角度看。
在此我要感谢学校为我们提供这次课程设计的机会，感谢老师对我的悉心指导，也感谢同学对我的帮助。
这次的课程设计让我理论联系实际，不仅巩固了我们的理论知识，还提高了我的动手能力，在这次课程设计中我所学到的知识是我的财富，让我终身受益。
参考文献
[1]张志涌．精通MATLAB6.5版[M]．北京：
北京航空航天大学出版社，2003．
[2]约翰·
Ｇ·
普罗克斯．数字信号处理[M]．西安：
西安交通大学出版社，2009．
[3]张小虹．信号系统与数字信号处理[M]．第1版．西安：
西安电子科技出版社，2002．
[4]谢德芳．数字信号处理[M]．北京：
科学出版社，2005．
[5]郝建军．数字通信[M]．第2版．北京：
北京邮电大学出版社，2010．
[6]张威．MATLAB基础与编程入门[M]．西安：
西安电子科技大学出版社，2010．
附录一：
语音信号滤波去噪——使用TUKEYWIN滤波器
%程序功能：
在Matlab中，用窗口设计法设计FIR滤波器
%程序作者：
%最后修改日期：
%程序一：
[x,fs,bits]=wavread（'
speech_dftwav'
）;
%读取音乐信号
sound（x,fs,bits）;
%播放音乐信号
N=length（x）;
%计算音乐信号的长度
fn=2136;
%单频噪声频率
t=0:
1/fs:
（N-1）/fs;
%计算时间范围，样本数除以采样频率
x=x（:
1）'
;
%将双声道转为单声道
y=x+sin（fn*2*pi*t）;
%添加噪声
sound（y,fs,bits）;
%播放干扰后的音乐信号，应该可以明显听出有尖锐的单频啸叫声
%程序二：
X=abs（fft（x））;
Y=abs（fft（y））;
%对原始信号和加噪信号进行fft变换,取幅度谱
X=X（1:
N/2）;
Y=Y（1:
%截取前半部分
figure
（1）
subplot（2,2,1）;
plot（t,x）;
xlabel（'
时间（t）'
ylabel（'
幅度'
title（'
原始语音信号'
axis（[0,2.5,-0.7,0.7]）;
subplot（2,2,2）;
plot（f,X）;
频率（f）'
幅度谱'
title（'
原始语音信号幅度谱'
axis（[0,4000,0,300]）;
subplot（2,2,3）;
plot（t,y）;
加噪声后的语音信号'
subplot（2,2,4）;
plot（f,Y）;
加噪声后的语音信号幅度谱'
%程序三：
fpd=2000;
fsd=2050;
fsu=2150;
fpu=2200;
%带阻滤波器设计指标
Rp=1;
As=22;
%通带波纹和阻带衰减
fcd=（fpd+fsd）/2;
fcu=（fpu+fsu）/2;
df=min（（fsd-fpd）,（fpu-fsu））;
%计算上下边带中心频率，和频率间隔
wcd=fcd/fs*2*pi;
wcu=fcu/fs*2*pi;
dw=df/fs*2*pi;
%将Hz为单位的模拟频率换算为rad为单位的数字频率
wsd=fsd/fs*2*pi;
wsu=fsu/fs*2*pi;
M=ceil（2.4*pi/dw）+1;
%计算Tukeywin窗设计该滤波器时需要的阶数
n=0:
M-1;
%定义时间范围
w_tuk=tukeywin（M）;
%产生M阶的Tukeywin窗
hd_bs=ideal_lp（wcd,M）+ideal_lp（pi,M）-ideal_lp（wcu,M）;
%调用自编函数计算理想带阻滤波器的脉冲响应
h_bs=w_tuk'
.*hd_bs;
%用窗口法计算实际滤波器脉冲响应
[db,mag,pha,grd,w]=freqz_m（h_bs,1）;
%调用自编函数计算滤波器的频率特性
figure
（2）
plot（w/pi,db）;
w/pi'
db'
FIR滤波器的幅度响应图'
axis（[0.30.7-7010]）;
line（[0,0.7],[-1,-1],'
Color'
'
r'
LineWidth'
2,'
LineStyle'
--'
line（[0,0.7],[-22,-22],'
line（[wsd/pi,wsd/pi],[-70,10],'
line（[wsu/pi,wsu/pi],[-70,10],'
gridon;
plot（w/pi,mag）;
幅度mag'
axis（[01-0.21.2]）;
plot（w/pi,pha）;
相位pha'
滤波器相位响应图'
axis（[01-44]）;
stem（n,h_bs）;
n'
h（n）'
滤波器脉冲响应图'
%程序四：
y_fil=fftfilt（h_bs,y）;
%用设计好的滤波器对y进行滤波
Y_fil=abs（fft（y_fil））;
Y_fil=Y_fil（1:
%计算频谱取前一半
sound（y_fil,fs,bits）;
wavwrite（y_fil,fs,bits,'
1111恢复.wav'
%保存恢复的语音信号
wavwrite（y,fs,bits,'
1111加噪.wav'
%保存加噪的语音信号
figure（3）
subplot（3,2,1）;
subplot（3,2,2）;
subplot（3,2,3）;
subplot（3,2,4）;
subplot（3,2,5）;
plot（t,y_fil）;
滤波后语音信号'
subplot（3,2,6）;
plot（f,Y_fil）;
滤波后语音信号幅度谱'
附录二：
函数FREQZ_M.M定义：
function[db,mag,pha,grd,w]=freqz_m（b,a）;
%freqz子程序的改进版本
%------------------------------------
%[db,mag,pha,grd,w]=freqz_m（b,a）;
%db=[0到pi弧度]区间内的相对振幅（db）
%mag=[0到pi弧度]区间内的绝对振幅
%pha=[0到pi弧度]区间内的相位响应
%grd=[0到pi弧度]区间内的群迟延
%w=[0到pi弧度]区间内的501个频率样本向量
%b=Ha（z）的分子多项式系数（对FIRb=h）
%a=Ha（z）的分母多项式系数（对FIR:
a=[1]）
%
[H,w]=freqz（b,a,1000,'
whole'
H=（H（1:
1:
501））'
w=（w（1:
mag=abs（H）;
db=20*log10（（mag+eps）/max（mag））;
pha=angle（H）;
%pha=unwrap（angle（H））;
grd=grpdelay（b,a,w）;
%grd=diff（pha）;
%grd=[grd
（1）grd];
%grd=[0grd（1:
500）;
grd;
grd（2:
501）0];
%grd=median（grd）*500/pi;
附录三：
函数IDEAL_LP.M定义：
functionhd=ideal_lp（wc,M）;
%理想低通滤波器计算
%--------------------------------
%[hd]=ideal_lp（wc,M）
%hd=0toM-1之间的理想脉冲响应
%wc=截止频率（弧度）
%M=理想滤波器的长度
alpha=（M-1）/2;
n=[0:
（M-1）];
m=n-alpha+eps;
hd=sin（wc*m）./（pi*m）;
升级会员 