基于matlab的音樂信號(hào)處理和分析

1、基于matlab的音樂信號(hào)處理和分析1音樂信號(hào)的音譜和頻譜的觀察使用windows下的錄音機(jī)錄制一段音樂信號(hào)或者采用其他軟件截取一段音樂信號(hào)（要求：時(shí)間不超過5s，文件格式為wav） 使用wavread語句讀取音樂信號(hào)獲取抽樣率；（注意：讀取的信號(hào)是雙聲道信號(hào)，即為雙列向量，需要分列處理）； 輸出音樂信號(hào)的波形和頻譜，觀察現(xiàn)象； 使用sound語句播放音樂信號(hào)，注意不同抽樣率下的音調(diào)變化，解釋現(xiàn)象程序：clear all;close all;clcy,fs,bit=wavread('e:music11'); %讀取音樂信號(hào)size(y) %看音樂信號(hào)是雙列還是單列y1=y(:,

2、1); %取單列fs=fs %獲取音樂信號(hào)的抽樣率n=length(y1); %音樂信號(hào)的長度fy1=fft(y1,n); %對(duì)信號(hào)做快速傅里葉變換w=2/n*0:n-1;figure %畫音樂信號(hào)的波形和頻譜subplot(2,1,1);plot(y1);grid ontitle('音樂信號(hào)的波形')xlabel('time/s')ylabel('magnitude')subplot(2,1,2);plot(w,abs(fy1);grid ontitle('音樂信號(hào)的頻譜')xlabel('frequency/pi

3、9;)ylabel('magnitude')sound(y1,fs) %以抽樣率fs播放音樂信號(hào)sound(y1,2*fs) %以抽樣率2fs播放音樂信號(hào)sound(y1,1/2*fs) %以抽樣率1/2fs播放音樂信號(hào)程序運(yùn)行結(jié)果如下：音樂信號(hào)的抽樣率fs=44100音樂信號(hào)的波形和頻譜如下圖：抽樣率較高時(shí)，聽到的音樂音調(diào)較高，速度較快；抽樣率較低時(shí)，聽到的音樂音調(diào)低，速度慢。由此可知，音樂信號(hào)的音調(diào)和抽樣率有關(guān)。2音樂信號(hào)的抽?。p抽樣） 觀察音樂信號(hào)頻率的上限，選擇適當(dāng)?shù)某槿￠g隔對(duì)信號(hào)進(jìn)行減抽樣（給出兩種抽取間隔，代表混疊與非混疊）； 輸出減抽樣音樂信號(hào)的波形和頻譜，觀察

4、現(xiàn)象，給出理論解釋； 播放減抽樣音樂信號(hào)，注意抽樣率的改變，比較不同抽取間隔下的聲音，解釋現(xiàn)象。音樂信號(hào)的頻率上限為0.2pi，這里分別采用3和60作為非混疊與混疊的抽樣間隔進(jìn)行抽樣，程序如下：clear all;close all;clcy,fs,bit=wavread('e:music22');%讀取音樂信號(hào)y1=y(:,1); %取單列x1=y1(1:180000); %截取音樂信號(hào)d=3;j=0;for i=1:d:length(x1) j=j+1; x2(j)=x1(i);end %對(duì)音樂信號(hào)以3為間隔做減抽樣n1=length(x1); %截?cái)嗟囊魳沸盘?hào)的長度n2=

