A high-precision heart rate estimation method and system based on millimeter-wave radar

Abstract

This invention proposes a high-precision heart rate estimation method and system based on millimeter-wave radar, belonging to the field of radar detection technology. The method includes: processing the target echo signal received by the millimeter-wave radar to extract a phase-time signal reflecting chest cavity vibration; extracting the respiratory fundamental frequency from the phase-time signal and generating a respiratory reference signal containing its harmonic components based on the respiratory fundamental frequency; and using an adaptive filter, with the respiratory reference signal as a reference input, filtering the phase-time signal to suppress the harmonic components generated by the respiratory fundamental frequency and separate the initial heartbeat signal. This invention enhances the high-frequency components of the heartbeat signal while suppressing noise amplification, significantly improving the quality of the heartbeat signal in low signal-to-noise ratio environments. Employing a high-resolution spectrum estimation method, it can clearly extract the heartbeat frequency spectrum peaks from the enhanced signal, thereby achieving high-precision and robust heart rate estimation in complex environments.

Description

Technical Field

[0001] This invention relates to the field of radar detection technology, and in particular to a high-precision heart rate estimation method and system based on millimeter-wave radar. Background Technology

[0002] Non-contact vital sign monitoring technology, particularly methods based on frequency-modulated continuous wave millimeter-wave radar, can continuously measure respiration and heart rate without the user's sensation by capturing and analyzing the phase information of radar echoes modulated by the subtle movements of the human chest cavity. This demonstrates significant potential in clinical monitoring, health tracking, and smart home applications. Compared to traditional contact sensors, this method avoids the discomfort and infection risks associated with skin contact, providing users with a more comfortable and convenient monitoring experience.

[0003] However, existing millimeter-wave radar heart rate estimation schemes still face severe challenges in terms of measurement accuracy and reliability in complex application scenarios. Specifically, the strong harmonic components generated by respiratory movements extend and cover the typical frequency band of the heartbeat signal, forming strong co-band interference. This makes methods based on traditional bandpass filtering or simple spectral peak search prone to misjudgment. At the same time, the heartbeat signal itself is extremely weak, and in low signal-to-noise ratio environments, its spectral peaks are easily submerged by noise. Conventional spectral analysis methods have limited frequency resolution, making it difficult to accurately separate and identify the spectral features characterizing heart rate from signals containing strong respiratory harmonic remnants and background noise. Although some adaptive signal decomposition techniques, such as variational mode decomposition, can be used in this scenario, their performance is heavily dependent on preset parameters. Parameter mismatch can easily lead to mode aliasing or signal distortion, lacking robustness and universality. Summary of the Invention

[0004] In view of this, this invention proposes a high-precision heart rate estimation method and system based on millimeter-wave radar. By constructing a reference signal based on the respiratory fundamental frequency and using an adaptive filter with adjustable parameters, it can achieve accurate tracking and cancellation of the respiratory fundamental frequency and harmonic components, fundamentally solving the problem of spectral interference of respiratory harmonics on the heartbeat frequency band. At the same time, through the synergistic processing of time-domain differential operation and multi-scale decomposition denoising, the high-frequency components of the heartbeat signal are enhanced while the noise amplification effect is suppressed, significantly improving the quality of the heartbeat signal in low signal-to-noise ratio environments. By adopting a high-resolution spectrum estimation method, the heartbeat frequency spectrum peaks can be clearly extracted from the enhanced signal, thereby achieving high-precision and high-robustness heart rate estimation in complex environments.

[0005] The technical solution of this invention is implemented as follows: On one hand, the present invention provides a high-precision heart rate estimation method based on millimeter-wave radar, comprising: The target echo signal received by the millimeter-wave radar is processed to extract the phase timing signal reflecting the vibration of the chest cavity; The respiratory fundamental frequency is extracted from the phase timing signal, and a respiratory reference signal containing its harmonic components is generated based on the respiratory fundamental frequency. An adaptive filter is used, with the respiratory reference signal as the reference input, to filter the phase timing signal to suppress the harmonic components generated by the respiratory fundamental frequency and separate the initial heartbeat signal. The parameters of the adaptive filter are adaptively adjusted according to the filtering error. The initial heartbeat signal is subjected to time-domain difference operation, multi-scale decomposition and denoising to obtain the enhanced target heartbeat signal; High-resolution spectral estimation of the target heartbeat signal is performed to obtain a pseudo-spectral map; Identify the main spectral peak within the heart rate frequency band of the pseudo-spectrum, and obtain the target heart rate value based on the frequency corresponding to the main spectral peak.

[0006] Based on the above technical solutions, preferably, the processing of the target echo signal received by the millimeter-wave radar to extract the phase timing signal reflecting the chest cavity vibration includes: The phasor mean elimination method is used to calculate the average value of the target echo signal in the slow time dimension, and then subtract the average value from each slow time pulse to remove static clutter and obtain the dynamic signal matrix. Perform a Fast Fourier Transform along the fast time dimension on the dynamic signal matrix to obtain a range profile containing target distance information; Within a preset distance range that includes the human body, the signal energy of each distance unit in the distance image is accumulated in the slow time dimension, and the distance unit with the largest accumulated energy is determined as the target distance unit. The phase is extracted from the signal of the target distance unit in the slow time dimension using arctangent operation, and then the phase jump caused by the periodicity of the arctangent function is corrected by the phase dewinding algorithm, thereby recovering the phase time sequence signal that is proportional to the displacement of the thoracic cavity.

[0007] Based on the above technical solution, preferably, after extracting the phase-time signal reflecting thoracic cavity vibration and before extracting the respiratory fundamental frequency from the phase-time signal, the method further includes adaptive baseline drift correction of the phase-time signal: A median filter with adaptively varying window length is used to filter the phase-time signal. The window length of the median filter is adaptively adjusted according to the local dynamic characteristics of the phase-time signal. The adaptive adjustment includes using a larger window length for filtering when the local segment of the phase-time signal is stable, and using a smaller window length for filtering when the local segment of the phase-time signal fluctuates violently. The output signal obtained by median filtering is used as the baseline trend line of the phase timing signal; Subtracting the original phase timing signal from the baseline trend line yields the baseline-corrected thoracic vibration signal.

