CN110276271A - Non-contact Heart Rate Estimation Method Fusion IPPG and Depth Information Anti-Noise Interference - Google Patents

Abstract

The invention relates to a non-contact heart rate estimation method integrating IPPG and depth information to resist noise interference, which includes the following steps: 1) collecting the face image of the subject and its corresponding depth information through an RGBD camera; 2) adopting a face detection algorithm Use template matching to detect and locate the collected face images to form a surveillance area; 3) Use face tracking algorithm to track the surveillance area corresponding to multiple continuously collected face images to capture the movement of the face, and obtain more 4) Use the weighted mean method to eliminate the motion noise interference of the obtained multiple face images with motion changes, and realize dynamic face recognition; 5) The face image after face recognition The depth information is used to adaptively select the ROI area; 6) The source signal is extracted from the ROI area by IPPG technology, the signal denoising is performed by the empirical decomposition mode EMD method, and the ARMA model is used for energy spectrum analysis to obtain the heart rate of the detected person.

Description

Non-contact Heart Rate Estimation Method Fusion IPPG and Depth Information Anti-Noise Interference

technical field

The invention relates to the technical field of heart rate estimation, in particular to a non-contact heart rate estimation method that combines image photoplethysmography (IPPG) and depth information to resist motion noise interference.

Background technique

According to the World Health Organization, more than 17.5 million people die each year due to cardiovascular disease, which is the leading cause of non-communicable disease mortality. Heart rate is an important physiological parameter that can reflect the state of the heart, and is an important indicator for cardiovascular disease prevention and clinical diagnosis.

The traditional method of heart rate measurement, the electrocardiogram is the standard method for heart rate measurement in clinical practice. It requires medical staff to place the lead electrodes according to the correct lead positions. The process is cumbersome and requires professional operation. During use, it is necessary to fill conductive media or apply pressure between the electrodes and the skin to fix it, so it brings a great deal to the subject. discomfort, and the complicated connecting lines are inconvenient to operate. The existing heart rate detection devices and related signal processing technologies mainly include the following solutions:

(1) Patch-type heart rate detector: paste electrode patches on some specific positions on the subject's body to extract the heart rate from the human body signal. However, the patch-type heart rate detector will cause some discomfort and stimulation to the subject when the electrode patch is attached during detection, and the test environment is limited, so it cannot be detected in real time for a long time, and the process is complicated and inconvenient to operate.

(2) Capacitive heart rate measurement: Based on the capacitive coupling principle, non-contact heart rate detection is completed; the capacitive non-contact ECG sensor has good ECG measurement performance and high reliability, and can be used for long-term measurement of ECG signals. However, the heart rate measurement result of the capacitive heart rate measurement is easily affected by motion noise, so the effect obtained by this method is not ideal.

(3) Photoelectric pulse measurement: Using the high opacity of blood and the characteristics that the penetration of light in general tissues is dozens of times larger than that in blood, the pulse signal is checked. However, whether the photoelectric pulse measurement adopts the reflection or projection method, in addition to being easily affected by motion noise, it also requires the sensor to be close to the skin in order to reduce the interference of ambient light noise.

(4) Heart rate measurement by photoplethysmography: The heart rate is calculated by recording the amount of light absorbed by the skin tissue after the light passes through the skin and tissue. Since this technology does not directly contact the measured object, it is a non-contact measurement. Has a better human-machine experience. However, the photoplethysmography method for heart rate detection is limited by the ambient light. Once the light changes or is not strong enough, the heart rate test will become inaccurate; and the tester needs to sit in a fixed position, and the change of position will also make the detection result. Inaccurate; in the moving state, the ROI area selection will also be inaccurate; and it is a technical difficulty to eliminate the noise in the signal due to the relative movement between the measured object and the image acquisition device.

(5) Image pulse measurement: The deformation of blood vessels according to the heartbeat can be observed by the resulting deformation of the skin around the blood vessels. However, this image pulse measurement is susceptible to artifacts and vibrations perturbing the measurement data and to the experimental light conditions.

(6) For signal processing, the current main processing technologies are: finite impulse response (finite impulse response, FIR) filter, band-pass filter and moving average filter.

In recent years, with the development of electronic technology and measurement methods. Non-contact heart rate measurement is proposed and implemented. Photoplethysmography is a non-invasive detection method that uses photoelectric sensors to detect changes in blood volume in living tissue. Due to the different absorption of light by different tissues of the human body, the received light intensity will change regularly after transmission or reflection. The regular contraction and relaxation of the heart is the main cause of the regular changes, so these pulsatile signals are directly related to the beating of the heart. However, the above measurement method has the problem that the collected signal is easily disturbed by motion noise, resulting in low accuracy.

