Heart rate extraction method based on adaptive comb filter

Abstract

The application discloses a heart rate extraction method based on an adaptive comb filter. The method comprises the following steps: 1, obtaining a life signal by processing a radar echo; 2, extracting a breathing rate and other features therefrom, and initializing parameters of the filter; 3, designing an ALE module for evaluating a filtering effect, and initializing parameters of the module; 4, passing the life signal through the adaptive comb filter, and iterating until the parameters to be optimized are basically converged; and 5, selecting a heartbeat frequency from a spectrum of an output signal of the life signal through the filter. The heart rate extraction method based on the adaptive comb filter considers the intermodulation characteristics of breathing and heart rate components of the life signal in a Doppler radar system, and effectively reduces the interference of breathing and higher harmonics on the selection of the heart rate. The experimental results show that when the higher harmonics of breathing in the life signal affect the measurement of the heart rate, the adaptive comb filter can accurately measure the heart rate of a person to be measured.

Description

Heart Rate Extraction Method Based on Adaptive Comb Filter Technical Field

[0001] This invention belongs to the field of signal processing technology, and in particular relates to a heart rate extraction method based on an adaptive comb filter. Background Technology

[0002] As people's living standards improve, the demand for quick and comfortable access to vital signs information such as respiration and heart rate is increasing. Low-cost Doppler radar, which utilizes the scattering characteristics of electromagnetic waves to acquire vital signs of individual targets, has been gradually applied in non-contact vital sign detection scenarios. However, in reality, accurately extracting vital sign information from radar echoes still faces significant challenges.

[0003] In existing extraction methods, due to the intermodulation effect of respiration and heart rate, a large number of harmonic spectra exist within the detection range of vital sign characteristic line spectra. For example, high-order harmonics of the respiratory frequency may interfere with the identification of the heartbeat frequency in the baseband spectrum of vital sign radar, making it difficult to accurately determine its frequency position. To combat the interference effect of high-order respiratory harmonics on heart rate measurement, some existing wavelet transform-based respiratory and heartbeat signal detection algorithms have successfully separated the respiratory and heartbeat components. However, the performance of these algorithms is highly dependent on the selection of the wavelet basis and cannot provide target distance information. There are also methods that use approximation function methods to separate respiratory and heartbeat signals, approximating the respiratory signal with an estimated time-domain approximation function. However, the extracted heartbeat signal waveform is poor, and the algorithm depends on the selection of the approximation function. In recent years, an adaptive heartbeat detection algorithm using the averaging removal method to remove clutter has been proposed. This method is simple and fast, but the introduced harmonic interference has a significant impact on the target scattering signal waveform.

[0004] This invention proposes an adaptive digital comb filter with strong respiratory harmonic suppression, specifically addressing the characteristics of vital signals. To address the common issues of frequency shift, diffusion, and unpredictable amplitude of higher respiratory harmonics, it employs an adaptive filter as a notch filter gain control module and uses an adaptive iterative algorithm to continuously update its weighting coefficients to optimize harmonic filtering. By adaptively adjusting the filter bandwidth and gain control parameters, it achieves better filtering of higher respiratory harmonics, thereby enabling more effective and accurate estimation of heart rate. Summary of the Invention

[0005] Technical Problem: The purpose of this invention is to provide a heart rate extraction method based on an adaptive comb filter, in order to solve the technical problem that it is difficult to accurately extract vital signs features using existing processing methods when the amplitude and position of interference items such as respiratory high-order harmonics in the detected vital signals are unstable and cause unpredictable effects on adjacent heart rate measurements.

[0006] Technical Solution: To solve the above-mentioned technical problems, the present invention provides a heart rate extraction method based on an adaptive comb filter, as follows:

[0007] Heart rate extraction involves separating the heart rate component from the vital signs and the high-harmonic components of the respiratory rate, which can cause interference. The specific steps include:

[0008] Step 1: Obtain the baseband signal of the quadrature sampling received by the Doppler radar, and process the baseband signal to obtain the vital signals containing breathing and heartbeat information;

[0009] Step 2: Based on the characteristics of the life signal mentioned above, adopt an adaptive comb filter with adjustable notch width and gain, and initialize the parameters in the comb filter according to the signal conditions;

[0010] Step 3: Adopt the structure of the adaptive line enhancer (ALE) for adaptively adjusting the parameters of the comb filter, and initialize the relevant parameters of the ALE;

[0011] Step 4: Obtain the update expressions for the notch bandwidth control parameter b, the tap coefficient g(n) of the notch gain control module, and the tap coefficient h(n) of the signal periodic component extraction module. Then, pass the demodulated vital sign signal waveform through an adaptive comb filter and iterate until the parameters to be optimized are basically converged.

[0012] Step 5: Select the heartbeat frequency from the spectrum of the output signal after passing through the adaptive comb filter.

[0013] Step 1 specifically includes the following steps:

[0014] Step 1.1: The Doppler radar transmitting antenna emits a single-frequency continuous wave electromagnetic wave signal to illuminate the human target. The signal returns to the receiving antenna, and the echo signal reaching the receiving antenna is written as:

