语音信号短时分析Matlab仿真：能量、过零率与基频估计

简介面向语音信号处理初学者与电子信息类专业学生这份Matlab仿真资源围绕语音信号的短时分析展开完整覆盖短时能量、短时过零率、短时自相关等经典特征提取方法也涉及短时幅度差、短时平均幅度等处理可帮助理解语音分帧、加窗与特征计算的核心流程。包内共10个文件以8个m脚本为主干包含主程序与分帧、特征提取等模块另附1个wav测试语音和1个avi操作录像压缩包仅约590KB结构精简便于对照源码逐段学习。当前已有617人浏览学习。代码均提供中文注释录像演示了在Matlab 2022A中设置当前文件夹路径并运行得到结果的完整过程能有效降低上手门槛尤其适合课程设计、实验报告或自学语音信号短时分析技术时参考使用。1. 语音信号短时分析matlab仿真要解决的核心问题语音信号的短时分析是MATLAB仿真里出现频率最高的题目之一从wav文件读入一段说话声依次计算短时能量、短时过零率、短时自相关再把三条特征曲线画在波形下方。它经典是因为三个特征恰好回答了语音处理入门时的三个问题语音从哪里开始和结束、每一帧是清音还是浊音、浊音部分的基频是多少。整个仿真做下来得到的不是孤立算法而是一套能直接接进端点检测和基频估计的流程。宏观上看一段语音由声带振动、声道摩擦、静音三部分组成波形特性在几十毫秒内就会切换所以它是不平稳的但把窗缩短到10到30毫秒每个片段内部又能看作相对平稳。短时分析正是基于这个前提把波形切成帧逐帧计算特征。仿真难度不在公式本身而在帧长换算、窗型选择、自相关搜索范围和坐标轴对齐这四件事上它们正好对应后面每章要处理的难点。2. 语音信号短时分析第一步分帧加窗与短时能量计算2.1 语音为什么不是平稳信号一段包含几个字的语音如果直接对整段做FFT得到的是把元音、辅音、静音全部混在一起的频谱既看不出字与字的边界也估不准基频。语音的生成机理决定了它不是平稳过程声带振动与声道的张合在几十毫秒内改变浊音段有准周期性、清音段接近噪声、静音段则几乎为零。宏观上不平稳微观上又存在相对稳定的时间段这正是短时分析成立的前提。但“相对稳定”要有一个尺度。人说话语速不同音节过渡快慢也不同用1毫秒的窗去看窗内采样点太少无法体现周期性用500毫秒的窗去看把多个音素混在一起又回到了整段FFT的老问题。工程上把稳定性尺度约定在10到30毫秒这是短时分析后面所有参数的基本依据。有了尺度下一步就是决定怎么切帧。直接截取相当于给信号加矩形窗矩形窗主瓣窄、分辨率好但旁瓣泄漏大帧边界处的突变会污染后续频谱分析多数仿真代码改用汉明窗旁瓣衰减明显更好代价是主瓣稍宽。对短时能量和过零率这类时域特征窗型影响没有对频谱那么大但一旦后面要接FFT做共振峰估计窗型就需要从矩形窗切到汉明窗或汉宁窗。2.2 帧长、帧移的参数表与采样率换算帧长设多少毫秒不直接写进程序而是先换算成采样点数。这是短时分析仿真中最容易出错的地方采样率fs为16000 Hz时一帧25毫秒对应round(16000 * 25 / 1000) 400个采样点如果把这段代码原样搬到8 kHz数据集上帧长就会从25毫秒变成50毫秒。常用参数换算关系如下采样率帧长10ms帧长20ms帧长25ms帧移10ms8000 Hz80点160点200点80点16000 Hz160点320点400点160点44100 Hz441点882点1103点441点帧移决定相邻帧的重叠程度。帧长400点时取帧移160点意味着相邻两帧有240点重叠。语音的基频和共振峰在时间上是连续变化的没有重叠的话帧与帧之间的过渡会被截断重叠比例一般取30%到50%既保留时间分辨率又不至于让计算量成倍增加。窗函数在这一步的作用不是滤波而是给每帧两端的采样点较低权重削弱边界突变。使用hamming时注意MATLAB的两种形态hamming(L)是对称窗hamming(L, periodic)是周期窗。分帧做频谱特征时用periodic它更利于后续FFT的循环延拓假设设计FIR滤波器时才用对称窗。这个区别在仿真验收时经常被问到下面代码统一用periodic。2.3 用矩阵化分帧计算短时能量短时能量的定义是把一帧内所有采样点平方后求和度量的是这一帧的响度。实现时分帧和计算分开分帧写成独立函数后面短时过零率、短时自相关都要复用同一个帧矩阵。function [frames, tAxis] frameSignal(x, fs, frameLenMs, frameShiftMs) % 语音信号短时分析的基础步骤: 分帧加窗 % x 单声道语音波形 % fs 采样率, 单位Hz % frameLenMs 帧长, 单位毫秒 % frameShiftMs 帧移, 单位毫秒 % frames 输出, nFrames x frameLen 矩阵, 每行是一帧 % tAxis 输出, 每帧中心时刻, 单位秒, 画图时与波形时间轴对齐 frameLen round(fs * frameLenMs / 1000); % 帧长换算为采样点 frameShift round(fs * frameShiftMs / 1000); % 帧移换算为采样点 nFrames floor((length(x) - frameLen) / frameShift) 1; win hamming(frameLen, periodic); % 行向量窗函数 frames zeros(nFrames, frameLen); for i 1:nFrames startIdx (i - 1) * frameShift 1; % 当前帧起点 frames(i, :) x(startIdx : startIdx frameLen - 1) .* win; end tAxis ((0 : nFrames - 1) * frameShift frameLen / 2) / fs; end分帧时用floor不用ceil结尾不够一帧的残段直接丢弃。丢弃残段看似浪费但语音尾部往往是渐弱的静音对后续特征影响很小强行补零反而会在帧尾制造阶跃抬高过零率。帧矩阵的行是帧索引列是帧内采样点位置短时能量就是逐行平方求和。[x, fs] audioread(speech.wav); if size(x, 2) 1 x x(:, 1); % 多声道录音取第一声道 end x x / max(abs(x)); % 幅度归一化, 避免特征受录音音量影响 [frames, tAxis] frameSignal(x, fs, 25, 10); energy sum(frames .^ 2, 2); % 每帧短时能量, 单位是归一化幅度平方参数说明25和10分别是帧长和帧移的毫秒值换算成采样点后是400点和160点。短时能量对幅值非常敏感录同一句话时麦克风距离稍有变化energy整体就会成倍变化所以幅度归一化要放在分帧之前。sum(frames .^ 2, 2)里维度参数2表示对每一行求和结果是一个nFrames乘1的列向量每个元素对应一帧能量。2.4 从能量曲线能读出什么从energy曲线能直接读到的第一层信息是语音边界静音段能量接近零第一个字起音位置对应能量抬升的地方最后一个音结束对应能量回落的起点。第二层信息是轻重音差异重读音节能量明显高于非重读音节。第三层信息是谐波强度的总体趋势但时域能量只给总强度分不出频谱细节所以能量曲线最可靠的用途还是边界检测。需要警惕的是幅度归一化只能让整段峰值等于1如果文件末尾有长时间背景噪声能量曲线不会归零而是拖出一条水平的底噪尾巴。处理办法是先取音频前50毫秒估计底噪均值再把这部分从energy里减掉边界检测结果会干净很多。这个预处理在端点检测场景几乎是必做的。3. 短时过零率的MATLAB实现清浊音边界判别3.1 过零率在语音里的物理含义短时过零率统计一帧内波形穿越零轴的次数单位是“次/样本点”。它的物理背景是语音频谱分布浊音段低频能量集中声带振动产生准周期脉冲波形在一个周期内过零次数少清音段如/s/、/sh/的频谱能量主要分布在2 kHz以上波形振荡快过零次数多。这个对比直接在时域表现为过零率高低差异。但过零率不能单独使用。静音段通常带有宽带底噪过零率可能比清音还高只用过零率做端点检测会把安静时段的随机噪声识别成语音。常规做法是把短时能量当过零率的闸门只有能量超过第一层阈值过零率才有判别意义。直流偏移是另一个常见干扰源。录音设备的AD转换、低切滤波器不彻底都会给信号叠加固定直流分量。波形整体抬升后本来该穿越零轴的波谷碰不到零轴过零率会被系统性压低。所以在计算过零率之前先减去帧内均值是比选窗型更优先的预处理。3.2 基于sign和diff的矩阵实现短时过零率直接做在帧矩阵上不需要逐帧循环。framesDC frames - mean(frames, 2); % 每帧去直流, 否则过零率被压低 signMat sign(framesDC); % 正数取1, 负数取-1 crossFlag abs(diff(signMat, 1, 2)); % 相邻样本符号不同记为2 zcr sum(crossFlag, 2) / (2 * size(framesDC, 2)); % 除以帧长的两倍diff(signMat, 1, 2)在每一行内部做一阶差分。数学上过零率定义为相邻样本符号改变的次数符号从1到-1时差分绝对值为2除以2后计数为1从0到正负符号时差分绝对值为1会产生半次伪计数。MATLAB的sign(0)返回0当一个采样点恰好落在零轴时上面代码会把它左右两侧各算半次对400点的帧来说误差极小通常不需要额外处理。分母里的2来自差分计数只和过零率定义有关不需要按采样率换算。如果你希望结果以Hz为单位把zcr乘上fs即可得到每秒过零次数便于与频率尺度对照。由于分母已经按每帧样本点归一化zcr本身是每样本点的过零概率在8 kHz和44.1 kHz采样率下阈值可以直接复用前提是帧长保持接近25毫秒。3.3 用能量和过零率区分静音、清音、浊音把短时能量和短时过零率放在一起看每一帧可以分成三类这个判别关系是前面两段代码的直接应用。帧类型短时能量短时过零率典型示例静音/底噪低中或波动句首句尾无语音段浊音高低, 通常小于0.1元音/a/、/i/、鼻音/m/清音低到中高, 通常大于0.2摩擦音/s/、/sh/、塞音除阻段0.1和0.2两个阈值以16 kHz采样率、25毫秒帧长换算后为参考因为zcr已经做过样本点归一化换采样率时阈值不必改。这个性质常被忽略仿真报告里写成“次/样本点”比写“次/秒”更稳妥。只依赖这三个特征的人声识别精度有限因为它们本质上都是宽带特征。汉语声母多为清音、韵母多为浊音单看全带过零率会把弱摩擦音漏掉更稳的做法是把信号拆成高低两个子带分别计算过零率高子带捕捉清音、低子带捕捉浊音。仿真只要算三个特征时用全带过零率足够到了端点检测环节再考虑分带。4. 短时自相关的基频估计实现延迟搜索范围与三特征对齐4.1 自相关峰与基频周期的关系短时自相关度量一帧波形与它自己平移后的重合程度。对浊音帧来说波形近似周期信号时移等于基波周期时重合度最高自相关出现强峰时移不等于周期时重合度明显下降。这个强峰对应的延迟反过来就能算出基频。但它有一个非常隐蔽的坑延迟等于0时信号与自己完全重合自相关必然取最大值。如果直接在整段自相关序列里找峰找到的一定是零延迟峰基频永远算不出来。所以扫峰时必须排除零延迟附近的区域只搜索与可能基频范围对应的延迟区间。基频范围由语音本身的声学特性决定多数成人语音的基频在60到500 Hz之间。换算成采样延迟时500 Hz对应2毫秒周期在16 kHz采样率下是32个采样点60 Hz对应16.7毫秒是267个采样点。扫描区间取32到267刚好落在帧长400点以内。4.2 带延迟搜索范围的归一化自相关实现代码先把每一帧在0到lagMax范围内的自相关值全部算出来再除以零延迟处的R0做归一化。归一化后每帧峰值都在0到1之间弱信号帧和强信号帧可以用同一套判断逻辑比较。function [f0, Rnorm] pitchByAutocorr(frames, fs, f0Min, f0Max) % 基于短时自相关的基频估计 % frames 分帧加窗后的帧矩阵 % fs 采样率 % f0Min, f0Max 基频搜索范围, 默认60和500 Hz % f0 输出每帧基频估计值, 静音帧可能出现NaN [~, frameLen] size(frames); lagMin max(2, floor(fs / f0Max)); % 最高基频对应的最小延迟 lagMax ceil(fs / f0Min); % 最低基频对应的最大延迟 if lagMax frameLen / 2 lagMax floor(frameLen / 2); % 延迟太接近帧长时估计不稳 end R zeros(size(frames, 1), lagMax 1); for lag 0 : lagMax % 列1到frameLen-lag 与 列1lag到frameLen 逐点相乘求和 R(:, lag 1) sum(frames(:, 1 : frameLen - lag) .* frames(:, 1 lag : frameLen), 2); end Rnorm R ./ R(:, 1); % 除以零延迟自相关, 幅值归一化 [~, peakIdx] max(Rnorm(:, lagMin 1 : lagMax 1), [], 2); lagEst lagMin peakIdx; % 实际峰值延迟 f0 fs ./ lagEst; endmax函数带着[], 2在行方向扫描peakIdx是子区间内的相对位置加上lagMin才是绝对延迟。搜索区间从lagMin1开始自然跳过了零延迟处必然出现的最大值。归一化除以R0在静音帧全零时会产生NaN这类帧在画图时通常直接跳过或者置为0再平滑。注意如果f0Min换算出的lagMax超过帧长一半自相关用到的有效样本太少峰形会变钝。优先加大帧长而不是缩小f0Min否则最低基频附近的周期峰可能被直接排除。4.3 用合成周期信号验证基频估计真实语音的基频曲线是未知的算法对错很难直接判断。合成信号的好处是基频已知生成一段120 Hz正弦波加噪声如果估计结果落在119到121 Hz之外代码里一定有错。fs 16000; t (0 : 16000) / fs; xSyn sin(2 * pi * 120 * t) 0.3 * randn(size(t)); % 120Hz正弦噪声 [frames, ~] frameSignal(xSyn, fs, 25, 10); [f0, ~] pitchByAutocorr(frames, fs, 60, 500); fprintf(基频估计均值: %.2f Hz\n, mean(f0));帧长25毫秒对应400点120 Hz的周期约133点一帧内包含约3个周期自相关搜索范围内有多个完整峰。噪声幅度0.3相对正弦幅值1已经不小若f0仍能稳定恢复说明分帧和自相关实现基本可靠。帧长如果缩短到10毫秒一帧内只剩1.2个周期自相关峰会被窗函数碾平这也是为什么基频估计场景不建议用太短的帧。4.4 把三个特征和波形画在同一个时间轴上仿真最后一步是画图这里最容易出问题的是时间轴错位。energy、zcr、f0都是按帧计算的每帧对应的时间点应取帧中心而不是帧起始位置。直接从帧序号除以fs会把整条特征曲线往左偏移半帧。tWave (0 : length(x) - 1) / fs; figure subplot(4,1,1); plot(tWave, x); ylabel(波形); subplot(4,1,2); plot(tAxis, energy); ylabel(短时能量); subplot(4,1,3); plot(tAxis, zcr); ylabel(短时过零率); subplot(4,1,4); plot(tAxis, f0); ylabel(基频/Hz); xlabel(时间/s);tAxis由frameSignal函数按“帧长一半 累计帧移”计算正是帧中心时刻。对比一条明显的浊音段波形对应区间能量抬高、过零率压低、基频曲线落在平直或略倾斜范围内。三个特征同时满足这三个表现才说明分帧、计算、绘图三个环节都串联正确。5. 短时特征用于端点检测与验证的三个实用技巧5.1 用能量粗定位、过零率找回清音边界端点检测最常见的问题是把以清音开头的字——比如“三”“七”“书”——的声母截掉。原因是清音段能量低仅靠能量阈值会把声母判成静音。常规做法是双门限先用能量找到稳定浊音段再向两侧扩展只要过零率仍高于清音特征就继续走直到连续若干帧低于阈值。energyThresh max(energy) * 0.1; voicedMask energy energyThresh; idx find(voicedMask); if isempty(idx) warning(未检测到语音段); startFrame 1; else startFrame max(idx(1) - 2, 1); while startFrame 1 zcr(startFrame) 0.15 startFrame startFrame - 1; % 向起始方向扩展找回清音 end end0.15是清音与浊音过零率的经验分界比3.3节里的0.1略高留了一点余量。while循环依赖前面frameSignal、短时过零率计算的结果串联执行而不是在一个独立函数里重算这样调试时可以直接看zcr曲线形状来判断阈值是否合理。5.2 自查三特征的两条快速规则第一浊音段过零率应当低于0.1清音段高于0.2如果整段zcr曲线都在0.05附近先查直流偏移是不是没扣干净再看分帧时是否把残段补零了。第二f0曲线在静音段出现随机大跳属于正常现象静音帧没有周期性看基频只关注浊音帧区间如果浊音段f0也剧烈抖动优先检查帧长低于10毫秒的帧长无法支撑低频基频估计。5.3 调参数时先动帧长再动帧移换采样率时最忌讳直接搬用旧的采样点数值。8000 Hz数据用25毫秒帧长是200点不是旧工程的400点。把换算逻辑写进frameSignal函数后改动只发生在参数行代码结构不需要变帧长动完再按30%到50%比例调帧移两个参数都落在合理区间的概率就大多了。本文还有配套的精品资源点击获取