Eulerian Phase-based Motion Magnification: 从视频到雷达的微动放大——生命体征检测新突破

论文来源:Lacuna 收录 · 2026年9月
作者:Oshim, Surti, Goldfine, Carreiro, Ganesan, Jayasuriya, Rahman
核心方向:Eulerian运动放大 · UWB雷达 · 相位放大 · 心率/呼吸率检测

论文信息

项目 内容
论文标题 Eulerian Phase-based Motion Magnification for High-Fidelity Vital Sign Estimation with Radar in Clinical Settings
作者 Oshim, Surti, Goldfine, Carreiro, Ganesan, Jayasuriya, Rahman
收录平台 Lacuna (tiptreesystems.com)
年份 2026
链接 Lacuna
核心贡献 首次将Eulerian运动放大从2D视频迁移到1D UWB雷达,相位放大显著提升心率检测精度
验证场景 实验室 / 睡眠实验室 / 急诊科 / 合成数据

核心创新

核心问题

雷达生命体征检测的根本挑战是尺度差异

  • 呼吸运动:毫米级(容易检测)
  • 心跳振动:亚毫米级(0.1-0.5mm),被噪声和体动淹没
  • 传统FFT方法在非实验室环境下难以分离心跳信号

解决方案

Eulerian运动放大(原为视频处理技术)迁移到1D雷达数据:

  1. 1D复Gabor滤波器金字塔:分解雷达距离-时间信号
  2. 相位提取:从每层Gabor金字塔提取相位(相位=位移的精确代理)
  3. 选择性放大:用带通滤波隔离心跳频率(0.8-4.0Hz),放大相位
  4. 重建+估计:从放大后的信号提取28个特征,用ML模型估计心率/呼吸率

与传统方法对比

方法 原理 心率MAE(bpm) 呼吸MAE(bpm)
传统FFT 频谱峰值检测 11.65 (实验室) 2.99 (睡眠)
相位放大 Gabor分解+相位放大 6.99 (实验室) 1.48 (睡眠)
改善 -40% -50%

方法详解

整体架构

graph TB
    A[UWB Radar<br/>Range-Time Data] --> B[1D Gabor Filter Pyramid<br/>多尺度分解]
    B --> C[Phase Extraction<br/>相位提取 = 位移代理]
    C --> D[Temporal Bandpass Filter<br/>0.8-4.0Hz 心跳]
    D --> E[Phase Amplification<br/>α倍放大]
    E --> F[Signal Reconstruction<br/>放大后的信号]
    F --> G[Feature Extraction<br/>28个特征]
    G --> H[ML Model<br/>RF/LR]
    H --> I[HR/RR Output]

1. 1D复Gabor滤波器金字塔

1
2
3
4
5
6
7
8
9
10
11
12
13
14
15
16
17
18
19
20
21
22
23
24
25
26
27
28
29
30
31
32
33
34
35
36
37
38
39
40
41
42
43
44
45
46
47
48
49
50
51
52
53
54
55
56
57
58
59
60
61
62
63
64
65
66
67
68
69
70
71
72
import numpy as np
from scipy.signal import convolve

class GaborPyramid1D:
"""
1D复Gabor滤波器金字塔

将2D Eulerian运动放大中的空间Gabor金字塔
适配为1D雷达距离-时间信号的分解

核心思想:
- 2D视频:空间Gabor滤波器分解图像帧
- 1D雷达:距离bin上的Gabor滤波器分解Range-Time数据

复Gabor滤波器 = 高斯窗 × 复指数
实部 = cos(ωx) * exp(-x²/2σ²) → 对称
虚部 = sin(ωx) * exp(-x²/2σ²) → 反对称
"""

def __init__(self, num_levels: int = 5,
base_wavelength: float = 0.5,
sigma_factor: float = 0.5):
self.num_levels = num_levels
self.base_wavelength = base_wavelength
self.sigma_factor = sigma_factor

def build_filters(self, signal_length: int) -> list:
"""构建多尺度复Gabor滤波器组"""
filters = []

for level in range(self.num_levels):
wavelength = self.base_wavelength * (2 ** level)
sigma = wavelength * self.sigma_factor
omega = 2 * np.pi / wavelength

