電子科技大學數(shù)字信號處理DSP課程設計-鋼琴音符識別

1、2014級數(shù)字信號處理課程設計報告題目：鋼琴音符識別姓名：邱晨曦學號：2014010909008答辯時間：2016/12/9推薦精選一. 題目要求：（1） 播放和記錄一段鋼琴音樂中的音符；（2） 記錄到音符以后，找到音符所對應的現(xiàn)代標準鋼琴的鋼琴鍵，并分析結(jié)果。二. 課程設計思路：（1） 涉及到的知識點：快速傅里葉變換、鋼琴音頻信號的時域和頻域的特性、能熵比的概念、頻率校正、頻率與音符的轉(zhuǎn)換關系。（2） 方案分析：A. 預處理部分：1. 直接用audioread函數(shù)讀出來的原始數(shù)據(jù)。優(yōu)點：準確率較高；缺點：數(shù)據(jù)量較大，采樣頻率為44kHz，遠大于奈奎斯特采樣率。2. 以11kHz的采樣率重新采

2、樣，并轉(zhuǎn)換為單聲道。優(yōu)點：數(shù)據(jù)量小了很多，易于處理；缺點：犧牲了部分的準確率，但對于音符的判斷影響可以忽略。 B. 端點檢測算法： . 雙門限法：1. 計算短時能量（高門限）和過零率（低門限）；2. 選取一個較高的門限，語音信號的能量包絡大部分都在此門限之上，進行一次初判，語音起止點位于該門限與短時能量包絡交點所對應的時間間隔之外；3. 根據(jù)噪聲能量，確定一個較低的門限，并從初判起點往左，從初判終點往右搜索，分別找到能零比曲線第一次與門限相交的兩個點，兩點之間段就是用雙門限方法所判定的語音段；4. 以短時平均過零率為準，從低門限點往左右搜索，找到短時平均過零率低于某閾值的兩點，為語音的起止點；

3、 圖1：雙門限法示意圖推薦精選說明：算法中的閥值是根據(jù)實驗過程調(diào)節(jié)的。該算法在實際應用的過程中發(fā)現(xiàn)：在語音信號頻率分布較為集中的時候，端點檢測出來的結(jié)果比較準確，但當語音信號頻率分布比較分散的時候，很難通過控制固定的閥值來檢測到每個音符；. 自相關法：由于兩種信號的自相關函數(shù)存在極大的差異，可以利用這種差別來提取語音端點。根據(jù)噪聲的情況，設置兩個閾值和，當相關函數(shù)最大值大于時，便判定是語音；當相關函數(shù)最大值大于或小于時，則判定為語音信號的端點。 該算法同樣存在當語音信號頻率分布較廣的時候，閥值比較難控制的問題。. 基于譜熵的端點檢測：基于譜熵語音端點檢測方法是通過檢測譜的平坦程度，來進行語音端

4、點檢測的，為了更好進行語音端點檢測，采用語音信號的短時功率譜構(gòu)造語音信息譜熵，從而對語音段和噪聲段進行區(qū)分。檢測思路：1. 對語音信號進行分幀加窗；2. 計算每一幀的譜能量；3. 計算出每一幀中每個樣本點的概率密度函數(shù)；4. 計算出每一幀的譜熵值（由信息論知識知道，熵值在自變量服從均勻分布的時候，熵值達到最大值，所以噪聲的熵值是比較大的，而鋼琴音符的熵值是比較小的，由此區(qū)別了噪聲和音符）；5. 設置判決門限；6. 根據(jù)各幀的譜熵值進行端點檢測。 在實驗過程中發(fā)現(xiàn)：依然存在當語音信號頻率分布較廣時，閥值不太好控制的問題。因此對該方法進行改進，引入，能熵比的概念：譜熵值類似于過零率，能熵比的表示為

5、。由于噪聲和信號的能熵比差別很大。因此在能熵比的圖像中，每一個“尖刺”就代表了一個特定頻率的語言信號。推薦精選圖2：能熵比圖中的“尖刺”在檢測過程中，依然不能通過簡單的設置閥值的辦法來進行端點檢測，原因是語音頻率分布較廣時，每個音符的能熵比變化范圍差別較大，如下圖所示，有的“尖刺”完全在門限之上，而有的則完全在門限之下。圖3：88階全音的能熵比圖因此，采用檢測能熵比中的“低谷點”（該點比左右兩邊的一定數(shù)目的點的能熵比都?。┑姆椒?。語音信號一定位于兩個低谷點之間的部分，再對低谷點進行適當?shù)淖笥乙苿幼鳛檎Z音信號的起止點。如下圖所示：圖4：標記起止點的能熵比圖（綠色為起始點，紅色為截止點）（3） 設

