Millimeter wave radar vital sign detection method based on HHO-CEEMDAM algorithm

Abstract

The invention discloses a millimeter wave radar vital sign detection method based on an HHO-CEEMDAM algorithm, and the method comprises the steps: building a millimeter wave radar experiment system, and collecting an intermediate frequency signal of a human body echo; performing data reading and recombination, extracting phase features, and enhancing a target signal through non-coherent accumulation to determine the distance between the chest of the human body and the radar; recovering the phase of the vital sign signal from the incoherent accumulation FFT result by using the linear characteristic of arc tangent demodulation, and performing phase unwrapping and phase difference to obtain optimized phase information; setting a CEEMDAN parameter initialization range and designing a fitness function; a CEEMDAN parameter is optimized by using an HHO algorithm; performing CEEMDAN decomposition by using the optimized parameters, and screening breathing and heartbeat IMF components; and according to the IMF component, obtaining estimated values of the respiratory rate and the heart rate. The method improves the estimation precision of the respiratory rate and the heart rate of the millimeter wave radar in fatigue driving detection.

Description

Technical Field

[0001] This invention belongs to the interdisciplinary field of biomedical signal processing and radar technology, and relates to a method for detecting human vital signs based on frequency modulated continuous wave (FMCW) radar, particularly a non-contact, high-precision real-time detection method for heart rate and respiratory rate. Background Technology

[0002] In vital sign detection methods, compared to traditional methods, millimeter-wave radar-based methods are unaffected by environmental conditions such as light and temperature, do not rely on visual information acquisition, and avoid privacy issues. Furthermore, its strong penetrating power allows it to penetrate obstacles such as clothing and seats to directly detect physiological signals such as the driver's heart rate and respiratory rate, providing a more accurate assessment of fatigue status. These characteristics give millimeter-wave radar a unique advantage in fatigue driving detection, providing a reliable and efficient solution for intelligent driving assistance systems.

[0003] However, respiratory and heartbeat signals are typically weak and easily affected by environmental noise and subtle bodily movements. Furthermore, current methods for estimating respiratory rate and heart rate heavily rely on signal decomposition algorithms. Existing methods, due to fixed parameters, suffer from insufficient mode aliasing and noise suppression, leading to increased errors in heart rate and respiratory rate estimation, especially under vehicle vibration or sudden bodily movement disturbances. Therefore, dynamically optimizing the decomposition algorithm parameters to improve the accuracy of vital sign estimation and thus enhance the overall accuracy of fatigue driving detection is crucial to overcoming the technological bottleneck of millimeter-wave radar fatigue detection. Summary of the Invention

[0004] To optimize the parameters of the signal decomposition algorithm and obtain accurate respiratory rate and heart rate estimation results, this invention provides a millimeter-wave radar vital sign detection method based on the HHO-CEEMDAM algorithm. This method uses frequency-modulated continuous wave radar for non-contact and non-invasive vital sign detection. Specifically, the acquired human echo intermediate frequency signal is read and reconstructed. The distance between the human chest cavity and the radar is determined through range-dimensional FFT and sliding cancellation. The phase of the vital sign signal is recovered and optimized using arctangent demodulation. Based on the HHO algorithm, the input parameters of CEEMDAM are optimized with minimizing the signal envelope entropy as the objective function. The optimized CEEMDAM algorithm is then used to decompose the preprocessed signal, filter effective components, and estimate respiratory rate and heart rate. Experimental results show that this method maintains the accuracy of heart rate and respiratory rate detection while exhibiting good consistency in estimation results under different environmental conditions.

[0005] The technical solution adopted by this invention to solve its technical problem is:

[0006] A millimeter-wave radar vital sign detection method based on the HHO-CEEMDAM algorithm includes the following steps:

[0007] Step 1: Set up a millimeter-wave radar experimental system, with the human body positioned in front of the radar, and collect the intermediate frequency signal of the human body echo;

[0008] Step 2: Read and reassemble the received intermediate frequency signal, obtain frequency domain information through distance-dimensional FFT and suppress noise through sliding cancellation technology, extract phase features, and finally enhance the target signal through non-coherent accumulation to determine the distance between the human chest cavity and the radar.

[0009] Step 3: Utilize the linearity of arctangent demodulation to recover the phase of the vital signs signal from the incoherently accumulated FFT results, and perform phase unwrapping and phase difference to obtain optimized phase information;

[0010] Step 4: Initialize the range of CEEMDAN parameters (noise amplitude, number of noise additions, and maximum number of iterations) and design the fitness function;

[0011] Step 5: Optimize CEEMDAN parameters using the HHO algorithm;

[0012] Step 6: Perform CEEMDAN decomposition using optimized parameters and filter the respiratory and cardiac IMF components;

[0013] Step 7: Obtain estimates of respiratory rate and heart rate based on the IMF components.

[0014] Furthermore, in step 1, the millimeter-wave radar experimental system uses an FMCW radar with an operating range from 76 GHz to 81 GHz and N transmitting antennas. TX =2 and the number of receiving antennas N RX =4, equivalent bandwidth B eq for:

[0015]

[0016] Where S is the linear frequency modulation slope, N is the number of effective ADC samples in linear frequency modulation, and f s T is the ADC sampling rate. c The distance resolution is calculated based on the duration of a single chirp:

[0017]

[0018] Where c is the speed of light, and V is the maximum speed. max Calculated as

[0019]

