python语音信号处理(一)「建议收藏」

大家好，又见面了，我是你们的朋友全栈君。

博客园地址：
https://www.cnblogs.com/LXP-Never/p/10078200.html

1.引言

python已经支持WAV格式的书写，而实时的声音输入输出需要安装pyAudio。最后我们还将使用pyMedia进行Mp3的解码和播放。

音频信号是模拟信号，我们需要将其保存为数字信号，才能对语音进行算法操作，WAV是Microsoft开发的一种声音文件格式，通常被用来保存未压缩的声音数据。

语音信号有四个重要的参数：声道数、采样频率、量化位数(位深)和比特率
声道数：可以是单声道、双声道 …

• 采样频率（Sample rate）：每秒内对声音信号采样样本的总数目，44100Hz采样频率意味着每秒钟信号被分解成44100份。换句话说，每隔1/144100秒就会存储一次，如果采样率高，那么媒体播放音频时会感觉信号是连续的。
• 量化位数（Bit depth）：也称为“位深”，每个采样点中信息的比特(bit)数。1 byte等于8 bit。通常有8bit、16bit、24bit、32bit…
• 比特率（Bit rate）: 每秒处理多少个Bit。比如一个单声道，用44.1KHz/16Bit的配置来说，它的比特率就为44100161=705600，单位是bit/s(或者bps)，因为通常计算出来的数字都比较大，大家就用kbit/s了，也就是705.6kbit/s。在对音频进行压缩时，比特率就成为了我们的一个要选的选项了，越高的比特率，其音质也就越好。一些常用的比特率有：

       - 32kbit/s： 一般只适用于语音
       - 96kbit/s： 一般用于语音或低质量流媒体
	   - 128或160kbit/s： 中等比特率质量
	   -  192kbit/s： 中等质量比特率
	   -  256kbit/s： 常用的高质量比特率
	   -  320kbit/s： MP3标准支持的最高水平

如果你需要自己录制和编辑声音文件，推荐使用Audacity，它是一款开源的、跨平台、多声道的录音编辑软件。在我的工作中经常使用Audacity进行声音信号的录制，然后再输出成WAV文件供Python程序处理。
3
　　如果想要快速看语音波形和语谱图，推荐使用Adobe Audition，他是Adobe公司开发专门处理音频的专业软件，微博关注vposy，下载地址见置顶。他破解了很多adobe公司的软件，包括PS、PR…

2. 读取音频文件

• wave库

wave库是python的标准库，对于python来说相对底层，wave不支持压缩/解压，但支持单声道/立体声语音的读取。

wave_read = wave.open(file,mode="rb")

参数：

f：语音文件名或文件路径
mode：读或写
"rb"：只读模式
"wb"：只写模式
返回：读取的文件流

该open()函数可用于with声明中。当with块完成时，wave_read.close()或wave_write.close()方法被调用

文件路径：
例如voice.wav文件在路径C:\Users\Never\Desktop\code for the speech的文件夹里
　　则file有以下三种填写格式：
　　　　r”C:\Users\Never\Desktop\code for the speech\voice.wav”
　　　　“C:/Users/Never/Desktop/code for the speech/voice.wav”
　　　　“C:\Users\Never\Desktop\code for the speech\voice.wav”
　　三者等价，右划线\为转意字符，如果要表达\则需要\，引号前面加r表示原始字符串。

wave_read.getparams()：一次性返回所有的音频参数，返回的是一个元组(声道数，量化位数(byte单位)，采样频率，采样点数，压缩类型，压缩类型的描述)。(nchannels, sampwidth, framerate, nframes, comptype, compname)wave模块只支持非压缩的数据，因此可以忽略最后两个信息。

str_data = wave_read.readframes(nframes)：读取的长度(以取样点为单位)，返回的是字符串类型的数据

wave_data = np.fromstring(str_data, dtype=np.float16)：将上面字符串类型数据转换为一维float16类型的数组。

现在的wave_data是一个一维的short类型的数组，但是因为我们的声音文件是双声道的，因此它由左右两个声道的取样交替构成：LR

wave_data.shape = (-1, 2)　　# -1的意思就是没有指定,根据另一个维度的数量进行分割，得到n行2列的数组。

wave_read.close（）　　关闭文件流wave
wave_read.getnchannels（）　　返回音频通道的数量（1对于单声道，2对于立体声）。
wave_read.getsampwidth（）　　以字节为单位返回样本宽度
wave_read.getframerate（）　　返回采样频率。
wave_read.getnframes（）　　　返回音频帧数。
wave_read.rewind（）　　　　　　将文件指针倒回到音频流的开头。
wave_read.tell（）　　　　　　返回当前文件指针位置。 

读Wave文件并且绘制波形

