本节内容我们来看下如何用Matlab和Python计算声音的声压级和响度。

声压级

1. 声压级定义

首先来看声压级，这个就是指的我们平时所说的声音有多少分贝。声压定义为声波在某一点产生的逾量瞬时压强的均方根值。由于声压容易被人耳感知，也易于测量，因此，通常使用声压作为描述声波大小的物理量。

声压级以符号SPL(sound pressure level)表示，其定义为将待测声压有效值p(e)与参考声压p(ref)的比值取常用对数，再乘以20，即：

在空气中参考声压p(ref)一般取为2e-5帕，这个数值是正常人耳对800赫声音刚刚能觉察其存在的声压值，也就是800赫声音的可听阈声压。一般讲，低于这一声压值，人耳就再也不能觉察出这个声音的存在了。显然该可听阈声压的声压级即为零分贝。

指的是声压有效值，就是一段声音信号的均方根(RMS)。设语音长度为T, 离散点数为N, 则有效声压的计算公式为：

常见的声音分贝值：

声音种类

声压级(分贝)

烈性炸药爆炸声

170

火箭、导弹发射场

150

飞机发动机

120

发电机工作

100

卡车

90

正常谈话

50

轻声耳语

30

农村静夜

10

Matlab代码

由前面的定义可很容易写出SPL的代码，需要注意的是，程序中的输入信号是数字信号，与实际的模拟信号大小成倍数关系。

function spl = SPLCal(x)

len = length(x);

%%

% 有效声压计算，即求RMS

pa = sqrt(sum(x.^2)/len);

%% 声压级计算

% 声压级值spl=20*log10(pa/p0),单位为dB

p0 = 2e-5;

spl = 20*log10(pa/p0);

end

完整代码如下：

clear all;clc;close all;

%%

[x,fs]=audioread('mySpeech.wav');

len=length(x);

%% 语音分帧

% 每帧大小为frameLen，当语音长度不是帧长的整数倍时：

% (1)若剩余长度大于等于帧长的二分之一，则补零至帧长

% (2)若剩余长度小于帧长的二分之一，则舍弃

% 常用的语音帧长：20ms、50ms、100ms、200ms

framTime = 100; % 单位：ms

% 每帧信号点数

frameLen=fs*framTime/1000;

% m为Length/frameLen后得到的余数

m = mod(len,frameLen);

if m >= frameLen/2 % 补零

x = [x;zeros(frameLen-m,1)];

len = length(x);

else   % 即m

nframe = floor(len/frameLen);

x = x(1:nframe*frameLen);

len = length(x);

end

% 最终的语音分帧总帧数

N = len/frameLen;

%% 计算声压级

s = zeros(1,frameLen);

% N帧信号的声压级值存储在spl向量里

spl = zeros(1,N);

for k = 1:N

s = x((k-1)*frameLen + 1:k*frameLen);

spl(k) = SPLCal(s);

end

%% 绘图

t = 1:len;

spl_rep = repmat(spl, frameLen, 1);

figure();

subplot(211);plot(t/fs,x);grid on;xlabel('时间(s)');title('输入语音波形');

subplot(212);stairs(t/fs,spl_rep(:),'r');grid on;xlabel('时间(s)');ylabel('声压级(dB)');

title('语音信号的声压级(dB)');

image-20210220185851447

Python代码

Python代码如下：

import pyaudio

import wave

import numpy as np

import matplotlib.pyplot as plt

def load_wav(wave_input_path):

wf = wave.open(wave_input_path, 'rb')  # 读 wav 文件

fs = wf.getframerate()

nframes = wf.getnframes()

str_data = wf.readframes(nframes)

wf.close()

wave_data = np.fromstring(str_data, dtype=np.short)

return wave_data.astype(np.float64), fs

def SPLCal(x):

Leng = len(x)

pa = np.sqrt(np.sum(np.power(x, 2))/Leng)

p0 = 2e-5

spl = 20 * np.log10(pa / p0)

return spl

if __name__ == '__main__':

x, fs = load_wav('audio.wav')

Leng = len(x)

frameTime = 100

frameLen = fs * frameTime // 1000

m = np.mod(Leng, frameLen)

if m>=frameLen/2:

x = np.append(x, np.zeros(int(frameLen-m)))

Leng = len(x)

else:

nframe = np.floor(Leng/frameLen)

x = x[0:nframe * frameLen + 1]

Leng = len(x)

N = Leng // frameLen

spl = np.array([])

for k in range(N):

