---

声音的向量表示
原理

1、向量x∈R^N表示时间区间上的音频信号 xi表示t=hi时的声压xi=αp(hi),i=1,…,N
2、每个xi称为样本
3、h(>0)为采样时间
4、1/h为采样率，典型的采样率为1/h=44100/sec或48000/sec
5、α被称为比率因子

（一 ~ 四）librosa学习点此处

（五 ~ 十七） librosa学习点此处

十八、写音频

（1）soundfile.write

在0.8.0以后的版本，librosa删除了librosa.output.write_wav(path, y, sr, norm=False)函数。因此写音频（将NumPy数组保存到WAV文件）应使用soundfile.write函数：

import soundfile
soundfile.write(file, data, samplerate)
#soundfile也可以用来读取音频：soundfile.read('D:\My life\music\some music/sweeter.mp3')

参数：

file：保存输出wav文件的路径
data：音频时间序列
sr：data的采样率

应用

#创建音频信号（numpy数组），以220Hz信号为例，将其传递给音频函数。
import numpy as np
import IPython.display as ipd

T = 5.0    # seconds，音频时长为5秒
t = np.linspace(0, T, int(T*sr), endpoint=False)   # time variable
x = 0.5*np.sin(2*np.pi*220*t)  # pure sine wave at 220 Hz

sr = 22050 # sample rate

#Saving the audio
soundfile.write('D:/My life/tone_220.wav', x, sr)

#Playing the audio
ipd.Audio(x, rate=sr) # load a NumPy array

#如果是0.8.0之前的版本，则用
①librosa.output.write_wav(path, y, sr, norm=False) #将NumPy数组保存到WAV文件。
参数：

路径：保存输出wav文件的路径
y：音频时间序列。
sr：y的采样率
norm：bool，是否启用幅度归一化。将数据缩放到[-1，+1]范围。

②librosa.ouput.annotation() 保存3列格式的注释。
③librosa.ouput.time_csv() 按CSV格式保存时间步骤。这可以用来存储打击跟踪器或分割算法的输出。

（2）缩放音频

音频信号每个样本的数值反应了振动离开平衡点的距离，即振幅。振幅越大，声音具有的能量越大，听起来响度就越大。对样本数值进行倍增，若倍增系数绝对值大于1则声音响度增加，若倍增系数绝对值小于1则声音响度下降。由于振动离开平衡点的距离是一个绝对值，因此倍增系数为负数时效果与倍增系数是其相反数时一致。
1、y=2y #增加一倍振幅，得到的音频响度略大于原音频
2、y=0.5y #减小振幅为原来一半，得到的音频响度略小于原音频
3、y=-y #翻转振幅，得到的音频响度与原音频一致
4、y=10*y #大幅增加振幅，得到的音频响度远大于原来音频

应用：

import matplotlib.pyplot as plt
plt.figure(figsize=(8, 6))#figsize用于定义子图大小
plt.subplots_adjust(left=None, bottom=None, right=None, top=None,wspace=0.2, hspace=1)#调整子图间距。参数：left：the left side of the subplots of the figure；right： the right side of the subplots of the figure；bottom：子图的底部 the bottom of the subplots of the figure；top：子图的高度；wspace：子图间空白部分的宽度；hspace：子图间空白部分的高度

y,sr = librosa.load("D:\My life\music\some music/sweeter.mp3",sr=None)    #y为长度等于采样率sr*时间的音频向量
plt.subplot(3, 2, 1)#subplot(行，列，第几副图)
librosa.display.waveplot(y, sr)    #创建波形图
plt.title('y')

y1=2*y
plt.subplot(3, 2, 3)
librosa.display.waveplot(y1, sr)    #创建波形图
plt.title('y1')

y2=0.5*y
plt.subplot(3, 2, 4)
librosa.display.waveplot(y2, sr)    #创建波形图
plt.title('y2')

y3=-y
plt.subplot(3, 2, 5)
librosa.display.waveplot(y3, sr)    #创建波形图
plt.title('y3')

y4=10*y
plt.subplot(3, 2, 6)
librosa.display.waveplot(y4, sr)    #创建波形图
plt.title('y4')

plt.show()    #显示波形图

附：

subplot

plt.subplots_adjust(left=None, bottom=None, right=None, top=None,wspace=0.2, hspace=1)#调整子图间距。
参数：