# -*- coding: utf-8 -*-
# 读Wave文件并且绘制波形
import wave
import matplotlib.pyplot as plt
import numpy as np

plt.rcParams['font.sans-serif'] = ['SimHei']  # 用来正常显示中文标签
plt.rcParams['axes.unicode_minus'] = False  # 用来正常显示符号

# 打开WAV音频
f = wave.open(r"C:\Windows\media\Windows Background.wav", "rb")

# 读取格式信息
# (声道数、量化位数、采样频率、采样点数、压缩类型、压缩类型的描述)
# (nchannels, sampwidth, framerate, nframes, comptype, compname)
params = f.getparams()
nchannels, sampwidth, framerate, nframes = params[:4]
# nchannels通道数 = 2
# sampwidth量化位数 = 2
# framerate采样频率 = 22050
# nframes采样点数 = 53395

# 读取nframes个数据，返回字符串格式
str_data = f.readframes(nframes)

f.close()

# 将字符串转换为数组，得到一维的short类型的数组
wave_data = np.fromstring(str_data, dtype=np.short)

# 赋值的归一化
wave_data = wave_data * 1.0 / (max(abs(wave_data)))

# 整合左声道和右声道的数据
wave_data = np.reshape(wave_data, [nframes, nchannels])
# wave_data.shape = (-1, 2) # -1的意思就是没有指定,根据另一个维度的数量进行分割

# 最后通过采样点数和取样频率计算出每个取样的时间
time = np.arange(0, nframes) * (1.0 / framerate)

plt.figure()
# 左声道波形
plt.subplot(2, 1, 1)
plt.plot(time, wave_data[:, 0])
plt.xlabel("时间/s",fontsize=14)
plt.ylabel("幅度",fontsize=14)
plt.title("左声道",fontsize=14)
plt.grid()  # 标尺

plt.subplot(2, 1, 2)
# 右声道波形
plt.plot(time, wave_data[:, 1], c="g")
plt.xlabel("时间/s",fontsize=14)
plt.ylabel("幅度",fontsize=14)
plt.title("右声道",fontsize=14)

plt.tight_layout()  # 紧密布局
plt.show()

读取通道数为2的音频信号

如果出现错误：np.fromstring string size must be a multiple of element size
参考解决方法：https://blog.csdn.net/veritasalice/article/details/104807415

• librosa库 (推荐)

这是我最常用也是最喜欢的语音库，librosa是python第三方库，我们在使用前需要在cmd终端运行： pip install librosa ，关于librosa的介绍我专门写了一篇博客librosa语音信号处理。

import librosa
y, sr = librosa.load(path, sr=fs)

该函数是会改变声音的采样频率的。如果 sr 缺省，librosa.load()会默认以22050的采样率读取音频文件，高于该采样率的音频文件会被下采样，低于该采样率的文件会被上采样。因此，如果希望以原始采样率读取音频文件，sr 应当设为 None。具体做法为 y, sr = librosa(filename, sr=None)。
音频数据 y 是直接经过归一化的数组

• scipy库

from scipy.io import wavfile
sampling_freq, audio = wavfile.read("***.wav")

audio 是直接经过归一化的数组

3.写音频文件

• wave库

在写入第一帧数据时，先通过调用setnframes()设置好帧数，setnchannels()设置好声道数，setsampwidth()设置量化位数，setframerate()设置好采样频率，然后writeframes(wave.tostring())用于写入帧数据。

wave_write = wave.open(file,mode="wb") 

wave_write是写文件流

wave_write.setnchannels（n）　　设置通道数。
wave_write.setsampwidth（n）　　将样本宽度设置为n个字节，量化位数
wave_write.setframerate（n）　　将采样频率设置为n。
wave_write.setnframes（n）　　将帧数设置为n
wave_write.setparams（tuple）　　以元组形式设置所有参数(nchannels, sampwidth, framerate, nframes,comptype, compname)
wave_write.writeframes（data）　　写入data个长度的音频，以采样点为单位
wave_write.tell（）　　返回文件中的当前位置

• 写wav文件
  写文件方法1

# Author:凌逆战
# -*- coding:utf-8 -*-
import wave
import numpy as np
import scipy.signal as signal

framerate = 44100   # 采样频率
time = 10           # 持续时间

t = np.arange(0, time, 1.0/framerate)

# 调用scipy.signal库中的chrip函数，
# 产生长度为10秒、取样频率为44.1kHz、100Hz到1kHz的频率扫描波
wave_data = signal.chirp(t, 100, time, 1000, method='linear') * 10000

# 由于chrip函数返回的数组为float64型，
# 需要调用数组的astype方法将其转换为short型。
wave_data = wave_data.astype(np.short)

# 打开WAV音频用来写操作
f = wave.open(r"sweep.wav", "wb")