Image photoplenthysmography (IPPG) is a non-contact physiological signal measurement technology developed on the basis of traditional photoplethysmography in recent years. The characteristics of non-contact measurement, low cost, and easy operation of IPPG technology, especially the non-contact measurement method, enable it to realize clinical and daily detection in some specific cases. The proposed image photoplethysmography (IPPG) generally uses a fixed ROI area and requires the subject to sit still in front of the image acquisition device without considering the subject's movement. The fixed ROI area set in the method cannot overcome the problem of inaccurate estimated heart rate signal caused by the noise introduced by the subject's movement.

SUMMARY OF THE INVENTION

The purpose of the present invention is to provide a non-contact heart rate estimation method integrating image photoplethysmography (IPPG) and depth information with reasonable design, high accuracy, and anti-motion noise interference, aiming at the deficiencies of the prior art.

To achieve the above object, the present invention adopts the following technical solutions:

A non-contact heart rate estimation method that combines IPPG and depth information to resist noise interference, which includes the following steps:

1) Information collection: collect the subject's face image and its corresponding depth information through the RGBD camera;

2) Face image detection: the face detection algorithm is used to detect and locate the collected face image using template matching to form a monitoring area;

3) face tracking: adopt the face tracking algorithm to track and capture the movement of the face in the monitoring area corresponding to the multiple face images collected continuously, and obtain a plurality of face images with motion changes;

4) Elimination of motion noise interference: use the weighted mean method to eliminate the motion noise interference of a plurality of acquired face images with motion changes, and realize dynamic face recognition;

5) ROI area extraction: adaptively select the ROI area by using the depth information for the face image after face recognition;

6) Heart rate calculation using IPPG technology: using IPPG technology to extract the source signal in the ROI area, using the empirical decomposition mode EMD method to denoise the signal, and using the ARMA model to perform energy spectrum analysis to obtain the heart rate value of the subject.

Preferably, the RGBD camera described in step 1) includes an RGB camera and a depth camera, the RGB camera is used to capture the RGB image of the human body, and the depth camera is used to capture the depth image of the human face. After the RGB camera and the depth camera are calibrated , the depth image captured by the depth camera and the RGB image captured by the RGB camera overlap.

Preferably, the method of step 2) face image detection is: taking the collected face image as an input image, searching for a rectangular area of possible scales and positions in the input image, that is, a candidate window to be detected, for each candidate window according to The following steps are used for processing: first use binocular template matching for rough screening, and then standardize the mean square error of the images in the window to standardize the distribution area of the face area, eliminate the influence of illumination changes, and then perform face template matching. The threshold of , then it is output as a candidate face.

Preferably, step 3) the method of face tracking is: assuming that only one face appears in each frame in the sequence, set the initial monitoring area as R ₀ , first detect the face in the monitoring area of each frame, if detected face, according to its positioning result, calculate the new monitoring area in the next frame according to the formula for solving the monitoring area, and process each frame in this way; if the detection fails, keep the monitoring area unchanged, and continue to detect in subsequent frames to prevent Tracking failure due to accidental missed detection; when the number of missed frames exceeds a certain number, that is, the missed detection times out, and the face is considered to have disappeared, it will restart to detect new faces that may appear in subsequent frames in the initial monitoring area;

Suppose, the surveillance area is represented by a six-tuple as LR=(x _min , x _max , y _min , y _max , ch _min , ch _max ), x _min , x _max are the range of the center position of the face in the x direction, y _min , y _max is the range of the center position of the face in the y direction, ch _min , ch _max is the range of the face scale; the face area is approximated as a square, which is represented by a triple as L=(x _center , x _center ,ch _face ), each element describes the center position and size of the face in turn; then, the formula for solving the monitoring area is as follows:

Among them, a and b represent the x and y positions of the face between two frames, respectively, and c represents the maximum change of the face scale.

Preferably, step 4) the method for eliminating motion noise interference is: according to the change of the left and right swing angle of the human face, define the weight calculation method of each face change, propose a construction strategy of the weighted average face, and then based on the pitch of the face The angle changes, it is divided into three levels of looking down, looking up and looking up, and constructing a weighted mean face in each level to form a weighted mean face matrix;

Assuming that the given motion change grayscale face image is I _j (x, y), the formula of the weighted mean face is as follows:

Among them, ZM represents the total number of different moving images of the same person, j represents the j-th changing image, x, y represent the two-dimensional image plane coordinates, and ω _j represents the weight value corresponding to the j-th face image;