[0015] x R (t)=σcos(2πf c (t-τ))

[0016] The Doppler radar hardware system uses the transmitted signal to mix the echo signal and then performs low-pass filtering to obtain the mixed signal, which can be written as:

[0017]

[0018] Among them: A T f is the signal amplitude. c Where is the carrier frequency, t is time, σ is the target's scattering coefficient, τ = 2d(t) / c is the target's echo delay, d(t) is the distance the radar reaches the target, and c is the electromagnetic wave propagation speed.

[0019] The mixed signal is subjected to quadrature sampling. Considering the quadrature sampling in the actual hardware system, i.e., the DC bias introduced by the inconsistency of the IQ sampling channel levels, the baseband signals I(t) and Q(t) obtained by the IQ quadrature sampling channels are as follows:

[0020]

[0021]

[0022] DC I DC Q These are the DC components of the in-phase and quadrature components in the baseband signal, respectively. A I A Q Let λ = c / f be the amplitudes of the in-phase and quadrature components in the baseband signal, respectively. c Indicates the carrier frequency wavelength;

[0023] If the sampling frequency f of the Doppler radar system s If the observation time is T, then the discrete sampled signal of the IQ channel is obtained as follows:

[0024]

[0025]

[0026] Where n = 1, ..., N, N = f s T represents the number of signal sampling points;

[0027] Doppler radar extracts respiratory and heart rate characteristics by sensing the relatively weak periodic distance changes in the chest wall of a human target. The distance d(n) can be expressed as:

[0028]

[0029] Where d0 is the nominal distance of the radar to the target, H i and R j f represents the amplitude of the periodic displacement of the chest wall caused by a certain harmonic of the target heartbeat and respiration, respectively. h and f r These represent the target heart rate and respiratory rate in the human body, respectively.

[0030] Step 1.2: Estimate the mean of the baseband signal and

[0031] The size of the estimated baseband signal mean processing block is N.

[0032] The method for processing the estimated baseband signal mean is as follows:

[0033]

[0034]

[0035] Step 1.3: Use mean filtering to remove the DC component of the baseband signal:

[0036] The mean filtering method is as follows:

[0037]

[0038]

[0039] Considering that amplitude distortion caused by inconsistency between the I and Q channels is usually small, and assuming that the signal amplitudes of the I and Q channels are basically the same after mean filtering, the above formula can be rewritten as:

[0040]

[0041]

[0042] Where A is the signal amplitude after mean filtering;

[0043] Step 1.4: Using the mean-filtered signal as the observable, construct the complex signal C(n):

[0044]

[0045] Where j is an imaginary number.

[0046] Step 2 includes the following steps:

[0047] Step 2.1: Initialize the filter order parameter N. The system function of the adaptive comb filter is:

[0048]

[0049] In the above formula, N represents the order of the comb filter, b represents the notch bandwidth control parameter of the filter, z is the z-transform operator, G(z) is the gain control function of the adaptive comb filter, and P is its order; for a comb filter of order N, its m-th notch frequency is

[0050]

[0051] It is evident that the first notch frequency w1 is entirely determined by N. If w1 is the respiratory frequency, the comb filter can filter out respiratory frequencies and their integer multiples of harmonic frequencies. However, doing so directly may also affect heart rate components that are close to a certain integer multiple of the respiratory harmonic frequency. Performing FFT processing on the complex data obtained in the previous step reveals that the respiratory frequency and its higher harmonics are more significant in the spectrum. The respiratory frequency, denoted as f, is selected from frequencies falling within [0.15Hz, 0.5Hz] and exhibiting a significant harmonic effect. r Then N is initialized as:

[0052] N = f s / f r

[0053] Where f s The sampling frequency;

[0054] Step 2.2: Initialize the filter gain function G(z)

[0055] The gain function G(z) ensures that the comb filter has appropriate fading at different frequencies. It is an FIR filter of order 2P+1. A lower order can reduce the computational cost of subsequent iterations, but the accuracy will be compromised. Considering system stability, P≤N must be satisfied. The gain of the traditional comb filter is a constant. Replacing it with a filter G(z) with tap weights g(n), the amplitude-frequency response of the final comb filter at different frequencies can be adjusted by adjusting g(n). It is desired that the amplitude-frequency response of the final comb filter has fading suitable for the actual situation of complex data at each breathing frequency harmonic point. The initialization of G(z) can be accomplished using the window function design method. Let the height of the harmonic point w of a certain breathing frequency of the complex data be |W(e jw )|, hoping it will fade to |S(e) after passing through the filter. jw If |, then the filter's amplitude-frequency response at w is

[0056]

[0057] Then the i-th element g of the tap coefficient g(n) of G(z) i for

[0058]

[0059] Where -P≤i≤P; the relationship between G(z) and g(n) is as follows:

[0060]

[0061] Step 2.3: Initialize parameter b

[0062] The parameter b controls the notch width, which is the width of each "comb tooth". Its value range is [0,1]. As b approaches 1, the notch band will become narrower. For actual data, it is more appropriate to assign an initial value to b within [0.2,0.6].

