# 语音信号处理基础（十）—

文章目录

1.定义
2.Mel频率分析
3.编写Mel滤波器函数

1.定义

在语音识别（SpeechRecognition）和话者识别（SpeakerRecognition）方面，最常用到的语音特征就是梅尔倒谱系数（Mel-scaleFrequency Cepstral Coefficients，简称MFCC）。根据人耳听觉机理的研究发现，人耳对不同频率的声波有不同的听觉敏感度。从200Hz到5000Hz的语音信号对语音的清晰度影响对大。两个响度不等的声音作用于人耳时，则响度较高的频率成分的存在会影响到对响度较低的频率成分的感受，使其变得不易察觉，这种现象称为掩蔽效应。由于频率较低的声音在内耳蜗基底膜上行波传递的距离大于频率较高的声音，故一般来说，低音容易掩蔽高音，而高音很难掩蔽低音。在低频处的声音掩蔽的临界带宽较高频要小。所以，人们从低频到高频这一段频带内按临界带宽的大小由密到疏安排一组带通滤波器，对输入信号进行滤波。将每个带通滤波器输出的信号能量作为信号的基本特征，对此特征经过进一步处理后就可以作为语音的输入特征。由于这种特征不依赖于信号的性质，对输入信号不做任何的假设和限制，又利用了听觉模型的研究成果。因此，这种参数比基于声道模型的LPCC的鲁棒性(所谓“鲁棒性”，也是指控制系统在一定（结构，大小）的参数变动下，维持其它某些性能的特性。比如说，计算机软件在输入错误、磁盘故障、网络过载或有意攻击情况下，能否不死机、不崩溃，就是该软件的鲁棒性。)更好，相比更符合人耳的听觉特性，而且当信噪比降低时仍然具有较好的识别性能。

2.Mel频率分析

我们刚才做了倒谱分析。给我们一段语音，我们可以得到了它的频谱包络（连接所有共振峰值点的平滑曲线）。但是，对于人类听觉感知的实验表明，人类听觉的感知只聚焦在某些特定的区域，而不是整个频谱包络。

而Mel频率分析就是基于人类听觉感知实验的。实验观测发现人耳就像一个滤波器组一样，它只关注某些特定的频率分量（人的听觉对频率是有选择性的）。也就说，它只让某些频率的信号通过，而压根就直接无视它不想感知的某些频率信号。但是这些滤波器在频率坐标轴上却不是统一分布的，在低频区域有很多的滤波器，他们分布比较密集，但在高频区域，滤波器的数目就变得比较少，分布很稀疏。

人的听觉系统是一个特殊的非线性系统，它响应不同频率信号的灵敏度是不同的。在语音特征的提取上，人类听觉系统做得非常好，它不仅能提取出语义信息, 而且能提取出说话人的个人特征。如果在语音识别系统中能模拟人类听觉感知处理特点，就有可能提高语音的识别率。

梅尔频率倒谱系数（Mel Frequency Cepstrum Coefficient, MFCC）是在Mel标度频率域提取出来的倒谱参数，Mel标度描述了人耳频率的非线性特性。也就是MFCC考虑到了人类的听觉特征，先将线性频谱映射到基于听觉感知的Mel非线性频谱中，然后转换到倒谱上。

将普通频率转化到Mel频率的公式是：
Fmel(f)=2595∗lg⁡(1+f/700)F_{m e l}(f)={\color{Red}2595} * \lg(1+f / 700) Fmel(f)=2595∗lg(1+f/700)

Fmel(f)是以Mel为单位的感知频率，f是以Hz为单位的实际频率。其中2595是以10为底得到的。以e为底得到1125F_{m e l}(f)是以Mel为单位的感知频率，f是以Hz为单位的实际频率。其中2595是以10为底得到的。以e为底得到1125 Fmel(f)是以Mel为单位的感知频率，f是以Hz为单位的实际频率。其中2595是以10为底得到的。以e为底得到1125

从梅尔回到频率的公式为：
Fmel−1(m)=700(10(m/2595)−1)F_{m e l}^{-1}(m)=700(10^{(m / 2595)} -1) Fmel−1(m)=700(10(m/2595)−1)
下图展示了Mel频率（感知频率）与线性频率的关系：它可以将不统一的频率转化为统一的频率，也就是统一的滤波器组。

在Mel频域内，人对音调的感知度为线性关系。举例来说，如果两段语音的Mel频率相差两倍，则人耳听起来两者的音调也相差两倍。

3.编写Mel滤波器函数

