如何去除周期信号噪声以优化傅里叶变换的滤波效果

解决方案 #

1. 前置信号清洗，剔除无效值 #

核心问题是全局FFT对全信号统计特征敏感，尾部异常值会拉偏主频计算结果。首先可以通过统计阈值过滤异常值，同时可选直接截断尾部连续无效段：

def preprocess_signal(sig, sigma=3, min_valid_len=10):
    # 3σ原则过滤离群异常点
    mean_val = np.mean(sig)
    std_val = np.std(sig)
    abnormal_mask = np.abs(sig - mean_val) > sigma * std_val
    sig_clean = sig.copy()
    
    # 异常点用相邻有效值填充
    for idx in np.where(abnormal_mask)[0]:
        if 0 < idx < len(sig_clean)-1:
            sig_clean[idx] = (sig_clean[idx-1] + sig_clean[idx+1])/2
        elif idx == 0:
            sig_clean[idx] = sig_clean[idx+1]
        else:
            sig_clean[idx] = sig_clean[idx-1]
    
    # 截断尾部连续无效段
    valid_count = 0
    valid_end = len(sig_clean)
    for idx in range(len(sig_clean)-1, -1, -1):
        if not abnormal_mask[idx]:
            valid_count += 1
            if valid_count >= min_valid_len:
                valid_end = idx + 1
                break
        else:
            valid_count = 0
    return sig_clean[:valid_end]

调用时只需要在滤波前加一步：sig = preprocess_signal(cutting_sig(test).F2.values)即可。

2. 优化FFT滤波逻辑 #

现有滤波逻辑只保留单主频分量，灵活性差，可改为保留主频附近的小范围频段，兼顾平滑效果和信号真实波动：

def optimized_fft_filter(sig, bandwidth_ratio=0.1):
    time_step = 1
    sig_fft = fftpack.fft(sig)
    power = np.abs(sig_fft)**2
    sample_freq = fftpack.fftfreq(sig.size, d=time_step)
    
    pos_mask = sample_freq > 0
    freqs = sample_freq[pos_mask]
    peak_freq = freqs[power[pos_mask].argmax()]
    
    # 保留主频±10%范围内的所有频率分量
    low = peak_freq * (1 - bandwidth_ratio)
    high = peak_freq * (1 + bandwidth_ratio)
    filtered_fft = sig_fft.copy()
    filtered_fft[(np.abs(sample_freq) < low) | (np.abs(sample_freq) > high)] = 0
    # 取实部避免复数干扰后续计算
    return np.real(fftpack.ifft(filtered_fft))

3. 备选平滑方案（抗干扰性更强） #

如果FFT滤波仍然不稳定，可以直接用Savitzky-Golay滤波替代FFT步骤，这类局部滤波方法对全局异常值的抗干扰能力更强，更适合当前场景：

from scipy.signal import savgol_filter
# 窗口大小可根据实际周期长度调整，一般取周期长度的1/5~1/3，多项式阶数选2即可
sig_smoothed = savgol_filter(sig_clean, window_length=21, polyorder=2)

4. 突变点检测优化 #

差分后添加阈值过滤，避免小波动产生的误检：

deriv = np.diff(sig_smoothed)
# 取差分绝对值的75分位数作为阈值，可根据实际效果调整
threshold = np.percentile(np.abs(deriv), 75)
pos_peaks, _ = find_peaks(deriv, height=threshold)
neg_peaks, _ = find_peaks(-deriv, height=threshold)

内容的提问来源于stack exchange，提问作者Han Zhengzu