1. 项目背景与核心目标
在可穿戴设备和远程医疗监测领域，光电容积图(PPG)信号处理一直是个既基础又关键的课题。去年我在开发一款健康手环时，发现市面上的心率算法在运动场景下普遍存在15-20%的误差率。这促使我深入研究PPG信号的本质特征，最终构建出这套基于短时傅里叶变换(STFT)的滤波方案。
PPG信号是通过皮肤表面的光电传感器采集的微血管搏动信号，理论上应该包含准确的心率信息。但实际采集时会混杂三种主要噪声：
- 运动伪影（频率0.1-5Hz）
- 环境光干扰（50/60Hz工频及其谐波）
- 传感器基线漂移（<0.1Hz）
传统方案采用固定带宽的带通滤波器，但运动时心率频段（0.8-3Hz）会与运动伪影频段重叠。这就是为什么你在跑步时手环心率经常"卡住"或乱跳的根本原因。
2. 信号生成与预处理
2.1 PPG信号模拟
我们先构建理想的PPG信号模型。一个完整的心搏周期包含：
- 收缩期陡升沿（0.1-0.2秒）
- 舒张期缓降沿（0.3-0.4秒）
- 重搏切迹（降支中段的微小凸起）
用以下公式生成10秒的模拟信号（采样率100Hz）：
python复制import numpy as np
t = np.linspace(0, 10, 1000)
heart_rate = 1.2 # Hz (72 BPM)
cardiac_phase = 2 * np.pi * heart_rate * t
# 基础PPG波形
ppg_clean = 0.3 * np.sin(cardiac_phase) + \
0.1 * np.sin(2*cardiac_phase + 0.5) + \
0.05 * np.sin(3*cardiac_phase + 1.2)
# 添加运动伪影（步行频率约1-2Hz）
motion_artifact = 0.5 * np.sin(2*np.pi*1.5*t) * np.random.normal(1, 0.1, len(t))
# 添加工频干扰
powerline_noise = 0.2 * np.sin(2*np.pi*50*t)
# 合成信号
ppg_noisy = ppg_clean + motion_artifact + powerline_noise
2.2 噪声特性分析
通过功率谱密度分析可以看到：
- 心率基频在1.2Hz处
- 运动伪影在1.5Hz形成谱峰重叠
- 50Hz工频干扰形成尖峰
关键发现：当运动频率与心率频率差值小于0.3Hz时，传统滤波器无法有效分离
3. STFT时频分析实现
3.1 窗口参数选择
短时傅里叶变换的核心是平衡时域和频域分辨率：
- 窗口长度越长，频域分辨率越高但会模糊快速变化
- 汉宁窗比矩形窗有更好的旁瓣抑制
经过实测对比，选择：
- 窗口长度：5秒（平衡心率变化跟踪能力）
- 重叠率：75%（确保平滑过渡）
- 窗函数：汉宁窗
python复制from scipy import signal
f, t, Zxx = signal.stft(ppg_noisy, fs=100, window='hann',
nperseg=500, noverlap=375)
3.2 动态滤波算法
创新点在于根据信号特征动态构建掩模：
- 定位最大能量频点作为候选心率
- 检测其谐波关系（2倍频处应有次高峰）
- 排除不符合生理规律的频点（>3.5Hz或<0.7Hz）
python复制# 构建时频掩模
mask = np.zeros_like(Zxx)
for i in range(Zxx.shape[1]):
# 当前时间片的频谱
spectrum = np.abs(Zxx[:,i])
# 寻找前三个峰值
peaks, _ = signal.find_peaks(spectrum, height=0.2*np.max(spectrum))
freqs = f[peaks]
# 谐波验证
valid_peaks = []
for pk in peaks:
if 0.7 < f[pk] < 3.5: # 生理范围
# 检查是否存在近似整数倍频
has_harmonic = any(np.abs(f[pk]*n - freqs).min() < 0.1
for n in [2,3])
if has_harmonic or len(valid_peaks) < 1:
valid_peaks.append(pk)
# 构建掩模
for pk in valid_peaks:
mask[pk-2:pk+3, i] = 1 # 保留±0.2Hz带宽
4. 信号重建与验证
4.1 逆变换与心率提取
应用掩模后执行逆STFT：
python复制_, ppg_filtered = signal.istft(Zxx * mask, fs=100,
window='hann', nperseg=500,
noverlap=375)
# 峰值检测
peaks, _ = signal.find_peaks(ppg_filtered, distance=50) # 最小间隔0.5秒
hr_instant = 60 / np.diff(peaks) # 瞬时心率(BPM)
4.2 性能评估指标
引入三个关键指标：
- 均方根误差(RMSE)：<3 BPM为优秀
- 皮尔逊相关系数：>0.9为合格
- 异常值比例：>5%需优化
实测结果：
- 静态场景：RMSE=1.2 BPM
- 慢跑场景：RMSE=2.8 BPM
- 对比传统带通滤波：误差降低63%
5. 工程实现优化
5.1 实时处理技巧
在嵌入式设备上需注意：
- 采用环形缓冲区减少内存拷贝
- 预计算窗函数系数
- 定点数优化FFT计算
c复制// ARM Cortex-M4 优化示例
void arm_stft_f32(const float32_t *input,
float32_t *output,
uint16_t fftSize) {
arm_rfft_fast_instance_f32 S;
arm_rfft_fast_init_f32(&S, fftSize);
// 应用汉宁窗并执行FFT
arm_mult_f32(input, hann_window, buffer, fftSize);
arm_rfft_fast_f32(&S, buffer, output, 0);
}
5.2 运动补偿进阶方案
对于高强度运动场景，建议：
- 增加三轴加速度计数据融合
- 采用自适应滤波消除运动伪影
- 机器学习分类器识别信号质量
实测数据：结合加速度补偿后，跑步误差可从4.2BPM降至2.1BPM
6. 常见问题排查
6.1 信号完全失真
可能原因：
- 传感器接触不良（检查DC分量）
- 采样率过低（需≥75Hz）
- 环境光过强（增加光学隔离）
6.2 心率跳变异常
排查步骤：
- 检查STFT时间窗是否过短
- 验证谐波检测阈值
- 查看时频图中是否有突发干扰
6.3 嵌入式部署问题
内存优化方案：
- 降低FFT点数（256点足够）
- 使用16位定点运算
- 分块处理减少RAM占用
我在实际部署中发现，对于STM32F4系列芯片，将浮点运算替换为Q15定点数格式，可使内存占用减少60%，而精度损失仅0.3BPM。