以M=6（6个滤波器）为例，fl、fh、fs、N、M的值一旦确定就固定不变。
将F_{mel}(fl)至F_{mel}(fh)的Mel频率范围均分为M+1=7段，
产生M+2=8个Mel频率值。
Fmel(fl)+0Fmel(fh)−Fmel(fl)M+1第1个Mel频率值Fmal(fl)+1Fmel(fh)−Fmel(fl)M+1第2个Mel频率值Fmal(fl)+2Fmel(fh)−Fmel(fl)M+1第3个Mel频率值Fmel(fl)+3Fmel(fh)−Fmel(fl)M+1第4个Mel频率值Fmel(fl)+4Fmel(fh)−Fmel(fl)M+1第5个Mel频率值Fmel(fl)+5Fmel(fh)−Fmel(fl)M+1第6个Mel频率值Fmal(fl)+6Fmel(fh)−Fmel(fl)M+1第7个Mel频率值Fmel(fl)+7Fmel(fh)−Fmel(fl)M+1第8个Mel频率值\begin{array}{l}{F_{m e l}\left(f_{l}\right)+0 \frac{F_{m e l}\left(f_{h}\right)-F_{m e l}\left(f_{l}\right)}{M+1}}第1个Mel频率值 \\ {F_{m a l}\left(f_{l}\right)+1 \frac{F_{m e l}\left(f_{h}\right)-F_{m e l}\left(f_{l}\right)}{M+1}} 第2个Mel频率值\\ {F_{m a l}\left(f_{l}\right)+2 \frac{F_{m e l}\left(f_{h}\right)-F_{m e l}\left(f_{l}\right)}{M+1}}第3个Mel频率值 \\ {F_{m e l}\left(f_{l}\right)+3 \frac{F_{m e l}\left(f_{h}\right)-F_{m e l}\left(f_{l}\right)}{M+1}}第4个Mel频率值 \\ {F_{m e l}\left(f_{l}\right)+4 \frac{F_{m e l}\left(f_{h}\right)-F_{m e l}\left(f_{l}\right)}{M+1}}第5个Mel频率值\\ {F_{m e l}\left(f_{l}\right)+5 \frac{F_{m e l}\left(f_{h}\right)-F_{m e l}\left(f_{l}\right)}{M+1}}第6个Mel频率值 \\ {F_{m a l}\left(f_{l}\right)+6 \frac{F_{m e l}\left(f_{h}\right)-F_{m e l}\left(f_{l}\right)}{M+1}} 第7个Mel频率值\\ {\quad F_{m e l}\left(f_{l}\right)+7 \frac{F_{m e l}\left(f_{h}\right)-F_{m e l}\left(f_{l}\right)}{M+1}}第8个Mel频率值\end{array} Fmel(fl)+0M+1Fmel(fh)−Fmel(fl)第1个Mel频率值Fmal(fl)+1M+1Fmel(fh)−Fmel(fl)第2个Mel频率值Fmal(fl)+2M+1Fmel(fh)−Fmel(fl)第3个Mel频率值Fmel(fl)+3M+1Fmel(fh)−Fmel(fl)第4个Mel频率值Fmel(fl)+4M+1Fmel(fh)−Fmel(fl)第5个Mel频率值Fmel(fl)+5M+1Fmel(fh)−Fmel(fl)第6个Mel频率值Fmal(fl)+6M+1Fmel(fh)−Fmel(fl)第7个Mel频率值Fmel(fl)+7M+1Fmel(fh)−Fmel(fl)第8个Mel频率值
分别求出这M+2=8个Mel频率对应的实际频率值
Fmel−1(Fmel(fl)+iFmel(fh)−Fmel(fl)M+1)0≤i≤8F_{m e l}^{-1}\left(F_{m e l}\left(f_{l}\right)+i \frac{F_{m e l}\left(f_{h}\right)-F_{m e l}\left(f_{l}\right)}{M+1}\right) \quad 0 \leq i \leq 8 Fmel−1(Fmel(fl)+iM+1Fmel(fh)−Fmel(fl))0≤i≤8
分别求出这M+2个实际频率值对应的FFT点数（M以6为例）
f(0)=(Nfs)Fmel−1(Fmel(f1)+0Fmel(fh)−Fmel(fl)M+1)f(1)=(Nfs)Fmel−1(Fmel(fl)+1Fmel(fh)−Fmel(fl)M+1)f(2)=(Nfs)Fmel−1(Fmel(fl)+2Fmel(fh)−Fmel(fl)M+1)f(3)=(Nfs)Fmel−1(Fmel(fl)+3Fmel(fh)−Fmel(fl)M+1)f(4)=(Nfs)Fmel−1(Fmel(fl)+5Fmel(fh)−Fmel(fl)M+1)f(6)=(Nfs)Fmel−1(Fmel(fl)+6Fmel(fh)−Fmel(fl)M+1)f(7)=(Nfs)Fmel−1(Fmal(fl)+7Fmel(fh)−Fmel(fl)M+1)\begin{array}{l}{f(0)=\left(\frac{N}{f s}\right) F_{mel}^{-1}\left(F_{m e l}\left(f_{1}\right)+0 \frac{F_{m e l}\left(f_{h}\right)-F_{m e l}\left(f_{l}\right)}{M+1}\right)} \\ {f(1)=\left(\frac{N}{f s}\right) F_{mel}^{-1}\left(F_{m e l}\left(f_{l}\right)+1 \frac{F_{m e l}\left(f_{h}\right)-F_{m e l}\left(f_{l}\right)}{M+1}\right)} \\ {f(2)=\left(\frac{N}{f s}\right) F_{mel}^{-1}\left(F_{m e l}\left(f_{l}\right)+2 \frac{F_{m e l}\left(f_{h}\right)-F_{m e l}\left(f_{l}\right)}{M+1}\right)} \\ {f(3)=\left(\frac{N}{f s}\right) F_{mel}^{-1}\left(F_{m e l}\left(f_{l}\right)+3 \frac{F_{m e l}\left(f_{h}\right)-F_{m e l}\left(f_{l}\right)}{M+1}\right)} \\ {f(4)=\left(\frac{N}{f s}\right) F_{mel}^{-1}\left(F_{m e l}\left(f_{l}\right)+5 \frac{F_{m e l}\left(f_{h}\right)-F_{m e l}\left(f_{l}\right)}{M+1}\right)} \\ {f(6)=\left(\frac{N}{f s}\right) F_{mel}^{-1}\left(F_{m e l}\left(f_{l}\right)+6 \frac{F_{mel}\left(f_{h}\right)-F_{m e l}\left(f_{l}\right)}{M+1}\right)} \\ {f(7)=\left(\frac{N}{f s}\right) F_{mel}^{-1}\left(F_{m a l}\left(f_{l}\right)+7 \frac{F_{mel}\left(f_{h}\right)-F_{m e l}\left(f_{l}\right)}{M+1}\right)}\end{array} f(0)=(fsN)Fmel−1(Fmel(f1)+0M+1Fmel(fh)−Fmel(fl))f(1)=(fsN)Fmel−1(Fmel(fl)+1M+1Fmel(fh)−Fmel(fl))f(2)=(fsN)Fmel−1(Fmel(fl)+2M+1Fmel(fh)−Fmel(fl))f(3)=(fsN)Fmel−1(Fmel(fl)+3M+1Fmel(fh)−Fmel(fl))f(4)=(fsN)Fmel−1(Fmel(fl)+5M+1Fmel(fh)−Fmel(fl))f(6)=(fsN)Fmel−1(Fmel(fl)+6M+1Fmel(fh)−Fmel(fl))f(7)=(fsN)Fmel−1(Fmal(fl)+7M+1Fmel(fh)−Fmel(fl))
编写Mel滤波函数MATLAB实现代码