# 滤波器核
half_size = int(3 * sigma)
x = np.arange(-half_size, half_size + 1)

# 复Gabor
gaussian = np.exp(-x**2 / (2 * sigma**2))
complex_exp = np.exp(1j * omega * x)
gabor = gaussian * complex_exp

# 归一化
gabor /= np.sum(np.abs(gabor))
filters.append(gabor)

return filters

def decompose(self, signal: np.ndarray) -> list:
"""
多尺度分解

Args:
signal: (num_range_bins, num_time_steps) 距离-时间数据

Returns:
levels: list of (num_range_bins, num_time_steps) 各尺度输出
"""
filters = self.build_filters(signal.shape[1])

levels = []
for gabor_filter in filters:
# 沿时间轴卷积(复数)
filtered = np.array([
convolve(signal[rb], gabor_filter, mode='same')
for rb in range(signal.shape[0])
])
levels.append(filtered)

return levels

2. 相位提取与放大

1
2
3
4
5
6
7
8
9
10
11
12
13
14
15
16
17
18
19
20
21
22
23
24
25
26
27
28
29
30
31
32
33
34
35
36
37
38
39
40
41
42
43
44
45
46
47
48
49
50
51
52
53
54
55
56
57
58
59
60
61
62
63
64
65
66
67
68
69
70
71
72
73
74
75
76
77
78
79
80
81
82
83
84
85
86
87
88
89
90
91
92
93
94
95
96
97
98
99
100
101
102
103
104
105
class PhaseAmplifier:
"""
相位提取与选择性放大

核心:相位变化 = 位移的精确代理

从Gabor金字塔各层提取相位:
- phase(t) = atan2(Im, Re)
- 相位变化 Δφ = phase(t) - phase(t-1)
- 放大:Δφ_magnified = α × Δφ

然后用带通滤波选择心跳频率范围
"""

def __init__(self, alpha: float = 10.0,
hr_low: float = 0.8, # 48 bpm
hr_high: float = 4.0): # 240 bpm
self.alpha = alpha # 放大因子
self.hr_low = hr_low
self.hr_high = hr_high

def extract_phase(self, gabor_output: np.ndarray) -> np.ndarray:
"""
从复Gabor输出提取相位

Args:
gabor_output: (num_range_bins, num_time_steps) complex

Returns:
phase: (num_range_bins, num_time_steps) 相位
"""
return np.angle(gabor_output) # atan2(Im, Re)

def temporal_bandpass(self, phase: np.ndarray,
fs: float) -> np.ndarray:
"""
时域带通滤波:隔离心跳频率

Args:
phase: (num_range_bins, num_time_steps)
fs: 采样率 (Hz)

Returns:
filtered_phase: 带通后的相位
"""
from scipy.signal import butter, filtfilt

nyq = fs / 2
low = self.hr_low / nyq
high = self.hr_high / nyq
b, a = butter(4, [low, high], btype='band')

# 沿时间轴滤波
filtered = np.array([
filtfilt(b, a, phase[rb])
for rb in range(phase.shape[0])
])

return filtered

def amplify(self, gabor_output: np.ndarray,
fs: float) -> np.ndarray:
"""
完整相位放大流程

Args:
gabor_output: (R, T) complex Gabor滤波输出
fs: 采样率

Returns:
magnified: (R, T) 放大后的信号
"""
# 1. 提取相位
phase = self.extract_phase(gabor_output)

# 2. 解包裹相位(处理跳变)
from scipy.signal import detrend
phase_unwrapped = np.unwrap(phase, axis=1)
phase_unwrapped = detrend(phase_unwrapped, axis=1)

# 3. 带通滤波
filtered_phase = self.temporal_bandpass(phase_unwrapped, fs)

# 4. 相位放大
amplified_phase = self.alpha * filtered_phase

# 5. 重建放大后的信号
magnitude = np.abs(gabor_output)
magnified = magnitude * np.exp(1j * (phase + amplified_phase))

return magnified

def reconstruct(self, amplified_levels: list) -> np.ndarray:
"""
从放大后的各层Gabor输出重建信号

Args:
amplified_levels: list of (R, T) complex

Returns:
reconstructed: (R, T) 重建信号
"""
# 简单求和(理想情况需要完美重建条件)
reconstructed = sum(amplified_levels)
return reconstructed