s = x[k*frameLen: (k+1)*frameLen]

spl = np.append(spl, SPLCal(s))

spl_rep = np.repeat(spl, frameLen)

plt.figure()

plt.subplot(211)

plt.plot(x)

plt.subplot(212)

plt.plot(spl_rep)

plt.show()

绘图如下：

image-20210221155302000

响度

当用同样的力气讲话的时候，为什么我们总觉得女性的声音要比男性的响？这就是我们下面要讲的响度。响度是听觉判断声音强弱的属性，跟人主观感觉有关。人主观感觉判断的声音强弱，即声音响亮的程度，根据它可以把声音排成由轻到响的序列。

当外界声振动传入人耳内，人们在主观感觉上形成听觉上声音强弱的概念。人们习惯于地用“响”与“不响”来描述声波的强度，但这一描述与声波的强度又不完全等同。人耳对声波响度的感觉还与声波的频率有关，即使相同声压级但频率不同的声音，人耳听起来会不一样响。例如，同样是60dB的两种声音，但一个声音的频率为100Hz，而另一个声音为1000Hz，人耳听起来1000Hz的声音要比100Hz的声音响。要使频率为100Hz的声音听起来和频率为1000Hz、声压级为60dB的声音同样响，则其声压级要达到67dB。

下面介绍几个相关的概念：

响度级： 按人耳对声音的感觉特性，依据声压和频率定出人对声音的主观音响感觉量，称为响度级，单位为方。

方(Phon)：当某一频率的纯音和1000Hz的纯音听起来同样响时，这时1000Hz纯音的声压级就定义为该待定声音的响度级。因此在1kHz的频率上，声压级为60dBSPL信号的响度为60方。对各个频率的声音作这样的听音比较，得出达到同样响度级时频率与声压级的关系曲线，这就是我们人耳的听觉等响曲线。

image-20210220112716013

从等响曲线图中我们发现，人耳对高频的声音更加敏感，同样声压级下的高频声音响度级比低频的高。一般女性发声的高频成分较多，而男性发声的低频成分相对较多，这就是在同样力气讲话时(声压级相同)，女性的声音听上去更加响的原因。

由于这种客观单位只是非常有限地表达了人耳对于响度的反应，因此可以引入一个关于响度的主观概念——宋。

宋(Sone)：表示人耳在自然状态下，根据声压级的变化所表现出的对于响度听感的变化。

“宋”与“方”的关系表现为1宋等于40方(即在等响曲线图中，1kHz处代表40dBSPL)，并且以1宋为标准，在2宋时响度增加一倍，在0.5宋时响度减小一倍。

根据IS0226--2003标准H等响度曲线的定义，声压级LP 为：

其中，

其中，为听力阈值；为响度感知指数；为以1000Hz为标准所计算的线性传输函数的幅值。这三个参数都可以在ISO226中查到。

image-20210220130314703

计算响度的Matlab代码如下：

function [spl, freq] = loudnessCal(phon)

f = [20 25 31.5 40 50 63 80 100 125 160 200 250 315 400 500 630 800 ...

1000 1250 1600 2000 2500 3150 4000 5000 6300 8000 10000 12500];

af = [0.532 0.506 0.480 0.455 0.432 0.409 0.387 0.367 0.349 0.330 0.315 ...

0.301 0.288 0.276 0.267 0.259 0.253 0.250 0.246 0.244 0.243 0.243 ...

0.243 0.242 0.242 0.245 0.254 0.271 0.301];

Lu = [-31.6 -27.2 -23.0 -19.1 -15.9 -13.0 -10.3 -8.1 -6.2 -4.5 -3.1 ...

-2.0  -1.1  -0.4   0.0   0.3   0.5   0.0 -2.7 -4.1 -1.0  1.7 ...

2.5   1.2  -2.1  -7.1 -11.2 -10.7  -3.1];

Tf = [ 78.5  68.7  59.5  51.1  44.0  37.5  31.5  26.5  22.1  17.9  14.4 ...

11.4   8.6   6.2   4.4   3.0   2.2   2.4   3.5   1.7  -1.3  -4.2 ...

-6.0  -5.4  -1.5   6.0  12.6  13.9  12.3];

Ln = phon;

%从响度级计算声压级

Af=4.47E-3 * (10.^(0.025*Ln) - 1.15) + (0.4*10.^(((Tf+Lu)/10)-9 )).^af;

Lp=((10./af).*log10(Af)) - Lu + 94;

spl = Lp;

freq = f;

end