简介：面向环境噪声、工业设备噪声及声学测试数据的合规性评估，这套MATLAB脚本实现A计权处理与1/3倍频程声压级计算。脚本内置ISO 266和IEC 61672推荐的滤波器系数，支持自定义采样率与频带范围，即使没有额外工具箱也可在MATLAB R2015b及以上版本直接运行，适合声学测试人员与相关专业学生快速掌握频带级噪声分析方法。压缩包内共8个文件，除MATLAB主程序外，还附带了Python辅助脚本、txt说明文档以及png曲线图等，便于理解A计权曲线与计权网络特性，整体大小约250KB，结构精简、查用方便。已有102人学习参考。输出结果以结构化数据集呈现，逐一给出各1/3倍频程中心频率对应的A加权声压级，可直接用于后续绘图、导出或与限值比对，在环境噪声监测、工业设备噪声评价及声学测试数据复现中具有较高实用价值。 我最近在整理设备噪声测试数据时，正好需要把一批时域录音转换成 A 计权 1/3 倍频程声压级频谱。翻了一圈现成工具，要么是商业软件收费不低，要么是自带算法和声级计对不上，干脆自己写了一个 MATLAB 脚本。这个脚本的核心功能就一句话：输入一段声压时域信号，输出 1/3 倍频程各频带的 A 计权声压级，附带频谱图。如果你也在做环境噪声评估、产品噪声排查，或者想跟声级计实测结果互相对标，这篇文章里从算法原理到代码实现、从校准到避坑，应该都能省你不少时间。
这里先给结论：用 MATLAB 做 A 计权 1/3 倍频程声压级分析，完全可行，而且如果按标准公式正确实现，误差可以控制在零点几分贝以内。但前提是你得明白两条技术路线的区别，以及在低频段和滤波器设计上会踩哪些坑。下面按我实际开发的思路逐步拆开讲。
1. 为什么噪声分析要用 A 计权和 1/3 倍频程
1.1 人耳听觉特性与 A 计权的由来
噪声评价里最常出现的指标就是 dB(A)，也就是 A 计权声压级。人耳不是对所有频率都同样敏感，同样是 60 dB 的声压级，1 kHz 声音听起来往往比 100 Hz 的响得多。A 计权曲线就是模拟人耳在中等响度下的频率响应特性，对低频做了大幅度衰减，对高频略有提升，在 1 kHz 位置修正量为 0 dB。
工程上无论环境噪声标准、家电产品噪声限值，还是职业健康评估，绝大多数都直接采用 A 计权声压级作为最终评价量。所以我们在做声学测试时，不能只看宽带总声压级，更要看频率分布，而 A 计权恰好把“人耳听到的大小”这个维度放进了每一段频谱里。
1.2 1/3 倍频程的分辨率为什么够用
用 FFT 做频谱分析可以得到几百上千条谱线，但人耳对频率的分辨能力是有限的，更接近对数频率尺度下的恒定比例分辨率。1/3 倍频程把整个可听频段划分为一系列频带，每个频带的带宽是中心频率的 23% 左右，相邻中心频率之比是 2 的 1/3 次方。20 Hz 到 20 kHz 范围内大约分布 30 个频带，这已经能很好地识别噪声的主要特征，比如低频嗡嗡声、中频风扇噪声、高频气动噪声。
相比窄带 FFT，1/3 倍频程谱更贴近听觉感知，也方便和标准限值表中的频带声压级直接对比。如果你最后要出噪声检测报告，用的几乎都是 1/3 倍频程谱，而不是几百根的 FFT 谱线。
1.3 脚本要解决的实际问题
我需要的工具具体包括：读取一段校准后的声压信号（单位 Pa）、按标准 1/3 倍频程中心频率划分频带、在每个频带内计算 A 计权声压级、以 dB(A) 为单位输出频谱表和图形。这个脚本不是要替代专业声学软件，而是用来在开发阶段快速评估噪声源、对比不同方案的降噪效果，顺便验证第三方测试报告的数据。
2. 脚本整体架构：两条技术路线怎么选
2.1 路线一：时域滤波链路（标准但重）
经典声级计内部处理流程是这样的：原始信号先经过带通滤波器组分成 1/3 倍频程各频带信号，每一路分别经过 A 计权滤波器，再求 RMS 值得到声压级。这种方式的优点是严格符合 IEC 61672 标准，能处理非稳态、瞬态信号，还能做实时连续监测。
缺点也很明显：需要为每个频带设计带通滤波器，低频段滤波器阶数高了容易数值不稳，阶数低了衰减不陡峭，频率边界容易出问题。加上还要做 A 计权时域滤波，代码量和调试成本都不小。如果是长时连续噪声监测，这条路是必须的；但只是为了分析一段几分钟的录音，有点杀鸡用牛刀。
2.2 路线二：FFT 频带能量积分（轻快但有前提）
另一条思路是对整段信号做 FFT，得到功率谱密度，然后按 1/3 倍频程频带边界把功率谱积分。这种方法的原理是帕塞瓦尔定理：时域 RMS 平方等于频域能量积分。由于 FFT 本身已经包含了全部频率信息，就可以在频域直接做 A 计权加权，不需要设计任何滤波器。
它的局限在于要求信号在分析时长内近似平稳，同时低频段频率分辨率受制于数据长度。比如采样率 48 kHz、分析时长 1 秒，FFT 频率分辨率为 1 Hz，在 20 Hz 频带内只有 7 个频点，积分误差会比较大。但对于常规产品噪声测试，采样几秒甚至几十秒完全可行，实测结果跟声级计对比差值通常在 0.5 dB 以内。
2.3 我的最终选型
我的需求是快速出频谱图、方便批量跑不同工况数据，不需要实时，所以我选择 FFT 频带积分法做主体脚本，同时留一个滤波法接口用于小信号比对。再强调一句：如果要出正式报告或用来做型式试验，建议还是老老实实用符合 IEC 61260 的滤波链路，或者直接用专业声学前端。FFT 积分法更适合研发自测和前期排查。
3. A 计权滤波器工程实现：从模拟公式到 MATLAB 代码
3.1 标准 A 计权频率响应公式
A 计权在频域有明确的解析式。IEC 61672 给出的模拟频率响应是：
R_A(f) = (12194^2 * f^4) / ((f^2 + 20.6^2) * sqrt((f^2 + 107.7^2) * (f^2 + 737.9^2)) * (f^2 + 12194^2))
然后换算成以 dB 为单位的 A 计权修正量：
A_db(f) = 20 * log10(R_A(f)) + 2.00
这个 +2.00 dB 的常数是为了让 1 kHz 处的修正量正好等于 0 dB。注意有些老资料里没加这个 2 dB，直接套用会发现 1 kHz 附近总差 2 dB，这是一个很容易翻车的点。
MATLAB 里写一个简单的函数：
function A = A_weight(f) f1 = 20.6; f2 = 107.7; f3 = 737.9; f4 = 12194; num = f4^2 * f.^4; den = (f.^2 + f1^2) .* sqrt((f.^2 + f2^2) .* (f.^2 + f3^2)) .* (f.^2 + f4^2); A = 20 * log10(num ./ den) + 2.00; end3.2 频域加权实现方式
FFT 积分法中最直接的 A 计权应用，是把频带内每一条谱线的功率密度乘以对应频率的线性加权系数，再做频带积分。线性加权系数是 10^(A_db/10)，因为 A_db 本身是 20log10 表示的电压/声压比，功率比就是平方关系。
对应代码片段：
W_linear = 10.^(A_weight(f_band) / 10); energy_A = trapz(f_band, pxx_band .* W_linear);其中f_band是当前频带内的频率向量，pxx_band是功率谱密度。这样做的精度比单纯取中心频率加权高很多，尤其是 200 Hz 以下，A 计权曲线斜率变化很快，用中心频率一点代替整个频带会带来明显误差。
3.3 时域 IIR 滤波（需要实时运行时用）
如果你的应用是实时测量或信号很长，不能整段 FFT，那需要设计 A 计权时域滤波器。工程上常用二阶级联滤波器逼近模拟曲线。MATLAB 里可以用designfilt尝试：
d = designfilt('weighting', 'FilterOrder', 4, ... 'FrequencyWeights', 'A', 'SampleRate', fs); x_a = filtfilt(d, x);不过designfilt的加权滤波器功能依赖 DSP System Toolbox，如果没有这个工具箱，可以用频域法做替代，或者把模拟 A 计权传递函数通过双线性变换离散化，自己写双二阶滤波。我个人建议非实时场景优先用频域加权，简单且不容易出错。
4. 1/3 倍频程带通滤波与声压级计算
4.1 中心频率与频带边界生成
1/3 倍频程中心频率以 1000 Hz 为基准，按 2^(1/3) 等比排列。标准列出的常用频率从 1 Hz 到 20 kHz 共有几十个，但工程测试一般只关心 20 Hz 到 20 kHz。可以用一行代码生成：
fc = 1000 * 2.^((-30:30)/3); % 生成以 1000 Hz 为基准的系列 fc = fc(fc >= 20 & fc <= 20000); % 截取需要的范围 fc = [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 16000 20000];建议直接用上面这个列表，因为标准标称频率是规定好的，自己算容易把 31.5 写成 31.25 这类误差带进去。每个频带的下限和上限分别为：
fl = fc / 2^(1/6)
fh = fc * 2^(1/6)
4.2 如果采用滤波法，带通滤波器怎么设计
滤波法需要为每个频带设计带通滤波器。这里我推荐用 4 阶 Butterworth，双边零相位滤波。设计一个频带的例子：
fc_i = 1000; fl = fc_i / 2^(1/6); fh = fc_i * 2^(1/6); [b, a] = butter(4, [fl, fh] / (fs / 2), 'bandpass'); y = filtfilt(b, a, x); Lp_i = 20 * log10(rms(y) / 2e-5);注意butter的归一化频率以 Nyquist 频率为 1，所以截止频率要除以 fs/2。所有频带共用一段原始信号x，逐频带滤波后计算 RMS 即可。
但滤波法计算量较大，低频段滤波器阶数不足时通带边界精度差，需要小心验证。如果只是想让结果接近标准，建议把 Butterworth 阶数提高到 6 到 8，并用freqz观察通带纹波。
4.3 声压级计算、校准与总声压级合成
无论哪种方法，最后都回到声压级公式：
L_p = 20 * log10(p_rms / p_ref)，其中 p_ref = 20 μPa = 2e-5 Pa
如果用periodogram计算功率谱密度，频带内 A 计权均方声压为energy_A，那么该频带声压级就是：
L_pA = 10 * log10(energy_A / p_ref^2)
这里用 10log10 是因为 energy_A 本身是平方量（Pa^2）。
校准是另一个容易忽略的点。如果你的信号已经是校准好的 Pa 值，直接算没问题。但大多数采集系统导出的只是数字码值，需要先乘以校准系数。标准做法是用声校准器在 1 kHz 位置输出 94 dB 声压级，对应有效声压为 p_ref * 10^(94/20) ≈ 1 Pa。采集到这段信号后计算数字 RMS，校准系数就是目标 RMS 除以实测数字 RMS。有了这个系数，后面所有时域数据乘上去就能转成实际声压。
总 A 计权声压级可以用频带能量叠加：
L_total = 10 * log10(sum(10.^(L_pA / 10)))
但要注意，频带叠加得到的总级和直接对全频带 A 计权信号求 RMS 的结果应当非常接近。如果二者差超过 0.5 dB，多半是频带覆盖不完整或者某个频带积分有问题。
我最终写的核心计算循环大致是这样的：
p_ref = 2e-5; N = length(x); [pxx, f] = periodogram(x, hanning(N), N, fs); % 单边功率谱密度 LpA = zeros(size(fc)); for i = 1:length(fc) fl = fc(i) / 2^(1/6); fh = fc(i) * 2^(1/6); idx = f >= fl & f <= fh; f_band = f(idx); pxx_band = pxx(idx); W = 10.^(A_weight(f_band) / 10); energy_A = trapz(f_band, pxx_band .* W); LpA(i) = 10 * log10(energy_A / p_ref^2); endtrapz用梯形积分是为了适应频带内频率点不均匀的问题，虽然periodogram返回的频率是均匀分布的，但用trapz更通用，例如从pwelch输出非均匀频率时也能用。
绘图可以直接用bar：
figure; bar(1:length(fc), LpA); set(gca, 'XTick', 1:4:length(fc), 'XTickLabel', fc(1:4:end)); xlabel('中心频率 / Hz'); ylabel('A计权声压级 / dB(A)'); grid on;如果频带太多，横轴用对数坐标更清晰，可以用semilogx画折线，或者直接把fc映射到等距位置方便观察刻度。
5. 实测校验与避坑：我踩过的几个坑
5.1 低频段结果抖得厉害
第一次用 1 秒数据跑 31.5 Hz 频带，发现声压级比声级计低 3 个多 dB。原因是 1 秒数据在 31.5 Hz 频带内只有大约 31 个周期，而且 Hann 窗对低频分量的能量泄露影响很大。FFT 频率分辨率只有 1 Hz，31.5 Hz 频带带宽才 9 Hz 左右，频点太少，积分误差自然大。
解决办法有三个：一是把分析时长加到至少 5 到 10 秒，频率分辨率到 0.2 Hz 以下；二是用periodogram时做重叠平均，降低方差；三是对低频段单独使用滤波法，频带 RMS 更稳定。我后来统一把测试录音时长加长到 10 秒，低频结果就稳了，和声级计对上了。
5.2 和声级计对不上，排查出 2 dB 偏差
有个朋友跑我脚本，发现所有频段都比声级计高 2 dB，问我怎么回事。我一看他用的 A 计权公式没加 2.00，导致 1 kHz 处不是 0 dB，而是 -2 dB，等于给全部结果加了 2 dB。这类问题在网上一堆代码里都有，注意核对 1 kHz 处的修正量是否等于 0。
另一个常见偏差是校准系数没对，声级计显示 94 dB，采集系统数字 RMS 可能不是 1，必须用校准系数统一到 Pa 单位。否则所有频带整体平移，差值固定。
5.3 带通滤波器边界频点归属问题
FFT 积分法里，频带边界的频率点不能简单算到某一个频带里，否则带间会重叠或泄漏。我建议用严格的大小判断，例如idx = f >= fl & f <= fh，这样相邻频带在边界处只归入一边。为了更严谨，可以把边界频点对半拆给两个频带，但工程上影响很小，一般不做也行。
5.4 批量处理时的脚本组织建议
如果有一段多工况数据，我建议把核心处理封装成函数，输入是x、fs，输出是fc和LpA。再写一个外层循环读取不同文件，最后把结果汇总到表格里。这样可以避免每个文件都重跑一遍 FFT 和频带循环。另外，.mat文件记得同时保存采样率fs，很多丢失校准信息的问题都是采样率对不上导致的。
个人体会是，自己写的脚本做研发分析非常够用，但不要盲目挑战标准合规。如果项目需要正式计量报告，还是找专业声学仪器对照验证，用手持声级计在实验室里打一遍同样信号，结果在误差范围内再大规模使用脚本。这也是我目前的工作流：先用脚本快速扫一遍所有数据，选出超标的几个工况，再用标准声学系统复测确认。这样既省时间，又不至于在关键数据上翻车。
本文还有配套的精品资源，点击获取