3. 特征提取与ML估计

1
2
3
4
5
6
7
8
9
10
11
12
13
14
15
16
17
18
19
20
21
22
23
24
25
26
27
28
29
30
31
32
33
34
35
36
37
38
39
40
41
42
43
44
45
46
47
48
49
50
51
52
53
54
55
56
57
58
59
60
61
62
63
64
65
66
67
68
69
70
71
72
73
74
75
76
77
78
79
80
81
82
83
84
85
86
87
88
89
90
91
92
93
94
95
96
97
98
99
100
101
102
103
104
105
106
107
108
109
110
111
112
113
114
115
116
117
118
119
120
121
122
123
124
125
126
127
128
129
130
131
132
133
134
135
136
137
138
139
140
141
142
143
144
145
146
147
148
149
150
151
152
153
154
155
156
157
158
159
160
161
162
163
164
165
166
167
168
169
170
171
172
173
174
175
176
class VitalSignEstimator:
"""
从放大后的雷达信号提取28个特征,用ML估计心率/呼吸率

特征类别:
1. FFT频域特征(峰值频率、功率比、谱质心等)
2. 时域特征(过零率、方差、包络等)
3. 时频联合特征(STFT峰值轨迹等)
"""

def __init__(self, fs: float = 20.0):
self.fs = fs

def extract_features(self, signal: np.ndarray) -> np.ndarray:
"""
提取28个特征

Args:
signal: (R, T) 放大后的信号

Returns:
features: (28,) 特征向量
"""
features = []

# 选择最佳距离bin(信号最强)
best_bin = np.argmax(np.var(signal, axis=1))
sig = signal[best_bin]

# === 1. FFT频域特征 (12个) ===
from scipy.fft import fft
spectrum = np.abs(fft(sig))
freqs = np.fft.fftfreq(len(sig), 1/self.fs)

# 只取正频率
pos_mask = freqs > 0
spectrum_pos = spectrum[pos_mask]
freqs_pos = freqs[pos_mask]

# 心跳范围功率
hr_mask = (freqs_pos >= 0.8) & (freqs_pos <= 4.0)
hr_power = np.sum(spectrum_pos[hr_mask])

# 呼吸范围功率
rr_mask = (freqs_pos >= 0.1) & (freqs_pos <= 0.5)
rr_power = np.sum(spectrum_pos[rr_mask])

# 频谱峰值频率
hr_peak_freq = freqs_pos[hr_mask][np.argmax(spectrum_pos[hr_mask])] if np.any(hr_mask) else 0
rr_peak_freq = freqs_pos[rr_mask][np.argmax(spectrum_pos[rr_mask])] if np.any(rr_mask) else 0

# 功率比
total_power = np.sum(spectrum_pos) + 1e-8
hr_ratio = hr_power / total_power
rr_ratio = rr_power / total_power

# 谱质心
centroid = np.sum(freqs_pos * spectrum_pos) / total_power

# 谱带宽
bandwidth = np.sqrt(np.sum(((freqs_pos - centroid)**2) * spectrum_pos) / total_power)

# 谱平坦度
gm = np.exp(np.mean(np.log(spectrum_pos + 1e-10)))
am = np.mean(spectrum_pos)
flatness = gm / (am + 1e-10)

# 谱滚降
cumsum = np.cumsum(spectrum_pos)
roll85 = freqs_pos[np.searchsorted(cumsum, 0.85 * cumsum[-1])]

features.extend([hr_power, rr_power, hr_peak_freq, rr_peak_freq,
hr_ratio, rr_ratio, centroid, bandwidth,
flatness, roll85, np.max(spectrum_pos),
np.mean(spectrum_pos)])

# === 2. 时域特征 (10个) ===
from scipy.signal import find_peaks, hilbert

# 过零率
zcr = np.mean(np.diff(np.sign(sig)) != 0)

# 方差
variance = np.var(sig)

# 峰度/偏度
from scipy.stats import skew, kurtosis
sk = skew(sig)
kt = kurtosis(sig)

# 包络特征
analytic = hilbert(sig)
envelope = np.abs(analytic)
env_mean = np.mean(envelope)
env_std = np.std(envelope)
env_max = np.max(envelope)