[0008] Based on the above technical solution, preferably, the adaptive filter adopts a recursive least squares criterion and uses the respiratory reference signal as a reference input. The adaptive filter minimizes the error between the phase timing signal and the filtered output signal, generates a harmonic component estimate related to the respiratory reference signal, and achieves suppression by subtracting the harmonic component estimate from the phase timing signal, thereby obtaining the initial heartbeat signal. The forgetting factor of the adaptive filter is adaptively adjusted according to the filtering error. During the initial convergence phase of the filter or when the instantaneous error is large, the forgetting factor is reduced to accelerate convergence and improve the tracking ability of the harmonic components. When the error is small, the forgetting factor is increased to reduce the steady-state error and maintain stable suppression of the respiratory harmonic components.

[0009] Based on the above technical solutions, preferably, the step of performing time-domain difference operations, multi-scale decomposition, and denoising processing on the initial heartbeat signal includes: The initial heartbeat signal is subjected to a first-order difference operation to suppress residual low-frequency components and enhance the high-frequency energy of the heartbeat signal, thereby obtaining the intermediate heartbeat signal. Discrete wavelet transform is performed on the intermediate heartbeat signal to obtain multi-level wavelet coefficients; Thresholding is performed on the wavelet coefficients of each layer based on a preset threshold function to suppress high-frequency coefficients corresponding to noise and retain coefficients corresponding to the main scale of heartbeat signal energy. Wavelet reconstruction is performed on the processed coefficients to obtain the enhanced target heartbeat signal.

[0010] Based on the above technical solutions, preferably, high-resolution spectral estimation is performed on the target heartbeat signal to obtain a pseudo-spectral map, including: Construct the autocorrelation matrix of the target heartbeat signal; The autocorrelation matrix is ​​subjected to eigenvalue decomposition to obtain the signal subspace and noise subspace; By utilizing the orthogonality between the signal subspace and the noise subspace, a pseudo-spectrum with high frequency resolution is generated.

[0011] Based on the above technical solutions, preferably, the step of generating a respiratory reference signal containing its harmonic components according to the respiratory fundamental frequency is wherein the generated respiratory reference signal is composed of a linear combination of the fundamental component corresponding to the respiratory fundamental frequency and the second and higher harmonic components of the respiratory fundamental frequency, wherein the fundamental component and the harmonic components are both periodic signals in the form of sine waves or cosine waves.

[0012] On the other hand, the present invention provides a high-precision heart rate estimation system based on millimeter-wave radar, comprising: The signal preprocessing and phase extraction module is configured to process the target echo signal received by the millimeter-wave radar to extract the phase timing signal reflecting the thoracic cavity vibration. The respiratory reference signal generation module is configured to extract the respiratory fundamental frequency from the phase timing signal and generate a respiratory reference signal containing its harmonic components based on the respiratory fundamental frequency. The adaptive harmonic suppression module is configured to use an adaptive filter, with the respiratory reference signal as the reference input, to filter the phase timing signal to suppress the harmonic components generated by the respiratory fundamental frequency and separate the initial heartbeat signal; wherein, the parameters of the adaptive filter are adaptively adjusted according to the filtering error; The heartbeat signal enhancement module is configured to perform time-domain difference operations, multi-scale decomposition, and denoising on the initial heartbeat signal to obtain the enhanced target heartbeat signal. The high-resolution heart rate calculation module is configured to perform high-resolution spectrum estimation on the target heartbeat signal to obtain a pseudo-spectral map, identify the main spectral peak within the heart rate frequency band of the pseudo-spectral map, and obtain the target heart rate value based on the frequency corresponding to the main spectral peak.

[0013] On the other hand, the present invention provides an electronic device including a processor, a memory, a user interface, and a network interface. The memory is used to store instructions, the user interface and the network interface are used to communicate with other devices, and the processor is used to execute the instructions stored in the memory to cause the electronic device to perform the above-described method.

[0014] On the other hand, the present invention provides a non-transitory computer-readable storage medium having a computer program stored thereon, which, when executed by a processor, implements the steps of the above-described method.

[0015] The high-precision heart rate estimation method and system based on millimeter-wave radar of the present invention have the following advantages over the prior art: 1. By constructing a reference signal based on the respiratory fundamental frequency and using an adaptive filter with adjustable parameters, it is possible to accurately track and cancel the respiratory fundamental frequency and harmonic components, fundamentally solving the problem of spectral interference of respiratory harmonics on the heartbeat frequency band. At the same time, through the synergistic processing of time-domain differential operation and multi-scale decomposition denoising, the high-frequency components of the heartbeat signal are enhanced while the noise amplification effect is suppressed, significantly improving the quality of the heartbeat signal in low signal-to-noise ratio environments. Using a high-resolution spectrum estimation method, the heartbeat frequency spectrum peaks can be clearly extracted from the enhanced signal, thereby achieving high-precision and high-robustness heart rate estimation in complex environments. 2. By performing baseline drift correction on the phase timing signal based on adaptive window length median filtering, it is possible to distinguish between low-frequency trends caused by human micro-movements and vital sign signals. While completely removing baseline drift, it preserves the high-frequency transient characteristics of the heartbeat signal to the maximum extent, providing a clean input for subsequent processing. 3. By adopting a variable step size RLS adaptive filter, the system can dynamically adjust the convergence characteristics according to the instantaneous error, quickly track when the respiratory waveform changes abruptly, and finely suppress it in steady state, thereby achieving faster and more thorough cancellation of the strong harmonic components generated by the respiratory fundamental wave, and significantly improving the separation purity of the heartbeat signal. 4. By cascading first-order difference and wavelet threshold denoising, the high-frequency energy of the heartbeat signal is effectively enhanced while accurately eliminating high-frequency noise introduced by the difference and inherent in the environment, thus substantially improving the signal-to-noise ratio of the heartbeat signal. At the same time, high-resolution spectrum estimation is performed through the MUSIC algorithm. Utilizing the orthogonality of the signal and noise subspaces, a pseudo-spectrum with sharp spectral peaks can be generated under conditions of extremely low signal-to-noise ratio and short data, thereby achieving extremely high accuracy and stability in the identification and estimation of the heartbeat frequency. Attached Figure Description