fs = 8000;    %采样频率
fl = 0;       %fl是设计滤波器的最低频率
fh = fs / 2;  %fh是设计滤波器的最高频率
bl = 1125 * log(1 + fl/700);    %最低频率对应的Mel频率
bh = 1125 * log(1 + fh/700);    %最高频率对应的Mel频率
p = 24;       %在fl和fh之间设计Mel滤波器的个数
N = 256;      %FFT点数,N越大,频率分辨率越高,误差越小
MelF = linspace(0, bh-bl, p+2);    %在0至bh-bl的Mel频率范围内产生p+2个Mel频率值
F = 700 * (exp(MelF/1125) - 1);    %将上一步产生的p+2个Mel频率值转化为p+2个实际频率值
df = fs / N;      %计算频率分辨率
n = N/2 + 1;      %计算fs/2内对应的FFT点数
f = (0:n-1) * df; %计算频率序列,共有n个频率点
bank = zeros(24, n);    %生成24行n列的全零矩阵
for m = 2:p+1     %循环起始数为2，因为下方语句中用到m-1的值,其中m-1必须为正整数,所以m为大于等于2的正整数%F为实际频率,分别找到三角滤波器左中右对应的频率F_left = F(m-1);F_mid = F(m); F_right = F(m+1);%使用ceil（）函数向上取整，计算实际频率对应的FFT点数n_left = ceil(F_left / df);   %频率F_left对应的FFT点数n_mid = ceil(F_mid / df);     %频率F_mid对应的FFT点数n_right = ceil(F_right / df); %频率F_right对应的FFT点数%用三角窗的频率响应函数作图for k = 1:n    %根据FFT点数画三角窗if k>=n_left & k<=n_midbank(m-1,k) = (k-n_left)/(n_mid-n_left);    %画第m-1个滤波器的左半部分elseif k>n_mid & k<=n_rightbank(m-1,k) = (n_right-k)/(n_right-n_mid);  %画第m-1个滤波器的右半部分endendplot(f,bank(m-1,:),'r');xlabel('频率/Hz');ylabel('相对幅值');title('Mel滤波器组的频率响应曲线')hold on;
end

