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
| import numpy as np
from scipy import signal as sig
class FatigueIndicatorAnalyzer:
"""
疲劳指标分析器
按论文方法分析各传感器指标的疲劳敏感性
"""
def analyze_hrv(self, ecg: np.ndarray,
fs: int = 500) -> dict:
"""
ECG→HRV指标
"""
r_peaks = self._detect_r(ecg, fs)
rr_intervals = np.diff(r_peaks) / fs
mean_hr = 60 / np.mean(rr_intervals)
sdnn = np.std(rr_intervals)
rmssd = np.sqrt(np.mean(np.diff(rr_intervals)**2))
freqs, power = self._welch_psd(rr_intervals, fs)
lf = np.sum(power[(freqs >= 0.04) & (freqs < 0.15)])
hf = np.sum(power[(freqs >= 0.15) & (freqs <= 0.4)])
return {
'mean_hr': mean_hr,
'sdnn': sdnn,
'rmssd': rmssd,
'lf_hf_ratio': lf / (hf + 1e-8),
'lf_power': lf,
'hf_power': hf,
'breathing_rate': self._extract_breathing(rr_intervals, fs)
}
def analyze_eeg(self, eeg: np.ndarray,
fs: int = 500) -> dict:
"""
EEG频域指标
"""
channels, n = eeg.shape
freqs, power = self._welch_psd(eeg, fs, axis=1)
theta = np.sum(power[(freqs >= 4) & (freqs < 8)], axis=1)
alpha = np.sum(power[(freqs >= 8) & (freqs < 13)], axis=1)
beta = np.sum(power[(freqs >= 13) & (freqs <= 30)], axis=1)
return {
'theta_power': theta,
'alpha_power': alpha,
'beta_power': beta,
'theta_alpha_ratio': theta / (alpha + 1e-8),
}
def analyze_eye(self, eye_data: np.ndarray,
fs: int = 100) -> dict:
"""
眼动指标
"""
blink_dur = self._detect_blinks(eye_data, fs)
pupil_diam = eye_data[:, 2]
return {
'blink_duration': np.mean(blink_dur),
'blink_rate': len(blink_dur) / (len(eye_data)/fs/60),
'pupil_diameter': np.mean(pupil_diam),
'perclos': self._compute_perclos(eye_data, fs)
}
def _detect_r(self, ecg, fs):
"""Pan-Tompkins R峰检测（简化版）"""
from scipy.signal import find_peaks
b, a = sig.butter(2, [5/(fs/2), 15/(fs/2)], btype='band')
filtered = sig.filtfilt(b, a, ecg)
diff = np.diff(filtered) ** 2
window = int(0.15 * fs)
smoothed = np.convolve(diff, np.ones(window)/window, mode='same')
peaks, _ = find_peaks(smoothed, distance=0.3*fs,
height=0.3*np.max(smoothed))
return peaks
def _welch_psd(self, data, fs, axis=-1):
"""Welch法功率谱估计"""
from scipy.signal import welch
if data.ndim == 1:
return welch(data, fs=fs, nperseg=min(256, len(data)))
else:
return welch(data, fs=fs, nperseg=min(256, data.shape[1]),
axis=axis)
def _detect_blinks(self, eye_data, fs):
"""眨眼检测"""
ear = eye_data[:, 3]
threshold = 0.2
below = ear < threshold
durations = []
i = 0
while i < len(below):
if below[i]:
start = i
while i < len(below) and below[i]:
i += 1
durations.append((i - start) / fs * 1000)
else:
i += 1
return durations
def _compute_perclos(self, eye_data, fs):
"""PERCLOS计算"""
ear = eye_data[:, 3]
threshold = 0.25
window = int(60 * fs)
perclos_values = []
for i in range(0, len(ear) - window, window):
segment = ear[i:i+window]
closed_ratio = np.mean(segment < threshold)
perclos_values.append(closed_ratio * 100)
return np.mean(perclos_values)
def _extract_breathing(self, rr_intervals, fs):
"""从RR间期提取呼吸率"""
rr = np.array(rr_intervals)
rr_interp = np.interp(
np.arange(0, len(rr), 0.01),
np.arange(len(rr)),
rr
)
freqs, power = self._welch_psd(rr_interp, 100)
resp_mask = (freqs >= 0.1) & (freqs <= 0.5)
if np.any(resp_mask):
peak_freq = freqs[resp_mask][np.argmax(power[resp_mask])]
return peak_freq * 60
return 0
|