CN114469037B - 一种基于毫米波雷达的心率测量方法 - Google Patents

一种基于毫米波雷达的心率测量方法 Download PDF

Info

Publication number
CN114469037B
CN114469037B CN202210111994.9A CN202210111994A CN114469037B CN 114469037 B CN114469037 B CN 114469037B CN 202210111994 A CN202210111994 A CN 202210111994A CN 114469037 B CN114469037 B CN 114469037B
Authority
CN
China
Prior art keywords
vector data
line vector
radar
signal line
phase information
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
CN202210111994.9A
Other languages
English (en)
Other versions
CN114469037A (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.)
Wuhan University WHU
Original Assignee
Wuhan University WHU
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 Wuhan University WHU filed Critical Wuhan University WHU
Priority to CN202210111994.9A priority Critical patent/CN114469037B/zh
Publication of CN114469037A publication Critical patent/CN114469037A/zh
Application granted granted Critical
Publication of CN114469037B publication Critical patent/CN114469037B/zh
Active legal-status Critical Current
Anticipated expiration legal-status Critical

Links

Classifications

    • AHUMAN NECESSITIES
    • A61MEDICAL OR VETERINARY SCIENCE; HYGIENE
    • A61BDIAGNOSIS; SURGERY; IDENTIFICATION
    • A61B5/00Measuring for diagnostic purposes; Identification of persons
    • A61B5/02Detecting, measuring or recording pulse, heart rate, blood pressure or blood flow; Combined pulse/heart-rate/blood pressure determination; Evaluating a cardiovascular condition not otherwise provided for, e.g. using combinations of techniques provided for in this group with electrocardiography or electroauscultation; Heart catheters for measuring blood pressure
    • A61B5/024Detecting, measuring or recording pulse rate or heart rate
    • AHUMAN NECESSITIES
    • A61MEDICAL OR VETERINARY SCIENCE; HYGIENE
    • A61BDIAGNOSIS; SURGERY; IDENTIFICATION
    • A61B5/00Measuring for diagnostic purposes; Identification of persons
    • A61B5/05Detecting, measuring or recording for diagnosis by means of electric currents or magnetic fields; Measuring using microwaves or radio waves 
    • A61B5/0507Detecting, measuring or recording for diagnosis by means of electric currents or magnetic fields; Measuring using microwaves or radio waves  using microwaves or terahertz waves
    • AHUMAN NECESSITIES
    • A61MEDICAL OR VETERINARY SCIENCE; HYGIENE
    • A61BDIAGNOSIS; SURGERY; IDENTIFICATION
    • A61B5/00Measuring for diagnostic purposes; Identification of persons
    • A61B5/72Signal processing specially adapted for physiological signals or for diagnostic purposes
    • A61B5/7203Signal processing specially adapted for physiological signals or for diagnostic purposes for noise prevention, reduction or removal
    • AHUMAN NECESSITIES
    • A61MEDICAL OR VETERINARY SCIENCE; HYGIENE
    • A61BDIAGNOSIS; SURGERY; IDENTIFICATION
    • A61B5/00Measuring for diagnostic purposes; Identification of persons
    • A61B5/72Signal processing specially adapted for physiological signals or for diagnostic purposes
    • A61B5/7235Details of waveform analysis
    • A61B5/725Details of waveform analysis using specific filters therefor, e.g. Kalman or adaptive filters
    • AHUMAN NECESSITIES
    • A61MEDICAL OR VETERINARY SCIENCE; HYGIENE
    • A61BDIAGNOSIS; SURGERY; IDENTIFICATION
    • A61B5/00Measuring for diagnostic purposes; Identification of persons
    • A61B5/72Signal processing specially adapted for physiological signals or for diagnostic purposes
    • A61B5/7235Details of waveform analysis
    • A61B5/7253Details of waveform analysis characterised by using transforms
    • A61B5/726Details of waveform analysis characterised by using transforms using Wavelet transforms

Landscapes

  • Health & Medical Sciences (AREA)
  • Life Sciences & Earth Sciences (AREA)
  • Engineering & Computer Science (AREA)
  • Surgery (AREA)
  • General Health & Medical Sciences (AREA)
  • Biophysics (AREA)
  • Pathology (AREA)
  • Veterinary Medicine (AREA)
  • Biomedical Technology (AREA)
  • Heart & Thoracic Surgery (AREA)
  • Medical Informatics (AREA)
  • Molecular Biology (AREA)
  • Public Health (AREA)
  • Animal Behavior & Ethology (AREA)
  • Physics & Mathematics (AREA)
  • Physiology (AREA)
  • Signal Processing (AREA)
  • Artificial Intelligence (AREA)
  • Computer Vision & Pattern Recognition (AREA)
  • Psychiatry (AREA)
  • Cardiology (AREA)
  • Nuclear Medicine, Radiotherapy & Molecular Imaging (AREA)
  • Radiology & Medical Imaging (AREA)
  • Radar Systems Or Details Thereof (AREA)
  • Measuring Pulse, Heart Rate, Blood Pressure Or Blood Flow (AREA)