The determination of ω _j is related to the left and right swing angle of the face movement. The calculation method of ω _j is as follows: Assuming that the coordinates of the eyes before the face swings left and right are A(x _A , y _A ) and B (x _B , y _B ), the human The coordinates of the eyes after the face is swayed to the right are A'(x _A' ,y _A' ) and B'(x _B' ,y _B' ) respectively, and the coordinates of the center of the circle where the cross-sectional view of the face is located is O(x _o , y _o ), the swing angle of the face to the right is ∠δ, then the calculation formula of the weight value ω _K of the K-th face image is as follows:

Among them, δ _k represents the left and right swing angle corresponding to the kth face image, then the solution of the weight ω _k is converted into the calculation of the swing angle δ _k ;

Eye coordinates A(x _A , y _A ), B(x _B , y _B ) and their eye coordinates A'(x _A' , y _A' ), B'(x _B' , y _B' ) after rocking left and right The positioning has been determined; assuming that the radius of the circle where the eyes rotate is R, then

The coordinates of the center point O(x _o , y _o ) are obtained through the following calculation process:

Combining the above two formulas (6) and (7), we can get:

The same can be obtained:

The coordinates of the center point o(x _o , y _o ) can be obtained by combining the two equations (8) and (9) above, and the result is brought into the equation (6) to obtain the radius R, then the solution can be obtained by the equation (5). The rocking angle δ; the rocking angle δ is a function of the coordinate value of the eyes, namely

δ=f(x _A ,x _B ,x _A' ,x _B' ) (10)

Assuming that the left and right sway angle δ, the pitch angle is γ, the relationship table between the division of the movement range of the face and its corresponding weighted mean face is obtained; given the gray-scale face image I _i (x, y) of the movement change, then The weighted mean face with a left-right swing angle of 0°-90° is defined as shown in formula (11):

Among them, ZM represents the total number of faces within a certain pitch angle and a certain left and right swing angle, j represents the jth face, x and y represent the two-dimensional image plane coordinates, and p and q are the mean values corresponding to the variation range. Face number, F' _{p, q} represents the weighted mean face with left and right swing angles between -90°-0°;

The calculated weighted mean face finally forms the weighted mean face matrix when the face has various motion changes:

For the face image after face recognition, frame it with a rectangle, and use h to represent the distance from the depth camera to the surface of the object, N to represent the measured value of the number of pixels in the face area framed in the image, and S to represent the actual size of the object. ; M represents the number of pixels in the entire field of view of the depth camera, S ₁ represents the area of the entire field of view of the depth camera; α and β represent the horizontal and vertical fields of view of the depth camera, respectively. According to the ratio of the field of view range of the depth camera and the number of pixels in the field of view range of the depth camera to the area size S of the actual object and the number of pixels M in the entire field of view of the depth camera at this distance is equal to the ratio of the entire field of view area s ₁ , that is

Among them, N is the number of pixels in the area where the face is located according to the statistics of the image, and S ₁ can be calculated from the geometric relationship according to the field of view of the depth camera, as shown in the following formula:

And because of the total number of pixels M in the field of view S ₁ , using this principle, for the built image acquisition system, the face area area S can be obtained when S ₁ , M, and N are known:

For the determination of the area of the ROI area, the present invention chooses to draw parallel lines from the upper edge to the lower 45%-75% area of the face image after face recognition, and takes 30% of the overall height of the image as the height of the ROI. The width is taken as 70% of the width of the image after recognition, and the ROI area containing only the nose and left and right cheeks is obtained.

Preferably, step 6) adopts IPPG technology to calculate the heart rate as follows: extract the source signal from the ROI area, divide the ROI area image into three channels: R, G, and B, and use the G channel with the least noise as the channel of the source signal , perform spatial pixel averaging on the separated G channel image:

where k is the number of image frames, K is the total number of frames of the image, Z(k) is the one-dimensional source signal of the R channel, zi _,j (n) is the color intensity value of the pixel (i,j) in the G channel, g, d are the height and width of the image, respectively;

Using the empirical decomposition mode EMD method to denoise the signal, the steps are as follows:

(1) Take the average value of the original signal u(t) to obtain the signal m ₁ (t);

(2) Calculate the first-order margin hp ₁ (t)=u(t)-m ₁ (t), check whether hp ₁ (t) satisfies the IMF condition, if not, return to step (1) and use hp ₁ (t) as the raw signal for the second screening, namely:

hp ₂ (t)=hp ₁ (t)-m ₂ (t) (17)

Repeat screening k times:

hp _k (t)=hp _k-1 (t)-m _k (t) (18)

Before hp _k (t) satisfies the IMF condition, the first component IMF ₁ of the IMF is obtained, namely

IMF ₁ = hp _k (t) (19)

(3) Subtract IMF ₁ from the original signal u(t) to obtain the residual amount r ₁ (t), namely

