面向人体呼吸心跳同时监测的雷达信号处理方法综述*

2023-12-26薛帅康王小月许致火刘微雪黄同辉

电讯技术 2023年12期

薛帅康,王小月,许致火,刘微雪,黄同辉

(南通大学 a.信息科学技术学院;b.交通与土木工程学院,江苏南通 226019)

0 引言

2017年联合国《世界人口老龄化发展报告》显示,全球60岁及以上的人口占总人口的13%,并以每年3%的速度增长。据2016年第四次中国城乡老人生活状况抽样调查数据显示,中国失能、半失能老年人大致有4 063万人,占老年人口的18.3%,健康老龄化成为中国人口发展目标[1]。人口老龄化带来老人健康问题,而大多数子女由于房贷和生活的压力,并不能总是在老人身边照顾,对老年人健康状况的监测变得十分必要。

呼吸、心跳等人体生命体征是现代卫生保健和医疗应用的重要指标[2],表征着人体的健康状况。传统的人体生命体征检测通常采用穿戴式传感器、粘贴式电极、呼吸带等仪器,如心电图、指夹式脉搏血氧仪、电子血压测量机等[3],虽然这些接触式心跳检测方法测量精度较高,但是电极的刺激会引起人体不适[4]。

与上述检测方法相比,雷达对于呼吸和心跳的监测具有很大优势。雷达具有较强的信号穿透性,信号波长越长,穿透能力越强,能在墙壁遮挡情况下对人体呼吸和心跳进行检测,有利于灾害应急救援[5]。同时雷达的检测范围更为广阔,可以实现一定区域的多人检测,且检测不易受外界环境条件的影响,稳定性更好。研究表明,无接触技术能够持续地监测用户的健康状况,这使其能够在发生晕厥、心脏骤停等紧急情况时及时示警,降低死亡或致残的风险,从而减轻社会医疗负担[6],这使得雷达检测技术在人体呼吸心跳监测领域拥有较为广泛的应用。

脉冲超宽带雷达[7]、连续波(Continuous Wave,CW)多普勒雷达[8-9]、调频连续波(Frequency Modulated Continuous Wave,FMCW)雷达[10-11]等都可实现对人体呼吸心跳的检测。雷达通过获取人体的胸腔起伏信号,实现对人体呼吸心跳的检测。人体胸腔的起伏很微弱,属于毫米级范畴[12],位于毫米波段的雷达可实现小位移的精确测量,对呼吸心跳的检测更准确[13]。

本文简要描述了不同雷达系统对人体呼吸和心跳信号检测的基本原理,总结了不同雷达系统在人体呼吸心跳检测方面的优缺点;从数据预处理、呼吸心跳分离、呼吸率和心率估计三个方面介绍了雷达信号处理算法,系统归纳了每一方面的最新研究进展;分析了现有研究中存在的局限性,并对未来的研究趋势进行了展望。

1 雷达检测原理

1.1 CW雷达检测原理

根据雷达射频回波信号处理结构的不同,CW雷达系统分为基于混频器的CW雷达和基于六端口接收机的CW雷达。CW雷达发射复信号为

STX(t)=ATXej[2πfct+φ0]。

(1)

式中:fc为载波频率;ATX为发射信号的振幅;φ0为初始相位。

假设人体的径向距离为d0,呼吸心跳产生的胸腔位置变化为x(t),由雷达作用人体后向散射得到的接收信号为

SRX(t)=ARXej[2πfc(t-td)+φ1]。

(2)

对于基于混频器的CW雷达,接收信号经过混频和低通滤波后,其基带复信号表示为

(3)

式中:Δφ=φ1-φ0为相位差;AIF为基带信号振幅。利用反正切法,可得到人体的胸腔位移信号x(t)为

(4)

图1为六端口接收机结构示意图,包括2个输入端口和4个输出端口。

图1 六端口接收机结构[14]

