Octave heartrate found

README.md

@@ -1,3 +1,11 @@
 # heatrateestimation
 
-This file contains our work of the heart rate estimation of a person from a video.
+This file contains our work of the heart rate estimation of a person from a video.
+
+## Python
+
+The python will deal with the image extraction and will have to be run each time for new videos.
+
+## Octave
+
+The octave file will work on the bulk of the project to process the data.

frame_extraction.py

@@ -1,4 +1,3 @@
-############
 # frame extraction
 # 
 # Task:
@@ -16,4 +15,21 @@
 # Author: Thomas Périn
 # Date: 17/02/2023
 #
+# Note: check video name when extracting a new one
 ############
+
+import cv2
+
+#Select video by name
+vidcap = cv2.VideoCapture('face_video.mp4')
+success,image = vidcap.read()
+
+#Counter to number the images
+count = 0
+
+#Extraction loop for all images
+while success:
+	cv2.imwrite("img%d.png" % count, image)     # save frame as JPEG file      
+	success,image = vidcap.read()
+	print('Read a new frame: ', success)
+	count += 1

heartRateEstimation.m

@@ -1,53 +1,51 @@
 %hr = heartRateEstimation(imgDirectory, fps, windowDuration, windowShift)
 %%%%%%%%%%%%%%%%%
 % hr = heartRateEstimation(imgDirectory, fps, windowDuration, windowShift)
-% 
+%
-% Task:
+% Task: Finding the heartrate of person from selcted images
 %
 % Inputs:
 %	-imgDirectory
-% 	-fps
+% -fps
-%	-windowDuration
-%	-windowShift
 %
 % Output:
-%	-hr:
+%	-heartrate:
 %
-% Author: Guillaume Gibert
+% Author: Thomas Périn
-% Date: 06/02/2023
+% Date: 17/02/2023
 %
 % Note: images were extracted with ffmpeg -i "video.mkv" ../img/gg_%04d.png
 %%%%%%%%%%%%%%%%%
 
-%temp var
+clc
-imgDirectory = 'img';
+clear all
-fps = 60;
-windowDuration = 30;
-windowShift = 1;
-roi = [225 750 825 1100];
 
-% global 
+% temp var
+fps = 15;
+windowDuration = 30;
+roi = [280 580 600 850];
+
+% global
 redChannel = [];
 greenChannel = [];
 blueChannel = [];
 
-%list the images
+%list the images with i the same directory
-listImg = dir([imgDirectory '/*.png']);
+listImg = dir(['*.png']);
 
+% iteration for every images available
 for l_img=1:windowDuration*fps
 	% load the current image
-	image_original = imread([imgDirectory '/' listImg(l_img).name]);
+  image_original = imread([listImg(l_img).name]);
-	%image(image_original);
+
-	
 	% crop the current image around the face
 	image_face = image_original(roi(1):roi(2), roi(3):roi(4), :);
-	%image(image_face);
+
-	
 	% spatial average for r, g, b channels
 	r_mean = mean(mean(image_face(:,:,1)));
 	g_mean = mean(mean(image_face(:,:,2)));
 	b_mean = mean(mean(image_face(:,:,3)));
-	
+
 	% store the current "average" pixel in global arrays
 	redChannel = [redChannel r_mean];
 	greenChannel = [greenChannel g_mean];
@@ -55,26 +53,39 @@ for l_img=1:windowDuration*fps
 end
 
 % estimate temporal average and standard deviation
+redChannel_avg = mean(redChannel);
+redChannel_std = std(redChannel);
 greenChannel_avg = mean(greenChannel);
 greenChannel_std = std(greenChannel);
+blueChannel_avg = mean(blueChannel);
+blueChannel_std = std(blueChannel);
 
-% normalize your data 
+% normalize your data
-% greenChannel_normalized(i) = (greenChannel(i) - greenChannel_avg)/greenChannel_std
+% first define length of our data
+L = columns(greenChannel);
+for i=1:L
+  greenChannel_normalized(i) = (greenChannel(i) - greenChannel_avg)/greenChannel_std;
+end
 
-% fft(greenChannel_normalized)
+% select only greenChannel for simplicity (no ICA)
+greenChannel_fft = fft(greenChannel_normalized);
 
-% power spectrum (https://www.mathworks.com/help/matlab/ref/fft.html)
+% power spectrum using link(https://www.mathworks.com/help/matlab/ref/fft.html)
+P2 = abs(greenChannel_fft/L);
+P1 = P2(1:L/2+1);
+P1(2:end-1) = 2*P1(2:end-1);
+f = 15*(0:(L/2))/L;
+plot(f,P1)
 
 % find the peak in the range ([0.75 4] Hz)
-% max -> value, index
+% define range and round to get intergers
-
+low_index = round(0.75*L/15);
-% convert the index from Hz to bpm
+high_index = round(2*L/15);
+% find max value and index inn new range
+[a, ia] = max(P1(low_index:high_index));
+% adapt the index for full range
+heart_frequency = f(ia+low_index)
 
 % determine heart rate
-
+% convert the index from Hz to bpm
-
+heart_rate = heart_frequency*60
-
-
-
-
-