简介面向语音信号处理初学者的MATLAB仿真资源围绕语音信号短时分析中的短时过零率、短时能量、短时自相关等经典特征提取方法展开适合通信工程、电子信息等专业学生进行课程设计或实验入门。资源基于MATLAB 2022A版本编写包内共含10个文件以8个.m源码文件为主体配套1个wav语音样本用于实验测试以及1个avi操作录像演示具体运行流程压缩包大小仅590KB轻量便捷。已有617人学习下载程序代码均配有中文注释关键函数如分帧、加窗、特征计算等均单独封装便于读者逐段理解算法实现同时录像中特别强调了MATLAB当前文件夹路径的设置要点可有效避免因路径问题导致的运行报错。整体结构清晰能帮助使用者快速掌握语音短时分析的编程实现与调试思路也可作为相关课程实验报告的参考素材。1. 从一段语音到“帧”短时分析为什么先把信号切碎接手语音信号处理的第一周很多人都会卡在同一个地方一段 5 秒的语音直接用傅里叶变换做频谱分析得到的是一条几乎看不出规律的曲线加窗、改点数、换函数都救不回来。原因不在算法而在你看待信号的方式。语音是典型的非平稳信号发音时声带振动频率、声道形状、音量都在变化用全局视角去分析整体特征结果自然发散。短时分析的核心做法是先承认这种非平稳性再把信号切成 20 到 30 毫秒的短帧在每一帧内默认信号近似平稳然后逐帧计算短时能量、短时过零率、短时自相关等特征将这些特征值连成曲线整段语音的变化规律才显现出来。这篇文章要解决的问题很具体如何用 MATLAB 从零搭建一套短时分析仿真覆盖短时能量、短时过零率、短时自相关三个经典参数并让代码达到“每行都有中文注释、能直接跑通、能验证结果”的程度。之所以强调这一点是因为网上大量带程序操作录像的语音信号 MATLAB 仿真包要么代码是直接拼接的要么注释少到没法改参数要么录像里演示的数据和脚本里预设的数据根本不是同一份。与其去下载一个来路不明的安装包不如按照工程上最常用的路径自己把分帧、加窗、特征计算、结果绘图这套流程完整实现一遍。这篇文章适合正在做语音处理课设、需要提取特征做分类或端点检测、以及想搞清楚短时分析参数到底怎么影响结果的人。2. 短时分析的三个基础参数能量、过零率、自相关各自能做什么2.1 为什么是这三兄弟一起出现短时能量、短时过零率、短时自相关经常被放在同一个文件里不是因为它们是语音处理教材的第一章而是因为它们恰好从三个互补的角度描述一帧信号能量反映幅度强弱过零率反映频率高低自相关反映周期性强弱。这三者组合起来就能在时域直接完成很多实际任务比如静音段和语音段的区分、清音和浊音的判断、基音周期的粗略估计。在实际使用中三者的适用范围有明确区别。短时能量对幅度变化最敏感环境噪声大时能量曲线会被抬高清音段和噪声段的区分度下降短时过零率对噪声和直流偏移同样敏感但对频率变化的分辨能力强短时自相关则专门用来找周期性但计算量明显高于前两者。因此短时分析的典型做法不是只看某一个参数而是把三个参数画在同一张图上用能量判断起止边界用过零率辅助判断清浊再用短时自相关在浊音段内估计基音周期。2.2 短时能量的定义和窗长的矛盾短时能量的计算逻辑非常直接把一帧信号逐点平方后求和再除以帧长做归一化。用公式表示就是$$E_n \sum_{m0}^{L-1} x_m^2 \cdot w_m$$其中 $L$ 是帧长$w_m$ 是窗函数序列。这里的核心参数是窗长。窗长取得太长时间分辨率就差两个相邻音节之间的能量变化会被抹平窗长取得太短帧内包含的样本点太少能量估计的方差会变得很大曲线抖动明显。工程上最常见的选择是 10 到 30 毫秒语音信号采样率为 16 kHz 时对应 160 到 480 个采样点。为什么是这个量级因为人在正常语速下一个音节的持续时间大约在 80 到 200 毫秒一个音素内部的平稳段约为 20 到 40 毫秒帧长介于两者之间才能既保留音节变化特征又不破坏帧内平稳假设。短时能量的直接用途是端点检测。细看能量曲线会发现语音段的能量明显高于静音段所以常设一个阈值来区分有声和无声。但阈值不能取固定值因为不同录音设备的增益不同同一段语音开头和结尾的音量也不同。一种做法是按整段信号最大能量的比例取阈值常见是最大能量的 10% 到 20% 作为门限同时配合一个最短语音时长限制来避免把短促的噪声尖峰误判成语音。2.3 短时过零率在清浊音区分中的角色过零率统计的是一帧内信号符号变化的次数。清音段比如“s”“f”这类辅音的频谱能量集中在较高频率波形在时间轴上往返穿越零点的次数多浊音段比如元音频谱能量集中在低频波形比较平缓过零率相对低。所以短时过零率常被用来做清浊音分类的辅助特征。计算时要注意一个实际问题数值截断和直流偏移会造成伪过零。如果语音信号在采集时混入了直流分量那么整段波形会被整体抬高或压低原本符号不变化的连续段也可能不断穿越零点导致过零率虚高。在计算过零率之前先对每一帧做去均值处理也就是把该帧的均值减掉再统计符号变化这是我在实际处理中一定会加的一步。过零率的第二个坑是符号函数本身对信号幅值不敏感。即使信号幅度极小只要正负交替过零率照样很高。因此单独的过零率不能区分噪声和语音要与短时能量联合判断。如果能量低但过零率很高大概率是背景噪声如果能量高且过零率低基本可以判定为浊音段。2.4 短时自相关周期性的定量化表达短时自相关计算的是同一帧信号与其延迟版本之间的相似度。语音信号中浊音段具有准周期性自相关函数会在基音周期的整数倍位置出现峰值通过寻找第一个明显峰值的位置就能反过来推算基音频率。数学定义是$$R_n(k) \sum_{m0}^{L-1-k} x_m \cdot x_{mk}$$计算时有几个要点。延迟点数 k 的范围决定了基音周期搜索范围对女声基频大约 200 到 350 Hz对男声大约 80 到 200 Hz采样率 16 kHz 时对应的周期范围分别是约 46 到 200 个采样点。所以延迟 k 的上限至少要覆盖最低基频对应的周期。其次自相关计算量随帧长和延迟范围增加而增加如果对每一帧都做完整循环MATLAB 跑起来会明显变慢。改进方法是先用短时能量判断这一帧是否属于浊音段只有能量高且周期明显的帧才做自相关计算能省掉不少时间。自相关结果中需要注意“半周期峰值”现象。当信号波形不是标准正弦而是接近对称形状时自相关函数可能在半个周期处出现一个幅度相当的峰值导致基音周期估计减半。解决的办法是在搜索峰值时设置一个峰值幅度比例阈值比如只有峰值高于全局最大值 80% 的候选点才算有效或者结合过零率先判断信号周期性程度周期性不强时直接跳过基音估计。3. 用 MATLAB 写出带中文注释的短时分析代码3.1 先确定输入信号和参数采样率、窗型、帧移写代码之前需要先把一组固定参数定义清楚。为了说明问题这里用一段合成语音来演示8 kHz 采样率包含一个 200 Hz 的浊音段和一个 1200 Hz 的清音段中间穿插静音。这样构造的好处是每个参数算出来的理论值都已知后面验证时可直接对照。如果用真实语音也可以但要额外注意录音设备引入的直流偏移这在后面的去均值流程里会专门处理。参数设置遵循语音分析最常用的组合窗长取 25 毫秒帧移取 10 毫秒窗型选汉明窗。8 kHz 采样率下25 毫秒就是 200 个采样点帧移 10 毫秒就是 80 个采样点。这两个参数决定了相邻两帧的重叠率为 60%既保证帧间特征的平滑过渡又不会引入过大的计算冗余。窗型选汉明窗是因为它在频域的主瓣宽度和旁瓣衰减之间平衡较好矩形窗在帧边缘会引起频谱泄漏汉宁窗衰减更彻底但主瓣更宽实际语音端点检测的工程代码里汉明窗出现频率最高。具体代码如下clear; close all; clc; %% 1. 生成测试信号浊音段 清音段 静音段 fs 8000; % 采样率 8 kHz t 0:1/fs:2; % 总时长 2 秒 x 0.6 * sin(2*pi*200*t); % 200 Hz 浊音信号幅度 0.6 % 在 0.6~0.9 秒之间叠加一段 1200 Hz 的类清音信号 idx_v round(0.6*fs):round(0.9*fs); x(idx_v) x(idx_v) 0.3 * sin(2*pi*1200*t(idx_v)); % 保存原始信号后面做能量和过零率分析 x_raw x; %% 2. 短时分析参数设置 win_len round(0.025 * fs); % 窗长 25ms200 点 hop_len round(0.010 * fs); % 帧移 10ms80 点 win hamming(win_len, periodic); % 汉明窗periodic 适合频谱分析这段代码先用正弦信号合成一个 2 秒的测试序列。之所以没有直接读 WAV 文件是因为合成信号的理论特征可精确计算比如 200 Hz 信号在 8 kHz 采样率下的周期是 40 个采样点后面短时自相关算出基音周期后可以直接对照。hamming函数中的periodic选项会在窗尾补一个零值使窗函数在离散傅里叶变换中满足周期延拓条件频谱分析场景下比symmetric更准确。3.2 手动分帧加窗不依赖 Signal Processing ToolboxMATLAB 自带buffer函数和spectrogram函数可以快速分帧加窗但它们的参数语义不直观对新手来说容易搞混。这里采用手动索引分帧逻辑更透明也方便你随时插入调试代码。分帧的核心是计算帧数公式为帧数 1 floor((信号长度 - 窗长) / 帧移)用这个公式能保证最后一帧不会超出信号边界。如果最后一帧剩余样本不足一个窗长常见做法是直接丢弃或者用零填充。端点检测场景下通常选择丢弃因为末尾多出的静音帧对特征无关紧要而在语音合成场景中则用零填充保持总帧数不变。%% 3. 分帧加窗 sig_len length(x_raw); num_frames 1 floor((sig_len - win_len) / hop_len); frames zeros(num_frames, win_len); for i 1:num_frames start_idx (i-1) * hop_len 1; end_idx start_idx win_len - 1; frame x_raw(start_idx:end_idx); frames(i, :) frame .* win; % 逐帧加汉明窗 end代码里frames是一个二维矩阵行数为帧数列数为窗长。这种存储方式的优点是后续计算短时能量、过零率、自相关时可以直接对frames做行运算不需要再写分帧循环。win是窗函数转置让维度能和 1×200 的帧向量做逐元素乘法。如果后面想要换矩形窗或汉宁窗只需替换win hamming(...)这一行。换窗之后过零率的数值会因为窗边缘被压低而略微改变但整体曲线形态不变这是正常现象。3.3 短时能量和短时过零率的 MATLAB 实现短时能量的代码只需要两行先对加窗后的帧内样本求平方再对每一帧的所有样本求和。由于frames已经是矩阵直接沿列方向求和即可。%% 4. 计算短时能量 frame_energy sum(frames .^ 2, 2); % 逐帧平方求和 frame_energy frame_energy / win_len; % 归一化到每样本能量 %% 5. 计算短时过零率先做帧内去均值 zcr zeros(num_frames, 1); for i 1:num_frames frame_centered frames(i, :) - mean(frames(i, :)); % 去直流分量 sign_changes diff(sign(frame_centered)); % 符号变化检测 zcr(i) sum(sign_changes ~ 0) / (2 * win_len); % 除以帧长归一化 end过零率计数用了sign函数将每帧映射成 1 和 -1 组成的序列再用diff求相邻差值。符号相同的位置差值为 0符号变化的位置差值为 ±2因此sign_changes ~ 0的个数就是真正的过零次数。之所以要除以2 * win_len而不是win_len是因为每两个样本之间最多只计一次过零把次数换算成“每样本平均过零率”时分子除以总样本间隔数才是正确的量纲。这个细节在不少参考代码里被忽略直接除win_len的值虽然只会差两倍但如果你用它对清浊音分类设置阈值就可能把阈值范围算错。去均值的步骤看似多余在合成信号场景下影响也不大但换成真实麦克风录音时声卡直流偏移在几百个样本内的均值可能达到信号幅度的 5% 到 10%不去掉的话静音段的过零率可能超过浊音段导致后续阈值判断全部失效。3.4 短时自相关用中心限幅法减少计算量自相关计算的直接实现是双重循环外层遍历延迟 k内层遍历样本求和。这种写法逻辑清晰但效率低对 2 秒信号、帧长 200 点的情况还能接受但换成 10 秒录音后 MATLAB 的运行时间会成倍增加。常用的优化方法是利用 MATLAB 的自带xcorr函数通过快速傅里叶变换来计算互相关对大帧长和大延迟范围都有明显加速效果。%% 6. 计算短时自相关并估计基音周期 max_lag round(fs / 60); % 最低基频 60 Hz对应最大延迟约 133 点 pitch_period zeros(num_frames, 1); for i 1:num_frames frame frames(i, :) - mean(frames(i, :)); [autocorr, lags] xcorr(frame, max_lag, coeff); % 取正半轴部分去掉零延迟点 pos_idx find(lags 0); ac autocorr(pos_idx); ac(1) 0; % 去掉零延迟主峰避免自相关最大值恒为 1 % 搜索第一个大于 0.3 倍主峰最大值的峰 peak_vals ac; threshold max(ac) * 0.3; [~, max_peak_idx] findpeaks(peak_vals, MinPeakHeight, threshold, ... MinPeakDistance, round(fs/400)); if ~isempty(max_peak_idx) pitch_period(i) max_peak_idx(1) - 1; % lags 从 0 开始需要减 1 end end参数max_lag fs/60设定了基音周期搜索的上限60 Hz 是最低基频的常见取值。coeff选项将自相关归一化到 [-1, 1] 范围便于使用统一的阈值判断峰值有效性。ac(1) 0这一行是重要的工程处理因为零延迟处的自相关值恒为 1不剔除它的话findpeaks找到的最大峰永远是延迟零基音周期估计就失去了意义。MinPeakDistance设置为 20 个采样点对应 400 Hz防止在一个基音周期内出现多个紧挨着的假峰。代码中把不做基音估计的帧的pitch_period保留为零这是有意的选择。清音段和静音段的自相关曲线没有明显峰值强制在这些帧上输出一个基音频率反而会干扰后续统计。所以最后的检索条件是“帧内存在超过阈值的峰才算有效”算是一种保守的做法。3.5 参数速查表与调用方式把上述代码封装成函数后调用方式变成[energy, zcr, pitch, times] short_time_analysis(x_raw, fs, ... WindowLen, 0.025, HopLen, 0.010, WindowType, hamming);函数内部的参数对应关系整理如下表参数名常用取值对应效果采样率 fs8000 / 16000 Hz决定窗长换算、延迟搜索范围窗长 WindowLen20–30 ms过短则能量抖动大过长则时间分辨率差帧移 HopLen5–15 ms越小曲线越平滑计算量越大窗型hamming / hann / rectwin影响旁瓣泄漏和过零率数值最大延迟 max_lagfs/60小于真实基音周期时峰值搜索失败峰值阈值0.3–0.5过小产生半周期假峰过大漏检弱浊音段这组参数是语音信号短时分析中最常见的一组适合端点检测、清浊音区分、基音周期估计等场景。如果是音乐信号分析窗长可以增大到 50 毫秒因为乐音的准平稳段比语音更长如果是爆破音或瞬时噪声分析窗长要缩短到 5 毫秒左右否则突发能量会被相邻帧平均掉。4. 用合成信号和操作录像思路验证仿真结果4.1 三种信号的对照实验设计上一节生成了包含浊音、清音、静音的 2 秒测试信号。现在用这段信号跑一遍完整分析并检查三个参数的输出是否与理论一致。期望结果如下短时能量在 0 到 0.6 秒和 0.9 到 2 秒区间稳定在 0.18 附近0.6 的平方除以 2在 0.6 到 0.9 秒区间增加到约 0.2250.6 和 0.3 两个分量平方和的均值静音段几乎为 0短时过零率在浊音段约为 0.05200 Hz 正弦每 40 个采样点两次过零折算每样本 0.05清音段明显增大短时自相关在浊音段能检测到 40 个采样点的峰值。%% 7. 绘制分析结果 frame_times (0:num_frames-1) * hop_len / fs; figure(Name, 短时分析结果, Position, [100 100 1000 800]); subplot(4,1,1); plot(t, x_raw); title(原始语音信号); ylabel(幅度); xlim([0 2]); subplot(4,1,2); plot(frame_times, frame_energy); title(短时能量); ylabel(E_n); xlabel(时间/s); xlim([0 2]); subplot(4,1,3); plot(frame_times, zcr); title(短时过零率); ylabel(ZCR); xlabel(时间/s); xlim([0 2]); subplot(4,1,4); plot(frame_times, pitch_period); title(短时自相关基音周期); xlabel(时间/s); ylabel(周期/采样点); ylim([0 80]); xlim([0 2]);运行后查看图像重点看两点一是过零率曲线在 0.6 到 0.9 秒区间的数值是否明显高于两侧二是基音周期曲线在 0.2 到 0.6 秒区间是否稳定在 40 附近。如果基音周期曲线上出现接近 80 的值且持续多个帧说明峰值搜索时跳过了第一峰值而选中了第二峰值如果出现接近 20 的值说明遇到半周期假峰问题需要调大MinPeakDistance或提高MinPeakHeight比例。4.2 查找包中自带的仿真波形与运行代码对于一份打包好的语音信号短时分析 MATLAB 仿真资源常见的目录结构是根目录下有一个主脚本文件比如short_time_analysis_demo.m一个或多个音频文件以及一个录像或操作演示目录。拿到资源后做的第一件事不是运行而是先核对主脚本里的音频文件名和音频文件列表是否一致。很多时候录像里演示的是speech.wav脚本里却写的是voice.wav运行必然报错。如果脚本运行后报“文件不存在”优先检查 MATLAB 当前路径。MATLAB 不会自动切换到脚本所在目录直接双击脚本运行时当前工作目录往往是之前打开过的路径。用cd(fileparts(mfilename(fullpath)))加在脚本开头或者在 MATLAB 编辑器的“运行”按钮下拉菜单中选择“运行并更改文件夹”就能避免路径问题。运行成功后对比脚本输出和录像中的曲线形状。如果波形整体一致但数值相差一个数量级先检查输入音频是否做了归一化如果曲线形状完全不同优先怀疑音频文件被替换或采样率不匹配。所有使用 8 kHz 采样率采样的 WAV 文件用 16 kHz 读入再做短时分析时间轴和频率轴都会正确一半错位这是最常见的隐蔽错误。4.3 替换成真实语音观察参数值的合理范围合成信号验证通过后把输入换成一段真实录音。需要注意的是真实语音的平均能量远低于合成信号因为人说话的动态范围大响度高的部分集中在元音段辅音和吸气声能量很低。调整方法有两种一种是先求整段信号的均方根值然后把信号乘以一个缩放系数使均方根值保持在 0.1 到 0.3 之间另一种是直接对每一帧做归一化但这样会破坏帧间的能量对比关系不适合端点检测场景。真实语音的过零率会比合成信号整体偏高原因是背景噪声的过零率通常比浊音段高。如果静音段的过零率与语音段的过零率区间重叠明显可以设置一个过零率的阈值下限比如只有能量超过最大能量的 15% 且过零率低于 0.15 时才判定为浊音段。这个组合判据在大多数安静环境下都能做到 90% 以上的端点识别准确率。在真实语音上运行前把frames的存储方式快速检查一遍。MATLAB 中如果输入音频是行向量frames的维度是帧数 × 帧长如果输入是列向量在分帧时出现了start_idx:end_idx不够长度的情况。稳妥做法是在读取音频后统一转成行向量if size(x_raw, 1) 1 x_raw x_raw; end这段代码放在读取音频之后、分帧之前可以避免 80% 由维度问题引起的报错。4.4 操作录像演示的常见惯用流程一套完整的语音信号短时分析 MATLAB 仿真操作录像通常包含四步运行主脚本、显示分帧示意、展示三条曲线、调整参数对比。注意录像的人会在运行前先用whos查看变量列表确认fs、x_raw、frames三个关键变量存在。你也可以用同样的方式自查whos fs x_raw frames frame_energy zcr pitch_period如果frames变量不存在而主脚本能正常出图说明脚本用的是实时分帧和逐帧绘图此时要重点检查绘图循环中frame变量是否在每次迭代后清空。如果把frame用来存储当前帧数据画能量图时又在同一循环里画了上一帧的数据曲线的第一个点和最后一个点就会出现明显偏差。5. 短时分析结果的验证方法从曲线形状反推帧长和窗型5.1 通过曲线形状判断帧长和帧移是否匹配有时拿到的代码不提供操作录像只有输出图。这时可以根据曲线的锯齿程度反推帧长和帧移是否合理。短时能量曲线相邻帧之间如果出现剧烈的上下跳动幅度超过整段曲线量程的 30%优先怀疑帧移过大导致帧间重叠率不足。把帧移从 10 毫秒改成 5 毫秒后曲线会明显平滑因为相邻两帧共享了 75% 的数据点能量值自然连续。反过来如果能量曲线的上升沿和下降沿被拉得很宽原本从静音到语音的变化过程跨了 4 到 5 帧说明窗长太长。此时压缩窗长到 10 毫秒曲线会在 1 到 2 帧内完成跳变端点检测用阈值定位时的误差也会随之缩小。经验法则是窗长与帧移的比值控制在 2 到 3 之间最合适对应 60% 到 70% 的重叠率低于 2 则曲线毛糙高于 3 则相邻帧几乎等于没移动等于浪费算力。5.2 用时间分辨率检查基音周期估计是否出错短时自相关输出的基音周期曲线是否可信同样可以通过肉眼判断。如果基音周期在相邻帧之间跳变超过 10 个采样点而实际说话者的发音在这 10 毫秒内不可能出现这么大的频率跳跃说明峰值检测选错了峰。典型的错误特征是曲线在一串 40 左右的值中插入几个 20 附近的点这是半个周期被误检的典型形态如果插入点是 80 附近的点则说明某个帧的起始相位导致第一周期被窗函数严重衰减搜到的第一个有效峰已经是第二周期。处理这类问题的参数优先级是先调MinPeakDistance到 20 个采样点再调阈值比例到 0.4。不要一开始就把阈值调高因为弱浊音段的周期峰值本身就只有最大值的 40% 到 60%调太高会把有效帧全部滤掉。另外一个隐含条件是findpeaks默认只检测峰值点如果某帧的自相关曲线在基音周期附近呈现一个宽平台而没有明显的单点尖峰findpeaks可能识别不出有效峰此时基音周期留空是合理的。5.3 峰值距离换算成频率的快速检查验证基音周期数值是否合理可以用一行命令快速换算pitch_hz fs ./ pitch_period(pitch_period 0); fprintf(基音频率范围: %.2f Hz ~ %.2f Hz\n, min(pitch_hz), max(pitch_hz));对男声录音输出范围应在 80 到 250 Hz对女声录音应在 150 到 350 Hz。如果最低值低于 60 Hz 或最高值超过 500 Hz都要怀疑是否出现了峰值选择的系统性偏移。实际操作中还有一种常见情况无声段被检测出 100 Hz 左右的“基音”这是因为背景噪声中低频部分的自相关虽然弱但仍能超过较低阈值。处理办法是增加一个能量条件只有当frame_energy大于全局能量最大值的 1/5 时才把基音周期写入结果数组否则保持为零。5.4 短时自相关与 FFT 频域对照法最后介绍一个排查短时分析实现错误的通用技巧把自相关峰值对应的频率和 FFT 频谱峰值对应的频率做对照。对浊音帧通过 R 得到的基音频率应与通过频谱分析得到的第一个明显谱峰一致。MATLAB 中可以这样验证i find(pitch_period 0, 1, first); N length(frames(i, :)); spec abs(fft(frames(i, :))); f (0:N/2-1) * fs / N; [~, spec_peak_idx] max(spec(2:round(N/2))); fprintf(自相关基音频率: %.2f Hz, FFT峰值频率: %.2f Hz\n, ... fs / pitch_period(i), f(spec_peak_idx));如果自相关测出的频率和 FFT 峰值频率相差很大问题大概率不在后期搜索参数上而是窗长选得太短导致频域分辨率不足以分辨基音频率。8 kHz 采样率、200 点窗长的频域分辨率是 40 Hz200 Hz 和 240 Hz 的基音差异在频谱图上无法区分但自相关在时域上可以分辨到 1 个采样点。所以自相关适合精细测量基音周期FFT 适合快速观察频谱包络两者互相校验比只看一个指标可靠得多。本文还有配套的精品资源点击获取 MATLAB语音短时分析实战：能量、过零率与自相关提取
简介面向语音信号处理初学者的MATLAB仿真资源围绕语音信号短时分析中的短时过零率、短时能量、短时自相关等经典特征提取方法展开适合通信工程、电子信息等专业学生进行课程设计或实验入门。资源基于MATLAB 2022A版本编写包内共含10个文件以8个.m源码文件为主体配套1个wav语音样本用于实验测试以及1个avi操作录像演示具体运行流程压缩包大小仅590KB轻量便捷。已有617人学习下载程序代码均配有中文注释关键函数如分帧、加窗、特征计算等均单独封装便于读者逐段理解算法实现同时录像中特别强调了MATLAB当前文件夹路径的设置要点可有效避免因路径问题导致的运行报错。整体结构清晰能帮助使用者快速掌握语音短时分析的编程实现与调试思路也可作为相关课程实验报告的参考素材。1. 从一段语音到“帧”短时分析为什么先把信号切碎接手语音信号处理的第一周很多人都会卡在同一个地方一段 5 秒的语音直接用傅里叶变换做频谱分析得到的是一条几乎看不出规律的曲线加窗、改点数、换函数都救不回来。原因不在算法而在你看待信号的方式。语音是典型的非平稳信号发音时声带振动频率、声道形状、音量都在变化用全局视角去分析整体特征结果自然发散。短时分析的核心做法是先承认这种非平稳性再把信号切成 20 到 30 毫秒的短帧在每一帧内默认信号近似平稳然后逐帧计算短时能量、短时过零率、短时自相关等特征将这些特征值连成曲线整段语音的变化规律才显现出来。这篇文章要解决的问题很具体如何用 MATLAB 从零搭建一套短时分析仿真覆盖短时能量、短时过零率、短时自相关三个经典参数并让代码达到“每行都有中文注释、能直接跑通、能验证结果”的程度。之所以强调这一点是因为网上大量带程序操作录像的语音信号 MATLAB 仿真包要么代码是直接拼接的要么注释少到没法改参数要么录像里演示的数据和脚本里预设的数据根本不是同一份。与其去下载一个来路不明的安装包不如按照工程上最常用的路径自己把分帧、加窗、特征计算、结果绘图这套流程完整实现一遍。这篇文章适合正在做语音处理课设、需要提取特征做分类或端点检测、以及想搞清楚短时分析参数到底怎么影响结果的人。2. 短时分析的三个基础参数能量、过零率、自相关各自能做什么2.1 为什么是这三兄弟一起出现短时能量、短时过零率、短时自相关经常被放在同一个文件里不是因为它们是语音处理教材的第一章而是因为它们恰好从三个互补的角度描述一帧信号能量反映幅度强弱过零率反映频率高低自相关反映周期性强弱。这三者组合起来就能在时域直接完成很多实际任务比如静音段和语音段的区分、清音和浊音的判断、基音周期的粗略估计。在实际使用中三者的适用范围有明确区别。短时能量对幅度变化最敏感环境噪声大时能量曲线会被抬高清音段和噪声段的区分度下降短时过零率对噪声和直流偏移同样敏感但对频率变化的分辨能力强短时自相关则专门用来找周期性但计算量明显高于前两者。因此短时分析的典型做法不是只看某一个参数而是把三个参数画在同一张图上用能量判断起止边界用过零率辅助判断清浊再用短时自相关在浊音段内估计基音周期。2.2 短时能量的定义和窗长的矛盾短时能量的计算逻辑非常直接把一帧信号逐点平方后求和再除以帧长做归一化。用公式表示就是$$E_n \sum_{m0}^{L-1} x_m^2 \cdot w_m$$其中 $L$ 是帧长$w_m$ 是窗函数序列。这里的核心参数是窗长。窗长取得太长时间分辨率就差两个相邻音节之间的能量变化会被抹平窗长取得太短帧内包含的样本点太少能量估计的方差会变得很大曲线抖动明显。工程上最常见的选择是 10 到 30 毫秒语音信号采样率为 16 kHz 时对应 160 到 480 个采样点。为什么是这个量级因为人在正常语速下一个音节的持续时间大约在 80 到 200 毫秒一个音素内部的平稳段约为 20 到 40 毫秒帧长介于两者之间才能既保留音节变化特征又不破坏帧内平稳假设。短时能量的直接用途是端点检测。细看能量曲线会发现语音段的能量明显高于静音段所以常设一个阈值来区分有声和无声。但阈值不能取固定值因为不同录音设备的增益不同同一段语音开头和结尾的音量也不同。一种做法是按整段信号最大能量的比例取阈值常见是最大能量的 10% 到 20% 作为门限同时配合一个最短语音时长限制来避免把短促的噪声尖峰误判成语音。2.3 短时过零率在清浊音区分中的角色过零率统计的是一帧内信号符号变化的次数。清音段比如“s”“f”这类辅音的频谱能量集中在较高频率波形在时间轴上往返穿越零点的次数多浊音段比如元音频谱能量集中在低频波形比较平缓过零率相对低。所以短时过零率常被用来做清浊音分类的辅助特征。计算时要注意一个实际问题数值截断和直流偏移会造成伪过零。如果语音信号在采集时混入了直流分量那么整段波形会被整体抬高或压低原本符号不变化的连续段也可能不断穿越零点导致过零率虚高。在计算过零率之前先对每一帧做去均值处理也就是把该帧的均值减掉再统计符号变化这是我在实际处理中一定会加的一步。过零率的第二个坑是符号函数本身对信号幅值不敏感。即使信号幅度极小只要正负交替过零率照样很高。因此单独的过零率不能区分噪声和语音要与短时能量联合判断。如果能量低但过零率很高大概率是背景噪声如果能量高且过零率低基本可以判定为浊音段。2.4 短时自相关周期性的定量化表达短时自相关计算的是同一帧信号与其延迟版本之间的相似度。语音信号中浊音段具有准周期性自相关函数会在基音周期的整数倍位置出现峰值通过寻找第一个明显峰值的位置就能反过来推算基音频率。数学定义是$$R_n(k) \sum_{m0}^{L-1-k} x_m \cdot x_{mk}$$计算时有几个要点。延迟点数 k 的范围决定了基音周期搜索范围对女声基频大约 200 到 350 Hz对男声大约 80 到 200 Hz采样率 16 kHz 时对应的周期范围分别是约 46 到 200 个采样点。所以延迟 k 的上限至少要覆盖最低基频对应的周期。其次自相关计算量随帧长和延迟范围增加而增加如果对每一帧都做完整循环MATLAB 跑起来会明显变慢。改进方法是先用短时能量判断这一帧是否属于浊音段只有能量高且周期明显的帧才做自相关计算能省掉不少时间。自相关结果中需要注意“半周期峰值”现象。当信号波形不是标准正弦而是接近对称形状时自相关函数可能在半个周期处出现一个幅度相当的峰值导致基音周期估计减半。解决的办法是在搜索峰值时设置一个峰值幅度比例阈值比如只有峰值高于全局最大值 80% 的候选点才算有效或者结合过零率先判断信号周期性程度周期性不强时直接跳过基音估计。3. 用 MATLAB 写出带中文注释的短时分析代码3.1 先确定输入信号和参数采样率、窗型、帧移写代码之前需要先把一组固定参数定义清楚。为了说明问题这里用一段合成语音来演示8 kHz 采样率包含一个 200 Hz 的浊音段和一个 1200 Hz 的清音段中间穿插静音。这样构造的好处是每个参数算出来的理论值都已知后面验证时可直接对照。如果用真实语音也可以但要额外注意录音设备引入的直流偏移这在后面的去均值流程里会专门处理。参数设置遵循语音分析最常用的组合窗长取 25 毫秒帧移取 10 毫秒窗型选汉明窗。8 kHz 采样率下25 毫秒就是 200 个采样点帧移 10 毫秒就是 80 个采样点。这两个参数决定了相邻两帧的重叠率为 60%既保证帧间特征的平滑过渡又不会引入过大的计算冗余。窗型选汉明窗是因为它在频域的主瓣宽度和旁瓣衰减之间平衡较好矩形窗在帧边缘会引起频谱泄漏汉宁窗衰减更彻底但主瓣更宽实际语音端点检测的工程代码里汉明窗出现频率最高。具体代码如下clear; close all; clc; %% 1. 生成测试信号浊音段 清音段 静音段 fs 8000; % 采样率 8 kHz t 0:1/fs:2; % 总时长 2 秒 x 0.6 * sin(2*pi*200*t); % 200 Hz 浊音信号幅度 0.6 % 在 0.6~0.9 秒之间叠加一段 1200 Hz 的类清音信号 idx_v round(0.6*fs):round(0.9*fs); x(idx_v) x(idx_v) 0.3 * sin(2*pi*1200*t(idx_v)); % 保存原始信号后面做能量和过零率分析 x_raw x; %% 2. 短时分析参数设置 win_len round(0.025 * fs); % 窗长 25ms200 点 hop_len round(0.010 * fs); % 帧移 10ms80 点 win hamming(win_len, periodic); % 汉明窗periodic 适合频谱分析这段代码先用正弦信号合成一个 2 秒的测试序列。之所以没有直接读 WAV 文件是因为合成信号的理论特征可精确计算比如 200 Hz 信号在 8 kHz 采样率下的周期是 40 个采样点后面短时自相关算出基音周期后可以直接对照。hamming函数中的periodic选项会在窗尾补一个零值使窗函数在离散傅里叶变换中满足周期延拓条件频谱分析场景下比symmetric更准确。3.2 手动分帧加窗不依赖 Signal Processing ToolboxMATLAB 自带buffer函数和spectrogram函数可以快速分帧加窗但它们的参数语义不直观对新手来说容易搞混。这里采用手动索引分帧逻辑更透明也方便你随时插入调试代码。分帧的核心是计算帧数公式为帧数 1 floor((信号长度 - 窗长) / 帧移)用这个公式能保证最后一帧不会超出信号边界。如果最后一帧剩余样本不足一个窗长常见做法是直接丢弃或者用零填充。端点检测场景下通常选择丢弃因为末尾多出的静音帧对特征无关紧要而在语音合成场景中则用零填充保持总帧数不变。%% 3. 分帧加窗 sig_len length(x_raw); num_frames 1 floor((sig_len - win_len) / hop_len); frames zeros(num_frames, win_len); for i 1:num_frames start_idx (i-1) * hop_len 1; end_idx start_idx win_len - 1; frame x_raw(start_idx:end_idx); frames(i, :) frame .* win; % 逐帧加汉明窗 end代码里frames是一个二维矩阵行数为帧数列数为窗长。这种存储方式的优点是后续计算短时能量、过零率、自相关时可以直接对frames做行运算不需要再写分帧循环。win是窗函数转置让维度能和 1×200 的帧向量做逐元素乘法。如果后面想要换矩形窗或汉宁窗只需替换win hamming(...)这一行。换窗之后过零率的数值会因为窗边缘被压低而略微改变但整体曲线形态不变这是正常现象。3.3 短时能量和短时过零率的 MATLAB 实现短时能量的代码只需要两行先对加窗后的帧内样本求平方再对每一帧的所有样本求和。由于frames已经是矩阵直接沿列方向求和即可。%% 4. 计算短时能量 frame_energy sum(frames .^ 2, 2); % 逐帧平方求和 frame_energy frame_energy / win_len; % 归一化到每样本能量 %% 5. 计算短时过零率先做帧内去均值 zcr zeros(num_frames, 1); for i 1:num_frames frame_centered frames(i, :) - mean(frames(i, :)); % 去直流分量 sign_changes diff(sign(frame_centered)); % 符号变化检测 zcr(i) sum(sign_changes ~ 0) / (2 * win_len); % 除以帧长归一化 end过零率计数用了sign函数将每帧映射成 1 和 -1 组成的序列再用diff求相邻差值。符号相同的位置差值为 0符号变化的位置差值为 ±2因此sign_changes ~ 0的个数就是真正的过零次数。之所以要除以2 * win_len而不是win_len是因为每两个样本之间最多只计一次过零把次数换算成“每样本平均过零率”时分子除以总样本间隔数才是正确的量纲。这个细节在不少参考代码里被忽略直接除win_len的值虽然只会差两倍但如果你用它对清浊音分类设置阈值就可能把阈值范围算错。去均值的步骤看似多余在合成信号场景下影响也不大但换成真实麦克风录音时声卡直流偏移在几百个样本内的均值可能达到信号幅度的 5% 到 10%不去掉的话静音段的过零率可能超过浊音段导致后续阈值判断全部失效。3.4 短时自相关用中心限幅法减少计算量自相关计算的直接实现是双重循环外层遍历延迟 k内层遍历样本求和。这种写法逻辑清晰但效率低对 2 秒信号、帧长 200 点的情况还能接受但换成 10 秒录音后 MATLAB 的运行时间会成倍增加。常用的优化方法是利用 MATLAB 的自带xcorr函数通过快速傅里叶变换来计算互相关对大帧长和大延迟范围都有明显加速效果。%% 6. 计算短时自相关并估计基音周期 max_lag round(fs / 60); % 最低基频 60 Hz对应最大延迟约 133 点 pitch_period zeros(num_frames, 1); for i 1:num_frames frame frames(i, :) - mean(frames(i, :)); [autocorr, lags] xcorr(frame, max_lag, coeff); % 取正半轴部分去掉零延迟点 pos_idx find(lags 0); ac autocorr(pos_idx); ac(1) 0; % 去掉零延迟主峰避免自相关最大值恒为 1 % 搜索第一个大于 0.3 倍主峰最大值的峰 peak_vals ac; threshold max(ac) * 0.3; [~, max_peak_idx] findpeaks(peak_vals, MinPeakHeight, threshold, ... MinPeakDistance, round(fs/400)); if ~isempty(max_peak_idx) pitch_period(i) max_peak_idx(1) - 1; % lags 从 0 开始需要减 1 end end参数max_lag fs/60设定了基音周期搜索的上限60 Hz 是最低基频的常见取值。coeff选项将自相关归一化到 [-1, 1] 范围便于使用统一的阈值判断峰值有效性。ac(1) 0这一行是重要的工程处理因为零延迟处的自相关值恒为 1不剔除它的话findpeaks找到的最大峰永远是延迟零基音周期估计就失去了意义。MinPeakDistance设置为 20 个采样点对应 400 Hz防止在一个基音周期内出现多个紧挨着的假峰。代码中把不做基音估计的帧的pitch_period保留为零这是有意的选择。清音段和静音段的自相关曲线没有明显峰值强制在这些帧上输出一个基音频率反而会干扰后续统计。所以最后的检索条件是“帧内存在超过阈值的峰才算有效”算是一种保守的做法。3.5 参数速查表与调用方式把上述代码封装成函数后调用方式变成[energy, zcr, pitch, times] short_time_analysis(x_raw, fs, ... WindowLen, 0.025, HopLen, 0.010, WindowType, hamming);函数内部的参数对应关系整理如下表参数名常用取值对应效果采样率 fs8000 / 16000 Hz决定窗长换算、延迟搜索范围窗长 WindowLen20–30 ms过短则能量抖动大过长则时间分辨率差帧移 HopLen5–15 ms越小曲线越平滑计算量越大窗型hamming / hann / rectwin影响旁瓣泄漏和过零率数值最大延迟 max_lagfs/60小于真实基音周期时峰值搜索失败峰值阈值0.3–0.5过小产生半周期假峰过大漏检弱浊音段这组参数是语音信号短时分析中最常见的一组适合端点检测、清浊音区分、基音周期估计等场景。如果是音乐信号分析窗长可以增大到 50 毫秒因为乐音的准平稳段比语音更长如果是爆破音或瞬时噪声分析窗长要缩短到 5 毫秒左右否则突发能量会被相邻帧平均掉。4. 用合成信号和操作录像思路验证仿真结果4.1 三种信号的对照实验设计上一节生成了包含浊音、清音、静音的 2 秒测试信号。现在用这段信号跑一遍完整分析并检查三个参数的输出是否与理论一致。期望结果如下短时能量在 0 到 0.6 秒和 0.9 到 2 秒区间稳定在 0.18 附近0.6 的平方除以 2在 0.6 到 0.9 秒区间增加到约 0.2250.6 和 0.3 两个分量平方和的均值静音段几乎为 0短时过零率在浊音段约为 0.05200 Hz 正弦每 40 个采样点两次过零折算每样本 0.05清音段明显增大短时自相关在浊音段能检测到 40 个采样点的峰值。%% 7. 绘制分析结果 frame_times (0:num_frames-1) * hop_len / fs; figure(Name, 短时分析结果, Position, [100 100 1000 800]); subplot(4,1,1); plot(t, x_raw); title(原始语音信号); ylabel(幅度); xlim([0 2]); subplot(4,1,2); plot(frame_times, frame_energy); title(短时能量); ylabel(E_n); xlabel(时间/s); xlim([0 2]); subplot(4,1,3); plot(frame_times, zcr); title(短时过零率); ylabel(ZCR); xlabel(时间/s); xlim([0 2]); subplot(4,1,4); plot(frame_times, pitch_period); title(短时自相关基音周期); xlabel(时间/s); ylabel(周期/采样点); ylim([0 80]); xlim([0 2]);运行后查看图像重点看两点一是过零率曲线在 0.6 到 0.9 秒区间的数值是否明显高于两侧二是基音周期曲线在 0.2 到 0.6 秒区间是否稳定在 40 附近。如果基音周期曲线上出现接近 80 的值且持续多个帧说明峰值搜索时跳过了第一峰值而选中了第二峰值如果出现接近 20 的值说明遇到半周期假峰问题需要调大MinPeakDistance或提高MinPeakHeight比例。4.2 查找包中自带的仿真波形与运行代码对于一份打包好的语音信号短时分析 MATLAB 仿真资源常见的目录结构是根目录下有一个主脚本文件比如short_time_analysis_demo.m一个或多个音频文件以及一个录像或操作演示目录。拿到资源后做的第一件事不是运行而是先核对主脚本里的音频文件名和音频文件列表是否一致。很多时候录像里演示的是speech.wav脚本里却写的是voice.wav运行必然报错。如果脚本运行后报“文件不存在”优先检查 MATLAB 当前路径。MATLAB 不会自动切换到脚本所在目录直接双击脚本运行时当前工作目录往往是之前打开过的路径。用cd(fileparts(mfilename(fullpath)))加在脚本开头或者在 MATLAB 编辑器的“运行”按钮下拉菜单中选择“运行并更改文件夹”就能避免路径问题。运行成功后对比脚本输出和录像中的曲线形状。如果波形整体一致但数值相差一个数量级先检查输入音频是否做了归一化如果曲线形状完全不同优先怀疑音频文件被替换或采样率不匹配。所有使用 8 kHz 采样率采样的 WAV 文件用 16 kHz 读入再做短时分析时间轴和频率轴都会正确一半错位这是最常见的隐蔽错误。4.3 替换成真实语音观察参数值的合理范围合成信号验证通过后把输入换成一段真实录音。需要注意的是真实语音的平均能量远低于合成信号因为人说话的动态范围大响度高的部分集中在元音段辅音和吸气声能量很低。调整方法有两种一种是先求整段信号的均方根值然后把信号乘以一个缩放系数使均方根值保持在 0.1 到 0.3 之间另一种是直接对每一帧做归一化但这样会破坏帧间的能量对比关系不适合端点检测场景。真实语音的过零率会比合成信号整体偏高原因是背景噪声的过零率通常比浊音段高。如果静音段的过零率与语音段的过零率区间重叠明显可以设置一个过零率的阈值下限比如只有能量超过最大能量的 15% 且过零率低于 0.15 时才判定为浊音段。这个组合判据在大多数安静环境下都能做到 90% 以上的端点识别准确率。在真实语音上运行前把frames的存储方式快速检查一遍。MATLAB 中如果输入音频是行向量frames的维度是帧数 × 帧长如果输入是列向量在分帧时出现了start_idx:end_idx不够长度的情况。稳妥做法是在读取音频后统一转成行向量if size(x_raw, 1) 1 x_raw x_raw; end这段代码放在读取音频之后、分帧之前可以避免 80% 由维度问题引起的报错。4.4 操作录像演示的常见惯用流程一套完整的语音信号短时分析 MATLAB 仿真操作录像通常包含四步运行主脚本、显示分帧示意、展示三条曲线、调整参数对比。注意录像的人会在运行前先用whos查看变量列表确认fs、x_raw、frames三个关键变量存在。你也可以用同样的方式自查whos fs x_raw frames frame_energy zcr pitch_period如果frames变量不存在而主脚本能正常出图说明脚本用的是实时分帧和逐帧绘图此时要重点检查绘图循环中frame变量是否在每次迭代后清空。如果把frame用来存储当前帧数据画能量图时又在同一循环里画了上一帧的数据曲线的第一个点和最后一个点就会出现明显偏差。5. 短时分析结果的验证方法从曲线形状反推帧长和窗型5.1 通过曲线形状判断帧长和帧移是否匹配有时拿到的代码不提供操作录像只有输出图。这时可以根据曲线的锯齿程度反推帧长和帧移是否合理。短时能量曲线相邻帧之间如果出现剧烈的上下跳动幅度超过整段曲线量程的 30%优先怀疑帧移过大导致帧间重叠率不足。把帧移从 10 毫秒改成 5 毫秒后曲线会明显平滑因为相邻两帧共享了 75% 的数据点能量值自然连续。反过来如果能量曲线的上升沿和下降沿被拉得很宽原本从静音到语音的变化过程跨了 4 到 5 帧说明窗长太长。此时压缩窗长到 10 毫秒曲线会在 1 到 2 帧内完成跳变端点检测用阈值定位时的误差也会随之缩小。经验法则是窗长与帧移的比值控制在 2 到 3 之间最合适对应 60% 到 70% 的重叠率低于 2 则曲线毛糙高于 3 则相邻帧几乎等于没移动等于浪费算力。5.2 用时间分辨率检查基音周期估计是否出错短时自相关输出的基音周期曲线是否可信同样可以通过肉眼判断。如果基音周期在相邻帧之间跳变超过 10 个采样点而实际说话者的发音在这 10 毫秒内不可能出现这么大的频率跳跃说明峰值检测选错了峰。典型的错误特征是曲线在一串 40 左右的值中插入几个 20 附近的点这是半个周期被误检的典型形态如果插入点是 80 附近的点则说明某个帧的起始相位导致第一周期被窗函数严重衰减搜到的第一个有效峰已经是第二周期。处理这类问题的参数优先级是先调MinPeakDistance到 20 个采样点再调阈值比例到 0.4。不要一开始就把阈值调高因为弱浊音段的周期峰值本身就只有最大值的 40% 到 60%调太高会把有效帧全部滤掉。另外一个隐含条件是findpeaks默认只检测峰值点如果某帧的自相关曲线在基音周期附近呈现一个宽平台而没有明显的单点尖峰findpeaks可能识别不出有效峰此时基音周期留空是合理的。5.3 峰值距离换算成频率的快速检查验证基音周期数值是否合理可以用一行命令快速换算pitch_hz fs ./ pitch_period(pitch_period 0); fprintf(基音频率范围: %.2f Hz ~ %.2f Hz\n, min(pitch_hz), max(pitch_hz));对男声录音输出范围应在 80 到 250 Hz对女声录音应在 150 到 350 Hz。如果最低值低于 60 Hz 或最高值超过 500 Hz都要怀疑是否出现了峰值选择的系统性偏移。实际操作中还有一种常见情况无声段被检测出 100 Hz 左右的“基音”这是因为背景噪声中低频部分的自相关虽然弱但仍能超过较低阈值。处理办法是增加一个能量条件只有当frame_energy大于全局能量最大值的 1/5 时才把基音周期写入结果数组否则保持为零。5.4 短时自相关与 FFT 频域对照法最后介绍一个排查短时分析实现错误的通用技巧把自相关峰值对应的频率和 FFT 频谱峰值对应的频率做对照。对浊音帧通过 R 得到的基音频率应与通过频谱分析得到的第一个明显谱峰一致。MATLAB 中可以这样验证i find(pitch_period 0, 1, first); N length(frames(i, :)); spec abs(fft(frames(i, :))); f (0:N/2-1) * fs / N; [~, spec_peak_idx] max(spec(2:round(N/2))); fprintf(自相关基音频率: %.2f Hz, FFT峰值频率: %.2f Hz\n, ... fs / pitch_period(i), f(spec_peak_idx));如果自相关测出的频率和 FFT 峰值频率相差很大问题大概率不在后期搜索参数上而是窗长选得太短导致频域分辨率不足以分辨基音频率。8 kHz 采样率、200 点窗长的频域分辨率是 40 Hz200 Hz 和 240 Hz 的基音差异在频谱图上无法区分但自相关在时域上可以分辨到 1 个采样点。所以自相关适合精细测量基音周期FFT 适合快速观察频谱包络两者互相校验比只看一个指标可靠得多。本文还有配套的精品资源点击获取