输入端口P1和P2对应发射信号STX和接收信号SRX,即

P1=ATXej[2πfct+φ0]=A+jB,

(5)

P2=ARXej[2πfc(t-td)+φ1]=C+jD。

(6)

式中:A=ATXcos(2πfct+φ0);B=ATXsin(2πfct+φ0);C=ARXcos(2πfc(t-td)+φ1);D=ARXsin(2πfc(t-td)+φ1)。P1和P2在输出端口的相位差为0,π/2,π和3π/2,通过计算得到输出信号B3,B4,B5,B6为

B3=0.25|P1+jP2|2,

(7)

B4=0.25|jP1+P2|2,

(8)

B5=0.25|jP1+jP2|2,

(9)

B6=0.25|P1-P2|2。

(10)

将式(5)、(6)代入到式(7)～(10)中得到

B3=0.25|(A-D)+j(B+C)|2,

(11)

B4=0.25|(C-B)+j(A+D)|2,

(12)

B5=0.25|-(B+D)+j(A+C)|2,

(13)

B6=0.25|(A-C)+j(B-D)|2。

(14)

而P1和P2的相位差可表示为

(15)

即由六端口接收机的输出,可得到人体的胸腔位移

(16)

1.2 FMCW雷达检测原理

FMCW雷达具有多种频率调制方式,在人体呼吸心跳检测中,大多采用线性调频方式。假设fc为载波频率,Tc为信号扫描周期,B为雷达信号的带宽,FMCW雷达发射复信号为

(17)

式中:ATX为发射信号的振幅;φ0为初始相位。

接收信号为发射信号延时td,

(18)

式中:ARX为接收信号的振幅;φ1为接收信号相位。

接收信号经过混频和低通滤波后,其基带复信号表示为

(19)

式中:Δφ=φ1-φ0相位差;AIF为基带信号振幅。为得到表征胸腔位移的相位信号,需要消除式(19)指数中的第一项。对基带复信号采样得到

(20)

式中:Ts为快时间采样周期;m,n为采样数。将采样信号组成快慢时间矩阵,如图2所示,先对快时间行进行快速傅里叶变换(Fast Fourier Transform,FFT)求得人体位置坐标,再提取位置坐标处的慢时间列,从而消除式(19)指数中的第一项,最后将慢时间列代入式(4)得到人体的胸腔位移信号x(t)。

图2 FMCW雷达基带信号快慢时间矩阵

1.3 脉冲超宽带雷达系统

典型的脉冲超宽带雷达由波形产生器、发射机、接收机、收发天线和信号处理器等部件组成。脉冲超宽带雷达通过测量从人体胸部反射的窄脉冲的飞行时间变化,实现对呼吸心跳的检测[15]。

超宽带雷达的发射复信号为

STX(t,nT)=p(t,nT)ej2πfct。

(21)

接收信号为

SRX(t,nT)=p(t-td,nT)ej2πfc(t-td)。

(22)

式中:p(t,nT)表示第n个持续时间为τ的窄脉冲,t为时间,T为脉冲重复周期;fc为载频。

接收信号经过混频和低通滤波后,其基带复信号表示为

SIF(t,nT)=p(t-td,nT)e-j2πfctd。

(23)

利用式(4)得到人体的胸腔位移信号x(t)。

1.4 不同雷达系统比较

雷达系统通过检测人体胸腔位移变化,实现人体呼吸和心跳信号检测,如图3所示。表1总结了不同雷达系统检测人体呼吸心跳信号的优缺点。

表1 不同雷达系统检测人体呼吸心跳信号优缺点

图3 雷达检测人体呼吸和心跳信号

1.5 雷达搭载平台

雷达搭载平台分为两种:固定式平台和移动式平台。在固定式平台中,雷达安装在桌面或墙壁上,如图4,整个检测过程中雷达静止不动。在移动式平台中,雷达安装在无人车或无人机上,如图5,雷达移动扫描某区域,并根据有效的微多普勒信号检测人体生命体征。