[0063] Step 3 includes the following steps:

[0064] Step 3.1: Design the adaptive spectral line enhancer module

[0065] The Adaptive Line Enhancer (ALE) is used to evaluate the filtering effect of the previous comb filter module and obtains the update expressions for several parameters to be updated based on its output e(n).

[0066] The input signal x(n) of the comb filter can be expressed as:

[0067] x(n) = d(n) + i(n)

[0068] Where d(n) is the desired useful signal and i(n) is the periodic interference signal, the filtered output y(n) after passing through the previous comb filter module is expressed as:

[0069]

[0070] in As the residual of the periodic interference signal, y(n) is used as the input of the ALE module. The ALE module will use past y(n) to approximate the current y(n). When y(n) is unpredictable, the ALE output power e(n) is approximately zero; conversely, when y(n) is unpredictable, the ALE output power is approximately zero. At this time, the ALE can predict y(n), therefore its output power e(n) increases; using the output power of the ALE as... The standard for measuring power is H(z), which is a U-order FIR filter with tap weights h(n).

[0071] The output e(n) of the ALE module can be expressed as:

[0072]

[0073] The 2-norm of e(n) is used as the cost function, and the update expressions for each parameter to be updated are derived based on this function.

[0074] The delay module parameter D and the impulse response of the filter function H(z) in the ALE section need to be initialized. The impulse response can be a zero sequence, and the U sampling points of the signal should cover at least one period of periodic interference, i.e., f r The reciprocal of , where D is the additional delay.

[0075] Step 4 includes the following steps:

[0076] Step 4.1: Update the notch bandwidth control parameter b. Use the SPRE algorithm to obtain the update expressions for each structural parameter of the adaptive comb filter. For b, use the 2-norm |e(n)| of the signal e(n). 2 As the cost function, the gradient with respect to b is:

[0077]

[0078] Where u(n) = x(n) - bu(nN), Filter[X(z)] represents the signal multiplied by it being passed through a filter with transfer function X(z), N is the order of the comb filter, and the superscript H indicates taking the conjugate transpose;

[0079] Considering system stability, b needs to be within the range (-1, 1), and the final update formula for b is:

[0080]

[0081]

[0082] Where the subscript T represents the iteration number, μ b Let b be the learning rate, and β be the upper limit of b, which can generally be taken as a value slightly less than 1.

[0083] Step 4.2: Update the tap coefficients g(n) of the notch filter gain control module.

[0084] Similarly, |e(n)| 2 Taking the gradient of g(n) yields

[0085]

[0086] P is the order of the gain function G(z), where

[0087] M u =[u(n-i+P)+u(n-N+P),u(n-i+P-1)+u(n-N+P-1),…,u(niP)+u(nNP)]#(4)

[0088] Ultimately there are

[0089]

[0090] Where μ g Let g(n) be the learning rate;

[0091] Step 4.3: Update the tap coefficients h(n) of the signal periodic component extraction module.

[0092] The ALE filter function H(z) also needs to be updated to obtain the optimal y(n) approximation tap coefficients, |e(n)| 2 Taking the gradient of h(n) yields

[0093]

[0094] in

[0095] M y =[y(nD),y(nD-1),…,y(nD-(U-1))]#(7)

[0096] Here, D and U are the additional delay of the ALE module and the order of the filter function H(z), respectively, ultimately resulting in:

[0097]

[0098] Where μ h Let h(n) be the learning rate;

[0099] Step 4.4: Pass the demodulated vital sign signal waveform through an adaptive comb filter and iterate according to the formulas in steps 4.1-4.3 until the parameters to be optimized, b, g(n), and h(n), are basically converged. The standard can be set to stop the iteration when the change of the three parameters is less than 0.05% after each iteration.

[0100] Step 5 includes the following steps:

[0101] Step 5.1: Perform FFT processing on the output signal after passing through the basically converged adaptive comb filter in the previous step to obtain its corresponding spectrum;

[0102] Step 5.2: Select the heartbeat frequency from the remaining spectral peaks.

[0103] Beneficial effects: The heart rate extraction method based on adaptive comb filter of the present invention has the following advantages:

[0104] The heart rate extraction method based on an adaptive comb filter proposed in this invention fully considers the intermodulation characteristics of heart rate and respiratory rate in Doppler radar systems. It utilizes adaptive comb filter technology to adaptively process orthogonally sampled complex signals, effectively reducing the interference of respiratory and its higher harmonic components on heart rate selection. Experimental results show that when the amplitude and position of respiratory higher harmonics in the detected vital signs affect heart rate measurement, the adaptive comb filter proposed in this invention for removing respiratory harmonics can measure the heart rate of the subject more accurately than existing methods. Attached Figure Description

[0105] Figure 1 is a flowchart of a heart rate extraction method based on an adaptive comb filter according to the present invention;

[0106] Figure 2 shows the spectrum of measured complex data in the embodiment;

[0107] Figure 3 shows a block diagram of the adaptive comb filter system in the embodiment;

