南郵通達語音時頻參數(shù)的提取和應用和語音編碼實驗.doc

資料收集于網(wǎng)絡 如有侵權請聯(lián)系網(wǎng)站 刪除 謝謝 南京郵電大學通達學院實 驗 報 告實驗名稱 語音時頻參數(shù)的提取和應用實驗 語音編碼實驗 課程名稱 信息處理技術專業(yè)綜合實驗 班級學號 姓 名 開課時間 2016/2017 學年, 第 二 學期 實驗一 語音時頻參數(shù)的提取和應用一、實驗目的熟悉基本的語音時頻參數(shù)提取方法，觀察并比較各類參數(shù)采用不同提取方法的效果，了解這些時頻參數(shù)在語音處理中的應用。二、實驗內容1、編寫并調試語音時域參數(shù)提取程序，包括短時功率、短時過零率等函數(shù)。2、編寫并調試語音頻域參數(shù)提取程序，包括頻譜、語譜圖、基音頻率、共振峰等。3、利用提取出的時頻參數(shù)對語音信號進行靜/清/濁音的劃分。三、實驗原理語音信號具有短時平穩(wěn)性，需對其進行分段（幀）處理（幀長一般取10至30ms之間）。根據(jù)激勵方式的不同，語音信號可以分成濁音和清音兩大類：濁音：激勵源可模擬成準周期性的脈沖串。該周期稱為基音周期，其倒數(shù)就是基音頻率，基音頻率一般在幾百赫茲以內。清音：激勵源可模擬成隨機白噪聲。簡化的語音生成數(shù)學模型如圖所示：對語音時頻波形進行繪制與觀察，有助于了解語音的清/濁音等特性。對基音周期等聲門參數(shù)和共振峰等聲道參數(shù)進行提取與分析，是語音處理與應用的重要內容。 對于時變語音信號，可以計算它的短時功率和短時過零率。第m幀（幀長為N個樣點）的功率計算公式為：過零率計算公式為：對于濁音語音，可以利用其頻譜具有豐富的諧波分量的特點，求出其諧波乘積譜，計算公式為：式中，R一般取為5。在諧波乘積譜中，基頻分量變得很大，更易于估計基音周期。四、實驗方法及程序1. 調用matlab中的wavread和wavplay讀入并播放語音文件。2. 調用fft命令計算語音的頻譜。3. 調用specgram命令得到語音的語譜圖。4. 調用plot命令繪制并觀察語音的時頻域波形。5. 根據(jù)短時功率、短時過零率、諧波乘積譜的計算公式，編寫分別用以計算短時功率、短時過零率、基音頻率的函數(shù)：function Px=stpower、function Zx = stzerocross(x, N) 和function HPSx = hpspectrum(x,N,R)。6. 利用短時功率和短時過零率這兩個參數(shù)可以對語音信號進行濁音/清音分類。編寫實現(xiàn)此分類功能的函數(shù)voiunvoi.m。五、實驗結果與分析1. 清音和濁音的短時功率、短時過零率各有何特點？清音段的能量一般比濁音段的小得多。濁音時能量集中于較低頻率段內，具有較低的過零率，而清音時能量集中于較高頻率段內，具有較高的過零率。2. 如何利用短時功率和短時過零率這兩個參數(shù)對語音信號進行濁音/清音的分類？ 發(fā)濁音時，由于聲門波引起了譜的高頻跌落。所以其語音能量集中在3kHz以下；發(fā)清音時，多數(shù)能量出現(xiàn)在高頻。高頻率意味著高的平均過零率，低頻率意味著低的平均過零率，一般可認為濁音具有較低的平均過零率，清音具有較高的平均過零率。3. 如何根據(jù)時域波形估計各幀元音的基音周期？ 看圖可知四幀語音的周期均為0.01s，第一幀0100Hz有10個峰值點，基頻約為10Hz，第二幀0100Hz有5個峰值點，基頻約為20Hz，第三幀0100Hz有5個峰值點，基頻約為20Hz，第四幀0100Hz有4個峰值點，基頻約為25Hz。4. 如何根據(jù)基于DFT的對數(shù)幅度譜估計某幀濁音的基音周期？01172為一個基因周期。5. 如何根據(jù)基于DFT的對數(shù)幅度譜，估計出共振峰頻率？有7個峰點，因此基頻約為167。6. 時域對語音信號進行加窗，反映在頻域，其窗譜對基于DFT的對數(shù)幅度譜有何影響？如何估計出窗譜的主瓣寬度？因為時間窗幅度的傅立葉變換為脈沖狀態(tài),則頻譜偏差很小,當定常過程x(n)的平均值為零且N很大時 ETN()=f()(12) 所以加窗的周期圖TN()是頻譜密度函數(shù)的無偏估計。固有TN()=1+()f()。如果原始信號的頻譜成份與FFT中的譜線完全一致，這種情況下采樣數(shù)據(jù)的長度為信號周期的整數(shù)倍，頻譜中只有主瓣。7. 如何根據(jù)窄帶、寬帶語譜圖，提取出語音的基音和共振峰軌跡？在基音提取中，廣泛采用語音波形或誤差信號波形的低通濾波。提取共振峰特性最簡便的手段是使用語譜儀。提取共振峰還有倒譜法、LPC分析法等更為有效、準確的方法。（做完試驗后，記住要回答這些問題）實驗過程記錄與結果分析1、進入matlab1.1 在Windows環(huán)境中，創(chuàng)建一個名為speech的文件夾，存放語音數(shù)據(jù)和與實驗相關的Matlab文件。如：E:speech。用于實驗的語音數(shù)據(jù)（.mat）及相應說明（.txt）包括：digits.mat % 英文數(shù)字 “0” 到 “9” 的發(fā)音digits.txt gliss.mat % 兩個包含/i/的滑音gliss.txt letters.mat % 英文字母表中26個字母的發(fā)音letters.txt ma1.txt ma1-1.mat % 一句連續(xù)語音ma1-2.mat % 一句連續(xù)語音ma1-3.mat % 一句連續(xù)語音ma1-4.mat % 一句連續(xù)語音timit.txttimit1.mat % 一句連續(xù)語音timit2.mat % 一句連續(xù)語音timit3.mat % 一句連續(xù)語音timit4.mat % 一句連續(xù)語音vowels.mat % 元音/a/, /i/, /o/, /u/的發(fā)音vowels.txtwords.mat % 十個孤立字的發(fā)音words.txt1.2 在Windows桌面上，雙擊Matlab的圖標，進入Matlab的工作環(huán)境。在命令窗口提示符后，鍵入cd E:speech指令，進入實驗文件夾路徑。1.3 在命令窗口鍵入clear all指令，清除工作空間的所有變量。2、加載語音數(shù)據(jù)，熟悉matlab命令2.1加載“timit1.mat”語音數(shù)據(jù)，命令為：load timit12.2用命令“who”和“whos”列出當前工作空間中的變量。變量“timit1”會以包含55911個元素的矩陣（列向量）形式出現(xiàn)Your variables are:timit1 Name Size Bytes Class Attributes timit1 55911x1 447288 double 2.3利用命令“l(fā)ength”和“size”查看語音信號timit1的長度和維數(shù)：m = length(timit1)m,n = size(timit1)（粘貼實驗結果）m = 55911m = 55911n = 12.4通過以下命令，分別查看語音信號timit1中間的1個、2個和20個數(shù)據(jù)值： timit1(5001)timit1(5001,5003)timit1(5001:5020)（粘貼實驗結果）ans = -0.0646ans = -0.0646 0.0111ans = -0.0646 0.0307 0.0111 -0.0200 0.2021 0.4012 0.1869 -0.1756 -0.1902 0.0823 0.2363 0.1712 0.0637 -0.0209 -0.0683 -0.0394 0.0081 0.0105 0.0311 0.06962.5 Matlab中，符號“”可以實現(xiàn)轉置運算，運行以下命令進行觀察：size (timit1)（粘貼實驗結果）ans = 1 559112.6使用函數(shù)max和min可以分別得到信號的最大值和最小值：max(timit1)min(timit1)（粘貼實驗結果）ans = 1ans = -0.80683、語音信號的時域分析3.1使用以下命令，畫出語音信號timit1（采樣頻率為16kHz，長度為3.5秒）的時域波形。plot (timit1)title (timit1)xlabel (Time Index,it n)ylabel (Amplitude)axis (0 length(timit1) -1 1)grid3.2分別畫出一幀濁音和一幀清音的語音時域波形（采樣頻率為16kHz，幀長為25ms，每幀有400個樣點）。subplot(211), plot(timit1(14501:14900) % /a/ in dark.subplot(212), plot(timit1(35501:35900) % /s/ in wash.3.3編寫個一用以計算信號短時功率的函數(shù)“stpower.m”，其Matlab代碼如下：function Px=stpower (x, N)M = length (x);Px = zeros (M, 1);Px (N)=x (1:N)*x (1:N)/N;For (m = (N+1): M)Px (m) = Px (m-1)+ (x (m) 2-x (m-N) 2)/N;3.4編寫一個用以計算信號短時過零率的函數(shù)“stzerocross.m”，其Matlab代碼如下：function Zx = stzerocross(x, N)M = length(x);Zx = zeros(M,1);Zx (N+1)=sum (abs (sign (x (2:N+1)sign (x (1:N) ) ) )/ (2*N);for (m= (N+2):M) Zx (m) = Zx (m-1) + (abs (sign (x (m) ) sign (x (m-1) ) ) abs(sign (x (m-N) ) sign (x (m-N-1 ) ) ) )/ (2*N);end3.5加載“digits.mat”語音數(shù)據(jù)，該數(shù)據(jù)包含英文數(shù)字 “0” 到 “9”的發(fā)音。計算其中的單詞“four”的短時功率和過零率（采樣頻率為10kHz，幀長為30ms，每幀有300個樣點）。load digits;N = 300; x=digits.four1;Px = stpower (x,N);Zx = stzerocross (x,N);plot(Px*1e -5 Zx x/2000)3.6加載“vowels.mat”語音數(shù)據(jù)，該數(shù)據(jù)包含元音/a/, /i/, /o/, /u/的發(fā)音。分別畫出一幀/a/、一幀/i/、一幀/o/ 和一幀/u/的時域波形（采樣頻率為10kHz，幀長為30ms，每幀有300個樣點）。load vowelssubplot(221)plot (vowels.a_1(2001:2300)subplot(222)plot (vowels.i_1(2001:2300)subplot(223)plot (vowels.o_1(2001:2300)subplot(224)plot (vowels.u_1(2001:2300)3.7編寫一個利用短時功率和短時過零率這兩個參數(shù)對語音信號進行濁音/清音分類的函數(shù)“voiunvoi.m”，其Matlab代碼如下：function voi = voiunvoi(x,N,Pth,Zth)% Short-time power Px and zero crossing Zx measures.Px = stpower(x,N);Zx = stzerocross(x,N);% Compare estimates with threshold values.voi = (PxPth*max(Px) & (ZxZth);% Shift the voi-flag N/2 samples to the left (middle sample in window).voi = voi(fix(N/2)+1:length(voi);voi(length(voi)*ones(fix(N/2),1);4、語音信號的頻域分析4.1加載“ma1_1”語音數(shù)據(jù)。基于DFT變換，畫出其中一幀數(shù)據(jù)（采樣頻率為8kHz，幀長為37.5ms，每幀有300個樣點）的頻域波形（對數(shù)幅度譜）。load ma1_1;x = ma1_1 (4161:4460); subplot(121), plot(x) N = 1024; k = 0:N/2-1;X = fftshift (fft (x.*hann (length (x),N);subplot(122), plot (k,20*log10 (abs(X(N/2:-1:1), axis (0 N/2-1 -inf inf )4.2編寫一個用以計算信號諧波乘積譜的函數(shù)“hpspectrum.m”，其Matlab代碼如下：function HPSx = hpspectrum(x,N,R)k = 1:R:N/2; K = length(k);X = fft (x.*hann(length(x), N);HPSx = X(k);for (r = R-1:-1:1) HPSx = HPSx.*X (1:r:r*K);end4.3利用函數(shù)“hpspectrum.m”，來估計一幀濁音的基音周期。代碼如下：X = ma1_1(4161:4460);N = 1024; R = 5;HPSx = hpspectrum (x, N, R);plot (20*log10 (abs (HPSx)4.4加載“vowels.mat”語音數(shù)據(jù)，分別畫出一幀/i/和一幀/u/（采樣頻率為10kHz，幀長為30ms，每幀有300個樣點）的基于DFT的對數(shù)幅度譜。其Matlab代碼如下：load vowelsx = vowels.i_1(2001:2300);N = 1024; k= -N/2:N/2-1;X = fftshift (fft (x.*hann (length(x),N);plot (k,20*log10(abs(X), axis(0 fix (N/2) 0 100)x = vowels.u_1(2001:2300);N= 1024; k = - N/2:N/2-1;X = fftshift (fft (x.*hann(length(x),N);plot (k,20*log10(abs(X), axis(0 fix(N/2) 0 100)4.5畫出一幀清音語音的基于DFT的對數(shù)幅度譜。語音數(shù)據(jù)為ma1_1中的第1570115860個樣點（采樣頻率為8kHz，幀長為20ms，每幀有160個樣點）。load ma1_1;x = ma1-1 (4161:4460); plot (x)N = 1024; k = - N/2;N/2-1;X = fftshift (fft (x.*hann (length (x),N);plot (k,20*log10 (abs(X), axis (0 fix (N/2) - inf inf )4.6 加載“timit1”語音數(shù)據(jù)。利用函數(shù)specgram，畫出該句語音的語譜圖。其Matlab代碼如下：load timit1;NFFT = 256 ; Fs = 16000; Win = 256; Noverlap = 128;specgram(timit1, NFFT, Fs, Win, Noverlap);4.7加載“gliss”語音數(shù)據(jù)。畫出該句語音的語譜圖。load gliss;NFFT = 256 ; Fs = 10000; Win = 256; Noverlap = 128;specgram(gliss.i_2, NFFT, Fs, Win, Noverlap);4.8 基于濁音/清音分類函數(shù)voiunvoi.m和諧波乘積譜函數(shù)hpspectrum.m，編寫一個用以計算濁音短時基音周期的函數(shù)“stpitch.m”。其Matlab代碼如下：function Fp = stpitch(x,N,Pth,Zth,NFFT,R,Fs)M = length(x);N = 2*fix(N/2);win = hann(N);Fp = zeros(M,1); voi = voiunvoi(x,N,Pth,Zth);invoi = 0; for (m=N:N/2:M) n1 = m-N+1:m-N/2; n2 = m-N/2+1:m; n = n1 n2; if (any(voi(n) if (invoi) Fp(n1) = Fmax; invoi = 0; end else HPSmax,Fmax = max(abs(hpspectrum(x(n),NFFT,R); if (invoi) Fp(n) = Fp(n) + win*Fmax; else Fp(n) = Fmax*ones(N/2,1) win(N/2+1:N); invoi = 1; end endend wsave = warning; warning(off);Fp = (Fp*Fs/NFFT).*(Fp./Fp);warning(wsave);1. 調用matlab中的wavread和wavplay讀入并播放語音文件。2. 調用fft命令計算語音的頻譜。3. 調用specgram命令得到語音的語譜圖。4. 調用plot命令繪制并觀察語音的時頻域波形。5. 根據(jù)短時功率、短時過零率、諧波乘積譜的計算公式，編寫分別用以計算短時功率、短時過零率、基音頻率的函數(shù)：function Px=stpower、function Zx = stzerocross(x, N) 和function HPSx = hpspectrum(x,N,R)。6. 利用短時功率和短時過零率這兩個參數(shù)可以對語音信號進行濁音/清音分類。編寫實現(xiàn)此分類功能的函數(shù)voiunvoi.m。4.9利用函數(shù)stpitch.m，分別求出兩句語音“timit1”和“timit2”的基音周期。實驗二 語音編碼一、實驗目的熟悉語音基本壓縮編碼的方法，觀察語音壓縮效果，加深對語音線性預測編碼（LPC）的理解。二、實驗內容1、編寫并調試語音LPC參數(shù)提取程序。2、編寫并調試語音基音周期提取程序。3、編寫并調試語音LPC合成程序。三、實驗原理語音信號中含有大量的冗余信息，采用各種信源編碼技術減除語音信號的冗余度，并充分利用人耳的聽覺掩蔽效應，就可以將其編碼速率壓縮很多倍，而仍能提供可懂語音。LPC聲碼器是一種比較簡單實用的語音壓縮方法，其基本原理是：根據(jù)語音生成模型，將語音看作激勵源通過一個線性時不變系統(tǒng)產(chǎn)生的輸出，利用線性預測分析對聲道參數(shù)進行估值，將求得的線性預測系數(shù)，結合基音周期等少量參數(shù)進行傳輸，就可以在接收端利用合成濾波器重構語音信號。線性預測系數(shù)的估計方法為：假設語音的當前樣值可以用過去的個語音樣值來進行預測 式中即為線性預測系數(shù)。實際值和預測值之間的均方誤差可表示為 要求均方誤差總和最小，將關于的偏導數(shù)設置為零，可以得到通過采用自相關法、協(xié)方差法或格形法求解該方程，即可得到最優(yōu)的。四、實驗方法及程序1. 調用xcorr命令計算一幀語音的自相關函數(shù)。2. 調用toeplitz命令形成該幀語音的自相關矩陣。3. 調用durbin命令，采用杜賓遞推算法計算該幀語音的線性預測系數(shù)。4. 編寫lpcauto.m函數(shù)，求取一句語音信號的線性預測系數(shù)及預測殘差。選擇設當?shù)拇昂瘮?shù)對語音信號進行分幀。5. 編寫lpcpitch函數(shù)，由殘差信號計算該句語音的基音周期。6. 編寫lpcgain函數(shù)，由預測殘差能量，求出該句語音的增益。7. 編寫lpcsyn函數(shù)，由該句語音的基音周期、預測殘差能量和增益進行LPC合成。五、實驗結果與分析1. 如何對全極點模型的線性預測參數(shù)進行提??？加窗等預處理對提取結果有何影響？自相關法、協(xié)方差法求解線性預測方程組、格型法。協(xié)方差法是先把計算均方誤差的間隔N固定下來，再研究它對計算n（i,1）的其中N是基音周期值。添加漢明窗和漢寧窗可以得到效果不錯的波形，很明顯兩個波形相似，矩形窗的波形相對不是很平穩(wěn)。通過改變幀長、幀移后波形更加的清晰，特征更加的明顯，更方便觀察。2. 如何實現(xiàn)線性預測參數(shù)的各種表現(xiàn)方式之間的轉換？對于語音信號，確定了各線性預測系數(shù)后，根據(jù)H（z）可得其頻率響應的估值即LPC譜。3. 使用預測殘差求基音周期有何優(yōu)點？ 該方法能進一步去除聲道共振峰的影響，簡化清濁音判決過程，提高判決精度，有效降低發(fā)生基音倍頻和半頻的錯誤。4. 如何根據(jù)線性預測系數(shù)求得LPC頻譜？對于語音信號，確定了各線性預測系數(shù)后，根據(jù)H（z）可得其頻率響應的估值即LPC譜。5. LPC頻譜與基于DFT求得的語音對數(shù)幅度譜相比有何不同？LPC譜優(yōu)點在于可以很好地表示共振峰結構而不出現(xiàn)額外的蜂起和起伏，基本不含有噪聲，其頻率靈敏度與人耳不相匹配；DET得到的對數(shù)幅度譜受基頻諧波的影響，最大值只出現(xiàn)在諧波頻率上，共振峰測定誤差較大。6. LPC合成語音的質量如何？有何改進措施？LPC合成技術的優(yōu)點是簡單直觀。其合成過程實質上只是一種簡單的解碼和拼接過程。另外，由于波形拼接技術的合成基元是語音的波形數(shù)據(jù)，保存了語音的全部信息，因而對于單個合成基元來說能夠獲得很高的自然度。但是，由于自然語流中的語音和孤立狀況下的語音有著極大的區(qū)別，如果只是簡單地把各個孤立的語音生硬地拼接在一起，其整個語流的質量勢必是不太理想的。而LPC技術從本質上來說只是一種錄音+重放，對于合成整個連續(xù)語流LPC合成技術的效果是不理想的。因此，LPC合成技術必須和其他技術相結合，才能明顯改善LPC合成的質量。實驗過程記錄與結果分析1、計算短時自相關函數(shù)1.1 分別計算正弦信號和白噪聲的短時自相關函數(shù)，估計正弦信號的基音周期。其Matlab代碼如下：x=sin(2*pi*0.01*(0:499);r,eta=xcorr(x,100,unbiased);stem(eta,r);w=randn(500,1);r,eta=xcorr(w,100,unbiased);stem(eta,r);1.2 計算正弦信號疊加白噪聲的短時自相關函數(shù)，試估計正弦信號的基音周期。其Matlab代碼如下：x=sin(2*pi*0.01*(0:499);w=randn(500,1);x1=x+w;r,eta=xcorr(x1,100,unbiased);stem(eta,r);1.3 分別畫出一幀濁音和一幀清音的語音自相關波形（采樣頻率為10kHz，幀長為25.6ms，每幀有256個樣點），估計濁音的基音周期。load digits; x=digits.three1;m=2756; N=256; n=m-N+1:m;r,eta=xcorr(x(n),250,unbiased);plot(eta,r);m=500; N=256; n=m-N+1:m;r,eta=xcorr(x(n),250,unbiased);plot(eta,r);2、求取LPC系數(shù)2.1加載“digits.three1”語音數(shù)據(jù)?；谧韵嚓P法，求出其中一幀數(shù)據(jù)（采樣頻率為10kHz，幀長為25.6ms，每幀有256個樣點）的14階LPC系數(shù)。畫出相應的LPC譜。load digits; x=digits.three1;m=2756; N=256; n=m-N+1:m; M=14;r,eta=xcorr(x(n),250,unbiased);Rx=toeplitz(r(M+1:2*M);rx=r(M+2:2*M+1);a=Rxrx;NFFT=1024; k=1:NFFT/2;X=fft(x(n).*hann(N),NFFT);Theta=1./fft(1;-a,NFFT);plot(k,20*log10(abs(353*Theta(k) X(k)axis(0 NFFT/2 inf inf)2.2編寫一個用以實現(xiàn)杜賓遞推算法的函數(shù)“durbin.m”，其Matlab代碼如下：function a,xi,kappa=durbin(r,M)kappa = zeros (M,1);a = zeros (M, 1);xi=r(1); zeros(M,1);for (j=1:M)kappa(j)=r(j+1)-a(1:j-1)*r(j:-1:2)/xi(j);a(j)=kappa(j);a(1:j-1)=a(1:j-1)-kappa(j)*a(j-1:-1:1);xi(j+1)=a(1:j-1)-kappa(j)*a(j-1:-1:1);xi(j+1)=xi(j)*(1-kappa(j)2);end2.3利用函數(shù)“durbin.m”，計算2.1中語音數(shù)據(jù)幀的14階LPC系數(shù)，并與2.1中的結果進行比較：x=digits.three1;m=2756; N=256; n=m-N+1:m; M=14;r,eta=xcorr(x(n),250,unbiased);aLD,xi,kappa=durbin(r(M+1:2*M+1),M);a, aLD, norm(a-aLD)a =Johannesburg n. 約翰內斯堡 0.3512 0.2072 0.3808 0.3163 0.0522 -0.1639 -0.2903 -0.2898 -0.1222 0.1960 0.3492 0.3024 0.1713 -0.8258aLD = 0.3512 0.2072 0.3808 0.3163 0.0522 -0.1639 -0.2903 -0.2898 -0.1222 0.1960 0.3492 0.3024 0.1713 -0.8258ans = 1.0047e-142.4編寫一個用以實現(xiàn)反射系數(shù)轉換為LPC系數(shù)的函數(shù)“rf2lpc.m”，其Matlab代碼如下：function a = rf2lpc(kappa)M = length(kappa);a = zeros(M,1);for (j=1:M) a(j) = kappa(j); a(1:j-1) = a(1:j-1) - kappa(j)*a(j-1:-1:1);end2.5 編寫一個用以實現(xiàn)LPC系數(shù)轉換為反射系數(shù)的函數(shù)“l(fā)pc2rf.m”，其Matlab代碼如下：function kappa = lpc2rf(a)M = length(a);kappa = zeros(M,1);for (j=M:-1:1) kappa(j) = a(j); a(1:j-1) = (a(1:j-1) + a(j)*a(j-1:-1:1)/(1 - kappa(j)2);end2.6使用函數(shù)“rf2lpc.m”和“l(fā)pc2rf.m”，檢驗反射系數(shù)和LPC系數(shù)相互轉換結果。norm(kappa-lpc2rf(aLD)norm(aLD-rf2lpc(kappa)ans = 1.4760e-14ans = 03、語音信號的逆濾波3.1 利用2.3中求出的一幀語音數(shù)據(jù)的LPC系數(shù)，構造逆濾波器，并畫出該幀語音信號的殘差波形。ehat=filter(1;-1,1,x(n);plot(x(n) ehat)3.2編寫一個用以實現(xiàn)語音信號LPC分析的函數(shù)“l(fā)pcauto”，其Matlab代碼如下：function ar,xi,e,m = lpcauto(x,M,win,Olap)Nx = length(x); N = length(win); if (N = 1) N = win; win = ones(N,1); endF = fix(Nx-Olap)/(N-Olap);ar = zeros(M+1,F);xi = zeros(M+1,F);e = zeros(Nx,1);m = zeros(F,1); n = 1:N; n1 = 1:Olap; n2 = N-Olap+1:N; n3 = Olap+1:N; win1 = win(n1)./(win(n1)+win(n2)+eps);win2 = win(n2)./(win(n1)+win(n2)+eps); for (f=1:F) r,eta = xcorr(x(n).*w