多通道方波噪声去除咨询：基于Scipy处理传感器电机干扰问题

原FFT方案失效的核心原因 #

• 采样率限制：单通道采样率仅1Hz，奈奎斯特频率为0.5Hz，电机产生的高频方波噪声已经发生频谱混叠，无法通过常规FFT识别真实噪声频率，自然无法通过频点滤除实现消噪
• 噪声非平稳特性：方波噪声为不定期出现的干扰，不属于全段信号的固定频率分量，全局60秒FFT提取的共有频率准确率极低，强行滤除会同时损伤有效信号
• 方波信号特征：方波包含多次谐波分量，单一或少数频点滤除无法完全消除噪声，反而会引入信号畸变

适配场景的可行处理方案 #

方案1：多通道共模噪声抵消（优先选择，实现最简单） #

电机干扰属于所有传感器共有的同步噪声，而各通道有效信号独立，直接提取共模分量扣除即可，无需频域处理，适配不定期出现的噪声场景，代码示例如下：

import numpy as np
from scipy import signal

# 输入数据形状约定为 (6, n_samples)，即6个通道、每个通道对应n个采样点
def remove_common_mode_noise(data, corr_threshold=0.8):
    # 用中位数计算共模分量，比均值更抗单通道异常值影响
    common_noise = np.median(data, axis=0)
    for ch_idx in range(data.shape[0]):
        # 计算单通道和共模分量的相关系数，确认存在共模噪声再扣除
        corr = np.corrcoef(data[ch_idx], common_noise)[0, 1]
        if corr > corr_threshold:
            data[ch_idx] = data[ch_idx] - common_noise
    return data

方案2：时域形态学滤波（适配方波特征匹配消噪） #

如果部分通道不存在共模噪声，可利用方波陡峭上升/下降沿、平顶的时域特征识别噪声段，用有效信号插值替换，代码示例如下：

def remove_square_noise(data, edge_threshold=0.1, flat_threshold=0.05):
    # 阈值可根据你的实际信号幅值范围等比例调整
    cleaned_data = data.copy()
    for ch_idx in range(data.shape[0]):
        sig = cleaned_data[ch_idx]
        # 一阶差分识别陡峭跳变沿
        diff = np.diff(sig, prepend=sig[0])
        edge_pos = np.where(np.abs(diff) > edge_threshold)[0]
        # 配对上升沿、下降沿，识别方波区间
        for i in range(0, len(edge_pos)-1, 2):
            start = edge_pos[i]
            end = edge_pos[i+1] if i+1 < len(edge_pos) else len(sig)-1
            # 验证区间内信号平坦，确认是方波噪声
            segment = sig[start:end+1]
            if np.max(segment) - np.min(segment) < flat_threshold:
                # 用前后有效信号线性插值替换噪声段
                if start > 0 and end < len(sig)-1:
                    cleaned_data[ch_idx, start:end+1] = np.linspace(sig[start-1], sig[end+1], end-start+1)
    return cleaned_data

方案3：改进频域滤波（仅适配固定频率周期性方波噪声） #

如果确认噪声为0.5Hz以下的周期性干扰，可采用短时傅里叶变换做滑动窗口的自适应滤波，避免全局FFT的缺陷，代码示例如下：

def stft_based_filter(data, fs=1, nperseg=30):
    cleaned_data = np.zeros_like(data)
    for ch_idx in range(data.shape[0]):
        f, t, Zxx = signal.stft(data[ch_idx], fs=fs, nperseg=nperseg)
        mag = np.abs(Zxx)
        # 取幅度前5%的高能量分量作为噪声置零，阈值可根据实际效果调整
        noise_threshold = np.percentile(mag, 95)
        Zxx[mag > noise_threshold] = 0
        # 逆变换回时域
        _, cleaned_sig = signal.istft(Zxx, fs=fs)
        cleaned_data[ch_idx] = cleaned_sig[:len(data[ch_idx])]
    return cleaned_data

实操建议 #

• 优先测试方案1，绝大多数共模干扰场景下该方案效果最优，计算量最小，不会损伤有效信号
• 可组合使用方案1+方案2，先消除共模噪声，再对残留的零星方波做形态学滤除
• 所有阈值参数可根据你实际采集的信号幅值、噪声特征做微调

内容的提问来源于stack exchange，提问作者Xu Chang