Abstract

本发明公开了一种基于毫米波雷达的高可靠心率测量方法。本发明对毫米波雷达采集的多组距离单元的雷达行向量数据进行归一化和累计能量计算,得到累计能量最大的雷达行向量数据;利用反正切法提取累计能量最大的行向量数据的相位信息;根据相位信息行向量进行解缠绕,一阶差分处理和距离转换得到振动信号行向量;对振动信号行向量进行MAD平滑和小波降噪,通过带通滤波器后得心跳信号行向量;对心跳信号行向量使用分帧加窗预处理后计算平均周期图功率谱,并以功率谱的重心频率作为人体心率的测量结果。其优点是:较现有技术提高了测量精度并具有强抗干扰能力,避免了由于高频噪声和现有技术的峰值检测法导致的测量结果波动大,可信度低的问题。

Description

一种基于毫米波雷达的心率测量方法
技术领域
本发明属于毫米波雷达数据处理技术领域,具体涉及一种基于毫米波雷达的心率测量方法。
背景技术
心率测量是人体生命体征检测的一个重要指标,目前在医疗领域,主要的测量人体心率的方法还是接触式测量,但是由于接触式测量必须要让设备与人体接触才能进行测量,准备工作繁琐且操作步骤复杂,在救援、探测等场合并不适用。近年来出现了许多非接触式心率测量技术,主要基于视频成像技术、雷达技术以及热成像技术等。其中,雷达技术由于具有可以穿透衣物以及障碍物的性质,在很多场合的心率测量中得到了广泛的应用。
目前应用于心率测量领域的雷达技术主要是毫米波生物雷达。心跳会引起人体胸腔的振动,通过雷达技术测量人体胸腔的振动信息,便可从中提取出心跳信息。目前已经设计出可以满足测量胸腔微小振动的高精度毫米波雷达,并且雷达结构也逐渐简化,便携性也得到了大大提高。
近年来,在毫米波生物雷达的应用领域已经实现了从发射波与反射波的合成波中提取出生命体征信号并分离出心跳信号,并通过短时傅里叶变换,峰值检测法等方法对心率进行估计,能够实现一般情况下对人体心率进行非接触式测量。而在实际应用中,由于环境噪声的影响,导致心率测量结果误差很大并且数据数据波动严重,测量结果可信度较低。
发明内容:
为了解决毫米波雷达测量心率时由于中高频噪声带来较大测量误差以及测量数据波动很大的问题,基于小波降噪技术与功率谱重心频率估计方法,本发明提出了一种基于毫米波雷达的心率测量方法。
本发明的技术方案为一种基于毫米波雷达的心率测量方法,其特征在于,包括以下步骤:
步骤1,通过毫米波雷达采集K组距离单元的雷达信号行向量数据,每组距离单元的雷达信号行向量数据的维度为N,将每组距离单元的雷达信号行向量数据进行归一化处理得到每组归一化处理后距离单元的雷达信号行向量数据,计算每组归一化处理后距离单元的雷达信号行向量数据的累计能量;在K组归一化处理后距离单元的雷达信号行向量数据的累计能量中选择累计能量最大对应的归一化处理后距离单元的雷达信号行向量数据,将累计能量最大对应的归一化处理后距离单元的雷达信号行向量数据作为人体所在部位的雷达数据,将累计能量最大对应的归一化处理后距离单元的雷达信号行向量数据进一步通过反正切法相位信号提取方法处理得到雷达相位信息行向量数据;
步骤2,将雷达相位信息行向量数据通过相位解缠绕处理得到解缠绕后的相位信息行向量数据;
步骤3,将解缠绕后的相位信息行向量数据进行一阶差分处理得到解缠绕后的相位行向差行向量数据,将到解缠绕后的相位行向差行向量数据进行距离转换得到距离转换后振动信号行向量数据;
步骤4,将距离转换后振动信号行向量数据通过MAD平滑方法进行平滑处理得到去除离群值的振动信号行向量数据;
步骤5,将去除离群值的振动信号行向量数据进行小波降噪处理得到小波降噪后振动信号行向量数据;
步骤6,将小波降噪后振动信号行向量数据通过带通滤波器处理进行滤波,得到心跳信号行向量数据;
步骤7,将心跳信号行向量数据首先进行分帧加窗得到加窗后心跳信号,根据加窗后心跳信号计算通过心跳信号的功率谱计算平均周期图功率谱,根据平均周期图功率谱计算重心频率,将重心频率作为人体估计心率;
作为优选,步骤1所述计算每组归一化处理后距离单元的雷达信号行向量数据的累计能量,具体方法为:
其中,xk(n)表示第k组归一化处理后距离单元的雷达信号行向量数据中第n个雷达数据,N为每组归一化处理后距离单元的雷达信号行向量数据的行向量数据维度,E(k)为第k组归一化处理后距离单元的雷达信号行向量数据的累计能量;
步骤1所述累计能量最大对应的归一化处理后距离单元的雷达信号行向量数据为:
xc(n),n∈[1,N]
其中,xc(n)表示累计能量最大对应的归一化处理后距离单元的雷达信号行向量数据中第n个雷达数据,即第c组归一化处理后距离单元的雷达信号行向量数据中第n个雷达数据,(xc(1),xc(2),...,xc(N))表示累计能量最大对应的归一化处理后距离单元的雷达信号行向量数,N为累计能量最大对应的归一化处理后距离单元的雷达信号行向量数据的行向量数据维;
步骤1所述反正切方法提取相位信息的处理方法,具体方法为:
其中,Ic(n)表示累计能量最大对应的归一化处理后距离单元的雷达信号行向量数据中第n个雷达数据的虚部,Rc(n)表示累计能量最大对应的归一化处理后距离单元的雷达信号行向量数据中第n个雷达数据的的实部,N为雷达相位信息行向量数据的行向量数据维,φc(n)表示雷达相位信息行向量数据中第n个雷达相位信息,(φc(1),φc(2),...,φc(N))表示雷达相位信息行向量数据;
作为优选,步骤2所述通过相位解缠绕处理得到解缠绕后的相位信息行向量数据,具体为:
其中,φc(n)表示雷达相位信息行向量数据中第n个雷达相位信息,
c(1),φc(2),...,φc(N))表示雷达相位信息行向量数据,表示解缠绕后的相位信息行向量数据中第n个解缠绕后相位信息,/>表示解缠绕后的相位信息行向量数据,N为解缠绕后的相位信息行向量数据的行向量数据维数;
作为优选,步骤3所述进行一阶差分处理得到解缠绕后的相位行向差行向量数据,具体为:
其中,表示解缠绕后的相位信息行向量数据中第n个解缠绕后相位信息,N为解缠绕后的相位行向差行向量数据的行向量数据维数,Δφc(n)为解缠绕后的相位行向差行向量数据中第n个相位数据的相位差,
表示解缠绕后的相位行向差行向量数据;
步骤3所述进行距离转换得到距离转换后振动信号行向量数据,具体为:
其中,λ为毫米波雷达发射的线性调频波的波长,Δxc(n)为距离转换后振动信号行向量数据中第n个距离转换后振动信号,(Δxc(1),Δxc(2),...,Δxc(N))表示距离转换后振动信号行向量数据,N表示距离转换后振动信号行向量数据的行向量数据维数;
作为优选,步骤5所述将去除离群值的振动信号行向量数据进行小波降噪处理为:
将去除离群值的振动信号行向量数据依次进行小波分解、阈值量化处理、信号重构;
所述小波分解是对去除离群值的振动信号行向量数据进行Q层分解,得到分解后小波信号;
所述阈值量化处理是对分解后小波信号的高频系数通过软阈值处理方式进行阈值量化处理,得到阈值量化处理之后的小波信号;
所述信号重构是对阈值量化处理之后的小波信号进行重构,得到小波降噪后振动信号行向量数据。
作为优选,步骤6所述带通滤波器为IIR带通滤波器或FIR带通滤波器,带通频率为0.8Hz-2Hz;
作为优选,步骤7所述心跳信号行向量数据首先进行分帧加窗得到加窗后心跳信号为:
将心跳信号行向量数据分为M帧分帧后的心跳信号,每帧分帧后的心跳信号的帧长为L;
所述分帧后的心跳信号表示为:
xce((m-1)*L+i),m∈[1,M],i∈[1,L]
其中,xce((m-1)*L+i)表示第m帧分帧后的心跳信号中第i个心跳信号,即心跳信号行向量数据中第(m-1)*L+i个心跳信号;
步骤7所述所述分帧加窗为:
xce1(m*L+i)=xce(m*L+i)Fw(m*L+i),m∈[1,L],i∈[1,M]
其中,xce(m*L+i),i∈[1,M]为第m帧分帧后的心跳信号,Fw(m*L+i)为第m帧分帧后的心跳信号的窗函数,xce1(m*L+i),i∈[1,M]为第m帧分帧加窗后的心跳信号;
步骤7所述根据加窗后心跳信号计算通过心跳信号的功率谱计算平均周期图功率谱,具体为:
其中,是归一化因子,Fw(n)为窗函数,ωg,g∈[1,G]为第g个离散角频率,G表示离散角频率的数量,M为心跳信号的帧数,L为每帧的帧长,xce1(m*L+n)为第m帧分帧加窗后的心跳信号中第n个心跳信号,/>为对应第g个离散角频率ωg的功率谱值;
步骤7所述通过心跳信号的功率谱计算得到重心频率,具体为:
其中,ωg为第g个离散角频率,为对应对应第g个离散角频率ωg的功率谱值,/>为重心频率;
本发明解决现有心率测量方法中由于中高频噪声导致的心率测量误差大以及心率测量结果波动大的问题。
本发明对毫米波雷达采集到的振动信号进行了基于MAD法的数据平滑处理,消除了经过一阶差分处理放大的粗大误差,避免了测量数据中的离群值对后续信号处理产生影响。
本发明对毫米波雷达采集到的振动信号进行了小波降噪处理,很好地抑制了采集到的信号中的中高频噪声,大大提升了心率测量的精度与抗干扰性,对信号进行Q层小波分解,既减小了数据处理模块的计算量,同时可以得到很好的降噪效果,采用软阈值处理方式,避免了去噪结果出现局部抖动现象。
本发明对毫米波雷达采集到的心跳信号进行了功率谱估计,采用重心频率估计的方法,能更好对一个测量周期内的平均心率进行估计,提升了测量结果的稳定性与可靠性。
附图说明:
图1:为本发明方法流程图;
具体实施方式:
为使本发明的目的、技术方案和优点更加清楚,下面将结合本发明实施例的附图,对本发明实施例的技术方案进行清楚、完整地描述。显然,所描述的实施例是本发明的一部分实施例,而不是全部的实施例。基于所描述的本发明的实施例,本领域普通技术人员在无需创造性劳动的前提下所获得的所有其他实施例,都属于本发明保护的范围。
一种基于毫米波雷达的心率测量方法,其特征在于,包括以下步骤:
步骤1,通过毫米波雷达采集K=19组距离单元的雷达信号行向量数据,每组距离单元的雷达信号行向量数据的维度为N=2048,将每组距离单元的雷达信号行向量数据进行归一化处理得到每组归一化处理后距离单元的雷达信号行向量数据,计算每组归一化处理后距离单元的雷达信号行向量数据的累计能量;在K组归一化处理后距离单元的雷达信号行向量数据的累计能量中选择累计能量最大对应的归一化处理后距离单元的雷达信号行向量数据,将累计能量最大对应的归一化处理后距离单元的雷达信号行向量数据作为人体所在部位的雷达数据,将累计能量最大对应的归一化处理后距离单元的雷达信号行向量数据进一步通过反正切法相位信号提取方法处理得到雷达相位信息行向量数据;
步骤1所述计算每组归一化处理后距离单元的雷达信号行向量数据的累计能量,具体方法为:
其中,xk(n)表示第k组归一化处理后距离单元的雷达信号行向量数据中第n个雷达数据,N=2048为每组归一化处理后距离单元的雷达信号行向量数据的行向量数据维度,E(k)为第k组归一化处理后距离单元的雷达信号行向量数据的累计能量;
步骤1所述累计能量最大对应的归一化处理后距离单元的雷达信号行向量数据为:
xc(n),n∈[1,N]
其中,xc(n)表示累计能量最大对应的归一化处理后距离单元的雷达信号行向量数据中第n个雷达数据,即第c组归一化处理后距离单元的雷达信号行向量数据中第n个雷达数据,(xc(1),xc(2),...,xc(N))表示累计能量最大对应的归一化处理后距离单元的雷达信号行向量数,N为累计能量最大对应的归一化处理后距离单元的雷达信号行向量数据的行向量数据维;
步骤1所述反正切方法提取相位信息的处理方法,具体方法为:
其中,Ic(n)表示累计能量最大对应的归一化处理后距离单元的雷达信号行向量数据中第n个雷达数据的虚部,Rc(n)表示累计能量最大对应的归一化处理后距离单元的雷达信号行向量数据中第n个雷达数据的的实部,N为雷达相位信息行向量数据的行向量数据维,φc(n)表示雷达相位信息行向量数据中第n个雷达相位信息,(φc(1),φc(2),...,φc(N))表示雷达相位信息行向量数据;
步骤2,将雷达相位信息行向量数据通过相位解缠绕处理得到解缠绕后的相位信息行向量数据;
步骤2所述通过相位解缠绕处理得到解缠绕后的相位信息行向量数据,具体为:
其中,φc(n)表示雷达相位信息行向量数据中第n个雷达相位信息,
c(1),φc(2),...,φc(N))表示雷达相位信息行向量数据,表示解缠绕后的相位信息行向量数据中第n个解缠绕后相位信息,/>表示解缠绕后的相位信息行向量数据,N为解缠绕后的相位信息行向量数据的行向量数据维数;
步骤3,将解缠绕后的相位信息行向量数据进行一阶差分处理得到解缠绕后的相位行向差行向量数据,将到解缠绕后的相位行向差行向量数据进行距离转换得到距离转换后振动信号行向量数据;
步骤3所述进行一阶差分处理得到解缠绕后的相位行向差行向量数据,具体为:
其中,表示解缠绕后的相位信息行向量数据中第n个解缠绕后相位信息,N为解缠绕后的相位行向差行向量数据的行向量数据维数,Δφc(n)为解缠绕后的相位行向差行向量数据中第n个相位数据的相位差,
表示解缠绕后的相位行向差行向量数据;
步骤3所述进行距离转换得到距离转换后振动信号行向量数据,具体为:
其中,λ为毫米波雷达发射的线性调频波的波长,Δxc(n)为距离转换后振动信号行向量数据中第n个距离转换后振动信号,(Δxc(1),Δxc(2),...,Δxc(N))表示距离转换后振动信号行向量数据,N表示距离转换后振动信号行向量数据的行向量数据维数;
步骤4,将距离转换后振动信号行向量数据通过MAD平滑方法进行平滑处理得到去除离群值的振动信号行向量数据;
步骤5,将去除离群值的振动信号行向量数据进行小波降噪处理得到小波降噪后振动信号行向量数据;
步骤5所述将去除离群值的振动信号行向量数据进行小波降噪处理为:
将去除离群值的振动信号行向量数据依次进行小波分解、阈值量化处理、信号重构;
所述小波分解是对去除离群值的振动信号行向量数据进行Q=2层分解,得到分解后小波信号;
所述阈值量化处理是对分解后小波信号的高频系数通过软阈值处理方式进行阈值量化处理,得到阈值量化处理之后的小波信号;
所述信号重构是对阈值量化处理之后的小波信号进行重构,得到小波降噪后振动信号行向量数据。
步骤6,将小波降噪后振动信号行向量数据通过带通滤波器处理进行滤波,得到心跳信号行向量数据;
步骤6所述带通滤波器为IIR带通滤波器或FIR带通滤波器,带通频率为0.8Hz-2Hz;
步骤7,将心跳信号行向量数据首先进行分帧加窗得到加窗后心跳信号,根据加窗后心跳信号计算通过心跳信号的功率谱计算平均周期图功率谱,根据平均周期图功率谱计算重心频率,将重心频率作为人体估计心率;
步骤7所述心跳信号行向量数据首先进行分帧加窗得到加窗后心跳信号为:
将心跳信号行向量数据分为M帧分帧后的心跳信号,每帧分帧后的心跳信号的帧长为L;
所述分帧后的心跳信号表示为:
xce((m-1)*L+i),m∈[1,M],i∈[1,L]
xce((m-1)*L+i),m∈[1,M],i∈[1,L]
其中,xce((m-1)*L+i)表示第m帧分帧后的心跳信号中第i个心跳信号,即心跳信号行向量数据中第(m-1)*L+i个心跳信号;
步骤7所述所述分帧加窗为:
xce1(m*L+i)=xce(m*L+i)Fw(m*L+i),m∈[1,L],i∈[1,M]
其中,xce(m*L+i),i∈[1,M]为第m帧分帧后的心跳信号,Fw(m*L+i)为第m帧分帧后的心跳信号的窗函数,xce1(m*L+i),i∈[1,M]为第m帧分帧加窗后的心跳信号;
步骤7所述根据加窗后心跳信号计算通过心跳信号的功率谱计算平均周期图功率谱,具体为:
其中,是归一化因子,Fw(n)为Hamming窗窗函数,ωg,g∈[1,G]为第g个离散角频率,G表示离散角频率的数量,M为心跳信号的帧数,L为每帧的帧长,xce1(m*L+n)为第m帧分帧加窗后的心跳信号中第n个心跳信号,/>为对应第g个离散角频率ωg的功率谱值;
步骤7所述通过心跳信号的功率谱计算得到重心频率,具体为:
其中,ωg为第g个离散角频率,为对应对应第g个离散角频率ωg的功率谱值,/>为重心频率。
应当理解的是,本说明书未详细阐述的部分均属于现有技术。
应当理解的是,上述针对实施例的描述较为详细,并不能因此而认为是对本发明专利保护范围的限制,本领域的普通技术人员在本发明的启示下,在不脱离本发明权利要求所保护的范围情况下,还可以做出替换或变形,均落入本发明的保护范围之内,本发明的请求保护范围应以所附权利要求为准。

