語(yǔ)音信號(hào)的數(shù)字濾波處理(五)

語(yǔ)音信號(hào)的數(shù)字濾波處理圖1所示，程序如下：n=0:50;A=10;a=10*sqrt(2)*pi;w0=10*sqrt(2)*pi;t0=0.001;t1=0.020;t2=0.005;x=A*exp(-a*n*t0).*sin(w0*n*t0);y1=A*exp(-a*n*t1).*sin(w0*n*t1);y2=A*exp(-a*n*t2).*sin(w0*n*t2);subplot(321);stem(n,x,'r'),gridontitle('離散時(shí)間信號(hào)')xlabel('n');ylabel('幅度')subplot(322);plot(n,x,'g'),gridontitle('連續(xù)時(shí)間信號(hào)')xlabel('t');ylabel('幅度')subplot(323);stem(n,y1),gridontitle('50Hz連續(xù)時(shí)間信號(hào)')xlabel('n');ylabel('幅度')subplot(324);stem(n,y2,'m'),gridontitle('200Hz理想采樣信號(hào)序列');xlabel('n');ylabel('幅度')subplot(325);stem(n,x,'k'),gridontitle('1000Hz連續(xù)時(shí)間信號(hào)')xlabel('n');ylabel('幅度')t=0:0.001:1;xx=A*exp(-a*n*t0).*sin(w0*n*t0);X=fft(xx);subplot(326)plot(abs(X));xlabel('f/Hz');ylabel('|X(ejw)|')gridon；圖2.4采樣定理圖示3課程設(shè)計(jì)的基本原理3.1IIR濾波器的設(shè)計(jì)原理IIR數(shù)字濾波器的設(shè)計(jì)一般是利用目前已經(jīng)很成熟的模擬濾波器的設(shè)計(jì)方法來(lái)進(jìn)行設(shè)計(jì)，通常采用模擬濾波器原型有巴特沃斯濾波器、切比雪夫?yàn)V波器、橢圓濾波器等。IIR數(shù)字濾波器單位取樣響應(yīng)是時(shí)寬無(wú)限的h(n)，即n0＜n＜∞。數(shù)字濾波器按照實(shí)現(xiàn)的方法和結(jié)構(gòu)形式分為遞歸型或非遞歸型兩類。遞歸型數(shù)字濾波器的當(dāng)前輸出y(n)是輸入x(n)的當(dāng)前值和以前各輸入值x(n)，x(n–1)，….，及以前各輸出值y(n)，y(n–1)，….的函數(shù)。一個(gè)N階遞歸型數(shù)字濾波器(IIR濾波器)的差分方程如公式2-1所示：2-1其中，式中的系數(shù)至少有一項(xiàng)不為零。說(shuō)明必須將延時(shí)的輸出序列進(jìn)行反饋。遞歸系統(tǒng)的傳統(tǒng)函數(shù)定義如公式2-2所示：2-2遞歸系統(tǒng)的傳遞函數(shù)H(z)在Z平面上不僅有零點(diǎn)，而且有極點(diǎn)。非遞歸型數(shù)字濾波器當(dāng)前的輸出值y(n)僅為當(dāng)前及以前的輸入序列的函數(shù)，而與以前的各個(gè)輸出值無(wú)關(guān)，因此從結(jié)構(gòu)上看非遞歸系統(tǒng)沒(méi)有反饋環(huán)路。3.2FIR濾波器的設(shè)計(jì)原理FIR數(shù)字濾波器的單位取樣響應(yīng)是時(shí)寬有限的h(n)，即n1≤n＜n2。一個(gè)N階的非遞歸型數(shù)字濾波器(FIR濾波器)的差分方程如公式2-3所示：2-3差分方程式中的系數(shù)等于單位取樣響應(yīng)的序列值h(n)，其系統(tǒng)函數(shù)H(z)可以表示為公式2-4所示的形式：2-4H(z)是的多項(xiàng)式，因此它的極點(diǎn)只能在Z平面的原點(diǎn)上。FIR濾波器通常采用窗函數(shù)方法來(lái)設(shè)計(jì)。窗設(shè)計(jì)的基本思想是，首先選擇一個(gè)適當(dāng)?shù)睦硐脒x頻濾波器（它總是具有一個(gè)非因果，無(wú)限持續(xù)時(shí)間脈沖響應(yīng)），然后街區(qū)（加窗）它的脈沖響應(yīng)得到線性相位和因果FIR濾波器。我們用表示理想的選頻濾波器，它在通帶上具有單位增益和線性相位，在阻帶上具有零響應(yīng)。一個(gè)帶寬wc<pi的低通濾波器由下式給定：為了從hd(n)得到一個(gè)FIR濾波器，必須同時(shí)在兩邊截取hd（n）。而要得到一個(gè)因果的線性相位濾波器，它的h(n)長(zhǎng)度為N，必須滿足公式2-5：2-5h(n)可以看做是hd(n)與窗函數(shù)w（n）的乘積，如公式2-6所示：2-6這種操作叫做加窗，其中根據(jù)w（n）的不同定義，可以得到不同的窗結(jié)構(gòu)。在頻域中，因果FIR濾波器響應(yīng)由和窗響應(yīng)的周期卷積得到，如公式2-7所示:2-7常用的窗函數(shù)有矩形窗、漢寧（HANNING）窗、海明（矩形）窗、布萊克曼（BLACKMAN）窗、凱澤（KAISER）窗等。4設(shè)計(jì)程序的調(diào)試及運(yùn)行結(jié)果4.1IIR濾波器的設(shè)計(jì)(1)選擇來(lái)自于D盤(pán)的“ding.wav”聲音作為語(yǔ)音信號(hào)(2)根據(jù)不同的濾波器選擇不同的噪聲信號(hào)加到語(yǔ)音信號(hào)中，得到被污染的語(yǔ)音信號(hào)。(3)分別設(shè)計(jì)巴特沃思濾波器（低通、帶通、高通）以及Kaiser窗濾波器（低通、帶通、高通）對(duì)被污染的語(yǔ)音信號(hào)濾波，濾掉相應(yīng)的噪音信號(hào)，得到符合要求的語(yǔ)音信號(hào)。方案設(shè)計(jì)流程圖如圖4.1所示。（圖為巴特沃思濾波器，Kaiser窗濾波器與巴特沃思類似）4.1總體方案設(shè)計(jì)流圖4.11巴特沃斯低通濾波器的設(shè)計(jì)(1)添加高頻的噪聲信號(hào)zs=0.05*cos(2*pi*10000*t/22050)，MATLAB程序如下。噪聲信號(hào)的時(shí)域及幅頻圖如圖3.2所示。%噪聲信號(hào)zt=0:length(x)-1;zs=0.05*cos(2*pi*10000*t/22050);zs0=0.05*cos(2*pi*10000*t/22050000);figure(2);subplot(2,1,1)plot(zs0)title('噪聲信號(hào)波形');zs1=fft(zs,1200);%sound(zs,FS,bits);%回放噪音subplot(2,1,2)plot(f(1:600),abs(zs1(1:600)));title('噪聲信號(hào)頻譜');圖4.2噪聲信號(hào)時(shí)域和頻譜圖(2)設(shè)計(jì)巴特沃思低通濾波器，MATLAB程序如下。濾波器圖如圖3.3所示。%低通濾波fp=3000;fs=3500;Fs=22050;rp=1;rs=10;wp=2*pi*fp/Fs;ws=2*pi*fs/Fs;Fs1=1;wap=2*tan(wp/2);was=2*tan(ws/2);[N,wc]=buttord(wap,was,rp,rs,'s');[B,A]=butter(N,wc,'s');[Bz,Az]=bilinear(B,A,Fs1);figure(4);[h,w]=freqz(Bz,Az,512,Fs1*22050);plot(w,abs(h));title('巴特沃斯低通濾波器');xlabel('頻率（HZ）');ylabel('耗損（dB）');gridon;圖4.3巴特沃思低通濾波器圖4.12巴特沃斯高通濾波器的設(shè)計(jì)(1)添加低頻噪聲信號(hào)zs0=0.05*cos(2*pi*100*t/22050)，MATLAB程序如下。噪聲信號(hào)的時(shí)域及幅頻圖如圖3.4所示。t=0:length(x)-1;zs0=0.05*cos(2*pi*100*t/22050);figure(2);subplot(2,1,1)plot(zs0)title('噪聲信號(hào)波形');zs1=fft(zs0,1200);%sound(zs,FS,bits);%回放噪音subplot(2,1,2)plot(f(1:600),abs(zs1(1:600)));title('噪聲信號(hào)頻譜');圖4.4噪聲信號(hào)時(shí)域和頻譜圖(2)設(shè)計(jì)巴特沃思高通濾波器，MATLAB程序如下。濾波器圖如圖3.5所示。%高通濾波fp=600;fs=400;Fs=22050;rp=1;rs=10;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]=buttord(wap,was,rp,rs,'s');[B,A]=butter(N,wc,'high','s');[Bz,Az]=bilinear(B,A,Fs1);figure(4);[h,w]=freqz(Bz,Az,512,Fs1*22050);plot(w,abs(h));title('巴特沃斯高通濾波器');xlabel('頻率（HZ）');ylabel('耗損（dB）');gridon;圖4.5巴特沃思高通濾波器圖4.13巴特沃斯高通濾波器的設(shè)計(jì)(1)添加低頻噪聲信號(hào)zs1=0.05*cos(2*pi*100*t/22050)，MATLAB程序如下。噪聲信號(hào)的時(shí)域及幅頻圖如圖3.6所示。t=0:length(x)-1;zs0=0.05*cos(2*pi*100*t/22050);figure(2);subplot(2,1,1)plot(zs0)title('噪聲信號(hào)波形');zs1=fft(zs0,1200);%sound(zs,FS,bits);%回放噪音subplot(2,1,2)plot(f(1:600),abs(zs1(1:600)));title('噪聲信號(hào)頻譜');圖4.6噪聲信號(hào)時(shí)域和頻譜圖(2)設(shè)計(jì)巴特沃思帶通濾波器，MATLAB程序如下。濾波器圖如圖3.7所示。%帶通濾波fp=[600,6000];fs=[400,7000];Fs=22050;rp=1;rs=10;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]=buttord(wap,was,rp,rs,'s');[B,A]=butter(N,wc,'s');[Bz,Az]=bilinear(B,A,Fs1);figure(4);[h,w]=freqz(Bz,Az,512,Fs1*22050);plot(w,abs(h));title('巴特沃斯帶通濾波器');xlabel('頻率（HZ）');ylabel('耗損（dB）');gridon;圖4.7巴特沃思帶通濾波器圖4.2FIR數(shù)字濾波器的設(shè)計(jì)4.21Kaiser窗低通程序的設(shè)計(jì)技術(shù)指標(biāo)：fr=1000。源程序如下。仿真波形見(jiàn)圖4.21。[y,fs,bit]=wavread('1');fr=1000;wr=2*pi*fr;Window=kaiser(8,beta);b=fir1(7,fr/(fs/2),Window);[h,w]=freqz(b,1);subplot(221);plot(w*fs/(2*pi),abs(h));title('低通濾波');y1=filter(b,1,y);y2=fftfilt(b,y);subplot(222);plot(y);title('原始語(yǔ)音信號(hào)');subplot(223);plot(y1);title('IIR濾波后語(yǔ)音信號(hào)');subplot(224);plot(y2);title('FIR濾波后語(yǔ)音信號(hào)');sound(y1,fs,bit);圖4.21凱塞窗低通濾波器仿真圖4.22Kaiser窗帶通濾波器技術(shù)指標(biāo)：fr=1000,fp=3000。源程序如下。仿真波形見(jiàn)圖4.22。[y,fs,bit]=wavread('1');fr=1000;fp=3000wr=2*pi*fr;Window=kaiser(16,beta);b=fir1(15,[fr/(fs/2),fp/(fs/2)],Window);[h,w]=freqz(b,1);subplot(221);plot(w*fs/(2*pi),abs(h));title('帶通濾波');y1=filter(b,1,y);y2=fftfilt(b,y);subplot(222);plot(y);title('原始語(yǔ)音信號(hào)');subplot(223);plot(y1);title('IIR濾波后語(yǔ)音信號(hào)');subplot(224);plot(y2);title('FIR濾波后語(yǔ)音信號(hào)');sound(y1,fs,bit);圖4.22凱塞窗帶通濾波器仿真波形4.23Kaiser窗高通濾波器指標(biāo)：fr=3000。源程序如下。仿真波形見(jiàn)圖4.23。[y,fs,bit]=wavread('1');fr=3000;wr=2*pi*fr;Window=kaiser(8,beta);b=fir1(7,fr/(fs/2),'high');[h,w]=freqz(b,1);subplot(221);plot(w*fs/(2*pi),abs(h));title('低通濾波');y1=filter(b,1,y);y2=fftfilt(b,y);subplot(222);plot(y);title('原始語(yǔ)音信號(hào)');subplot(223);plot(y1);title('IIR濾波后語(yǔ)音信號(hào)');subplot(224);plot(y2);title('FIR濾波后語(yǔ)音信號(hào)');sound(y1,fs,bit);圖4.23凱塞窗高通濾波器仿真波形4.3回放語(yǔ)音信號(hào)經(jīng)過(guò)以上的加噪聲處理后，可在Matlab中用函數(shù)sound對(duì)聲音進(jìn)行回放。其調(diào)用格式：sound(y,Fs)，sound(y)和sound(y,Fs,bits)。可以察覺(jué)濾波前后的聲音有明顯的變化。5心得與體會(huì)在本次課程設(shè)計(jì)的過(guò)程中，主要遇到的問(wèn)題是原始音樂(lè)信號(hào)的選取。有些音樂(lè)信號(hào)在MATLAB中運(yùn)行時(shí)出錯(cuò)，開(kāi)始以為自己的操作問(wèn)題，在跟同學(xué)商討多次更換音樂(lè)信號(hào)后得到了解決。其次就是在濾波器參數(shù)的選取上，開(kāi)始時(shí)不清楚怎么選取，又在課本中找了好久才找到。學(xué)習(xí)的過(guò)程是相互討論共同進(jìn)步的，多多討論課題中遇到的問(wèn)題，可以鞏固我們的知識(shí)掌握能力，增加熟練運(yùn)用度。從本次課程設(shè)計(jì)的中心來(lái)看，課題是希望將數(shù)字信號(hào)處理技術(shù)應(yīng)用于某一實(shí)際