4.3 利用滤波器处理加有噪声的音频波形
（1） 利用窗函数法设计的低通滤波器处理加有噪声的音频波形
程序参见附录二3-（1）利用窗函数法设计的低通滤波器处理加噪声的音频波形
加噪前音频语音波形的时域图
10.5
幅度(Y )
0-0.5-10.05
0.1
0.20.25时间(t)
加噪前音频波形的频域图0.15
0.3
0.35
200150
幅度(F Y )
1005000
100
200
300
400
500600频率(f)
700
800
900
1000
图4-16 原始音频的时域与频域
10.5
幅度(Y 1)
0-0.5-10.05
0.1
0.20.25时间(t)
加噪声后音频波形的频域图0.15
0.3
0.35
200幅度(F Y 1)
150
100500400
500
600
700
800频率(f)
900
1000
1100
1200
图4-17 加噪声后音频的时域与频域波形
滤波器的增益响应
200
-20-40增益/分贝
-60-80-100-120-140
0200400
600800频率/赫兹
[1**********]0
图4-18 滤波器的增益响应
10.5
幅度(s f )
0-0.5-10.05
0.1
0.20.25时间(t)
滤波后音频波形的频域图0.15
0.3
0.35
200150
100500400
500
600
700
800频率(f)
900
1000
1100
1200
幅度(F s f )
图4-19 滤波后音频的时域与频域波形
从参考程序及以上的四个图像中可以得到如下结论：
①从原始信号波形的频域图可以看到其频率分量主要在500到900Hz 之间，噪声的频率分量主要集中在950Hz ，利用通带截频为800Hz 的低通滤波器可以滤除噪声。对比图4-16和图4-19滤波前后的波形和频谱，可以看到波形得到了重现 ②滤波器的采样频率为22050Hz ，滤波器的阶数为266
③滤波器的通带截频0.8π，阻带截频0.82π，过渡带宽0.02π ④通带衰减为0.019dB ，阻带衰减约为53dB
（2） 利用频率抽样法设计的高通滤波器处理加有噪声的音频波形 程序参见附录二3-（2）
加噪前信号波形的时域图
10.5
幅度(Y )
0-0.5-10
0.05
0.1
0.20.25时间(t)
加噪前信号波形的频域图
0.15
0.3
0.35
0.4
200150
幅度(F Y )
1005000
100
200
300
400
500600频率(f)
700
800
900
1000
图4-20 加噪前信号的时域与频域波形
加噪声后信号波形的时域图
10.5
幅度(Y 1)
0-0.5-10
0.05
0.1
0.20.250.3时间(t)
加噪声后信号波形的频域图
0.15
0.35
0.4
200幅度(F Y 1)
150
1005000
100
200
300
400
500600频率(f)
700
800
900
1000
图4-21 加噪后信号的时域与频域波形
滤波器的增益响应
20
-20增益/分贝
-40
-60
-80
-100
00.10.20.3
0.40.50.6归一化频率
0.70.80.91
图4-22 滤波器的增益响应
信号滤波后时域图
10.5
幅度
0-0.5-10
0.05
0.1
0.15
0.20.25时间/秒
信号滤波后频域图
0.3
0.35
0.4
200150
幅度
1005000
100
200
300
400500频率/赫兹
600
700
800
900
图4-23 信号滤波后的时域图和频域图
从参考程序及以上的四个图像中可以得到如下结论：
①从原始信号波形的频域图可以看到其频率分量主要在500到900Hz 之间，噪声的频率分量主要集中在250Hz ，利用通带截频为300Hz 的低通滤波器可以滤除噪声。对比图4-20和图4-23滤波前后的波形和频谱，可以看到波形得到了重现 ②滤波器的采样频率为22050Hz ，滤波器的阶数为266 ③滤波器的通带截频0.4π，阻带截频0.3π，过渡带宽0.1π
3. 利用滤波器处理加有噪声的音频波形
(1)利用窗函数法设计的低通滤波器处理加有噪声的音频的程序:
f3=950;%所加噪声正弦函数的频率 [Y,fs,bits]=wavread('E:\马世超.wav'); %利用wavread 产生音频的函数及采样频率
L=length(Y);t=0:1/fs:(L-1)/fs;%定义时间的范围及步长 y=0.005*sin(2*pi*f3*t); n1=floor(L/2);%所加噪声 f1=(0:n1)*fs/L;Y=Y(:,1);sound(Y,fs);%输出加噪前音频 Y1=y+Y';%给音频加噪声
FY1=abs(fft(Y1,L));FY=abs(fft(Y,L));sound(Y1,fs);%输出加噪后的音频 figure(1) subplot(211)
plot(t(1:1000),Y(1:1000)); grid on;%加噪前音频的时域图 xlabel('时间(t)');ylabel('幅度(Y)'); title('加噪前音频语音波形的时域图'); axis([0.05 0.35 -1 1]); subplot(212)
plot(f1,FY(1:n1+1)); grid on;%加噪前音频的频域图 xlabel('频率(f)');ylabel('幅度(FY)'); title('加噪前音频波形的频域图'); axis([400 1200 0 200]);
figure(2) subplot(211)
plot(t(1:1000),Y1(1:1000)); grid on;%加噪后音频的时域图 xlabel('时间(t)');ylabel('幅度(Y1)'); title('加噪声后音频波形的时域图'); axis([0.05 0.35 -1 1]); subplot(212)
plot(f1,FY1(1:n1+1)); grid on;%加噪后音频的频域图 xlabel('频率(f)');ylabel('幅度(FY1)'); title('加噪声后音频波形的频域图'); axis([400 1200 0 200]); m=0.03; M=round(8/m); N=M-1;%定义滤波器的阶数 b=fir1(N,0.6); figure(3)
[h,f]=freqz(b,1,512);%滤波器的幅频特性图
plot(f*fs/(2*pi),20*log10(abs(h)))%参数分别是频率与幅值
xlabel('频率/赫兹');ylabel('增益/分贝');title('滤波器的增益响应'); figure(4)
sf=filter(b,1,Y1);%使用filter 函数对信号进行滤波 Fsf=abs(fft(sf,L)); subplot(211)
plot(t(1:1000),sf(1:1000)); grid on;%滤波后音频的时域图 xlabel('时间(t)');ylabel('幅度(sf)'); title('滤波后音频波形的时域图'); axis([0.05 0.35 -1 1]); subplot(212)
plot(f1,Fsf(1:n1+1)); grid on;%滤波后音频的频域图 xlabel('频率(f)');ylabel('幅度(Fsf)');
title('滤波后音频波形的频域图'); axis([400 1200 0 200]);sound(sf,fs);
(2)利用频率抽样法设计的高通滤波器处理加有噪声的信号的程序:
f3=250;%所加噪声正弦函数的频率 [Y,fs,bits]=wavread('E:\马世超.wav'); %利用wavread 产生信号的函数及采样频率
L=length(Y);t=0:1/fs:(L-1)/fs;%定义时间的范围及步长 y=0.025*sin(2*pi*f3*t); n1=floor(L/2);%所加噪声 f1=(0:n1)*fs/L;Y=Y(:,1);sound(Y,fs);%输出加噪前信号 Y1=y+Y';%给信号加噪声
FY1=abs(fft(Y1,L));FY=abs(fft(Y,L));sound(Y1,fs);%输出加噪后的信号 figure(1) subplot(211)
plot(t(1:1000),Y(1:1000)); grid on;%加噪前信号的时域图 xlabel('时间(t)');ylabel('幅度(Y)'); title('加噪前信号波形的时域图'); subplot(212)
plot(f1,FY(1:n1+1)); grid on;%加噪前信号的频域图 xlabel('频率(f)');ylabel('幅度(FY)'); title('加噪前信号波形的频域图'); axis([0 1000 0 200]); figure(2) subplot(211)
plot(t(1:1000),Y1(1:1000)); grid on;%加噪后信号的时域图 xlabel('时间(t)');ylabel('幅度(Y1)'); title('加噪声后信号波形的时域图'); subplot(212)
plot(f1,FY1(1:n1+1)); grid on;%加噪后信号的频域图
xlabel('频率(f)');ylabel('幅度(FY1)'); title('加噪声后信号波形的频域图'); axis([0 1000 0 200]); M=32;%所需频率采样点个数 Wp=0.6*pi;%通带截止频率 m=0:M/2;%阻频带上的采样点 Wm=2*pi*m./(M+1);%阻带截止频率
mtr=ceil(Wp*(M+1)/(2*pi));%向正方向舍入ceil(3.5)=4;ceil(-3.2)=-3; Ad=[Wm>=Wp]; Ad(mtr)=0.28;
Hd=Ad.*exp(-j*0.5*M*Wm);%构造频域采样向量H(k) Hd=[Hd conj(fliplr(Hd(2:M/2+1)))];
%fliplr函数实现矩阵的左右翻转conj 是求复数的共轭 h=real(ifft(Hd));%h（n ）=IDFT[H(k)]
w=linspace(0,pi,1000);%用于产生0,pi 之间的1000点行矢量 H=freqz(h,[1],w);%滤波器的幅频特性图 figure(3)
plot(w/pi,20*log10(abs(H)));grid on;%参数分别是归一化频率与幅值 xlabel('归一化频率');ylabel('增益/分贝');title('滤波器的增益响应'); axis([0 1 -100 20]); figure(4)
sf=filter(h,1,Y1);%使用filter 函数对信号进行滤波 subplot(211)
plot(t,sf)%滤波后的信号图像
xlabel('时间/秒');ylabel('幅度');title('信号滤波后时域图'); axis([0 0.4 -1 1]);%限定图像坐标范围 subplot(212)
plot(f1,Fsf(1:n1+1)); grid on;%滤波后的信号频域图
xlabel('频率/赫兹');ylabel('幅度');title('信号滤波后频域图');
axis([0 1000 0 200]);
4.3 利用滤波器处理加有噪声的音频波形
（1） 利用窗函数法设计的低通滤波器处理加有噪声的音频波形
程序参见附录二3-（1）利用窗函数法设计的低通滤波器处理加噪声的音频波形
加噪前音频语音波形的时域图
10.5
幅度(Y )
0-0.5-10.05
0.1
0.20.25时间(t)
加噪前音频波形的频域图0.15
0.3
0.35
200150
幅度(F Y )
1005000
100
200
300
400
500600频率(f)
700
800
900
1000
图4-16 原始音频的时域与频域
10.5
幅度(Y 1)
0-0.5-10.05
0.1
0.20.25时间(t)
加噪声后音频波形的频域图0.15
0.3
0.35
200幅度(F Y 1)
150
100500400
500
600
700
800频率(f)
900
1000
1100
1200
图4-17 加噪声后音频的时域与频域波形
滤波器的增益响应
200
-20-40增益/分贝
-60-80-100-120-140
0200400
600800频率/赫兹
[1**********]0
图4-18 滤波器的增益响应
10.5
幅度(s f )
0-0.5-10.05
0.1
0.20.25时间(t)
滤波后音频波形的频域图0.15
0.3
0.35
200150
100500400
500
600
700
800频率(f)
900
1000
1100
1200
幅度(F s f )
图4-19 滤波后音频的时域与频域波形
从参考程序及以上的四个图像中可以得到如下结论：
①从原始信号波形的频域图可以看到其频率分量主要在500到900Hz 之间，噪声的频率分量主要集中在950Hz ，利用通带截频为800Hz 的低通滤波器可以滤除噪声。对比图4-16和图4-19滤波前后的波形和频谱，可以看到波形得到了重现 ②滤波器的采样频率为22050Hz ，滤波器的阶数为266
③滤波器的通带截频0.8π，阻带截频0.82π，过渡带宽0.02π ④通带衰减为0.019dB ，阻带衰减约为53dB
（2） 利用频率抽样法设计的高通滤波器处理加有噪声的音频波形 程序参见附录二3-（2）
加噪前信号波形的时域图
10.5
幅度(Y )
0-0.5-10
0.05
0.1
0.20.25时间(t)
加噪前信号波形的频域图
0.15
0.3
0.35
0.4
200150
幅度(F Y )
1005000
100
200
300
400
500600频率(f)
700
800
900
1000
图4-20 加噪前信号的时域与频域波形
加噪声后信号波形的时域图
10.5
幅度(Y 1)
0-0.5-10
0.05
0.1
0.20.250.3时间(t)
加噪声后信号波形的频域图
0.15
0.35
0.4
200幅度(F Y 1)
150
1005000
100
200
300
400
500600频率(f)
700
800
900
1000
图4-21 加噪后信号的时域与频域波形
滤波器的增益响应
20
-20增益/分贝
-40
-60
-80
-100
00.10.20.3
0.40.50.6归一化频率
0.70.80.91
图4-22 滤波器的增益响应
信号滤波后时域图
10.5
幅度
0-0.5-10
0.05
0.1
0.15
0.20.25时间/秒
信号滤波后频域图
0.3
0.35
0.4
200150
幅度
1005000
100
200
300
400500频率/赫兹
600
700
800
900
图4-23 信号滤波后的时域图和频域图
从参考程序及以上的四个图像中可以得到如下结论：
①从原始信号波形的频域图可以看到其频率分量主要在500到900Hz 之间，噪声的频率分量主要集中在250Hz ，利用通带截频为300Hz 的低通滤波器可以滤除噪声。对比图4-20和图4-23滤波前后的波形和频谱，可以看到波形得到了重现 ②滤波器的采样频率为22050Hz ，滤波器的阶数为266 ③滤波器的通带截频0.4π，阻带截频0.3π，过渡带宽0.1π
3. 利用滤波器处理加有噪声的音频波形
(1)利用窗函数法设计的低通滤波器处理加有噪声的音频的程序:
f3=950;%所加噪声正弦函数的频率 [Y,fs,bits]=wavread('E:\马世超.wav'); %利用wavread 产生音频的函数及采样频率
L=length(Y);t=0:1/fs:(L-1)/fs;%定义时间的范围及步长 y=0.005*sin(2*pi*f3*t); n1=floor(L/2);%所加噪声 f1=(0:n1)*fs/L;Y=Y(:,1);sound(Y,fs);%输出加噪前音频 Y1=y+Y';%给音频加噪声
FY1=abs(fft(Y1,L));FY=abs(fft(Y,L));sound(Y1,fs);%输出加噪后的音频 figure(1) subplot(211)
plot(t(1:1000),Y(1:1000)); grid on;%加噪前音频的时域图 xlabel('时间(t)');ylabel('幅度(Y)'); title('加噪前音频语音波形的时域图'); axis([0.05 0.35 -1 1]); subplot(212)
plot(f1,FY(1:n1+1)); grid on;%加噪前音频的频域图 xlabel('频率(f)');ylabel('幅度(FY)'); title('加噪前音频波形的频域图'); axis([400 1200 0 200]);
figure(2) subplot(211)
plot(t(1:1000),Y1(1:1000)); grid on;%加噪后音频的时域图 xlabel('时间(t)');ylabel('幅度(Y1)'); title('加噪声后音频波形的时域图'); axis([0.05 0.35 -1 1]); subplot(212)
plot(f1,FY1(1:n1+1)); grid on;%加噪后音频的频域图 xlabel('频率(f)');ylabel('幅度(FY1)'); title('加噪声后音频波形的频域图'); axis([400 1200 0 200]); m=0.03; M=round(8/m); N=M-1;%定义滤波器的阶数 b=fir1(N,0.6); figure(3)
[h,f]=freqz(b,1,512);%滤波器的幅频特性图
plot(f*fs/(2*pi),20*log10(abs(h)))%参数分别是频率与幅值
xlabel('频率/赫兹');ylabel('增益/分贝');title('滤波器的增益响应'); figure(4)
sf=filter(b,1,Y1);%使用filter 函数对信号进行滤波 Fsf=abs(fft(sf,L)); subplot(211)
plot(t(1:1000),sf(1:1000)); grid on;%滤波后音频的时域图 xlabel('时间(t)');ylabel('幅度(sf)'); title('滤波后音频波形的时域图'); axis([0.05 0.35 -1 1]); subplot(212)
plot(f1,Fsf(1:n1+1)); grid on;%滤波后音频的频域图 xlabel('频率(f)');ylabel('幅度(Fsf)');
title('滤波后音频波形的频域图'); axis([400 1200 0 200]);sound(sf,fs);
(2)利用频率抽样法设计的高通滤波器处理加有噪声的信号的程序:
f3=250;%所加噪声正弦函数的频率 [Y,fs,bits]=wavread('E:\马世超.wav'); %利用wavread 产生信号的函数及采样频率
L=length(Y);t=0:1/fs:(L-1)/fs;%定义时间的范围及步长 y=0.025*sin(2*pi*f3*t); n1=floor(L/2);%所加噪声 f1=(0:n1)*fs/L;Y=Y(:,1);sound(Y,fs);%输出加噪前信号 Y1=y+Y';%给信号加噪声
FY1=abs(fft(Y1,L));FY=abs(fft(Y,L));sound(Y1,fs);%输出加噪后的信号 figure(1) subplot(211)
plot(t(1:1000),Y(1:1000)); grid on;%加噪前信号的时域图 xlabel('时间(t)');ylabel('幅度(Y)'); title('加噪前信号波形的时域图'); subplot(212)
plot(f1,FY(1:n1+1)); grid on;%加噪前信号的频域图 xlabel('频率(f)');ylabel('幅度(FY)'); title('加噪前信号波形的频域图'); axis([0 1000 0 200]); figure(2) subplot(211)
plot(t(1:1000),Y1(1:1000)); grid on;%加噪后信号的时域图 xlabel('时间(t)');ylabel('幅度(Y1)'); title('加噪声后信号波形的时域图'); subplot(212)
plot(f1,FY1(1:n1+1)); grid on;%加噪后信号的频域图
xlabel('频率(f)');ylabel('幅度(FY1)'); title('加噪声后信号波形的频域图'); axis([0 1000 0 200]); M=32;%所需频率采样点个数 Wp=0.6*pi;%通带截止频率 m=0:M/2;%阻频带上的采样点 Wm=2*pi*m./(M+1);%阻带截止频率
mtr=ceil(Wp*(M+1)/(2*pi));%向正方向舍入ceil(3.5)=4;ceil(-3.2)=-3; Ad=[Wm>=Wp]; Ad(mtr)=0.28;
Hd=Ad.*exp(-j*0.5*M*Wm);%构造频域采样向量H(k) Hd=[Hd conj(fliplr(Hd(2:M/2+1)))];
%fliplr函数实现矩阵的左右翻转conj 是求复数的共轭 h=real(ifft(Hd));%h（n ）=IDFT[H(k)]
w=linspace(0,pi,1000);%用于产生0,pi 之间的1000点行矢量 H=freqz(h,[1],w);%滤波器的幅频特性图 figure(3)
plot(w/pi,20*log10(abs(H)));grid on;%参数分别是归一化频率与幅值 xlabel('归一化频率');ylabel('增益/分贝');title('滤波器的增益响应'); axis([0 1 -100 20]); figure(4)
sf=filter(h,1,Y1);%使用filter 函数对信号进行滤波 subplot(211)
plot(t,sf)%滤波后的信号图像
xlabel('时间/秒');ylabel('幅度');title('信号滤波后时域图'); axis([0 0.4 -1 1]);%限定图像坐标范围 subplot(212)
plot(f1,Fsf(1:n1+1)); grid on;%滤波后的信号频域图
xlabel('频率/赫兹');ylabel('幅度');title('信号滤波后频域图');
axis([0 1000 0 200]);