f.setnchannels(1)           # 配置声道数
f.setsampwidth(2)           # 配置量化位数
f.setframerate(framerate)   # 配置取样频率
comptype = "NONE"
compname = "not compressed"

# 也可以用setparams一次性配置所有参数
# outwave.setparams((1, 2, framerate, nframes,comptype, compname))

# 将wav_data转换为二进制数据写入文件
f.writeframes(wave_data.tobytes())
f.close()

写WAV文件方法2

# Author:凌逆战
# -*- coding:utf-8 -*-
import wave
import numpy as np
import struct

f = wave.open(r"C:\Windows\media\Windows Background.wav", "rb")
params = f.getparams()
nchannels, sampwidth, framerate, nframes = params[:4]
strData = f.readframes(nframes)
waveData = np.fromstring(strData,dtype=np.int16)
f.close()
waveData = waveData*1.0/(max(abs(waveData)))

# wav文件写入
# 待写入wav的数据，这里仍然取waveData数据
outData = waveData
outwave = wave.open("write.wav", 'wb')
nchannels = 1   # 通道数设置为1
sampwidth = 2   # 量化位数设置为2
framerate = 8000    # 采样频率8000
nframes = len(outData)    # 采样点数

comptype = "NONE"
compname = "not compressed"
outwave.setparams((nchannels, sampwidth, framerate, nframes,
    comptype, compname))

for i in outData:
        outwave.writeframes(struct.pack('h', int(i * 64000 / 2)))

        # struct.pack(FMT, V1)将V1的值转换为FMT格式字符串
outwave.close()

• librosa库 (推荐)

librosa.output.write_wav(path, y, sr, norm=False)

参数：

path：str，保存输出wav文件的路径
y：np.ndarry 音频时间序列
sr：y的采样率
norm：True/False，是否启动幅值归一化

在0.8.0以后的版本，librosa都会将这个函数删除，推荐用下面的函数：

import soundfile
soundfile.write(file, data, samplerate)

参数：

file：保存输出wav文件的路径
data：音频数据
samplerate：采样率

• scipy库

from scipy.io.wavfile import write
write(output_filename, freq, audio)

import numpy as np
import matplotlib.pyplot as plt
from scipy.io.wavfile import write

# 定义存储音频的输出文件
output_file = 'output_generated.wav'

# 指定音频生成的参数
duration = 3            # 单位秒
sampling_freq = 44100   # 单位Hz
tone_freq = 587         # 音调的频率
min_val = -2 * np.pi
max_val = 2 * np.pi

# 生成音频信号
t = np.linspace(min_val, max_val, duration * sampling_freq)
audio = np.sin(2 * np.pi * tone_freq * t)

# 添加噪声(duration * sampling_freq个(0,1]之间的随机值)
noise = 0.4 * np.random.rand(duration * sampling_freq)
audio += noise

scaling_factor = pow(2,15) - 1  # 转换为16位整型数
audio_normalized = audio / np.max(np.abs(audio))    # 归一化
audio_scaled = np.int16(audio_normalized * scaling_factor)  # 这句话什么意思

write(output_file, sampling_freq, audio_scaled) # 写入输出文件

audio = audio[:300] # 取前300个音频信号

x_values = np.arange(0, len(audio), 1) / float(sampling_freq)
x_values *= 1000    # 将时间轴单位转换为秒

plt.plot(x_values, audio, color='blue')
plt.xlabel('Time (ms)')
plt.ylabel('Amplitude')
plt.title('Audio signal')
plt.show()

• 合成有音调的音乐

import json
import numpy as np
from scipy.io.wavfile import write
import matplotlib.pyplot as plt          

# 定义合成音调
def Synthetic_tone(freq, duration, amp=1.0, sampling_freq=44100):
    # 建立时间轴
    t = np.linspace(0, duration, int(duration * sampling_freq))
    # 构建音频信号
    audio = amp * np.sin(2 * np.pi * freq * t)
    return audio.astype(np.int16)


# json文件中包含一些音阶以及他们的频率
tone_map_file = 'tone_freq_map.json'

# 读取频率映射文件
with open(tone_map_file, 'r',encoding='UTF-8') as f:
    tone_freq_map = json.loads(f.read())
    print(tone_freq_map)
# {'A': 440, 'Asharp': 466, 'B': 494, 'C': 523, 'Csharp': 554, 'D': 587, 'Dsharp': 622, 'E': 659, 'F': 698, 'Fsharp': 740, 'G': 784, 'Gsharp': 831}

# 设置生成G调的输入参数
input_tone = 'G'
duration = 2             # seconds
amplitude = 10000        # 振幅
sampling_freq = 44100    # Hz
# 生成音阶
synthesized_tone = Synthetic_tone(tone_freq_map[input_tone], duration, amplitude, sampling_freq)