5、length(x2); %以3為間隔減抽樣后音樂信號(hào)的長度fx1=fft(x1,n1); %對(duì)截?cái)嘈盘?hào)做快速傅里葉變換fx2=fft(x2,n2); %對(duì)減抽樣后的信號(hào)做快速傅里葉變換w1=2/n1*0:n1-1;w2=2/n2*0:n2-1;figure %畫出截?cái)嗪蟮男盘?hào)的波形和頻譜subplot(2,1,1);plot(x1);grid ontitle('截?cái)嗪蟮牟ㄐ?#39;)xlabel('time/s')ylabel('magnitude')subplot(2,1,2);plot(w1,abs(fx1);grid ontitle('截

6、斷信號(hào)的頻譜')xlabel('frequency/pi')ylabel('magnitude')figure %畫出以3為間隔減抽樣后的信號(hào)的波形和頻譜subplot(2,1,1);plot(x2);grid ontitle('對(duì)截?cái)嘈盘?hào)以3為間隔減抽樣后的波形')xlabel('time/s')ylabel('magnitude')subplot(2,1,2);plot(w2,abs(fx2);grid ontitle('對(duì)截?cái)嘈盘?hào)以3為間隔減抽樣后的頻譜')xlabel('fr

7、equency/pi')ylabel('magnitude')d=60;j=0;for i=1:d:length(x1) j=j+1; x3(j)=x1(i);end %以60為間隔對(duì)截?cái)嘈盘?hào)進(jìn)行減抽樣n3=length(x3); fx3=fft(x3,n3);w3=2/n3*0:n3-1;figure %畫以60為間隔減抽樣后的音樂信號(hào)的波形和頻譜subplot(2,1,1);plot(x3);grid ontitle('對(duì)截?cái)嘈盘?hào)以60為間隔減抽樣后的波形')xlabel('time/s')ylabel('magnitude&

8、#39;)subplot(2,1,2);plot(w3,abs(fx3);grid ontitle('對(duì)截?cái)嘈盘?hào)以60為間隔減抽樣后的頻譜')xlabel('frequency/pi')ylabel('magnitude')sound(x1,fs) %播放截?cái)嘈盘?hào)sound(x2,1/3*fs) %播放以3為間隔做減抽樣后的音樂信號(hào)sound(x3,1/60*fs) %播放以60為間隔做減抽樣后的音樂信號(hào)實(shí)驗(yàn)結(jié)果如圖：使用兩種抽樣間隔進(jìn)行實(shí)驗(yàn)： 以3為間隔進(jìn)行減抽樣，得到的音樂信號(hào)時(shí)域長度變?yōu)樵瓉淼?/3，頻譜寬度變?yōu)樵瓉淼?倍，這是時(shí)域截?cái)嘣斐?/p>

9、的頻譜泄漏。由于原音樂信號(hào)的頻率上限為0.2，故變?yōu)樵瓉淼?倍后仍不混疊，依然能夠聽到音樂； 以60為間隔進(jìn)行減抽樣，得到的音樂信號(hào)時(shí)域點(diǎn)數(shù)變?yōu)樵瓉淼?/60，相應(yīng)的，頻譜寬度變?yōu)樵瓉淼?0倍，此時(shí)頻譜已經(jīng)混疊混疊，播放音樂信號(hào)時(shí)什么都聽不到。3音樂信號(hào)的am調(diào)制 觀察音樂信號(hào)的頻率上限，選擇適當(dāng)調(diào)制頻率對(duì)信號(hào)進(jìn)行調(diào)制（給出高、低兩種調(diào)制頻率）； 輸出調(diào)制信號(hào)的波形和頻譜，觀察現(xiàn)象，給出理論解釋； 播放調(diào)制音樂信號(hào)，注意不同調(diào)制頻率下的聲音，解釋現(xiàn)象。音樂信號(hào)的頻率上限為0.2pi，這里采用0.5pi與0.95pi作為調(diào)制頻率對(duì)音樂信號(hào)進(jìn)行調(diào)制，程序如下：close all;clear all

10、;clcy,fs,bit=wavread('e:music22');y1=y(:,1);x1=y1(1:180000);n1=1:180000;c1=cos(n1*pi/2);%余弦序列 x4=x1.*c1' %對(duì)信號(hào)進(jìn)行調(diào)制n4=length(x4); %調(diào)制后信號(hào)的長度fx4=fft(x4,n4); %對(duì)調(diào)制信號(hào)做快速傅里葉變換w4=2/n4*0:n4-1;figure %畫出調(diào)制后的信號(hào)的波形和頻譜subplot(2,1,1);plot(x4);grid ontitle('對(duì)信號(hào)以頻率pi/2調(diào)制后的波形')xlabel('time/s&#