[0016] To more clearly illustrate the technical solutions in the embodiments of the present invention or the prior art, the drawings used in the description of the embodiments or the prior art will be briefly introduced below. Obviously, the drawings described below are only some embodiments of the present invention. For those skilled in the art, other drawings can be obtained based on these drawings without creative effort.

[0017] Figure 1 This is a flowchart illustrating the high-precision heart rate estimation method based on millimeter-wave radar of the present invention. Figure 2 This is an example diagram showing the target heart rate value result of the high-precision heart rate estimation method based on millimeter-wave radar of the present invention. Detailed Implementation

[0018] To make the objectives, technical solutions, and advantages of this invention clearer, the technical solutions of this invention will be clearly and completely described below with reference to the accompanying drawings. Obviously, the described embodiments are only some, not all, of the embodiments of this invention. All other embodiments obtained by those skilled in the art based on the embodiments of this invention without creative effort are within the scope of protection of this invention.

[0019] like Figure 1-2As shown, the high-precision heart rate estimation method based on millimeter-wave radar of the present invention consists of four core stages in sequence: signal preprocessing and phase extraction, adaptive harmonic suppression and signal separation, heartbeat signal enhancement and purification, and high-resolution heart rate calculation.

[0020] The first stage involves static clutter suppression and target range cell localization of the raw radar echo, followed by phase demodulation and dewinding to obtain a continuous phase-series signal characterizing thoracic vibration. The second stage estimates the respiratory fundamental frequency based on this phase signal and constructs a respiratory reference signal containing harmonic components. An adaptively adjustable filter is used to efficiently cancel respiratory harmonics, thus obtaining the initial heartbeat signal. The third stage sequentially performs time-domain differential enhancement and multi-scale denoising on the initial heartbeat signal to improve its signal-to-noise ratio and suppress noise interference. The fourth stage finally performs high-resolution spectral analysis on the enhanced heartbeat signal, outputting an accurate heart rate estimate through precise peak searching and frequency conversion. The following will describe each step in detail with reference to specific examples.

[0021] The target echo signal received by the millimeter-wave radar is processed to extract the phase timing signal reflecting the vibration of the chest cavity.

[0022] First, the raw intermediate frequency signal acquired by the millimeter-wave radar is preprocessed to extract the high signal-to-noise ratio chest vibration phase timing signal. This process accurately locates the human target from the raw data containing environmental noise, static clutter, and various interferences, and calculates the phase change information modulated by chest micro-movements caused purely by vital signs, including breathing and heartbeat. This stage of processing is the foundation for all subsequent advanced signal processing and can be further subdivided into the following four sub-steps.

[0023] The phasor mean elimination method is used to calculate the average value of the target echo signal in the slow time dimension, and then subtract the average value from each slow time pulse to remove static clutter and obtain the dynamic signal matrix.

[0024] Static clutter primarily originates from stationary objects within the radar's field of view, as well as the DC bias and fixed-phase noise inherent in the radar hardware itself. Its intensity is typically much higher than the micro-motion signals caused by vital signs. To eliminate this interference, phasor averaging can be employed. Specifically, the average value of the echo signals from all pulses in the slow-time dimension is taken; this average value represents the complex value of the static background. Subsequently, this average static background value is subtracted from the complex signal of each slow-time pulse. This operation is equivalent to a high-pass filter with an extremely low cutoff frequency, which can almost completely remove static reflection components while preserving the dynamic signals of vital signs, thereby significantly improving the signal-to-noise ratio of the dynamic signal and laying the foundation for subsequent processing.

[0025] In one specific embodiment, the raw ADC data acquired by the radar is assumed to be a two-dimensional matrix. Where m represents the fast-time dimension sampling point index, corresponding to distance information, and n represents the slow-time dimension Chirp sequence index, corresponding to time information. The calculation process for static clutter suppression is shown in the following formula:

[0026] in, This represents the original signal at the nth chirp and the mth fast-time sampling point, where N is the total number of chirs in the slow-time dimension. The calculation involves taking the average value across all N chirps for each fixed fast-time sampling point m. This average value is the estimate of the static clutter. This refers to the target echo data matrix after removing static clutter.

[0027] After suppressing static clutter, it is necessary to accurately locate the signal region containing the reflection from the human chest cavity in the distance dimension. First, a Fast Fourier Transform is performed on the preprocessed signal matrix along the fast time dimension to transform the signal from the time domain to the distance domain, obtaining a one-dimensional distance image. This distance image reflects the reflection intensity of the target at different distance units. Within a preset distance range containing the human body, the signal energy of each distance unit in the distance image in the slow time dimension is accumulated, and the distance unit with the largest accumulated energy is determined as the target distance unit. From the signal of the target distance unit in the slow time dimension, the phase is extracted using arctangent operation, and then the phase jump caused by the periodicity of the arctangent function is corrected by the phase dewinding algorithm, thereby recovering a phase time sequence signal that is proportional to and continuous with the chest cavity displacement.

[0028] In one specific embodiment, the data matrix after removing static clutter... For each chirp, perform a fast Fourier transform along the fast time dimension m to convert the signal to the range-frequency domain, obtaining the range spectrum matrix. .

[0029] Then, for each distance gate index m, calculate its accumulated signal energy over all N chirs. Its calculation expression is: .

[0030] Iterate through all distance gate indices m to find the one that maximizes the accumulated energy value. The index that reaches the maximum value :