[0108] Figure 4 shows the variation curves of the amplitude-frequency response of the traditional comb filter in the embodiment under different parameters b.

[0109] Figure 5. Block diagram of the adaptive spectral line enhancer;

[0110] Figure 6 shows the amplitude-frequency response of the comb filter after initializing each parameter in the embodiment.

[0111] Figure 7 shows the sum error norm variation curves of the three parameters to be updated after basic convergence in the embodiment;

[0112] Figure 8 shows the amplitude-frequency response of the comb filter after basic convergence in the embodiment;

[0113] Figure 9 is a comparison of the spectra of the input signal, the output signal filtered by a traditional comb filter, and the output signal filtered by a comb filter after all parameters have basically converged in the embodiment. Detailed Implementation

[0114] To better understand the purpose, structure, and function of this invention, the following detailed description of a heart rate extraction method based on an adaptive comb filter, in conjunction with the accompanying drawings, is provided.

[0115] Example 1:

[0116] This method addresses the technical challenge of accurately extracting vital signs using existing processing techniques when the detected spectral heartbeat component is weak or affected by higher harmonics of respiration. To verify the effectiveness of this method, this implementation uses the following scenario as an example to measure the vital signs of a subject: radar frequency 77 GHz, subject's reference respiratory rate 0.3 Hz, and heart rate 1.3 Hz. The digital sampling frequency is 20 Hz, and the input signal-to-noise ratio is 5 dB.

[0117] This example discloses a heart rate extraction method based on an adaptive comb filter. The vital signs include respiratory rate and heart rate features. The specific implementation block diagram is shown in Figure 1, and includes the following steps:

[0118] Step 1: Acquire the baseband signal of the quadrature sampling received by the Doppler radar, process the baseband signal to finally obtain the life signal containing breathing and heartbeat information;

[0119] Step 1.1: When the Doppler radar transmitting antenna emits a single-frequency continuous wave electromagnetic wave signal that illuminates a human target, the echo signal returning to the receiving antenna can be written as:

[0120] x R (t)=σcos(2πf c (t-τ))

[0121] Where A T f is the signal amplitude. c Let σ be the carrier frequency, t be the time, σ be the scattering coefficient of the target, τ = 2d(t) / c be the echo delay of the target, d(t) be the distance of the radar to the target, and c be the electromagnetic wave propagation speed.

[0122] The Doppler radar hardware system uses the transmitted signal to mix the echo signal and then performs low-pass filtering to obtain the mixed signal, which can be written as:

[0123]

[0124] The mixed signal is sampled quadraturely. Considering the DC bias introduced by the inconsistency in the IQ sampling channel levels in the actual hardware system, the baseband signals I(t) and Q(t) obtained by the IQ quadrature sampling channels are as follows:

[0125]

[0126]

[0127] DC I DC QThese are the DC components of the in-phase and quadrature components in the baseband signal, respectively. A I A Q Let λ = c / f be the amplitudes of the in-phase and quadrature components in the baseband signal, respectively. c Indicates the carrier frequency wavelength.

[0128] Assuming the sampling frequency f of the Doppler radar system s If the observation time is T, then the discrete sampled signal of the IQ channel is obtained as follows:

[0129]

[0130]

[0131] Where n = 1, ..., N, N = f s T represents the number of signal sampling points.

[0132] Doppler radar extracts respiratory and heart rate characteristics by sensing the relatively weak periodic distance changes in the chest wall of a human target. The distance d(n) can be expressed as:

[0133]

[0134] Where d0 is the nominal distance of the radar to the target, H i and R j f represents the amplitude of the periodic displacement of the chest wall caused by a certain harmonic of the target heartbeat and respiration, respectively. h and f r These represent the target heart rate and respiratory rate, respectively.

[0135] Step 1.2: Estimate the mean of the baseband signal and

[0136] The size of the estimated baseband signal mean processing block is N.

[0137] The method for processing the estimated baseband signal mean is as follows:

[0138]

[0139]

[0140] Step 1.3: Use mean filtering to remove the DC component of the baseband signal:

[0141] The mean filtering method is as follows:

[0142]

[0143]

[0144] Considering that the amplitude distortion caused by the inconsistency between the I and Q channels is usually small, we can assume that the signal amplitudes of the I and Q channels are basically the same after mean filtering. The above formula can be rewritten as:

[0145]

[0146]

[0147] Where A is the signal amplitude after mean filtering.

[0148] Step 1.4: Using the mean-filtered signal as the observable, construct the complex signal C(n):

[0149]

[0150] Where j is an imaginary number. Since d0 is the nominal distance from the radar to the target, it remains essentially constant during the observation time and can therefore be considered a time-independent quantity. The above equation can be rewritten as:

[0151]

[0152] in J n (a) represents the nth-order Bessel function of the first kind, where a is an argument, and Div(x) represents all divisors of the integer x. In this formula, term (a) represents heart rate and its harmonic components, and term (b) represents respiratory rate and its harmonic components. According to this formula, the spectrum will contain not only peaks corresponding to respiratory rate and heart rate, but also numerous peaks of their higher harmonics. This may cause the true respiratory rate / heart rate peaks to be affected by these terms. Because the heartbeat has a relatively small amplitude and weak peaks, and the amplitude of the higher harmonic components of respiratory rate is difficult to determine and there is a possibility that a certain harmonic may deviate from the predetermined frequency, all of these factors make heart rate measurement more difficult.