图4 桌面固定式平台

图5 机载移动平台[20]

2 雷达信号处理算法

现有雷达系统使用I/Q解调器,其具有提高信噪比、避免零点等优点[21]。I/Q解调器、六端口接收机都可将雷达复信号转换为I/Q信号,本文主要讨论利用的雷达采集的I/Q数据进行人体呼吸和心跳信号检测的算法流程。

以正常人的数据为例,整个人体呼吸和心跳信号检测信号处理流程如图6所示,主要分成数据预处理、人体呼吸心跳信号分离、呼吸率心率估计三部分,其中数据预处理包含去直流偏置、I/Q不平衡消除和相位解缠。

图6 人体呼吸和心跳信号检测处理流程

2.1 数据预处理

2.1.1 去直流偏置与I/Q不平衡消除

由于呼吸心跳信号比较微弱,由误差干扰引起的人体呼吸和心跳信号失真,会导致对呼吸率和心率的错误估计[22],因而需要消除数据采集过程中混入的误差。来自静止物体的反射(杂波)和硬件缺陷产生直流偏移[23],混频器之间的幅度和相位失配引起I/Q通道不平衡[24],这两部分是误差的主要来源。误差消除方法分为模拟域消除方法和数字域消除方法。在模拟域消除方面,优化硬件电路可保证I路和Q路的平衡,包括选用更好的混频器、放大器和低通滤波器[25],但这必定会提高成本,且降低电路的通用性。与之相比,数字域消除方法即利用算法消除误差,无需改变硬件电路,通用性更强。本节主要介绍现有的数字域误差消除方法。

在直流偏置补偿方面,将数据的均值作为直流偏置的估计是最常用的算法[26],但此种算法精度较低。为提高估计精度,文献[27]提出了基于压缩感知的直流偏置估计方法,将直流偏置补偿问题化为圆拟合问题,通过l1范数最小化,估计直流偏置,但在提高估计精度的同时也提高了其复杂度。为了降低复杂度,文献[28]提出使用圆心动态直流偏移跟踪方法,使用高效的梯度下降算法来实现动态直流偏移跟踪与校正,有效提高直流偏置估计的实时性。在I/Q不平衡消除方面,文献[24,29]中提出了基于椭圆拟合I/Q不平衡补偿,利用类似于Gram-Schmidt正交化的方法校正I/Q不平衡。为了提高估计的鲁棒性和精度,文献[25]进一步提出了Levenberg-Marquardt方法,通过最小化与椭圆上投影的正交距离,迭代来实现椭圆的拟合,消除I/Q不平衡。数字域误差消除方法通过建立 I/Q不平衡的校正与直流偏置补偿的模型,将其转化为优化问题,再利用相应的算法求解,实现数据补偿。通过这种方式的误差消除,会降低系统实时性,且存在残留误差,不适合实时性和精度要求较高的应用场景。

2.1.2 相位解缠

利用式(4)返回的相位值区间为[-π,π],超过区间范围的相位会被折叠到[-π,π]内,产生相位缠绕。为了恢复实际的相位,需要对相位解缠绕。

现有的相位解缠方法有差分法、微分与交叉乘法(Differential and Cross-Multiply,DACM)、扩展微分与交叉乘法。差分法通过对缠绕相位整体±π实现相位解缠。DACM解缠法对式(4)中的相位进行微分,再整体积分直接得到真实的相位。扩展DACM解缠法与DACM相比,求解表达式更简单,计算量更小。各方法的比较如表2所示。

表2 雷达相位解缠方法比较

2.2 人体呼吸心跳信号分离方法

人体呼吸和心跳信号分离方法大致可分为频域滤波信号分离方法、小波信号分解重构方法、模态分解重构方法、最优化信号处理算法。

2.2.1 频域滤波信号分离法

