使用matlab进行简单音乐合成

1、信号与系统 综合实验之音乐合成(1)(1) 请根据东方红片断的简谱和“十二平均律”计算出该片断中各个乐音的 频率，在 MATLABMATLAB生成幅度为 1 1、抽样频率为 8kHz8kHz 的正弦信号表示这些乐 音。请用 soundsound函数播放每个乐音，听一听音调是否正确。最后用这一系列 乐音信号拼出东方红片断，注意控制每个乐音持续的时间要符合节拍， 用 sousou ndnd 播放你合成的音乐，听起来感觉如何代码如下：f =8000;t2=0:1/f:1;t4=0:1/f:;t8=0:1/f:;omg5=;omg6=;omg2=392;omg1=; omg6l=; 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)(2) 你一定注意到 (1)(1) 的乐曲中相邻乐音之间有“啪”的杂声，这是由于相位不 连续产生了高频分量。这种噪声严重影响合成音

3、乐的质量，丧失真实感。为 了消除它，我们可以用图 所示包络修正每个乐音，以保证在乐音的邻接处 信号幅度为零。此外建议用指数衰减的包络来表示。我采用的是指数衰减的包络。代码如下：f =8000;t2=0:1/f:1;t4=0:1/f:;t8=0:1/f:;omg5=;omg6=;omg2=392;omg1=;omg6l=;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(

4、-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),exp(t), 没有系数； 后来根据不同节拍，更改了不同的衰减系数，这样声音听起来感觉更加圆润。(3)(3) 请用最简单的方法将 (2)(2) 中的音乐分别升高和降低一个八度。(提示：音乐 播放的时间可以

5、变化)再难一些，请用 resampleresample 函数(也可以用 interpinterp 和 decimatedecimate 函数)将上述音乐升高半个音阶。(提示：视计算复杂度，不必 特别精确)答：最简单的方法是直接更改抽样频率 f f。将 f f 从 8K8K 改为 4K4K,则升高一个八度，并 且播放速度增快了一倍；将 f f 从 8k8k 改为 16k,16k,则降低一个八度，速度也变慢了一倍。 升高半个音阶，只须在(2 2)代码最后加一句 resampleresample ( m,1000,1059m,1000,1059)即可。f =8000;%改为4000或者16000t2

6、=0:1/f:1;t4=0:1/f:;t8=0:1/f:;omg5=;omg6=;omg2=392;omg1=;omg6l=;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=

7、exp(-1*t2).*sin(2*pi*omg2*t2);m=m1 m2 m3 m4 m6 m7 m8 m9; resample(m,1000,1059);sound(m);(4)(4) 试着在 (2)(2) 的音乐中增加一些谐波分量，听一听音乐是否更有“厚度”了注 意谐波分量的能量要小，否则掩盖住基音反而听不清音调了。(如果选择基 波幅度为 1 1 ，二次谐波幅度 0:20:2 ，三次谐波幅度 0:0: 3 3 ，听起来像不像象风琴) 代码如下：f =8000;t2=0:1/f:1;t4=0:1/f:;t8=0:1/f:;omg5=;omg6=;omg2=392;omg1=; omg6l=

8、;m1=exp(-2*t4).*sin(2*pi*omg5*t4)+*exp(-2*t4).*sin(2*pi*2*omg5*t4)+*exp(-2*t4).*sin( 2*pi*3*omg5*t4);m2=exp(-4*t8).*sin(2*pi*omg5*t8)+*exp(-4*t8).*sin(2*pi*2*omg5*t8)+*exp(-4*t8).*sin( 2*pi*3*omg5*t8);m3=exp(-4*t8).*sin(2*pi*omg6*t8)+*exp(-4*t8).*sin(2*pi*2*omg6*t8)+*exp(-4*t8).*sin( 2*pi*3*omg6*t8)

9、;m4=exp(-1*t2).*sin(2*pi*omg2*t2)+*exp(-1*t2).*sin(2*pi*2*omg2*t2)+*exp(-1*t2).*sin( 2*pi*3*omg2*t2);m6=exp(-2*t4).*sin(2*pi*omg1*t4)+*exp(-2*t4).*sin(2*pi*2*omg1*t4)+*exp(-2*t4).*sin( 2*pi*3*omg1*t4);m7=exp(-4*t8).*sin(2*pi*omg1*t8)+*exp(-4*t8).*sin(2*pi*2*omg1*t8)+*exp(-4*t8).*sin( 2*pi*3*omg1*t8)

