MATLAB音乐合成

一、实验目的1.熟悉MATLAB的软件和语言指令的使用；2.学习利用MATLAB进行连续信号的时域、频域分析；3.通过电子音乐合成方面的练习增进对傅里叶级数的理解。二、实验内容1.请根据《东方红》片断的简谱和“十二平均律”计算出该片断中各个乐音的频率，在MATLAB中生成幅度为1、抽样频率为8kHz的正弦信号表示这些乐音。请用sound函数播放每个乐音，听一听音调是否正确。最后用这一系列乐音信号拼出《东方红》片断，注意控制每个乐音持续的时间要符合节拍，用sound播放你合成的音乐，听起来感觉如何？俗话说万事开头难，由于这一题是我自己编写的第一个MATLAB程序，之前都是看老师的演示程序，对MATLAB指令并不熟练，所以这样一个很简单程序我也花了不少功夫。编写过程中，我主要还是对照着课本例题，一步一步，从生成抽样时间到编辑音调。最后在help的帮助下，参看了help中的Hallelujah范例，学会了sound函数的用法。代码如下：clearall,closeall,clc;t=linspace(0,4-1/8000,4*8000)';%4小节每小节1s,以4分音符为1拍，每小节2拍y=0*t;y(t0.5)=sin(392*2*pi*t(t0.5));%第一个音'5'，下依次y(0.5t&t0.75)=sin(392*2*pi*t(0.5t&t0.75));y(0.75=t&t1)=sin(440*2*pi*t(0.75=t&t1));y(1t&t1.5)=sin(293.66*2*pi*t(1t&t1.5));y(2=t&t2.5)=sin(261.63*2*pi*t(2=t&t2.5));y(2.5t&t2.75)=sin(261.63*2*pi*t(2.5t&t2.75));y(2.75=t&t3)=sin(220*2*pi*t(2.75=t&t3));y(3t&t3.5)=sin(293.66*2*pi*t(3t&t3.5));sound(y,8000)%播放音乐运行程序可以正确播放东方红前4小节，并且符合节拍。不过，确实存在第二题中所说的‘啪’声，并且两个同样的音之间没有分开。在之后的编程过程中我发现，这第一题我的方法其实很笨，而且代码利用率不高。于是在后边题中我对程序进行了改进，具体的改进方法我将在后边题中说明。2.你一定注意到(1)的乐曲中相邻乐音之间有“啪”的杂声，这是由于相位不连续产生了高频分量。这种噪声严重影响合成音乐的质量，丧失真实感。为了消除它，我们可以用包络修正每个乐音，以保证在乐音的邻接处信号幅度为零。此外建议用指数衰减的包络来表示1。消除‘啪’声的原理是让乐音的邻接处信号幅度为零，因此我在每个音波形开始和结束时与一个从0到1递增或从1到0递减的序列相点乘，实现渐变，保证了乐音邻接处信号幅度为0。clearall,closeall,clc;%优化代码后的第二问t1=linspace(0,0.25-1/8000,2000)';%修正8分音符包络set1=ones(2000,1);set1(t10.05)=20*t1(t10.05);set1(t10.15)=10*(0.25-t1(t10.15));t2=linspace(0,0.5-1/8000,4000)';%修正4分音符包络set2=ones(4000,1);set2(t20.05)=20*t1(t20.05);set2(t20.4)=10*(0.5-t2(t20.4));y0=zeros(4000,1);y1=set2.*note4(392);%G1y2=set1.*note8(392);%G1y3=set1.*note8(440);%A1y4=set2.*note4(293.66);%D1y5=set2.*note4(261.63);%C1y6=set1.*note8(261.63);%C1y7=set1.*note8(220);%A0y8=set2.*note4(293.66);%D1y=[y1;y2;y3;y4;y0;y5;y6;y7;y8;y0];sound(y,8000)functionnote=note4(x)%4分音符函数x为频率t=linspace(0,0.5-1/4000,4000)';note=sin(x*2*pi*t);functionnote=note8(x)%8分音符函数x为频率t=linspace(0,0.25-1/2000,2000)';note=sin(x*2*pi*t);从代码可以看出，我在完成此题程序的编写时，对题（1）中的代码进行了优化。我先构造两个函数note4和note8，根据各自的频率参数分别产生对应音调的4分音符和8分音符波形序列，然后再对各自的包络修正。最后将所有音的串起来，即得到了最终音乐。下为包络修正序列：3.请用最简单的方法将(2)中的音乐分别升高和降低一个八度。（提示：音乐播放的时间可以变化）再难一些，请用resample函数（也可以用interp和decimate函数）将上述音乐升高半个音阶。（提示：视计算复杂度，不必特别精确）此题在（2）的基础上添加一些代码即可实现，相同部分我就不重复附上，只说明后边附加的部分。开始时我采用的是将原波形序列点数增加或减少一倍来实现的。后发现，其实只需要调整sound函数的采样频率即可间接实现将频率翻倍或减半的效果。z1=ones(16000,1);%升高一个八度form=1:16000z1(m)=y(2*m);endsound(y,8000*2);%通过提高采样频率升高一个八度sound(y,8000/2);%通过降低采样频率降低一个八度利用resample函数可以很方便地变更原波形的频率。它的原理是通过插值P和抽取Q，以Q/P为倍数改变w=resample(y,10000,10595);%p/q=1/1.0595即升半音sound(w,8000);4.试着在(2)的音乐中增加一些谐波分量，听一听音乐是否更有“厚度”了？注意谐波分量的能量要小，否则掩盖住基音反而听不清音调了。（如果选择基波幅度为1，二次谐波幅度0:2，三次谐波幅度0:3，听起来像不像象风琴？）题（2）代码不用修改，直接在note4和note8两个频率波形产生函数中加入谐波分量即可，函数名为note4a和note8a。代码如下functionnote=note4a(x)%加谐波后的4分音符函数，x为频率t=linspace(0,0.5-1/4000,4000)';note=1*sin(x*2*pi*t)+0.2*sin(2*x*2*pi*t)+0.3*sin(3*x*2*pi*t);functionnote=note8a(x)%加谐波后的8分音符函数，x为频率t=linspace(0,0.25-1/2000,2000)';note=1*sin(x*2*pi*t)+0.2*sin(2*x*2*pi*t)+0.3*sin(3*x*2*pi*t);加入谐波分量后声音明显感觉有立体感了。5.自选其它音乐合成，例如贝多芬第五交响乐的开头两小节我选的是国歌的前奏部分，这部分既有8分音符，也有4分音符，适合做测试合成的音乐。代码如下：clearall,closeall,clc;%国歌前奏部分t1=linspace(0,0.25-1/8000,2000)';%修正8分音符set1=ones(2000,1);set1(t10.05)=20*t1(t10.05);set1(t10.15)=10*(0.25-t1(t10.15));t2=linspace(0,0.5-1/8000,4000)';%修正4分音符set2=ones(4000,1);set2(t20.05)=20*t1(t20.05);set2(t20.4)=10*(0.5-t2(t20.4));y0=zeros(4000,1);y1=set1.*note8(261.63);%C1y2=set1.*note8(329.63);%E1y3=set1.*note8(392);%G1y4=set2.*note4(440);%A1y5=set2.*note4(392);%G1y6=set2.*note4(329.63);%E1y7=set2.*note4(261.63);%C1y8=set1.*note8(392);%G1y9=set1.*note8(196);%G0y=[y1;y2;y3;y3;y4;y5;y2;y1;y8;y8;y6;y7;y9;y9;y9;y9;y7;y0];sound(y,8000);这一题比较简单，完成得比较顺利。6.先用wavread函数载入光盘中的fmt.wav文件，播放出来听听效果如何？是否比刚才的合成音乐真实多了？在阅读了helpwavread后，很轻松地实现了载入并播放。代码如下：[Y,FS]=wavread('fmt.wav');sound(Y,FS);与真实的吉他音乐相比，之前所做的合成音乐很不真实，因为它缺乏包络的变化（即真实乐器的振幅特点）以及频率丰富的和弦、揉弦等，而这些是真实的乐器很容易做到的。7.你知道待处理的wave2proc是如何从真实值realwave中得到的么？这个预处理过程可以去除真实乐曲中的非线性谐波和噪声，对于正确分析音调是非常重要的。提示：从时域做，可以继续使用resample函数。这一题可以说是我在做大作业中遇到的第一个坎。开始我没想明白提示中所说的resample函数在此题中如何使用。提示说从时域做，我的第一反应是用自相关，但后来想到这样只能检出波形而不能滤掉非线性谐波和噪声。在仔细观察了realwave和wave2proc的波形后，我发现其处理后恰好是10个周期，波形近乎相同。于是我想到如果对于每个周期的波形来说，非线性谐波和噪声都是随机的。如此一来，如果将很多周期的波形叠加，最后再除以周期数，那么所得的取“平均”后的周期波形，非线性谐波和噪声的影响因抵消而减小，而有价值的波形本身不受影响。这样，resample函数就派上了用场，可以利用它将realwave抽样点变为10倍，即10个周期，每个周期的抽样点数量都为整数243。如此只要将realwave每次向左错开1个周期然后相加，一共进行10次，即可得到取“平均”后的波形。再利用resample函数将抽样频率还原即可。代码如下：clearall,closeall,clc;load('guitar.mat');y=resample(realwave,10,1);%用10倍采样频率采样z=y;x=zeros(2430,1);form=1:10%一共10个周期z=[z(244:2430);z(1:243)];%左移一个周期x=x+z;%累加endx=x/10;%求平均y=resample(x,1,10);%还原回采样频率运行后发现所得波形y与wave2proc一模一样，分毫不差，当时成就感油然而生。可能我采用的方法与老师的方法相同，才能得到如此一致的结果。y和wave2proc的波形如下：8.这段音乐的基频是多少？是哪个音调？请用傅里叶级数或者变换的方法分析它的谐波分量分别是什么。提示：简单的方法是近似取出一个周期求傅里叶级数但这样明显不准确，因为你应该已经发现基音周期不是整数（这里不允许使用resample函数）。复杂些的方法是对整个信号求傅里叶变换（回忆周期性信号的傅里叶变换），但你可能发现无论你如何提高频域的分辨率，也得不到精确的包络（应该近似于冲激函数而不是sinc函数），可选的方法是增加时域的数据量，即再把时域信号重复若干次，看看这样是否效果好多了？请解释之。首先在不增加时域数据量的情况下，采用教材的矩阵计算法，对整个信号求傅里叶变换，代码如下：clearall,closeall,clc;%未增加时域数据量load('guitar.mat');N=243;t=linspace(0,(N-1)/8000,N);f=wave2p