[0153] Figure 2 shows the spectrum of a measured complex vital sign data. RR represents the respiratory rate component, and HR represents the heart rate component. nRR represents the nth harmonic component of respiration (corresponding to the component in term (b)), and RR-HR represents the cross term of respiration and heart rate (corresponding to the component where p = -1 and l = 1 in the above formula). As can be seen from the figure, the spectrum of vital signs in actual situations is quite complex. For example, the usually weaker third harmonic of respiration is stronger than the second harmonic, and all harmonic components of respiration exhibit diffuse diffusion. This is because the target chest wall movement of the human body includes complex higher harmonic components and the gradual changes in respiratory and heartbeat movements themselves. In this situation, the true heart rate is affected by nearby higher harmonics of respiration, and since the higher harmonics of respiration themselves deviate somewhat from the predetermined frequency, accurately extracting the subject's vital sign information becomes increasingly difficult.

[0154] Traditional comb filters can filter out periodic interference, but when applied to complex data, they perform poorly in cases of frequency shift in higher harmonics, failing to completely filter out diffuse respiratory harmonics, and residual components still affect heart rate measurement. Therefore, designing a comb filter that can automatically adjust its parameters based on the intensity and location of each harmonic would make heart rate extraction more feasible after the higher harmonics of the adjacent respiratory rate have faded. Step 2: Based on the above characteristics of vital signals, design an adaptive comb filter with adjustable notch width and gain, and initialize the parameters of the comb filter according to the signal conditions;

[0155] Step 2.1: Initialize the filter order parameter N. The system function of the adaptive comb filter proposed in this patent is:

[0156]

[0157] The system block diagram is shown in Figure 3. The red dashed line represents the main part of the comb filter. The blue dashed line represents the ALE module, which will be introduced in the next section. In the system functions, N represents the order of the comb filter. For a comb filter of order N, its m-th notch frequency is...

[0158]

[0159] It is evident that the first notch frequency, w1, is entirely determined by N. If w1 is the respiratory rate, the comb filter can filter out respiratory frequencies and their integer multiples of harmonic frequencies. However, doing so directly may also affect heart rate components that are close to certain integer multiples of respiratory harmonic frequencies. Performing FFT processing on the complex data obtained in the previous step reveals that the respiratory rate and its higher harmonics are more significant in the spectrum. The respiratory rate, denoted as f, is selected from frequencies falling within [0.15Hz, 0.5Hz] and exhibiting a significant harmonic effect.r Then N is initialized as:

[0160] N = f s / f r

[0161] Where f s The sampling frequency.

[0162] Step 2.2: Initialize the filter gain function G(z).