10、;m8=exp(-4*t8).*sin(2*pi*omg6l*t8)+*exp(-4*t8).*sin(2*pi*2*omg6l*t8)+*exp(-4*t8).*si n(2*pi*3*omg6l*t8);m9=exp(-1*t2).*sin(2*pi*omg2*t2)+*exp(-1*t2).*sin(2*pi*2*omg2*t2)+*exp(-1*t2).*sin( 2*pi*3*omg2*t2);m=m1 m2 m3 m4 m6 m7 m8 m9;sound(m);加入谐波分量后，音色有所变化，感觉更加清脆一些。(5)(5) 自选其它音乐合成，例如贝多芬第五交响乐的开头两小节。 我选取

11、的是晴天的第一句代码如下：f=8000;t2=0:1/f:1;t4=0:1/f:;t8=0:1/f:;t=0:1/f:;omg1=392;omg2=440;omg3=;omg4=;omg5=;omg6=;omg7=; omg5l=;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(-1*t2).*sin(2*pi*omg1*t2); m5=exp(-2*t4).*sin(2*pi*omg2*t4);m6=exp(-2*t

12、2).*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=exp(-4*t8).*sin(2*pi*omg2*t8);m14=exp(-4*t8).*sin(2*pi*omg1*t8); m15=

13、exp(-2*t2).*sin(2*pi*omg5l*t2); m=m0 m1 m2 m3 m4 m5 m6 m7 m8m9 m10 m11 m12 m13 m14 m15; sound(m);(6)(6) 先用 wavreadwavread 函数载入光盘中的 文件，播放出来听听效果如何是否比刚才 的合成音乐真实多了x=wavread( );sound(x);( 7 7)你知道待处理的 wave2procwave2proc 是如何从真实值 realwaverealwave 中得到的么这个预处 理过程可以去除真实乐曲中的非线性谐波和噪声，对于正确分析音调是非常重要 的。提示：从时域做，可以继续使

14、用 resampleresample 函数。realwaverealwave 中的波形有十个周期，要除去其中的噪声可以采用时域求均值的方 法。步骤如下：首先用 resampleresample 函数将其采样率增大为十倍，再等分十份取平均， 然后重复这个平均后的波形十次，还原其长度，最后在用 resampleresample 函数还原采样率到 1/101/10 。代码如下：load( );wavetemp = zeros(length(realwave), 1); waveresampled = resample(realwave,10,1);for n = 1 : 10wavetemp = w

15、avetemp + waveresampled(n - 1) * length(realwave) + 1 : n * length(realwave) / 10;end mywave2proc = repmat(wavetemp, 10, 1); mywave2proc = resample(mywave2proc, 1, 10);(8)(8) 这段音乐的基频是多少是哪个音调请用傅里叶级数或者变换的方法分析它 的谐波分量分别是什么。提示：简单的方法是近似取出一个周期求傅里叶级数但 这样明显不准确，因为你应该已经发现基音周期不是整数(这里不允许使用 resampleresample 函数)。复

16、杂些的方法是对整个信号求傅里叶变换(回忆周期性信号 的傅里叶变换)，但你可能发现无论你如何提高频域的分辨率， 也得不到精确的 包络(应该近似于冲激函数而不是 sinsine e函数)，可选的方法是增加时域的数据 量，即再把时域信号重复若干次，看看这样是否效果好多了请解释之。答：共计 221-2=219221-2=219 个周期；基因频率约为 80008000 / / (219(219 / / 9)=9)=. .用 FFTFFT 寸该信号做 DFDF变换：load();my_wave2proc = repmat(wave2proc, 25, 1);Fs = 8000;Len gth = len

