1. 声品质与粗糙度先搞懂你在评估什么1.1 声品质不只是“有多响”做噪声振动工程的人都会有这么个经历某个部件的A声级已经压得很低了但客户试驾、试听、试用之后依然皱眉头抱怨“声音不好听”“有点吵”“不够高级”。单纯用分贝数评价声音在低频成分多、调制特征明显的场景下经常失灵。这就引出“声品质”这个概念——它描述的是人耳对声音主观感受的综合评价包括响度、尖锐度、粗糙度、波动度、清晰度等多个维度。其中粗糙度Roughness专门捕捉那些让人烦躁的快速调制成分是评价机械噪声、风噪、电机噪声时绕不开的指标。我在实际项目里接触声品质分析最早就是从粗糙度切入的。原因是某型电机的电磁噪声在频谱上不算突出但人的主观感受就是“嗡嗡的、刺刺的、不安逸”拿响度、尖锐度解释半天也没法量化最后用粗糙度参数一比问题立刻清晰了。从那以后我就确信声品质分析里粗糙度绝不是凑数的附带参数而是能真正揭示“听感问题”的关键工具。这篇文章把我用MATLAB做粗糙度计算、以及围绕声品质做应用分析的经验整理一遍给同样要跟声音主观评价打交道的工程师一点参考。1.2 粗糙度的心理声学机理为什么调制让人烦躁粗糙度对应的主观感受是声音在幅度上产生中速调制时带来的“粗糙、沙哑、刺耳”的感觉。举一个比较直观的例子两个喇叭同时播放频率接近的正弦波人耳听到的不是两个音而是一个音在“打拍子”音量随着时间快速起伏。这种起伏叫“拍频”振幅包络以两个频率差值作为调制频率周期性变化。当调制频率落在20Hz到300Hz区间时人就明显感受到声音“变粗”了调制频率太高或太低粗糙感都会减弱。峰值感知大致出现在调制频率60Hz到80Hz附近。为什么会有这种感知特性跟人耳听觉外周的信号处理机制有关。人耳基底膜在接收声音后会把宽带信号分解成不同临界频带Bark尺度上的窄带成分每个频带相当于一个带通滤波器。听觉系统对这些频带内信号的幅度包络变化非常敏感包络的起伏率、起伏深度共同决定粗糙度的大小。换句话说粗糙度本质上是对“听觉谱-时间”域里调制能量的一种加权统计。而调制频率落在响度感知的积分时间常数范围内时感知会达到峰值超过这个范围听觉系统来不及“跟踪”包络变化粗糙感就下降。这个机理也是后面设计MATLAB算法时要遵循的主线。1.3 粗糙度与其他声品质参数怎么分工声品质参数各有分工响度对应“大不大”尖锐度对应“尖不尖、亮不亮”粗糙度对应“粗不粗、杂不杂”波动度对应“晃不晃、抖不抖”。实际工程里这些参数经常一起算、一起看才能解释复杂听感。比如一个声音又响又尖消费者可能说“吵得头疼”但如果还带明显粗糙度就会说“听起来很低级、很没质感”——粗糙度提升给人“工艺粗糙”的心理暗示这种影响用A声级完全解释不了。了解这个分工有个实际价值拿到一个噪声录音先算响度、尖锐度等多个参数再用粗糙度做补充定位往往能更准确地锁定噪声源。比如某设备噪声既有轴承缺陷引起的周期性调制也有齿轮啮合的高频成分单看窄带频谱可能要猜半天算完粗糙度和尖锐度之后缺陷特征、啮合特征的贡献就分开了。这也是我在应用分析方法里把多参数联立作为“基本功”的原因。2. 读懂算法粗糙度计算的核心逻辑2.1 时变响度模型是粗糙度计算的地基粗糙度计算在数学上并不复杂但它的每一步都依赖一个听觉模型。业界比较成熟的做法是先计算“时变响度谱”也就是把声音信号映射到Bark尺度上输出随时间变化的响度密度值然后对特定频带范围内的响度信号做包络提取分析包络的调制深度和调制率最后经过权重求和得到粗糙度值。这个框架听起来清晰实际工程里细节非常多。2.2 Bark滤波器组与等响度校正人耳对频率的分辨是非线性的频率越高可分辨的带宽越大。模拟这种特性最常用的是Bark尺度将0Hz到约16kHz分为24个临界频带。计算时常使用Gammatone滤波器或双二阶滤波器组模拟听觉外围滤波效果。MATLAB里可以用filter函数配合不同中心频率和带宽的系数来实现也可以直接调用Audio Toolbox中的滤波器组工具不过Audio Toolbox的滤波器组是针对音频效果的直接用于声品质计算通常还需要自己校准带宽和增益我建议还是按照经典文献里的系数表来搭。等响度校正常被新手忽略但它很关键。人耳对3kHz附近的声波最敏感对低频和高频感知偏弱同一物理声压在不同频率上产生的响度感知完全不同。粗糙度计算需要反映真实感知输入到模型里的信号不能直接使用线性声压值而是要做频率计权处理。我的做法是在进入滤波器组之前对信号做FFT按1/3倍频程或Bark带统计声压谱然后叠加等响度曲线的修正量再反变换到时域。这个步骤能够显著减少低频成分对粗糙度结果的虚假放大。2.3 包络提取与调制分析在得到各Bark带的时变响度信号后下一步是提取每个频带内响度包络的调制深度。经典算法通常对时变响度信号进行高通滤波隔离掉慢变化的成分保留20Hz到300Hz的调制成分。我参考的模型在这里用了一个系数可调的滤波器中心频率范围对应粗糙度感知最灵敏的调制频段。调制深度的计算可以用Hilbert变换求瞬时包络也可以对信号做短时傅里叶变换来追踪包络的变化幅度。Hilbert变换的方法在MATLAB里实现很简洁直接用hilbert函数求解析信号再对包络做标准差统计就能得到调制深度。但要注意Hilbert变换对信号平稳性比较敏感实际采集的噪声信号往往有转速波动直接globally处理容易失真。我的经验是先做分帧处理每帧时长取100ms到200ms帧间重叠50%逐帧算调制深度再按帧做能量加权平均这样能把非平稳因素的影响控制住。2.4 粗糙度合成与权重设计各频带的调制深度算出来后需要按贡献权重叠加成单一粗糙度值。权重设计参考的是心理声学实验结论中频段的调制贡献最大低频和高频的调制贡献相对较低调制率接近70Hz时权重最高偏离后权重递减。经典模型会用相对调制深度、调制率权重函数等因子对乘积做进一步处理最终输出单位为asper1 asper定义为1kHz、60dB、100%调制深度、70Hz调制率引起的粗糙感。我实现的版本还额外加入一个针对Bark带的带宽修正因子防止窄带信号被过度放大。整个合成公式看起来不长但每个系数都是从文献实验数据里拟合出的不要随手乱改。如果只是为了做算法对比可以按我下文给出的参数搭一套默认值再拿标准信号验证确认结果与参考值误差在可接受范围内。3. 基于MATLAB的粗糙度计算工具实现3.1 数据准备与信号预处理写代码之前先处理好音频文件。粗糙度计算对采样率不敏感但建议用44.1kHz或48kHz标准采样率过低采样率会丢失高频Bark带的信息过高则徒增计算量。时长方面理想情况是计算2秒以上的稳态段对瞬态声音至少保留0.5秒才能获得稳定包络统计。数据准备的几个关键点确认声道。粗糙度计算针对单声道信号多声道录音要合成单声道后再处理。混音时注意别把相位抵消的成分直接平均掉建议取能量较大的通道。去除直流偏置。采集设备经常混入DC分量对包络提取造成干扰用detrend函数或直接减去均值即可。单位统一。MATLAB读取wav文件后得到的是归一化幅值需要根据校准系数换算为声压值或声压级才能得到准确的绝对粗糙度。如果只是相对比较保持线性幅值也能计算但结果只能做趋势比较。3.2 核心代码结构与关键函数下面是我常用的一个粗糙度计算核心函数框架不贴完整工程代码把关键步骤和参数说明写清楚。function R roughness_calc(x, fs) % x: 单声道信号物理声压Pa % fs: 采样率Hz % 1. 分帧参数 frameLen round(fs * 0.2); % 200ms 分帧 hopLen round(frameLen * 0.5); % 50%重叠 % 2. Bark滤波器组简化示例实际系数查文献表 nBands 24; [filterCoeffs, centerFreqs] getBarkFilters(fs, nBands); % 3. 逐帧处理 R_sum 0; weight_sum 0; for bi 1:nBands y filter(filterCoeffs{bi}, 1, x); % 提取时变响度对每一帧求RMS得到响度包络 env zeros(nFrames, 1); for fi 1:nFrames seg y(idx(fi):idx(fi)frameLen-1); env(fi) rms(seg); end % 高通滤波包络保留调制成分 env_f highpass(env, 20, fs/frameLen); % 计算调制深度归一化标准差 modDepth std(env_f) / mean(env); % 调制率估计用FFT找包络主频 [modRate, modMag] estimateModRate(env, fs/frameLen); % 权重函数 w roughnessWeight(modRate, bi); R_sum R_sum w * modDepth * modMag; weight_sum weight_sum w; end R R_sum / weight_sum; end不瞒你说这个简化版本在算法细节上还不够专业真正的实现需要整合时变响度模型、响度密度谱成形、相对调制深度计算等功能。但对于想快速跑通“粗糙度≈包络调制程度”逻辑的人来说这个框架是很好的起点。有几个细节我在实践里反复踩过坑时变响度的计算不能直接用线性RMS而要经过Bark带响度密度转换。简单做法是先将滤波器输出做窄带声压级分析再套用ISO 532B中的响度密度近似公式。高通滤波器的截止频率和阶数会影响调制成分的保留。默认高通20Hz、IIR滤波、2阶就够了太高阶会引入相位失真影响包络形状。highpass函数在较新版本MATLAB中依赖Signal Processing Toolbox如果没有该工具箱可以用代码手动实现两阶Butterworth高通效果相近。3.3 GUI工具设计让非算法工程师也能用算法写好后纯命令行调用对日常试验分析已经够用但工程团队的结构通常比较复杂。试验人员、NVH工程师、设计工程师各司其职不是每个人都有MATLAB使用经验。为了方便团队使用、减少沟通成本我后来把粗糙度计算封装成了一个小型GUI工具界面包含音频加载、参数设置、结果曲线和报告导出几个区域。界面实现采用App Designer底层调用上面提到的算法函数。关键交互设计有几点经验文件和参数分离。音频文件通过文件选择按钮加载参数设置面板单独放采样率、信号类型、频率权重开关等选项。图形展示时把时域波形、Bark谱图、时变响度和粗糙度结果放在同一个选项卡里方便对照查看。导出功能直接生成Excel报告内容包括信号基本参数、各Bark带调制深度表、最终粗糙度值以及判断结论可选。算法运行前先做一次信号完整性检查提示是否存在削波、过短、空文件等情况避免后续计算出无意义结果。GUI工具上线之后团队内部对声品质的关注度明显提高。数据显示过去一个噪声问题从发现到定位平均需要一周现在有了量化指标通常一天内就能锁定关键频带效率提升很直观。3.4 计算结果的验证方法算法搭完之后必须验证结果可靠直接用matlab自带的audiowrite生成已知参数的声音来校验是比较快的做法。我的验证流程生成1kHz正弦波幅度调制频率70Hz、调制深度100%理论粗糙度接近1 asper计算值与理论值误差应在5%以内。生成没有调制信号的稳态白噪声粗糙度应该很低接近0。生成不同调制频率的信号30Hz、70Hz、150Hz粗糙度应呈现峰值在70Hz附近的趋势。用已知声品质分析软件如Head Acoustics Artemis导出的测试音频做交叉验证如果计算结果趋势一致说明算法可信。不少人在第一步就容易卡住因为不同模型对“标准声”的定义略有不同。我的经验是拿文献里公开的标准测试信号做输入如果计算结果和文献引用值偏差在合理范围内就果断相信算法不要为了数值精确到小数点后两位而过度调参粗糙度是主观感知的量化本身带有统计不确定性差异控制在15%以内对工程分析完全够用。4. 声品质应用分析的方法与实践4.1 多参数联动的完整分析框架粗糙度单独算出来只是第一步声品质应用于工程决策靠的是多个参数联立对比。我在实际项目里常用一套“三阶分析法”第一阶物理信号分析。对噪声信号做FFT、阶次分析、小波分析定位主要频率成分和调制特征建立“声音问题在哪”的初步判断。第二阶声品质参数计算。对同一组信号批量计算响度、尖锐度、粗糙度、波动度形成参数矩阵看哪个参数随工况变化最明显哪个参数与主观评价最能对上号。第三阶统计相关分析。把声品质参数与主观评价打分做相关分析或对多个工况的参数做聚类找到影响听感的关键自变量。这三步走完通常能形成一个因果关系链结构共振→产生强烈调制成分→粗糙度上升→主观听感变差→用户投诉。链路清晰下一步做结构优化或声学包装就有了明确靶点。多参数联立时要特别注意参数之间的相关性。比如尖锐度升高往往伴随粗糙度升高因为它们都受高频成分影响但这不代表可以相互替代。以我的经验粗糙度对“周期性冲击声”更敏感尖锐度对“宽带高频嘶声”更敏感两者结合才能区分出问题类型。比如设备出现轴承早期损伤时域波形里有周期性冲击粗糙度显著上升而尖锐度变化不大如果出现气流啸叫尖锐度显著上升而粗糙度变化较小。4.2 工况对比与设计迭代声品质分析最重要的应用场景之一是不同设计方案的对比评估。举个例子某产品有两个结构设计方案方案A的声压级略高但主观试听反馈“更柔和”方案B的声压级更低但听感“有点烦躁”。这个时候单看声压级会得出“方案B更好”的结论算完粗糙度和波动度之后却能发现方案B的粗糙度比方案A高出一大截。再进一步分析是结构件的模态频率与激励源的调制频率接近产生了共振放大的调制成分。这样设计方向就清楚了不是简单降低声压级而是要避开某个频段的调制放大。工况对比时有个操作细节值得提一下确保所有工况的转速、负载、环境噪声本底一致计算窗口尽量取自稳态段。如果条件允许每个工况重复采集3次以上取粗糙度平均值和标准差能够有效避免偶发噪声对结果的干扰。我们曾遇到过同型号两台设备粗糙度差异超过20%的情况开始以为是装配差异后来反复测量才发现其中一台旁边刚好有个变频器在间歇工作采集窗里混入了外部噪声源。4.3 从实验室数据到产品定义声品质参数最终是为产品开发服务的。工程上用于建立噪声指标体系里常用声品质参数作为“主观验收指标”配合传统声压级限值一起控制产品质量。粗糙度作为其中一项常常对应客户“听起来很细腻”的需求。在实际制定指标时建议采用“分级阈值”的方法而不是单一限值。比如某类产品优秀等级的粗糙度小于0.5 asper合格等级小于0.8 asper超出则判定为问题项。分级阈值配合主观评价数据回归得到保证标准不会过于严格或过于宽松。我记得某个项目中产品在粗糙度指标上卡住了很久后来发现是标准采样窗口取的时长太短导致包络统计不稳定延长到3秒后数据自然就稳定了。这类细节直接影响指标能否落地。4.4 面向报告与决策的呈现技巧工程分析的最终交付不只是给一个数字而是要让决策者看得懂、信得过。我处理粗糙度结果时通常会画三类图第一类是频谱图与时变包络图展示声音信号的物理基础第二类是Bark带调制深度热力图显示粗糙度贡献的主要频带第三类是粗糙度随工况变化的曲线图直接呈现不同转速、负载下的趋势。三类图配合既能说明“为什么得出这个结论”也让非声学专业的人直观理解。报告里尽量不要用“听起来不舒服”之类的主观措辞而是用“40Hz调制频率附近包络调制深度达到0.35对应听感粗糙感明显”这样的客观描述。对决策者来说一个量化的异常仓比十句定性描述更有说服力。5. 实操中遇到的常见问题与排查技巧5.1 信号采集端的坑粗糙度计算对采集质量要求比普通FFT分析更高因为算法依赖时域包络任何幅值异常都会直接影响包络提取结果。信号削波是头号问题。输入超过ADC量程时波形被削平包络失真严重粗糙度计算值会异常偏高。排查方法是加载信号后检查max(abs(x))是否接近满量程如果削波比例超1%建议重新采集。量化噪声也会干扰高频Bark带的包络提取设备底噪大时计算前必须先做时域降噪处理但注意低通滤波不能把20Hz到300Hz的调制信号滤掉。还有采样率和文件格式问题。有的采集系统默认输出16bit PCM量化信噪比有限如果做精细的声品质分析建议至少24bit采样。另外多通道设备记录时注意各通道相位一致性不同通道之间的微小时延会显著影响合成后的单声道信号包络分析前要确认通道对齐。5.2 参数设置不当导致的异常结果滤波器组设计参数不匹配是常见问题。比如Bark带滤波器阶数过低频带之间泄漏严重包络信号相互污染粗糙度计算结果偏大。Bark滤波器组通常至少4阶带外衰减要达到20dB以上才可用。分帧长度设置也很关键。帧太短小于50ms包络统计不充分调制深度估计随机波动大帧太长大于500ms非平稳成分被平均掉真实的粗糙度变化被掩盖。以我的经验稳态噪声分析用200ms帧长、50%重叠最合适瞬态冲击声分析建议缩短到50ms到100ms否则冲击特征会被抹平。高通滤波器截止频率设置过高会滤掉一部分有效调制粗糙度偏小设置过低又无法隔离开慢变化成分。建议固定为20Hz不要随意下调。我对不同设备噪声做过敏感度分析20Hz到30Hz区间对结果影响尚可接受超过50Hz之后明显抑制粗糙度需要注意。5.3 计算结果的合理性判断拿到一个粗糙度数值之后先别急着下结论先做合理性检查。不同声音类型的典型粗糙度范围如下声音类型典型粗糙度范围纯正弦稳态音0.1 asper以下普通空调室内噪声0.2-0.5 asper电机电磁噪声0.3-0.8 asper齿轮啮合冲击声0.5-1.5 asper严重异响1.5 asper以上如果计算结果明显超出对应声音类型的典型范围先检查信号是否存在削波、采集是否包含异常瞬态、滤波器参数是否设置正确再考虑算法实现是否有误。切忌硬套经验值不同模型、不同计权模式下的绝对数值会有差异横向对比时务必保持同一套参数。5.4 计算性能优化粗糙度算法涉及滤波器组逐带处理加上分帧统计数据量大的时候计算速度确实感人。一段60秒的音频24个Bark带全跑一遍纯MATLAB代码可能要几分钟。我在实际项目里做了几个优化滤波器组的处理可以改为矩阵运算用filter一次处理多通道数据减少循环开销。分帧时可以先用buffer函数一次性生成帧矩阵避免循环切帧。如果只是批量对比多个文件可以先用动态时间规划或VAD检测截取关键稳态段再针对截断数据计算能省一半以上时间还不影响判断结论。如果项目需要实时性可以考虑将核心循环改为C代码生成或MEX编译速度提升明显。写在最后的几点体会声品质粗糙度的MATLAB实现表面上是写算法实际上是三步工程问题理解听感机理、搭对听觉模型、跑通应用链路。很多人在第一步就没站稳以为粗糙度就是“先算响度再求标准差”结果换一个场景就失灵。我的建议是动手码代码之前先画一张信号处理流程图从采集到最终输出每一步处理的原因都写清楚这样做出来的模型才是真正可解释、可迁移的。另外一个小技巧共享一下在处理项目数据时保留一组“基准音频”并常态化跑一遍这样即使算法改进了也能快速确认新旧结果之间的差异是算法层面的还是数据层面的。我在团队里就是因为坚持维护这条标准曲线几次排查问题都能快速定位原因节省了大量时间。这个内容后续还可以扩展的方向包括对时变响度模型做更精细的适配、把粗糙度算法拓展到波动度联合计算、用深度学习做调制模式的自动分类等等。我目前正在做的是把多参数声品质计算整合到一个统一的MATLAB工具箱中让整个声品质分析流程更加自动化。如果你也在做相关工作欢迎多交流踩坑经验。