CN-122004823-A - Thoracic cavity displacement waveform extraction method based on self-adaptive time-frequency analysis and two-dimensional filtering
Abstract
The application belongs to the technical field of radar signal processing, and provides a thoracic cavity displacement waveform extraction method based on self-adaptive time-frequency analysis and two-dimensional filtering, which comprises the following steps: firstly, acquiring a radar complex baseband sampling sequence containing human vital sign information on a target range gate, initializing self-adaptive time-frequency analysis related parameters, performing self-adaptive time-frequency analysis on the complex baseband sampling sequence, recording an optimal window length sequence, calculating a representative width parameter, calculating an initial time-frequency mask by using SO-CFAR detection, constructing a ridge line neighborhood time-frequency mask, calculating a final time-frequency mask, acquiring a denoising complex time spectrum by using the initial time-frequency mask, calculating a reconstructed time-domain denoising complex baseband signal, and finally extracting a chest displacement waveform. According to the application, through self-adaptive time-frequency analysis and two-dimensional SO-CFAR filtering, noise is effectively suppressed, chest displacement waveforms are accurately extracted, and robustness and accuracy of vital sign monitoring under the condition of low signal-to-noise ratio are improved.
Inventors
- YAO SHUAI
- YANG FENGSHENG
- WANG ZIXU
- Nie Yunteng
- TAN QIANNING
- LIU JIEYIN
- LUO HAORAN
- Dai Xuejie
- SUN YAN
- WU QISONG
Assignees
- 东南大学
Dates
- Publication Date
- 20260512
- Application Date
- 20260325
Claims (9)
- 1. The thoracic cavity displacement waveform extraction method based on self-adaptive time-frequency analysis and two-dimensional filtering is characterized by comprising the following steps of: s1, acquiring a radar complex baseband sampling sequence z (N), n=1, 2, N, containing vital sign information of a detected human body on a target distance gate to be processed; S2, initializing a plurality of preset parameters of self-adaptive time-frequency analysis; S3, calculating a complex time spectrum Z (L, k), a time-frequency power spectrum P (L, k) and an optimal window length sequence W opt (L) of self-adaptive time-frequency analysis based on a complex baseband sampling sequence Z (n), wherein l=1, 2; s4, estimating a representative width parameter B ref of the target to be processed based on a time-frequency power spectrum P (l, k) of the self-adaptive time-frequency analysis; S5, calculating an initial time-frequency mask M 0 (l, k) by utilizing minimum-selection constant false alarm SO-CFAR detection according to the time-frequency power spectrum P (l, k) and the representative width parameter B ref ; S6, constructing a ridge line neighborhood time-frequency mask M r (l, k) based on P (l, k), calculating a final time-frequency mask M (l, k) by combining the initial time-frequency mask M 0 (l, k), and obtaining a denoising complex time-frequency spectrum according to the final time-frequency mask M (l, k) and the complex time-frequency spectrum Z (l, k) ; S7, denoising complex time spectrum And an optimal window length sequence W opt (l), calculating a reconstructed time-domain denoised complex baseband signal z d (n d ); S8, extracting instantaneous phase of the reconstructed time domain denoising complex baseband signal z d (n d ), performing phase unwrapping, trend term removal and scaling coefficient conversion, and extracting chest cavity displacement waveform 。
- 2. The method according to claim 1, wherein the method for acquiring the Lei Dafu baseband sampling sequence z (n) in S1 specifically comprises the following steps: Acquiring a Lei Dafu baseband sampling sequence z (N) containing vital sign information on a range gate of a target to be processed by a frequency modulation continuous wave FMCW radar data acquisition device, wherein n=1, 2, & gt, N is a time index of the Lei Dafu baseband sampling sequence z (N), N is the total number of sampling points of a complex baseband sequence, and acquiring a slow time sampling frequency f s and a radar working wavelength lambda corresponding to the Lei Dafu baseband sampling sequence z (N), wherein the slow time sampling frequency and the radar working wavelength respectively meet that f s is more than or equal to 10Hz and less than or equal to 50 Hz, lambda is more than or equal to 3mm and less than or equal to 13 mm, and N is the sampling point number of not less than 30 seconds, namely 。
- 3. The method according to claim 2, wherein the step of S2 initializing a plurality of preset parameters of the adaptive time-frequency analysis specifically comprises the steps of: S201, initializing a minimum candidate window length W min and a maximum candidate window length W max of self-adaptive time-frequency analysis: , w min and W max are positive integers, and the units are sampling points; s202, initializing window length search steps delta W, wherein delta W=1, and the unit is a sampling point; S203, initializing a window length candidate set omega of self-adaptive time-frequency analysis: ; Wherein { · } is a set; s204, initializing the number Q of elements in the set omega: ; S205, initializing the q-th candidate window length W q in a window length candidate set omega of self-adaptive time-frequency analysis: ; wherein q is the window length index in the set Ω; s206, initializing frame shift H, wherein H=1, and the unit is a sampling point; S207, initializing total time frame number L: ; S208, initializing a time frame index L, namely sliding and framing Lei Dafu baseband sampling sequences z (n) according to frame shift H to obtain l=1, 2, and L; S209, initializing the number N fft :W max ≤N fft of points of Fourier transform to be less than or equal to N, wherein N fft is a positive integer; s210, initializing a discrete frequency resolution Δf: ; s211, the structure is centered at 0 Discrete frequency set with uniform sampling in range : ; Where k is the sampling frequency index; S212, initializing a trade-off coefficient for balancing the energy focusing index and the window length continuity penalty term : ; S213, initialize the optimal window length sequence W opt (l):W opt (L) = 0,l =1, 2.
- 4. A method according to claim 3, wherein S3, based on the complex baseband sampling sequence Z (n), calculates a complex time spectrum Z (l, k), a time-frequency power spectrum P (l, k) and an optimal window length sequence W opt (l) for the adaptive time-frequency analysis, specifically as follows: s301, adopting a Hanning window, sequentially constructing a short-time analysis window function with the length of W q , multiplying W q sampling points corresponding to the first time frame by the analysis window function point by point, adding zero to N fft points at the tail end of a multiplied sequence, carrying out N fft point discrete Fourier transform to obtain a discrete Fourier transform result, carrying out spectrum centering on the discrete Fourier transform result to obtain a candidate complex time spectrum corresponding to the first time frame and the q candidate window length And calculates the corresponding candidate time-frequency power spectrum : ; Wherein, the Is an absolute function, k is a bin index, k=1, 2,..n fft ; s302, candidate time-frequency power spectrum based on the ith time frame and the qth candidate window length Calculating a normalization factor C q (l): ; Wherein, the Candidate time-frequency power spectrum representing the length of the q-th candidate window in the l-th time frame Searching for a maximum in the frequency dimension; According to candidate time-frequency power spectrum And a normalization factor C q (l), calculating a normalized candidate time-frequency power spectrum : ; S303, calculating an energy focusing index J E,q (l) corresponding to the first time frame and the q candidate window length according to the fourth moment of the power in the frequency dimension: ; S304, calculating a window length continuity penalty term J C,q (l) corresponding to the first time frame and the q candidate window length: ; S305, carrying out weighted synthesis on the energy focusing index and the window length continuity penalty term to obtain a comprehensive evaluation function J q (l) corresponding to the ith time frame and the (q) th candidate window length: ; s306, for the first time frame, at all candidate window lengths Selecting a candidate window length for maximizing the comprehensive evaluation function J q (l) as an optimal window length W opt (l) of the time frame, and recording an index q opt (l) corresponding to the optimal window length; S307, spectrum at candidate complex Selecting complex time spectrum corresponding to the optimal window length as complex time spectrum of self-adaptive time-frequency analysis under the first time frame : ; And a corresponding time-frequency power spectrum P (l, k): 。
- 5. The method according to claim 4, wherein the step of estimating the representative width parameter B ref of the target to be processed based on the adaptive time-frequency analysis time-frequency power spectrum P (i, k) is performed as follows: S401, initializing a local energy threshold coefficient eta, wherein eta is a positive real number between 0 and 1; S402, for each time frame index l, calculating a first frame time-frequency power spectrum maximum value P max (l): ; Wherein, the Represents the maximum value of P (l, k) in the range of 1≤k≤N fft on the first frame; Calculating a ridge line center index I max (l) when the time-frequency power spectrum takes the maximum value on the first frame: ; Wherein, the The value of k when P (l, k) takes the maximum value in the range of 1-k-N fft on the first frame; S403, searching in the directions N fft and 1 respectively with the ridge line center index I max (l) as a starting point in the current time frame until the power value is smaller than eta P max (l) for the first time or reaches the frequency index boundary, obtaining a signal frequency upper boundary index I H,η (l) and a signal frequency lower boundary index I L,η (l), and meeting the following conditions: 1≤I L,η (l)≤I H,η (l)≤N fft ; s404, calculating a frequency bandwidth sequence B (l): ; Wherein max (·) represents the maximum value taken between them; S405, calculating a representative width parameter B ref : 。
- 6. The method according to claim 5, wherein S5, based on the time-frequency power spectrum P (l, k) and the representative width parameter B ref , calculates an initial time-frequency mask M 0 (l, k) by using SO-CFAR detection, specifically as follows: S501, initializing parameters required by two-dimensional SO-CFAR detection, wherein the method specifically comprises the following steps: initializing a time dimension protection unit half-width G t and a time dimension training unit half-width T t : ; ; Initializing a frequency dimension protection unit half-width G f and a frequency dimension training unit half-width T f : ; ; Wherein, the Is an upward rounding function; Initializing SO-CFAR false alarm probability P fa : ; Calculating the number N tr,t of time dimension training units and the number N tr,f of frequency dimension training units: ; ; Calculating a time-dimensional threshold amplification coefficient alpha t and a frequency-dimensional threshold amplification coefficient alpha f : ; ; S502, constructing a time dimension cross window training area in a time dimension by using a half-width G t of a time dimension protection unit and a half-width T t of a time dimension training unit, and calculating a time dimension noise power average value : ; Wherein S t (l, k) is a time-dimensional reference window power sum, and S t (l, k) is specifically expressed as: ; wherein i is an index variable summed by the time dimension training unit; S503, constructing a frequency dimension cross window training area in the frequency dimension by using the half width G f of the frequency dimension protection unit and the half width T f of the frequency dimension training unit, and calculating the frequency dimension noise power average value : ; Wherein S f (l, k) is the frequency-dimensional reference window power sum, and S f (l, k) is specifically expressed as: ; wherein j is an index variable summed by the frequency dimension training unit; S504, calculating a time dimension detection threshold Thr t (l, k) and a frequency dimension detection threshold Thr f (l, k) according to the time dimension noise power average value and the frequency dimension noise power average value: ; ; s505, adopting a selectivity minimum criterion, and synthesizing a time-dimensional noise power average value and a frequency-dimensional detection threshold into a two-dimensional SO-CFAR detection threshold Thr SO (l, k): ; wherein, min (·) represents taking the minimum value of the two; S506, calculating an initial time-frequency mask M 0 (l, k): 。
- 7. the method of claim 6, wherein S6 constructs a ridge neighborhood time-frequency mask M r (l, k) based on P (l, k) and calculates a final time-frequency mask in combination with the initial time-frequency mask M 0 (l, k) Based on the final time-frequency mask M (l, k) and complex time spectrum Obtaining denoised complex time spectrum The method is characterized by comprising the following steps: S601, constructing a ridge line neighborhood time-frequency mask M r (l, k) on each time frame l: ; S602, multiplying the ridge line neighborhood time-frequency mask M r (l, k) and the initial time-frequency mask M 0 (l, k) point by point to obtain a final time-frequency mask M (l, k): ; S603, calculating the denoising complex time spectrum by using the final time-frequency mask M (l, k) : 。
- 8. The method of claim 7, wherein S7 is based on denoising complex time spectrum And an optimal window length sequence W opt (l), calculating a reconstructed time-domain denoised complex baseband signal z d (n d ), specifically as follows: S701, calculating a reconstructed signal length: ; S702, sequences U num (n d ) and U den (n d ) of construction length N rec ): U num (n d )=0,U den (n d )=0, ; Wherein n d is the time index of the reconstructed signal, U num (n d ) is the superimposed molecular accumulation sequence, U den (n d ) is the superimposed denominator accumulation sequence; s703, for each time frame index l=1, 2, L performs the following processing: for denoising complex time spectrum on the first frame Performing N fft -point inverse discrete Fourier transform, and intercepting the previous W opt (l) sampling points to obtain a time domain frame signal M=1, 2,..w opt (l), where m is the intra discrete-time sample point index; Constructing a hanning window of length W opt (l) over the first frame as a function of the synthesis window ; Time domain frame signal of the first frame And synthesizing window functions Point-by-point multiplication for calculating windowed frame signal : ; Based on windowed frame signals And synthesizing window functions Accumulated computations U num (n d ) and U den (n d ): ; s704, calculating a reconstructed time-domain denoised complex baseband signal z d (n d for each sampling point): 。
- 9. The method of claim 8, wherein S8, the reconstructed time-domain denoised complex baseband signal z d (n d ) extracts the instantaneous phase and performs phase unwrapping, trend term removal and scaling factor conversion to extract the chest displacement waveform The method is characterized by comprising the following steps: s801, calculating an instantaneous phase sequence from sampling points according to the reconstructed time domain denoising complex baseband signal z d (n d : ; Wherein, the To extract Is a radial angle of (2); S802, pair Performing phase unwrapping and trend removal term processing to obtain a purified phase sequence , ; S803 to purge the phase sequence Mapping to thoracic cavity displacement, and calculating to obtain thoracic cavity displacement waveform : 。
Description
Thoracic cavity displacement waveform extraction method based on self-adaptive time-frequency analysis and two-dimensional filtering Technical Field The application belongs to the technical field of radar signal processing, and particularly relates to a thoracic cavity displacement waveform extraction method based on self-adaptive time-frequency analysis and two-dimensional filtering. Background The non-contact vital sign monitoring utilizes radar echo phase modulation to extract the chest cavity displacement of the human body, has the advantages of comfort, portability and the like, and is of great interest in the field of medical health. Current methods typically directly perform phase extraction, unwrapping, and trend term removal on the complex baseband signal of the selected range gate in the time domain to obtain the respiratory displacement waveform. However, the current method has the problem that the time domain phase extraction is easy to be interfered by noise, and the resolution of the fixed time frequency analysis method is limited. Disclosure of Invention The embodiment of the application provides a thoracic cavity displacement waveform extraction method based on self-adaptive time-frequency analysis and two-dimensional filtering, which can solve the problems that the time domain phase extraction is easy to be interfered by noise and the resolution of a fixed time-frequency analysis method is limited in the current method. In a first aspect, the embodiment of the application provides a thoracic cavity displacement waveform extraction method based on self-adaptive time-frequency analysis and two-dimensional filtering, which comprises the steps of S1, acquiring a radar complex baseband sampling sequence Z (N), n=1, 2, N, on a target distance gate to be processed, wherein the radar complex baseband sampling sequence Z (N) contains measured human body vital sign information, S2, initializing a plurality of preset parameters of the self-adaptive time-frequency analysis, S3, calculating a complex time spectrum Z (l, k) of the self-adaptive time-frequency analysis based on the complex baseband sampling sequence Z (N) Time-frequency power spectrum P (L, k) and optimal window length sequence W opt (L), l=1, 2,.. L, S4, estimating representative width parameter B ref of the target to be processed based on self-adaptive time-frequency analysis time-frequency power spectrum P (L, k), S5, Based on the time-frequency power spectrum P (l, k) and the representative width parameter B ref, calculating an initial time-frequency mask M 0 (l, k) by using minimum-selection constant false alarm SO-CFAR detection, S6, Constructing a ridge neighborhood time-frequency mask M r (l, k) based on P (l, k), calculating a final time-frequency mask M (l, k) by combining the initial time-frequency mask M 0 (l, k), and obtaining a denoised complex time-frequency spectrum according to the final time-frequency mask M (l, k) and the complex time-frequency spectrum Z (l, k)S7, denoising complex time spectrumAnd the optimal window length sequence W opt (l), calculating a reconstructed time domain denoising complex baseband signal z d (nd), S8, extracting instantaneous phase and carrying out phase unwrapping, trend term removal and scaling coefficient conversion on the reconstructed time domain denoising complex baseband signal z d (nd), and extracting chest cavity displacement waveforms。 In a possible implementation manner of the first aspect, the method for acquiring the Lei Dafu baseband sampling sequence z (n) in S1 specifically includes the following steps: Acquiring a Lei Dafu baseband sampling sequence z (N) containing vital sign information on a range gate of a target to be processed by using a frequency modulation continuous wave FMCW radar data acquisition device, wherein n=1, 2, & gt, N is a time index of the Lei Dafu baseband sampling sequence z (N), N is the total number of sampling points of a complex baseband sequence, and acquiring a slow time sampling frequency f s and a radar working wavelength lambda corresponding to the Lei Dafu baseband sampling sequence z (N), wherein the slow time sampling frequency and the radar working wavelength respectively satisfy that f s ≤50 Hz,3 mm≤ λ is less than or equal to 10 Hz and less than or equal to 13 mm, and N is sampling points of not less than 30 seconds, namely 。 Optionally, in another possible implementation manner of the first aspect, S2, initializing a plurality of preset parameters of the adaptive time-frequency analysis specifically includes the following steps: S201, initializing a minimum candidate window length W min and a maximum candidate window length W max of self-adaptive time-frequency analysis: , w min and W max are positive integers, and the units are sampling points; s202, initializing window length search steps delta W, wherein delta W=1, and the unit is a sampling point; S203, initializing a window length candidate set omega of self-ada