[0031] The index The corresponding distance gate is determined to be the optimal distance gate required for subsequent processing in this invention, i.e., the target distance unit.

[0032] After accurately locking the target distance unit, the phase signal reflecting the micro-motion information of the thoracic cavity needs to be extracted and corrected. Specifically, the instantaneous phase needs to be calculated by the arctangent demodulation method. However, the output range of the standard arctangent function arctan is limited to a certain range. When the actual phase change exceeds this range, a jump will occur. This phenomenon is called phase winding. The wound phase cannot truly reflect the continuous displacement change of the thoracic cavity and must be corrected, i.e., phase dewinding.

[0033] During phase unwinding, the periodic displacement information of the thoracic cavity caused by vital signs is modulated onto the phase of the radar echo signal. For each complex signal sampling point s(n) corresponding to a slow time index n, its wrapping phase is calculated using the four-quadrant arctangent function. :

[0034] Among them, the wrapping phase The range of values ​​is restricted to Within the radian interval, s(n) = I(n) + j·Q(n) represents the complex signal extracted from the target distance cell at the slow time index n. I(n) represents the in-phase component of the complex signal, Q(n) represents the quadrature component of the complex signal, and atan2(Q,I) is the arctangent function in the four quadrants.

[0035] For the package phase sequence Sequential processing is performed, and the difference between adjacent phase points is detected to identify and correct the issue. Jump, initialize the first value of the unwound phase sequence. For n=2 to N, iterative calculation is performed according to the following rules:

[0036] Among them, phase increment The determination is made by comparing the phase difference between adjacent packages with π, and the calculation rules are as follows: .

[0037] The signal obtained after phase dewinding It still includes low-frequency, nonlinear baseline drift caused by involuntary micro-movements of the human body (such as slight body swaying), which overlaps with the frequency band of vital signs signals. Therefore, after extracting the phase-time signal reflecting chest cavity vibration and before extracting the respiratory fundamental frequency from the phase-time signal, it also includes adaptive baseline drift correction of the phase-time signal.

[0038] In some specific embodiments, a median filter with adaptively varying window length is used to filter the phase timing signal. The window length of the median filter is adaptively adjusted according to the local dynamic characteristics of the phase timing signal. The adaptive adjustment includes using a larger window length for filtering when the local segment of the phase timing signal is stable, and using a smaller window length for filtering when the local segment of the phase timing signal fluctuates violently.

[0039] Set the minimum value L for the median filter window length. min and maximum value L max The window length is measured in slow-time index points, and can be set according to the signal sampling rate. ,in, For the slow time dimension sampling rate, the calculated values ​​need to be rounded to ensure that the window length is an odd number.

[0040] For phase signals For each data point n, take W points before and after it as the center, forming a local analysis window of length 2W+1. Calculate the gradient of the signal within this window. For example, using a first-order difference approximation:

[0041]

[0042] Subsequently, the variance of the gradient sequence is calculated. :

[0043] in, It is the average gradient within this local window, and the gradient variance. The magnitude of the signal directly reflects the degree of drastic change in the signal at the current moment.

[0044] For dynamically calculating the window size, it can be based on the calculated local gradient variance. The median filter window length L(n) applicable to the current point n is dynamically determined through a preset mapping function, for example, using a linear mapping function:

[0045] Where 'a' is a preset scaling factor used to control the degree to which the gradient variance affects the window size.

[0046] Finally, L(n) is rounded down to ensure it is odd and restricted to a certain value. Within the range.

[0047] After setting the adaptive window length median filter, for each data point in the phase signal Use its corresponding dynamically calculated window length To select neighboring data points, sort all data points within the neighborhood by numerical value, and take the median value as the filtered output value for the current point. This value represents the estimated baseline drift trend term.

[0048] The output signal obtained through median filtering is used as the baseline trend line for the phase time series signal. After setting the adaptive window length median filter, for each data point in the phase signal... Use its corresponding dynamically calculated window length To select neighboring data points, sort all data points within the neighborhood by numerical value, and take the median value as the filtered output value for the current point. This value represents the estimated baseline drift trend term.

[0049] Subtracting the baseline trend line from the original phase timing signal yields the baseline-corrected thoracic vibration signal. Subtracting this trend term from the original phase signal gives the pure thoracic vibration signal after removing baseline drift. : .

[0050] The respiratory fundamental frequency is extracted from the phase-time signal, and a respiratory reference signal containing its harmonic components is generated based on the respiratory fundamental frequency. The purpose of this step is to first separate and determine the fundamental frequency of respiratory motion from the extracted phase-time signal, and then actively construct a theoretical signal containing its harmonics based on this fundamental frequency as a template for filtering out respiratory interference.

[0051] Specifically, the respiratory reference signal containing its harmonic components is generated based on the respiratory fundamental frequency. The generated respiratory reference signal is formed by linearly combining the fundamental component corresponding to the respiratory fundamental frequency and the second and higher harmonic components of the respiratory fundamental frequency. The fundamental component and the harmonic components are periodic signals in the form of sine waves or cosine waves.

[0052] In one specific embodiment, the baseline-corrected signal Perform short-time Fourier transform or sliding window FFT analysis to dynamically estimate the fundamental frequency of the respiratory signal. Subsequently, a reference input vector is constructed based on this fundamental frequency. The reference vector typically contains the sine and cosine pairs of the first K harmonics of the respiratory fundamental frequency to ensure coverage of the frequency bands that the heartbeat signal may exist. Usually, K=2 or 3 can effectively cover frequencies above 2.0 Hz.

[0053] Reference vector The specific form is as follows:

[0054] Among them, T s The sampling interval is in the slow time dimension, K is the harmonic order under consideration, and the length M of the vector is 2K.

[0055] An adaptive filter is used, with the respiratory reference signal as the reference input, to filter the phase timing signal to suppress the harmonic components generated by the respiratory fundamental frequency and separate the initial heartbeat signal. The parameters of the adaptive filter are adaptively adjusted according to the filtering error.