r ₁ (t)=u(t)-IMF ₁ (20)

(4) Let u ₁ (t)=r ₁ (t), take u ₁ (t) as the new original signal, repeat the above steps to obtain the second IMF component IMF ₂ , repeat this n times:

(5) When the nth component r _n (t) has become a monotonic function and can no longer decompose the IMF, the entire EMD decomposition process is completed; the original signal u(t) can be expressed as n IMF components and an average trend component r The combination of _n (t), namely:

The ARMA model is used to analyze the energy spectrum of the signal with a frequency of 0.75-2.0Hz in the heartbeat frequency band after EMD decomposition, which corresponds to the normal human heart rate range of 45-120 beats/min, and the frequency corresponding to the highest energy point is the heartbeat frequency f _h , the heart rate is:

XLV= _60fh (22).

The present invention adopts the above technical solutions, in order to further improve the accuracy of image photoplethysmography (IPPG) heart rate estimation, the present invention proposes a non-contact heart rate estimation method that combines image photoplethysmography (IPPG) and depth information. The invention uses an image that can be measured in real time by an RGBD camera. By adding depth information and fusing it into the IPPG signal extraction, an adaptive ROI area extraction method is designed to reduce the interference of the IPPG source signal due to the introduction of noise due to inappropriate ROI area selection. At the same time, considering that the subject's head swing will cause the up and down pitch change and left and right swing change of the face image, the weighted average method is used to reduce the motion noise interference caused by the subject's head movement, and the image-based photoplethysmography ( IPPG) accuracy of heart rate estimation. Compared with the existing method, the heart rate acquisition method proposed by the present invention adopts an RGBD camera to acquire face image information and depth information non-contact, and does not need to attach electrodes during heart rate monitoring, which is more convenient to use, and is not limited by ambient light and motion. The influence of noise improves the accuracy of heart rate monitoring.

Description of drawings

Now in conjunction with accompanying drawing, the present invention is further elaborated:

Fig. 1 is the flow chart of the non-contact heart rate estimation method of fusion IPPG technology and depth information of the present invention;

2 is a schematic flowchart of a face detection algorithm of the present invention;

3 is a schematic flowchart of a face tracking algorithm of the present invention;

4 is a cross-sectional view of a face image of the present invention;

Fig. 5 is the relation table of the motion variation range division of the human face of the present invention and its corresponding weighted mean human face;

FIG. 6 is a view of the field of view of the depth camera of the present invention.

Detailed ways

As shown in one of Figures 1-6, the non-contact heart rate estimation method that integrates IPPG technology and depth information of the present invention, and the non-contact heart rate estimation method that integrates IPPG and depth information to resist noise interference, includes the following steps:

1) Information collection: collect the subject's face image and its corresponding depth information through the RGBD camera;

2) Face image detection: the face detection algorithm is used to detect and locate the collected face image using template matching to form a monitoring area;

3) face tracking: adopt the face tracking algorithm to track and capture the movement of the face in the monitoring area corresponding to the multiple face images collected continuously, and obtain a plurality of face images with motion changes;

4) Elimination of motion noise interference: use the weighted mean method to eliminate the motion noise interference of a plurality of acquired face images with motion changes, and realize dynamic face recognition;

5) ROI area extraction: adaptively select the ROI area by using the depth information for the face image after face recognition;

6) Heart rate calculation using IPPG technology: using IPPG technology to extract the source signal in the ROI area, using the empirical decomposition mode EMD method to denoise the signal, and using the ARMA model to perform energy spectrum analysis to obtain the heart rate value of the subject.

Preferably, the RGBD camera described in step 1) includes an RGB camera and a depth camera, the RGB camera is used to capture the RGB image of the human body, and the depth camera is used to capture the depth image of the human face. After the RGB camera and the depth camera are calibrated , the depth image captured by the depth camera and the RGB image captured by the RGB camera overlap.

As shown in Figure 2, step 2) the method of face image detection is: take the collected face image as the input image, search for the rectangular area of possible scale and position in the input image, that is, the candidate window to be detected, for each The candidate window is processed according to the following steps: first use binocular template matching for rough screening, and then standardize the mean square error of the images in the window to standardize the distribution area of the face area and eliminate the influence of illumination changes, and then perform face template matching. If it exceeds the set threshold, it will be output as a candidate face.

As shown in Figure 3, step 3) the method of face tracking is: assuming that only one face appears in each frame in the sequence, set the initial monitoring area as R ₀ , first detect the face in the monitoring area of each frame, If a human face is detected, the new monitoring area in the next frame is calculated according to the formula for solving the monitoring area according to the positioning result, and each frame is processed in this way; if the detection fails, the monitoring area is kept unchanged, and the detection is continued in subsequent frames. , to prevent the tracking failure caused by accidental missed detection; when the number of missed frames exceeds a certain number, that is, the missed detection times out, and the face is considered to have disappeared, it will restart to detect new people that may appear in subsequent frames in the initial monitoring area. Face;