11、39;)ylabel('magnitude')subplot(2,1,2);plot(w4,abs(fx4);grid ontitle('對(duì)信號(hào)以頻率pi/2調(diào)制后的頻譜')xlabel('frequency/pi')ylabel('magnitude')n1=1:180000;c1=cos(n1*pi*0.95);x5=x1.*c1' n5=length(x5);fx5=fft(x5,n5);w5=2/n5*0:n5-1;figuresubplot(2,1,1);plot(x5);grid ontitle('對(duì)截

12、斷信號(hào)以頻率pi*0.95調(diào)制后的波形')xlabel('time/s')ylabel('magnitude')subplot(2,1,2);plot(w5,abs(fx5);grid ontitle('對(duì)截?cái)嘈盘?hào)以頻率pi*0.95調(diào)制后的頻譜')xlabel('frequency/pi')ylabel('magnitude')sound(x1,fs)sound(x4,fs)sound(x5,fs)對(duì)音樂信號(hào)用兩種頻率進(jìn)行調(diào)制：調(diào)制頻率為pi/2時(shí)，調(diào)制后頻譜搬移到了0.5pi處，不混疊，可以聽到聲音，但

13、由于在高頻處，音調(diào)很高，是刺耳的一段聲音；調(diào)制頻率為0.95pi時(shí)，調(diào)制后頻譜搬移到0.95pi處，頻譜混疊，不能聽到聲音。4 am調(diào)制音樂信號(hào)的同步解調(diào) 計(jì)巴特沃斯濾波器完成同步解調(diào)；觀察濾波器頻率響應(yīng)曲線； 窗函數(shù)法設(shè)計(jì)fir濾波器完成同步解調(diào)，觀察濾波器頻率響應(yīng)曲線（要求：分別使用矩形窗和布萊克曼窗，進(jìn)行比較）； 輸出解調(diào)信號(hào)的波形和頻譜，觀察現(xiàn)象，給出理論解釋； 播放解調(diào)音樂信號(hào)，比較不同濾波器下的聲音，解釋現(xiàn)象。設(shè)計(jì)巴特沃斯濾波器對(duì)解調(diào)信號(hào)進(jìn)行濾波：clear all;close all;clcy,fs,bit=wavread('e:music22');y1=y(:,

14、1);x1=y1(1:180000);n1=1:180000;c1=cos(n1*pi/2);%生成余弦序列x3=x1.*c1' %調(diào)制x4=x3.*c1' %解調(diào)n4=length(x4);fx4=fft(x4,n4);w4=2/n4*0:n4-1;figure %畫解調(diào)后信號(hào)的波形和頻譜subplot(2,1,1);plot(x4);grid ontitle('同步解調(diào)后音樂信號(hào)的波形')xlabel('time/s')ylabel('magnitude')subplot(2,1,2);plot(w4,abs(fx4);gri

15、d ontitle('同步解調(diào)后音樂信號(hào)的頻譜')xlabel('frequency/pi')ylabel('magnitude')wp=0.2;ws=0.3;rp=1;rs=15; n,wc=buttord(wp,ws,rp,rs);b,a=butter(n,wc);hd,w=freqz(b,a); %設(shè)計(jì)數(shù)字巴特沃斯濾波器figureplot(w/pi,abs(hd)grid ontitle('數(shù)字巴特沃斯濾波器')xlabel('digital frequency/pi')ylabel('magnit

16、ude')y2=filter(b,a,x4); %用巴特沃斯濾波器進(jìn)行濾波n=length(y2);fy2=fft(y1,n);w=2/n*0:n-1;figuresubplot(2,1,1);plot(y2);grid ontitle('巴特沃斯濾波后音樂信號(hào)的波形')xlabel('time/s')ylabel('magnitude')subplot(2,1,2);plot(w,abs(fy2);grid ontitle('巴特沃斯濾波后音樂信號(hào)的頻譜')xlabel('frequency/pi')yl

17、abel('magnitude')sound(x4,fs)sound(y2,fs)用窗函數(shù)法設(shè)計(jì)矩形窗和布萊克曼窗對(duì)解調(diào)信號(hào)進(jìn)行濾波：clear all;close all;clcy,fs,bit=wavread('e:music22');size(y)y1=y(:,1);fs=fsx1=y1(1:90000);n1=1:90000;c1=cos(n1*pi/2);x3=x1.*c1'x4=x3.*c1'n4=length(x4);fx4=fft(x4,n4);w4=2/n4*0:n4-1;figuresubplot(2,1,1);plot(x4

18、);grid ontitle('同步解調(diào)后音樂信號(hào)的波形')xlabel('time/s')ylabel('magnitude')subplot(2,1,2);plot(w4,abs(fx4);grid ontitle('同步解調(diào)后音樂信號(hào)的頻譜')xlabel('frequency/pi')ylabel('magnitude')n=33;wc=pi/4;hd=ideal(n,wc);w1=boxcar(n);w2=blackman(n);h1=hd.*w1'h2=hd.*w2'm=

19、512;fh1=fft(h1,m);fh2=fft(h2,m);w=2/m*0:m-1;figuresubplot(2,1,1);plot(w,abs(fh1);grid ontitle('矩形窗的頻譜')xlabel('頻率')ylabel('幅度')subplot(2,1,2);plot(w,abs(fh2);grid ontitle('布萊克曼窗的頻譜')xlabel('頻率')ylabel('幅度')y3=conv(x4,h1);fy3=fft(y3,n4);figuresubplot(2,

20、1,1);plot(y3);grid ontitle('矩形窗濾波后的波形')xlabel('time/s')ylabel('magnitude')subplot(2,1,2);plot(w4,abs(fy3);grid ontitle('矩形窗濾波后的頻譜')xlabel('frequency/pi')ylabel('magnitude')y4=conv(x4,h2);fy4=fft(y4,n4);figuresubplot(2,1,1);plot(y4);grid ontitle('布萊

21、克曼窗濾波后的波形')xlabel('time/s')ylabel('magnitude')subplot(2,1,2);plot(w4,abs(fy4);grid ontitle('布萊克曼窗濾波后的頻譜')xlabel('frequency/pi')ylabel('magnitude') 布萊克曼窗濾波效果更好，更平滑。但是聽上去不同濾波器的效果差別不大。5 音樂信號(hào)的濾波去噪 原始信號(hào)疊加幅度為0.05，頻率為3khz,5khz,8khz的三余弦混合噪聲，觀察噪聲頻譜以及加噪后音樂信號(hào)的音頻和頻譜，并

22、播放音樂，感受噪聲對(duì)音樂信號(hào)的影響； 給原始音樂信號(hào)疊加幅度為0.5的隨機(jī)白噪聲（可用rand語句產(chǎn)生），觀察噪聲頻譜以及加噪后音樂信號(hào)的音譜和頻譜，并播放音樂，感受噪聲對(duì)音樂信號(hào)的影響； 根據(jù)步驟、觀察到的頻譜，選擇合適指標(biāo)設(shè)計(jì)濾波器進(jìn)行濾波去噪，關(guān)注去噪后信號(hào)音譜和頻譜，并播放音樂，解釋現(xiàn)象。clear all;close all;clcy,fs,bit=wavread('e:music22');y0=y(:,1);y1=y0(1:180000) %截取音樂信號(hào)n=length(y1);%取截?cái)嘈盘?hào)的長度n=0:n-1;x1=0.05*cos(n*2*pi*3000/fs)

23、+0.05*cos(n*2*pi*5000/fs)+0.05*cos(n*2*pi*8000/fs);%三余弦噪聲y2=y1+x1'%音樂信號(hào)和三余弦噪聲疊加fx1=fft(x1,n);fy2=fft(y2,n);w1=2/n*0:n-1;figuresubplot(3,1,1);plot(w1,abs(fx1);grid ontitle('三余弦混合噪聲的頻譜')xlabel('頻率')ylabel('幅度')subplot(3,1,2);plot(w1,abs(fy2);grid ontitle('加三余弦混合噪聲后音樂信號(hào)的

24、頻譜')xlabel('頻率')ylabel('幅度')subplot(3,1,3);plot(y2);grid ontitle('加三余弦混合噪聲后音樂信號(hào)的音譜')xlabel('時(shí)間')ylabel('幅度')x2=0.5*(rand(n,1)-0.5);%產(chǎn)生白噪聲y3=y1+x2;%音樂信號(hào)和白噪聲疊加fx2=fft(x2,n);fy3=fft(y3,n);figuresubplot(3,1,1);plot(w1,abs(fx2);grid ontitle('白噪聲的頻譜')xla

25、bel('頻率')ylabel('幅度')subplot(3,1,2);plot(w1,abs(fy3);grid ontitle('加白噪聲后音樂信號(hào)的頻譜')xlabel('頻率')ylabel('幅度')subplot(3,1,3);plot(y3);grid ontitle('加白噪聲后音樂信號(hào)的音譜')xlabel('時(shí)間')ylabel('幅度')wp=0.09;ws=0.11;rp=1;rs=15;n,wc=buttord(wp,ws,rp,rs);b,

26、a=butter(n,wc);hd,w=freqz(b,a);%設(shè)計(jì)巴特沃斯濾波器y4=filter(b,a,y2);%對(duì)加三余弦噪聲的音樂信號(hào)進(jìn)行濾波n2=length(y4);fy4=fft(y4,n2);w2=2/n2*0:n2-1;figuresubplot(2,1,1);plot(y4);grid ontitle('對(duì)加三余弦混合噪聲的音樂信號(hào)進(jìn)行濾波之后的音譜')xlabel('時(shí)間')ylabel('幅度')subplot(2,1,2);plot(w2,abs(fy4);grid ontitle('對(duì)加三余弦混合噪聲的音樂信

27、號(hào)進(jìn)行濾波之后的頻譜')xlabel('頻率')ylabel('幅度')wp=0.15;ws=0.17;rp=1;rs=15;n,wc=buttord(wp,ws,rp,rs);b,a=butter(n,wc);hd,w=freqz(b,a);%設(shè)計(jì)巴特沃斯濾波器y5=filter(b,a,y3);%對(duì)加白噪聲的音樂信號(hào)進(jìn)行濾波n3=length(y5);fy5=fft(y5,n3);w3=2/n3*0:n3-1;figuresubplot(2,1,1);plot(y5);grid ontitle('對(duì)加白噪聲的音樂信號(hào)進(jìn)行濾波之后的音譜'

28、;)xlabel('時(shí)間')ylabel('幅度')subplot(2,1,2);plot(w2,abs(fy5);grid ontitle('對(duì)加白噪聲的音樂信號(hào)進(jìn)行濾波之后的頻譜')xlabel('頻率')ylabel('幅度')sound(y2,fs)sound(y3,fs)sound(y4,fs)sound(y5,fs)對(duì)于加三余弦混合噪聲的音樂信號(hào)，播放的過程中，背后一直有一個(gè)音調(diào)很高的“滴”的聲音，這是由于余弦信號(hào)頻率較高；對(duì)于加白噪聲的音樂信號(hào)，播放過程中，背后一直有“絲啦絲啦”的干擾噪聲，就像信號(hào)不

29、好時(shí)聽收音機(jī)的聲音；根據(jù)濾波后的音樂播放效果可知，三余弦混合噪聲可以完全濾去，這是由于余弦信號(hào)不在音樂信號(hào)所在的頻率處；白噪聲卻不能完全濾去，這是由于白噪聲在低頻處也有分布，音樂信號(hào)也存在低頻處，因此低頻處的白噪聲無法濾除，所以聽到的音樂信號(hào)仍然會(huì)有干擾的聲音，但是效果比濾波前要好。6 音樂信號(hào)的幅頻濾波及相頻分析 設(shè)計(jì)低通濾波器（可自行選擇不同的截止頻率），濾除原始音樂信號(hào)的高頻信息，觀察濾波前后的幅度頻譜，并比較濾波前后的音樂效果，感受高頻信息對(duì)音樂信號(hào)的影響； 設(shè)計(jì)高通濾波器（可自行選擇不同的截止頻率），濾除原始音樂信號(hào)的低頻信息，觀察濾波前后的幅度頻譜，并比較濾波前后的音樂效果，感受低

30、頻信息對(duì)音樂信號(hào)的影響； 選取兩端不同的音樂信號(hào)，分別將其幅度譜與相位譜交叉組合構(gòu)成信的音樂信號(hào)，播放比較組合后的音樂與原始音樂，感受相頻信息對(duì)音樂信號(hào)的影響。clear all;close all;clcy,fs,bit=wavread('e:music22');y1=y(:,1);wp=0.06;ws=0.09;rp=1;rs=15;n,wc=buttord(wp,ws,rp,rs);b,a=butter(n,wc);hd,w=freqz(b,a);%設(shè)計(jì)巴特沃斯低通濾波器y2=filter(b,a,y1);%濾除音樂信號(hào)的高頻信息n1=length(y1);n2=leng

31、th(y2);fy1=fft(y1,n1);fy2=fft(y2,n2);w1=2/n1*0:n1-1;w2=2/n2*0:n2-1;figure%繪制低通濾波前后音樂信號(hào)的幅頻特性subplot(2,1,1);plot(w1,abs(fy1);grid ontitle('低通濾波前音樂信號(hào)的頻譜')xlabel('頻率')ylabel('幅度')subplot(2,1,2);plot(w2,abs(fy2);grid ontitle('低通濾波后音樂信號(hào)的頻譜')xlabel('頻率')ylabel('幅

32、度')wp=0.06;ws=0.09;rp=1;rs=15;n,wc=buttord(wp,ws,rp,rs);b,a=butter(n,wc,'high');hd,w=freqz(b,a);%設(shè)計(jì)巴特沃斯高通濾波器y3=filter(b,a,y1);%濾除音樂信號(hào)的低頻信息n3=length(y3);w3=2/n3*0:n3-1;fy3=fft(y3,n3);figure%繪制高通濾波前后音樂信號(hào)的幅頻特性subplot(2,1,1);plot(w1,abs(fy1);grid ontitle('高通濾波前音樂信號(hào)的頻譜')xlabel('頻率

33、')ylabel('幅度')subplot(2,1,2);plot(w3,abs(fy3);grid ontitle('高通濾波后音樂信號(hào)的頻譜')xlabel('頻率')ylabel('幅度')sound(y1,fs)sound(y2,fs)sound(y3,fs)低通濾波后的音樂信號(hào)音調(diào)較原音樂信號(hào)低，高通濾波后的音樂信號(hào)音調(diào)較原音樂信號(hào)高。clear all;close all;clcx1,fs1,bit1=wavread('e:music44');x2,fs2,bit2=wavread('e

34、:music22');y11=x1(:,1);y1=y11(1:300000);y22=x2(:,1);y2=y22(1:300000);n1=length(y1);fy1=fft(y1,n1);w1=2/n1*0:n1-1;n2=length(y2);fy2=fft(y2,n2);w2=2/n2*0:n2-1;figure %畫兩信號(hào)的幅度譜subplot(2,2,1);plot(w1,abs(fy1);grid ontitle('音樂信號(hào)44的幅度譜')xlabel('frequency/hz')ylabel('magnitude')subplot(2,2,2);plot(w2,abs(fy2);grid ontitle('音樂信號(hào)22的幅度譜')xlabel('frequency/hz')ylabel('magnitude')subplot(2,2,3);plot(w1,angle(fy1);grid on %畫兩信