[0056] In some embodiments, the adaptive filter employs a recursive least squares criterion and uses the respiratory reference signal as a reference input. The adaptive filter minimizes the error between the phase timing signal and the filtered output signal, generates a harmonic component estimate related to the respiratory reference signal, and achieves suppression by subtracting this harmonic component estimate from the phase timing signal, thereby obtaining the initial heartbeat signal. The forgetting factor of the adaptive filter is adaptively adjusted according to the filtering error. During the initial convergence phase of the filter or when the instantaneous error is large, the forgetting factor is reduced to accelerate convergence and improve the tracking ability of harmonic components. When the error is small, the forgetting factor is increased to reduce steady-state error and maintain stable suppression of respiratory harmonic components.

[0057] In some more specific embodiments, the parameters of the adaptive filter are first set, and the initial value of the weight vector can be set to the zero vector, i.e. Initialize the inverse correlation matrix: ,in I is a small regularization constant, typically between 0.01 and 0.1, used to ensure the initial stability of the algorithm. The identity matrix.

[0058] Then, based on empirical values ​​or other methods, the boundary value of the forgetting factor can be set, for example... Set to 0.95, Setting these boundaries to 0.998 ensures that the forgetting factor always varies within a reasonable and effective range.

[0059] Set sensitivity parameters This parameter controls the sensitivity of the error to the adjustment of the forgetting factor. The larger the value, the more sensitive the forgetting factor is to changes in error. A suitable value can be determined through simulation experiments, for example, by setting it to 10.

[0060] For each moment Perform the following recursive operation: a. Update forgetting factors: The forgetting factor at the current time is dynamically calculated based on the squared filtering error from the previous time step. :

[0061] in, It is the prior error at time n-1. This formula ensures that: when the system initially converges or a sudden change occurs, the error... When it is large, Automatically decreases and approaches This assigns higher weights to new data, thereby accelerating convergence and improving the ability to track changes; as the algorithm stabilizes and the error... When smaller, Automatically increases and approaches This enhances the memory of past data, thereby achieving smaller steady-state errors and better harmonic suppression.

[0062] b. Calculate the gain vector :

[0063] Where k(n) is the gain vector at the current time (n), P(n-1) is the estimate of the inverse correlation matrix at the previous time, and u H (n) is the conjugate transpose of vector u(n).

[0064] c. Calculate the prior error e(n):

[0065] Among them, expected response That is, the current input signal. This error reflects the difference between the current filter output and the desired output, w H (n-1) is the conjugate transpose of the weight vector w(n-1).

[0066] d. Update the filter weight vector w(n):

[0067] For real signal systems, complex conjugate That is This step itself adjusts the filter coefficients based on the error and gain.

[0068] e. Update the inverse correlation matrix P(n):

[0069] This step is the core of the RLS algorithm, which recursively updates the relevant information of the input signal.

[0070] In the above iteration process, the prior error This is the output of the filter at time n, which represents the output of the mixed signal. In order to cancel out the reference signal as much as possible The simulated respiration and its harmonic components are followed by the remaining signal; therefore, the sequence calculated throughout the iteration process is used as the initially purified heartbeat signal. : .

[0071] The initial heartbeat signal is subjected to time-domain difference operation, multi-scale decomposition and denoising processing to obtain the enhanced target heartbeat signal.

[0072] The initial heartbeat signal obtained after adaptive harmonic suppression usually still contains residual extremely low frequency trend terms, residual respiratory harmonics that are not completely canceled, and amplified measurement noise. This stage effectively solves the above problems through the cascade processing of time-domain differential enhancement and multi-scale denoising, providing high-quality input for the final high-precision frequency estimation.

[0073] In some specific embodiments, the initial heartbeat signal is first subjected to a first-order differential operation to suppress residual low-frequency components and enhance the high-frequency energy of the heartbeat signal, thereby obtaining an intermediate heartbeat signal.

[0074] The first-order difference operation is essentially a high-pass filter. Its core function is to further suppress the extremely low-frequency components remaining in the initial heartbeat signal, while simultaneously highlighting the relatively high-frequency characteristics of the heartbeat signal itself. Its mathematical expression is:

[0075] in, It is the sampled value of the preliminarily purified heartbeat signal output from the aforementioned steps at time n. It is the sampled value of the signal at the previous moment. It is the new signal obtained after differential operation, representing the change in signal between adjacent time points, and is also the intermediate heartbeat signal in this scheme.

[0076] This operation is equivalent to a high-pass filter in the frequency domain, and its frequency response function is: The amplitude-frequency characteristic shows that the gain increases with increasing frequency. This characteristic brings two key benefits: First, the heartbeat signal is a high-frequency component compared to the respiratory signal, typically with a heartbeat frequency > 0.8 Hz and a respiratory frequency < 0.5 Hz. Differential operations amplify high-frequency components, thus effectively enhancing the amplitude of the heartbeat signal and making it easier to detect in subsequent processing. Second, the operation also further attenuates any residual extremely low-frequency trend terms and near-DC components in the signal. These components may originate from slow baseline drift that was not completely eliminated by the previous filtering step or leakage of the respiratory fundamental frequency.

[0077] Through phase difference operations, the time-domain variation characteristics of the heartbeat signal become more apparent, and the signal-to-noise ratio is further improved. This creates more favorable conditions for subsequent noise suppression steps and the final high-resolution frequency estimation, resulting in a higher output intermediate heartbeat signal. It will be used for subsequent wavelet transform denoising.

