CN103169449B - 呼吸信号检测方法和装置 - Google Patents

呼吸信号检测方法和装置 Download PDF

Info

Publication number
CN103169449B
CN103169449B CN201310066728.XA CN201310066728A CN103169449B CN 103169449 B CN103169449 B CN 103169449B CN 201310066728 A CN201310066728 A CN 201310066728A CN 103169449 B CN103169449 B CN 103169449B
Authority
CN
China
Prior art keywords
signal
slow
temporal filtering
time domain
breath signal
Prior art date
Legal status (The legal status is an assumption and is not a legal conclusion. Google has not performed a legal analysis and makes no representation as to the accuracy of the status listed.)
Active
Application number
CN201310066728.XA
Other languages
English (en)
Other versions
CN103169449A (zh
Inventor
乔登宇
李鑫
李烨
胡波平
孙佳平
孙东芳
张钦
Current Assignee (The listed assignees may be inaccurate. Google has not performed a legal analysis and makes no representation or warranty as to the accuracy of the list.)
Zhuhai Institute Of Advanced Technology Chinese Academy Of Sciences Co ltd
Original Assignee
Shenzhen Institute of Advanced Technology of CAS
Priority date (The priority date is an assumption and is not a legal conclusion. Google has not performed a legal analysis and makes no representation as to the accuracy of the date listed.)
Filing date
Publication date
Application filed by Shenzhen Institute of Advanced Technology of CAS filed Critical Shenzhen Institute of Advanced Technology of CAS
Priority to CN201310066728.XA priority Critical patent/CN103169449B/zh
Publication of CN103169449A publication Critical patent/CN103169449A/zh
Application granted granted Critical
Publication of CN103169449B publication Critical patent/CN103169449B/zh
Active legal-status Critical Current
Anticipated expiration legal-status Critical

Links

Landscapes

  • Radar Systems Or Details Thereof (AREA)

Abstract

本发明涉及一种呼吸信号检测方法和装置,包括以下步骤:采集原始信号;对所述原始信号进行滤波处理,且所述滤波处理中采用的滤波参数为根据预先分析得到的慢时间域呼吸信号的谐波结构确定的;计算滤波处理后的原始信号的慢时间域信号的周期图;估算滤波处理后的原始信号的慢时间域信号的噪声功率;根据所述慢时间域呼吸信号的谐波结构、慢时间域信号的周期图及慢时间域信号的噪声功率计算谐波图;判断所述谐波图是否大于预设门限,若是,则表示存在呼吸信号,若否,则表示不存在呼吸信号。上述呼吸信号检测方法和装置,采用呼吸信号的谐波结构确定滤波参数,进行滤波处理,较好的滤除环境噪声且尽可能的保留呼吸信号,提高了检测的性能。

Description

呼吸信号检测方法和装置
技术领域
本发明涉及生理信号检测领域,特别是涉及一种呼吸信号检测方法和装置。
背景技术
非接触式生命探测主要用于探测人体的呼吸运动。为了在灾后救援探测废墟中是否存在活的生命体,往往采集非接触式生命探测方式进行探测。雷达生命探测仪是一种融合雷达技术和生物医学工程技术的可穿透非金属介质的非接触、远距离的探测生命体的一种探测系统。UWB(Ultra Wide Band,超宽带)脉冲雷达的高穿透性能和高定位性能使得它在应急救援、反恐监测等领域具有广泛的应用前景。
传统的UWB脉冲雷达生命探测技术往往将人体的呼吸信号看作某个单一频点的正弦信号,当在空间的某一位置存在某一频率的呼吸信号时,相应的空间频域谱点相比于不存在呼吸信号的谱点具有较高的能量,基于该能量的显著性,完成对呼吸信号的检测。然而,在实际的应急救援、反恐监测等应用场景下,呼吸信号相对于环境噪声干扰等很微弱,采用传统的UWB脉冲雷达生命探测检测呼吸信号时,因仅将呼吸信号视作单一频点的正弦信号,忽略了环境噪声干扰,导致检测的性能较差。
发明内容
基于此,有必要针对传统的呼吸信号检测性能差的问题,提供一种能提高检测性能的呼吸信号检测方法。
此外,还有必要针对传统的呼吸信号检测性能差的问题,提供一种能提高检测性能的呼吸信号检测装置。
一种呼吸信号检测方法,包括以下步骤:
采集原始信号;
对所述原始信号进行滤波处理,且所述滤波处理中采用的滤波参数为根据预先分析得到的慢时间域呼吸信号的谐波结构确定的;
计算滤波处理后的原始信号的慢时间域信号的周期图;
估算滤波处理后的原始信号的慢时间域信号的噪声功率;
根据所述慢时间域呼吸信号的谐波结构、慢时间域信号的周期图及慢时间域信号的噪声功率计算谐波图;
判断所述谐波图是否大于预设门限,若是,则表示存在呼吸信号,若否,则表示不存在呼吸信号。
一种呼吸信号检测装置,包括:
采集模块,用于采集原始信号;
噪声抑制模块,用于对所述原始信号进行滤波处理,且所述滤波处理中采用的滤波参数为根据预先分析得到的慢时间域呼吸信号的谐波结构确定的;
周期图计算模块,用于计算滤波处理后的原始信号的慢时间域信号的周期图;
噪声功率估算模块,用于计算滤波处理后的原始信号的慢时间域信号的噪声功率;
谐波图计算模块,用于根据所述慢时间域呼吸信号的谐波结构、慢时间域信号的周期图及慢时间域信号的噪声功率计算谐波图;
判断模块,用于判断所述谐波图是否大于预设门限,若是,则表示存在呼吸信号,若否,则表示不存在呼吸信号。
上述呼吸信号检测方法和装置,采用慢时间域呼吸信号的谐波结构确定滤波参数,根据确定后的滤波参数滤波器对采集的原始信号进行滤波处理,若原始信号中存在呼吸信号,则通过此滤波处理可较好的滤除环境噪声且尽可能的保留呼吸信号,根据呼吸信号谐波结构及计算得到的慢时间域信号周期图和噪声功率计算谐波图,通过谐波图与预设门限进行判断,提高了检测的性能,更加有效的检测呼吸信号是否存在。
附图说明
图1为一个实施例中呼吸信号检测方法的流程示意图;
图2为一个实施例中预先分析得到慢时间域呼吸信号的谐波结构的流程示意图;
图3为胸腔表面与雷达间距离随时间变化的示意图;
图4为仿真得到的回波矩阵实验图;
图5为采用离散傅立叶级数分析呼吸信号谐波结构的结果图;
图6为慢时间域呼吸信号的谐波结构图;
图7A为快时间域滤波器的幅度谱响应图;
图7B为快时间域滤波器的冲击响应图;
图8为存在呼吸信号的雷达回波实验数据图;
图9为快时间域滤波后的输出矩阵图;
图10A为带通慢时间滤波器的幅度谱图;
图10B为带通慢时间域滤波器的冲击响应图;
图11为带通慢时间域滤波处理后的输出矩阵图;
图12A为高通慢时间滤波器的幅度谱图;
图12B为高通慢时间域滤波器的冲击响应图;
图13为高通慢时间域滤波处理后的输出矩阵图;
图14为带通慢时间域滤波器输出矩阵的周期图;
图15为谐波图;
图16为两种检测方法的检测性能对比示意图;
图17为一个实施例中呼吸信号检测装置的结构示意图;
图18为另一个实施例中呼吸信号检测装置的结构示意图。
具体实施方式
下面结合具体的实施例及附图对呼吸信号检测方法和装置的技术方案进行详细的描述,以使其更加清楚。
如图1所示,在一个实施例中,一种呼吸信号检测方法,包括以下步骤:
步骤S110,采集原始信号。
具体的,可通过雷达探测仪在灾后救援现场或反恐场景中采集原始信号,该原始信号中可能包括呼吸信号。该原始信号组成雷达回波矩阵,该雷达回波矩阵中通常包含大量噪声,杂波和干扰。雷达回波矩阵可由式(1)表示:
式(1)中,表示雷达回波矩阵;r[m,n]表示由于胸腔表面运动产生的信号;c[n]表示由背景环境产生的杂波;w[m,n]表示系统噪声,通常符合高斯白噪声模型;d[m]表示由采样引入的不稳定的直流基线;l[m,n]表示由于采样幅值不稳定性引入的线性趋势,在慢时间域上通常显示为基线漂移。
步骤S120,对该原始信号进行滤波处理,且该滤波处理中采用的滤波参数为根据预先分析得到的慢时间域呼吸信号的谐波结构确定的。
对原始信号进行滤波处理的目的是为了抑制c[n]、w[m,n]、d[m]和l[m,n],保留r[m,n]。
在一个实施例中,在步骤S110之前,还包括步骤:预先分析得到慢时间域呼吸信号的谐波结构。
如图2所示,在一个实施例中,预先分析得到慢时间域呼吸信号的谐波结构包括如下步骤:
步骤S210,获取并根据脉冲回波延迟时间、慢时间域采样周期和快时间域采样周期构建离散时间回波矩阵的仿真模型,并设置该离散时间回波矩阵的仿真参数。
由于UWB脉冲雷达回波数据存在两个维度,即快时间维度和慢时间维度,因而通常被称为回波矩阵。UWB雷达对人体呼吸运动的检测主要是对人体胸腔周期起伏运动的检测。当UWB雷达对人体探测时,人体胸腔的起伏运动会导致胸腔表面与雷达之间的径向距离的周期性变化。这一变化导致回波脉冲的延迟时间的周期波动。当人体吸气时,胸腔膨胀,胸腔表面与雷达间的距离变短,经胸腔表面反射的脉冲回波的延迟时间缩短;反之,当人体呼气,胸腔收缩,胸腔表面与雷达间的距离变长,经胸腔表面反射的脉冲回波的延迟时间延长。回波脉冲延迟时间的周期波动反应在二维回波数据(快时间域和慢时间域)上是在慢时间域上信号具有周期性。为此,呼吸信号是指因胸腔运动产生回波脉冲延迟时间的周期波动,再由该周期波动导致二维回波数据在慢时间域上形成的周期信号。其中,快时间域是指单个脉冲发射及回波时间形成的时间域,慢时间域是指每个脉冲发射时间点形成的时间域。
本实施例中,采用高阶余弦信号对人体呼吸导致胸腔表面的起伏运动进行建模。设被测目标与雷达间的距离为d0,人体胸腔起伏运动的幅度为B,人体的呼吸率为fr,则人体胸腔表面与雷达间的距离d(t)随时间t的变化可用式(2)表示:
d(t)=d0-B×(cosπfrt)u    (2)
式(2)中,u表示高阶余弦信号的阶数。通常u取6时,(2)式可以较好的近似人体胸腔表面与雷达间距离的变化过程。图3给出了当d0=4m,B=0.01m,fr=0.25Hz,u=6的情况下,根据该模型给出的胸腔表面与雷达间距离随时间变化的示意图,其中,m为米,Hz为赫兹。图3中横坐标为时间,单位为秒,纵坐标为人体胸腔表面与雷达间的距离,单位为米。
根据d(t)计算反射的脉冲回波的延迟时间td(t),如式(3)所示:
t d ( t ) = 2 × d ( t ) c - - - ( 3 )
式(3)中,c表示光速。
采用p(t)表示雷达发射机产生的脉冲波形。设雷达慢时间域采样周期为Tst,雷达快时间域采样周期为Tft。根据脉冲回波延迟时间td(t)、慢时间域采样周期Tst和快时间域采样周期Tft构建离散时间回波矩阵的仿真模型,可得到离散时间回波矩阵的仿真模块如式(4)所示:
r[m,n]=p(nTft-td(mTst))    (4)
本实施案例中,雷达发射机采用的是一阶高斯脉冲,脉冲波形p(t)表达式如下:
p ( t ) = - A × 2 t σ 2 × e - t 2 σ 2 - - - ( 5 )
式(5)中,A为脉冲幅度因子,σ为时间因子。因脉冲幅度因子的取值并不影响后续分析,故在此不作特殊设置。时间因子σ在本实施案例中取值128ps,ps为皮秒,即10-12秒。把式(5)代入式(4)可得在本实施案例下离散时间回波矩阵表达式:
r [ m , n ] = - A 2 × ( d 0 - B ( cos π f r m T st ) u ) σ 2 × e - ( n T ft - 2 ( d 0 - B ( cos π f r m T st ) u ) c ) 2 σ 2
其中,Tst表示慢时间域采样周期,Tft表示雷达快时间域采样周期。
设置该离散时间回波矩阵的仿真参数,Tst和Tft的取值应该确保采样在慢时间域和快时间域上均满足采样定理。为了确保离散时间回波矩阵在慢时间域上是周期信号,在仿真参数设置时应确保为fr的整数倍。通常,为fr的100倍。该值较为合理,一方面既保证为fr的整数倍,另一方面,较大,可确保在慢时间域上采样满足采样定理。Tft的取值可通过脉冲波形p(t)计算脉冲频谱,根据频谱确定脉冲的频带宽度,根据频带宽度确定Tft,保证在快时间域上采样满足采样定理。
步骤S220,根据该仿真模型和仿真参数获取离散时间回波矩阵的仿真结果。
具体的,采用matlab软件作为仿真平台。根据仿真模型,设置仿真参数,经matlab进行仿真,以获得图4所示的回波矩阵,图4中,d0=4m,B=10mm,fr=0.25Hz,A=10-9,σ=128ps,Tst=0.04s,Tft=10ps,其中,m为米,mm为毫米,Hz为赫兹,ps为皮秒,横坐标为快时间域变量n,纵坐标为慢时间域变量m,且右边为与左边对应的灰度刻度条(即图中的-6至6)。
步骤S230,根据该仿真结果提取慢时间域呼吸信号。
具体的,对于回波矩阵r[m,n],rm[n]用于表示单路快时间域信号,其中,快时间域离散时间变量n可看作该函数的自变量,慢时间域离散时间变量m可看作固定参数;类似的,rn[m]用于表示单路慢时间域信号,其中,慢时间域离散时间变量m可看作该函数的自变量,快时间域离散时间变量n可看作固定参数。
在雷达信号处理中,为了抑制杂波,任何信号的直流分量都会被移除,因此,对各路rn[m]取均值,获得不含直流分量的慢时间域信号计算均值的方法是:由于rn[m]的周期已知(在本实施案例中为100),首先仿真出一个周期的rn[m],然后对一个周期的仿真数据求和再除以周期即得到rn[m]的均值。典型的慢时间域呼吸信号指的是在回波矩阵中具有最大功率的单路慢时间域呼吸信号其中,各路慢时间域信号功率可通过对信号平方再取平均即可获得。
在本实施案例中,按照上述方法进行计算,得到本案例典型的慢时间域呼吸信号n=2663。
步骤S240,利用离散傅立叶级数分析该慢时间域呼吸信号,得到该慢时间域呼吸信号的谐波结构向量。
首先定义呼吸信号的谐波结构为一向量,该向量的第n个元素对应第n次谐波,该向量的第n个元素的取值为第n次谐波的功率在信号总功率中的比重。由于典型的慢时间域呼吸信号是离散周期信号,可以采用离散傅立叶级数对其进行分析,获得其谐波结构。
在本实施案例中,由图4中所给参数设置知,fst=100×fr,fst为慢时间域采样频率,且这意味着,这些呼吸信号经过慢时间域采样后都是周期的,且离散周期TM=100。
对于离散周期信号,采用离散傅立叶级数进行分析。离散傅立叶级数的计算公式为:
R n k [ k ] = Σ m = 0 T M - 1 r h t [ m ] W T M km , W T M = e - j ( 2 π T M ) , k=0,1,...,TM-1   (6)
在本实施案例中,对做离散傅立叶级数的计算,结果如图5所示,图5中离散傅立叶级数k=0,1,...,9。
由于直流分量已被移除, 对应第k次谐波,且基于此,式(7)给出了计算第k次谐波的功率在总功率中所占比重Ratio(k)的计算方法:
Ratio ( k ) = 2 × R h t [ k ] 2 Σ k = 1 T M | R h t [ k ] | 2 - - - ( 7 )
根据图5给出的离散傅立叶级数计算呼吸信号的谐波结构向量。图6用画图的方式显示了典型的慢时间域呼吸信号的谐波结构向量,横坐标为谐波阶次,纵坐标为功率比重。从图6可以看出,在本实施案例中,1次谐波和2次谐波占据了信号的绝大部分功率,因此,在谐波结构向量中只保留1次和2次谐波,得到典型的慢时间域呼吸信号的谐波结构向量:
上述通过离散傅立叶级数分析得到了呼吸信号的谐波结构,方便结合呼吸信号的谐波结构设置相应的滤波参数,以提高滤波的质量,避免滤除存在的呼吸信号。
进一步的,在一个实施例中,步骤S120包括以下步骤:
(1)对原始信号中每一路信号进行快时间域滤波处理得到相应路的快时间域滤波输出信号,并组成快时间域滤波输出矩阵。
具体的,快时间域滤波处理采用快时间域滤波器进行滤波处理。该快时间域滤波器为带通线性相位有限脉冲响应数字滤波器。原始信号为包含噪声的回波矩阵,回波矩阵的每一路快时间信号经过该带通线性相位有限脉冲响应数字滤波器获得相应的一路输出信号,多路输出信号组成快时间域滤波输出矩阵。该带通线性相位有限脉冲响应数字滤波器的频率响应的带宽应与雷达脉冲频带相匹配。例如,雷达脉冲的带宽是从1GHz到3GHz,则为了确保雷达脉冲在幅度谱上不出现失真,滤波器的频率响应在1GHz到3GHz上应尽可能保持近似为1,而在其它频段,为了抑制噪声,应尽可能接近为0。此外,快时间域滤波处理采用带通线性相位有限脉冲响应数字滤波器是为了保持回波脉冲的时域波形。
快时间域滤波输出矩阵rftF[m,n]的表达式如下:
式(8)中,为雷达回波矩阵(即原始信号组成的回波矩阵),hftF[n]为带通线性相位有限脉冲响应数字滤波器的单位脉冲响应;*n表示快时间域的时域卷积计算。
在本实施案例中,雷达脉冲的带宽为0.45GHz到3.555GHz,而由于使用的天线的带宽为0.9GHz到5GHz,综合考虑,回波脉冲的带宽为0.9GHz到3.555GHz。采用经典的线性相位有限脉冲响应数字滤波器设计技术,得到本实施案例中所使用的快时间域滤波参数如图7A和图7B所示,图7A表示快时间域滤波器的幅度谱响应,图7A中横坐标为频率,纵坐标为幅值;图7B表示快时间域滤波器的冲击响应,图7B中横坐标为时间,纵坐标为幅值。该快时间域滤波器的阶数为83,群延时为42。
通过该快时间域滤波,回波矩阵r中的d[m]能被有效抑制。
在本实施案例中,通过实验,采集了一套实验数据,如图8所示,为存在呼吸信号的雷达回波实验数据示意图,且右边为与左边对应的灰度刻度条(即图8中的50.3至51)。在图8中,各点的数值都在50.3至51.1之间变化,d[m]>50,其值远远大于其它信号的数值,为了抑制d[m],对实验数据在快时间域上滤波,输出如图9所示。图9为快时间域滤波后的输出矩阵图,且右边为与左边对应的灰度刻度条(图9中的-0.15至0.2)。从图9中可以看到,各点数值在-0.2至0.2之间变化,d[m]得到有效抑制。
(2)对该快时间域滤波输出矩形的每一路快时间域滤波输出信号进行带通慢时间域滤波处理得到相应路的带通慢时间域滤波输出信号,并组成慢时间域滤波输出矩阵,其中,该带通慢时间域滤波处理中采用的滤波参数为根据该慢时间域呼吸信号的谐波结构确定的。
具体的,慢时间域滤波处理包括带通慢时间域滤波处理和高通慢时间域滤波处理,慢时间域滤波处理采用慢时间域滤波器进行,该慢时间域滤波器由2个线性相位有限脉冲响应数字滤波器组成,即带通慢时间域滤波器和高通慢时间域滤波器。快时间域滤波输出矩阵的每一路快时间域滤波输出信号经过慢时间域滤波器获得相应的输出信号。
带通慢时间域滤波器作用是保留回波矩阵中的r[m,n],抑制c[n]、w[m,n]和l[m,n]。根据呼吸信号的谐波结构确定带通慢时间域滤波器的参数,即根据呼吸信号的谐波结构确定带通慢时间域滤波处理中采用的滤波参数。在本实施案例中,在呼吸信号的谐波结构分析中,已知呼吸信号主要包含1次谐波和2次谐波。基于此,考虑到人体呼吸运动的频率主要在0.2Hz到0.5Hz之间浮动,呼吸信号的频带主要为0.2Hz至1Hz。保留r[m,n]对于慢时间域滤波来说就是保留呼吸信号,因而,带通慢时间域滤波处理中的滤波参数包括滤波频率,该滤波频率为0.2赫兹到1赫兹。
为此,带通慢时间域滤波器的频率响应在0.2Hz到1Hz上应尽可能保持近似为1,而在其它频段,为了抑制噪声,应尽可能接近为0。这里需要指出,c[n]和l[m,n]在慢时间域上通常集中在极低频段,因而通过该慢时间域滤波后均会被有效抑制。
此外,慢时间域滤波处理采用带通线性相位有限脉冲响应数字滤波器是为了保持呼吸信号的时域波形。
带通慢时间域滤波输出矩阵rstF[m,n]的计算如下式所示:
rstF[m,n]=rftF[m,n]*mhstF[m]    (9)
式(9)中,hstF[m]为带通慢时间域滤波器的单位脉冲响应;*m表示慢时间域的时域卷积计算。
本实施案例中所使用的带通慢时间域滤波器参数如图10A和图10B所示。图10A表示带通慢时间滤波器的幅度谱,图10A中横坐标为频率,纵坐标为幅值;图10B表示带通慢时间域滤波器的冲击响应,图10B中横坐标为时间,纵坐标为幅值。该带通慢时间域滤波器的阶数为401,群延时为200。
本实施案例中,经带通慢时间域滤波器处理后的回波矩阵(即实验输出数据)如图11所示,且右边为与左边对应的灰度刻度条(即图11中的-0.015至0.015)。
(3)对该快时间域滤波输出矩阵的每一路快时间域滤波输出信号进行高通慢时间域滤波处理得到相应路的高通慢时间域滤波输出信号,并组成高通慢时间域滤波输出矩阵。
高通慢时间域滤波器主要是为了获取一个只包含w[m,n],用于接下来估计噪声功率。高通慢时间域滤波器的性能要求就是下截止频率足够高,能够有效抑制r[m,n]、c[n]和l[m,n]。因此根据前述分析,在本实施案例中,高通慢时间域滤波器的下截止频率至少应大于1Hz。
高通慢时间域滤波输出矩阵的计算如下式所示:
式(10)中,为高通慢时间域滤波器的单位脉冲响应;*m表示慢时间域的时域卷积计算。
本实施案例中所使用的带通慢时间域滤波器参数如图12A和图12B所示。图12A表示高通慢时间滤波器的幅度谱,图12A中横坐标为频率,纵坐标为幅值;图12B表示高通慢时间域滤波器的冲击响应,图12B中横坐标为时间,纵坐标为幅值。该高通慢时间域滤波器的阶数为201,群延时为100。
本实施案例中,经高通慢时间域滤波器处理后的回波矩阵(即实验输出数据)如图13所示,且右边为与左边对应的灰度刻度条(即图13中的-1.5至1.5)。
步骤S130,计算滤波处理后的原始信号的慢时间域信号的周期图。
具体的,步骤130具体包括:根据带通慢时间域滤波输出矩阵的每一路带通慢时间域滤波输出信号进行周期图计算,得到各路慢时间域信号的周期图。
在本实施案例中,基于带通慢时间域滤波器实验输出矩阵,计算各路慢时间域信号的周期图In[k],计算公式如下:
I n [ k ] = 1 M | R stF n [ k ] | 2 , k=0,1,2,...,M-1    (11)
式(11)中,M为带通慢时间域滤波器输出矩阵的单路慢时间域信号的数据长度。在本实施案例中M=2000,为慢时间域信号的快速傅立叶变换(FFT)。图14中给出了带通慢时间域滤波器输出矩阵的周期图,且右边为与左边对应的灰度刻度条(即图14中的0.005至0.02)。图14中,在n=175附近,周期图In[k]存在较显著的一次谐波和二次谐波信号。
步骤S140,估算滤波处理后的原始信号的慢时间域信号的噪声功率。
具体的,步骤S140具体包括:根据该高通慢时间域滤波输出矩阵的每一路高通慢时间域滤波输出信号估算各路慢时间域信号的噪声功率。
根据高通慢时间域滤波输出矩阵的每一路慢时间域信号估计各路慢时间域信号的噪声功率,计算公式如下:
式(12)中,Pn表示第n路慢时间域信号的噪声功率估计值,为高通慢时间域滤波输出矩阵,为高通慢时间域滤波器的单位脉冲响应。
步骤S150,根据该慢时间域呼吸信号的谐波结构、慢时间域信号的周期图及慢时间域信号的噪声功率计算谐波图。
具体的,谐波图Hn[k]计算公式如下:
H n [ k ] = Σ j = 1 J h ( a j × I n [ j × k ] ) P n - - - ( 13 )
其中,Jh表示周期信号所含谐波阶次的最大值;aj等于谐波结构向量的第j个分量的取值,In[k]为各路慢时间域信号的周期图。
在本实施案例中,由谐波结构分析知,Jh=2,a1=0.8562,a2=0.1404。
根据谐波图定义,由实验数据计算所得谐波图如图15所示,横坐标为快时间域,纵坐标为频率,且右边为与左边对应的灰度刻度条(即图15中的5至40)。
步骤S160,判断该谐波图是否大于预设门限,若是,则执行步骤S170,若否,则执行步骤S180。
基于谐波图,由呼吸信号的二元假设判决,判决公式如下:
H n [ k ] > H 0 < H 1 r - - - ( 14 )
式(14)中,H1表示存在呼吸信号,H0表示不存在呼吸信号。当Hn[k]大于预设门限r时,H1为判决结果,反之,Hn[k]小于等于预设门限r时,H0为判决结果。当预设门限r过高,会降低检测算法的检测概率,而预设门限r过低,又会增大检测算法的虚警率。
在图15所示的本实施案例的谐波图中,通过适当的选取门限值,例如r=10,将能成功检测到呼吸信号。
步骤S170,表示存在呼吸信号。
步骤S180,表示不存在呼吸信号。
上述呼吸信号检测方法,采用慢时间域呼吸信号的谐波结构确定滤波参数,根据确定后的滤波参数滤波器对采集的原始信号进行滤波处理,若原始信号中存在呼吸信号,则通过此滤波处理可较好的滤除环境噪声且尽可能的保留呼吸信号,根据呼吸信号谐波结构及计算得到的慢时间域信号周期图和噪声功率计算谐波图,通过谐波图与预设门限进行判断,提高了检测的性能,更加有效的检测呼吸信号是否存在。
进一步的,在一个实施例中,当存在呼吸信号时,上述呼吸信号检测方法还包括:根据该谐波图的最大值估算目标距离和/或目标呼吸率。
具体的,在实际系统中,由于各种电路元器件和数字滤波器的非理想特性,脉冲往往会被展宽,同时,考虑到在应用场景中不可避免的多经效应,单个人体的呼吸运动在回波矩阵的快时间域中会被展宽和延伸。这一点在图15所示的谐波图中很好的得到了验证。但是,考虑到人体胸腔表面的直接反射信号往往强于其它多径信号,当判决存在呼吸信号时,对目标距离和呼吸率进行估算的步骤包括:
(1)获取Hn[k]的最大值,其对应的快时间域变量n用表示,其对应的离散时间频率用表示。
(2)将作为呼吸率的估计值。
(3)估计目标距离d。
具体的,按传统的雷达脉冲信号计算距离的方式结合快时间域滤波器对距离估计的影响可得,目标距离计算公式如下:
d = c &times; T ft &times; ( n ^ - delay ftF ) 2 - - - ( 15 )
式(15)中,delayftF为快时间域滤波器的群延时,c表示光速,Tft为雷达快时间域采样周期,表示快时间域变量n。目标距离是指被测目标与检测装置之间的距离。
在本实施案例中,图15中小箭头“”所指向的Hn[k]即为谐波图的最大值,根据该最大值即可按式(15)估计呼吸率和目标距离。
上述呼吸信号检测方法对呼吸信号的检测性能明显优于传统的UWB脉冲探测方法(即将呼吸信号视作单一频点的正弦信号的检测方法)对呼吸信号检测。如图16所示为两种检测方法的检测性能对比示意图,横坐标为信噪比(单位dB),纵坐标为检测概率。图16所示,在高斯白噪声信道模型下,通过MonteCarlo仿真,获得两种方法在虚警率Pfa为0.01和0.001的情况下对于不同的信噪比(SNR)所能达到的检测概率的对比图,图中黑色圆点(即图16中“●”)连线为传统UWB探测方法在Pfa为0.01时不同信噪比所达到的检测概率,黑色矩形点(即图16中“■”)连线为本发明中呼吸信号检测方法在Pfa为0.01时不同信噪比所达到的检测概率,无色圆点(即图16中“○”)连线为传统UWB探测方法在Pfa为0.001时不同信噪比所达到的检测概率,无色矩形点(即图16中“□”)连线为本发明中呼吸信号检测方法在Pfa为0.001时不同信噪比所达到的检测概率,从图16中可以看出本发明中呼吸信号检测方法的检测性能明显优于传统的UWB探测方法。
如图17所示,在一个实施例中,一种呼吸信号检测装置,包括采集模块110、噪声抑制模块120、周期图计算模块130、噪声功率估算模块140、谐波图计算模块150和判断模块160。其中:
采集模块110用于采集原始信号。
具体的,采样模块110为雷达探测仪或其他信号检测装置,本实施例中,通过雷达探测仪在灾后救援现场或反恐场景中采集原始信号,该原始信号中可能包括呼吸信号。该原始信号组成雷达回波矩阵,该雷达回波矩阵中通常包含大量噪声,杂波和干扰。雷达回波矩阵可由式(1)表示:
式(1)中,表示雷达回波矩阵;r[m,n]表示由于胸腔表面运动产生的信号;c[n]表示由背景环境产生的杂波;w[m,n]表示系统噪声,通常符合高斯白噪声模型;d[m]表示由采样引入的不稳定的直流基线;l[m,n]表示由于采样幅值不稳定性引入的线性趋势,在慢时间域上通常显示为基线漂移。
噪声抑制模块120,用于对该原始信号进行滤波处理,且该滤波处理中采用的滤波参数为根据预先分析得到的慢时间域呼吸信号的谐波结构确定的。
对原始信号进行滤波处理的目的是为了抑制c[n]、w[m,n]、d[m]和l[m,n],保留r[m,n]。
进一步的,上述呼吸信号检测装置还包括呼吸信号谐波结构分析模块170,用于预先分析得到慢时间域呼吸信号的谐波结构。
如图17中,呼吸信号谐波结构分析模块170包括建模单元171、仿真单元173、提取单元175和分析单元177。其中:
建模单元171,用于获取并根据脉冲回波延迟时间、慢时间域采样周期和快时间域采样周期构建离散时间回波矩阵的仿真模型,并设置该离散时间回波矩阵的仿真参数。
本实施例中,采用高阶余弦信号对人体呼吸导致胸腔表面的起伏运动进行建模。设被测目标与雷达间的距离为d0,人体胸腔起伏运动的幅度为B,人体的呼吸率为fr,则人体胸腔表面与雷达间的距离d(t)随时间t的变化可用式(2)表示:
d(t)=d0-B×(cosπfrt)u    (2)
式(2)中,u表示高阶余弦信号的阶数。通常u取6时,(2)式可以较好的近似人体胸腔表面与雷达间距离的变化过程。图3给出了当d0=4m,B=0.01m,fr=0.25Hz,u=6的情况下,根据该模型给出的胸腔表面与雷达间距离随时间变化的示意图,其中,m为米,Hz为赫兹。图3中横坐标为时间,单位为秒,纵坐标为人体胸腔表面与雷达间的距离,单位为米。
根据d(t)计算反射的脉冲回波的延迟时间td(t),如式(3)所示:
t d ( t ) = 2 &times; d ( t ) c - - - ( 3 )
式(3)中,c表示光速。
采用p(t)表示雷达发射机产生的脉冲波形。设雷达慢时间域采样周期为Tst,雷达快时间域采样周期为Tft。根据脉冲回波延迟时间td(t)、慢时间域采样周期Tst和快时间域采样周期Tft构建离散时间回波矩阵的仿真模型,可得到离散时间回波矩阵的仿真模块如式(4)所示:
r[m,n]=p(nTft-td(mTst))    (4)
本实施案例中,雷达发射机采用的是一阶高斯脉冲,脉冲波形p(t)表达式如下:
p ( t ) = - A &times; 2 t &sigma; 2 &times; e - t 2 &sigma; 2 - - - ( 5 )
式(5)中,A为脉冲幅度因子,σ为时间因子。因脉冲幅度因子的取值并不影响后续分析,故在此不作特殊设置。时间因子σ在本实施案例中取值128ps,ps为皮秒,即10-12秒。把式(5)代入式(4)可得在本实施案例下离散时间回波矩阵表达式:
r [ m , n ] = - A 2 &times; ( d 0 - B ( cos &pi; f r m T st ) u ) &sigma; 2 &times; e - ( n T ft - 2 ( d 0 - B ( cos &pi; f r m T st ) u ) c ) 2 &sigma; 2
其中,Tst表示慢时间域采样周期,Tft表示雷达快时间域采样周期。
设置该离散时间回波矩阵的仿真参数,Tst和Tft的取值应该确保采样在慢时间域和快时间域上均满足采样定理。为了确保离散时间回波矩阵在慢时间域上是周期信号,在仿真参数设置时应确保为fr的整数倍。通常,为fr的100倍。该值较为合理,一方面既保证为fr的整数倍,另一方面,较大,可确保在慢时间域上采样满足采样定理。Tft的取值可通过脉冲波形p(t)计算脉冲频谱,根据频谱确定脉冲的频带宽度,根据频带宽度确定Tft,保证在快时间域上采样满足采样定理。
仿真单元173,用于根据该仿真模型和仿真参数获取离散时间回波矩阵的仿真结果。
具体的,采用matlab软件作为仿真平台。根据仿真模型,设置仿真参数,经matlab进行仿真,以获得图4所示的回波矩阵,图4中,d0=4m,B=10mm,fr=0.25Hz,A=10-9,σ=128ps,Tst=0.04s,Tft=10ps,其中,m为米,mm为毫米,Hz为赫兹,ps为皮秒,横坐标为快时间域变量n,纵坐标为慢时间域变量m,且右边为与左边对应的灰度刻度条(即图中的-6至6)。
提取单元175,用于根据该仿真结果提取慢时间域呼吸信号。
具体的,对于回波矩阵r[m,n],rm[n]用于表示单路快时间域信号,其中,快时间域离散时间变量n可看作该函数的自变量,慢时间域离散时间变量m可看作固定参数;类似的,rn[m]用于表示单路慢时间域信号,其中,慢时间域离散时间变量m可看作该函数的自变量,快时间域离散时间变量n可看作固定参数。
为了抑制杂波,任何信号的直流分量都会被移除,因此,对各路rn[m]取均值,获得不含直流分量的慢时间域信号计算均值的方法是:由于rn[m]的周期已知(在本实施案例中为100),首先仿真出一个周期的rn[m],然后对一个周期的仿真数据求和再除以周期即得到rn[m]的均值。典型的慢时间域呼吸信号指的是在回波矩阵中具有最大功率的单路慢时间域呼吸信号其中,各路慢时间域信号功率可通过对信号平方再取平均即可获得。
在本实施案例中,按照上述方法进行计算,得到本案例典型的慢时间域呼吸信号n=2663。
分析单元177,用于利用离散傅立叶级数分析该慢时间域呼吸信号,得到该慢时间域呼吸信号的谐波结构向量。
首先定义呼吸信号的谐波结构为一向量,该向量的第n个元素对应第n次谐波,该向量的第n个元素的取值为第n次谐波的功率在信号总功率中的比重。由于典型的慢时间域呼吸信号是离散周期信号,可以采用离散傅立叶级数对其进行分析,获得其谐波结构。
在本实施案例中,由图4中所给参数设置知,fst=100×fr,fst为慢时间域采样频率,且这意味着,这些呼吸信号经过慢时间域采样后都是周期的,且离散周期TM=100。
对于离散周期信号,采用离散傅立叶级数进行分析。离散傅立叶级数的计算公式为:
R n k [ k ] = &Sigma; m = 0 T M - 1 r h t [ m ] W T M km , W T M = e - j ( 2 &pi; T M ) , k=0,1,...,TM-1    (6)
在本实施案例中,对做离散傅立叶级数的计算,结果如图5所示,图5中离散傅立叶级数k=0,1,...,9。
由于直流分量已被移除, 对应第k次谐波,且基于此,式(7)给出了计算第k次谐波的功率在总功率中所占比重Ratio(k)的计算方法:
Ratio ( k ) = 2 &times; R h t [ k ] 2 &Sigma; k = 1 T M | R h t [ k ] | 2 - - - ( 7 )
根据图5给出的离散傅立叶级数计算呼吸信号的谐波结构向量。图6用画图的方式显示了典型的慢时间域呼吸信号的谐波结构向量,横坐标为谐波阶次,纵坐标为功率比重。从图6可以看出,在本实施案例中,1次谐波和2次谐波占据了信号的绝大部分功率,因此,在谐波结构向量中只保留1次和2次谐波,得到典型的慢时间域呼吸信号的谐波结构向量:
上述通过离散傅立叶级数分析得到了呼吸信号的谐波结构,方便结合呼吸信号的谐波结构设置相应的滤波参数,以提高滤波的质量,避免滤除存在的呼吸信号。
进一步的,噪声抑制模块120包括:
快时间域滤波模块121,用于对该原始信号中每一路信号进行快时间域滤波处理得到相应路的快时间域滤波输出信号,并组成快时间域滤波输出矩阵。
具体的,快时间域滤波模块121为快时间域滤波器。该快时间域滤波器为带通线性相位有限脉冲响应数字滤波器。原始信号为包含噪声的回波矩阵,回波矩阵的每一路快时间信号经过该带通线性相位有限脉冲响应数字滤波器获得相应的一路输出信号,多路输出信号组成快时间域滤波输出矩阵。该带通线性相位有限脉冲响应数字滤波器的频率响应的带宽应与雷达脉冲频带相匹配。例如,雷达脉冲的带宽是从1GHz到3GHz,则为了确保雷达脉冲在幅度谱上不出现失真,滤波器的频率响应在1GHz到3GHz上应尽可能保持近似为1,而在其它频段,为了抑制噪声,应尽可能接近为0。此外,快时间域滤波器采用带通线性相位的滤波器是为了保持回波脉冲的时域波形。
快时间域滤波输出矩阵rftF[m,n]的表达式如下:
式(8)中,为雷达回波矩阵(即原始信号组成的回波矩阵),hftF[n]为带通线性相位有限脉冲响应数字滤波器的单位脉冲响应;*n表示快时间域的时域卷积计算。
在本实施案例中,雷达脉冲的带宽为0.45GHz到3.555GHz,而由于使用的天线的带宽为0.9GHz到5GHz,综合考虑,回波脉冲的带宽为0.9GHz到3.555GHz。采用经典的线性相位有限脉冲响应数字滤波器设计技术,得到本实施案例中所使用的快时间域滤波参数如图7A和图7B所示,图7A表示快时间域滤波器的幅度谱响应,图7A中横坐标为频率,纵坐标为幅值;图7B表示快时间域滤波器的冲击响应,图7B中横坐标为时间,纵坐标为幅值。该快时间域滤波器的阶数为83,群延时为42。
通过该快时间域滤波模块121滤波后回波矩阵中的d[m]能被有效抑制。
在本实施案例中,通过实验,采集了一套实验数据,如图8所示,为存在人体呼吸运动的雷达回波实验数据示意图。在图8中,各点的数值都在50.3至51.1之间变化,d[m]>50,其值远远大于其它信号的数值,为了抑制d[m],对实验数据在快时间域上滤波,输出如图9所示。图9为快时间域滤波后的实验输出数据示意图。从图9中可以看到,各点数值在-0.2至0.2之间变化,d[m]得到有效抑制。
慢时间域滤波模块123,用于对快时间域滤波输出矩阵的每一路快时间域滤波输出信号滤波获得相应的输出信号。慢时间域滤波模块123由2个线性相位有限脉冲响应数字滤波器组成,即带通慢时间域滤波器1231和高通慢时间域滤波器1233。
带通慢时间域滤波器1231,用于对快时间域滤波输出矩形的每一路快时间域滤波输出信号进行带通慢时间域滤波处理得到相应路的带通慢时间域滤波输出信号,并组成慢时间域滤波输出矩阵,其中,该带通慢时间域滤波处理中采用的滤波参数为根据该慢时间域呼吸信号的谐波结构确定的。
带通慢时间域滤波器作用是保留回波矩阵中的r[m,n],抑制c[n]、w[m,n]和l[m,n]。根据呼吸信号的谐波结构确定带通慢时间域滤波器的参数,即根据呼吸信号的谐波结构确定带通慢时间域滤波处理中采用的滤波参数。在本实施案例中,在呼吸信号的谐波结构分析中,已知呼吸信号主要包含1次谐波和2次谐波。基于此,考虑到人体呼吸运动的频率主要在0.2Hz到0.5Hz之间浮动,呼吸信号的频带主要为0.2Hz至1Hz。保留r[m,n]对于慢时间域滤波来说就是保留呼吸信号,因而,带通慢时间域滤波处理中的滤波参数包括滤波频率,该滤波频率为0.2赫兹到1赫兹。
为此,带通慢时间域滤波器的频率响应在0.2Hz到1Hz上应尽可能保持近似为1,而在其它频段,为了抑制噪声,应尽可能接近为0。这里需要指出,c[n]和l[m,n]在慢时间域上通常集中在极低频段,因而通过该慢时间域滤波后均会被有效抑制。
此外,慢时间域滤波处理采用带通线性相位有限脉冲响应数字滤波器是为了保持呼吸信号的时域波形。
带通慢时间域滤波输出矩阵rstF[m,n]的计算如下式所示:
rstF[m,n]=rftF[m,n]*mhstF[m]    (9)
式(9)中,hstF[m]为带通慢时间域滤波器的单位脉冲响应;*m表示慢时间域的时域卷积计算。
本实施案例中所使用的带通慢时间域滤波器参数如图10A和图10B所示。图10A表示带通慢时间滤波器的幅度谱,图10A中横坐标为频率,纵坐标为幅值;图10B表示带通慢时间域滤波器的冲击响应,图10B中横坐标为时间,纵坐标为幅值。该带通慢时间域滤波器的阶数为401,群延时为200。
本实施案例中,经带通慢时间域滤波器处理后的回波矩阵如图11所示。
高通慢时间域滤波器1233,用于对该快时间域滤波输出矩阵的每一路快时间域滤波输出信号进行高通慢时间域滤波处理得到相应路的高通慢时间域滤波输出信号,并组成高通慢时间域滤波输出矩阵。
高通慢时间域滤波器主要是为了获取一个只包含w[m,n],用于接下来估计噪声功率。高通慢时间域滤波器的性能要求就是下截止频率足够高,能够有效抑制r[m,n]、c[n]和l[m,n]。因此根据前述分析,在本实施案例中,高通慢时间域滤波器的下截止频率至少应大于1Hz。
高通慢时间域滤波输出矩阵的计算如下式所示:
式(10)中,为高通慢时间域滤波器的单位脉冲响应;*m表示慢时间域的时域卷积计算。
本实施案例中所使用的带通慢时间域滤波器参数如图12A和图12B所示。图12A表示高通慢时间滤波器的幅度谱,图12A中横坐标为频率,纵坐标为幅值;图12B表示高通慢时间域滤波器的冲击响应,图12B中横坐标为时间,纵坐标为幅值。该高通慢时间域滤波器的阶数为201,群延时为100。
本实施案例中,经高通慢时间域滤波器处理后的回波矩阵如图13所示。
周期图计算模块130,用于计算滤波处理后的原始信号的慢时间域信号的周期图。
具体的,周期图计算模块130还用于根据带通慢时间域滤波输出矩阵的每一路带通慢时间域滤波输出信号进行周期图计算,得到各路慢时间域信号的周期图。
在本实施案例中,基于带通慢时间域滤波器实验输出矩阵,计算各路慢时间域信号的周期图In[k],计算公式如下:
I n [ k ] = 1 M | R stF n [ k ] | 2 , k=0,1,2,...,M-1   (11)
式(11)中,M为带通慢时间域滤波器实验输出矩阵的单路慢时间域信号的数据长度。在本实施案例中M=2000,为慢时间域信号的快速傅立叶变换(FFT)。图14中给出了带通慢时间域滤波器输出矩阵的周期图。图14中,在n=175附近,周期图In[k]存在较显著的一次谐波和二次谐波信号。
噪声功率估算模块140,用于计算滤波处理后的原始信号的慢时间域信号的噪声功率。
具体的,噪声功率估算模块140还用于根据该高通慢时间域滤波输出矩阵的每一路高通慢时间域滤波输出信号估算各路慢时间域信号的噪声功率。
根据高通慢时间域滤波输出矩阵的每一路慢时间域信号估计各路慢时间域信号的噪声功率,计算公式如下:
式(12)中,Pn表示第n路慢时间域信号的噪声功率估计值,为高通慢时间域滤波输出矩阵,为高通慢时间域滤波器的单位脉冲响应。
谐波图计算模块150,用于根据所述慢时间域信号的周期图及慢时间域信号的噪声功率计算谐波图。
具体的,谐波图Hn[k]计算公式如下:
H n [ k ] = &Sigma; j = 1 J h ( a j &times; I n [ j &times; k ] ) P n - - - ( 13 )
其中,Jh表示周期信号所含谐波阶次的最大值;aj等于谐波结构向量的第j个分量的取值,In[k]为各路慢时间域信号的周期图。
在本实施案例中,由谐波结构分析知,Jh=2,a1=0.8562,a2=0.1404。
根据谐波图定义,由实验数据计算所得谐波图如图15所示。
判断模块160,用于判断所述谐波图是否大于预设门限,若是,则表示存在呼吸信号,若否,则表示不存在呼吸信号。
基于谐波图,由呼吸信号的二元假设判决,判决公式如下:
H n [ k ] > H 0 < H 1 r - - - ( 14 )
式(14)中,H1表示存在呼吸信号,H0表示不存在呼吸信号。当Hn[k]大于预设门限r时,H1为判决结果,反之,Hn[k]小于等于预设门限r时,H0为判决结果。当预设门限r过高,会降低检测算法的检测概率,而预设门限r过低,又会增大检测算法的虚警率。
在图15所示的本实施案例的谐波图中,通过适当的选取门限值,例如r=10,将能成功检测到呼吸信号。
上述呼吸信号检测装置,采用慢时间域呼吸信号的谐波结构确定滤波参数,根据确定后的滤波参数滤波器对采集的原始信号进行滤波处理,若原始信号中存在呼吸信号,则通过此滤波处理可较好的滤除环境噪声且尽可能的保留呼吸信号,根据呼吸信号谐波结构及计算得到的慢时间域信号周期图和噪声功率计算谐波图,通过谐波图与预设门限进行判断,提高了检测的性能,更加有效的检测呼吸信号是否存在。
进一步的,如图18所示,在一个实施例中,上述呼吸信号检测装置,除了包括采集模块110、噪声抑制模块120、周期图计算模块130、噪声功率估算模块140、谐波图计算模块150和判断模块160,还包括距离和呼吸率估算模块180。其中:
距离和呼吸率估算模块180,用于在判断存在呼吸信号时,根据所述谐波图的最大值估算目标距离和/或目标呼吸率。
具体的,当判决存在呼吸信号时,距离和呼吸率估算模块180,对目标距离和呼吸率进行估算包括以下过程:
(1)获取Hn[k]的最大值,其对应的快时间域变量n用表示,其对应的离散时间频率用表示。
(2)将作为呼吸率的估计值。
(3)估计目标距离d。
具体的,按传统的雷达脉冲信号计算距离的方式结合快时间域滤波器对距离估计的影响可得,目标距离计算公式如下:
d = c &times; T ft &times; ( n ^ - delay ftF ) 2 - - - ( 15 )
式(15)中,,delayftF为快时间域滤波器的群延时,c表示光速,Tft为雷达快时间域采样周期,表示快时间域变量n。目标距离是指被测目标与检测装置之间的距离。
在本实施案例中,图15中小箭头“”所指向的Hn[k]即为谐波图的最大值,根据该最大值即可按式(15)估计呼吸率和目标距离。
以上所述实施例仅表达了本发明的几种实施方式,其描述较为具体和详细,但并不能因此而理解为对本发明专利范围的限制。应当指出的是,对于本领域的普通技术人员来说,在不脱离本发明构思的前提下,还可以做出若干变形和改进,这些都属于本发明的保护范围。因此,本发明专利的保护范围应以所附权利要求为准。

Claims (8)

1.一种呼吸信号检测方法,包括以下步骤:
采集原始信号;
对所述原始信号进行滤波处理,且所述滤波处理中采用的滤波参数为根据预先分析得到的慢时间域呼吸信号的谐波结构确定的;
计算滤波处理后的原始信号的慢时间域信号的周期图;
估算滤波处理后的原始信号的慢时间域信号的噪声功率;
根据所述慢时间域呼吸信号的谐波结构、慢时间域信号的周期图及慢时间域信号的噪声功率计算谐波图;
判断所述谐波图是否大于预设门限,若是,则表示存在呼吸信号,若否,则表示不存在呼吸信号;
所述呼吸信号检测方法还包括预先分析得到慢时间域呼吸信号的谐波结构的步骤,所述预先分析得到慢时间域呼吸信号的谐波结构包括:
获取并根据脉冲回波延迟时间、慢时间域采样周期和快时间域采样周期构建离散时间回波矩阵的仿真模型,并设置所述离散时间回波矩阵的仿真参数;
根据所述仿真模型和仿真参数获取离散时间回波矩阵的仿真结果;
根据所述仿真结果提取慢时间域呼吸信号;
利用离散傅立叶级数分析所述慢时间域呼吸信号,得到所述慢时间域呼吸信号的谐波结构向量;
所述对原始信号进行滤波处理的步骤包括:
对所述原始信号中每一路信号进行快时间域滤波处理得到相应路的快时间域滤波输出信号,并组成快时间域滤波输出矩阵;
对所述快时间域滤波输出矩形的每一路快时间域滤波输出信号进行带通慢时间域滤波处理得到相应路的带通慢时间域滤波输出信号,并组成慢时间域滤波输出矩阵,其中,所述带通慢时间域滤波处理中采用的滤波参数为根据所述慢时间域呼吸信号的谐波结构确定的;
对所述快时间域滤波输出矩阵的每一路快时间域滤波输出信号进行高通慢时间域滤波处理得到相应路的高通慢时间域滤波输出信号,并组成高通慢时间域滤波输出矩阵。
2.根据权利要求1所述的呼吸信号检测方法,其特征在于,所述带通慢时间域滤波处理中的滤波参数包括滤波频率,所述滤波频率为0.2赫兹到1赫兹。
3.根据权利要求1所述的呼吸信号检测方法,其特征在于,所述计算滤波处理后的原始信号的慢时间域信号的周期图的步骤包括:
根据所述带通慢时间域滤波输出矩阵的每一路带通慢时间域滤波输出信号进行周期图计算,得到各路慢时间域信号的周期图;
所述估算滤波处理后的原始信号的慢时间域信号的噪声功率的步骤包括:
根据所述高通慢时间域滤波输出矩阵的每一路高通慢时间域滤波输出信号估算各路慢时间域信号的噪声功率。
4.根据权利要求1所述的呼吸信号检测方法,其特征在于,当存在呼吸信号时,所述呼吸信号检测方法还包括:
根据所述谐波图的最大值估算目标距离和/或目标呼吸率。
5.一种呼吸信号检测装置,其特征在于,包括:
采集模块,用于采集原始信号;
噪声抑制模块,用于对所述原始信号进行滤波处理,且所述滤波处理中采用的滤波参数为根据预先分析得到的慢时间域呼吸信号的谐波结构确定的;
周期图计算模块,用于计算滤波处理后的原始信号的慢时间域信号的周期图;
噪声功率估算模块,用于计算滤波处理后的原始信号的慢时间域信号的噪声功率;
谐波图计算模块,用于根据所述慢时间域呼吸信号的谐波结构、慢时间域信号的周期图及慢时间域信号的噪声功率计算谐波图;
判断模块,用于判断所述谐波图是否大于预设门限,若是,则表示存在呼吸信号,若否,则表示不存在呼吸信号;
所述装置还包括:
呼吸信号谐波结构分析模块,用于预先分析得到慢时间域呼吸信号的谐波结构;所述呼吸信号谐波结构分析模块包括:
建模单元,用于获取并根据脉冲回波延迟时间、慢时间域采样周期和快时间域采样周期构建离散时间回波矩阵的仿真模型,并设置所述离散时间回波矩阵的仿真参数;
仿真单元,用于根据所述仿真模型和仿真参数获取离散时间回波矩阵的仿真结果;
提取单元,用于根据所述仿真结果提取慢时间域呼吸信号;
分析单元,用于利用离散傅立叶级数分析所述慢时间域呼吸信号,得到所述慢时间域呼吸信号的谐波结构向量;
所述噪声抑制模块包括:
快时间域滤波模块,用于对所述原始信号中每一路信号进行快时间域滤波处理得到相应路的快时间域滤波输出信号,并组成快时间域滤波输出矩阵;
慢时间域滤波模块,包括:
带通慢时间域滤波器,用于对所述快时间域滤波输出矩形的每一路快时间域滤波输出信号进行带通慢时间域滤波处理得到相应路的带通慢时间域滤波输出信号,并组成慢时间域滤波输出矩阵,其中,所述带通慢时间域滤波处理中采用的滤波参数为根据所述慢时间域呼吸信号的谐波结构确定的;
高通慢时间域滤波器,用于对所述快时间域滤波输出矩阵的每一路快时间域滤波输出信号进行高通慢时间域滤波处理得到相应路的高通慢时间域滤波输出信号,并组成高通慢时间域滤波输出矩阵。
6.根据权利要求5所述的呼吸信号检测装置,其特征在于,所述带通慢时间域滤波单元的滤波频率为0.2赫兹到1赫兹。
7.根据权利要求5所述的呼吸信号检测装置,其特征在于,所述周期图计算模块还用于根据所述带通慢时间域滤波输出矩阵的每一路带通慢时间域滤波输出信号进行周期图计算,得到各路慢时间域信号的周期图;
所述噪声功率估算模块还用于根据所述高通慢时间域滤波输出矩阵的每一路高通慢时间域滤波输出信号估算各路慢时间域信号的噪声功率。
8.根据权利要求5所述的呼吸信号检测装置,其特征在于,当存在呼吸信号时,所述呼吸信号检测装置还包括:
距离和呼吸率估算模块,用于根据所述谐波图的最大值估算目标距离和/或目标呼吸率。
CN201310066728.XA 2013-03-01 2013-03-01 呼吸信号检测方法和装置 Active CN103169449B (zh)

Priority Applications (1)

Application Number Priority Date Filing Date Title
CN201310066728.XA CN103169449B (zh) 2013-03-01 2013-03-01 呼吸信号检测方法和装置

Applications Claiming Priority (1)

Application Number Priority Date Filing Date Title
CN201310066728.XA CN103169449B (zh) 2013-03-01 2013-03-01 呼吸信号检测方法和装置

Publications (2)

Publication Number Publication Date
CN103169449A CN103169449A (zh) 2013-06-26
CN103169449B true CN103169449B (zh) 2014-12-10

Family

ID=48629968

Family Applications (1)

Application Number Title Priority Date Filing Date
CN201310066728.XA Active CN103169449B (zh) 2013-03-01 2013-03-01 呼吸信号检测方法和装置

Country Status (1)

Country Link
CN (1) CN103169449B (zh)

Families Citing this family (5)

* Cited by examiner, † Cited by third party
Publication number Priority date Publication date Assignee Title
CN105640502B (zh) * 2015-12-29 2018-08-14 深圳先进技术研究院 一种呼吸信号检测方法和装置
CN106991708B (zh) * 2017-04-27 2020-04-14 飞依诺科技(苏州)有限公司 超声多普勒血流成像的处理方法及处理系统
WO2020042444A1 (zh) * 2018-08-31 2020-03-05 深圳迈睿智能科技有限公司 人体存在探测器及其人体存在探测方法
CN109521422B (zh) * 2018-10-15 2020-06-09 中国人民解放军第四军医大学 一种基于雷达信号的多目标生命探测方法及探测雷达
CN113893033B (zh) * 2021-07-01 2023-05-12 中国科学院苏州生物医学工程技术研究所 一种肺部经皮穿刺导航方法及系统

Family Cites Families (5)

* Cited by examiner, † Cited by third party
Publication number Priority date Publication date Assignee Title
CN101324666A (zh) * 2007-06-16 2008-12-17 电子科技大学 隐蔽目标生命迹象探测方法及隐蔽目标探测装置
CN102018503B (zh) * 2010-10-21 2012-12-12 中国科学院深圳先进技术研究院 生命探测雷达中的呼吸与心跳信号的提取方法及装置
EP2517621A1 (en) * 2011-04-29 2012-10-31 Nederlandse Organisatie voor toegepast -natuurwetenschappelijk onderzoek TNO A radar apparatus for detecting multiple life-signs of a subject, a method and a computer program product
CN102298822A (zh) * 2011-08-15 2011-12-28 东北大学秦皇岛分校 自助式无源生命呼救方法及系统
CN102841341B (zh) * 2012-09-03 2014-08-27 深圳先进技术研究院 一种脉冲雷达动目标检测方法

Also Published As

Publication number Publication date
CN103169449A (zh) 2013-06-26

Similar Documents

Publication Publication Date Title
CN108226892B (zh) 一种基于深度学习的复杂噪声环境下的雷达信号恢复方法
CN103169449B (zh) 呼吸信号检测方法和装置
CN106093868B (zh) 一种基于双源ir-uwb生物雷达的强反射杂波消除方法
US10401479B2 (en) Remote sensing of human breathing at a distance
Li et al. Advanced signal processing for vital sign extraction with applications in UWB radar detection of trapped victims in complex environments
CN102512157B (zh) 基于模型的动态心电图t波交替定量分析方法
CN106468770B (zh) K分布杂波加噪声下的近最优雷达目标检测方法
CN113440120B (zh) 一种基于毫米波雷达的人员呼吸心跳检测方法
CN106019254B (zh) 一种uwb冲击生物雷达多人体目标距离向分离辨识方法
CN107167802A (zh) 一种基于超宽带雷达的呼吸信号检测算法
CN102488517A (zh) 一种检测脑电信号中爆发抑制状态的方法以及装置
Yao et al. An adaptive seismic signal denoising method based on variational mode decomposition
CN113640792A (zh) 一种基于机器学习的车内生命体毫米波雷达检测方法
Nejadgholi et al. Time-frequency based contactless estimation of vital signs of human while walking using PMCW radar
CN205758555U (zh) 一种人体呼吸的检测装置
CN106199539A (zh) 基于白化滤波器的地杂波抑制方法
Sarkar et al. Accurate sensing of multiple humans buried under rubble using IR-UWB SISO radar during search and rescue
CN114027805B (zh) 一种基于微波雷达的生命体征测量方法及系统
CN110269642A (zh) 基于分数阶傅里叶变换和小波变换的多普勒心率估计方法
Abdulatif et al. Power-based real-time respiration monitoring using FMCW radar
CN115345216A (zh) 一种融合先验信息的fmcw雷达干扰消除方法
Li et al. Macro-motion detection using ultra-wideband impulse radar
Kakouche et al. Respiration Rate Estimation Using IR-UWB Radar Based on the Second Maximum Eigenvalues
CN110850488A (zh) 一种冲激脉冲超声波的生命探测方法
CN102353943A (zh) 一种抑制射频干扰的方法以及雷达式生命探测仪

Legal Events

Date Code Title Description
C06 Publication
PB01 Publication
C10 Entry into substantive examination
SE01 Entry into force of request for substantive examination
C14 Grant of patent or utility model
GR01 Patent grant
TR01 Transfer of patent right
TR01 Transfer of patent right

Effective date of registration: 20240201

Address after: 519085 101, Building 5, Longyuan Smart Industrial Park, No. 2, Hagongda Road, Tangjiawan Town, High-tech Zone, Zhuhai City, Guangdong Province

Patentee after: ZHUHAI INSTITUTE OF ADVANCED TECHNOLOGY CHINESE ACADEMY OF SCIENCES Co.,Ltd.

Country or region after: China

Address before: 1068 No. 518055 Guangdong city in Shenzhen Province, Nanshan District City Xili University School Avenue

Patentee before: SHENZHEN INSTITUTES OF ADVANCED TECHNOLOGY CHINESE ACADEMY OF SCIENCES

Country or region before: China