Claims (1)

1.一种基于毫米波雷达的心率测量方法,其特征在于,包括以下步骤:
步骤1,通过毫米波雷达采集K组距离单元的雷达信号行向量数据,每组距离单元的雷达信号行向量数据的维度为N,将每组距离单元的雷达信号行向量数据进行归一化处理得到每组归一化处理后距离单元的雷达信号行向量数据,计算每组归一化处理后距离单元的雷达信号行向量数据的累计能量;在K组归一化处理后距离单元的雷达信号行向量数据的累计能量中选择累计能量最大对应的归一化处理后距离单元的雷达信号行向量数据,将累计能量最大对应的归一化处理后距离单元的雷达信号行向量数据作为人体所在部位的雷达数据,将累计能量最大对应的归一化处理后距离单元的雷达信号行向量数据进一步通过反正切法相位信号提取方法处理得到雷达相位信息行向量数据;
步骤2,将雷达相位信息行向量数据通过相位解缠绕处理得到解缠绕后的相位信息行向量数据;
步骤3,将解缠绕后的相位信息行向量数据进行一阶差分处理得到解缠绕后的相位行向差行向量数据,将到解缠绕后的相位行向差行向量数据进行距离转换得到距离转换后振动信号行向量数据;
步骤4,将距离转换后振动信号行向量数据通过MAD平滑方法进行平滑处理得到去除离群值的振动信号行向量数据;
步骤5,将去除离群值的振动信号行向量数据进行小波降噪处理得到小波降噪后振动信号行向量数据;
步骤6,将小波降噪后振动信号行向量数据通过带通滤波器处理进行滤波,得到心跳信号行向量数据;
步骤7,将心跳信号行向量数据首先进行分帧加窗得到加窗后心跳信号,根据加窗后心跳信号计算通过心跳信号的功率谱计算平均周期图功率谱,根据平均周期图功率谱计算重心频率,将重心频率作为人体估计心率;
步骤1所述计算每组归一化处理后距离单元的雷达信号行向量数据的累计能量,具体方法为:
其中,xk(n)表示第k组归一化处理后距离单元的雷达信号行向量数据中第n个雷达数据,N为每组归一化处理后距离单元的雷达信号行向量数据的行向量数据维度,E(k)为第k组归一化处理后距离单元的雷达信号行向量数据的累计能量;
步骤1所述累计能量最大对应的归一化处理后距离单元的雷达信号行向量数据为:
xc(n),n∈[1,N]
其中,xc(n)表示累计能量最大对应的归一化处理后距离单元的雷达信号行向量数据中第n个雷达数据,即第c组归一化处理后距离单元的雷达信号行向量数据中第n个雷达数据,(xc(1),xc(2),...,xc(N))表示累计能量最大对应的归一化处理后距离单元的雷达信号行向量数,N为累计能量最大对应的归一化处理后距离单元的雷达信号行向量数据的行向量数据维;
步骤1所述反正切法提取相位信息的处理方法,具体方法为:
其中,Ic(n)表示累计能量最大对应的归一化处理后距离单元的雷达信号行向量数据中第n个雷达数据的虚部,Rc(n)表示累计能量最大对应的归一化处理后距离单元的雷达信号行向量数据中第n个雷达数据的的实部,N为雷达相位信息行向量数据的行向量数据维,φc(n)表示雷达相位信息行向量数据中第n个雷达相位信息,(φc(1),φc(2),...,φc(N))表示雷达相位信息行向量数据;
步骤2所述通过相位解缠绕处理得到解缠绕后的相位信息行向量数据,具体为:
其中,φc(n)表示雷达相位信息行向量数据中第n个雷达相位信息,
c(1),φc(2),...,φc(N))表示雷达相位信息行向量数据,表示解缠绕后的相位信息行向量数据中第n个解缠绕后相位信息,/>表示解缠绕后的相位信息行向量数据,N为解缠绕后的相位信息行向量数据的行向量数据维数;
步骤3所述进行一阶差分处理得到解缠绕后的相位行向差行向量数据,具体为:
其中,表示解缠绕后的相位信息行向量数据中第n个解缠绕后相位信息,N为解缠绕后的相位行向差行向量数据的行向量数据维数,Δφc(n)为解缠绕后的相位行向差行向量数据中第n个相位数据的相位差,/>为解缠绕后的相位信息行向量数据中第n-1个解缠绕后相位信息;/>表示解缠绕后的相位行向差行向量数据;
步骤3所述进行距离转换得到距离转换后振动信号行向量数据,具体为:
其中,λ为毫米波雷达发射的线性调频波的波长,Δxc(n)为距离转换后振动信号行向量数据中第n个距离转换后振动信号,(Δxc(1),Δxc(2),...,Δxc(N))表示距离转换后振动信号行向量数据,N表示距离转换后振动信号行向量数据的行向量数据维数;
步骤5所述将去除离群值的振动信号行向量数据进行小波降噪处理为:
将去除离群值的振动信号行向量数据依次进行小波分解、阈值量化处理、信号重构;
所述小波分解是对去除离群值的振动信号行向量数据进行Q层分解,得到分解后小波信号;
所述阈值量化处理是对分解后小波信号的高频系数通过软阈值处理方式进行阈值量化处理,得到阈值量化处理之后的小波信号;
所述信号重构是对阈值量化处理之后的小波信号进行重构,得到小波降噪后振动信号行向量数据;
步骤6所述带通滤波器为IIR带通滤波器或FIR带通滤波器,带通频率为0.8Hz-2Hz;
步骤7所述心跳信号行向量数据首先进行分帧加窗得到加窗后心跳信号为:
将心跳信号行向量数据分为M帧分帧后的心跳信号,每帧分帧后的心跳信号的帧长为L;
所述分帧后的心跳信号表示为:
xce((m-1)*L+i),m∈[1,M],i∈[1,L]
其中,xce((m-1)*L+i)表示第m帧分帧后的心跳信号中第i个心跳信号,即心跳信号行向量数据中第(m-1)*L+i个心跳信号,*表示乘法符号;
步骤7所述分帧加窗为:
xce1(m*L+i)=xce(m*L+i)Fw(m*L+i),
其中,xce(m*L+i)为第m+1帧分帧后的心跳信号中第i个心跳信号,Fw(m*L+i)为第m+1帧分帧后的心跳信号的窗函数,xce1(m*L+i),为第m+1帧分帧加窗后的心跳信号中第i个心跳信号,*表示乘法符号;
步骤7所述根据加窗后心跳信号计算通过心跳信号的功率谱计算平均周期图功率谱,具体为:
其中,是归一化因子,Fw(p)为窗函数,ωg,g∈[1,G]为第g个离散角频率,G表示离散角频率的数量,M为心跳信号的帧数,L为每帧的帧长,xce1(m*L+p)为第m+1帧分帧加窗后的心跳信号中第p个心跳信号,/>为对应第g个离散角频率ωg的功率谱值,*表示乘法符号;
步骤7所述通过心跳信号的功率谱计算得到重心频率,具体为:
其中,ωg为第g个离散角频率,为对应对应第g个离散角频率ωg的功率谱值,/>为重心频率。
CN202210111994.9A 2022-01-29 2022-01-29 一种基于毫米波雷达的心率测量方法 Active CN114469037B (zh)