[0078] Discrete wavelet transform is performed on the intermediate heartbeat signal to obtain multi-level wavelet coefficients. The signal is then enhanced through phase difference operation. While the heartbeat component is amplified, the differential operation, while boosting the high-frequency signal, inevitably amplifies environmental and system noise present in the high-frequency band. To further optimize the signal-to-noise ratio after signal enhancement and accurately extract heartbeat features, this embodiment introduces wavelet transform denoising as the final post-processing step. Wavelet transform has the characteristics of multi-resolution analysis, enabling simultaneous localized analysis of signals in both the time and frequency domains, making it highly suitable for denoising non-stationary signals such as vital signs. This invention preferably uses the Daubechies 4 (db4) wavelet as the wavelet basis function because its waveform characteristics are relatively close to the shape of the heartbeat signal, which is beneficial for better preserving signal features during the decomposition process.

[0079] In some specific embodiments, the differentially enhanced signal Perform multi-resolution wavelet decomposition. Set the decomposition level to J=5, and use wavelet transform to decompose the signal into sub-band signals within different frequency ranges:

[0080]

[0081] in, These are the approximation coefficients of the Jth layer, representing the coarsest low-frequency trend in the signal. arrive It is the first The detail coefficients from layer 1 to layer 1 represent different frequency band components of the signal from low frequency to high frequency. The main energy of the heartbeat signal is usually concentrated in a specific intermediate frequency detail layer, while the higher frequency detail layers mainly contain noise.

[0082] Thresholding is applied to the wavelet coefficients of each layer based on a preset threshold function to suppress high-frequency coefficients corresponding to noise, while retaining coefficients corresponding to the main energy scale of the heartbeat signal. This step aims to suppress coefficients representing noise while preserving coefficients representing the true heartbeat signal. For each layer of detail coefficients... The threshold for this layer is adaptively determined using an unbiased risk estimation threshold rule ('rigrsure'). This rule is based on the principle of Stan's unbiased risk estimation. It obtains the optimal threshold by minimizing a risk function, and is particularly suitable for situations where noise characteristics are unknown, exhibiting strong adaptability.

[0083] A soft thresholding function is used to process detail coefficients. The soft thresholding function not only selects coefficients whose absolute values ​​are less than the threshold... Setting the coefficient to zero will also shrink coefficients with absolute values ​​greater than the threshold towards zero. The quantity, its mathematical expression is:

[0084] in, It is the k-th wavelet coefficient. These are the coefficients after thresholding. Soft thresholding functions provide smoother signal estimates than hard thresholding functions, which may introduce artificial jitter. sign() is responsible for preserving the sign of the original input.

[0085] Wavelet reconstruction is performed on the processed coefficients to obtain the enhanced target heartbeat signal. The thresholded detail coefficients from each layer are then used. and the original J-th layer approximation coefficients Perform inverse wavelet transform to reconstruct the denoised time-domain signal. : .

[0086] High-resolution spectral estimation of the target heartbeat signal is performed to obtain a pseudo-spectral map.

[0087] In order to enhance and denoise the heartbeat signal after the aforementioned steps Accurately estimating heart rate requires high-resolution spectral analysis methods. Traditional periodogram methods based on fast Fourier transform have frequency resolution limited by the length of data records, and their performance degrades under low signal-to-noise ratio and short data conditions. Therefore, this embodiment adopts the Multiple Signal Classification (MUSIC) algorithm, a high-resolution spectral estimation method based on subspace decomposition, which can still produce sharp pseudo-spectral peaks under the above-mentioned unfavorable conditions, thereby significantly improving the accuracy and reliability of frequency estimation.

[0088] In one specific embodiment, the autocorrelation matrix of the target heartbeat signal is first constructed, and then derived from the final heartbeat signal. Estimate its autocorrelation function, and construct a [structure / system] using this autocorrelation function. Hermitian autocorrelation matrix R:

[0089] Here, U is a unitary matrix whose column vectors are the eigenvectors of R. It is a diagonal matrix whose diagonal elements are the eigenvalues ​​of R, arranged in descending order, i.e. .

[0090] Eigenvalue decomposition is performed on the autocorrelation matrix to obtain the signal and noise subspaces. The eigenvalues ​​are then sorted and analyzed to determine the dimensions of the signal and noise subspaces. Assuming the signal contains P complex sinusoidal signals (i.e., the frequency components to be estimated; in this application, we mainly focus on the principal component of heart rate, but it may contain a small amount of harmonics or residual interference; P is usually a small value, such as 1~3), then the eigenvectors corresponding to the first P largest eigenvalues ​​are... The remaining MP smaller and similarly sized eigenvectors correspond to the eigenvectors. Zhang Cheng noise subspace .

[0091] Right now: .

[0092] By leveraging the orthogonality between the signal and noise subspaces, a pseudo-spectral map with high frequency resolution is generated. The core idea of ​​the MUSIC algorithm is to utilize the orthogonality between the steering vectors corresponding to the signal frequency components and the noise subspace. Within the frequency range to be detected (e.g., 0.8 Hz to 2.5 Hz), the frequency is scanned at a certain step size f. For each candidate frequency, its corresponding complex sinusoidal steering vector is constructed. :

[0093]

[0094] Where j is the unit of description, T s Let M be the sampling interval of the signal, and M be the order of the signal model, e j2πfTs It is the basic complex sine factor.

[0095] Calculate the MUSIC pseudospectral value at this frequency point. :

[0096] Where a(f) is the steering vector corresponding to frequency f, a H (f) is the conjugate transpose of the guiding vector a(f), U N To estimate the noise subspace matrix, U N H For the noise subspace matrix U N The conjugate transpose of U N U N H Let be the projection matrix of the noise subspace.

[0097] Since the signal subspace and noise subspace are orthogonal, when the scanning frequency f equals the actual frequency components present in the signal, the steering vector... It will be approximately located in the signal subspace, thus intersecting with the noise subspace. Approximately orthogonal, such that the denominator is... Approaching zero, ultimately leading to A sharp peak is generated at this frequency.

[0098] The main spectral peak is identified within the heart rate frequency band of the pseudo-spectrum, and the target heart rate value is obtained based on the frequency corresponding to the main spectral peak. This is done after obtaining the MUSIC pseudo-spectrum. Afterwards, a spectral peak search is needed to determine the heart rate.

