《(完整word版)使用matlab进行简单音乐合成.doc》由会员分享，可在线阅读，更多相关《(完整word版)使用matlab进行简单音乐合成.doc（15页珍藏版）》请在三一文库上搜索。
1、信号与系统 综合实验之音乐合成(1) 请根据东方红片断的简谱和“十二平均律”计算出该片断中各个乐音的频率，在MATLAB 中生成幅度为 1 、抽样频率为 8kHz 的正弦信号表示这些乐音。请用 sound 函数播放每个乐音，听一听音调是否正确。最后用这一系列乐音信号拼出东方红片断，注意控制每个乐音持续的时间要符合节拍，用 sound 播放你合成的音乐，听起来感觉如何？代码如下：f =8000;t2=0:1/f:1;t4=0:1/f:0.5;t8=0:1/f:0.25;omg5=523.35;omg6=587.33;omg2=392;omg1=349.23;omg6l=293.66;m1=sin
2、(2*pi*omg5*t4);m2=sin(2*pi*omg5*t8);m3=sin(2*pi*omg6*t8);m4=sin(2*pi*omg2*t2);m6=sin(2*pi*omg1*t4);m7=sin(2*pi*omg1*t8);m8=sin(2*pi*omg6l*t8);m9=sin(2*pi*omg2*t2);m=m1 m2 m3 m4 m6 m7 m8 m9;sound(m);听的时候发现在相邻乐音之间有杂音，这是由于相位不连续造成的。(2) 你一定注意到 (1) 的乐曲中相邻乐音之间有“啪”的杂声，这是由于相位不连续产生了高频分量。这种噪声严重影响合成音乐的质量，丧失真实感。
3、为了消除它，我们可以用图 1.5 所示包络修正每个乐音，以保证在乐音的邻接处信号幅度为零。此外建议用指数衰减的包络来表示 。我采用的是指数衰减的包络。代码如下：f =8000;t2=0:1/f:1;t4=0:1/f:0.5;t8=0:1/f:0.25;omg5=523.35;omg6=587.33;omg2=392;omg1=349.23;omg6l=293.66;m1=exp(-2*t4).*sin(2*pi*omg5*t4);m2=exp(-4*t8).*sin(2*pi*omg5*t8);m3=exp(-4*t8).*sin(2*pi*omg6*t8);m4=exp(-1*t2).*si
4、n(2*pi*omg2*t2);m6=exp(-2*t4).*sin(2*pi*omg1*t4);m7=exp(-4*t8).*sin(2*pi*omg1*t8);m8=exp(-4*t8).*sin(2*pi*omg6l*t8);m9=exp(-1*t2).*sin(2*pi*omg2*t2);m=m1 m2 m3 m4 m6 m7 m8 m9;sound(m);第一次我采用的指数衰减没有时间前面的系数， 即每个都只乘 exp(t), 没有系数；后来根据不同节拍，更改了不同的衰减系数，这样声音听起来感觉更加圆润。(3) 请用最简单的方法将 (2) 中的音乐分别升高和降低一个八度。（提示：音乐
5、播放的时间可以变化）再难一些，请用 resample 函数（也可以用 interp 和 decimate 函数）将上述音乐升高半个音阶。（提示：视计算复杂度，不必特别精确）答：最简单的方法是直接更改抽样频率 f 。将f 从 8K改为 4K，则升高一个八度，并且播放速度增快了一倍； 将 f 从 8k改为 16k，则降低一个八度， 速度也变慢了一倍。升高半个音阶，只须在（ 2）代码最后加一句 resample （m,1000,1059 ）即可。f =8000;%改为 4000 或者 16000t2=0:1/f:1;t4=0:1/f:0.5;t8=0:1/f:0.25;omg5=523.35;omg
6、6=587.33;omg2=392;omg1=349.23;omg6l=293.66;m1=exp(-2*t4).*sin(2*pi*omg5*t4);m2=exp(-4*t8).*sin(2*pi*omg5*t8);m3=exp(-4*t8).*sin(2*pi*omg6*t8);m4=exp(-1*t2).*sin(2*pi*omg2*t2);m6=exp(-2*t4).*sin(2*pi*omg1*t4);m7=exp(-4*t8).*sin(2*pi*omg1*t8);m8=exp(-4*t8).*sin(2*pi*omg6l*t8);m9=exp(-1*t2).*sin(2*pi*o
7、mg2*t2);m=m1 m2 m3 m4 m6 m7 m8 m9;resample(m,1000,1059);sound(m);(4) 试着在 (2) 的音乐中增加一些谐波分量， 听一听音乐是否更有“厚度”了？注意谐波分量的能量要小，否则掩盖住基音反而听不清音调了。（如果选择基波幅度为 1 ，二次谐波幅度 0:2 ，三次谐波幅度 0:3 ，听起来像不像象风琴？）代码如下：f =8000;t2=0:1/f:1;t4=0:1/f:0.5;t8=0:1/f:0.25;omg5=523.35;omg6=587.33;omg2=392;omg1=349.23;omg6l=293.66;m1=exp(-
8、2*t4).*sin(2*pi*omg5*t4)+0.2*exp(-2*t4).*sin(2*pi*2*omg5*t4)+0.3*exp(-2*t4).*sin(2*pi*3*omg5*t4);m2=exp(-4*t8).*sin(2*pi*omg5*t8)+0.2*exp(-4*t8).*sin(2*pi*2*omg5*t8)+0.3*exp(-4*t8).*sin(2*pi*3*omg5*t8);m3=exp(-4*t8).*sin(2*pi*omg6*t8)+0.2*exp(-4*t8).*sin(2*pi*2*omg6*t8)+0.3*exp(-4*t8).*sin(2*pi*3*om
9、g6*t8);m4=exp(-1*t2).*sin(2*pi*omg2*t2)+0.2*exp(-1*t2).*sin(2*pi*2*omg2*t2)+0.3*exp(-1*t2).*sin(2*pi*3*omg2*t2);m6=exp(-2*t4).*sin(2*pi*omg1*t4)+0.2*exp(-2*t4).*sin(2*pi*2*omg1*t4)+0.3*exp(-2*t4).*sin(2*pi*3*omg1*t4);m7=exp(-4*t8).*sin(2*pi*omg1*t8)+0.2*exp(-4*t8).*sin(2*pi*2*omg1*t8)+0.3*exp(-4*t8)
10、.*sin(2*pi*3*omg1*t8);m8=exp(-4*t8).*sin(2*pi*omg6l*t8)+0.2*exp(-4*t8).*sin(2*pi*2*omg6l*t8)+0.3*exp(-4*t8).*sin(2*pi*3*omg6l*t8);m9=exp(-1*t2).*sin(2*pi*omg2*t2)+0.2*exp(-1*t2).*sin(2*pi*2*omg2*t2)+0.3*exp(-1*t2).*sin(2*pi*3*omg2*t2);m=m1 m2 m3 m4 m6 m7 m8 m9;sound(m);加入谐波分量后，音色有所变化，感觉更加清脆一些。(5) 自选
11、其它音乐合成，例如贝多芬第五交响乐的开头两小节。我选取的是晴天的第一句代码如下：f=8000;t2=0:1/f:1;t4=0:1/f:0.5;t8=0:1/f:0.25;t=0:1/f:0.125;omg1=392;omg2=440;omg3=493.88;omg4=523.25;omg5=587.33;omg6=659.25;omg7=698.45;omg5l=293.66;m0=0;m1=exp(-2*t4).*sin(2*pi*omg5*t4);m2=exp(-2*t4).*sin(2*pi*omg5*t4);m3=exp(-2*t4).*sin(2*pi*omg1*t4);m4=exp
12、(-1*t2).*sin(2*pi*omg1*t2);m5=exp(-2*t4).*sin(2*pi*omg2*t4);m6=exp(-2*t2).*sin(2*pi*omg3*t2);m7=exp(-2*t4).*sin(2*pi*omg5*t4);m8=exp(-2*t4).*sin(2*pi*omg5*t4);m9=exp(-2*t4).*sin(2*pi*omg1*t4);m10=exp(-2*t4).*sin(2*pi*omg1*t4);m11=exp(-4*t8).*sin(2*pi*omg2*t8);m12=exp(-4*t8).*sin(2*pi*omg3*t8);m13=ex
13、p(-4*t8).*sin(2*pi*omg2*t8);m14=exp(-4*t8).*sin(2*pi*omg1*t8);m15=exp(-2*t2).*sin(2*pi*omg5l*t2);m=m0 m1 m2 m3 m4 m5 m6 m7 m8 m9 m10 m11 m12 m13 m14 m15; sound(m);(6) 先用 wavread函数载入光盘中的 fmt.wav文件，播放出来听听效果如何？是否比刚才的合成音乐真实多了？x=wavread('fmt.wav');sound(x);（ 7）你知道待处理的 wave2proc 是如何从真实值 realwave 中
14、得到的么？这个预处理过程可以去除真实乐曲中的非线性谐波和噪声， 对于正确分析音调是非常重要的。提示：从时域做，可以继续使用resample函数。realwave 中的波形有十个周期， 要除去其中的噪声可以采用时域求均值的方法。步骤如下：首先用 resample 函数将其采样率增大为十倍， 再等分十份取平均，然后重复这个平均后的波形十次， 还原其长度，最后在用 resample 函数还原采样率到 1/10 。代码如下：load('Guitar.MAT');wavetemp = zeros(length(realwave), 1);waveresampled = resample(
15、realwave,10,1);forn = 1 : 10wavetemp = wavetemp + waveresampled(n - 1) * length(realwave) + 1 : n* length(realwave) / 10; endmywave2proc = repmat(wavetemp, 10, 1);mywave2proc = resample(mywave2proc, 1, 10);figure;subplot(3, 1, 1);plot(realwave);subplot(3, 1, 2);plot(wave2proc);subplot(3, 1, 3);plot(
16、mywave2proc);0.50-0.55010015020025000.50-0.55010015020025000.50-0.5050100150200250(8) 这段音乐的基频是多少？是哪个音调？请用傅里叶级数或者变换的方法分析它的谐波分量分别是什么。 提示：简单的方法是近似取出一个周期求傅里叶级数但这样明显不准确， 因为你应该已经发现基音周期不是整数 （这里不允许使用resample函数）。复杂些的方法是对整个信号求傅里叶变换（回忆周期性信号的傅里叶变换） ，但你可能发现无论你如何提高频域的分辨率，也得不到精确的包络（应该近似于冲激函数而不是sinc 函数），可选的方法是增加时域的
17、数据量，即再把时域信号重复若干次，看看这样是否效果好多了？请解释之。答：共计 221-2=219个周期；基因频率约为 8000 / (219 / 9) = 328.7671Hz.用 FFT对该信号做 DFT变换：load('Guitar.MAT');my_wave2proc = repmat(wave2proc, 25, 1);Fs = 8000;Length = length(my_wave2proc);NFFT = 2 nextpow2(Length);Y = fft(my_wave2proc, NFFT) / Length;amplitude = 2 * abs(Y(1
18、: NFFT / 2 + 1);frequency = Fs / 2 * linspace(0, 1 ,NFFT / 2 + 1);plot(frequency, amplitude);max_rate_of_grade, max_position = max(amplitude(1 : 100);frequency(max_position)0.080.070.060.050.040.030.020.01005001000150020002500300035004000算得基频为 329.1016Hz, 和前面通过周期估算得到的很相近。对照可得，为 C大调 mi。(9) 再次载入 fmt.w
19、av ，现在要求你写一段程序，自动分析出这段乐曲的音调和节拍！如果你觉得太难就允许手工标定出每个音调的起止时间， 再不行你就把每个音调的数据都单独保存成一个文件， 然后让 MATLAB 对这些文件进行批处理。注意：不允许逐一地手工分析音调。编辑音乐文件，推荐使用 CoolEdit" 编辑软件。此问参考了学长的版本。看了版本才明白了算法。主要分2 个步骤：1、分割出单个音符，计算每个音符的时间，2、计算每个音符的频率，把频率转化为音名。分割音符即将一个一个的音符片段从整个曲子中切下。音符的起始都伴随时域上瞬时能量的激增，时域上波形幅度突然变大。利用这一点特征对音符进行提取。选择小步幅对
20、fmt.wav 进行扫描，步幅选为0.01 秒最为适宜。具体方法如下：在即0.01秒的范围内找极大点， 而且是后续能量连续下降的极大点， 这个极大点如果满足以下三个条件，则为有效的音符的起始点：1、 这个极大点比前一段的极大点能量高出80%；Pw>1.8*Pwp2、 这个极大点比整个乐曲的平均能量高；St>Stp+10003、 距前面最近的一个音符的起始点时间差1/8 秒以上。Pw>Avg分割的结果：音符分割0.50.4langi0.3S0.20.10024681012Time4x 10一共分割得30 个音，和CoolEdit 的分析基本一致。过度部分每个音符的持续时间Tim
21、e=diff(StartPoint Len)/8000;用于傅立叶变换的区间EndPoint=StartPoint(2:Num) Len-200;能量修正（归一）Power=Power/max(Power);用于傅立叶变换的区间长度STime=Time-200/8000;计算频率并转化对每个音符内部的一个区间（我选这个音符的起始到结束前的0.025 秒，即StartPoint到 EndPoint ）做快速傅立叶变换，然后在频谱中选出基频。关于自动寻找基频：一般的，基波不是幅度最大的频率，因为有很多音调的谐波，甚至非线性谐波的幅度大于基波。一般来说，一次谐波的幅度不太大，不会大于基波幅度的2倍，
22、这样， 在有可能出现基波的地方，利用每个点的幅度值除以该点的频率，也就是每个点对于原点的斜率做为判定基波的标准，经过这样的简化处理，找到斜率最大的点，基本上就能够判定基波出现在什么位置。FO=F_fft./omg;%每个点的幅度值除以该点的频率M1 M2=max(FO);%找到斜率最大的点，判定为基波Freq(p)=omg(M2)/2/pi;%记下频率值Name(p)=Freq2Name(Freq(p);%把频率值转化为音名（关于函数 Freq2Name(Freq)：把频率与十几个标准音的频率作比较，选出最接近的一个，作为这个频率对应的音名。 ）display('Score ：
23、9;);display(Name);display(Time);%输出结果：#音名时长#音名时长#音名时长#音名时长1-60.16769-40.504417-60.38902510.47182-61.486610-40.74841850.268726-40.53363-70.476411-71.23731940.230527-60.46274-60.45911230.78692030.201028-70.5771520.458413-60.22362120.261029-60.6410630.450014-60.47132210.473130-51.39747-50.466015-60.472
24、823-70.44378-40.44691660.49792420.5484如此，便完成了乐曲的音调和节拍的自动分析。把此结果用函数 CreateMusic 放，听到的音乐是很接近原来的 fmt.wav 的，仅仅有两个左右的音符有偏差。代码如下：处理并回fmt=wavread('fmt.wav');%sound(fmt, 8000)Len=length(fmt);Avg=norm(fmt)/sqrt(Len);Num=0;Step=80;Pwp=0;Stp=-1000;St=1;forn=1:Step:Len-StepM=0;St=n;whileM=1 & St+Ste
25、p-1<=LenPw M=max(fmt(St:St+Step-1);St=St+M-1;endifPw>1.8*Pwp & St>Stp+1000 & Pw>AvgNum=Num+1;StartPoint(Num)=St;Power(Num)=Pw;Stp=St;endPwp=Pw;endTime=diff(StartPoint Len)/8000;STime=Time-200/8000;EndPoint=StartPoint(2:Num) Len-200;Power=Power/max(Power);figure;holdon ;plot(fmt);
26、axis(0 length(fmt) 0 0.6);title(' ò ?· ?· ?');xlabel('Time');ylabel('Signal'forp=1:Numplot(StartPoint(p) StartPoint(p),0 0.6,end);'r*-') ;FreqLimit=170;forp=1:Numf=fmt(StartPoint(p):EndPoint(p)-1);F_fft,t,omg=FastFourier(f);omg(floor(STime(p)*(4000-Freq
27、Limit):ceil(STime(p)*(4000+FreqLimit)=1e5;FO=F_fft./omg;M1 M2=max(FO);Freq(p)=omg(M2)/2/pi;Name(p)=Freq2Name(Freq(p);enddisplay('Score o');display(Name);display(Time);UnknownMusic=Name;Time;Power;Music=CreateMusic('C',UnknownMusic,'Exp', 'Guitar+');PlayAndPlot(Music);
28、(10) 用 (7) 计算出来的傅里叶级数再次完成第 (4) 题，听一听是否像演奏fmt.wav的吉他演奏出来的？答：在（ 7）中计算出来的傅里叶级数中取前四项，后面的谐波并不显著将其略去。基波和 2、3、4次谐波的系数为如下，将其归一后再合成即可。合成的方法参考了学长的报告。music_score = 5562116 2;000000-10;10.50.5210.50.52;time = 5;harmonic_coefficient = 0.0538 0.0753 0.0464 0.0482;harmonic_coefficient = harmonic_coefficient /norm(
29、harmonic_coefficient);diationic_scale = 220 * 2.(3 5 7 8 10 12 14 + 5) / 12);t = 0 : 1 / 8000 : time;music_wave = zeros(1, length(t);time_start = 0;speed = 0.5;forn = 1 : length(music_score)frequency = diationic_scale(music_score(1, n) * 2 music_score(2, n);time_end = time_start + music_score(3, n)
30、* speed ;envelope = GetEnvelope3(time, time_start, time_end);form = 1 : 4music_wave = music_wave + harmonic_coefficient(m) * sin(2 * pi * m * frequency * t) .* (t >= time_start & t < time_end) .* envelope;endtime_start = time_end;endmusic_wave = music_wave / max(music_wave);sound(music_wav
31、e, 8000);(11) 也许 (9) 还不是很像，因为对于一把泛音丰富的吉他而言， 不可能每个音调对应的泛音数量和幅度都相同。但是通过完成第 (8) 题，你已经提取出 fmt.wav中的很多音调， 或者说，掌握了每个音调对应的傅里叶级数，大致了解了这把吉他的特征。现在就来演奏一曲东方红吧。提示：如果还是音调信息不够，那就利用相邻音调的信息近似好了，毕竟可以假设吉他的频响是连续变化的。答：和上题类似。music_score = 55621162;000000-10;10.5 0.5 210.5 0.5 2;time = 5;load('HarmonicCoefficients.mat
32、');diationic_scale = 220 * 2.(3 5 7 8 10 12 14 + 5) / 12);t = 0 : 1 / 8000 : time;music_wave = zeros(1, length(t);time_start = 0;speed = 0.5;forn = 1 : length(music_score)frequency = diationic_scale(music_score(1, n) * 2 music_score(2,n);time_end = time_start + music_score(3, n) * speed ;envelop
33、e = GetEnvelope3(time, time_start, time_end);musical_name_index= round(12*log(frequency/220)/log(2)+ 13;form = 1 : 7music_wave = music_wave +harmonic_coefficient(musical_name_index,m) * sin(2* pi*m * frequency* t) .* (t >= time_start & t < time_end) .* envelope; endtime_start = time_end;endmusic_wave = music_wave / max(music_wave);plot(music_wave); sound(music_wave, 8000);两种吉他版本听起来都更加清脆悦耳一些。实验总结：这次试验总体感觉要比第一次难写。 一开始我不太会对于音乐信号进行傅里叶变化并且进行后续分析。 后来才考了学长的报告， 以及和同学讨论后才懂得对于这些函数的运用方法。最后由于时间以及能力所限，并没有做图形界面。