[0020] Where λ is the wavelength, the operating frequency of the FMCW radar changes within the linear frequency modulation period, and the transmitted signal is represented as:

[0021] s t (t)=A t cos(2πf0t+πSt 2 (5)

[0022] Among them, A t f0 is the amplitude of the transmitted wave, and f0 is the starting frequency of the chirp.

[0023] The human body echo signal is delayed by τ, and the received signal is:

[0024] s r (t)=A r cos(2πf0(t-τ+πS(t-τ) 2 (6)

[0025] After mixing, the difference frequency signal frequency f b The distance R between the target and the target is:

[0026] f b =Sτ (7)

[0027]

[0028] The subject collected data at any location within the platform for 51.2 seconds. The DCA1000 high-speed data acquisition card was used to capture the intermediate frequency signal of the IWR1843 in real time, and the received intermediate frequency radar human respiratory and heart rate data (I / Q data) was transmitted to the computer.

[0029] Furthermore, step 2 is as follows:

[0030] (2.1) The binary data is parsed into a complex signal matrix by alternating I / Q operations, and the data is reassembled to obtain a dimension of (N). RX N ADC ×N chirp ×N TX ×N Frame A matrix of N; where N chirp N ADC and N Frame These represent the number of chirps, the number of ADC sampling points, and the number of frames in the original echo data, respectively.

[0031] (2.2) Perform a distance-dimensional FFT on the reconstructed result to obtain a complex matrix, where each row corresponds to the frequency domain information of a distance cell;

[0032] (2.3) Arrange the FFT results of each chirp row by row to obtain complete frequency domain data. Introduce a sliding cancellation technique: for each pair of consecutive chirps, differencing their FFT results removes background noise and clutter, making the target signal clearer. The differencing FFT results are used for subsequent target detection and analysis.

[0033] (2.4) Using the real and imaginary parts of the Fourier transform result, the phase value at each frequency point is calculated using the arctangent function. The calculation formula is as follows:

[0034] θ=arctan2(Imag(S),Real(S)) (9)

[0035] Where S is the FFT result after difference, and Real(S) and Imag(S) are the real and imaginary parts of S, respectively.

[0036] (2.5) By analyzing the amplitude and phase of each distance point, the distance and motion state information of the target are obtained. Non-coherent accumulation technology is used to accumulate the amplitude of the signal to enhance the strength of the target signal and reduce noise.

[0037] Furthermore, step 3 is as follows:

[0038] (3.1) Utilizing the linearity of arctangent demodulation, the phase of the point with maximum energy is extracted from the accumulated FFT results by traversal.

[0039] (3.2) The phase unwrapping algorithm is used to detect and track changes in the signal phase, adjusting the phase value to a continuous range to eliminate discontinuities caused by periodic signal changes. The processing flow is as follows: calculate the phase difference Δθ between two adjacent phase points. The phase unwrapping expression is as follows:

[0040] θ i =θ i-1 +Δθ+2kπ (10)

[0041] Where, when |Δθ|<π, k=0; when Δθ>π, k=-1; when Δθ<-π, k=1;

[0042] (3.3) Using phase difference, the phase difference Δθ = θ between two adjacent points is calculated. i -θ i-1 This allows us to obtain the phase change of the target, which can then be used for subsequent target trajectory analysis and motion prediction.

[0043] The process of step 4) is as follows:

[0044] (4.1) Before using HHO to optimize the key parameters of CEEMDAN, the range of values ​​for each parameter must be given. For example, the range of noise amplitude ε is 0.1 to 1.0, the range of noise addition N is 50 to 200, and the range of maximum iteration number MaxIter is 100 to 500.

[0045] (4.2) Using the minimum envelope entropy as the fitness function, for a given signal component x(t), its analytic signal is first obtained through Hilbert transform:

[0046] z(t)=x(t)+jH[x(t)] (11)

[0047] To eliminate the influence of amplitude, the envelope signal e(t) = |z(t)| needs to be normalized to obtain the probability distribution sequence:

[0048]

[0049] Choose the minimum envelope entropy value E e As the fitness function value F:

[0050]

[0051] The process of step 5) is as follows:

[0052] (5.1) Within the range of parameter values, the HHO algorithm randomly generates a set of CEEMDAN parameters (ε, N, MaxIter) as the initial population;

[0053] (5.2) Calculate the current energy factor

[0054]

[0055] Where E0 is the initial energy value, t is the current iteration number, and T is the maximum iteration number. The working mode is selected based on the |E| value.

[0056] (5.3) When |E|≥1, the algorithm enters the exploration phase, updates the position, and randomly searches for potential optimal regions in the parameter space;

[0057]

[0058] Where X(t) represents the current position of the individual, X rand (t) represents the randomly selected individual location, X rabbit (t) represents the location of the prey (the current optimal solution), X mean (t) represents the average position of the population, r1, r2, r3, r4 and q are random numbers, and LB and UB are the lower and upper bounds of the search space, respectively.

[0059] (5.4) When |E|<1, the algorithm enters the development phase, adopting different hunting strategies based on the escape probability r and the energy state. r is a random number, and the strategies include the following four:

[0060] When |E|≥0.5 and r≥0.5, a soft-encirclement strategy is adopted, and its position update formula is:

[0061] X(t+1)=ΔX(t)-E|JX rabbit (t)-X(t)| (16)

[0062] When |E| < 0.5 and r ≥ 0.5, a hard-bound strategy is used to update its position:

[0063] X(t+1)=X rabbit (t)-E|Δx(t)| (17)

[0064] When |E|≥0.5 and r<0.5, a gradual rapid dive strategy is used to update its position:

[0065] X(t+1)=X rabbit (t)-E|JX rabbit (t)-X(t)| (18)

[0066] When |E| < 0.5 and r < 0.5, a random raid strategy is used to update its position:

[0067] X(t+1)=X rabbit (t)-E|JX rabbit (t)-X(t)|+Levy(D) (19)

[0068] Where, ΔX(t)=X rabbit (t)-X(t) represents the position difference between the prey and the current individual, J is the random jump intensity, Levy(D) represents the random step size of Levy's flight, and D is the dimension of the problem;

[0069] (5.5) Update the position of individuals in the population according to the hunting strategy and calculate the fitness value. Finally, stop the search and output the optimal parameter combination when the maximum number of iterations is reached or the fitness value converges.

[0070] The process of step 6) is as follows:

[0071] (6.1) Set m = 1 and set the initial residual r0(t) to the original signal x(t). Based on the optimized noise amplitude ε and the number of noise additions N, add Gaussian white noise ε0·ω to the original signal. i (t), where ω i (t) represents the noise added for the i-th time, i = 1, 2, ..., N;

[0072] (6.2) Perform EMD decomposition on the noisy signal to obtain the first IMF component. Calculate the mean of all noisy signal decomposition results to obtain the final value of the first IMF component:

[0073]

[0074] The updated residual r1(t) is the original signal minus the first IMF component:

[0075] r1(t)=x(t)-c1(t) (21)

[0076] (6.3) For the m-th IMF component (m≥2), for the current residual r m-1 (t) Add noise ε m-1 ·E m-1 (ω i (t)), where E m-1 (·) represents the (m-1)th IMF component after EMD decomposition of the noise. EMD decomposition of the noisy signal yields the m-th IMF component. By calculating the mean of all noisy signal decomposition results, the final value of the m-th IMF component is obtained:

[0077]

[0078] Then, update the residual r. m (t) represents the current residual minus the m-th IMF component:

[0079] r m (t)=r m-1 (t)-c m (t) (23)

[0080] (6.4) Residual r m The algorithm terminates when (t) becomes a monotonic function or can no longer be decomposed into IMF components. Based on the decomposed IMF components, the IMF component corresponding to 0.1–0.6 Hz for respiratory signals (such as IMF9) is selected; and the IMF component corresponding to 0.8–2 Hz for heartbeat signals (such as IMF8) is selected.

[0081] The process of step 7) is as follows:

[0082] (7.1) Perform spectral analysis on the selected IMF components and extract the dominant frequency components;

[0083] (7.2) Obtain the respiratory rate R b and heart rate H b The estimated values ​​are all in BPM (times per minute).

[0084] R b =f r ×60 (24)

[0085] H b =f h ×60 (25)

[0086] Among them, f r f is the dominant frequency of the detected respiratory signal. h The main frequency of the detected heartbeat signal is expressed in Hz.

[0087] The technical concept of this invention is as follows: Data reading and reconstruction, range-dimensional FFT, and sliding cancellation are performed on the acquired human echo intermediate frequency signal to determine the distance between the human chest cavity and the radar. The phase of the vital signs signal is recovered and optimized using arctangent demodulation. The input parameters of CEEMDAN are optimized based on the HHO algorithm with the objective function of minimizing the signal envelope entropy. The optimized CEEMDAN algorithm is then used to decompose the preprocessed signal, filter effective components, and estimate respiratory rate and heart rate. This millimeter-wave radar-based method for detecting human vital signs is an innovative monitoring technology with broad application prospects and market demand. With continuous technological development and improvement, this method will play an important role in multiple fields such as medicine, home security, and safety, improving people's quality of life and safety levels.

[0088] The beneficial effects of this invention are mainly reflected in:

[0089] 1. By leveraging the high resolution and high precision of FMCW radar, this invention can accurately detect minute vital signs in the human body, such as breathing and heartbeat.

[0090] 2. This invention utilizes the long-range detection capability of FMCW radar to achieve non-contact monitoring of human vital signs. This avoids the inconvenience and potential infection risks associated with sensor attachment in traditional methods, thus improving user comfort and safety.

[0091] 3. By introducing the sliding cancellation technique, the FFT results of the differential adjacent chirps are used to effectively suppress static clutter (such as environmental reflections) and low-frequency noise (such as equipment vibration).

[0092] 4. The present invention performs phase unwinding and phase difference on the extracted phase, thereby improving the effectiveness and recognition accuracy of obtaining phase feature information.

[0093] 5. CEEMDAN optimized using the HHO algorithm achieves adaptive parameter optimization, significantly improving the separation accuracy of respiratory and heartbeat signals, effectively suppressing spectral aliasing interference, and exhibiting good consistency in different environments. Attached Figure Description

[0094] Figure 1 This is a flowchart illustrating the implementation process of the vital signs detection method of the present invention;

[0095] Figure 2 This is a schematic diagram of the distance dimension FFT result of a single chirp according to the present invention;

[0096] Figure 3 This is a schematic diagram of the original signal after phase extraction according to the present invention;

[0097] Figure 4 This is a schematic diagram of the phase difference result of the present invention;

[0098] Figure 5 This is a flowchart of the HHO algorithm of the present invention;