[0099] In one specific embodiment, the search range is first limited to a reasonable heart rate range based on common physiological knowledge, for example... The corresponding heart rate is 48 to 150 bpm, which can effectively avoid misinterpreting other non-heartbeat peaks in the pseudo-spectrum as heart rate.

[0100] Within the set heart rate search range Inside, searching for pseudo-music spectra. The global maximum value, i.e., the highest peak:

[0101] Here, argmax is a function used to find the index of the maximum element, and this peak value... The corresponding frequency is then identified as the final estimated heart rate.

[0102] The estimated heart rate Unit: Hz, converted to clinically accepted heart rate values ​​HR (beats per minute, bpm): .

[0103] The high-precision heart rate estimation system based on millimeter-wave radar of the present invention includes a signal preprocessing and phase extraction module, a respiratory reference signal generation module, an adaptive harmonic suppression module, a heartbeat signal enhancement module, and a high-resolution heart rate calculation module.

[0104] The signal preprocessing and phase extraction module is configured to process the target echo signal received by the millimeter-wave radar to extract the phase timing signal reflecting the thoracic cavity vibration.

[0105] The respiratory reference signal generation module is configured to extract the respiratory fundamental frequency from the phase timing signal and generate a respiratory reference signal containing its harmonic components based on the respiratory fundamental frequency.

[0106] The adaptive harmonic suppression module is configured to use an adaptive filter, with the respiratory reference signal as the reference input, to filter the phase timing signal to suppress the respiratory fundamental frequency and its harmonic components, and to separate the initial heartbeat signal; wherein, the parameters of the adaptive filter are adaptively adjusted according to the filtering error.

[0107] The heartbeat signal enhancement module is configured to perform time-domain difference operations, multi-scale decomposition, and denoising on the initial heartbeat signal to obtain the enhanced target heartbeat signal.

[0108] The high-resolution heart rate calculation module is configured to perform high-resolution spectrum estimation on the target heartbeat signal to obtain a pseudo-spectral map, identify the main spectral peak within the heart rate frequency band of the pseudo-spectral map, and obtain the target heart rate value based on the frequency corresponding to the main spectral peak.

[0109] The electronic device of the present invention includes a processor, a memory, a user interface, and a network interface. The memory is used to store instructions, the user interface and the network interface are used to communicate with other devices, and the processor is used to execute the instructions stored in the memory to cause the electronic device to perform the above-described method.

[0110] The present invention provides a non-transitory computer-readable storage medium having a computer program stored thereon, which, when executed by a processor, implements the steps of the above-described method.

[0111] It should be noted that, for the sake of simplicity, the foregoing method embodiments are all described as a series of actions. However, those skilled in the art should understand that this application is not limited to the described order of actions, as some steps may be performed in other orders or simultaneously according to this application. Furthermore, those skilled in the art should also understand that the embodiments described in the specification are preferred embodiments, and the actions and modules involved are not necessarily essential to this application.

[0112] In the above embodiments, the descriptions of each embodiment have different focuses. For parts not described in detail in a certain embodiment, please refer to the relevant descriptions in other embodiments.

[0113] In the several embodiments provided in this application, it should be understood that the disclosed apparatus can be implemented in other ways. For example, the apparatus embodiments described above are merely illustrative; for instance, the division of units is only a logical functional division, and in actual implementation, there may be other division methods. For example, multiple units or components may be combined or integrated into another system, or some features may be ignored or not executed. Furthermore, the shown or discussed mutual couplings or direct couplings or communication connections may be through some service interfaces; indirect couplings or communication connections between apparatuses or units may be electrical or other forms.

[0114] The units described as separate components may or may not be physically separate. The components shown as units may or may not be physical units; that is, they may be located in one place or distributed across multiple network units. Some or all of the units can be selected to achieve the purpose of this embodiment according to actual needs.

[0115] Furthermore, the functional units in the various embodiments of this application can be integrated into one processing unit, or each unit can exist physically separately, or two or more units can be integrated into one unit. The integrated unit can be implemented in hardware or as a software functional unit.

[0116] If the integrated unit is implemented as a software functional unit and sold or used as an independent product, it can be stored in a computer-readable storage device (CMD). Based on this understanding, the technical solution of this application, in essence, or the part that contributes to the prior art, or all or part of the technical solution, can be embodied in the form of a software product. This computer software product is stored in a memory and includes several instructions to cause a computer device (which may be a personal computer, server, or network device, etc.) to execute all or part of the steps of the methods of the various embodiments of this application. The aforementioned memory includes various media capable of storing program code, such as USB flash drives, portable hard drives, magnetic disks, or optical disks.

[0117] The above description is only a preferred embodiment of the present invention and is not intended to limit the present invention. Any modifications, equivalent substitutions, improvements, etc., made within the spirit and principles of the present invention should be included within the protection scope of the present invention.

Claims

1. A high-precision heart rate estimation method based on millimeter-wave radar, characterized in that, include: The target echo signal received by the millimeter-wave radar is processed to extract the phase timing signal reflecting the vibration of the chest cavity; The respiratory fundamental frequency is extracted from the phase timing signal, and a respiratory reference signal containing its harmonic components is generated based on the respiratory fundamental frequency. An adaptive filter is used, with the respiratory reference signal as the reference input, to filter the phase timing signal to suppress the harmonic components generated by the respiratory fundamental frequency and separate the initial heartbeat signal. The parameters of the adaptive filter are adaptively adjusted according to the filtering error. The initial heartbeat signal is subjected to time-domain difference operation, multi-scale decomposition and denoising to obtain the enhanced target heartbeat signal; High-resolution spectral estimation of the target heartbeat signal is performed to obtain a pseudo-spectral map; Identify the main spectral peak within the heart rate frequency band of the pseudo-spectrum, and obtain the target heart rate value based on the frequency corresponding to the main spectral peak.