运行结果

MATLAb自带Mel滤波函数实现代码

clear all; clc; close all;
% 调用melbankm函数，在0-0.5区间设计24个Mel滤波器，用三角形窗函数
% melbank(p, n, fs, fl, fh, w);
% fs为采样频率
% fl是设计滤波器的最低频率（用fs进行归一，一般取0）
% fh是设计滤波器的最高频率（用fs进行归一，一般取0.5）
% p是在fl和fh之间设计Mel滤波器的个数
% n是一帧FFT后数据的长度
% w是窗函数，'t'代表三角窗，'m'代表汉明窗，'n'代表汉宁窗
bank = melbankm(24, 256, 8000, 0, 0.5, 't');
bank = full(bank);
bank = bank / max(bank(:)); % 幅值归一化
df = 8000 / 256;    % 计算分辨率
ff = (0 : 128) * df;
for k = 1 : 24plot(ff, bank(k,:), 'k');hold on;
end
hold off;
grid on;
xlabel('频率/Hz');
ylabel('相对幅值');
title('Mel滤波器组的频率响应曲线')

运行结果

C语言结合gnuplot、tcc实现代码

#include <stdio.h>
#include <math.h>int main ()
{int fs,fl,fh,P,N,f,i,j,M,k,m;double bl,bh,df;fs = 8000;         //采样频率fl = 0;           //fl是设计滤波器的最低频率fh = fs / 2;      //fh是设计滤波器的最高频率bl = 1125 * log(1 + fl/700);    //最低频率对应的Mel频率bh = 1125 * log(1 + fh/700);    //最高频率对应的Mel频率P = 24;       //在fl和fh之间设计Mel滤波器的个数N = 256;      //FFT点数,N越大,频率分辨率越高,误差越小double MelF[24+2],F[24+2];int n[24+2];for(i=0; i<=P+1; i++)MelF[i]=bl + 1.0*i*((bh-bl)/(P+1));//在0至bh-bl的Mel频率范围内分成P+1段，产生p+2个Mel频率值for(i=0; i<=P+1; i++)F[i] = 700 * (exp(MelF[i]/1125.0) - 1);    //将上一步产生的p+2个Mel频率值转化为p+2个实际频率值df = fs*1.0 / N;           //计算频率分辨率for(i=0; i<=P+1; i++)n[i]= ceil(1.0*F[i]/df);      //计算fs/2内对应的FFT点数int bank[fh+1];    for (m=1;m<=P;m++)     //循环起始数为1，因为下方语句中用到m-1的值,其中m-1必须为整数,所以m为大于等于1的整数{//用三角窗的频率响应函数作图for (k=0;k<=fh;k++){//根据FFT点数画三角窗if (k>=n[m-1] && k<n[m])bank[k] = 1.0*(k-n[m-1])/(n[m]-n[m-1]);    //画第m个滤波器的左半部分else if (k>=n[m] && k<=n[m+1])bank[k] = 1.0*(n[m+1]-k)/(n[m+1]-n[m]);  //画第m个滤波器的右半部分elsebank[k] =0;}   for(k=0;k<=N/2;k++)printf("%f\t%f\n", k*df, bank[k]);}}

tcc命令

%编译生成可执行文件Mel_filtering.exe
tcc Mel_filtering.c
%通过重定向的方法另存数据
Mel_filtering.exe > Mel_filtering.txt

gnuplot画图命令

%画图
gnuplot
%作出满足均匀分布的随机数的线状图(线的类型采用红色 linetype 7)
plot [0:3999] [0:1] "Mel_filtering.txt"u 1:2 w l lt 7

运行结果