[0099] Figure 6 This is a time-domain diagram of the IMF5-IMF8 components of the present invention;

[0100] Figure 7 This is a spectrum diagram of the MF5-IMF8 components of the present invention;

[0101] Figure 8 This is a schematic diagram of the respiratory time-domain waveform of the present invention;

[0102] Figure 9 This is a schematic diagram of the heartbeat time-domain waveform of the present invention;

[0103] Figure 10 This is a schematic diagram of the average error rate of breathing frequency in a vehicle (left) and indoor environment (right) according to the present invention;

[0104] Figure 11 This is a schematic diagram of the average error rate of heart rate in the vehicle (left) and indoor environment (right) according to the present invention; Specific implementation methods

[0105] The present invention will now be further described with reference to the accompanying drawings.

[0106] Reference Figures 1 to 11 A method for detecting vital signs using millimeter-wave radar based on the HHO-CEEMDAM algorithm includes the following steps:

[0107] Step 1: Set up a millimeter-wave radar experimental system, with the human body positioned in front of the radar, and collect the intermediate frequency signal of the human body echo;

[0108] The millimeter-wave radar experimental system platform described above is deployed using the IWR1843 millimeter-wave frequency modulated continuous wave (FMCW) radar sensor developed by Texas Instruments, operating in the 76 GHz to 81 GHz frequency band. This example uses N transmitting antennas. TX =2 and the number of receiving antennas N RX =4, equivalent bandwidth B eqfor:

[0109]

[0110] Where S is the linear frequency modulation slope, N is the number of effective ADC samples in linear frequency modulation, and f s T is the ADC sampling rate. c The distance resolution is calculated based on the duration of a single chirp:

[0111]

[0112] Where c is the speed of light, and V is the maximum speed. max Calculated as

[0113]

[0114] Where λ is the wavelength, the operating frequency of the FMCW radar changes within the linear frequency modulation period, and the transmitted signal is represented as:

