清华电子系大二小学期MATLAB音乐合成大作业（完整项目代码+文档说明）

# MATLAB 音乐合成大作业 无15 向明义 2021012760 > 源代码及资源文件均放在 src/ 目录下，代码以 hw_x_x.m 命名，还有一些自定义函数。 > > 生成的音频文件放在 results/ 目录下。 > > report.md 中的图片放在 images/ 目录下。 > > 软件版本为 MATLAB R2022b。 ## 实验目的 - 增进对傅里叶级数的理解。 - 熟练运用 MATLAB 基本指令。 - 了解乐理和电子音乐的基本知识。 ## 实验内容 ### 0. 前置知识 首先要知道怎么合成do re mi。 一段时长1s的523Hz频率的正弦波，发出的就是小字二组c的do音，如果用8192Hz的采样频率，1s要采8192个点，用 t 表示时间采样点，则每个点上的音的幅度是 $\sin(2\pi ft)$。 这样一段代码可以生成一个do音： ```matlab Fs = 8192; T = 0.5; t = linspace(0, T, Fs*T); % Fs*T是采样点数 f = 523.25; y = sin(2*pi*f*t); sound(y); ``` 而幅度的大小就是音量的大小，如果给 $y$ 套上一个包络，就可以得到音量随包络渐变的效果。 谐波就是把另外频率的波叠加上去。 ### 1. 简单的合成音乐 #### 1.1 东方红第一小节 考虑到乐曲中用到的一般不会超过选定调子内的低一组和高一组音，所以计算3*7=21个频率备用： ```matlab baseA = 220; realBase = baseA * 2^(8/12); % 1=F调 % 乐曲中用到的基本是低一组和高一组 diffs = [-12,-10,-8,-7,-5,-3,-1,0,2,4,5,7,9,11,12,14,16,17,19,21,23]; baseFreq = realBase*2.^(diffs/12); mid = @(x) baseFreq(x+7); % mid(1) = F low = @(x) baseFreq(x); high = @(x) baseFreq(x+14); ``` 然后按照简谱写出需要的音调（第二列表示拍数）： ```matlab song = [ % 按照2拍分节 mid(5),1; mid(5),0.5; mid(6),0.5; mid(2),2; mid(1),1; mid(1),0.5; low(6),0.5; mid(2),2; ]; ``` 最后枚举每一行，将对应的频率和时间按照前置知识部分所述转换为采样数据拼接起来即可。 代码详见 `src/1_1.mlx`。对应音频文件在 `results/1_1.wav`。 这样合成的音乐，音量大小处处相同，同时在两个音调切换处有明显的“啪”音。 #### 1.2 加入包络 相邻乐音间的“啪”杂音是由于两个音频率不一致，在转换时产生了高频分量，为了消除它，需要给前一个音的末尾加上渐弱消音：  我们知道 $xe^{-x}$ 的图像可以很好的满足这个特征曲线： <img src=".\images\image-20230729203231164.png" alt="image-20230729203231164" style="zoom: 50%;" /> 经过简单的调参，确定使用 $6xe^{-8x}$ 作为包络，使用的方法为将 $t/T$ 作为 $x$ 传入函数然后相乘。 此时啪声已经可以较好地消去，为了使乐音更有渐进感，将相邻两乐音做迭接。具体方法为将原乐音延长一定比例，然后将后一个乐音与延长部分叠加。 关键改动如下： ```matlab shiftLen = 0; for k = 1:rows f = song(k,1); T = beatTime * song(k,2) * shiftRatio; % 这里对单音做延长 t = linspace(0, T, Fs*T)'; res = sin(2*pi*f*t) .* envelope(t/T); % 这里乘上包络 tunes = [ tunes(1:end-shiftLen); tunes(end-shiftLen+1:end) + res(1:shiftLen); res(shiftLen+1:end); ]; % 拼接 shiftLen = round(Fs * beatTime * song(k, 2) * (shiftRatio - 1)); % 取到延长部分的点数 end ``` 最终波形： <img src=".\images\image-20230729210914376.png" alt="image-20230729210914376" style="zoom: 67%;" /> #### 1.3 改变音调 升降八度最简单的方法是在 `sound` 输出时将采样频率加倍或减半： ```matlab sound(tunes, Fs*2); % 升八度（频率加倍，时间减半） ``` ```matlab sound(tunes, Fs/2); % 降八度（频率减半，时间加倍） ``` 对原音乐进行 `resample` 重新采样，如果采样频率上升（采样点数会增多），但播放仍按原采样率播放，则相当于原音乐的频率降低了。 因此要升高半个音阶，则采样频率应降低 $2^{1/12}$： ```matlab tunes = resample(tunes, Fs, round(Fs*2^(1/12))); ``` #### 1.4 加入谐波 谐波频率是基波的整数倍。假设基波幅度为1，设置二次谐波相对幅度为0.3，三次谐波相对幅度为0.15： ```matlab harmoConfs = [1, 0.3, 0.15]; res = sin(2*pi*f*t*(1:length(harmoConfs))) * harmoConfs' .* envelope(t/T) ; ``` <img src=".\images\image-20230729215936949.png" alt="image-20230729215936949" style="zoom: 50%;" /> 将波形放大查看，可以发现不再是完全的单频波。 #### 1.5 自选音乐合成 选用《传说的世界》前半部分。为了更好听，把重叠部分扩到了2倍，因此 shiftLen 部分做了一定修改。 参见 `src/hw_1_5.m` 及 `results/1_5.wav` ### 2. 用傅里叶级数分析音乐 #### 2.1 播放fmt.wav 使用 `audioread('fmt.wav')` 载入，再以 8K Hz采样率播放。效果明显更贴近真实音效。 #### 2.2 消除噪声  这是一段周期性的音频，由于噪声是偏随机的，因此多次取平均可以消除噪声。 realwave中有约10个周期，将这10个周期取平均，采用 resample 的方法，设原音频采样长度为 len，将原音频以10倍频率重新采样，这样新采样数据的每 len 个点就是原音频的一个周期，叠加之后复制10次，再以1/10频率重新采样即得到 wave2proc。 ```matlab wav = resample(realwave, 10, 1); % 10倍频重新采样 len = length(realwave); res = zeros([len 1]); for k = 1:10 res = res + wav((k-1)*len+1:k*len); % 周期叠加 end res = repmat(res, 10, 1); res = resample(res, 1, 10); ```  #### 2.3 傅里叶分析谐波分量 离散傅里叶变换与分量的关系： ```matlab X1 = fft(res(1:N)); plot([0:N-1]*Fs/N, abs(X1)); ``` 取一个周期做快速傅里叶变换，得到： <img src=".\images\image-20230730163659663.png" alt="image-20230730163659663" style="zoom:50%;" /> 由傅里叶变换的知识可知，这样得到的是各频率分量幅度的包络。为了得到更准确的变换，我们把信号重复若干次再做变换，以下是分别取1、10、1000个周期得到的变换图：  随着时域周期数增加，频谱越来越近似冲激函数。 接下来根据变换提取各分量的频率以及幅度，由于采样频率为8000Hz，所以谐波频率不超过4000Hz。 首先观察一下，基波分量在 200Hz 到 400Hz 之间，取范围内幅度最大值为基波频率： ```matlab Rg = round([200,400]/Fs*N); [~, base] = max(abs(X(Rg(1):Rg(2)))); base = base + Rg(1) - 2; % 注意matlab下标是从1开始的 freq = base/N*Fs; ``` 得到 $freq = 329.22~Hz$ 。 然后按照基频的整数倍取点，得到各次谐波分量的相对强度： ``` 0 0.0739 1 1.0000 2 1.4572 3 0.9587 4 1.0999 5 0.0523 6 0.1099 7 0.3589 8 0.1240 9 0.1351 10 0.0643 11 0.0019 12 0.0058 ``` #### 2.4 自动分析乐曲音调和节拍 <img src=".\images\image-20230730192819258.png" alt="image-20230730192819258" style="zoom: 50%;" /> 给出这样一段音频，根据波形自动划分它的节拍，从而得到每一个音调片段，对片段做傅里叶分析得到它的基频、谐波幅度等信息保存下来，从而可以合成与该乐器相似的曲调。 ##### 2.4.1 划分节拍 划分节拍的方式参考了谷源涛老师2021年信号与系统大作业中的节奏点划分方法： 总体分为五步：平方、加窗、差分、半波整流、提取峰值。下面逐步阐释： - 平方 平方得到能量（这一步换成绝对值好像也行） ```matlab y1 = wav.^2; ``` - 加窗 这里使用matlab的barthannwin窗，窗的长度取采样率的1/10，时域做卷积，即频域相乘。