# 写入输出文件
write('output_tone.wav', sampling_freq, synthesized_tone)

# 音阶及其连续时间
tone_seq = [('D', 0.3), ('G', 0.6), ('C', 0.5), ('A', 0.3), ('Asharp', 0.7)]

# 构建基于和弦序列的音频信号
output = np.array([])
for item in tone_seq:
    input_tone = item[0]
    duration = item[1]
    synthesized_tone = Synthetic_tone(tone_freq_map[input_tone], duration, amplitude, sampling_freq)
    output = np.append(output, synthesized_tone, axis=0)

# 写入输出文件
write('output_tone_seq.wav', sampling_freq, output)

tone_freq_map.json

{
    "A": 440,
    "Asharp": 466,
    "B": 494,
    "C": 523,
    "Csharp": 554,
    "D": 587,
    "Dsharp": 622,
    "E": 659,
    "F": 698,
    "Fsharp": 740,
    "G": 784,
    "Gsharp": 831
}

– 音频播放
wav文件的播放用到的是pyaudio库
p = pyaudio.PyAudio() stream = p.open(format = p.get_format_from_width(sampwidth) , channels ,rate ,output = True) stream.write(data)　　# 播放data数据

以下列出pyaudio对象的open()方法的主要参数：

rate：取样频率
channels：声道数
format：取样值的量化格式 (paFloat32, paInt32, paInt24, paInt16, paInt8 …)。在上面的例子中，使用get_format_from_width方法将wf.sampwidth()的返回值2转换为paInt16
input：输入流标志，如果为True的话则开启输入流
output：输出流标志，如果为True的话则开启输出流
input_device_index：输入流所使用的设备的编号，如果不指定的话，则使用系统的缺省设备
output_device_index：输出流所使用的设备的编号，如果不指定的话，则使用系统的缺省设备
frames_per_buffer：底层的缓存的块的大小，底层的缓存由N个同样大小的块组成
start：指定是否立即开启输入输出流，缺省值为True

播放wav音频

# -*- coding: utf-8 -*-
import pyaudio
import wave

chunk = 1024
wf = wave.open(r"c:\WINDOWS\Media\Windows Background.wav", 'rb')
p = pyaudio.PyAudio()

# 打开声音输出流
stream = p.open(format = p.get_format_from_width(wf.getsampwidth()),
                channels = wf.getnchannels(),
                rate = wf.getframerate(),
                output = True)

# 写声音输出流到声卡进行播放
while True:
    data = wf.readframes(chunk)
    if data == "":
        break
    stream.write(data)

stream.stop_stream()
stream.close()
p.terminate()   # 关闭PyAudio

录音

以SAMPLING_RATE为采样频率，每次读入一块有NUM_SAMPLES个采样的数据块，当读入的采样数据中有COUNT_NUM个值大于LEVEL的取样的时候，将数据保存进WAV文件，一旦开始保存数据，所保存的数据长度最短为SAVE_LENGTH个块。WAV文件以保存时的时刻作为文件名。

从声卡读入的数据和从WAV文件读入的类似，都是二进制数据，由于我们用paInt16格式(16bit的short类型)保存采样值，因此将它自己转换为dtype为np.short的数组。

''' 以SAMPLING_RATE为采样频率， 每次读入一块有NUM_SAMPLES个采样点的数据块， 当读入的采样数据中有COUNT_NUM个值大于LEVEL的取样的时候， 将采样数据保存进WAV文件， 一旦开始保存数据，所保存的数据长度最短为SAVE_LENGTH个数据块。 从声卡读入的数据和从WAV文件读入的类似，都是二进制数据， 由于我们用paInt16格式(16bit的short类型)保存采样值， 因此将它自己转换为dtype为np.short的数组。 '''


from pyaudio import PyAudio, paInt16
import numpy as np
import wave

# 将data中的数据保存到名为filename的WAV文件中
def save_wave_file(filename, data):
    wf = wave.open(filename, 'wb')
    wf.setnchannels(1)          # 单通道
    wf.setsampwidth(2)          # 量化位数
    wf.setframerate(SAMPLING_RATE)  # 设置采样频率
    wf.writeframes(b"".join(data))  # 写入语音帧
    wf.close()


NUM_SAMPLES = 2000      # pyAudio内部缓存块的大小
SAMPLING_RATE = 8000    # 取样频率
LEVEL = 1500           # 声音保存的阈值，小于这个阈值不录
COUNT_NUM = 20 # 缓存快类如果有20个大于阈值的取样则记录声音
SAVE_LENGTH = 8 # 声音记录的最小长度：SAVE_LENGTH * NUM_SAMPLES 个取样