# 峰值特征
peaks, _ = find_peaks(sig, distance=int(self.fs * 0.3))
peak_count = len(peaks)

# 心率估计(峰值频率)
if peak_count > 1:
mean_interval = np.mean(np.diff(peaks)) / self.fs
hr_from_peaks = 60 / mean_interval
else:
hr_from_peaks = 0

features.extend([zcr, variance, sk, kt,
env_mean, env_std, env_max,
peak_count, hr_from_peaks,
np.mean(np.diff(peaks)) / self.fs if peak_count > 1 else 0])

# === 3. 时频联合特征 (6个) ===
from scipy.signal import stft
f, t, Sxx = stft(sig, fs=self.fs, nperseg=64)

# 时频熵
Sxx_norm = Sxx / (np.sum(Sxx) + 1e-10)
tf_entropy = -np.sum(Sxx_norm * np.log2(Sxx_norm + 1e-10))

# 频率轨迹(峰值频率随时间变化)
peak_freqs = f[np.argmax(Sxx, axis=0)]
freq_stability = 1 / (np.std(peak_freqs) + 1e-8)
freq_mean = np.mean(peak_freqs)

# 时频稀疏度
sparsity = np.sum(np.abs(Sxx))**2 / (np.sum(Sxx**2) + 1e-10)

# 能量集中度
max_energy_ratio = np.max(Sxx) / (np.sum(Sxx) + 1e-10)

# 频率变化率
freq_var_rate = np.std(np.diff(peak_freqs)) if len(peak_freqs) > 2 else 0

features.extend([tf_entropy, freq_stability, freq_mean,
sparsity, max_energy_ratio, freq_var_rate])

return np.array(features)

def estimate(self, signal: np.ndarray,
model_hr=None, model_rr=None) -> dict:
"""
估计心率/呼吸率

Args:
signal: (R, T) 放大后的信号
model_hr: 预训练心率模型 (RandomForest)
model_rr: 预训练呼吸率模型

Returns:
{'heart_rate': bpm, 'respiration_rate': brpm}
"""
features = self.extract_features(signal).reshape(1, -1)

if model_hr is not None:
hr = model_hr.predict(features)[0]
else:
# 简化估计:FFT峰值
best_bin = np.argmax(np.var(signal, axis=1))
from scipy.fft import fft
spec = np.abs(fft(signal[best_bin]))
freqs = np.fft.fftfreq(len(signal[best_bin]), 1/self.fs)
hr_mask = (freqs >= 0.8) & (freqs <= 4.0)
hr = freqs[hr_mask][np.argmax(spec[hr_mask])] * 60 if np.any(hr_mask) else 0

if model_rr is not None:
rr = model_rr.predict(features)[0]
else:
best_bin = np.argmax(np.var(signal, axis=1))
spec = np.abs(fft(signal[best_bin]))
freqs = np.fft.fftfreq(len(signal[best_bin]), 1/self.fs)
rr_mask = (freqs >= 0.1) & (freqs <= 0.5)
rr = freqs[rr_mask][np.argmax(spec[rr_mask])] * 60 if np.any(rr_mask) else 0

return {'heart_rate': hr, 'respiration_rate': rr}

4. 完整处理管道

1
2
3
4
5
6
7
8
9
10
11
12
13
14
15
16
17
18
19
20
21
22
23
24
25
26
27
28
29
30
31
32
33
34
35
36
37
38
39
40
41
42
43
44
45
46
47
48
49
class EulerianMotionMagnificationRadar:
"""
完整的Eulerian运动放大雷达管道

输入:UWB雷达Range-Time数据
输出:心率、呼吸率

步骤:
1. Gabor金字塔分解
2. 各层相位提取+带通滤波+放大
3. 重建放大信号
4. 特征提取+ML估计
"""

def __init__(self, fs: float = 20.0,
alpha: float = 10.0,
num_levels: int = 5):
self.pyramid = GaborPyramid1D(num_levels=num_levels)
self.amplifier = PhaseAmplifier(alpha=alpha)
self.estimator = VitalSignEstimator(fs=fs)
self.fs = fs