left：the left side of the subplots of the figure；
right： the right side of the subplots of the figure；
bottom：子图的底部 the bottom of the subplots of the figure；
top：子图的高度；
wspace：子图间空白部分的宽度；
hspace：子图间空白部分的高度

绘图代码：

fig = plt.figure()  #头
...具体内容
plt.show()    #尾

（3）线性组合和混音

numpy数组提供的线性运算功能可以用于实现混音。原理：

1、对多个音频信号进行线性运算y = a1x1 + a2x2 + … + akxk可以实现自动混音。
2、此时每个音频信号xk称为音轨
3、混合后的结果y称为混合
4、每个系数ak是音轨在混合中的权重

由于相同频率的人声和伴奏中容易撞频（尽管其谐波分量基本不同），使混合的人声和伴奏难以区分。因此对于双声道音频，常常采用偏置的方法，为两个声道中人声和伴奏设不同的权重（通常一个声道人声大，一个声道伴奏声大）参考资料

应用：

#尝试人声与伴奏的混合，取一段人声轨和一段伴奏轨，分别赋予权重(0.25,0.75)、(0.5,0.5)、(0.4,0.6)进行线性混合
y=0.25x1+0.75x2
y=0.5x1+0.5x2
y=0.6x1+0.4x2
#得到的音频为在伴奏上加上人声的效果
#分析：音频是线性相加，得到的波形图中无法直观分出人声和伴奏，但由于人声和伴奏中各频率谐波分量与伴奏中不同，线性相加在频域被分配到了各频率谐波上，所以可以在频谱图上清晰区分出来。为了得到足够清晰的混音音频，需要不断更改每个混音音轨的权重以得到最适合的混合模式
#拓展：虽然人声和伴奏中谐波分量基本不同，但相同频率上的相加也容易导致撞频的发生，使得混合的人声和伴奏难以区分，因此对于双声道音频，常常采用偏置的方法，为两个声道中人声和伴奏设不同的权重（通常一个声道人声大，一个声道伴奏声大）。形状为(2,n)的numpy数组可以用于容纳双声道音频信号。两个声道分别取权重(0.4,0.6)、(0.6,0.4)进行混合

（4）变声

①变速（更改采样率）

import librosa
import soundfile
y,sr = librosa.load("D:/My life/music/some music/sweeter.mp3")
# 通过改变采样率来改变音速，相当于播放速度X2
soundfile.write("D:/My life/music/some music/sweeter_resample.wav",y,sr*2)
#Ps:两倍速的声音好可爱！

②变调

import librosa
import soundfile
y,sr = librosa.load("D:/My life/music/some music/sweeter.mp3")
# 通过移动音调变声 ，14是上移14个半步， 如果是 -14 下移14个半步
b = librosa.effects.pitch_shift(y, sr, n_steps=14)
soundfile.write("D:/My life/music/some music/sweeter_pitch_shift.wav",b,sr)

男女声变调方法为进行频谱搬移（音频信号乘一个余弦函数）:

男女声音域范围（HZ）:
男低音：82.4（E）~ 329.6（e1） 女低音：164.8（e）~ 698.46（f2）
男中音：110（A） ~ 440（a4） 女中音： 196（g）~ 880（a2）
男高音：130.8（c）~ 523.25（c2） 女高音： 261.63（c1）~1046.5（c3）

③傅里叶变换变声

import librosa 
import matplotlib.pyplot as plt
import numpy as np
import soundfile
y,sr = librosa.load("D:/My life/music/some music/sweeter.mp3")

# stft 短时傅立叶变换
a = librosa.stft(y)
length = len(a)

# 改变或去除某些值，可以改变声音
r_a = a[10:length-10]

# istft 逆短时傅立叶变换，变回去
b = librosa.istft(r_a)

soundfile.write("D:/My life/music/some music/sweeter_stft.wav",b,sr)

# 以下是显示频谱图
fig = plt.figure()
plt.subplots_adjust(left=None, bottom=None, right=None, top=None,wspace=0.2, hspace=1)#调整子图间距。
s1 = fig.add_subplot(3,1,1)
s2 = fig.add_subplot(3,1,2)
s3 = fig.add_subplot(3,1,3)

s1.plot(y)
s2.plot(a)
s3.plot(b)

plt.show()

变声实例

（5）加载、拼接、截取、叠加、静音、调整音量