# 开启声音输入
pa = PyAudio()
stream = pa.open(format=paInt16, channels=1, rate=SAMPLING_RATE, input=True,
                frames_per_buffer=NUM_SAMPLES)

save_count = 0  # 用来计数
save_buffer = []    #

while True:
    # 读入NUM_SAMPLES个取样
    string_audio_data = stream.read(NUM_SAMPLES)
    # 将读入的数据转换为数组
    audio_data = np.fromstring(string_audio_data, dtype=np.short)
    # 计算大于LEVEL的取样的个数
    large_sample_count = np.sum( audio_data > LEVEL )
    print(np.max(audio_data))
    # 如果个数大于COUNT_NUM，则至少保存SAVE_LENGTH个块
    if large_sample_count > COUNT_NUM:
        save_count = SAVE_LENGTH
    else:
        save_count -= 1

    if save_count < 0:
        save_count = 0

    if save_count > 0:
        # 将要保存的数据存放到save_buffer中
        save_buffer.append( string_audio_data )
    else:
        # 将save_buffer中的数据写入WAV文件，WAV文件的文件名是保存的时刻
        if len(save_buffer) > 0:
            filename = "recorde" + ".wav"
            save_wave_file(filename, save_buffer)
            print(filename, "saved")
            break

3.语音信号处理

3.1语音信号的产生和感知

我们要对语音进行分析，首先要提取能够表示该语音的特征参数，有了特征参数才可能利用这些参数进行有效的处理，在对语音信号处理的过程中，语音信号的质量不仅取决于处理方法，同时取决于时候选对了合适的特征参数。

语音信号是一个非平稳的时变信号，但语音信号是由声门的激励脉冲通过声道形成的，而声道(人的口腔、鼻腔)的肌肉运动是缓慢的，所以“短时间”(10~30ms)内可以认为语音信号是平稳时不变的。由此构成了语音信号的“短时分析技术”。

在短时分析中，将语音信号分为一段一段的语音帧，每一帧一般取10~30ms，我们的研究就建立在每一帧的语音特征分析上。

提取的不同的语音特征参数对应着不同的语音信号分析方法：时域分析、频域分析、倒谱域分析...由于语音信号最重要的感知特性反映在功率谱上，而相位变化只起到很小的作用，所有语音频域分析更加重要。

3.2 信号加窗

通常对信号截断、分帧需要加窗，因为截断都有频域能量泄露，而窗函数可以减少截断带来的影响。
窗函数在scipy.signal信号处理工具箱中，如hanning窗：

import matplotlib.pyplot as plt
import scipy.signal as signal
plt.figure(figsize=(6,2))
plt.plot(signal.hanning(512))
plt.show()

3.3 信号分帧

在分帧中，相邻两帧之间会有一部分重叠，帧长($wlen$) = 重叠($overlap$)+帧移($inc$)，如果相邻两帧之间不重叠，那么由于窗函数的形状，截取到的语音帧边缘会出现损失，所以要设置重叠部分。$inc$为帧移，表示后一帧对前一帧的偏移量，$f_s$表示采样率，$f_n$表示一段语音信号的分帧数。

$$\frac{N-overlap}{inc}=\frac{N-wlen+inc}{inc}$$
信号分帧的理论依据，其中$x$是语音信号，$w$是窗函数：
$$y(n)=\sum _{-n=\frac{N}{2}+1}^{\frac{N}{2}}x(m)w(n-m)$$

加窗截断类似采样，为了保证相邻帧不至于差别过大，通常帧与帧之间有帧移，其实就是插值平滑的作用。

给出示意图：

这里主要用到numpy工具包，涉及的指令有：

np.repeat：主要是直接重复
np.tile：主要是周期性重复

对比一下：
向量情况：

矩阵情况：
对于数据

repeat操作:

tile操作：

对应结果：

对应分帧的代码实现：
这是没有加窗的示例：

• 没有加窗的语音分帧

import numpy as np
import wave
import os
#import math

def enframe(signal, nw, inc):
   '''将音频信号转化为帧。 参数含义： signal:原始音频型号 nw:每一帧的长度(这里指采样点的长度，即采样频率乘以时间间隔) inc:相邻帧的间隔（同上定义） '''
   signal_length=len(signal) #信号总长度
   if signal_length<=nw: #若信号长度小于一个帧的长度，则帧数定义为1
       nf=1
   else: #否则，计算帧的总长度
       nf=int(np.ceil((1.0*signal_length-nw+inc)/inc))
   pad_length=int((nf-1)*inc+nw) #所有帧加起来总的铺平后的长度
   zeros=np.zeros((pad_length-signal_length,)) #不够的长度使用0填补，类似于FFT中的扩充数组操作
   pad_signal=np.concatenate((signal,zeros)) #填补后的信号记为pad_signal
   indices=np.tile(np.arange(0,nw),(nf,1))+np.tile(np.arange(0,nf*inc,inc),(nw,1)).T  #相当于对所有帧的时间点进行抽取，得到nf*nw长度的矩阵
   indices=np.array(indices,dtype=np.int32) #将indices转化为矩阵
   frames=pad_signal[indices] #得到帧信号