2. The high-precision heart rate estimation method based on millimeter-wave radar as described in claim 1, characterized in that, The process of processing the target echo signal received by the millimeter-wave radar to extract the phase timing signal reflecting the chest cavity vibration includes: The phasor mean elimination method is used to calculate the average value of the target echo signal in the slow time dimension, and then subtract the average value from each slow time pulse to remove static clutter and obtain the dynamic signal matrix. Perform a Fast Fourier Transform along the fast time dimension on the dynamic signal matrix to obtain a range profile containing target distance information; Within a preset distance range that includes the human body, the signal energy of each distance unit in the distance image is accumulated in the slow time dimension, and the distance unit with the largest accumulated energy is determined as the target distance unit. The phase is extracted from the signal of the target distance unit in the slow time dimension using arctangent operation, and then the phase jump caused by the periodicity of the arctangent function is corrected by the phase dewinding algorithm, thereby recovering the phase time sequence signal that is proportional to the displacement of the thoracic cavity.

3. The high-precision heart rate estimation method based on millimeter-wave radar as described in claim 1, characterized in that, After extracting the phase-time signal reflecting thoracic vibration, and before extracting the respiratory fundamental frequency from the phase-time signal, the method further includes adaptive baseline drift correction of the phase-time signal: A median filter with adaptively varying window length is used to filter the phase-time signal. The window length of the median filter is adaptively adjusted according to the local dynamic characteristics of the phase-time signal. The adaptive adjustment includes using a larger window length for filtering when the local segment of the phase-time signal is stable, and using a smaller window length for filtering when the local segment of the phase-time signal fluctuates violently. The output signal obtained by median filtering is used as the baseline trend line of the phase timing signal; Subtracting the original phase timing signal from the baseline trend line yields the baseline-corrected thoracic vibration signal.

4. The high-precision heart rate estimation method based on millimeter-wave radar as described in claim 1, characterized in that, The adaptive filter employs a recursive least squares criterion and uses the respiratory reference signal as a reference input. The adaptive filter minimizes the error between the phase-time signal and the filtered output signal, generating a harmonic component estimate related to the respiratory reference signal. This harmonic component estimate is subtracted from the phase-time signal to achieve suppression, thereby obtaining the initial heartbeat signal. The forgetting factor of the adaptive filter is adaptively adjusted based on the filtering error. During the initial convergence phase or when the instantaneous error is large, the forgetting factor is reduced to accelerate convergence and improve the tracking ability of harmonic components. When the error is small, the forgetting factor is increased to reduce steady-state error and maintain stable suppression of respiratory harmonic components.

5. The high-precision heart rate estimation method based on millimeter-wave radar as described in claim 1, characterized in that, The process of performing time-domain difference operations, multi-scale decomposition, and denoising on the initial heartbeat signal includes: The initial heartbeat signal is subjected to a first-order difference operation to suppress residual low-frequency components and enhance the high-frequency energy of the heartbeat signal, thereby obtaining the intermediate heartbeat signal. Discrete wavelet transform is performed on the intermediate heartbeat signal to obtain multi-level wavelet coefficients; Thresholding is performed on the wavelet coefficients of each layer based on a preset threshold function to suppress high-frequency coefficients corresponding to noise and retain coefficients corresponding to the main scale of heartbeat signal energy. Wavelet reconstruction is performed on the processed coefficients to obtain the enhanced target heartbeat signal.

6. The high-precision heart rate estimation method based on millimeter-wave radar as described in claim 1, characterized in that, High-resolution spectral estimation of the target heartbeat signal is performed to obtain a pseudo-spectral map, including: Construct the autocorrelation matrix of the target heartbeat signal; The autocorrelation matrix is ​​subjected to eigenvalue decomposition to obtain the signal subspace and noise subspace; By utilizing the orthogonality between the signal subspace and the noise subspace, a pseudo-spectrum with high frequency resolution is generated.

7. The high-precision heart rate estimation method based on millimeter-wave radar as described in claim 1, characterized in that, The method of generating a respiratory reference signal containing its harmonic components based on the respiratory fundamental frequency is wherein the generated respiratory reference signal is composed of a linear combination of the fundamental component corresponding to the respiratory fundamental frequency and the second and higher harmonic components of the respiratory fundamental frequency, wherein the fundamental component and the harmonic components are periodic signals in the form of sine waves or cosine waves.

8. A high-precision heart rate estimation system based on millimeter-wave radar, characterized in that, include: The signal preprocessing and phase extraction module is configured to process the target echo signal received by the millimeter-wave radar to extract the phase timing signal reflecting the thoracic cavity vibration. The respiratory reference signal generation module is configured to extract the respiratory fundamental frequency from the phase timing signal and generate a respiratory reference signal containing its harmonic components based on the respiratory fundamental frequency. The adaptive harmonic suppression module is configured to use an adaptive filter, with the respiratory reference signal as the reference input, to filter the phase timing signal to suppress the harmonic components generated by the respiratory fundamental frequency and separate the initial heartbeat signal; wherein, the parameters of the adaptive filter are adaptively adjusted according to the filtering error; The heartbeat signal enhancement module is configured to perform time-domain difference operations, multi-scale decomposition, and denoising on the initial heartbeat signal to obtain the enhanced target heartbeat signal. The high-resolution heart rate calculation module is configured to perform high-resolution spectrum estimation on the target heartbeat signal to obtain a pseudo-spectral map, identify the main spectral peak within the heart rate frequency band of the pseudo-spectral map, and obtain the target heart rate value based on the frequency corresponding to the main spectral peak.

9. An electronic device, characterized in that, The device includes a processor, a memory, a user interface, and a network interface. The memory is used to store instructions, the user interface and the network interface are used to communicate with other devices, and the processor is used to execute the instructions stored in the memory to cause the electronic device to perform the method as described in any one of claims 1-7.

10. A non-transitory computer-readable storage medium having a computer program stored thereon, characterized in that, When the computer program is executed by a processor, it implements the steps of the method described in any one of claims 1-7.

Images