从‘啊’到‘哦’：用Python+Codec2拆解语音编码，可视化理解基频、共振峰与LPC

语音编码技术听起来高深莫测，但它的核心原理其实就藏在我们每天发出的简单元音里。想象一下，当你发出"啊"和"哦"这两个音时，你的声带和口腔正在执行一场精密的物理交响乐——而这正是语音编码算法试图用数学语言描述的对象。本文将带你用Python和Codec2库，通过可视化手段拆解这段"交响乐"的乐谱，无需深奥的数学公式，只需几行代码和直观的图表，就能理解语音如何被压缩到每秒几百比特的秘密。

1. 准备工作：搭建语音分析实验室

在开始实验前，我们需要准备以下工具：

pip install librosa matplotlib numpy sounddevice scipy

这些Python库将构成我们的"数字听诊器"：

• Librosa ：专业的音频分析工具包
• Matplotlib ：绘制各种语音分析图表
• NumPy ：处理音频信号的数值计算
• Sounddevice ：实时录音和播放
• SciPy ：信号处理算法支持

提示：建议使用Jupyter Notebook进行实验，可以即时查看可视化结果并交互式调整参数

录制第一个样本很简单：

import sounddevice as sd
import numpy as np

fs = 8000  # 采样率8kHz，与Codec2标准一致
duration = 2  # 录音时长2秒

print("开始录制元音'a'...")
myrecording = sd.rec(int(duration * fs), samplerate=fs, channels=1)
sd.wait()  # 等待录音完成
np.save("a_original.npy", myrecording)  # 保存原始音频

2. 可视化语音的"指纹图谱"

2.1 波形图：声音的时域肖像

波形图是最直观的声音表现形式，它展示了声压随时间的变化。用以下代码绘制"啊"的波形：

import matplotlib.pyplot as plt
import librosa.display

plt.figure(figsize=(12, 4))
librosa.display.waveshow(myrecording, sr=fs)
plt.title("'啊'的波形图")
plt.xlabel("时间(s)")
plt.ylabel("振幅")
plt.show()

你会看到类似这样的特征：

• 周期性振动 ：浊音的典型特征
• 振幅包络 ：反映发音的强弱变化
• 过零率 ：波形穿过零点的频率，清浊音的重要区分指标

2.2 频谱图：声音的频域密码

将声音转换到频域，秘密就更加清晰了：

D = np.abs(librosa.stft(myrecording.flatten(), n_fft=1024))
plt.figure(figsize=(12, 6))
librosa.display.specshow(librosa.amplitude_to_db(D, ref=np.max),
                         y_axis='log', x_axis='time', sr=fs)
plt.colorbar(format='%+2.0f dB')
plt.title("'啊'的频谱图")
plt.show()

这张图会揭示三个关键要素：

1. 基频(F0) ：最低的明显峰，决定音高
2. 谐波 ：基频整数倍的峰，形成音色基础
3. 共振峰 ：几个突出的频带，决定元音特性

注意：不同元音的共振峰分布不同，这正是区分"啊"和"哦"的关键

2.3 语谱图：声音的时空演变

语谱图是频谱的时间序列，用热力图形式展示：

plt.figure(figsize=(12, 6))
S = librosa.feature.melspectrogram(y=myrecording.flatten(), sr=fs)
S_dB = librosa.power_to_db(S, ref=np.max)
librosa.display.specshow(S_dB, x_axis='time', y_axis='mel', sr=fs)
plt.colorbar(format='%+2.0f dB')
plt.title("'啊'的语谱图")
plt.show()

观察重点：

• 横杠模式 ：稳定的共振峰表现为水平条纹
• 能量集中区 ：前三个共振峰(F1-F3)的位置
• 时间变化 ：发音起始和结束时的瞬态特征

3. 解剖元音：基频、谐波与共振峰

3.1 测量基频：声音的"心跳"

基频是声带振动的频率，用自相关法可以准确测量：