# win=np.tile(winfunc(nw),(nf,1)) #window窗函数，这里默认取1
# return frames*win #返回帧信号矩阵
   return frames
def wavread(filename):
   f = wave.open(filename,'rb')
   params = f.getparams()
   nchannels, sampwidth, framerate, nframes = params[:4]
   strData = f.readframes(nframes)#读取音频，字符串格式
   waveData = np.fromstring(strData,dtype=np.int16)#将字符串转化为int
   f.close()
   waveData = waveData*1.0/(max(abs(waveData)))#wave幅值归一化
   waveData = np.reshape(waveData,[nframes,nchannels]).T
   return waveData

filepath = "./data/" #添加路径
dirname= os.listdir(filepath) #得到文件夹下的所有文件名称 
filename = filepath+dirname[0]
data = wavread(filename)
nw = 512
inc = 128
Frame = enframe(data[0], nw, inc) 

• 加窗的语音分帧

def enframe(signal, nw, inc, winfunc):
    '''将音频信号转化为帧。 参数含义： signal:原始音频型号 nw:每一帧的长度(这里指采样点的长度，即采样频率乘以时间间隔) inc:相邻帧的间隔（同上定义） '''
    signal_length=len(signal) #信号总长度
    if signal_length<=nw: #若信号长度小于一个帧的长度，则帧数定义为1
        nf=1
    else: #否则，计算帧的总长度
        nf=int(np.ceil((1.0*signal_length-nw+inc)/inc))
    pad_length=int((nf-1)*inc+nw) #所有帧加起来总的铺平后的长度
    zeros=np.zeros((pad_length-signal_length,)) #不够的长度使用0填补，类似于FFT中的扩充数组操作
    pad_signal=np.concatenate((signal,zeros)) #填补后的信号记为pad_signal
    indices=np.tile(np.arange(0,nw),(nf,1))+np.tile(np.arange(0,nf*inc,inc),(nw,1)).T  #相当于对所有帧的时间点进行抽取，得到nf*nw长度的矩阵
    indices=np.array(indices,dtype=np.int32) #将indices转化为矩阵
    frames=pad_signal[indices] #得到帧信号
    win=np.tile(winfunc,(nf,1))  #window窗函数，这里默认取1
    return frames*win   #返回帧信号矩阵

• overlap and add
  　　将分帧好的语音拼接回完整的语音

def overlap_add(x, window_size, hop_size):
    # x (frames, frame_length)
    frames, frame_length = x.shape
    wav_len = frames * hop_size + window_size - hop_size  # 帧长
    print("帧长", wav_len)
    wav = np.zeros((wav_len,), dtype=x.dtype)

    for frame in range(frames):
        if frame == frames - 1:
            # 最后一帧
            wav[hop_size * frame: hop_size * frame + window_size] += x[frame]  # 把整帧前加进去
        else:
            wav[hop_size * frame: hop_size * frame + hop_size] += x[frame][:hop_size]  # 把帧前加进去

    return wav

3.4 语音信号的短时时域处理

短时能量和短时平均幅度

短时能量和短时平均幅度的主要用途：

• 区分浊音和清音段，因为浊音的短时能量E(i)E(i)比清音大很多；
• 区分声母和韵母的分界和无话段和有话段的分界
  短时平均过零率

对于连续语音信号，过零率意味着时域波形通过时间轴，对于离散信号，如果相邻的取样值改变符号，则称为过零。

作用：

• 发浊音时由于声门波引起谱的高频跌落，所以语音信号能量约集中在3kHz以下
• 发清音时多数能量集中在较高的频率上，

因为高频意味着高的短时平均过零率，低频意味着低的短时平均过零率，所以浊音时具有较低的过零率，而清音时具有较高的过零率

1.利用短时平均过零率可以从背景噪声中找出语音信号，
2.可以用于判断寂静无话段与有话段的起点和终止位置。
3.在背景噪声较小的时候，用平均能量识别较为有效，在背景噪声较大的时候，用短时平均过零率识别较为有效。

短时自相关函数

短时自相关函数主要应用于端点检测和基音的提取，在韵母基因频率整数倍处将出现峰值特性，通常根据除R(0)外的第一峰值来估计基音，而在声母的短时自相关函数中看不到明显的峰值。

短时平均幅度差函数

用于检测基音周期，而且在计算上比短时自相关函数更加简单。