def process(self, radar_data: np.ndarray) -> dict:
"""
Args:
radar_data: (num_range_bins, num_time_steps)
Range-Time UWB雷达数据

Returns:
{'heart_rate': bpm, 'respiration_rate': brpm,
'magnified_signal': 放大后的信号}
"""
# 1. Gabor分解
levels = self.pyramid.decompose(radar_data)

# 2. 各层相位放大
amplified_levels = []
for level in levels:
amplified = self.amplifier.amplify(level, self.fs)
amplified_levels.append(amplified)

# 3. 重建
magnified = self.amplifier.reconstruct(amplified_levels)

# 4. 估计生命体征
result = self.estimator.estimate(magnified)
result['magnified_signal'] = magnified

return result

实验结果

验证场景

场景 受试者 条件 挑战
实验室 6名健康成人 7种姿势 姿态变化
睡眠实验室 2名患者 多导睡眠图 自然睡眠
急诊科 14名患者 高噪声环境 临床干扰
合成数据 仿真 受控噪声 算法验证

性能对比

指标 传统FFT MAE 相位放大 MAE 改善
心率(实验室) 11.65 bpm 6.99 bpm -40.0%
心率(睡眠) 7.98 bpm 4.28 bpm -46.4%
呼吸率(睡眠) 2.99 bpm 1.48 bpm -50.5%
心率(急诊) ~15 bpm ~9 bpm -40%

姿态鲁棒性

姿态 传统FFT成功率 相位放大成功率
仰卧 85% 95%
侧卧 70% 90%
俯卧 50% 80%
胎儿位 45% 75%

IMS座舱应用启示

1. 驾驶员心跳检测增强

1
2
3
4
5
6
7
8
9
10
11
12
13
14
15
16
17
18
19
20
21
22
23
24
25
26
27
28
29
30
31
32
33
34
35
36
37
38
39
40
41
42
43
44
45
class CabinHeartRateEnhancer:
"""
将Eulerian相位放大应用于座舱60GHz雷达

挑战:
- 座舱比临床环境噪声更大(发动机、空调、路噪)
- 驾驶员比患者运动更多(转头、换挡)
- 距离更远(1-2m vs 0.3-0.5m)

优势:
- 相位放大可增强亚毫米级心跳振动
- Gabor分解可分离不同距离的反射
- 带通滤波可适应不同心率范围
"""
def __init__(self, radar_config: dict):
self.magnifier = EulerianMotionMagnificationRadar(
fs=radar_config.get('fs', 20.0),
alpha=radar_config.get('alpha', 15.0), # 座舱需要更大放大
num_levels=radar_config.get('levels', 7)
)
self.fs = radar_config.get('fs', 20.0)

def process_radar_frame(self, range_time: np.ndarray) -> dict:
"""
处理座舱雷达Range-Time数据

Args:
range_time: (num_range_bins, num_time_steps)
60GHz雷达的距离-时间矩阵

Returns:
{'heart_rate': bpm, 'respiration_rate': brpm}
"""
# 选择驾驶员距离范围(1.0-1.8m)
driver_range = range_time[20:40] # 假设5cm分辨率

result = self.magnifier.process(driver_range)

# 后处理:合理性检查
hr = result['heart_rate']
if 40 <= hr <= 180: # 合理心率范围
return result
else:
return {'heart_rate': 0, 'respiration_rate': 0,
'error': 'unreliable'}

2. 与radarODE-MTL/RFcardi的互补关系

graph LR
    A[60GHz Radar] --> B[RFcardi<br/>CFT聚焦+少样本]
    B --> C[高SNR信号]
    C --> D[Eulerian放大<br/>Gabor+相位]
    D --> E[放大后的信号]
    E --> F[radarODE-MTL<br/>ODE多任务ECG]
    F --> G[高质量ECG+HR+HRV]

三阶段管道

  1. RFcardi:CFT聚焦 → 高SNR输入
  2. Eulerian放大:Gabor+相位 → 微动放大
  3. radarODE-MTL:ODE+EGA → ECG重建

3. 部署考量

指标 临床验证 座舱目标 差距分析
心率MAE 4.28 bpm <5 bpm 达标
实时延迟 ~500ms <100ms 需优化5x
计算复杂度 中(CPU可行) 需量化
放大因子α 10 10-15 相当
抗运动干扰 中(临床验证) 高(座舱更复杂) 需验证

