【自学记录】PyTorch语音识别实战

第三章：音频信号处理的理论与Python实战

3.4 梅尔频率倒谱系数的计算过程

第一步：数据准备

import numpy as np
from torchaudio import datasets
import sound_utils
import soundfile as sf
from matplotlib import pyplot as plt

"--------------第一步：数据准备-------------"
#signal：NumPy 数组，包含了音频文件中的音频数据。每个元素代表音频信号中的一个样本
#sample_rate：整数，表示音频文件的采样率（以赫兹为单位）。采样率是指每秒从音频信号中采样的次数。例如，如果采样率是 44100 Hz，那么每秒就有 44100 个样本被读取。
#signal_len=176400,sample_rate=44100,4秒的音频，采样率为44100HZ，采样点=176400=4*44100

signal , sample_rate = sf.read("./carsound.wav", dtype="float32")
signal_num = np.arange(len(signal))#生成[0,1,2,3...,len(signal)]
##signal_num/sample_rate =时间

plt.figure()
plt.subplot(2,2,1)

plt.plot(signal_num/sample_rate, signal, color='blue')

plt.rcParams['font.sans-serif']=['SimHei']
plt.rcParams['axes.unicode_minus'] = False
plt.xlabel("Time(s)")
plt.ylabel("Amplitude(波幅)")
plt.title("signal of Voice(原始) ")

第二步：信号的均值滤波处理

"--------------第二步：信号的均值滤波处理-------------"
pre_emphasis = 0.97
#np.append(arr, values, axis=None)如果axis没有给出，那么arr，values都将先展平成一维数组
#signal原始形状（176400，2）执行完上述代码后，emphasized_signal形状为（352800,）
emphasized_signal = np.append(signal[0], signal[1:] - pre_emphasis * signal[:-1])
#print(np.shape(emphasized_signal))
#print(emphasized_signal[:10])

emphasized_signal_num = np.arange(len(emphasized_signal))
plt.subplot(2,2,2)
plt.plot(emphasized_signal_num/sample_rate, emphasized_signal, color='black')

plt.xlabel("Time(s)")
plt.ylabel("Amplitude")
plt.title("signal of Voice(均值滤波)")

第三步:信号的分帧处理

"--------------第三步:信号的分帧处理-------------"
frame_size = 0.025 #单位：秒
frame_stride = 0.01
#frame_length：一帧对应多少采样点，frame_step：帧_步幅 对应多少采样点
frame_length, frame_step = frame_size * sample_rate, frame_stride * sample_rate  #将帧时长转化为采样率
#emphasized_signal形状为（352800,）
#signal_length = 352800
signal_length = len(emphasized_signal)

frame_length = int(round(frame_length))#round 对浮点数四舍五入

frame_step = int(round(frame_step))
#1、计算帧数量 num_frames
#用于计算在给定的信号长度（signal_length）和帧长度（frame_length）以及帧步长（frame_step）下，最少需要多少帧（num_frames）来覆盖或近似覆盖整个信号。
#np.ceil 向上取整
num_frames = int(np.ceil(float(np.abs(signal_length - frame_length)) / frame_step)) #确保最少有一个frame存在
#2、计算信号长度（总采样数） pad_signal_length
pad_signal_length = num_frames * frame_step + frame_length
#3、对 emphasized_signal 进行填充，以便其长度达到 pad_signal_length
pad_signal = np.append(emphasized_signal, np.zeros((pad_signal_length - signal_length)))
#4、创建帧索引
indices = np.tile(np.arange(0,frame_length),(num_frames,1))+np.tile(np.arange(0,num_frames*frame_step,frame_step), (frame_length, 1)).T
#5、提取帧
frames = pad_signal[indices.astype(np.int32)]  

pad_signal_num = np.arange(len(pad_signal))

plt.subplot(2,2,3)
plt.plot(pad_signal_num/sample_rate, pad_signal, color='blue')
print(pad_signal_num/sample_rate)
plt.xlabel("Time(s)")
plt.ylabel("Amplitude")
plt.title("signal of Voice-分帧 ")

第四步-加汉明窗

# "--------------第四步-加汉明窗-------------"
frames *= np.hamming(frame_length)

第五步-FFT

# "--------------第五步-FFT-------------"
NFFT = 512
#频谱
mag_frames = np.absolute(np.fft.rfft(frames, NFFT))  # Magnitude of the FFT
#功率谱:对频谱取模平方
pow_frames = (1.0 / NFFT) * (mag_frames ** 2)

第六步：基于梅尔滤波器组的梅尔变换

# "--------------第六步：基于梅尔滤波器组的梅尔变换-------------"
low_freq_mel = 0
nfilt = 40
# 将频率转化为梅尔系数，计算最高频率的梅尔值，通常基于采样率和人耳的听觉范围。
high_freq_mel = (2595 * np.log10(1 + (sample_rate / 2) / 700)) 
#使用np.linspace在梅尔频率范围内均匀分布nfilt + 2个点（包括起点和终点）
mel_points = np.linspace(low_freq_mel, high_freq_mel, nfilt + 2)
# 将梅尔系数转化为频率
hz_points = (700 * (10**(mel_points / 2595) - 1))  
#计算每个滤波器在FFT频谱中的边界
#计算每个梅尔频率点对应的FFT bin索引。
bin = np.floor((NFFT + 1) * hz_points / sample_rate)

#4、创建梅尔滤波器组
fbank = np.zeros((nfilt, int(np.floor(NFFT / 2 + 1))))#np.floor 函数用于对数组中的每个元素进行向下取整
for m in range(1, nfilt + 1):
    f_m_minus = int(bin[m - 1])
    f_m = int(bin[m])
    f_m_plus = int(bin[m + 1])

    for k in range(f_m_minus, f_m):
        fbank[m - 1, k] = (k - bin[m - 1]) / (bin[m] - bin[m - 1])
    for k in range(f_m, f_m_plus):
        fbank[m - 1, k] = (bin[m + 1] - k) / (bin[m + 1] - bin[m])
filter_banks = np.dot(pow_frames, fbank.T)

第七步：对数运算

# "--------------第七步：对数运算-------------"
filter_banks = np.where(filter_banks == 0, np.finfo(float).eps, filter_banks)  # Numerical Stability
filter_banks = 20 * np.log10(filter_banks)

第八步：离散余弦计算

# "--------------第八步：离散余弦计算-------------"
from scipy.fftpack import dct
num_ceps = 12
mfcc = dct(filter_banks, type=2, axis=1, norm='ortho')[:, 1 : (num_ceps + 1)]
(nframes, ncoeff) = mfcc.shape
n = np.arange(ncoeff)
cep_lifter =22
lift = 1 + (cep_lifter / 2) * np.sin(np.pi * n / cep_lifter)
mfcc *= lift


mfcc_features=mfcc.T
plt.subplot(2,2,4)
ax = plt.matshow(mfcc_features[:48,:48], cmap=plt.cm.Blues)
plt.title('MFCC')
plt.show()

ax = plt.matshow(mfcc_features[:48,:96], cmap=plt.cm.Blues)
plt.title('MFCC')
plt.show()