椭圆函数汉宁窗语音滤波.doc

语音信号的数字滤波处理（11）1 绪论1.1设计目的与要求(1)掌握数字信号处理的基本概念，基本理论和基本方法。(2)熟悉离散信号和系统的时域特性。(3)掌握序列快速傅里叶变换方法。(4)学会MATLAB的使用，掌握MATLAB的程序设计方法。(5)掌握利用MATLAB对语音信号进行频谱分析。(6)掌握滤波器的网络结构。(7)掌握MATLAB设计IIR、FIR数字滤波器的方法和对信号进行滤波的方法。1.2 设计内容1.2.1预习题部分(1)设计卷积运算的演示程序：可输入任意两个序列x1(n)、x2(n)，并指定x1(n)为自己的学号,本实验为x1(n)=2,0,0,8,8,4,2,5,0,1,1,0。x2(n)的内容和长度自选。例如x2(n)= 1, 2.43, 6.17,12.93,22.17,32.25,40.88, 45.87, 45.87, 40.88, 32.25,22.17, 12.93, 6.17, 2.43,1.0000。(2)编写程序演示采样定理（时域采样、频谱周期延拓），同时演示采样频率小于2fc时，产生的混叠效应： 对下面连续信号进行采样：， A为幅度因子，a为衰减因子，为模拟角频率，其中n为学号(例如，王墨同学n=23) 要求输入采样频率fs（根据程序处理需要指定范围）后，在时域演示信号波形、采样脉冲及采样后信号；在频域演示不同采样频率下对应信号的频谱。1.2.2设计题部分(1)利用Windows下的录音机或其他软件，选择Windows系统的“叮咚.”（tada.wav），并对该信号进行采样；(2)语音信号的频谱分析，画出采样后语音信号的时域波形和频谱图；(3)产生噪声信号加到语音信号中，得到被污染的语音信号，并回放语音信号；(4)污染信号的频谱分析，画出被污染的语音信号时域波形和频谱； (5)根据有关的频谱特性，采用间接法设计IIR数字滤波器，并画出相应滤波器的幅频图（设计3个IIR滤波器）a. 模拟滤波器类型：椭圆函数滤波器（低通、带通、高通）b总体要求：MATLAB原程序+仿真波形+技术指标 (6)根据有关的频谱特性，采用直接法设计FIR数字滤波器，并画出相应滤波器的幅频图（设计3个FIR滤波器）a. 滤波器类型：hanning窗（低通、带通、高通）b总体要求：MATLAB原程序+仿真波形+技术指标+窗函数 (7)用自己设计的这些滤波器分别对被不同噪声污染的信号进行滤波；(8)分析得到信号的频谱，画出滤波后信号的时域波形和频谱，并对滤波前后的信号进行对比，分析信号的变化；(9)回放语音信号。2 课程设计基本原理2.1设计思路(1)选择来自于C盘的“ir_begin.wav”声音作为语音信号(2)根据不同的滤波器选择不同的噪声信号加到语音信号中，得到被污染的语音信号。(3)分别设计椭圆函数滤波器（低通、带通、高通）以及hanning窗滤波器（低通、带通、高通）对被污染的语音信号滤波，滤掉相应的噪音信号，得到符合要求的语音信号。2.2 设计原理 语音信号的频谱范围主要为800HZ左右，并且在5000HZ左右有一个小信号，因此，在设计低通滤波器时，应把噪声频谱设定在5000HZ以上，这样，通过低通滤波器，即可滤除噪声信号从而还原语音信号；在设计高通滤波器时，应把噪声设定在800HZ以内，以通过高通滤波器滤除低频的噪声信号，从而还原相对频率较高的语音信号；在设计带通滤波器时，可把噪声设计在低于800HZ或高于5000HZ的频谱上，已通过带通滤波器还原带通范围内的语音信号。2.2.1设计IIR数字滤波器的基本思想设计IIR数字滤波器的方法主要有基于冲激响应不变法的IIR数字滤波器设计，基于双线性Z变换法的IIR数字滤波器设计，数字高通、带通及带阻IIR滤波器设计，基于MATLAB函数直接设计IIR数字滤波器。本实验中采用双线性变换法变换的巴特沃思数字滤波器。2.2.1.1椭圆函数低通数字滤波器 (1)选择来自于C盘的“ir_begin.wav”声音作为语音信号（用如下语句调用 x,FS,bits=wavread(C:WindowsMediair_begin.wav)）。 (2)给信号加一个大频率的噪声（取噪声频率远大于语音信号的最大频率），产生污染信号。 (3)设计一个椭圆函数低通滤波器，通带范围包括语音信号，阻带频率设定为小于噪声信号频率。 (4)将设计好的椭圆函数低通滤波器滤除被噪声污染后的语音信号。还原语音信号。2.2.1.2椭圆函数高通数字滤波器 (1)选择来自于C盘的“ir_begin.wav”声音作为语音信号（用如下语句调用x,FS,bits=wavread(C:WindowsMediair_begin.wav)）。 (2)给信号加一个小频率的噪声（取噪声频率远小于语音信号的最小频率），产生污染信号。 (3)设计一个椭圆函数高通滤波器，通带范围包括语音信号，阻带频率设定为大于噪声信号频率。 (4)将设计好的椭圆函数高通滤波器滤除被噪声污染后的语音信号。还原语音信号。2.2.1.3椭圆函数带通数字滤波器 (1)选择来自于C盘的“ir_begin.wav”声音作为语音信号（用如下语句调用x,FS,bits=wavread(C:WindowsMediair_begin.wav)）。(2)给信号加一个小频率或大频率的噪声（取噪声频率远小于语音信号的最小频率或大于语音信号的最大频率），产生污染信号。本实验取小频率的噪声信号。(3)设计一个椭圆函数带通滤波器，通带范围包括语音信号，阻带频率设定为不包括噪声信号频率。(4)将设计好的椭圆函数带通滤波器滤除被噪声污染后的语音信号。还原语音信号。2.2.2设计FIR滤波器的基本思想FIR滤波器通常采用窗函数方法来设计。正确地选择窗函数可以提高设计数字滤波器的性能，或者在满足设计要求的情况下，减小FIR数字滤波器的阶次。常用的窗函数有以下几种：矩形窗(Rectangular window)、三角窗(Triangular window)、汉宁窗(Hanning window)、海明窗(Hamming window)、布拉克曼窗(Bartlett window)、切比雪夫窗(Chebyshev window)、巴特里特窗(Bartlett window)及凯塞窗(Kaiser window)。本实验中选用布拉克曼窗(Bartlett window)设计滤波器。2.2.1.1 Hanning窗低通滤波器 (1)选择来自于C盘的“ir_begin.wav”声音作为语音信号（用如下语句调用x,FS,bits=wavread(C:WindowsMediair_begin.wav)）。(2)给信号加一个大频率的噪声（取噪声频率远大于语音信号的最大频率），产生污染信号。(3) 设计一个Hanning窗低通滤波器，通带范围包括语音信号，阻带频率设定为小于噪声信号频率。(4)将设计好的Hanning窗低通滤波器滤除被噪声污染后的语音信号。还原语音信号。2.2.1.2 Hanning窗高通滤波器 (1)选择来自于C盘的“ir_begin.wav”声音作为语音信号（用如下语句调用x,FS,bits=wavread(C:WindowsMediair_begin.wav)）。(2)给信号加一个小频率的噪声（取噪声频率远小于语音信号的最小频率），产生污染信号。(3)设计一个Hanning窗高通滤波器，通带范围包括语音信号，阻带频率设定为大于噪声信号频率。(4)将设计好的Hanning窗高通滤波器滤除被噪声污染后的语音信号。还原语音信号。2.2.1.3 Hanning窗带通滤波器 (1)选择来自于C盘的“ir_begin.wav”声音作为语音信号（用如下语句调用x,FS,bits=wavread(C:WindowsMediair_begin.wav)）。(2)给信号加一个小频率或大频率的噪声（取噪声频率远小于语音信号的最小频率或大于语音信号的最大频率），产生污染信号。本实验取小频率的噪声信号。(3)设计一个Hanning窗带通滤波器，通带范围包括语音信号，阻带频率设定为不包括噪声信号频率。(4)将设计好的Hanning窗带通滤波器滤除被噪声污染后的语音信号。还原语音信号。3 设计步骤和过程3.1 设计步骤(1)语音信号的采集要求利用windows下的录音机（开始程序附件娱乐录音机，文件属性立即转换8000Hz，8位，单声道）录制一段自己的话音，或者采用Windows自带的声音文件（默认为22050Hz），时间控制在几秒左右。然后在MATLAB软件平台下，利用函数wavread对语音信号进行采样，记住采样频率和采样点数。通过wavread函数的使用，要求理解采样频率、采样位数等概念。wavread函数调用格式：y=wavread(file)，读取file所规定的wav文件，返回采样值放在向量y中。y,fs,nbits=wavread(file)，采样值放在向量y中，fs表示采样频率（Hz），nbits表示采样位数。y=wavread(file,N)，读取前N点的采样值放在向量y中。y=wavread(file,N1,N2)，读取从N1点到N2点的采样值放在向量y中。(2)语音信号的频谱分析要求首先画出语音信号的时域波形；然后对语音信号进行频谱分析，在MATLAB中，可以利用函数fft对信号进行快速付立叶变换，得到信号的频谱特性；从而加深对频谱特性的理解。(3)设计数字滤波器和画出频率响应根据语音信号的特点给出有关滤波器的性能指标，例如：1）低通滤波器性能指标，fp=1000Hz，fc=1200 Hz， As=100dB，Ap=1dB；2）高通滤波器性能指标，fc=2800 Hz，fp=3000 Hz，As=100dB，Ap=1dB；3）带通滤波器性能指标，fp1=1200 Hz，fp2=3000 Hz，fc1=1000 Hz，fc2=3200 Hz，As=100dB，Ap=1dB。(4)回放语音信号在MATLAB中,函数sound可以对声音进行回放。其调用格式：sound(x,fs,bits);可以感觉滤波前后的声音有变化。这里可以通过两种方式进行信号采集，一是在自己电脑C:WINDOWSMedia中有很多系统声音，我选择了自己电脑里的“Windows XP 叮当声”这个声音文件（.wav）。二用电脑自带录音软件(开始所有程序附件娱乐录音机)，或者使用其它专业的录音软件，录制时需要配备录音硬件（如麦克风），为便于比较，需要在安静、干扰小的环境下录音。3.2设计详细过程3.2.1语音信号的采集与频谱分析将电脑C:ir_begin.wav这一声音文件导入MATLAB软件中，绘出采样后的语音信号的时域波形和频谱图。MATLAB程序如下：x,FS,bits=wavread(C:WindowsMediair_begin.wav);x=x(:,1);figure(1);subplot(2,1,1);plot(x);%sound(x,FS,bits); %回放语音title(语音信号时域波形图)y=fft(x,3260);f=(FS/3260)*(1:1630); subplot(2,1,2);plot(f(1:1630),abs(y(1:1630);title(语音信号频谱图);运行程序得到语音信号时域和频谱图形如图3.1所示。图3.1 语音信号时域和频谱图3.2.2 IIR数字滤波器的设计下面分别用椭圆函数数字滤波器（低通、带通、高通）对被噪声污染后的信号进行滤波。3.2.2.1椭圆函数低通滤波器(1)添加高频的噪声信号zs=0.05*cos(2*pi*12000*t/FS)，MATLAB程序如下。噪声信号的时域及幅频图如图3.2所示。 %产生噪声信号并加到语音信号t=0:length(x)-1;zs=0.05*cos(2*pi*12000*t/FS);% zs0=0.05*cos(2*pi*10000*t/FS);figure(2);subplot(2,1,1)plot(zs)title(噪声信号波形);zs1=fft(zs,3260);%sound(zs,FS,bits); %回放噪音 subplot(2,1,2)plot(f(1:1630),abs(zs1(1:1630);title(噪声信号频谱);图3.2噪声信号时域和频谱图(2) 设计椭圆函数低通滤波器，MATLAB程序如下。滤波器图如图3.3所示。%低通滤波fp=9000;fs=10000; rp=1;rs=100; wp=2*pi*fp/FS; ws=2*pi*fs/FS; T=1;Fs1=1; wap=2*tan(wp/2); was=2*tan(ws/2); N,wc=ellipord(wap,was,rp,rs,s); B,A=ellip(N,rp,rs,wc,s); Bz,Az=bilinear(B,A,Fs1); figure(4); h,w=freqz(Bz,Az,1630,Fs1*FS); plot(w,abs(h); title(椭圆低通滤波器); xlabel(频率HZ);ylabel(耗损dB); grid on; 图3.3椭圆函数低通滤波器图3.2.2.2椭圆函数高通滤波器(1)添加低频噪声信号zs0=0.05*cos(2*pi*100*t/FS)，MATLAB程序如下。噪声信号的时域及幅频图如图3.4所示。%产生噪声信号并加到语音信号 t=0:length(x)-1; zs0=0.05*cos(2*pi*100*t/FS); figure(2); subplot(2,1,1); plot(zs0) title(噪声信号波形); zs1=fft(zs0,3260); %sound(zs0,FS,bits); %回放噪音 subplot(2,1,2); plot(f(1:1630),abs(zs1(1:1630); title(噪声信号频谱); 图3.4 噪声信号时域和频谱图(2)设计椭圆函数高通滤波器，MATLAB程序如下。滤波器图如图3.5所示。%高通滤波 fp=300;fs=200;Fs=22050; rp=1;rs=100; wp=2*pi*fp/Fs; ws=2*pi*fs/Fs; T=1;Fs1=1; wap=2*tan(wp/2); was=2*tan(ws/2); N,wc=ellipord(wap,was,rp,rs,s); B,A=ellip(N,rp,rs,wc,high,s); Bz,Az=bilinear(B,A,Fs1); figure(4); h,w=freqz(Bz,Az,1630,Fs1*FS); plot(w,abs(h); title(椭圆高通滤波器); xlabel(频率HZ);ylabel(耗损dB); grid on; 3.5 椭圆函数高通滤波器图3.2.2.3椭圆函数带通滤波器(1)添加低频噪声信号0.02*cos(2*pi*12000*t/FS)，MATLAB程序如下。噪声信号的时域及幅频图如图3.6所示。%产生噪声信号并加到语音信号 t=0:length(x)-1; zs0=0.02*cos(2*pi*12000*t/FS); figure(2); subplot(2,1,1);plot(zs0) title(噪声信号波形); zs1=fft(zs0,3260); %sound(zs,FS,bits); %回放噪音 subplot(2,1,2);plot(f(1:1630),abs(zs1(1:1630); title(噪声信号频谱); 图3.6 噪声信号时域和频谱图(2)设计椭圆函数带通滤波器，MATLAB程序如下。滤波器图如图3.7所示。%带通滤波 fp1=300;fp2=8000;fs1=200;fs2=9000; rp=1;rs=50; wp=2*pi*fp1/FS,2*pi*fp2/FS; ws=2*pi*fs1/FS,2*pi*fs2/FS; T=1;Fs1=1; wap=2*tan(wp/2); was=2*tan(ws/2); N,wc=ellipord(wap,was,rp,rs,s); B,A=ellip(N,rp,rs,wc,s); Bz,Az=bilinear(B,A,Fs1); figure(4); h,w=freqz(Bz,Az,1630,Fs1*FS); plot(w,abs(h); title(椭圆带通滤波器);xlabel(频率HZ);ylabel(耗损dB); grid on; 图3.7 椭圆函数带通滤波器图3.2.3 FIR数字滤波器的设计3.2.3.1 语音处理hanning窗低通滤波器(1)添加噪声信号zs=0.05*cos(2*pi*10000*t/FS)，MATLAB程序如下。噪声信号的时域及幅频图如图3.8所示。%产生噪声信号并加到语音信号t=0:length(x)-1;zs=0.05*cos(2*pi*10000*t/FS);zs0=0.05*cos(2*pi*10*t/FS);figure(2);subplot(2,1,1)plot(zs0)title(噪声信号波形);zs1=fft(zs,3260);%sound(zs,FS,bits); %回放噪音 subplot(2,1,2)plot(f(1:1630),abs(zs1(1:1630);title(噪声信号频谱);图3.8 噪声信号时域和频谱图(2)设计Hanning窗低通滤波器，MATLAB程序如下。滤波器图如图3.9所示。%低通加窗滤波fp=3000,fc=3500; wp=2*pi*fp/FS;ws=2*pi*fc/FS;Bt=ws-wp; N0=ceil(11*pi/Bt); N=N0+mod(N0+1,2);wc=(wp+ws)/2/pi; hn=fir1(N-1,wc,hanning(N); figure(4);H,W=freqz(hn,1,512);plot(W*FS/(2*pi),20*log(abs(H); title(hanning低通窗);图3.9 hanning窗低通滤波器图3.2.3.2 语音处理hanning窗高通滤波器(1)添加噪声信号zs0=0.05*cos(2*pi*100*t/FS)，MATLAB程序如下。噪声信号的时域及幅频图如图3.10所示。%产生噪声信号并加到语音信号t=0:length(x)-1;zs0=0.05*cos(2*pi*100*t/FS);figure(2);subplot(2,1,1)plot(zs0)title(噪声信号波形);zs1=fft(zs0,3260);%sound(zs,FS,bits); %回放噪音 subplot(2,1,2)plot(f(1:1630),abs(zs1(1:1630);title(噪声信号频谱);图3.10 噪声信号时域和频谱图 (2)设计hannig窗高通滤波器，MATLAB程序如下。滤波器图如图3.11所示。%高通加窗滤波fp=600,fc=400; wp=2*pi*fp/FS;ws=2*pi*fc/FS;Bt=wp-ws; N0=ceil(11*pi/Bt); N=N0+mod(N0+1,2);wc=(wp+ws)/2/pi; hn=fir1(N-1,wc,high,hanning(N); figure(4);H,W=freqz(hn,1,512);plot(W*FS/(2*pi),20*log(abs(H); title(hanning高通窗);图3.11 hanning窗高通滤波器图3.2.3.3 语音处理hanning窗带通通滤波器(1)添加噪声信号zs0=0.05*cos(2*pi*100*t/FS)，MATLAB程序如下。噪声信号的时域及幅频图如图3.12所示。%产生噪声信号并加到语音信号t=0:length(x)-1;zs0=0.05*cos(2*pi*100*t/FS);figure(2);subplot(2,1,1)plot(zs0)title(噪声信号波形);zs1=fft(zs0,3260);%sound(zs,FS,bits); %回放噪音 subplot(2,1,2)plot(f(1:1630),abs(zs1(1:1630);title(噪声信号频谱);图3.12 噪声信号的时域及频谱图(2)设计hanning 窗带通滤波器，MATLAB程序如下。滤波器图如图3.13所示。%带通加窗滤波wlp=600*2*pi/FS,wup=6000*2*pi/FS;wls=400*2*pi/FS;wus=7000*2*pi/FS; Bt=wlp-wls; N0=ceil(11*pi/Bt); N=N0+mod(N0+1,2);wc=(wls+wlp)/2/pi,(wus+wup)/2/pi; hn=fir1(N-1,wc,hanning(N); figure(4);H,W=freqz(hn,1,512);plot(W*FS/(2*pi),20*log(abs(H); title(hanning带通窗);图3.13 hanning窗带通滤波器图4 程序调试及运行结果4.1 椭圆函数低通数字滤波器将语音信号和噪音信号叠加可以得到含噪声信号，分析其频谱特性之后用椭圆函数低通滤波器进行滤波。MATLAB程序如下，仿真处理如图4.1所示。x,FS,bits=wavread(C:WindowsMediair_begin.wav);x=x(:,1);figure(1);subplot(2,1,1);plot(x);%sound(x,FS,bits); %回放语音title(语音信号时域波形图)y=fft(x,3260);f=(FS/3260)*(1:1630); subplot(2,1,2);plot(f(1:1630),abs(y(1:1630);title(语音信号频谱图);%产生噪声信号并加到语音信号t=0:length(x)-1;zs=0.05*cos(2*pi*12000*t/FS);% zs0=0.05*cos(2*pi*10000*t/FS);figure(2);subplot(2,1,1)plot(zs)title(噪声信号波形);zs1=fft(zs,3260);%sound(zs,FS,bits); %回放噪音 subplot(2,1,2)plot(f(1:1630),abs(zs1(1:1630);title(噪声信号频谱);x1=x+zs;%sound(x1,FS,bits); %回放加入噪声后的语音y1=fft(x1,3260);figure(3);subplot(2,1,1);plot(x1);title(加入噪声后的信号波形);subplot(2,1,2);plot(f(1:1630),abs(y1(1:1630);title(加入噪声后的信号频谱);%低通滤波fp=9000;fs=10000; rp=1;rs=100; wp=2*pi*fp/FS; ws=2*pi*fs/FS; T=1;Fs1=1; wap=2*tan(wp/2); was=2*tan(ws/2); N,wc=ellipord(wap,was,rp,rs,s); B,A=ellip(N,rp,rs,wc,s); Bz,Az=bilinear(B,A,Fs1); figure(4); h,w=freqz(Bz,Az,1630,Fs1*FS); plot(w,abs(h); title(椭圆低通滤波器); xlabel(频率HZ);ylabel(耗损dB); grid on; yd=filter(Bz,Az,x1);figure(5);subplot(2,1,1);plot(yd); title(滤波后信号波形);ydd=fft(yd,3260);subplot(2,1,2);plot(f(1:1630),abs(ydd(1:1630);title(滤波后信号频谱);sound(yd,FS,bits);图4.1椭圆函数低通滤波器处理噪声图4.2 椭圆函数高通数字滤波器将语音信号和噪音信号叠加可以得到含噪声信号，分析其频谱特性之后用椭圆函数高通滤波器进行语音信号的滤波。MATLAB程序如下，仿真处理如图4.2所示。x,FS,bits=wavread(C:WINDOWSMediair_begin.wav); x=x(:,1); figure(1); subplot(2,1,1); plot(x);%sound(x,FS,bits); %回放语音 title(语音信号时域波形图) y=fft(x,3260); f=(FS/3260)*1:1630; subplot(2,1,2); plot(f(1:1630),abs(y(1:1630); title(语音信号频谱图); %产生噪声信号并加到语音信号 t=0:length(x)-1; zs0=0.05*cos(2*pi*100*t/FS); figure(2); subplot(2,1,1); plot(zs0) title(噪声信号波形); zs1=fft(zs0,3260); %sound(zs0,FS,bits); %回放噪音 subplot(2,1,2); plot(f(1:1630),abs(zs1(1:1630); title(噪声信号频谱); x1=x+zs0; %sound(x1,FS,bits); %回放加入噪声后的语音 y1=fft(x1,3260); figure(3); subplot(2,1,1);plot(x1); title(加入噪声后的信号波形); subplot(2,1,2); plot(f(1:1630),abs(y1(1:1630); title(加入噪声后的信号频谱); %高通滤波 fp=300;fs=200;Fs=22050; rp=1;rs=100; wp=2*pi*fp/Fs; ws=2*pi*fs/Fs; T=1;Fs1=1; wap=2*tan(wp/2); was=2*tan(ws/2); N,wc=ellipord(wap,was,rp,rs,s); B,A=ellip(N,rp,rs,wc,high,s); Bz,Az=bilinear(B,A,Fs1); figure(4); h,w=freqz(Bz,Az,1630,Fs1*FS); plot(w,abs(h); title(椭圆高通滤波器); xlabel(频率HZ);ylabel(耗损dB); grid on; yd=filter(Bz,Az,x1); figure(5); subplot(2,1,1);plot(yd); title(滤波后信号波形); ydd=fft(yd,3260); subplot(2,1,2);plot(f(1:1630),abs(ydd(1:1630); title(滤波后信号频谱); sound(yd,FS,bits); 图4.2 椭圆函数高通滤波器处理噪声图 4.3椭圆函数带通数字滤波器将语音信号和噪音信号叠加可以得到含噪声信号，分析其频谱特性之后用巴特沃思带通滤波器进行语音信号的滤波。MATLAB程序如下，仿真处理如图4.3所示。clc,clear,close all;x,FS,bits=wavread(C:WindowsMediair_begin.wav); x=5*x(:,1); figure(1); subplot(2,1,1); plot(x); %sound(x,FS,bits); %回放语音 title(语音信号时域波形图); y=fft(x,3260); f=(FS/3260)*(1:1630); subplot(2,1,2); plot(f(1:1630),abs(y(1:1630); title(语音信号频谱图);%产生噪声信号并加到语音信号 t=0:length(x)-1; zs0=0.02*cos(2*pi*12000*t/FS); figure(2); subplot(2,1,1);plot(zs0) title(噪声信号波形); zs1=fft(zs0,3260); %sound(zs,FS,bits); %回放噪音 subplot(2,1,2);plot(f(1:1630),abs(zs1(1:1630); title(噪声信号频谱); x1=x+zs0; %sound(x1,FS,bits); %回放加入噪声后的语音 y1=fft(x1,3260); figure(3); subplot(2,1,1);plot(x1); title(加入噪声后的信号波形); subplot(2,1,2); plot(f(1:1630),abs(y1(1:1630); title(加入噪声后的信号频谱); %带通滤波 fp1=300;fp2=8000;fs1=200;fs2=9000; rp=1;rs=50; wp=2*pi*fp1/FS,2*pi*fp2/FS; ws=2*pi*fs1/FS,2*pi*fs2/FS; T=1;Fs1=1; wap=2*tan(wp/2); was=2*tan(ws/2); N,wc=ellipord(wap,was,rp,rs,s); B,A=ellip(N,rp,rs,wc,s); Bz,Az=bilinear(B,A,Fs1); figure(4); h,w=freqz(Bz,Az,1630,Fs1*FS); plot(w,abs(h); title(椭圆带通滤波器);xlabel(频率HZ);ylabel(耗损dB); grid on; yd=filter(Bz,Az,x1); figure(5); subplot(2,1,1);plot(yd); ydd=fft(yd,3260); subplot(2,1,2);plot(f(1:1630),abs(ydd(1:1630);sound(yd,FS,bits);图4.3 椭圆函数带通滤波器处理噪声图4.4 hanning窗低通滤波器将语音信号和噪音信号叠加可以得到含噪声信号，分析其频谱特性之后用hanning 窗低通滤波器进行语音信号的滤波。MATLAB程序如下，仿真处理如图4.4所示。x,FS,bits=wavread(C:WINDOWSMediair_begin.wav);x=x(:,1);figure(1);subplot(2,1,1);plot(x);%sound(x,FS,bits); %回放语音title(语音信号时域波形图)y=fft(x,3260);f=(FS/32600)*1:1630; subplot(2,1,2);plot(f(1:1630),abs(y(1:1630);title(语音信号频谱图);%产生噪声信号并加到语音信号t=0:length(x)-1;zs=0.05*cos(2*pi*10000*t/FS);zs0=0.05*cos(2*pi*10*t/FS);figure(2);subplot(2,1,1)plot(zs0)title(噪声信号波形);zs1=fft(zs,3260);%sound(zs,FS,bits); %回放噪音 subplot(2,1,2)plot(f(1:1630),abs(zs1(1:1630);title(噪声信号频谱);x1=x+zs;%sound(x1,FS,bits); %回放加入噪声后的语音y1=fft(x1,3260);figure(3);subplot(2,1,1);plot(x1);title(加入噪声后的信号波形);subplot(2,1,2);plot(f(1:1630),abs(y1(1:1630);title(加入噪声后的信号频谱);%低通加窗滤波fp=3000,fc=3500; wp=2*pi*fp/FS;ws=2*pi*fc/FS;Bt=ws-wp; N0=ceil(11*pi/Bt); N=N0+mod(N0+1,2);wc=(wp+ws)/2/pi; hn=fir1(N-1,wc,hanning(N); figure(4);H,W=freqz(hn,1,512);plot(W*FS/(2*pi),20*log(abs(H); title(hanning低通窗);X=conv(hn,x); sound(X,FS,bits);X1=fft(X,3260);figure(5);subplot(211);plot(X);title(滤波后的信号波形);subplot(212);plot(f(1:1630),abs(X1(1:1630); title(滤波后的信号频谱)图4.4 hanning窗低通滤波器处理噪声图4.5 hanning窗高通滤波器将语音信号和噪音信号叠加可以得到含噪声信号，分析其频谱特性之后用hanning窗高通滤波器进行语音信号的滤波。MATLAB程序如下，仿真处理如图4.5所示。x,FS,bits=wavread(C:WINDOWSMediair_begin.wav);x=x(:,1);figure(1);subplot(2,1,1);plot(x);%sound(x,FS,bits); %回放语音title(语音信号时域波形图)y=fft(x,3260);f=(FS/3260)*1:1630; subplot(2,1,2);plot(f(1:1630),abs(y(1:1630);title(语音信号频谱图);%产生噪声信号并加到语音信号t=0:length(x)-1;zs0=0.05*cos(2*pi*100*t/FS);figure(2);subplot(2,1,1)plot(zs0)title(噪声信号波形);zs1=fft(zs0,3260);%sound(zs,FS,bits); %回放噪音 subplot(2,1,2)plot(f(1:1630),abs(zs1(1:1630);title(噪声信号频谱);x1=x+zs0;%sound(x1,FS,bits); %回放加入噪声后的语音y1=fft(x1,3260);figure(3);subplot(2,1,1);plot(x1);title(加入噪声后的信号波形);subplot(2,1,2);plot(f(1:1630),abs(y1(1:1630);title(加入噪声后的信号频谱);%高通加窗滤波fp=600,fc=400; wp=2*pi*fp/FS;ws=2*pi*fc/FS;Bt=wp-ws; N0=ceil(11*pi/Bt); N=N0+mod(N0+1,2);wc=(wp+ws)/2/pi; hn=fir1(N-1,wc,high,hanning(N); figure(4);H,W=freqz(hn,1,512);plot(W*FS/(2*pi),20*log(abs(H); title(hanning高通窗);X=conv(hn,x); sound(X,FS,bits);X1=fft(X,3260);figure(5);subplot(211);plot(X);title(滤波后的信号波形);subplot(212);plot(f(1:1630),abs(X1(1:1630); title(滤波后的信号频谱)图4.5 hanning窗高通滤波器处理噪声图 4.6 hanning窗带通滤波器将语音信号和噪音信号叠加可以得到含噪声信号，分析其频谱特性之后用hanning窗带通滤波器进行语音信号的滤波。MATLAB程序如下，仿真处理如图4.6所示。Fs=22050;x,FS,bits=wavread(C:WINDOWSMediair_begin.wav);x=x(:,1);figure(1);subplot(2,1,1);plot(x);%sound(x,FS,bits); %回放语音title(语音信号时域波形图)y=fft(x,3260);f=(FS/3260)*1:1630; subplot(2,1,2);plot(f(1:1630),abs(y(1:1630);title(语音信号频谱图);%产生噪声信号并加到语音信号t=0:length(x)-1;zs0=0.05*cos(2*pi*100*t/FS);figure(2);subplot(2,1,1)plot(zs0)title(噪声信号波形);zs1=fft(zs0,3260);%sound(zs,FS,bits); %回放噪音 subplot(2,1,2)plot(f(1:1630),abs(zs1(1:1630);title(噪声信号频谱);x1=x+zs0;%sound(x1,FS,bits); %回放加入噪声后的语音y1=fft(x1,3260);figure(3);subplot(2,1,1);plot(x1);title(加入噪声后的信号波形);subplot(2,1,2);plot(f(1:1630),abs(y1(1:1630);title(加入噪声后的信号频谱);%带通加窗滤波wlp=600*2*pi/FS,wup=6000*2*pi/FS;wl