The monitoring area refers to the possible position and scale range of a certain face in a certain frame. It is assumed that the monitoring area is represented by a six-tuple as LR=(x _min , x _max , y _min , y _max , ch _min , ch _max ), x _min , x _max is the range of the center position of the face in the x direction, y _min , y _max is the range of the center position of the face in the y direction, ch _min , ch _max is the range of the face scale; The area is approximately a square, which is represented by a triple as L=(x _center , x _center , ch _face ), and each element sequentially describes the center position and size of the face;

The task of face tracking is to track and capture the motion of a face in consecutive images. The tracking algorithm predicts the surveillance area corresponding to each face in a new frame according to the existing face localization results, and detects the face within its constraints. Due to the randomness of human movement (such as a sudden change of movement direction), it is difficult to predict a new surveillance area by using the existing human face movement trajectory without understanding the behavior of the person. Using a relatively simple method, only the tracking results of the previous frame are used, and the maximum change of the face between two frames is limited according to the actual application requirements, and the new monitoring area is further solved. In this way, although the processing speed will be affected to some extent, the robustness of the tracking can be guaranteed. The face tracking problem is aimed at the frontal face image sequence, and the movement of the face is mainly manifested in two forms: translation and scale change. Therefore, the formula for solving the surveillance area is as follows:

Among them, a and b represent the x and y positions of the face between two frames, respectively, and c represents the maximum change of the face scale, which is related to the movement speed of the face and the relative position relationship between the face and the camera.

The change of the face in three-dimensional space is divided into translation and rotation along the horizontal axis, vertical axis and number axis, and some of the translation and rotation can be effectively overcome by the method of set normalization. But for face images and up and down pitch changes and left and right roll changes, geometric normalization is often insurmountable. Aiming at the difficulty of face detection during motion, a face recognition algorithm based on weighted mean is adopted.

Preferably, step 4) the method for eliminating motion noise interference is: according to the change of the left and right swing angle of the human face, define the weight calculation method of each face change, propose a construction strategy of the weighted average face, and then based on the pitch of the face The angle change is divided into three levels: looking down, looking up and looking up, and constructing a weighted mean face in each level to form a weighted mean face matrix to realize dynamic face recognition; The face images are superimposed to form an average face. Such an image contains the comprehensive information of multiple moving faces and can reflect the changes of faces with different movements.

Assuming that the given motion change grayscale face image is I _j (x, y), the formula of weighted mean face (WMF) is as follows:

Among them, ZM represents the total number of different moving images of the same person, j represents the j-th changing image, x, y represent the two-dimensional image plane coordinates, and ω _j represents the weight value corresponding to the j-th face image;

As shown in Figure 4, the determination of ω _j is related to the left and right swing angle of the face movement. The calculation method of ω _j is: Assume that the coordinates of the eyes before the face swing is A(x _A , y _A ) and B(x _B respectively , y _B ), the coordinates of the eyes after the face swings to the right are A'(x _A' ,y _A' ) and B'(x _B' ,y _B' ), respectively, the center of the circle where the cross-sectional view of the face is located The coordinates are O(x _o , y _o ), and the swing angle of the face to the right is ∠δ, then the calculation formula of the weight value ω _K of the K-th face image is as follows:

Among them, δ _k represents the left and right swing angle corresponding to the kth face image, then the solution of the weight ω _k is converted into the calculation of the swing angle δ _k ;

Eye coordinates A(x _A , y _A ), B(x _B , y _B ) and their eye coordinates A'(x _A' , y _A' ), B'(x _B' , y _B' ) after rocking left and right The positioning has been determined; assuming that the radius of the circle where the eyes rotate is R, then

The coordinates of the center point O(x _o , y _o ) are obtained through the following calculation process:

Combining the above two formulas (6) and (7), we can get:

The same can be obtained:

The coordinates of the center point o(x _o , y _o ) can be obtained by combining the two equations (8) and (9) above, and the result is brought into the equation (6) to obtain the radius R, then the solution can be obtained by the equation (5). The rocking angle δ; the rocking angle δ is a function of the coordinate value of the eyes, namely

δ=f(x _A ,x _B ,x _A' ,x _B' ) (10)

Assuming the left and right swing angle δ, the pitch angle is γ, the motion variation range division of the face and the relation table of the corresponding weighted mean face are obtained, as shown in Figure 5; the grayscale face image I _i ( x, y), then the left and right swing angle is the weighted mean face definition within 0°-90° as shown in formula (11):