6、計框架和流程：1. 用audioread函數(shù)讀入鋼琴音樂，并用sound函數(shù)播放；2. 為了方便處理，對信號以11.025kHz的頻率進行重新采樣，并統(tǒng)一轉(zhuǎn)換成單聲道的信號；3. 因為語言信號可以在短時間內(nèi)認為是平穩(wěn)的，因此對語音信號進行分幀的處理，設置幀長320，為了減小誤差，兩幀之間設置重疊部分，因此幀移取80；4. 計算每一幀的能熵比；5. 找到能熵比中的“低谷點”（該點比左右兩邊的一定數(shù)目的點的能熵比都?。?；6. 如果兩個低谷點之間的距離大于miniL(認為持續(xù)長度超過一定長度的為音符，最小長度推薦精選miniL可自行設置)。則低谷點右移sr (即shift right,數(shù)值可自行調(diào)節(jié)

7、)幀作為一段信號的起始點，將低谷點左移sl（即shift left,數(shù)值可自行調(diào)節(jié)）幀作為截止點注：采用該方法的優(yōu)點是通過調(diào)節(jié)相關參數(shù)能適應多種情況，缺點是檢測環(huán)境發(fā)生較大變化時，需要重新設置參數(shù)；7. 將找到的語音段轉(zhuǎn)換成未分幀時對應坐標的語音段，并對每段做快速傅里葉變換；8. 找到每段快速傅里葉變換中的最大值以及最大值所對應的橫坐標（fft點），將橫坐標轉(zhuǎn)換成相應的頻率，得到的頻率即為該段音符的頻率；9. 利用比值法進行頻率的校正，窗函數(shù)選擇矩形窗；10. 根據(jù)檢測到的頻率確定音符，計算公式為：，為第幾個按鍵，再通過查表得到對應音符；11. 分析結(jié)果。三. 具體設計過程： （1） 部分代碼

8、（測試部分缺?。褐骱瘮?shù)部分：x,fs=audioread(鋼琴音頻.WAV); format short;wlen=320; inc=80; % 分幀的幀長和幀移overlap=wlen-inc; % 幀之間的重疊部分sound(x,fs); % 播放音樂x=calsample(x,fs); % 為了方便處理，重新以11025Hz的頻率采樣，并轉(zhuǎn)換成單聲道x=x-mean(x); % 消去直流分量x=x/max(abs(x); % 幅值歸一化y = Enframe(wlen,inc,x); % 分幀fn = size(y,2); % 取得幀數(shù)time = (0 : length(x)-1)/

9、11025; % 計算時間坐標frameTime = frame2time(fn, wlen, inc, 11025); % 計算各幀對應的時間坐標sr=2;sl=13;miniL=33; % 配置左右移動的幀數(shù)和要求的最短幀數(shù)voicesegment,vos,Ef=get_segment(y,fn,sr,sl,miniL); % 獲得語音段real_f,ft,ax=get_f(x,voicesegment,vos,wlen,inc); % 檢測頻率的結(jié)果for i=1:length(real_f) real_f(i)=roundn(real_f(i),-4);end for i=1:leng

10、th(real_f) real_node(i)=get_node(real_f(i)end%*繪圖部分*subplot 211;stem(real_f);title(頻率檢測結(jié)果);xlabel(音符/個);ylabel(頻率/Hz);推薦精選subplot 212;stem(real_node,r);title(音符檢測結(jié)果);xlabel(音符/個);ylabel(對應按鍵);figure(2);subplot 211;plot(Ef);title(能熵比圖及語音起止點);xlabel(幀數(shù)/個);ylabel(能熵比);for i=1:length(voicesegment) text(

11、voicesegment(i).begin,Ef(voicesegment(i).begin),o,color,g) text(voicesegment(i).end,Ef(voicesegment(i).end),o,color,r)end subplot 212, plot(time,x,k); title(語音信號端點檢測結(jié)果)axis(0 max(time) -1 1); ylabel(幅值);for k=1 : vos % 標出有話段 nx1=voicesegment(k).begin; nx2=voicesegment(k).end; nxl=voicesegment(k).dur

12、ation; fprintf(%4d %4d %4d %4dn,k,nx1,nx2,nxl); subplot 212 line(frameTime(nx1) frameTime(nx1),-1 1,color,r,linestyle,-); line(frameTime(nx2) frameTime(nx2),-1 1,color,r,linestyle,-);end 其中的用到的子函數(shù)：1.calsample.m （調(diào)整采樣率和聲道）function sample = calsample(sampledata,FS)temp_sample = resample(sampledata,1,F

13、S/11025); %調(diào)整采樣頻率,n = size(temp_sample);if (n = 2) %轉(zhuǎn)換成單聲道 sample = temp_sample(:,1);else sample = temp_sample;endend2. Enframe.m （分幀函數(shù)）function f=Enframe(len,inc,x) %對讀入的語音進行分幀，len為幀長， %inc為幀重疊樣點數(shù)，x為輸入語音數(shù)據(jù)fh=fix(size(x,1)-len)/inc)+1); %計算幀數(shù) f=zeros(fh,len); %設置一個零矩陣，行為幀數(shù)，列為幀長i=1; 推薦精選n=1;while i=f