Priority Applications (1)

Application Number Priority Date Filing Date Title
CN202210111994.9A CN114469037B (zh) 2022-01-29 2022-01-29 一种基于毫米波雷达的心率测量方法

Applications Claiming Priority (1)

Application Number Priority Date Filing Date Title
CN202210111994.9A CN114469037B (zh) 2022-01-29 2022-01-29 一种基于毫米波雷达的心率测量方法

Publications (2)

Publication Number Publication Date
CN114469037A CN114469037A (zh) 2022-05-13
CN114469037B true CN114469037B (zh) 2024-01-12

Family

ID=81477998

Family Applications (1)

Application Number Title Priority Date Filing Date
CN202210111994.9A Active CN114469037B (zh) 2022-01-29 2022-01-29 一种基于毫米波雷达的心率测量方法

Country Status (1)

Country Link
CN (1) CN114469037B (zh)

Citations (15)

* Cited by examiner, † Cited by third party
Publication number Priority date Publication date Assignee Title
JPH08275934A (ja) * 1995-04-04 1996-10-22 Agency Of Ind Science & Technol 相関調査システム
JP2010264095A (ja) * 2009-05-15 2010-11-25 Nissan Motor Co Ltd 心拍数測定装置および心拍数測定方法
CN102499664A (zh) * 2011-10-24 2012-06-20 西双版纳大渡云海生物科技发展有限公司 基于视频图像的非接触式生命体征的检测方法及检测系统
CN103536282A (zh) * 2013-11-06 2014-01-29 中国人民解放军第三军医大学 基于Fast-ICA法的磁感应心肺活动信号分离方法
CN108158573A (zh) * 2017-12-26 2018-06-15 智慧康源(厦门)科技有限公司 基于自适应阈值小波变换的心电信号降噪方法
CN111419208A (zh) * 2020-04-13 2020-07-17 河北工业大学 一种基于加速度传感器的无束缚实时心率监测方法与系统
JP2020130745A (ja) * 2019-02-22 2020-08-31 株式会社豊田中央研究所 心身管理システム
CN111685741A (zh) * 2020-06-11 2020-09-22 中山大学 基于正交解调脉冲超宽带雷达探测人体呼吸率心率的方法
CN111743524A (zh) * 2020-06-19 2020-10-09 联想(北京)有限公司 一种信息处理方法、终端和计算机可读存储介质
CN112438707A (zh) * 2019-08-16 2021-03-05 富士通株式会社 生命体征的检测装置、方法及系统
CN112754441A (zh) * 2021-01-08 2021-05-07 杭州环木信息科技有限责任公司 一种基于毫米波的非接触式心跳检测方法
JP2021146061A (ja) * 2020-03-23 2021-09-27 株式会社リコー 生体情報取得装置、生体情報取得方法及びプログラム
CN113729677A (zh) * 2021-10-12 2021-12-03 九州云合(山东)智能科技有限公司 一种智能生命体征监护方法
CN113854990A (zh) * 2021-10-27 2021-12-31 青岛海信日立空调系统有限公司 一种心跳检测方法及装置
CN113854981A (zh) * 2021-09-24 2021-12-31 南京六季光电技术研究院有限公司 一种基于毫米波雷达的运动多人生命体征监测方法