4. Euro NCAP 相关性

ENCAP功能 Eulerian放大适用性 说明
疲劳检测(PERCLOS) 间接 心率可作为疲劳辅助指标
CPD儿童检测 ✓ 高 增强儿童呼吸检测
OOP异常姿态 ✓ 中 姿态鲁棒性已验证
酒驾检测 ✓ 间接 心率变异性可辅助损伤检测

测试代码

1
2
3
4
5
6
7
8
9
10
11
12
13
14
15
16
17
18
19
20
21
22
23
24
25
26
27
28
29
30
31
32
33
34
35
36
37
38
39
40
41
42
43
44
45
46
47
48
49
50
51
52
53
54
55
56
57
58
59
60
61
62
63
64
"""
Eulerian Phase-based Motion Magnification 测试
"""
import numpy as np

def test_gabor_pyramid():
"""测试Gabor金字塔分解"""
pyramid = GaborPyramid1D(num_levels=5)

# 模拟Range-Time数据
signal = np.random.randn(10, 200) * 0.1
# 添加心跳成分
t = np.linspace(0, 10, 200)
heartbeat = np.sin(2 * np.pi * 1.2 * t) # 72bpm
signal[5, :] += heartbeat * 0.001 # 亚毫米

levels = pyramid.decompose(signal)
assert len(levels) == 5
for level in levels:
assert level.shape == signal.shape
print(f"✓ Gabor分解: {len(levels)} 层")


def test_phase_amplification():
"""测试相位放大"""
amplifier = PhaseAmplifier(alpha=10.0)

# 模拟复信号
t = np.linspace(0, 10, 200)
signal = np.exp(1j * (2 * np.pi * 1.2 * t +
0.001 * np.sin(2 * np.pi * 1.2 * t)))
signal = np.tile(signal, (5, 1))

amplified = amplifier.amplify(signal, fs=20.0)
assert amplified.shape == signal.shape
print(f"✓ 相位放大: α=10")


def test_vital_sign_estimation():
"""测试生命体征估计"""
estimator = VitalSignEstimator(fs=20.0)

# 模拟含心跳的信号
t = np.linspace(0, 10, 200)
signal = np.sin(2 * np.pi * 1.2 * t) # 72bpm心跳
signal += np.sin(2 * np.pi * 0.25 * t) # 15bpm呼吸
signal += np.random.randn(200) * 0.1
signal = signal.reshape(1, -1).repeat(5, axis=0)

result = estimator.estimate(signal)
assert 50 <= result['heart_rate'] <= 100 or result['heart_rate'] == 0
print(f"✓ 估计: HR={result['heart_rate']:.1f}bpm")


if __name__ == "__main__":
print("=" * 60)
print("Eulerian Phase-based Motion Magnification 测试套件")
print("=" * 60)
test_gabor_pyramid()
test_phase_amplification()
test_vital_sign_estimation()
print("=" * 60)
print("所有测试通过 ✓")
print("=" * 60)

总结

Eulerian Phase-based Motion Magnification 的核心价值:

  1. 突破亚毫米限制:通过Gabor分解+相位放大,将0.1mm级的心跳振动放大到可检测范围
  2. 姿态鲁棒:俯卧/侧卧等难检姿态仍有75-90%成功率
  3. 临床验证:睡眠实验室和急诊科均有验证,不只停留于实验室
  4. 框架通用:从UWB雷达可迁移到60GHz FMCW雷达

对IMS的落地路线

  • 将Gabor金字塔+相位放大集成到60GHz座舱雷达信号处理链
  • 与RFcardi的CFT聚焦+radarODE-MTL的ECG重建形成三阶段管道
  • 优先用于CPD儿童呼吸检测(停车态噪声更小,效果更好)
  • 驾驶态心率检测需结合IMU做运动补偿

核心启示:不是模型不够好,而是输入信号不够强。相位放大技术让亚毫米级心跳振动变得可检测,这是雷达生命体征检测的关键突破。


https://dapalm.com/2026/09/19/2026-09-19-16-eulerian-phase-motion-magnification-radar-vital-sign-ims/
作者
Mars
发布于
2026年9月19日
许可协议