语音信号的短时频域处理
　　在语音信号处理中，信号在频域或其他变换域上的分析处理占重要的位置，在频域上研究语音可以使信号在时域上无法表现出来的某些特征变得十分明显，一个音频信号的本质是由其频率内容决定的，

将时域信号转换为频域信号一般对语音进行短时傅里叶变换。

fft_audio = np.fft.fft(audio)

绘制语音信号的频谱图

import numpy as np
from scipy.io import wavfile
import matplotlib.pyplot as plt

sampling_freq, audio = wavfile.read(r"C:\Windows\media\Windows Background.wav")   # 读取文件

audio = audio / np.max(audio)   # 归一化，标准化

# 应用傅里叶变换
fft_signal = np.fft.fft(audio)
print(fft_signal)
# [-0.04022912+0.j -0.04068997-0.00052721j -0.03933007-0.00448355j
# ... -0.03947908+0.00298096j -0.03933007+0.00448355j -0.04068997+0.00052721j]

fft_signal = abs(fft_signal)
print(fft_signal)
# [0.04022912 0.04069339 0.0395848 ... 0.08001755 0.09203427 0.12889393]

# 建立时间轴
Freq = np.arange(0, len(fft_signal))

# 绘制语音信号的
plt.figure()
plt.plot(Freq, fft_signal, color='blue')
plt.xlabel('Freq (in kHz)')
plt.ylabel('Amplitude')
plt.show()

• 提取频域特征
  将信号转换为频域之后，还需要将其转换为有用的形式，梅尔频率倒谱系数(MFCC)，MFCC首先计算信号的功率谱，然后用滤波器组和离散余弦变换的组合来提取特征。

提取MFCC特征

import numpy as np
import matplotlib.pyplot as plt
from scipy.io import wavfile
from python_speech_features import mfcc, logfbank

# 读取输入音频文件
sampling_freq, audio = wavfile.read("input_freq.wav")

# 提取MFCC和滤波器组特征
mfcc_features = mfcc(audio, sampling_freq)
filterbank_features = logfbank(audio, sampling_freq)

print('\nMFCC:\n窗口数 =', mfcc_features.shape[0])
print('每个特征的长度 =', mfcc_features.shape[1])
print('\nFilter bank:\n窗口数 =', filterbank_features.shape[0])
print('每个特征的长度 =', filterbank_features.shape[1])

# 画出特征图，将MFCC可视化。转置矩阵，使得时域是水平的
mfcc_features = mfcc_features.T
plt.matshow(mfcc_features)
plt.title('MFCC')
# 将滤波器组特征可视化。转置矩阵，使得时域是水平的
filterbank_features = filterbank_features.T
plt.matshow(filterbank_features)
plt.title('Filter bank')

plt.show()

语谱图

绝大部分信号都可以分解为若干不同频率的正弦波。
这些正弦波中，频率最低的称为信号的基波，其余称为信号的谐波。
基波只有一个，可以称为一次谐波，谐波可以有很多个，每次谐波的频率是基波频率的整数倍。谐波的大小可能互不相同。
以谐波的频率为横坐标，幅值（大小）为纵坐标，绘制的系列条形图，称为频谱。频谱能够准确反映信号的内部构造。

语谱图综合了时域和频域的特点，明显的显示出来了语音频率随时间的变化情况，语谱图的横轴为时间，纵轴为频率任意给定频率成分在给定时刻的强弱用颜色深浅表示。颜色深表示频谱值大，颜色浅表示频谱值小，语谱图上不同的黑白程度形成不同的纹路，称为声纹，不用讲话者的声纹是不一样的，可以用做声纹识别。
其实得到了分帧信号，频域变换取幅值，就可以得到语谱图，如果仅仅是观察，matplotlib.pyplot有specgram指令：

import wave
import matplotlib.pyplot as plt
import numpy as np

f = wave.open(r"C:\Windows\media\Windows Background.wav", "rb")
params = f.getparams()
nchannels, sampwidth, framerate, nframes = params[:4]
strData = f.readframes(nframes)#读取音频，字符串格式
waveData = np.fromstring(strData,dtype=np.int16)#将字符串转化为int
waveData = waveData*1.0/(max(abs(waveData)))#wave幅值归一化
waveData = np.reshape(waveData,[nframes,nchannels]).T
f.close()

plt.specgram(waveData[0],Fs = framerate, scale_by_freq = True, sides = 'default')
plt.ylabel('Frequency(Hz)')
plt.xlabel('Time(s)')
plt.show()

matlab语谱图

[Y,FS]=audioread('p225_355_wb.wav');

% specgram(Y,2048,44100,2048,1536);
    %Y1为波形数据
    %FFT帧长2048点(在44100Hz频率时约为46ms)
    %采样频率44.1KHz
    %加窗长度，一般与帧长相等
    %帧重叠长度，此处取为帧长的3/4