[0163] The gain function G(z) ensures that the comb filter has appropriate fading at different frequencies; it is an FIR filter of order 2P+1. A lower order reduces the computational cost of subsequent iterations, but accuracy will suffer. Considering system stability, P ≤ N must be satisfied. The gain of a traditional comb filter is a constant. Replacing it with a filter G(z) with tap weights g(n), the amplitude-frequency response of the final comb filter at different frequencies can be adjusted by changing g(n). We want the final amplitude-frequency response of the comb filter to have fading suitable for the actual situation of complex data at each breathing frequency harmonic point. The initialization of G(z) can be accomplished using the window function design method. Let the height of the harmonic point w at a certain breathing frequency of the complex data be |W(e^(-1 / 2))|. jw We hope that it will fade to |S(e)| after passing through the filter. jw )|(Here S(e jw (|can be set by yourself), then the filter's amplitude frequency response at w is

[0164]

[0165] Then the i-th (-P≤i≤P) element g(n) of the tap coefficients of G(z) i for

[0166]

[0167] The relationship between G(z) and g(n) is as follows:

[0168]

[0169] Step 2.3: Initialize parameter b.

[0170] The parameter b controls the notch width, i.e., the width of each "comb tooth". Its value range is [0,1]. Figure 4 shows the system amplitude-frequency response when the gain function G(z) = 0 and N = 5, and b is 0, 0.4, and 0.8 respectively. It can be seen that as b approaches 1, the notch band becomes narrower. For practical data, assigning an initial value to b within [0.2,0.6] is more reasonable. Figure 6 shows the amplitude-frequency response of the comb filter after initializing the parameters in this example.

[0171] Step 3: Design the structure of the adaptive line enhancer (ALE) for adaptively adjusting the parameters of the comb filter, and initialize the relevant parameters of the ALE;

[0172] Step 3.1: Design the adaptive spectral line enhancer module.

[0173] The adaptive line enhancer (ALE) is used to evaluate the filtering effect of the previous comb filter module and obtain the update expressions for several parameters to be updated based on its output e(n). Its system block diagram is shown in Figure 5.

[0174] We express the input signal x(n) of the comb filter as follows:

[0175] x(n) = d(n) + i(n)

[0176] Where d(n) is the desired useful signal and i(n) is the periodic interference signal. After passing through the previous comb filter module, the filtered output y(n) can be expressed as:

[0177]

[0178] in This is the residue of a periodic interference signal. Using y(n) as input to the ALE module, the ALE module will approximate the current y(n) using past y(n). When When y(n) is unpredictable, ALE cannot predict it. In this case, the power of ALE's output e(n) is approximately zero. Conversely, when y(n) is unpredictable... At this time, the ALE can predict y(n), therefore its output power e(n) increases. Therefore, we can use the output power of the ALE as... The standard for measuring power. H(z) is a U-order FIR filter with tap weights h(n).

[0179] The output e(n) of the ALE module can be expressed as:

[0180]

[0181] The 2-norm of e(n) can be used as the cost function, and the update expressions for each parameter to be updated can be derived based on this.

[0182] We need to initialize the delay module parameters D and the impulse response of the filter function H(z) in the ALE section. The impulse response can be a zero sequence, and the U sampling points of the signal should cover at least one period of periodic interference, i.e., f r The reciprocal of . D represents the additional delay.

[0183] Step 4: Obtain the update expressions for the notch bandwidth control parameter b, the tap coefficient g(n) of the notch gain control module, and the tap coefficient h(n) of the signal periodic component extraction module. Then, pass the demodulated vital sign signal waveform through an adaptive comb filter and iterate until the parameters to be optimized are basically converged.

[0184] Step 4.1: Update the notch bandwidth control parameter b. We use the SPRE algorithm to obtain the update expressions for each structural parameter of the adaptive comb filter. For b, let |e(n)| 2 As the cost function, the gradient with respect to b is:

[0185]

[0186] Where u(n) = x(n) - bu(nN), and Filter[X(z)] represents the signal multiplied by it that is passed through a filter with transfer function X(z).

[0187] Considering system stability, b needs to be within the range (-1, 1). The final update formula for b is:

[0188]

[0189]

[0190] Where the subscript T represents the iteration number, μ b Let b be the learning rate, and β be the upper limit of b, which can generally be taken as a value slightly less than 1, such as 0.99.

[0191] Step 4.2: Update the tap coefficients g(n) of the notch gain control module.

[0192] Similarly, |e(n)| 2 Taking the gradient of g(n) yields

[0193]

[0194] in

[0195] M u =[u(n-i+P)+u(n-N+P),u(n-i+P-1)+u(n-N+P-1),…,u(niP)+u(nNP)]#(4)

[0196] Ultimately there are

[0197]

[0198] Where μ g Let g(n) be the learning rate.

[0199] Step 4.3: Update the tap coefficient h(n) of the signal periodic component extraction module.

[0200] The ALE filter function H(z) also needs to be updated to obtain the optimal y(n) approximation tap coefficients. |e(n)| 2 Taking the gradient of h(n) yields

[0201]

[0202] in

[0203] M y =[y(nD),y(nD-1),…,y(nD-(U-1))]#(7)

[0204] Ultimately there are

[0205]

[0206] Where μ h Let h(n) be the learning rate.

[0207] Step 4.4: Pass the demodulated vital sign signal waveform through an adaptive comb filter, and iterate according to steps 4.1-4.3 until the parameters to be optimized, b, g(n), and h(n), are basically converged. The criterion can be set to stop iterating when the change of the three parameters is less than 0.05% after each iteration.

[0208] As the iteration progresses, the notch bandwidth and gain at each notch point of the adaptive filter are dynamically adjusted according to the data conditions. The learning rate of the three parameters should be flexibly set according to factors such as signal amplitude, and should not be too large. Generally speaking, a larger learning rate will result in faster convergence speed for the corresponding parameters, and greater fluctuations in the parameters near the convergence value. Conversely, a larger learning rate will result in slower convergence speed and smaller fluctuations near the convergence value, but can avoid situations where oscillations are difficult to achieve convergence. The output power of the ALE module is |e(n)|. 2The parameters should exhibit a convergence trend, similar to those to be updated, and become smaller over time. Figure 7 shows an example of how b, g(n), and h(n) converge to near-convergence with increasing iteration count. Figure 8 shows the amplitude-frequency gain of the current gain filter G(z) and the amplitude-frequency gain of the entire adaptive comb filter.

[0209] Step 5: Select the heart rate from the spectrum of the output signal after the vital signs signal waveform has passed through the adaptive comb filter.

[0210] Step 5.1: Perform FFT processing on the output signal after passing through the basically converged adaptive comb filter in the previous step to obtain its corresponding spectrum;

[0211] Figure 9 shows a comparison of the spectra of the input signal in this example, the output signal filtered by a traditional comb filter, and the output signal filtered by an adaptive comb filter after the parameters proposed in this patent have basically converged. In the figure, RR represents respiratory rate, and HR represents heart rate. It can be seen that the respiratory harmonics are more effectively suppressed, especially the fourth and fifth harmonics near the heart rate component.

[0212] Step 5.2: Select the heartbeat frequency from the remaining spectral peaks.

[0213] In this example, within the range of [0.75, 2.0] Hz, there are two significant peaks at 1.3 Hz and 1.0 Hz. 1.0 Hz is the cross term of the true heart rate of 1.3 Hz and the respiratory rate of 0.3 Hz.

[0214] It is understood that the present invention has been described through some embodiments, and those skilled in the art will recognize that various changes or equivalent substitutions can be made to these features and embodiments without departing from the spirit and scope of the invention. Furthermore, under the teachings of the present invention, these features and embodiments can be modified to adapt to specific situations and materials without departing from the spirit and scope of the invention. Therefore, the present invention is not limited to the specific embodiments disclosed herein, and all embodiments falling within the scope of the claims of this application are within the protection scope of the present invention.

Claims

1. A heart rate extraction method based on an adaptive comb filter, characterized in that, Heart rate extraction involves separating the heart rate frequency component and the high-order harmonic components of the respiratory frequency, which can cause interference, from vital signs. The specific steps include: Step 1: Acquiring the baseband signal from orthogonal sampling received by Doppler radar, and processing the baseband signal to obtain a vital signal containing respiratory and heart rate information; Step 2: Using an adaptive comb filter with adjustable notch width and gain based on the characteristics of the vital signal, and initializing the parameters in the comb filter according to the signal conditions; Step 3: Employing an adaptive line enhancer (ALE) structure for adaptively adjusting the parameters of the comb filter, and initializing the relevant ALE parameters; Step 4: Obtaining the notch bandwidth control parameter b, the tap coefficient g(n) of the notch gain control module, and the tap coefficient of the signal periodic component extraction module. The update expression of the number h(n) is used to iterate the demodulated vital sign signal waveform through an adaptive comb filter until the parameters to be optimized are basically converged; Step 5: Select the heartbeat frequency from the spectrum of the output signal after passing through the adaptive comb filter; Wherein: Step 3 includes the following steps: Step 3.1: Design an adaptive spectral enhancer module. The adaptive spectral enhancer (ALE) is used to evaluate the filtering effect of the previous comb filter module and obtain the update expressions of several parameters to be updated based on its output e(n). The input signal x(n) of the previous comb filter is expressed as: x(n) = d(n) + i(n) where d(n) is the desired useful signal and i(n) is the periodic interference signal. After passing through the previous comb filter module, the filtered output y(n) is expressed as: ;in As the residual of the periodic interference signal, y(n) is used as the input of the ALE module. The ALE module will use past y(n) to approximate the current y(n). When y(n) is unpredictable, the ALE output power e(n) is approximately zero; conversely, when y(n) is unpredictable, the ALE output power is approximately zero. At this time, the ALE can predict y(n), therefore its output power e(n) increases; using the output power of the ALE as... The standard for measuring power is H(z), which is a U-order FIR filter with tap weights h(n); the output e(n) of the ALE module can be expressed as: The 2-norm of e(n) can be used as the cost function, and the update expressions for each parameter to be updated can be derived based on this. Here, it is necessary to initialize the delay module parameter D and the impulse response of the filter function H(z) in the ALE part. The impulse response can be taken as a 0 sequence, and the U sampling points of the signal should cover at least one periodic interference period, i.e., f r The reciprocal of , where D is the additional delay.

2. The heart rate extraction method based on an adaptive comb filter according to claim 1, characterized in that, Step 1 specifically includes the following steps: Step 1.1, the Doppler radar transmitting antenna transmits a single-frequency continuous wave signal electromagnetic wave to illuminate the human target, which returns to the receiving antenna. The echo signal reaching the receiving antenna is written as: The Doppler radar hardware system uses the transmitted signal to mix the echo signal and then performs low-pass filtering to obtain the mixed signal, which can be written as: Among them: A T f is the signal amplitude. c Let t be the carrier frequency and t be the time. Let τ = 2d(t) / c be the scattering coefficient of the target, d(t) be the distance the radar reaches the target, and c be the electromagnetic wave propagation speed. The mixed signal is then subjected to quadrature sampling. Considering the quadrature sampling in the actual hardware system, i.e., the DC bias introduced by the inconsistency in the IQ sampling channel levels, the baseband signals I(t) and Q(t) obtained by the IQ quadrature sampling channels are respectively: ； DC I DC Q These are the DC components of the in-phase and quadrature components in the baseband signal, respectively. A I A Q Let λ = c / f represent the amplitudes of the in-phase and quadrature components in the baseband signal, respectively. c Indicates the carrier frequency wavelength; if the sampling frequency f of the Doppler radar system... s If the observation time is T, then the discrete sampled signal of the IQ channel is obtained as follows: ； Where n=1,…,N, N=f s T represents the number of signal sampling points; Doppler radar extracts respiratory and heart rate characteristics by sensing the weak periodic distance changes in the chest wall of a human target, and its distance d(n) is represented as: Where d0 is the nominal distance of the radar to the target, and H i and R j f represents the amplitude of the periodic displacement of the chest wall caused by a certain harmonic of the target heartbeat and respiration, respectively. h and f r These represent the target heart rate and respiratory rate, respectively; Step 1.2: Estimate the mean of the baseband signal. and The size of the estimated baseband signal mean processing block is N, and the estimated baseband signal mean processing method is as follows: ； Step 1.3: Remove the DC component of the baseband signal using mean filtering: The mean filtering method is as follows: ； Considering that amplitude distortion caused by inconsistency between the I and Q channels is usually small, assuming that the signal amplitudes of the I and Q channels are basically the same after mean filtering, the above formula can be rewritten as: ； Where A is the signal amplitude after mean filtering; Step 1.4: Using the mean-filtered signal as the observation, construct the complex signal C(n): ; where j is an imaginary number.

3. The heart rate extraction method based on an adaptive comb filter according to claim 1, characterized in that... Step 2 includes the following steps: Step 2.1, Initialize the filter order parameter N. The system function of the adaptive comb filter is: In the above formula, N represents the order of the comb filter, b represents the notch bandwidth control parameter of the filter, z is the z-transform operator, G(z) is the gain control function of the adaptive comb filter, and P is its order. For a comb filter of order N, its m-th notch frequency is... It can be seen that its first notch frequency w1 is completely determined by N. If w1 is the respiratory frequency, the comb filter can filter out respiratory and respiratory integer multiple harmonic frequency components. However, doing so directly may also affect the heart rate component that is close to a certain integer multiple harmonic frequency of breathing. Performing FFT processing on the complex data obtained in the previous step, the respiratory frequency and its higher harmonics will be more significant in the spectrum. Select the frequency that falls within [0.15Hz, 0.5Hz] and has a significant harmonic effect as the respiratory frequency, denoted as f. r Then N is initialized as: N = f s / f r Where f s The sampling frequency is given by step 2.

2. Initialize the filter gain function G(z). The gain function G(z) ensures that the comb filter has appropriate fading at different frequencies. It is an FIR filter of order 2P+1. A lower order can reduce the computational load of subsequent iterations, but the accuracy will be lacking. Considering the system stability, P≤N must be satisfied. The gain of the traditional comb filter is a constant. Replacing it with a filter G(z) with tap weights g(n) allows the amplitude-frequency response of the final comb filter at different frequencies to be adjusted by adjusting g(n). It is hoped that the amplitude-frequency response of the final comb filter has fading suitable for the actual situation of complex data at each breathing frequency harmonic point. The initialization of G(z) can be accomplished by the window function design method. Let the height of the harmonic point w of a certain breathing frequency of complex data be given by step 2.

2. Hopefully, it will fade to [value] after passing through the filter. Then the filter's amplitude-frequency response at w is Then the i-th element of the tap coefficient g(n) of G(z) for Where -P≤i≤P; the relationship between G(z) and g(n) is: Step 2.3: Initialize the notch bandwidth control parameter b. Parameter b controls the notch width, that is, the width of each "comb tooth". Its value range is [0,1]. As b approaches 1, the notch band will become narrower. For actual data, b can be assigned an initial value within [0.2,0.6].

4. The heart rate extraction method based on an adaptive comb filter according to claim 1, characterized in that... Step 4 includes the following steps: Step 4.1, update the notch bandwidth control parameter b, and use the SPRE algorithm to obtain the update expressions for each structural parameter of the adaptive comb filter. For b, the signal... 2-norm As the cost function, the gradient with respect to b is: Where u(n) = x(n) - bu(nN), Filter[X(z)] represents the signal multiplied by it passed through a filter with transfer function X(z), N is the order of the comb filter, and the superscript H indicates taking the conjugate transpose; considering system stability, b needs to be between (-1, 1), and the final update formula for b is Where the subscript T represents the iteration number, μ b Let β be the learning rate for parameter b, and β be the upper limit of b, which can generally be taken as a value slightly less than 1; Step 4.2, update the tap coefficients g(n) of the notch gain control module. Similarly, Taking the gradient of g(n) yields ; P is the order of the gain function G(z), where Ultimately there are ;where μ g Let g(n) be the learning rate; Step 4.3: Update the tap coefficients h(n) of the signal periodic component extraction module. The ALE filtering function H(z) also needs to be updated to obtain the optimal y(n) approximation tap coefficients. Taking the gradient of h(n) yields ;in Here, D and U are the additional delay of the ALE module and the order of the filter function H(z), respectively, ultimately resulting in: ;where μ h Let h(n) be the learning rate; Step 4.4: Pass the demodulated vital sign signal waveform through an adaptive comb filter and iterate according to the formulas in Steps 4.1-4.3 until the parameters to be optimized, b, g(n), and h(n), are basically converged; The standard can be set to stop the iteration when the change of the three parameters is less than 0.05% after each iteration.

5. The heart rate extraction method based on an adaptive comb filter according to claim 1, characterized in that... Step 5 includes the following steps: Step 5.1, perform FFT processing on the output signal after passing through the basically converged adaptive comb filter in the previous step to obtain its corresponding spectrum; Step 5.2, select the heartbeat frequency from the remaining spectral peaks.

Images