人体呼吸的频域范围在0.1～0.5 Hz,心跳的频域范围在0.8～2 Hz[21],根据呼吸和心跳频域范围设计带通滤波器,分离呼吸和心跳信号。设计无限长脉冲响应(Infinite Impulse Response,IIR)带通滤波器,比如,巴特沃斯型和椭圆型,实现呼吸心跳分离。有限长脉冲响应(Finite Impulse Response,FIR)滤波器阶数达到百阶,使用较少。

2.2.2 小波分解重构信号处理方法

小波分解重构的实现人体呼吸心跳分离的基本原理如图7所示,φ(t)表示雷达相位信号,LP为低通滤波器,HP为高通滤波器,↓2表示下2采样。小波变换将雷达相位信号分解为不同频率分量,在呼吸频率范围进行信号重构得到呼吸信号,在心跳频率范围进行信号重构得到心跳信号[31-32]。小波分解重构信号处理方法中小波函数的选择直接影响小波分解重构的效果,一般情况下,小波函数阶数越高,分解越精细,但是计算量越大,应用时需要根据情况选择合适的小波函数。

图7 小波分解重构原理

2.2.3 模态分解重构方法

模态分解方法分为经验模态分解(Empirical Mode Decomposition,EMD)与变分模态分解(Variational Mode Decomposition,VMD)两类。对于EMD,图8给出了分解流程:雷达相位信号被分解成N个本征模态函数(Intrinsic Mode Function,IMF),对IMF分量进行FFT,将频谱峰值在0.1～0.5 Hz范围内的IMF分量归为呼吸信号,将频谱峰值在0.8～2 Hz范围内的IMF分量归为心跳信号,其余IMF分量则视为噪声干扰,用分类的IMF分量重构得到呼吸和心跳信号。相比于EMD算法,一些改进EMD算法,如集合经验模式分解(Ensemble Empirical Mode Decomposition,EEMD)、互补集合经验模态分解(Complete Ensemble Empirical Mode Decomposition,CEEMD),在人体呼吸和心跳信号分离的应用中具有更高的精度。而VMD是一种自适应、完全非递归的模态变分方法,其通过在频域中分割信号,得到信号的有效分解分量,进而重构出人体呼吸和心跳信号[34-35]。EMD、EEMD、VMD等模态分解算法种类众多,技术相对成熟,但是其计算量较大,实时性差,应用时需要根据经验选择合适IMF分解次数,减少迭代次数,以提高算法实时性。

图8 EMD算法的流程

2.2.4 最优化信号处理算法

最优化信号处理算法基于呼吸心跳信号的频率范围不同,建立优化模型,通过求解模型实现呼吸心跳的分离。文献[11]提出了一种基于正交匹配追踪的压缩感知方法,通过设置稀疏值,得到呼吸心跳的频率。当呼吸谐波频率位于心跳频率范围内时,会对心跳估计产生干扰,使得上述算法估计误差增大。为了抑制呼吸谐波产生的干扰,文献[36]中提出了一种基于双参数最小均方滤波器的信号恢复和分离方法,在基频处提取呼吸信号,并从混合信号中剔除呼吸信号得到心跳信号,消除了呼吸谐波的干扰。虽然上述算法消除了呼吸谐波干扰,但其步骤繁琐且实时性差。文献[6]从相位和距离两种信息的组合中提取生命信号,通过距离积分消除相互干扰,用自回归方法检索生命信号,算法实现简单且实时性较高。为了提高呼吸心跳信号估计精度,文献[37]提出了外推变换和多重信号分类(Multiple Signal Classification,MUSIC)的二维参数相结合的方法,用外推变换提高距离信息的分辨率,针对目标外推变换后的幅度和相位信息,应用MUSIC算法估计出呼吸和心跳频率。将优化方法应用到呼吸心跳分离中,相较于传统的频域分析,在研究思路上有很大创新,值得进一步研究。

2.3 人体呼吸速率与心率估计方法