def get_pitch(signal, sr):
    signal = signal - np.mean(signal)
    signal = signal / np.max(np.abs(signal))
    
    # 自相关计算
    corr = np.correlate(signal, signal, mode='full')
    corr = corr[len(corr)//2:]
    
    # 寻找峰值
    d = np.diff(corr)
    start = np.nonzero(d > 0)[0][0]
    peak = np.argmax(corr[start:]) + start
    return sr / peak

pitch = get_pitch(myrecording.flatten(), fs)
print(f"测量到的基频：{pitch:.1f} Hz")

成年男性的典型基频范围是85-180Hz，女性是165-255Hz。这个基础频率决定了我们感知到的音高。

3.2 谐波结构：声音的"泛音列"

谐波是基频的整数倍频率成分，它们决定了声音的丰富程度。提取前10个谐波：

n_harmonics = 10
harmonics = []
for i in range(1, n_harmonics+1):
    harmonics.append(pitch * i)
    
print("谐波频率(Hz)：", [f"{h:.1f}" for h in harmonics])

3.3 共振峰定位：元音的"身份证"

共振峰是声道形状的声学表现，定位前三个共振峰：

def get_formants(signal, sr):
    # 预加重
    signal = np.append(signal[0], signal[1:] - 0.97 * signal[:-1])
    
    # 计算LPC系数
    n_coeff = 12  # 通常取采样率(kHz)+4
    a = librosa.lpc(signal, order=n_coeff)
    
    # 求根并转换为频率
    roots = np.roots(a)
    roots = roots[np.imag(roots) >= 0]  # 保留上半平面根
    angles = np.arctan2(np.imag(roots), np.real(roots))
    freqs = angles * (sr / (2 * np.pi))
    
    # 按频率排序并返回前几个
    freqs = np.sort(freqs)
    return freqs[:3]  # 返回前三个共振峰

formants = get_formants(myrecording.flatten(), fs)
print(f"共振峰F1-F3(Hz): {formants[0]:.1f}, {formants[1]:.1f}, {formants[2]:.1f}")

典型元音的共振峰范围：

• "啊"(/a/)：F1≈700Hz, F2≈1100Hz
• "哦"(/o/)：F1≈500Hz, F2≈900Hz

4. LPC分析：数学建模声道形状

线性预测编码(LPC)的核心思想是：当前采样点可以用前几个采样点的线性组合来预测。这实际上是在模拟声道的滤波特性。

4.1 计算LPC系数

order = 12  # 预测器阶数
a = librosa.lpc(myrecording.flatten(), order=order)
print("LPC系数：", a)

这些系数描述了一个虚拟的"数字声道"，其频率响应会逼近原始语音的频谱包络。

4.2 可视化LPC频谱包络

# 计算原始频谱
fft_size = 1024
X = np.fft.fft(myrecording.flatten(), fft_size)
X = np.abs(X[:fft_size//2])
freqs = np.linspace(0, fs/2, fft_size//2)

# 计算LPC频谱
w, h = scipy.signal.freqz(1, a, worN=fft_size//2)
h = np.abs(h)

# 绘制对比
plt.figure(figsize=(12, 6))
plt.plot(freqs, 20 * np.log10(X / np.max(X)), label="原始频谱")
plt.plot(freqs, 20 * np.log10(h / np.max(h)), 'r', label="LPC包络")
plt.xlabel("频率(Hz)")
plt.ylabel("幅度(dB)")
plt.legend()
plt.title("原始频谱与LPC包络对比")
plt.show()

你会看到LPC曲线如何"拥抱"住频谱的峰值，这正是共振峰的位置。

4.3 合成重建语音

用LPC系数可以重建近似语音：

def synthesize_lpc(excitation, a):
    # excitation是激励信号，a是LPC系数
    return scipy.signal.lfilter([1], a, excitation)

# 创建脉冲串激励（模拟浊音）
pulse_train = np.zeros(8000)
pulse_train[::int(fs/pitch)] = 1  # 按基频间隔放置脉冲

# 合成语音
synthesized = synthesize_lpc(pulse_train, a)

# 保存并播放
sd.play(synthesized, fs)

虽然听起来机械感明显，但已经能辨认出元音特征，这就是低码率编码的基本原理。

5. Codec2实战：从原理到应用

5.1 安装与基本使用

Codec2有多种Python绑定，推荐使用pycodec2：

pip install pycodec2

编码解码一段语音：

import pycodec2

# 初始化编解码器(2400bps模式)
c2 = pycodec2.Codec2(2400)

# 编码
encoded = c2.encode(myrecording.flatten())

# 解码
decoded = c2.decode(encoded)

# 保存解码结果
sf.write('a_decoded.wav', decoded, fs)

5.2 编解码前后对比分析

让我们对比原始与解码后的语音特征：

# 波形对比
plt.figure(figsize=(12, 6))
plt.subplot(2,1,1)
librosa.display.waveshow(myrecording, sr=fs)
plt.title("原始语音波形")
plt.subplot(2,1,2)
librosa.display.waveshow(decoded, sr=fs)
plt.title("解码后语音波形")
plt.tight_layout()
plt.show()

# 频谱对比
plt.figure(figsize=(12, 6))
plt.subplot(2,1,1)
D_orig = np.abs(librosa.stft(myrecording.flatten()))
librosa.display.specshow(librosa.amplitude_to_db(D_orig, ref=np.max),
                         y_axis='linear', x_axis='time', sr=fs)
plt.title("原始语音频谱")
plt.subplot(2,1,2)
D_dec = np.abs(librosa.stft(decoded))
librosa.display.specshow(librosa.amplitude_to_db(D_dec, ref=np.max),
                         y_axis='linear', x_axis='time', sr=fs)
plt.title("解码后语音频谱")
plt.tight_layout()
plt.show()

虽然波形细节有所丢失，但关键的共振峰结构和基频信息都得到了保留，这正是低码率编码的智慧所在——只保留人耳最敏感的特征。

5.3 码率计算与优化

Codec2在2400bps模式下，每秒语音仅需300字节：

bits_per_second = 2400
bytes_per_second = bits_per_second / 8
print(f"每秒语音数据量：{bytes_per_second}字节")

相比之下，未压缩的8kHz 16-bit单声道音频需要16kB/s，压缩比达到53:1。这种极致的压缩效率来自于：

1. 参数化编码 ：只传输基频、LPC系数等关键参数
2. 感知优化 ：专注于保留人耳敏感的特征
3. 动态比特分配 ：对不同参数采用不同精度的量化

6. 进阶实验：探索语音编码的边界

6.1 不同元音的编码对比

录制"哦"、"咦"、"呜"等不同元音，比较它们的编码效果：

vowels = ['a', 'o', 'i', 'u']
formant_data = []

for vowel in vowels:
    print(f"请录制元音'{vowel}'...")
    recording = sd.rec(int(duration * fs), samplerate=fs, channels=1)
    sd.wait()
    
    # 分析共振峰
    f1, f2, f3 = get_formants(recording.flatten(), fs)
    formant_data.append([vowel, f1, f2, f3])
    
    # 编码解码
    encoded = c2.encode(recording.flatten())
    decoded = c2.decode(encoded)
    
    # 保存对比
    sf.write(f'{vowel}_original.wav', recording, fs)
    sf.write(f'{vowel}_decoded.wav', decoded, fs)

# 展示共振峰数据
formant_table = pd.DataFrame(formant_data, columns=['元音', 'F1', 'F2', 'F3'])
print(formant_table)

你会发现不同元音在F1-F2平面上形成清晰的聚类，这正是语音编码能区分它们的基础。

6.2 极端码率下的语音质量

Codec2支持从700bps到3200bps多种模式，对比不同码率的效果：

rates = [700, 1200, 2400, 3200]
recordings = {}

for rate in rates:
    c2 = pycodec2.Codec2(rate)
    encoded = c2.encode(myrecording.flatten())
    decoded = c2.decode(encoded)
    recordings[rate] = decoded

随着码率降低，你会注意到：

• 高频细节首先丢失
• 噪声成分增加
• 但元音特征仍保持可懂度

6.3 实际应用中的优化技巧

在实际语音通信系统中，还会采用以下技术提升质量：

1. 舒适噪声生成(CNG) ：在静默段插入人工噪声避免突兀感
2. 丢包隐藏(PLC) ：当数据包丢失时智能插值恢复
3. 动态码率调整 ：根据网络状况自动切换编码模式

实现一个简单的PLC算法：

def packet_loss_concealment(previous_frame, current_lost=True):
    if current_lost:
        # 简单重复上一帧
        return previous_frame
    else:
        return None  # 正常处理

# 模拟丢包场景
packets = [encoded[i:i+8] for i in range(0, len(encoded), 8)]
lost_indices = [3, 7]  # 假设第4和第8个包丢失

decoded_with_plc = []
prev_frame = None
for i, pkt in enumerate(packets):
    if i in lost_indices:
        frame = packet_loss_concealment(prev_frame)
    else:
        frame = c2.decode(pkt)
        prev_frame = frame
    if frame is not None:
        decoded_with_plc.extend(frame)