基于Numpy的谐波滤波器边缘不平滑问题排查与解决

问题描述 #

尝试编写谐波滤波器处理WAV文件，运行后输出信号边缘不平滑，当前输出效果与预期平滑效果存在差异。已知输入WAV文件的采样率sample_rate和数据长度len(data)均为10000，已尝试调整参数、核查基频索引但未解决问题。

原代码 #

from scipy.io import wavfile
import numpy as np

wave_file_path = 'in.wav'
num_harmonics = 10
base_frequency = 50
sample_rate, data = wavfile.read(wave_file_path)
#  after this - sample_rate = 10,000 & len(data) = 10,000

frequencies = np.fft.fftfreq(len(data), 1.0 / sample_rate)
spectrum = np.fft.fft(data, len(data))

base_index = int(base_frequency * len(data) / sample_rate)
# the above - base_index = 50 

# Remove harmonics
for harmonics in range(2, num_harmonics + 1):
# is this the right multiply? harmonics * base_index? it maximum=500, whenever the number of spectrum elements is much bigger (10,000)
    harmonic_index = int(harmonics * base_index)
    spectrum[harmonic_index] = 0

# Inverse Fourier Transform to get filtered signal
filtered_data = np.round(np.fft.ifft(np.fft.fftshift(spectrum)).real).astype(np.int16)

output_file_path = 'out.wav'

wavfile.write(output_file_path, len(filtered_data), filtered_data)

代码问题分析 #

• 频谱对称性处理缺失：np.fft.fft输出的频谱是共轭对称的（包含正、负频率），仅置零正频率的谐波位置会破坏对称性，逆变换后信号出现失真和边缘不平滑。
• 单点置零无法有效滤除谐波：实际信号的谐波能量并非集中在单一FFT bin，而是分布在相邻bin中，仅将单个索引置零无法彻底消除谐波。
• fftshift误用：逆变换前使用fftshift会打乱原始频谱的频率轴顺序，导致逆变换结果错误。
• WAV写入参数错误：wavfile.write的第二个参数要求是采样率，原代码传入len(filtered_data)会导致输出文件采样率异常，引发播放时的信号畸变。

修正后的谐波滤波器代码 #

from scipy.io import wavfile
import numpy as np

wave_file_path = 'in.wav'
num_harmonics = 10
base_frequency = 50
sample_rate, data = wavfile.read(wave_file_path)

# 处理立体声输入，转为单声道
if len(data.shape) > 1:
    data = data[:, 0]

# 计算FFT频谱与对应频率轴
spectrum = np.fft.fft(data)
freqs = np.fft.fftfreq(len(data), 1.0 / sample_rate)

# 定义谐波频率的滤波带宽（围绕目标频率的区间，可根据实际调整）
bin_width = sample_rate / len(data)
filter_bandwidth = bin_width * 1.5

# 逐一滤除各次谐波，同时处理正、负频率
for harmonic_order in range(2, num_harmonics + 1):
    target_freq = harmonic_order * base_frequency
    # 匹配目标谐波的正、负频率区间
    freq_mask = (np.abs(freqs - target_freq) < filter_bandwidth) | (np.abs(freqs + target_freq) < filter_bandwidth)
    spectrum[freq_mask] = 0

# 逆FFT得到滤波后的实信号
filtered_data = np.round(np.fft.ifft(spectrum).real).astype(np.int16)

# 写入正确采样率的WAV文件
output_file_path = 'out.wav'
wavfile.write(output_file_path, sample_rate, filtered_data)

关键修正说明 #

• 保证频谱对称性：同时对正、负频率对应的谐波区间置零，确保频谱共轭对称，逆变换后得到无失真的实信号。
• 区间滤波替代单点置零：通过设置带宽，覆盖谐波能量分布的相邻FFT bin，彻底滤除谐波成分。
• 移除错误的fftshift：直接对原始FFT频谱处理后做逆变换，避免频率轴错位。
• 修正WAV写入参数：传入正确的采样率sample_rate，保证输出文件格式正常。
• 兼容立体声输入：自动将立体声转为单声道处理，避免多通道数据引发的错误。

内容的提问来源于stack exchange，提问作者חדווה מזרחי