14、h %幀間循環(huán) j=1; while j=len %幀內(nèi)循環(huán) f(i,j)=x(n); j=j+1; n=n+1; end n=n-len+inc; %下一幀開始位置 i=i+1;end3. frame2time.m（坐標刻度轉(zhuǎn)換）function frameTime=frame2time(frameNum,framelen,inc,fs)frameTime=(1:frameNum)-1)*inc+framelen/2)/fs; % 求對應的時間坐標4.get_segment.m（端點檢測，確定音符段）function voicesegment,vos,Ef=get_segment(y,fn,

15、sr,sl,miniL)if size(y,2)=fn, y=y; end % 把y轉(zhuǎn)換為每列數(shù)據(jù)表示一幀語音信號wlen=size(y,1); % 取得幀長for i=1:fn Sp = abs(fft(y(:,i); % FFT取幅值 Sp = Sp(1:wlen/2+1); % 只取正頻率部分 Esum(i) = sum(Sp.*Sp); % 計算能量值 prob = Sp/(sum(Sp); % 計算概率 H(i) = -sum(prob.*log(prob+eps); % 求譜熵值endhindex=find(H=miniL % 剔除持續(xù)幀長度小于miniL的音符段推薦精選 j=j+

16、1; temp1=x2(k-1)+sr; % 將低谷點右移sr個幀，作為一個音符的起始點 voicesegment(j).begin=temp1; temp2=x2(k)-sl; % 將低谷點左移sl個幀，作為一個音符的截止點 voicesegment(j).end=temp2; voicesegment(j).duration=voicesegment(j).end-voicesegment(j).begin+1; %音符持續(xù)幀數(shù) endendvos=length(voicesegment); %返回音符個數(shù)end5.get_max.m（找到高峰點）function max_x,max=ge

17、t_max(x)l=length(x); %獲得數(shù)組的長度max=;max_x=;j=1;for i=20:l-15 %找到“峰值點” if (x(i)x(i-1)&(x(i)x(i-2)&(x(i)x(i-3)&(x(i)x(i-4)&(x(i)x(i-5). &(x(i)x(i-6)&(x(i)x(i-7)&(x(i)x(i-8)&(x(i)x(i-9)&(x(i)x(i-10). &(x(i)x(i-11)&(x(i)x(i-12)&(x(i)x(i-13)&(x(i)x(i-14). &(x(i)x(i-15)&(x(i)x(i-16)&(x(i)x(i-17). &(x(i)x(i+

18、1)&(x(i)x(i+2)&(x(i)x(i+3)&(x(i)x(i+4)&(x(i)x(i+5). &(x(i)x(i+6)&(x(i)x(i+7)&(x(i)x(i+8)&(x(i)x(i+9)&(x(i)x(i+10). &(x(i)x(i+11)&(x(i)x(i+11)&(x(i)x(i+12)&(x(i)x(i+13)&(x(i)x(i+14). &(x(i)x(i+15)&(x(i)x(i+16). &(x(i)0.1) max_x(j)=i; %找到后賦值給返回參數(shù) max(j)=x(i); i=i+1; j=j+1; else推薦精選 i=i+1; end endend6.

19、get_min.m（找到低谷點）function min_x,min=get_min(x)l=length(x); %獲得數(shù)組的長度min=;min_x=;j=1;for i=100:l-10 %尋找低谷點 if (x(i)x(i-1)&(x(i)x(i-2)&(x(i)x(i-3)&(x(i)x(i-4)&(x(i)x(i-5). &(x(i)x(i-6)&(x(i)x(i-7)&(x(i)x(i-8)&(x(i)x(i-9)&(x(i)x(i-10). &(x(i)x(i-11)&(x(i)x(i-12)&(x(i)x(i-13)&(x(i)x(i-14). &(x(i)x(i-15)&(

20、x(i)x(i-16)&(x(i)x(i-17)&(x(i)x(i-18)&(x(i)x(i-19). &(x(i)x(i-20)&(x(i)x(i-21)&(x(i)x(i-22)&(x(i)x(i-23)&(x(i)x(i-24). &(x(i)x(i-25)&(x(i)x(i-26)&(x(i)x(i-27)&(x(i)x(i-28)&(x(i)x(i-29). &(x(i)x(i-30)&(x(i)x(i-31)&(x(i)x(i-32)&(x(i)x(i-33)&(x(i)x(i-34). &(x(i)x(i+1)&(x(i)x(i+2)&(x(i)x(i+3)&(x(i)x(i+4

21、)&(x(i)x(i+5). &(x(i)x(i+6)&(x(i)x(i+7)&(x(i)x(i+8)&(x(i)x(i+9)&(x(i)nx1nx,mx=size(x);if mx=1, x=x;end %轉(zhuǎn)換成行矩陣nx,mx=size(x);M=fix(N/2)+mod(N,2); xf=fft(x);%xf=xf(1:M)*2/N;ddf=fs/N; % 頻率分辨率n1=fix(nx1/ddf); % 將頻率轉(zhuǎn)換成fft對應的點n2=round(nx2/ddf); A=abs(xf); % 取fft的模值Amax,index=max(A(n1:n2); %找到fft點n1到n2之間幅值

22、的最大值index=index+n1-1; % 移動到n1和n2中間%比值法%加矩形窗indsecL=A(index-1)A(index+1); %滿足條件則為1，不滿足則為0df=indsecL.*A(index-1)./(Amax+A(index-1)-(1-indsecL).*A(index+1)./(Amax+A(index+1);Z=(index-1-df)*ddf; %修正后的頻率end9.get_node.m（確定音符）function node=get_node(f)推薦精選node=12*log2(f/440)+49;node=round(node); %四舍五入確定對應按鍵

23、end （2） 生成的圖與數(shù)據(jù)等端點檢測結(jié)果序號起始幀截止幀持續(xù)幀數(shù)14097582112166553181235554250303545318372556387441557456510558525579559594648551066371452117297865812801854541386992355149389935615100810514416106611225717113711764018119112677719128213375620135213914021140614757022149015132423152815886124160316666425168117305026176

24、0180849271823187351281903194442291959201759302032209059312105213632322177222549332240227839342317234024352380242950推薦精選推薦精選對數(shù)據(jù)的分析：真實按鍵真實音符真實頻率/Hz檢測頻率/Hz對應按鍵檢測音符絕對誤差/Hz相對誤差52G7/A7523.251524.897352G7/A71.64630.31%53G7554.365555.043753G70.67870.12%54F7/G7587.33589.876854F7/G72.54680.43%55F7622.254623.1

25、58555F70.90450.15%56E7659.255659.707956E70.45290.07%57D7/E7698.456699.296857D7/E70.84080.12%58D7739.989741.850358D71.86130.25%59C7/D7783.991785.948559C7/D71.95750.25%60C7Double high C830.609832.41460C7Double high C1.8050.22%61B6880883.230361B63.23030.37%62A6/B6932.328933.856562A6/B61.52850.16%63A6987.767989.996463A62.22940.23%64G6/A61046.51047.0164G6/A60.510.05%65G61108.731112.724465G63.99440.36%66F6/G61174.661177.113766F6/G62.45370.21%67F61244.511248.123867F63.61380.29%68E61318.511322.60968E64.0990.31%69D6/E61396.911403.008469D6/E66.09840.44%70D61479.981486.316270D66.33620.43%71C6/D