【正文】 1. [9] 薛定宇 , 陳陽泉 . 基于 MATLAB/ Simulink 的系統(tǒng)仿真技術與應用 [M]. 北京：清華大學出版社， 2021. [10] 鄒理和 . 數(shù)字信號處理 [M]. 北京：國防工業(yè)出版社， 1985. [11] 王世一 . 數(shù)字信號處理 .修訂版 [M]. 北京： 北京理工大學出版社， 1997. [12] 侯朝煥等 . 使用 FFT 信號處理技術 [M]. 北京：海洋出版社， 1999. [13] Nussbaumer H J. FastFourier Transform and Convolution Algorithms2nded[M]. 德國： SpringerVerlag， 1982 . [14] An InfraRed Remote Control System Designed for Universal Control[J]. IEEE Transactions on Consumer Electronics， 1999. 31 致謝 本課題在選題及研究過程中得到羅朝明老師的悉心指導。這樣的軟件 雖然 有很 多種，其中最具有代表性的就是 MATLAB。不過只有通過自己不斷的遇到困難并且通過不斷的研究以及解決困難才是我們獲得知識的源泉。 本課題加入噪聲后的聲音與原始的聲音有著明顯的不同，有很大的干擾聲。 wavwrite(z2,Fs,39。信號經(jīng)過低通濾波器（頻域） 39。 %對信號進行低通濾波 figure(9)。 N1=ceil(*pi/Bt)。 %%%%%%%%%4000HZ 的低通濾波器對噪音信號進行濾波 %%%%%%%%% wp1=2*pi*3800/Fs。)。 figure(1)。 [b,a]=butter(Nd,wdc)。s39。 [n,wc]=buttord(wp,ws,Rp,As,39。 N=length(x1) %計算音頻信號的長度 wp=2*3800/Fs。%將帶通濾波后的音樂信號保存在語音播放中 如圖 所示 圖 經(jīng)過帶通濾波器的信號頻譜 將信號 通過 4000HZ 的巴特沃斯低通濾波器后，因為只剩下低頻段的信號了，我們可以明顯發(fā)現(xiàn)其聲音與原始信號的聲音相比聲音較低沉。 sound(z1,Fs)。 plot(ftemp,fftshift(abs(fft(z1))))。 %長度為 N1的矩形窗 Window b1=fir1(N31,wc,Window)。 Bt=wp1ws1。 wp1=2*pi*3205/Fs。其信號進過帶通濾波器程序： [x1,Fs,bits]=wavread(39。 wavwrite(z1,Fs,39。 title(39。,Window)。 N2=N1+mod(N1+1,2)。 ws1=2*pi*4000/Fs。C:\代碼 \39。%將低通濾波后的音樂信號保存在語音播放中。 sound(z1,Fs)。 plot(ftemp,fftshift(abs(fft(z1))))。 %長度為 N1 的矩形窗 Window b1=fir1(N1,wc1,Window)。 wc1=(wp1+ws1)/2/pi。)。 圖 巴特 沃斯低通濾波器的頻率響應 巴特沃斯低通濾波器零極點圖 24 第四章 語音信號的濾波 語音信號的濾波測試 本節(jié)是通過將采集的原始語音信號分別通過已設計好的低通、高通和帶通濾波器，并觀察和分析它們各自的頻譜圖，并利用函數(shù) sound 函數(shù)播放濾波后的語音號，仔細聆聽其中聲音的不同。 [Bz,Az]=bilinear(B,A,Fs)。)。 Rp=3。C:\代碼 \39。也就是說， S 平面的左半平面映射到 Z平面的單位圓內， S 平面的右半平面映射到 Z平面的單位圓外，S平面的虛軸映射到 Z 平面的單位圓上。 當 Ω 1 由 π /T 經(jīng)過 0 變化到 π /T 時， Ω 由 ∞經(jīng)過 0 變化到 +∞，也即映射了整個 jΩ 軸。 如 圖所示 21 圖 帶通濾波器頻率響應 用雙線性法設計 IIR 數(shù)字濾波器 采用非線性頻率壓縮方法，將整個頻率軸上的頻率范圍壓縮到 π/T～ π/T 之間，再用 z=esT轉換到 Z 平面上。 Freqz(b1,1,512)。 N3=ceil(*pi/Bt)。 %帶通濾波器的通帶范圍 ws1=2*pi*3000/fs。高通濾波器的頻率響應 39。,Window)。 N2=N1+mod(N1+1,2)。 ws1=2*pi*4000/Fs。 如 圖所示 圖 低通濾波器頻率響應 高通濾波器設計程序： [x1,Fs,bits]=wavread(39。 freqz(b1,1,512)。 N1=ceil(*pi/Bt)。 wp1=2*pi*3800/Fs。 以下為利用窗函數(shù)法設計低通、高通和帶通濾波器。 (3) 計算濾波器的單位脈沖響應 h（ n） ，它是理想單位脈沖響應和窗函數(shù)的乘積。但是這些公式都是近似的，得出的窗口長度 N 還需要在計算中不斷修正以達到最好的要求。通常在運用中，因為數(shù)字濾波器往往被用來實現(xiàn)選頻操作，所以指標形式通常為在頻域中以分貝值給出的相應響應和幅度響應。 IIR 數(shù)字濾波器可用一個 n 階差分方程表示： ( ) ( ) ( )rky n b x n r a y n k? ? ? ??? （ ） 或用它的 Z域系統(tǒng)函數(shù) : () 1rrkkbzHZ az???? ?? （ ） 對照模擬濾波器的傳遞函數(shù) : 1101()mmmmnnnnb S b S bHS a S a S a????? ? ?? ? ? ? （ ） 可以看出 ,設計數(shù)字濾波器的思路與模擬濾波器相仿 ,設計實質都是尋找一組系數(shù) {b,a},去逼近所求的頻率響應 ,以便在性能上滿足預定的技術要求；不同的是數(shù)字濾波器是在 Z 平面尋找合適的 H(Z), 模擬濾波器則是在 S 平面用數(shù)學逼近法尋找近似的所需特性 H(S)。 由于在運算中會出現(xiàn)舍入處理，會造成誤差不斷累積，偶爾會產(chǎn)生微弱的 寄生振蕩 。 根據(jù)這四類 FIR 濾波器 ,可得到 其 響應的頻率響應的特 征 : )()()( ???? ? ?? jrj eHeH 其中 )(?rH 為振幅響應 ,與之區(qū)別的是 它與幅值特性 )( ?jeH 不同 ,前者 的取值 可正可負 ,而后者只能為非負值。 于是 ,根據(jù)β值的不同 以及 m的奇偶性 ,就 分別 產(chǎn)生了四種類型的線性 FIR相位濾波器 .它們是 : 1型線性相位 FIR 濾波器 : 0,m?? 為奇數(shù) , )(nh 以中心點 ( 1)/2m? 對稱 。這 使各個領域的研究人員可以直觀方便地進行科學研究與工程應用。 Matlab 的出現(xiàn) 為數(shù)字濾波的研究和應用提供了一個很好的平臺。低通、高通、帶通和帶阻濾波器都稱作是經(jīng)典濾波器，它們每一種又有模擬濾波器和數(shù)字濾波器兩種形式。 按脈沖響應來分類 ，數(shù)字濾波器分為 IIR 和 FIR，即無限沖激響應濾波器和有限沖激響應濾波器；其中 IIR 網(wǎng)絡中有反饋回路， FIR 網(wǎng)絡中沒有 反 饋回路。疊加信號的頻譜 39。 z1=y1+fft(x1_high)。 plot(ftemp,fftshift(abs(fft(x1_high))))。,Window)。 N2=N1+mod(N1+1,2)。 ws1=2*pi*4000/Fs。 ylabel(39。)。) subplot(2,1,2)。 xlabel(39。 plot(x1_z)。 %將兩個信號疊加成一個新的信號 —— 加噪聲處理 n=length(x1)。 如圖 和 所示 圖 原信號時域圖 圖 原信號頻譜圖 語音信號的加噪處理及頻譜分析 在 MATLAB 中 ，運用 randn 函數(shù)產(chǎn)生與原信號等長度的隨機噪音信號。原始語音信號的頻譜 39。 ftemp=[N/2:N/21].*Fs/N。)。時間（ ms） 39。 %采樣間隔 ts=t*1000 %化成 ms t1=(0:n11)*ts。)。 首先畫出語音信號的時域波形，然后用函數(shù) FFT 對語音信號進行傅立葉變換，得到信號的頻譜特性。 信號頻譜分析是將信號源發(fā)出的信號強度按頻率順序展開，使其成為頻率的函數(shù)，并考察變化規(guī)律，稱為頻譜分析。 采樣頻率： Fs=44100Hz bits=16 采樣點數(shù)： N=1329408 對原始語音信號進行采集并播放： [x1,Fs,bits]=wavread(39。顯然采樣頻率越 高，計算機攝取的圖片越多，對于原始音頻的還原也越加精確。無論采樣頻率如何，理論上來說采樣的位數(shù)決定了音頻數(shù)據(jù)的最大力度范圍。如圖 所示 1()Xk )()( 21 kXWkX kN? kNW 2()Xk )()( 21 kXWkX kN? 蝶形信號流圖 語音信號的采集及頻譜分析 采樣頻 率 在進行模擬∕數(shù)字信號的轉換過程中，當采樣頻率 大于信號中最好頻率 fmax 的 2 倍時，即： =2fmax，則采樣之后的數(shù)字信號完整的保留了原始信號中的信息，一般實際應用中保證采樣頻率為信號最高頻率的5— 10 倍，采樣定理又稱為耐奎斯特定理。 12) ( ) ( )kNX k X k W X k??（ （ ） （ 5） 求出后半部的表示式 ： )2/( 2/2/ NkrNrkN WW ??? ????????????????? ?? ???????????12/012/022/2)2/(2/2212/012/012/1)2/(2/11)()()()2()()()()2(NrNrrkNkNrNNrNrrkNkNrNkXWrxWrxkNXkXWrxWrxkNX （ ） 看出：后半部的 k 值所對應的 1()Xk， 2()Xk則完全重復了前半部分的 k 值所對應的 1()Xk， 2()Xk的值。12/,0)。 1314x (2 ) x ( )x (2 1) x ( )ll??? ??? 141,0 ?? Nl ，? （ ） ?? ?? ??? ??? 14 0 k)12 2/114 0 k2 2/11 )12(x)2(x)(NllNNllN WlWlkX （ ?? ???? ?? 14 0 k 4/4k 2/14 0 k 4/3 )(x)(x Nl lNNNl lN WlWWl )k(k 42/3 XWX kN?? ）（ 14,1,0 ?? Nk ? （ ） 且： )k()k()k4N( 42/31 XWXX kN??? 14,1,0 ?? Nk ? 式中： ???? 14 0 4/33 )()k( Nl lkNWlxX （ ） ???? 14 0 4/44 )()k( Nl lkNWlxX （ ） 根據(jù)上面同樣的分析 可以得出： 利用四個 /4N 點的 DFT 及兩級蝶形組合運算來計算 N 點 DFT， 相 比只用一次分解蝶形組合方式的計算量又減少了一大半。把兩個 /2N 點 DFT 8 和成為 N 點 DFT 時， 其中 有 /2N 個蝶形運算還需要 /2N 次復數(shù)乘法及2 / 2NN??次 復 數(shù) 加 法 。 再考慮 kNW 的以下性質： kk2/k2 NNNNNN WWWW ???? ）（ 這 樣，把上述各式帶入，就可以將 ()Xk表達為前后兩部分： )k()k()k( 2k1 XWXX N?? k=0,1,…, 12?N 12,1,0k),k()k()2k()2k()2Nk( 212)2k(1 ????????? ? NXWXNXWNXX kNNN ?（ ） 因此，只要求出 0 到 ( /2 1)N ? 區(qū)間的所有 1()Xk， 2()Xk值，即可求出 0到( 1)N? 區(qū)間內的所有 ()Xk值，顯然節(jié)省了運算量。 設序列 ()xn 長度為 N ,且滿足 MN 2? (M 為正