呼吸速率与心率的估计分为时域估计和频域估计。时域估计根据呼吸和心跳的波形变化得到呼吸心跳速率。频域估计根据呼吸心跳的频域值得到频率估计,进而转换成呼吸率与心率。

在频域估计方面,最常用的方法为最大频率估计法[38],通过对呼吸与心跳时域信号进行FFT,将信号转换到频域,其呼吸心跳的频率范围内的最大频率作为呼吸与心跳的频率估计。但是最大频率估计方法对噪声抑制较差,Welch、MUSIC等谱估计方法可抑制噪声,提高估计精度。

在时域估计方面,文献[39]中使用峰值检测算法,检测呼吸和心跳波形的极大值(或极小值)数目,再结合时间长度,得到呼吸和心跳速率的估计。当噪声干扰较大时,上述算法估计精度下降。为此,文献[40]利用峰值检测算法对同一目标进行多次测量,取多次测量的平均值作为对速率的估计,从而提高估计精度。虽然采用多次测量的均值可以提高估计精度,但是上述算法不能识别信号中的错误峰值。文献[6]使用算法识别并保存呼吸心跳波形的峰值点,根据呼吸心跳频率范围去除错误的峰值点,得到呼吸和心跳的速率估计,其估计精度更高。文献[18]对呼吸心跳信号进行自相关运算,计算出自相关运算后的峰值数,作为呼吸和心跳速率的估计。

3 研究趋势展望

由于呼吸心跳信号引起的变化十分微弱,在系统误差和噪声干扰的情况下,尽管许多文献针对不同的问题提出了相应的解决方法,但是这些方法都存在不足,仍有改进的空间。

3.1 呼吸和心跳信号分离方面

呼吸和心跳的分离主要利用其频范围不同,频率范围基于正常人体数据划定。对呼吸心跳信号进行分离时,位于呼吸心跳频域范围的干扰噪声并不能消除,这将对心跳信号产生影响,目前尚未有相关的方法用于呼吸心跳频域范围噪声的消除。由于呼吸信号引起的位移远大于心跳信号引起的位移,使得呼吸谐波的能量与心跳相近,在进行呼吸心跳的分离时,由呼吸信号产生的位于心跳频率范围的谐波分量将会干扰对心跳信号的判别。有少量文献针对这一问题提出了相应的方法,例如文献[41]中利用对心跳的二次谐波进行估计来减少呼吸谐波的干扰,但并没有学者对呼吸谐波的抑制进行深入研究并提出相应的方法。

3.2 心率和呼吸速率估计方面

现有的呼吸心率估计方法基于假定呼吸心跳在一段时间内不会发生较大的波动。在实际应用中,当发生突发情况,例如受到惊吓、刺激等,呼吸心跳会发生急剧变化。在这种情况下,现有的估计方法对呼吸心跳的估计不准确,失去实际意义,并不能及时反映真实的呼吸心跳的速率。

3.3 随机运动消除方面

现有的很多呼吸心跳检测方法建立在人体静止不动且雷达离人体距离较近又对准胸口位置的理想条件下。在实际应用中,人处于随机运动状态,且雷达不能一直正对人体,对于随机运动的消除不可或缺。现有的一些随机运动消除的方法,实现起来较为复杂,且只对剧烈运动有较好的消除效果[42]。如何简化随机运动消除的系统,对慢运动带来的干扰实现准确的消除,是未来的一个研究方向。

3.4 多人呼吸心跳检测方面

现有的呼吸心跳检测方法主要基于单目标的检测,而在实际应用中需要实现对多人的检测。现有的许多文献都是使用超宽带雷达实现对多人的检测,基于其他雷达的多人呼吸心跳检测的相关研究较少,且利用雷达进行呼吸心跳检测的可靠性仍然有待提高。因此,如何实现可靠的多目标检测,让雷达系统人体呼吸心跳信号检测用于医疗健康服务,是未来的一个主要趋势。