Patent Citations (15)

* Cited by examiner, † Cited by third party
Publication number Priority date Publication date Assignee Title
JPH08275934A (ja) * 1995-04-04 1996-10-22 Agency Of Ind Science & Technol 相関調査システム
JP2010264095A (ja) * 2009-05-15 2010-11-25 Nissan Motor Co Ltd 心拍数測定装置および心拍数測定方法
CN102499664A (zh) * 2011-10-24 2012-06-20 西双版纳大渡云海生物科技发展有限公司 基于视频图像的非接触式生命体征的检测方法及检测系统
CN103536282A (zh) * 2013-11-06 2014-01-29 中国人民解放军第三军医大学 基于Fast-ICA法的磁感应心肺活动信号分离方法
CN108158573A (zh) * 2017-12-26 2018-06-15 智慧康源(厦门)科技有限公司 基于自适应阈值小波变换的心电信号降噪方法
JP2020130745A (ja) * 2019-02-22 2020-08-31 株式会社豊田中央研究所 心身管理システム
CN112438707A (zh) * 2019-08-16 2021-03-05 富士通株式会社 生命体征的检测装置、方法及系统
JP2021146061A (ja) * 2020-03-23 2021-09-27 株式会社リコー 生体情報取得装置、生体情報取得方法及びプログラム
CN111419208A (zh) * 2020-04-13 2020-07-17 河北工业大学 一种基于加速度传感器的无束缚实时心率监测方法与系统
CN111685741A (zh) * 2020-06-11 2020-09-22 中山大学 基于正交解调脉冲超宽带雷达探测人体呼吸率心率的方法
CN111743524A (zh) * 2020-06-19 2020-10-09 联想(北京)有限公司 一种信息处理方法、终端和计算机可读存储介质
CN112754441A (zh) * 2021-01-08 2021-05-07 杭州环木信息科技有限责任公司 一种基于毫米波的非接触式心跳检测方法
CN113854981A (zh) * 2021-09-24 2021-12-31 南京六季光电技术研究院有限公司 一种基于毫米波雷达的运动多人生命体征监测方法
CN113729677A (zh) * 2021-10-12 2021-12-03 九州云合(山东)智能科技有限公司 一种智能生命体征监护方法
CN113854990A (zh) * 2021-10-27 2021-12-31 青岛海信日立空调系统有限公司 一种心跳检测方法及装置