import numpy as np
import librosa
import soundfile

①加载音频文件

audio_path1 = 'D:/My life/music/some music/sweeter.mp3'
y1, sr1 = librosa.load(audio_path1)
dur1 = librosa.get_duration(y1, sr=sr1)
print('数据x类型和采样率sr类型', type(y1), type(sr1))
print('数据x尺寸和采样率', y1.shape, sr1)
print('该音频的时长为：', dur1)

audio_path2 = 'D:/My life/music/some music/糸.mp3'
y2, sr2 = librosa.load(audio_path2)
dur2 = librosa.get_duration(y2, sr=sr2)
print('数据x类型和采样率sr类型', type(y2), type(sr2))
print('数据x尺寸和采样率', y2.shape, sr2)
print('该音频的时长为：', dur2)

②音频拼接

audio_dst = np.hstack((y1, y2))  #在第一段音频后接上第二段音频
print('数据x尺寸和采样率', audio_dst.shape)
soundfile.write('D:/My life/music/some music/拼接.wav', audio_dst, sr2)

③音频截取

start = 12
duration = 20 #音频时长为20s
stop = start +duration
audio_dst = y2[start*sr2:stop*sr2] #第二段音频从12s处开始，截取到12+20s处。
soundfile.write('D:/My life/music/some music/截取-2.wav', audio_dst, sr2)

④音频叠加

num_min = min(len(y1), len(y2))
num_max = max(len(y1), len(y2))
audio_dst = np.zeros(num_max)
for i in range(num_max):
    if i<num_min:
        audio_dst[i] = y1[i] + y2[i]
    else:
        if len(y1)>=len(y2):
            audio_dst[i] = y1[i]
        else:
            audio_dst[i] = y2[i]
soundfile.write('D:/My life/music/some music/叠加.wav', audio_dst, sr2)#两段音乐同时播放

⑤生成静音音频

num = 3*sr2 
audio_dst = np.zeros(num)
soundfile.write('D:/My life/music/some music/静音.wav', audio_dst, sr2)#生成3s静音音频

⑥调整音频音量

soundfile.write('D:/My life/music/some music/调整音量.wav', y2*0.1, sr2)#减小音量：第二段音频音量*0.1,

十九、音乐

二十世纪八十年代，有专家研究巴赫《第一勃兰登堡协奏曲》的音乐信号时发现，音乐信号的功率谱与人类大脑生理信号的功率谱相似，符合1/f信号公式。还发现，音乐信号α越靠近数值1越好听，从科学上找到一个近似参数来判定音乐的悦耳度。2012年加拿大麦吉尔大学音乐系主任分析发现，音乐节奏也满足这个规律，α值为0.8。不同音乐体裁的α值不一样，所以也可以用这个数值反推音乐的风格体裁。不同作曲家风格音乐的α值不一样，但是作曲家们所作出来的各种风格体裁的音乐的参数是相似的。参考资料

（1）音高pitch and 曲调tuning

原理

1、f = 440Hz是中央A
2、一个八度在频率上翻了一倍
3、十二平均律中，每个半音在频率上翻了2^(1/12)倍
4、每个半音的距离就是黑白键上两个键之间的距离，从do到升do，各白键之间的差音为全全半全全全半，而每个全音中的半音被分到了黑键上
5、任意两个音之间的距离称为音程，音程的单位是度，相同单音之间音程为1度，之后每差一个音级增加1度

实验步骤：

1、设置基频f为440Hz，用numpy以440Hz生成正弦函数
2、每隔1/44100s设置一个采样点，每间隔1s生成1s的440Hz正弦波
3、在这基础上，第n次生成时再生成1s频率为440⋅2¹Hz的正弦波，与原正弦波线性相加进行和弦，从而模拟小二度到纯八度的所有音程

代码：

x=np.empty((0,))    #生成空数组
for i in range(0,13):
x1=np.linspace(0,1, num=44100, endpoint=True, dtype=float)    #生成采样点
x1=5*np.sin(2*np.pi*440*np.power(2,i/12)*x1)+5*np.sin(2*np.pi*440*x1)    #生成波形与和弦
x2=np.zeros((44100,))    #留空部分
x=np.concatenate((x,x1,x2),axis=0)    #连接数组

分析