Among them, ZM represents the total number of faces within a certain pitch angle and a certain left and right swing angle, j represents the jth face, x and y represent the two-dimensional image plane coordinates, and p and q are the mean values corresponding to the variation range. Face number, F' _{p, q} represents the weighted mean face with left and right swing angles between -90°-0°;

The calculated weighted mean face finally forms the weighted mean face matrix when the face has various motion changes:

As shown in Figure 6, the face image after face recognition is framed by a rectangle, and h represents the distance from the depth camera to the surface of the object, N represents the measured value of the number of pixels in the face area framed in the image, S represents the actual area of the object; M represents the number of pixels in the entire field of view of the depth camera, S ₁ represents the area of the entire field of view of the depth camera; α and β represent the horizontal and vertical fields of view of the depth camera, respectively.

According to the ratio of the field of view range of the depth camera and the number of pixels in the field of view range of the depth camera to the area size S of the actual object and the number of pixels M in the entire field of view of the depth camera at this distance is equal to the ratio of the entire field of view area s ₁ , that is

Among them, N is the number of pixels in the area where the face is located according to the statistics of the image, and S ₁ can be calculated from the geometric relationship according to the field of view of the depth camera, as shown in the following formula:

And because of the total number of pixels M in the field of view S ₁ , using this principle, for the built image acquisition system, the face area area S can be obtained when S ₁ , M, and N are known:

For the determination of the area of the ROI area, the present invention chooses to draw parallel lines from the upper edge to the lower 45%-75% area of the face image after face recognition, and takes 30% of the overall height of the image as the height of the ROI. The width is taken as 70% of the width of the image after recognition, and the ROI area containing only the nose and left and right cheeks is obtained.

Preferably, step 6) adopts IPPG technology to calculate the heart rate as follows: extract the source signal from the ROI area, divide the ROI area image into three channels: R, G, and B, and use the G channel with the least noise as the channel of the source signal , perform spatial pixel averaging on the separated G channel image:

where k is the number of image frames, K is the total number of frames of the image, Z(k) is the one-dimensional source signal of the R channel, zi _,j (n) is the color intensity value of the pixel (i,j) in the G channel, g, d are the height and width of the image, respectively;

Using the empirical decomposition mode EMD method to denoise the signal, the steps are as follows:

(1) Take the average value of the original signal u(t) to obtain the signal m ₁ (t);

(2) Calculate the first-order margin hp ₁ (t)=u(t)-m ₁ (t), check whether hp ₁ (t) satisfies the IMF condition, if not, return to step (1) and use hp ₁ (t) as the raw signal for the second screening, namely:

hp ₂ (t)=hp ₁ (t)-m ₂ (t) (17)

Repeat screening k times:

hp _k (t)=hp _k-1 (t)-m _k (t) (18)

Before hp _k (t) satisfies the IMF condition, the first component IMF ₁ of the IMF is obtained, namely

IMF ₁ = hp _k (t) (19)

(3) Subtract IMF ₁ from the original signal u(t) to obtain the residual amount r ₁ (t), namely

r ₁ (t)=u(t)-IMF ₁ (20)

(4) Let u ₁ (t)=r ₁ (t), take u ₁ (t) as the new original signal, repeat the above steps to obtain the second IMF component IMF ₂ , repeat this n times:

(5) When the nth component r _n (t) has become a monotonic function and can no longer decompose the IMF, the entire EMD decomposition process is completed; the original signal u(t) can be expressed as n IMF components and an average trend component r The combination of _n (t), namely:

The ARMA model is used to analyze the energy spectrum of the signal with a frequency of 0.75-2.0Hz in the heartbeat frequency band after EMD decomposition, which corresponds to the normal human heart rate range of 45-120 beats/min, and the frequency corresponding to the highest energy point is the heartbeat frequency f _h , the heart rate is:

XLV= _60fh (22).

The present invention adopts the above technical solutions, in order to further improve the accuracy of image photoplethysmography (IPPG) heart rate estimation, the present invention proposes a non-contact heart rate estimation method that combines image photoplethysmography (IPPG) and depth information. The invention uses an image that can be measured in real time by an RGBD camera. By adding depth information and fusing it into the IPPG signal extraction, an adaptive ROI area extraction method is designed to reduce the interference of the IPPG source signal due to the introduction of noise due to inappropriate ROI area selection. At the same time, considering that the subject's head swing will cause the up and down pitch change and left and right swing change of the face image, the weighted average method is used to reduce the motion noise interference caused by the subject's head movement, and the image-based photoplethysmography ( IPPG) accuracy of heart rate estimation. Compared with the existing method, the heart rate acquisition method proposed by the present invention adopts an RGBD camera to acquire face image information and depth information non-contact, and does not need to attach electrodes during heart rate monitoring, which is more convenient to use, and is not limited by ambient light and motion. The influence of noise improves the accuracy of heart rate monitoring.