Non-Patent Citations (5)

* Cited by examiner, † Cited by third party
Title
Doppler Radar System for Long Range Detection of Respiration and Heart Rate;Ook, Park Seong;The Journal Of Korean Institute of Electromagnetic Engineering and Science;第25卷(第4期);全文 *
Spectral characteristics of heart rate fluctuations during parabolic flight;Bart Verheyden, et al;European journal of applied physiology;第95卷;全文 *
何正嘉,张西宁.现代信号处理及工程应用.西安:西安交通大学出版社,2007,210-212. *
基于视频分析的非接触式心率检测系统研究;陈怡然;中国优秀硕士学位论文全文数据库 (工程科技Ⅱ辑)(第5期);全文 *
王展.现代数字信号处理.长沙:国防科技大学出版社 ,2016,212-218. *

Also Published As

Publication number Publication date
CN114469037A (zh) 2022-05-13

Similar Documents

Publication Publication Date Title
CN109965858B (zh) 基于超宽带雷达人体生命体征检测方法和装置
Wang et al. Blind source extraction of acoustic emission signals for rail cracks based on ensemble empirical mode decomposition and constrained independent component analysis
CN110208785B (zh) 基于稳健稀疏分数阶傅立叶变换的雷达机动目标快速检测方法
CN109805931A (zh) 基于太赫兹多普勒雷达的远距离生命微动信号检测方法
Aboy et al. Adaptive modeling and spectral estimation of nonstationary biomedical signals based on Kalman filtering
CN115644840A (zh) 基于毫米波雷达的生命体征检测方法
Tang et al. ECG de-noising based on empirical mode decomposition
US20090030319A1 (en) Doppler signal spectrum calculating method and apparatus
CN116153329A (zh) 一种基于cwt-lbp的声音信号时频纹理特征提取方法
RU2367970C2 (ru) Устройство обнаружения узкополосных шумовых гидроакустических сигналов на основе вычисления интегрального вейвлет-спектра
Yamamoto et al. Non-contact heartbeat detection by heartbeat signal reconstruction based on spectrogram analysis with convolutional LSTM
CN110742593B (zh) 一种基于线谱跟踪的生命信号特征提取方法
CN114469037B (zh) 一种基于毫米波雷达的心率测量方法
CN114820847A (zh) 一种用于透射衰减超声层析成像的幅值提取方法
CN117838083A (zh) 一种基于毫米波雷达的体征快速精确检测方法
CN114325740A (zh) 一种应用于激光多普勒测速仪的速度测量方法与系统
Hu et al. The effect of multi-directional on remote heart rate measurement using PA-LI joint ICEEMDAN method with mm-wave FMCW radar
CN110420022B (zh) 一种基于双密度小波变换的p波检测方法
CN113095113B (zh) 一种用于水下目标识别的小波线谱特征提取方法及系统
RU2464588C1 (ru) Устройство обнаружения шумовых гидроакустических сигналов в виде звукоряда на основе вычисления интегрального вейвлет-спектра
CN116058818A (zh) 基于多序列woa-vmd算法的超宽带雷达心率检测方法
CN114246571A (zh) 一种基于雷达和变分模态分解算法的心率测量方法
Amazeen et al. Phase detection of R waves in noisy electrocardiograms
Kang et al. Pulse Rate Estimation using PPG Affected with Motion Artifacts Base d on VMD and Hilbert Transform
CN113189570A (zh) 一种基于复域压缩感知的阵列信号处理方法及系统

Legal Events

Date Code Title Description
PB01 Publication
PB01 Publication
SE01 Entry into force of request for substantive examination
SE01 Entry into force of request for substantive examination
GR01 Patent grant
GR01 Patent grant