小二度音程（相差一个半音）、大七度音程（相差五个全音一个半音）是极不协和音程
大二度音程（相差一个全音）、小七度音程（相差五个全音）、三全音（三个全音）是不协和音程
小三度音程（相差一个全音一个半音）、大三度音程（相差两个全音）、小六度音程（相差四个全音）、大六度音程（相差四个全音一个半音）是不完全协和音程
纯四度音程（相差两个全音一个半音）、纯五度音程（相差三个全音一个半音）是完全协和音程
纯八度音程（相差六个全音）是极完全协和音程

Dynamic Time Warping
1、dtw
2、 fill off diagonal

①estimate_tuning()

估计音频序列的音调或者频谱输入

②pitch_tuning()

给定一个集合,估计其调谐偏移(一个bin的分数)相对于A440 = 440.0Hz。

③piptrack()

阈值抛物线插值STFT上的节距跟踪。

将音乐旋律转化为音符：

import matplotlib.pyplot as plt
import librosa
import librosa.display
import librosa.util
import numpy as np
import pandas as pd
 
#要转换的音频文件
input="D:/My life/music/some music/Fabrizio Paterlini - My Misty Mornings.mp3"
 
y,sr=librosa.load(input,sr=None,duration=None)
chroma=librosa.feature.chroma_cqt(y=y, sr=sr) 
 
c=pd.DataFrame(chroma)
c0=(c==1)
c1=c0.astype(int)
labels=np.array(range(1,13))
note_values=labels.dot(c1)
 
plt.figure(figsize=(15,20))
 
plt.subplots_adjust(wspace=1, hspace=0.2)
 
plt.subplot(311)
librosa.display.specshow(chroma, y_axis='chroma', x_axis='time')
plt.xlabel('note')
plt.ylabel('beat')
note_values=labels.dot(c1)
 
plt.subplot(312)
librosa.display.waveplot(y, sr=sr)     #可以(,color)定义颜色
plt.xlabel('second')
plt.ylabel('amplitude')


plt.subplot(313)
plt.grid(linewidth=0.5)
plt.xticks(range(0,600,50))
plt.yticks(range(1,13),["C","C#","D","D#","E","F","F#","G","G#","A","A#","B"])
plt.scatter(range(len(note_values)),note_values,marker="s",s=1,color="red")
#plt.plot(note_values,color="#008080",linewidth=0.8)
#plt.hlines(note_values, range(len(note_values)),np.array(range(len(note_values)))+1 ,color="red", linewidth=10)
 
plt.show()

完美案例欣赏。想要提取出完美音阶太难了，音频稍微有一点噪音就会使结果不准。

（2）回音

D = np.repeat(D, 2, axis=1)

（3）间断

D[:,::2] = 0

（4）音色

D = np.roll(D, 50, axis=0)

1、对于周期性信号p(t) = ∑k=1~ K ( akcos(2πfkt) + bksin(2πfkt))，每个分解信号为谐波或泛音
2、f为频率
3、a与b为谐波系数
4、在K足够大时，任意周期性信号都可以通过积分变换为这种形式
5、谐波振幅的配比决定了音色，对于每个频率由k决定的谐波，其谐波振幅为ck =sqrt(a_k^2+b_k2)，因此谐波振幅可以组成一个长度为k的向量c = ( 0.3,0.4,…)，足够多的谐波以不同振幅数值混合可以组成不同音色
在220-11000Hz之间分别生成1、5、10、20、50个频段的正弦、余弦信号，随机生成谐波系数

分析

使用不同足够多的谐波以不同谐波系数混合可以组成不同音色音频
现实中已知乐器也可以利用频谱分析提取音色特征，供计算机再现
利用谐波原理可以制作各种各样的电子合成器，供作曲家使用

应用：

#使用 Librosa 对音色可视化
import librosa 
import matplotlib.pyplot as plt
import numpy as np
import librosa.display
import os

class specshowplot:
    def __init__(self):
        self.plotlist = []

def addplot(self, path, index):
    self.plotlist.insert(index, path)

def show(self):
    for index, path in  enumerate(self.plotlist):
        self._show(path, index)
    plt.show()