The above description should not limit the protection scope of the present invention in any way.

Claims (7)

1. merging the non-contact heart rate estimation technique of IPPG and depth information anti-noise jamming, it is characterised in that: it includes following step It is rapid:

1) facial image and its corresponding depth information of measured information collection: are acquired by RGBD camera；

2) facial image detects: carrying out Face datection using template matching to collected facial image using Face datection algorithm And positioning, form monitor area；

3) face tracking: monitor area corresponding to the multiple facial images arrived using Face tracking algorithm to continuous acquisition is carried out The movement of tracking face, obtains multiple facial images with motion change；

4) motion artifacts interference is eliminated: eliminating multiple facial images with motion change obtained with weighted mean method Motion artifacts interference, realize dynamic human face identification；

5) ROI region is extracted: to the facial image after recognition of face using depth information come adaptive selection ROI region；

6) rate calculation is carried out using IPPG technology: source extraction is carried out to ROI region using IPPG technology, utilizes experience point Solution mode EMD method carries out signal denoising and carries out the heart rate value that energy spectrum analysis seeks detected person using arma modeling.

2. the non-contact heart rate estimation technique of fusion IPPG and depth information anti-noise jamming according to claim 1, special Sign is: RGBD camera described in step 1) includes RGB camera and depth camera, and the RGB camera is for shooting The RGB image of human body, depth camera are used to shoot the depth image of face, the RGB camera and depth camera calibration Afterwards, the depth image of depth camera shooting and the RGB image of RGB camera shooting are overlapped.

3. the non-contact heart rate estimation technique of fusion IPPG and depth information anti-noise jamming according to claim 1, special Sign is: the method for step 2) facial image detection are as follows: using collected facial image as input picture, searches for input picture The rectangular area of middle possibility scale and position, i.e., candidate window to be detected carry out each candidate window in accordance with the following steps Processing: coarse sizing first is carried out using eyes template matching, then mean square deviation standardization is carried out to image in window, makes human face region point Cloth area normalization eliminates the influence of illumination variation, then carries out face template matching, if being not above the threshold value of setting, It is exported as candidate face.

4. the non-contact heart rate estimation technique of fusion IPPG and depth information anti-noise jamming according to claim 1, special Sign is: the method for step 3) face tracking are as follows: assuming that only occurring a face in sequence in each frame, sets initial surveillance zone Domain is R₀, face first is detected in the monitor area of each frame, if detecting face, is monitored according to its positioning result according to solution The formula in region calculates monitor area new in next frame, so handles each frame；If detection failure, keeps monitor area not Become, and continue to test in subsequent frames, to prevent from tracking failure as caused by accidental missing inspection；When missing inspection frame is more than a fixed number Amount, i.e. missing inspection time-out, it is believed that face has disappeared, then restarts to be likely to occur in detection subsequent frame in initial monitor area New face；

It is assumed that monitor area is expressed as LR=(x with one hexa-atomic group_min,x_max,y_min,y_max,ch_min,ch_max), x_min,x_maxIt is Range of the face center in the direction x, y_min,y_maxIt is range of the face center in the direction y, ch_min,ch_maxIt is face The range of scale；It is approximately a square by human face region, is shown as L=(x with a triple table_center,x_center, ch_face), each element successively describes center and the size of face；Then, the formula for solving monitor area is as follows:

Wherein, a, b, which are respectively indicated, defines face x between two frames, y location, and c indicates the maximum variation of face scale.

5. the non-contact heart rate estimation technique of fusion IPPG and depth information anti-noise jamming according to claim 1, special Sign is: step 4) motion artifacts interference elimination method are as follows: is swung left and right the variation of angle according to face, defines every width face The weight calculation method of variation proposes the construction strategy of weighted average face, the pitch angle variation of face is then based on, by it It is divided into vertical view, look squarely and looks up three levels, and constructs weighted mean face in each level, forms weighted mean face Matrix；

It is assumed that given motion change Gray Face image is I_j(x, y), then, the formula of weighted mean face are as follows:

Wherein, ZM indicates the different motion total number of images of same people, and j represents jth width modified-image, and x, y indicate that two dimensional image is flat Areal coordinate, ω_jIndicate the weighted value of corresponding jth width facial image；