specgram(Y,2048,FS,2048,1536);
xlabel('时间(s)')
ylabel('频率(Hz)')
title('语谱图')

3.5 语音识别

import os
import numpy as np
import scipy.io.wavfile as wf
import python_speech_features as sf
import hmmlearn.hmm as hl


# 1. 读取training文件夹中的训练音频样本，每个音频对应一个mfcc矩阵，每个mfcc都有一个类别(apple...)
def search_file(directory):
    """ :param directory: 训练音频的路径 :return: 字典{'apple':[url, url, url ... ], 'banana':[...]} """
    # 使传过来的directory匹配当前操作系统
    directory = os.path.normpath(directory)
    objects = { 
   }
    # curdir：当前目录
    # subdirs: 当前目录下的所有子目录
    # files: 当前目录下的所有文件名
    for curdir, subdirs, files in os.walk(directory):
        for file in files:
            if file.endswith('.wav'):
                label = curdir.split(os.path.sep)[-1]  # os.path.sep为路径分隔符
                if label not in objects:
                    objects[label] = []
                # 把路径添加到label对应的列表中
                path = os.path.join(curdir, file)
                objects[label].append(path)
    return objects


# 读取训练集数据
train_samples = search_file('../machine_learning_date/speeches/training')

""" 2. 把所有类别为apple的mfcc合并在一起，形成训练集。 训练集: train_x：[mfcc1,mfcc2,mfcc3,...],[mfcc1,mfcc2,mfcc3,...]... train_y：[apple],[banana]... 由上述训练集样本可以训练一个用于匹配apple的HMM。"""

train_x, train_y = [], []
# 遍历字典
for label, filenames in train_samples.items():
    # [('apple', ['url1,,url2...'])
    # [("banana"),("url1,url2,url3...")]...
    mfccs = np.array([])
    for filename in filenames:
        sample_rate, sigs = wf.read(filename)
        mfcc = sf.mfcc(sigs, sample_rate)
        if len(mfccs) == 0:
            mfccs = mfcc
        else:
            mfccs = np.append(mfccs, mfcc, axis=0)
    train_x.append(mfccs)
    train_y.append(label)

# 3.训练模型，有7个句子，创建了7个模型
models = { 
   }
for mfccs, label in zip(train_x, train_y):
    model = hl.GaussianHMM(n_components=4, covariance_type='diag', n_iter=1000)
    models[label] = model.fit(mfccs)  # # {'apple':object, 'banana':object ...}

""" 4. 读取testing文件夹中的测试样本， 测试集数据： test_x [mfcc1, mfcc2, mfcc3...] test_y [apple, banana, lime] """
test_samples = search_file('../machine_learning_date/speeches/testing')

test_x, test_y = [], []
for label, filenames in test_samples.items():
    mfccs = np.array([])
    for filename in filenames:
        sample_rate, sigs = wf.read(filename)
        mfcc = sf.mfcc(sigs, sample_rate)
        if len(mfccs) == 0:
            mfccs = mfcc
        else:
            mfccs = np.append(mfccs, mfcc, axis=0)
    test_x.append(mfccs)
    test_y.append(label)


# 5.测试模型
# 1. 分别使用7个HMM模型，对测试样本计算score得分。
# 2. 取7个模型中得分最高的模型所属类别作为预测类别。
pred_test_y = []
for mfccs in test_x:
    # 判断mfccs与哪一个HMM模型更加匹配
    best_score, best_label = None, None
    # 遍历7个模型
    for label, model in models.items():
        score = model.score(mfccs)
        if (best_score is None) or (best_score < score):
            best_score = score
            best_label = label
    pred_test_y.append(best_label)

print(test_y)   # ['apple', 'banana', 'kiwi', 'lime', 'orange', 'peach', 'pineapple']
print(pred_test_y)  # ['apple', 'banana', 'kiwi', 'lime', 'orange', 'peach', 'pineapple']

我对上面这段代码专门写了一篇博客来进一步讲解和分析，想详细了解的读者可以移步https://www.cnblogs.com/LXP-Never/p/11415110.html，语音数据集在这里。

4 参考文献

网址：用python做科学计算 http://old.sebug.net/paper/books/scipydoc/index.html#

python标准库wave模块https://docs.python.org/3.6/library/wave.html

《python机器学习经典案例》美Prateek Joshi著

傅里叶变换的介绍：http://www.thefouriertransform.com/

各种音阶及其对应的频率 http://pages.mtu.edu/~suits/notefreqs.html

这篇博客的代码https://github.com/LXP-Neve/Speech-signal-processing

这个网站有很多numpy写的语音信号处理代码

【知乎文章】采样率，位深以及比特率