matlab处理音乐数字信号.docx

1、matlab处理音乐数字信号数字信号处理学院：物理工程学院 姓名：金义飞 一、读取音乐信号并观察波形与频谱程序w,fs,bit=wavread(你是我的眼.wav); %音频信号w。采样率fs 精度bit wav=(w(:,1); %取第一列信号 sound(w,fs) %听取抽样后的wav音乐sound(w,2*fs) %听取抽样后的wav音乐sound(w,fs/2) %听取抽样后的wav音乐figure; %创建一个窗口subplot(2,1,1);plot(wav); %画出时域图 fwav=fft(wav); %对音频信号做傅里叶变换 lwav=round(length(fwav)/

2、2); %音频信号正半轴长度 nwav=0:lwav-1; wwav=nwav/(lwav);f=wwav/2*fs; %x轴坐标 实际频率 subplot(2,1,2);plot(wwav,abs(fwav(1:lwav); %画出音乐频域图图像现象说明：从音谱即时域谱可以看出音乐信号随时间的幅度变化；从频谱可看到音乐信号主要分布在低频段，高频成分较少，在0.5pi 以后几乎无音乐信号二、音乐信号的抽取（减抽样）三倍减抽样1，程序clear all;close all;clc %清除存储空间，关闭所有窗口，清除命令w,fs,bit=wavread(你是我的眼.wav)wav=(w(:,1);

3、fwav=fft(wav);lwav=round(length(fwav)/2);j=0;d1=3; for i=1:d1:lwav; j=j+1; dwav1(j)=wav(i);end %本段每隔三个点取一个点，相当于三倍减抽样sound(dwav1,fs/d1);figure;subplot(2,1,1);plot(dwav1); fwav=fft(dwav1);lwav=round(length(fwav)/2);nwav=0:lwav-1;wwav=nwav/(lwav);f=wwav/fs;subplot(2,1,2);plot(wwav,abs(fwav(1:lwav);2图像2

4、0倍减抽样1，程序clear all;close all;clcw,fs,bit=wavread(你是我的眼.wav)wav=(w(:,1);fwav=fft(wav);lwav=round(length(fwav)/2);j=0;d1=20; for i=1:d1:lwav; j=j+1; dwav1(j)=wav(i);endsound(dwav1,fs/d1);figure;subplot(2,1,1);plot(dwav1); fwav=fft(dwav1);lwav=round(length(fwav)/2);nwav=0:lwav-1;wwav=nwav/(lwav);f=wwav

5、/fs;subplot(2,1,2);plot(wwav,abs(fwav(1:lwav); 分别选取抽样间隔为3倍频和20倍频进行减抽样，后者产生混叠，前者不混叠，当选取3倍频减抽样时，频率样间隔变大但未超过最小抽样间隔，尽管频谱有所延伸，但并不产生混叠，即减抽样后频谱上限未超过折叠频率；而选取20倍频减抽样时，由于时域抽样间隔变大超过最小抽样间隔，频谱延伸且产生混叠，即频率上限超过折叠频率，与3倍减抽样相比声音信号中多了很多杂音，这是由于抽样间隔的增大使信号频谱向高频搬移，且间隔越大，搬移越厉害。 三、音乐信号的AM调制1、调制频率wc=0.7*pi程序w,fs,bit=wavread(你

6、是我的眼.wav);wav=(w(:,1);lam=length(wav);fwav=fft(wav);lwav=length(fwav);nwav=0:lwav-1;amw=0.7*pi;cwav=cos(amw*0:lam-1);amwav=wav.*cwav;sound(amwav,fs);figure;subplot(2,1,1);plot(amwav);title(你是我的眼时域调制)famwav=fft(amwav);subplot(2,1,2);wwav=nwav/lwav;plot(2*wwav,abs(famwav(1:lwav);title(你是我的眼频域调制)图像2、调制

7、频率wc=pi/2程序w,fs,bit=wavread(你是我的眼.wav);wav=(w(:,1);lam=length(wav);fwav=fft(wav);lwav=length(fwav);nwav=0:lwav-1;amw=pi/2;cwav=cos(amw*0:lam-1);amwav=wav.*cwav;sound(amwav,fs);figure;subplot(2,1,1);plot(amwav);title(你是我的眼时域调制)famwav=fft(amwav);subplot(2,1,2);wwav=nwav/lwav;plot(2*wwav,abs(famwav(1:l

8、wav);title(你是我的眼频域调制)图像声音原信号频率上限大约为0.4pi,选取高频调制频率为0.7pi,低频调制频率为Pi/2; 高频调制后，信号频谱从原点处搬移至0.7pi处，且频谱左右各搬移0.7pi,由于频率上限为0.4pi,故搬移后信号最高频率超过折叠频率，从而产生了混叠，将会使信号失真;低频调制后，信号频谱从原点搬移至pi/20处，且频谱左右各搬移Pi/2,但搬移后信号频率上限并未超出折叠频率，故信号未发生混叠；播放调制后的声音信号，可以听出高频调制后声音信号有刺耳的噪音，而低频调制仍能清楚的听出其音调。因为高频调制产生混叠而低频没有产生混叠。 四、AM调制音乐信号的同步解调

9、1、巴特沃斯滤波程序clear all;close all;clcw,fs,bit=wavread(你是我的眼.wav);wav=(w(:,1);fwav=fft(wav);lwav=round(length(fwav)/2);nwav=0:lwav-1;wwav=nwav/(lwav);lam=length(wav);amw=pi/2;cwav=cos(amw*0:lam-1);amwav=wav.*cwav;pmwav=amwav.*cwav;fpmwav=fft(pmwav);figure;subplot(2,1,1);plot(pmwav);title(解调时域)fpmwav=fft(

10、pmwav);subplot(2,1,2);plot(wwav,abs(fpmwav(1:lwav);title(解调频域)T=1;wp=2/T*tan(0.3*pi/2);ws=2/T*tan(0.4*pi/2);Rp=1;Rs=15; N,wc=buttord(wp,ws,Rp,Rs,s);B,A=butter(N,wc,s);%求模拟滤波器系统函数的系数C,D=bilinear(B,A,1/T); %求数字滤波器系统函数的系数W=0:0.001*pi:0.5*pi;H=freqs(B,A,W);%模拟滤波器的频率响应Hd=freqz(C,D,W);%数字滤波器频率响应figuresubp

11、lot(3,1,1);plot(W/pi,abs(Hd);title(数字巴特沃斯滤波特图)grid on;y=filter(C,D,pmwav);%对信号滤波fy=fft(y);subplot(3,1,2);plot(abs(y);title(滤波后音频);subplot(3,1,3);plot(wwav,abs(fy(1:lwav);title(滤波后频谱);sound(y,fs);图像2、矩形窗滤波clear all;close all;clcw,fs,bit=wavread(你是我的眼.wav);wav=(w(:,1);fwav=fft(wav);lwav=round(length(f

12、wav)/2);nwav=0:lwav-1;wwav=nwav/(lwav);lam=length(wav);amw=pi/2;cwav=cos(amw*0:lam-1);amwav=wav.*cwav;pmwav=amwav.*cwav;fpmwav=fft(pmwav);N=39;wc=0.35*pi;%参数设置M=512;hd=ideal(N,wc);w2=boxcar(N);%设置矩形窗h2=hd.*w2;%求滤波器频率响应fh2=fft(h2,M);lx=length(fh2);wx1=(2/M)*0:M-1;db1=-20*log10(abs(fh2(1)./(abs(fh2)+e

13、ps);figure;subplot(3,1,1);plot(wx1,db1);title(矩形窗滤波特性);grid on;y0=pmwav;yy2=conv(y0,h2);subplot(3,1,2);plot(yy2);title(矩形窗滤波后音谱);grid on;fyy2=fft(yy2);l2=round(length(fyy2)/2);wx2=(0:l2-1)/l2;subplot(3,1,3);plot(wx2,abs(fyy2(1:l2);title(矩形窗滤波后频谱);grid on;sound(yy2,fs);定义函数文件function hd=ideal(N,wc)fo

14、r n=0:N-1 if n=(N-1)/2 hd(n+1)=wc/pi; else hd(n+1)=sin(wc*(n-(N-1)/2)/(pi*(n-(N-1)/2); endend图像4、布莱克曼窗滤波clear all;close all;clcw,fs,bit=wavread(你是我的眼.wav);wav=(w(:,1);fwav=fft(wav);lwav=round(length(fwav)/2);nwav=0:lwav-1;wwav=nwav/(lwav);lam=length(wav);amw=pi/2;cwav=cos(amw*0:lam-1);amwav=wav.*cwa

15、v;pmwav=amwav.*cwav;fpmwav=fft(pmwav);N=75;wc=0.15*pi;hd=ideal(N,wc);w2=blackman(N);%布莱克曼窗h2=hd.*w2;fh2=fft(h2);lx=length(h2);wx1=(0:lx-1)*2/lx;db1=-20*log10(abs(fh2(1)./(abs(fh2)+eps);figure;subplot(3,1,1);plot(wx1,db1);title(blackman窗滤波特性);grid on;y0=pmwav;yy2=conv(y0,h2);subplot(3,1,2);plot(yy2);

16、title(blackman窗滤波后音谱);grid on;fyy2=fft(yy2);l2=round(length(fyy2)/2);wx2=(0:l2-1)/l2;subplot(3,1,3);plot(wx2,abs(fyy2(1:l2);title(blackman窗滤波后频谱);grid on;sound(yy2,fs);函数文件function hd=ideal(N,wc)for n=0:N-1 if n=(N-1)/2 hd(n+1)=wc/pi; else hd(n+1)=sin(wc*(n-(N-1)/2)/(pi*(n-(N-1)/2); endend图像分析：巴特沃斯滤

17、波器的特点是通带内的频率响应曲线最大限度平坦，没有起伏，而在阻频带则逐渐下降为零。五、音乐信号的滤波去噪1、原始信号叠加余弦混合噪声clear all;close all;clcw,fs,bit=wavread(你是我的眼.wav)wav=(w(:,1);fwav=fft(wav);lwav=round(length(fwav)/2);nwav=0:lwav-1;wwav=nwav/(lwav);sf1=3000;sf2=5000;sf3=8000;l=length(wav);T=2/fs;t=0:T:(l-1)*T;s1=0.05*(cos(2*pi*sf1*t)+cos(2*pi*sf2*

18、t)+cos(2*pi*sf3*t); wav_s1=wav+s1;sound(wav_s1,fs)fwav_s1=fft(wav_s1);f_s1=fft(s1);figure;subplot(2,2,1);plot(s1);title(余弦噪声时域);subplot(2,2,3);plot(wwav,abs(f_s1(1:lwav);title(余弦噪声频域);subplot(2,2,2);plot(wav_s1);title(加噪信号时域);subplot(2,2,4);plot(wwav,abs(fwav_s1(1:lwav);title(加噪信号频域);用巴特沃斯滤波clear al

19、l;close all;clcw,fs,bit=wavread(你是我的眼.wav)wav=(w(:,1);fwav=fft(wav);lwav=round(length(fwav)/2);nwav=0:lwav-1;wwav=nwav/(lwav);sf1=3000;sf2=5000;sf3=8000;l=length(wav);T=2/fs;t=0:T:(l-1)*T;XX文库 - 让每个人平等地提升自我s1=0.05*(cos(2*pi*sf1*t)+cos(2*pi*sf2*t)+cos(2*pi*sf3*t);wav_s1=wav+s1;T=1;wp=2/T*tan(0.15*pi/

20、2);ws=2/T*tan(0.2*pi/2);Rp=1;Rs=10; N,wc=buttord(wp,ws,Rp,Rs,s);B,A=butter(N,wc,s);%求模拟滤波器系统函数的系数C,D=bilinear(B,A,1/T); %求数字滤波器系统函数的系数W=0:0.001*pi:0.5*pi;H=freqs(B,A,W);%模拟滤波器的频率响应Hd=freqz(C,D,W);%数字滤波器频率响应figuresubplot(3,1,1);plot(W/pi,abs(Hd);title(数字巴特沃斯滤波特图)grid on;y=filter(C,D,wav_s1);%对信号滤波fy=

21、fft(y);subplot(3,1,2);plot(abs(y);title(滤波后音频);subplot(3,1,3);plot(wwav,abs(fy(1:lwav);title(滤波后频谱);sound(y,fs);图像分析：三余弦噪声在频谱是三条独立的直线，且不处于低频段，播放时听到有刺耳的噪音，所以要把高频信号滤掉，巴特沃斯滤波器具有很好的这方面的特性。3、原始信号叠加随机白噪声并对其滤波。clear all;close all;clcw,fs,bit=wavread(你是我的眼.wav)wav=(w(:,1);fwav=fft(wav);lwav=round(length(fwa

22、v)/2);nwav=0:lwav-1;wwav=nwav/(lwav);l=length(wav);s2=(rand(1,l)-0.5)/5;%用rand语句设置幅度为0.5的随机白噪声f_s2=fft(s2);wav_s2=wav+s2;fwav_s2=fft(wav_s2);figure;subplot(2,2,1);plot(s2);title(随机白噪声时域);subplot(2,2,3);plot(wwav,abs(f_s2(1:lwav);title(随机白噪声频域);subplot(2,2,2);plot(wav_s2);title(加噪噪声时域);subplot(2,2,4)

23、;plot(wwav,abs(fwav_s2(1:lwav);title(加噪噪声频域);图片程序clear all;close all;clcw,fs,bit=wavread(你是我的眼.wav)wav=(w(:,1);fwav=fft(wav);lwav=round(length(fwav)/2);nwav=0:lwav-1;wwav=nwav/(lwav);l=length(wav);s2=(rand(1,l)-0.5)/5;%用rand语句设置幅度为0.5的随机白噪声f_s2=fft(s2);wav_s2=wav+s2;fwav_s2=fft(wav_s2);N=39;wc=0.35*

24、pi;%参数设置M=512;hd=ideal(N,wc);w2=boxcar(N);%设置矩形窗h2=hd.*w2;%求滤波器频率响应fh2=fft(h2,M);lx=length(fh2);wx1=(2/M)*0:M-1;db1=-20*log10(abs(fh2(1)./(abs(fh2)+eps);figure;subplot(3,1,1);plot(wx1,db1);title(矩形窗滤波特性);grid on;y0=wav_s2;yy2=conv(y0,h2);subplot(3,1,2);plot(yy2);title(矩形窗滤波后音谱);grid on;fyy2=fft(yy2)

25、;l2=round(length(fyy2)/2);wx2=(0:l2-1)/l2;subplot(3,1,3);plot(wx2,abs(fyy2(1:l2);title(矩形窗滤波后频谱);grid on;sound(yy2,fs);图像分析：白噪声分布均匀，声音不刺耳但是还是伴随着均匀幅度较小的噪音，加滤波器之后声音改善较多但是并不能完全滤除这是因为白噪声分布在所有频段内，想要滤完几乎是不可能的。思考：1，音乐信号的音调与信号的什么特征有关？ 答： 音调主要由声音的频率决定。音调随频率的升降而升降；听起来的特征：音调高：轻，短，细，音调低：重，长，粗。2音乐信号的音色与信号的什么特征有关

26、？ 答：音色是指声音的感觉特性，具体的来说，音色的类型是由振源的结构材质等特性有关，某些音色具有多种特性，例如人声的音色既柔软又暗淡，双簧管的音色既坚硬又明亮，圆号同时具有暗淡和明亮的音色。波形和音色是有密切关系的，确定的波形具有确定的音色。 3两种不同音色的音乐信号叠加混叠后，为何人耳还可以分辨？ 答：不同的发声体由于材料、结构不同，发出声音的音色也就不同，这样我们就可以通过音色的不同去分辨不同的发声体音色是声音的特色，根据不同的音色，即使在同一音高和同一声音强度的情况下，也能区分出是不同乐器或人发出的。同样的音量和音调上不同的音色就好比同样色度和亮度配上不同的色相的感觉一样。 4音乐信号的幅度与相位特征对信号有哪些影响？ 答：音乐信号的幅度与信号的响度有关，响度越大则幅度越大；由于相位中包含着信号传播过程中的反射，衍射，折射等信息，所以人耳可以根据相位信息判断发生物体的位置。