ω_jDetermination and face the movement angle that is swung left and right it is related, ω_jCalculation method be: assuming that face be swung left and right before it is double Eye coordinates are respectively A (x_A,y_A) and B (x_B,y_B), it is respectively A'(x that the eyes coordinate after waving to the right, which occurs, for face_A',y_A') and B' (x_B',y_B'), the central coordinate of circle of circle is O (x where face cross-sectional view_o,y_o), the angle of oscillation of face to the right is ∠ δ, then K width Facial image weighted value ω_KCalculation formula it is as follows:

Wherein, δ_kIndicate the angle that is swung left and right of corresponding kth width facial image, then weight ω_kSolution be converted to calculating angle of oscillation Spend δ_k；

Eyes coordinate A (x_A,y_A)、B(x_B,y_B) and its eyes coordinate A'(x after being swung left and right_A',y_A')、B'(x_B',y_B') Determine position；Assuming that the radius of circle is R where eyes rotation, then

Centre point O (the x_o,y_o) coordinate obtained by following calculating process:

(6) (7) two formula above simultaneous, abbreviation arrangement can obtain:

It can similarly obtain:

Centre point o (x can be obtained by (8) (9) two formula simultaneous above_o,y_o) coordinate, bring the result into formula (6), find out half Diameter R then can solve swing angle δ by formula (5)；Angle of oscillation δ is the function of eyes coordinate value, i.e.,

δ=f (x_A,x_B,x_A',x_B') (10)

Assuming that the angle δ that is swung left and right, pitch angle γ, the motion change range that face is made divides and its corresponding weighted mean face Relation table；The Gray Face image I of given motion change_i(x, y), the then angle that is swung left and right are the weighted mean in 0 ° -90 ° Face definition is as shown in formula (11):

Wherein, ZM represents a certain pitch angle, a certain face sum being swung left and right in angular range, and j indicates jth width face, X, y represent two dimensional image plane coordinate, and p, q correspond to the mean value face number of the variation range, F '_{P, q}Expression is swung left and right Weighted mean face of the angle between -90 ° -0 °；

The weighted mean face being calculated, ultimately forming face, there are the weighted mean face squares under multi-motion situation of change Battle array:

6. the non-contact heart rate estimation technique of fusion IPPG and depth information anti-noise jamming according to claim 1, special Sign is: the method that step 5) ROI region is extracted are as follows: to the facial image after recognition of face from top edge 45%-75% down Region draw parallel lines, take 30% height as ROI of its image whole height, the width of ROI takes picture traverse after identification 70%, obtain only include nose and left and right cheek ROI region.

7. the non-contact heart rate estimation technique of fusion IPPG and depth information anti-noise jamming according to claim 1, special Sign is: the method that step 6) carries out rate calculation using IPPG technology are as follows: source extraction is carried out to ROI region, by the area ROI Area image is divided into tri- channels R, G, B, using the smallest channel G of noise as the channel of source signal, to the G channel image after separation Carry out space pixel average treatment:

Wherein k is number of image frames, and K is the totalframes of image, and Z (k) is the one-dimensional source signal in the channel R, z_i,j(n) for pixel (i, J) in the channel G color intensity value, g, d are the height and width of image respectively；

Signal denoising is carried out using experience resolution model EMD method, its step are as follows:

(1) average value for taking original signal u (t) obtains signal m₁(t)；

(2) single order surplus hp is calculated₁(t)=u (t)-m₁(t), hp is checked₁(t) whether meet IMF condition, if it is not, then returning to Step (1), uses hp₁(t) as the original signal of programmed screening, it may be assumed that

hp₂(t)=hp₁(t)-m₂(t) (17)

Repeat screening k times:

hp_k(t)=hp_k-1(t)-m_k(t) (18)

In hp_k(t) before meeting IMF condition, the one-component IMF of IMF is obtained₁, i.e.,

IMF₁=hp_k(t) (19)

(3) original signal u (t) subtracts IMF₁Surplus r can be obtained₁(t), i.e.,

r₁(t)=u (t)-IMF₁ (20)

(4) u is enabled₁(t)=r₁(t), by u₁(t) as new original signal, above-mentioned steps is repeated and obtain second IMF points Measure IMF₂, such repeatedly n times:

(5) as n-th of component r_n(t) monotonic function is had become, when can not decompose IMF again, the decomposable process of entire EMD is completed；It is former Beginning signal u (t) can be expressed as n IMF component and an average tendency component r_n(t) combination, it may be assumed that

The signal that its frequency is located at 0.75-2.0 Hz in heartbeat frequency band after being decomposed with arma modeling to EMD does energy spectrum analysis, i.e., Corresponding 45-120 times/min of human normal heart rate range, the corresponding frequency in energy highest point is palmic rate f_h, heart rate are as follows:

XLV=60f_h (22)。