17、gth(my_wave2proc);-0.50-0.50100150200250figure;subplot(3, 1,1);plot(realwave);NFFT = 2人nextpow2(Le ngth);Y = fft(my_wave2proc, NFFT) / Len gth;amplitude = 2 * abs(Y(1 : NFFT / 2 + 1);freque ncy = Fs / 2 * lin space(0, 1 ,NFFT / 2 + 1);plot(freque ncy, amplitude);max_rate_of_grade, max_positi on = ma

18、x(amplitude(1 : 100);freque ncy(max_positi on)算得基频为，和前面通过周期估算得到的很相近。对照可得，为 C C 大调 mimi。(9)(9)再次载入，现在要求你写一段程序，自动分析出这段乐曲的音调和节拍！ 如果你觉得太难就允许手工标定出每个音调的起止时间， 再不行你就把每个音调 的数据都单独保存成一个文件， 然后让 MATLAMATLAB B对这些文件进行批处理。注意： 不允许逐一地手工分析音调。编辑音乐文件，推荐使用CoolEditCoolEdit 编辑软件。此问参考了学长的版本。看了版本才明白了算法。主要分2个步骤：1分割出单个音符，计算每个音

19、符的时间，2、计算每个音符的频率，把频率转化为音名。分割音符即将一个一个的音符片段从整个曲子中切下。音符的起始都伴随时域上瞬时能量的激增，即时域上波形幅度突然变大。利用这一点特征对音符进行提取。选择小步幅对进行扫描， 步幅选为秒最为适宜。 具体方法如下：在秒的范围内找极大点，而且是后续能量连续下降的极大点，这个极大点如果满足以下三个条件，则为有效的音符的起始点：1、 这个极大点比前一段的极大点能量咼出80%Pw*Pwp2、 这个极大点比整个乐曲的平均能量高；StStp+1OOO3、 距前面最近的一个音符的起始点时间差1/8秒以上。PwAvg分割的结果:EndPoint)做快速傅立叶变换，然后在

20、频谱中选出基频。关于自动寻找基频：一般的，基波不是幅度最大的频率，因为有很多音调的谐波，甚至非线性谐波的幅度大于基波。一般来说，一次谐波的幅度不太大，不会大于基波幅度的2倍，这样，在有可能出现基波的地方，禾U用每个点的幅度值除以该点的频率，也就是每个点 对于原点的斜率做为判定基波的标准，经过这样的简化处理，找到斜率最大的点， 基本上就能够判定基波出现在什么位置。FO=F_fft./omg;%!个点的幅度值除以该点的频率M1 M2=max(FO);%找到斜率最大的点，判定为基波Freq(p)=omg(M2)/2/pi;%己下频率值Name(p)=Freq2Name(Freq(p);喘巴频率值转化

21、为音名(关于函数Freq2Name(Freq):把频率与十几个标准音的频率作比较，选出最接近的一 个，作为这个频率对应的音名。)一共分割得30个音，和CoolEdit的分过度部分每个音符的持续时间 用于傅立叶变换的区间 能量修正(归一) 用于傅立叶变换的区间长度 计算频率并转化 对每个音符内部的一个区间Time=diff(StartPoint Len)/8000;En dPoi nt=StartPoi nt(2:Num) Len卜200; Power=Power/max(Power);STime=Time-200/8000;(我选这个音符StartPoi nt到音符分割5 5 4 4 3 3

22、2 2 1010o o o o o o o o o olalangsngs放，听到的音乐是很接近原来的的，仅仅有两个左右的音符有偏差。 代码如下：fmt=wavread( );%sound(fmt, 8000)Len=length(fmt);Avg=norm(fmt)/sqrt(Len);Num=0;Step=80;Pwp=0;Stp=-1000;St=1;for n=1:Step:Len-StepM=0;St=n;while M=1 & St+Step-1*Pwp & StStp+1000 & PwAvgNum=Num+1;StartPoint(Num)=St;Pow

23、er(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;hold on;plot(fmt);axis(0 length(fmt) 0 );title(o?);xlabel( Time );ylabel( Signal );for p=1:Num如此，便完成了乐曲的音调和节拍的自动分析。把此结果用函数CreateMusic处理并回%俞出结果:display(Sc

24、ore:);display(Name);display(Time);plot(StartPoint(p) StartPoint(p),0 ,r*-) ;endFreqLimit=170;for p=1:Numf=fmt(StartPoint(p):EndPoint(p)-1);F_fft,t,omg=FastFourier(f);omg(floor(STime(p)*(4000-FreqLimit):ceil(STime(p)*(4000+FreqLimit)=1e5;FO=F_fft./omg;M1 M2=max(FO);Freq(p)=omg(M2)/2/pi;Name(p)=Freq2N

25、ame(Freq(p);end display( Scoreo);display(Name);display(Time);UnknownMusic=Name;Time;Power;Music=CreateMusic( C ,UnknownMusic, Exp , Guitar+ );PlayAndPlot(Music);(10)(10) 用(7)(7) 计算出来的傅里叶级数再次完成第 (4)(4) 题，听一听是否像演奏 的吉 他演奏出来的答：在( 7 7)中计算出来的傅里叶级数中取前四项，后面的谐波并不显著将其略 去。基波和 2 2、3 3、4 4 次谐波的系数为如下，将其归一后再合成即可。合成的方法 参考了学长的报告。music_score = 5 5 6 2 1 1 6 2;0 0 0 0 0 0 -1 0;1 2 1 2;time = 5;harmo