[0115] s t (t)=A t cos(2πf0t+πSt 2 (5)

[0116] Among them, A t f0 is the amplitude of the transmitted wave, and f0 is the starting frequency of the chirp.

[0117] The human body echo signal is delayed by τ, and the received signal is:

[0118] s r (t)=A r cos(2πf0(t-τ+πS(t-τ) 2 (6)

[0119] After mixing, the difference frequency signal frequency f b The distance R between the target and the target is:

[0120] f b =Sτ (7)

[0121]

[0122] The subject collected data at any location within the platform for 51.2 seconds. The DCA1000 high-speed data acquisition card was used to capture the intermediate frequency signal of the IWR1843 in real time, and the received intermediate frequency radar human respiratory and heart rate data (I / Q data) was transmitted to the computer.

[0123] Step 2: Read and reassemble the received intermediate frequency signal, obtain frequency domain information through distance-dimensional FFT and suppress noise through sliding cancellation technology, extract phase features, and finally enhance the target signal through non-coherent accumulation to determine the distance between the human chest cavity and the radar.

[0124] The process of reading and reconstructing the received intermediate frequency signal, obtaining frequency domain information through distance-dimensional FFT and suppressing noise using sliding cancellation technology, extracting phase features, and finally enhancing the target signal through non-coherent accumulation to determine the distance between the human chest cavity and the radar is as follows:

[0125] (2.1) The binary data is parsed into a complex signal matrix by alternating I / Q operations, and the data is reassembled to obtain a dimension of (N). RX N ADC ×N chirp ×N TX ×N Frame A matrix of N; where N chirp N ADC and N Frame These represent the number of chirps, the number of ADC sampling points, and the number of frames in the original echo data, respectively.

[0126] (2.2) Perform a distance-dimensional FFT on the reconstructed result to obtain a complex matrix, where each row corresponds to the frequency domain information of a distance cell;

[0127] (2.3) Arrange the FFT results of each chirp in rows to obtain complete frequency domain data. Introduce the sliding cancellation technique. For each pair of consecutive chirps, the FFT results are differentially divided to remove the influence of background noise and clutter, making the target signal clearer. The differentially divided FFT results are used for subsequent target detection and analysis.

[0128] (2.4) Using the real and imaginary parts of the Fourier transform result, the phase value at each frequency point is calculated using the arctangent function. The calculation formula is as follows:

[0129] θ=arctan 2(Imag(S),Real(S)) (9)

[0130] Where S is the FFT result after difference, and Real(S) and Imag(S) are the real and imaginary parts of S, respectively.

[0131] (2.5) By analyzing the amplitude and phase of each distance point, the distance and motion state information of the target are obtained. Non-coherent accumulation technology is used to accumulate the amplitude of the signal to enhance the strength of the target signal and reduce noise.

[0132] Step 3: Utilize the linearity of arctangent demodulation to recover the phase of the vital signs signal from the incoherently accumulated FFT results, and perform phase unwrapping and phase difference to obtain optimized phase information;

[0133] The process of recovering the phase of vital signs signals from incoherently accumulated FFT results using the linear characteristics of arctangent demodulation, and then performing phase unwrapping and phase difference to obtain optimized phase information is as follows:

[0134] (3.1) Utilizing the linearity of arctangent demodulation, the phase of the point with maximum energy is extracted from the accumulated FFT results by traversal.

[0135] (3.2) The phase unwrapping algorithm is used to detect and track changes in the signal phase, adjusting the phase value to a continuous range to eliminate discontinuities caused by periodic signal variations. The processing flow is as follows: calculate the phase difference Δθ between two adjacent phase points, with the phase unwrapping expression as follows:

[0136] θ i =θ i-1 +Δθ+2kπ (10)

[0137] Where, when |Δθ|<π, k=0; when Δθ>π, k=-1; when Δθ<-π, k=1;

[0138] (3.3) Using phase difference, the phase difference Δθ = θ between two adjacent points is calculated. i -θ i-1 This allows us to obtain the phase change of the target, which can then be used for subsequent target trajectory analysis and motion prediction.

[0139] Step 4: Initialize the range of CEEMDAN parameters (noise amplitude, number of noise additions, and maximum number of iterations) and design the fitness function;

[0140] The process of initializing the CEEMDAN parameters (noise amplitude, number of noise additions, and maximum number of iterations) and designing the fitness function is as follows:

[0141] (4.1) Before using HHO to optimize the key parameters of CEEMDAN, the range of values ​​for each parameter must be given. For example, the range of noise amplitude ε is 0.1 to 1.0, the range of noise addition N is 50 to 200, and the range of maximum iteration number MaxIter is 100 to 500.

[0142] (4.2) Using the minimum envelope entropy as the fitness function, for a given signal component x(t), its analytic signal is first obtained through Hilbert transform:

[0143] z(t)=x(t)+jH[x(t)] (11)

[0144] To eliminate the influence of amplitude, the envelope signal e(t) = |z(t)| needs to be normalized to obtain the probability distribution sequence:

[0145]

[0146] Choose the minimum envelope entropy value E e As the fitness function value F:

[0147]

[0148] Step 5: Optimize CEEMDAN parameters using the HHO algorithm;

[0149] The process of optimizing CEEMDAN parameters using the HHO algorithm is as follows:

[0150] (5.1) Within the range of parameter values, the HHO algorithm randomly generates a set of CEEMDAN parameters (ε, N, MaxIter) as the initial population;

[0151] (5.2) Calculate the current energy factor

[0152]

[0153] Where E0 is the initial energy value, t is the current iteration number, and T is the maximum iteration number. The working mode is selected based on the |E| value;

[0154] (5.3) When |E|≥1, the algorithm enters the exploration phase, updates the position, and randomly searches for potential optimal regions in the parameter space;

[0155]

[0156] Where X(t) represents the current position of the individual, X rand (t) represents the randomly selected individual location, X rabbit (t) represents the location of the prey (the current optimal solution), X mean (t) represents the average position of the population, r1, r2, r3, r4 and q are random numbers, and LB and UB are the lower and upper bounds of the search space, respectively.

[0157] (5.4) When |E| < 1, the algorithm enters the development phase. Different hunting strategies are adopted based on the escape probability r and the energy state, where r is a random number. The strategies include the following four:

[0158] When |E|≥0.5 and r≥0.5, a soft-encirclement strategy is adopted, and its position update formula is:

[0159] X(t+1)=ΔX(t)-E|JX rabbit (t)-X(t)| (16)

[0160] When |E| < 0.5 and r ≥ 0.5, a hard-bound strategy is used to update its position:

[0161] X(t+1)=X rabbit (t)-E|Δx(t)| (17)

[0162] When |E|≥0.5 and r<0.5, a gradual rapid dive strategy is used to update its position:

[0163] X(t+1)=X rabbit (t)-E|JX rabbit (t)-X(t)| (18)

[0164] When |E| < 0.5 and r < 0.5, a random raid strategy is used to update its position:

[0165] X(t+1)=X rabbit (t)-E|JX rabbit (t)-X(t)|+Levy(D) (19)

[0166] Where, ΔX(t)=X rabbit (t)-X(t) represents the position difference between the prey and the current individual, J is the random jump intensity, Levy(D) represents the random step size of Levy's flight, and D is the dimension of the problem;

[0167] (5.5) Update the positions of individuals in the population according to the hunting strategy and calculate the fitness value. Finally, stop the search and output the optimal parameter combination when the maximum number of iterations is reached or the fitness value converges;

[0168] Step 6: Perform CEEMDAN decomposition using optimized parameters and filter the respiratory and cardiac IMF components;

[0169] The process of using optimized parameters to perform CEEMDAN decomposition and screen for respiratory and cardiac IMF components is as follows:

[0170] (6.1) Set m = 1 and set the initial residual r0(t) to the original signal x(t). Based on the optimized noise amplitude ε and the number of noise additions N, add Gaussian white noise ε0·ω to the original signal. i (t), where ω i (t) represents the noise added for the i-th time, i = 1, 2, ..., N;

[0171] (6.2) Perform EMD decomposition on the noisy signal to obtain the first IMF component. Calculate the mean of all noisy signal decomposition results to obtain the final value of the first IMF component:

[0172]

[0173] The updated residual r1(t) is the original signal minus the first IMF component:

[0174] r1(t)=x(t)-c1(t) (21)

[0175] (6.3) For the m-th IMF component (m≥2), for the current residual r m-1 (t) Add noise ε m-1 ·E m-1 (ω i (t)), where E m-1 (·) represents the (m-1)th IMF component after EMD decomposition of the noise. EMD decomposition of the noisy signal yields the m-th IMF component. The final value of the m-th IMF component is obtained by calculating the mean of all noisy signal decomposition results:

[0176]

[0177] Then, update the residual r. m (t) represents the current residual minus the m-th IMF component:

[0178] r m (t)=r m-1 (t)-c m (t) (23)

[0179] (6.4) Residual r m The algorithm terminates when (t) becomes a monotonic function or can no longer be decomposed into IMF components. Based on the decomposed IMF components, the IMF component corresponding to 0.1–0.6 Hz for respiratory signals (such as IMF9) is selected; and the IMF component corresponding to 0.8–2 Hz for heartbeat signals (such as IMF8) is selected.

[0180] Step 7: Obtain estimated values ​​for respiratory rate and heart rate based on the IMF components.

[0181] The process of obtaining estimated values ​​of respiratory rate and heart rate based on IMF components is as follows:

[0182] (7.1) Perform spectral analysis on the selected IMF components and extract the dominant frequency components;

[0183] (7.2) Obtain the respiratory rate R b and heart rate H b The estimated values ​​are all in BPM (bpm), times per minute.

[0184] R b =f r ×60 (24)

[0185] H b =f h ×60 (25)

[0186] Among them, f r f is the dominant frequency of the detected respiratory signal. h The main frequency of the detected heartbeat signal is expressed in Hz.

[0187] The experimental system deployed uses an IWR1843 radar sensor and a DCA 1000 data acquisition card. The IWR1843 radar sensor is equipped with two transmitting antennas and four receiving antennas. mmWave Studio, used in conjunction with a Windows graphical interface program, is responsible for configuring and controlling the TI radar sensor and receiving the data stream from the analog-to-digital converter. When the radar sensor is activated, its emitted signal is reflected by the human chest. The transmitted and received signals are mixed and filtered to form an intermediate frequency (IF) signal, which is then converted into an ADC data stream through an analog-to-digital converter. The ADC data stream is transmitted to the DCA 1000 via wires. The DCA 1000 is connected to a computer via an Ethernet interface, and the collected ADC data stream is transmitted to the computer for further analysis via the Ethernet interface. The collected ADC data stream is processed using the algorithm presented in this chapter to separate respiratory and heart rate data.

[0188] In this embodiment, the person stands still in front of the radar device, which has a distance resolution of 3.76 centimeters.

[0189] The echo signal reflected from the radar's transmitted signal after encountering the target is read and reconstructed. Range-dimensional FFT is applied to extract the target range gate information from the echo signal. A schematic diagram of the range dimension of a chirp is shown below. Figure 2 As shown.

[0190] A sliding cancellation technique is used to differencing the FFT results of each pair of consecutive chirps. Precise phase information at each frequency point is obtained through phase extraction. Simultaneously, non-coherent accumulation technology is combined to enhance the target signal strength. Phase extraction is performed on the points of maximum energy using an ergonomic approach. A schematic diagram of the original signal after phase extraction is shown below. Figure 3 As shown.

[0191] Phase unwrapping and phase differencing are performed to obtain optimized phase information. The signal after phase differencing is as follows: Figure 4 As shown.

[0192] The HHO algorithm initializes the range of CEEMDAN parameters (noise amplitude, number of noise additions, and maximum number of iterations) and uses the minimum envelope entropy as the fitness function. A set of CEEMDAN parameters is randomly generated as the initial population. A global search is performed in the parameter space to find the optimal parameter combination. Based on the escape probability and energy state, different hunting strategies (such as soft encirclement, hard encirclement, gradual rapid dive, and random raid) are selected to gradually approach the optimal solution. Then, the positions of individuals in the population are updated according to the hunting strategy, and the fitness value is calculated. Finally, when the maximum number of iterations is reached or the fitness value converges, the search stops and the optimal parameter combination is output. The flowchart of the HHO algorithm is as follows: Figure 5 As shown.

[0193] Based on the optimized noise amplitude and the number of noise additions, Gaussian white noise is added to the original signal, and the signal is decomposed. The results of multiple decompositions are averaged to obtain stable IMF components. The time-domain plots of the IMF5-IMF8 components are shown below. Figure 6 As shown. The spectrum of the IMF5-IMF8 components is as follows. Figure 7 As shown;

[0194] IMF components were selected based on the frequency ranges of respiration (0.1–0.6 Hz) and heart rate (0.8–2 Hz). Spectral analysis was performed on the selected IMF components to extract the dominant frequency components, and estimates of respiratory rate and heart rate were calculated. A schematic diagram of the respiratory time-domain waveform is shown below. Figure 8 As shown in the diagram. A schematic diagram of the heartbeat time-domain waveform is shown below. Figure 9 As shown.

[0195] To explore the practicality of this method, respiratory rate and heart rate were measured at a distance of 0.5 on 10 test subjects aged 22-25 years with normal physical indicators and normal daily routines on the test day. The measured values ​​of respiratory rate and heart rate of some test subjects (left) are shown in the comparison with the actual values ​​displayed by the electrocardiogram. The error rate of respiratory rate and heart rate of the test subjects was controlled between 0% and 6%.

[0196] To evaluate the impact of environmental changes on the accuracy of respiratory and heart rate estimation algorithms, this paper compares the estimation results of respiratory rate and heart rate in an in-vehicle environment and an open indoor environment. The comparison results of the average error rates of respiratory rate and heart rate in the two scenarios are as follows: Figure 10 and Figure 11 As shown:

[0197] according to Figure 10 Data shows that the average error rates for respiratory rate and heart rate of the 10 test subjects in the in-vehicle environment were 3.28% and 3.30%, respectively; according to Figure 11Data shows that in an open indoor environment, the average error rates for these two indicators are 2.97% and 2.96%, respectively. Although there is more noise and interference in the in-vehicle environment, the difference in error rates is small, reflecting that the estimation results under different environmental conditions have good consistency.

[0198] The embodiments described in this specification are merely examples of implementations of the inventive concept and are for illustrative purposes only. The scope of protection of this invention should not be considered limited to the specific forms described in these embodiments; rather, it extends to equivalent technical means conceived by those skilled in the art based on the inventive concept.

Claims

1. A millimeter-wave radar vital sign detection method based on the HHO-CEEMDAM algorithm, characterized in that, The method includes the following steps: Step 1: Set up a millimeter-wave radar experimental system, with the human body positioned in front of the radar, and collect the intermediate frequency signal of the human body echo; Step 2: Read and reassemble the received intermediate frequency signal, obtain frequency domain information through distance-dimensional FFT and suppress noise through sliding cancellation technology, extract phase features, and finally enhance the target signal through non-coherent accumulation to determine the distance between the human chest cavity and the radar. Step 3: Utilize the linearity of arctangent demodulation to recover the phase of the vital signs signal from the incoherently accumulated FFT results, and perform phase unwrapping and phase difference to obtain optimized phase information; Step 4: CEEMDAN parameter initialization range setting and fitness function design, wherein the CEEMDAN parameters include noise amplitude, number of noise additions, and maximum number of iterations; Step 5: Optimize CEEMDAN parameters using the HHO algorithm; Step 6: Perform CEEMDAN decomposition using optimized parameters and filter the respiratory and cardiac IMF components; Step 7: Obtain estimates of respiratory rate and heart rate based on the IMF components.

2. The millimeter-wave radar vital sign detection method based on the HHO-CEEMDAM algorithm as described in claim 1, characterized in that, In step 1, the millimeter-wave radar experimental system uses an FMCW radar with an operating range from 76 GHz to 81 GHz and N transmitting antennas. TX =2 and the number of receiving antennas N RX =4, equivalent bandwidth B eq for: Where S is the linear frequency modulation slope, N is the number of effective ADC samples in linear frequency modulation, and f s T is the ADC sampling rate. c The distance resolution is calculated based on the duration of a single chirp: Where c is the speed of light, and V is the maximum speed. max Calculated as Where λ is the wavelength, the operating frequency of the FMCW radar changes within the linear frequency modulation period, and the transmitted signal is represented as: s t (t)=A t cos(2πf0t+πSt 2 ) (5) Among them, A t f0 is the amplitude of the transmitted wave, and f0 is the starting frequency of the chirp. The human body echo signal is delayed by τ, and the received signal is: s r (t)=A r cos(2πf0(t-τ+πS(t-τ) 2 ) (6) After mixing, the difference frequency signal frequency f b The distance R between the target and the target is: f b =Sτ (7) The subject collected data at any location within the platform for 51.2 seconds. The DCA1000 high-speed data acquisition card was used to capture the intermediate frequency signal of the IWR1843 in real time, and the received intermediate frequency radar human respiratory and heart rate data, i.e., I / Q data, was transmitted to the computer.

3. The millimeter-wave radar vital sign detection method based on the HHO-CEEMDAM algorithm as described in claim 1 or 2, characterized in that, The process of step 2 is as follows: (2.1) The binary data is parsed into a complex signal matrix by alternating I / Q operations, and the data is reassembled to obtain a dimension of (N). RX N ADC ×N chirp ×N TX ×N Frame A matrix of N; where N chirp N ADC and N Frame These represent the number of chirps, the number of ADC sampling points, and the number of frames in the original echo data, respectively. (2.2) Perform a distance-dimensional FFT on the reconstructed result to obtain a complex matrix, where each row corresponds to the frequency domain information of a distance cell; (2.3) Arrange the FFT results of each chirp in rows to obtain complete frequency domain data. Introduce the sliding cancellation technique. For each pair of consecutive chirps, the FFT results are differentially divided to remove the influence of background noise and clutter, making the target signal clearer. The differential FFT results are used for subsequent target detection and analysis. (2.4) Using the real and imaginary parts of the Fourier transform result, the phase value at each frequency point is calculated using the arctangent function. The calculation formula is as follows: θ=arctan2(Imag(S),Real(S)) (9) Where S is the FFT result after difference, and Real(S) and Imag(S) are the real and imaginary parts of S, respectively; (2.5) By analyzing the amplitude and phase of each distance point, the distance and motion state information of the target can be obtained. Non-coherent accumulation technology is used to accumulate the amplitude of the signal to enhance the strength of the target signal and reduce noise.

4. The millimeter-wave radar vital sign detection method based on the HHO-CEEMDAM algorithm as described in claim 1 or 2, characterized in that, The process of step 3 is as follows: (3.1) Utilizing the linearity of arctangent demodulation, the phase of the point with maximum energy is extracted from the accumulated FFT results by traversal. (3.2) The phase unwrapping algorithm is used to detect and track changes in the signal phase, adjusting the phase value to a continuous range to eliminate discontinuities caused by periodic signal changes. The processing flow is as follows: calculate the phase difference Δθ between two adjacent phase points. The phase unwrapping expression is as follows: i i =θ i-1 +Δθ+2kπ (10) Where, when |Δθ|<π, k=0; when Δθ>π, k=-1; when Δθ<-π, k=1; (3.3) Using phase difference, the phase difference Δθ = θ between two adjacent points is calculated. i -θ i-1 This allows us to obtain the phase change of the target, which can then be used for subsequent target trajectory analysis and motion prediction.

5. The millimeter-wave radar vital sign detection method based on the HHO-CEEMDAM algorithm as described in claim 1 or 2, characterized in that, The process of step 4 is as follows: (4.1) Before using HHO to optimize the key parameters of CEEMDAN, the range of values ​​for each parameter must be given. The range of noise amplitude ε is 0.1 to 1.0, the range of noise addition number N is 50 to 200, and the range of maximum iteration number MaxIter is 100 to 500. (4.2) Using the minimum envelope entropy as the fitness function, for a given signal component x(t), its analytic signal is first obtained through Hilbert transform: z(t)=x(t)+jH[x(t)] (11) To eliminate the influence of amplitude, the envelope signal e(t) = |z(t)| needs to be normalized to obtain the probability distribution sequence: Choose the minimum envelope entropy value E e As the fitness function value F:

6. The millimeter-wave radar vital sign detection method based on the HHO-CEEMDAM algorithm as described in claim 1 or 2, characterized in that, The process of step 5 is as follows: (5.1) Within the range of parameter values, the HHO algorithm randomly generates a set of CEEMDAN parameters (ε, N, MaxIter) as the initial population; (5.2) Calculate the current energy factor Where E0 is the initial energy value, t is the current iteration number, and T is the maximum iteration number. The working mode is selected based on the |E| value. (5.3) When |E|≥1, the algorithm enters the exploration phase, updates the position, and randomly searches for potential optimal regions in the parameter space; Where X(t) represents the current position of the individual, X rand (t) represents the randomly selected individual location, X rabbit (t) represents the location of the prey, i.e., the current optimal solution, X. mean (t) represents the average position of the population, r1, r2, r3, r4 and q are random numbers, and LB and UB are the lower and upper bounds of the search space, respectively. (5.4) When |E|<1, the algorithm enters the development phase, adopting different hunting strategies based on the escape probability r and the energy state. r is a random number, and the strategies include the following four: When |E|≥0.5 and r≥0.5, a soft-encirclement strategy is adopted, and its position update formula is: X(t+1)=ΔX(t)-E|JX rabbit (t)-X(t)| (16) When |E| < 0.5 and r ≥ 0.5, a hard-bound strategy is used to update its position: X(t+1)=X rabbit (t)-E|Δx(t)| (17) When |E|≥0.5 and r<0.5, a gradual rapid dive strategy is used to update its position: X(t+1)=X rabbit (t)-E|JX rabbit (t)-X(t)| (18) When |E| < 0.5 and r < 0.5, a random raid strategy is used to update its position: X(t+1)=X rabbit (t)-E|JX rabbit (t)-X(t)|+Levy(D) (19) Where, ΔX(t)=X rabbit (t)-X(t) represents the position difference between the prey and the current individual, J is the random jump intensity, Levy(D) represents the random step size of Levy's flight, and D is the dimension of the problem; (5.5) Update the position of individuals in the population according to the hunting strategy and calculate the fitness value. Finally, stop the search and output the optimal parameter combination when the maximum number of iterations is reached or the fitness value converges.

7. The millimeter-wave radar vital sign detection method based on the HHO-CEEMDAM algorithm as described in claim 1 or 2, characterized in that, The process of step 6 is as follows: (6.1) Set m = 1 and set the initial residual r0(t) to the original signal x(t). Based on the optimized noise amplitude ε and the number of noise additions N, add Gaussian white noise ε0·ω to the original signal. i (t), where ω i (t) represents the noise added for the i-th time, i = 1, 2, ..., N; (6.2) Perform EMD decomposition on the noisy signal to obtain the first IMF component. Calculate the mean of all noisy signal decomposition results to obtain the final value of the first IMF component: The updated residual r1(t) is the original signal minus the first IMF component: r1(t)=x(t)-c1(t) (21) (6.3) For the m-th IMF component (m≥2), for the current residual r m-1 (t) Add noise ε m-1 ·E m-1 (ω i (t)), where E m-1 (·) represents the (m-1)th IMF component after EMD decomposition of the noise. EMD decomposition of the noisy signal yields the m-th IMF component. The final value of the m-th IMF component is obtained by calculating the mean of all noisy signal decomposition results: Then, update the residual r. m (t) represents the current residual minus the m-th IMF component: r m (t)=r m-1 (t)-c m (t) (23) (6.4) Residual r m The algorithm terminates when (t) becomes a monotonic function or can no longer be decomposed into IMF components. Based on the decomposed IMF components, the IMF components corresponding to 0.1–0.6 Hz are selected for respiratory signals, and the IMF components corresponding to 0.8–2 Hz are selected for heartbeat signals.

8. The millimeter-wave radar vital sign detection method based on the HHO-CEEMDAM algorithm as described in claim 1 or 2, characterized in that, The process of step 7 is as follows: (7.1) Perform spectral analysis on the selected IMF components and extract the dominant frequency components; (7.2) Obtain the respiratory rate R b and heart rate H b The estimated values ​​are all in BPM (times per minute). R b =f r ×60 (24) H b =f h ×60 (25) Among them, f r f is the dominant frequency of the detected respiratory signal. h The main frequency of the detected heartbeat signal is expressed in Hz.

Images