def _show(self, path, index):
    y,sr = librosa.load(path)
    count = len(self.plotlist)
    plt.subplot(count, 2, index*2+1)

    plt.title(os.path.basename(path))
    plt.subplots_adjust(left=None, bottom=None, right=None, top=None, wspace=0.4, hspace=1)
    librosa.display.waveplot(y)

    s = librosa.stft(y)
    s[abs(s)<s.max()/5] = 0
    y = librosa.istft(s)
    tone = librosa.tone(y,sr,len(y))
    d = librosa.feature.melspectrogram( y=tone )
    d = librosa.power_to_db(d,ref=np.max)
    s2 = plt.subplot(count, 2, index*2+2)
    librosa.display.specshow(d,x_axis='time', y_axis='mel')


if __name__ == '__main__':
outputpath = "D:/My life/music/some music/sweeter.mp3"
path = "D:/My life/music/some music/mi.mp3"
sp = specshowplot()
sp.addplot(path, 0)
sp.addplot(outputpath, 1)
sp.show()

#音色谱

import librosa
import librosa.display
import numpy as np
 
y,sr=librosa.load("D:/My life/music/some music/sweeter.mp3")
print(y.shape)
print(sr)
 
plt.plot(y)
 
#音色谱
chroma_stft = librosa.feature.chroma_stft(y=y, sr=sr,n_chroma=12, n_fft=4096)
#另一种常数Q音色谱
chroma_cq = librosa.feature.chroma_cqt(y=y, sr=sr)
#功率归一化音色谱
chroma_cens = librosa.feature.chroma_cens(y=y, sr=sr)
print(chroma_cens.shape)
 
plt.figure(figsize=(15,15))
plt.subplot(3,1,1)
librosa.display.specshow(chroma_stft, y_axis='chroma')
plt.title('chroma_stft')
plt.colorbar()
plt.subplot(3,1,2)
librosa.display.specshow(chroma_cq, y_axis='chroma', x_axis='time')
plt.title('chroma_cqt')
plt.colorbar()
plt.subplot(3,1,3)
librosa.display.specshow(chroma_cens, y_axis='chroma', x_axis='time')
plt.title('chroma_cens')
plt.colorbar()
plt.tight_layout()

（5）节拍

librosa.beat	用于估计节拍和检测节拍事件

应用：读取一段音频，判断节奏，并画出时频特性：

import  librosa
import  matplotlib.pyplot as plt
import  librosa.display
import numpy as np

y, sr  =  librosa.load('D:\My life\music\some music/sweeter.mp3')
tempo, beats  =  librosa.beat.beat_track(y = y, sr = sr)#tempo是节拍的个数，beat_frames 是节拍的帧数 

# Or use timing instead of frame indices
times  =  librosa.frames_to_time(beats, sr = sr)
y_beat_times  =  librosa.clicks(times = times, sr = sr)#beat_times 是将帧数转化为时间
#Or：y_beats  =  librosa.clicks(frames = beats, sr = sr)
#Or generate a signal of the same length as y： y_beats  =  librosa.clicks(frames = beats, sr = sr, length = len (y)) 
#Or with a click frequency of 880Hz and a 500ms sample：y_beat_times880  =  librosa.clicks(times = times, sr = sr, click_freq = 880 , click_duration = 0.5 )
 
# Display click waveform next to the spectrogram
plt.figure()
S  =  librosa.feature.melspectrogram(y = y, sr = sr)
ax  =  plt.subplot( 2 , 1 , 2 )
librosa.display.specshow(librosa.power_to_db(S, ref = np.max ),
                          x_axis = 'time' , y_axis = 'mel' )
plt.subplot( 2 , 1 , 1 , sharex = ax)
librosa.display.waveplot(y_beat_times, sr = sr, label = 'Beat clicks' )
plt.legend()
plt.xlim( 15 ,  30 )
plt.tight_layout()

卡点视频如何制作？
DCT：DCT变换实际上是做了Fourier的逆变换，这里是通过离散余弦变换来实现，做完DCT后会得到一些系数向量，我们称之为倒谱系数；

应用：为音频加节拍器

tempo, beats = librosa.beat.beat_track(onset_envelope = onset_envelope_y)
print(tempo)

注：接着要用到mir_eval，所以需要安装pip install mir_eval
官网链接为：https://pypi.org/project/mir_eval/0.4/

import mir_eval
beat_times = librosa.frames_to_time(beats)
y_click = mir_eval.sonify.clicks(beat_times, sr, length = len(y))
ipd.Audio(data=y+y_click, rate=sr)

可以得到加了节拍器的音频。试试运行一下听听吧！

---

1. (n−1)/12 ↩︎