Revision 11:58a443cf3f7c, committed 2019-03-16

Comitter:

jeannie9809

Date:

Sat Mar 16 07:46:28 2019 +0000

Parent:

10:3ddf92687e34

Child:

12:95d141421e18

Commit message:

FINALLY DONE WITH HRV

Changed in this revision

--- a/8cee5929f4d8/main.cpp	Sat Mar 16 06:13:01 2019 +0000
+++ b/8cee5929f4d8/main.cpp	Sat Mar 16 07:46:28 2019 +0000
@@ -7,7 +7,11 @@
 #include <math.h>
 #include <vector>
 using namespace std;
+
 #define FIFO_DATA_MAX 288
+//#define M_PI 3.14159265358979323846
+#define M_PI 3.14159
+
 void UpdateSensorData(void);
 void StartHaptic(void);
 void StopHaptic(void const *n);
@@ -37,11 +41,9 @@
 uint8_t testsignal = 60;
 
 // I added this
-const int num_samples = 1200; // for 30 sec
-double TIME[num_samples];
+const int num_samples = 500; // for 30 sec
 int ppg[num_samples];
 int SDNN_n, SDNN;
-bool first_sample_set = true;
 double arousal, valence, HF_LF, HF_LF_n;
 bool ready = false;
 
@@ -91,8 +93,8 @@
         /*Send Temperature at 25 degrees Celsius
         kw40z_device.SendTemperature(temperature);
 
-        /*Send Pressure at 100kPA */
-        //kw40z_device.SendPressure(pressure);
+        //Send Pressure at 100kPA 
+        //kw40z_device.SendPressure(pressure); */
 
         /*Send Mag,Accel,Gyro Data.
         kw40z_device.SendGyro(x,y,z);
@@ -197,7 +199,7 @@
 
                 num = *(uint32_t *) sample;
                 
-                if (num < 310000) {
+                if (num < 0) {
                     ppg_single_sample = 0;
 //                    printf("keep closer to your hand \r\n");
                 } else {
@@ -221,7 +223,7 @@
         }
 
         interruptStatus.all = 0xFF;
-        
+
         if (mask_ppg == 1) {
             interruptStatus.bits.ppg_rdy = 0;
         }
@@ -233,159 +235,188 @@
     }
     */
     //ppg = temp_ppg;
-    
+
 }
 
 void interruptHandler()
 {
     evqueue.call(interruptHandlerQueued);
 //    interruptHandlerQueued();
+//    printf("\nLezzgo\n");
     if(ready) {
-    printf("\nStarting...  ");
-    
-    int i = 0;
-    int j = 0;
-    
-    // moving average
-    double movave[num_samples];
-    int avecap = 25;
-    for(i = 0; i < (avecap-1)/2; i++) {
-        movave[i] = (double)ppg[0];
-        for(j = 1; j < 1+(avecap-1)/2; j++) {
-            movave[i] = movave[i] + (double)ppg[j];
+        printf("\nStarting...  ");
+
+        int i = 0;
+        int j = 0;
+        
+        double TIME[num_samples];
+        for(int i = 0; i < num_samples; i++) {
+            TIME[i] = (double)i*(5.0/(double)num_samples); // change to 30.0 later
+        }
+        
+        // moving average
+        double movave[num_samples];
+        int avecap = 25;
+        for(i = 0; i < (avecap-1)/2; i++) {
+            movave[i] = (double)ppg[0];
+            for(j = 1; j < 1+(avecap-1)/2; j++) {
+                movave[i] = movave[i] + (double)ppg[j];
+            }
+            movave[i] = (double)(movave[i]/(i+(avecap-1)/2-1));
+        }
+        for(i = (num_samples-(avecap-1)/2); i < num_samples; i++) {
+            movave[i] = (double)ppg[i-(avecap-1)/2-1];
+            for(j = (i-(avecap-1)/2); j < num_samples; j++) {
+                movave[i] = movave[i] + (double)ppg[j];
+            }
+            movave[i] = movave[i]/(double)(num_samples-i+(avecap-1)/2);
+        }
+        for(i = (avecap-1)/2; i < (num_samples-(avecap-1)/2); i++) {
+            movave[i] = (double)ppg[i-(avecap-1)/2-1];
+            for(j = i-(avecap-1)/2; j < i+(avecap-1)/2; j++) {
+                movave[i] = movave[i] + (double)ppg[j];
+            }
+            movave[i] = movave[i]/(double)avecap;
         }
-        movave[i] = (double)(movave[i]/(i+(avecap-1)/2-1));
-    }
-    for(i = (num_samples-(avecap-1)/2); i < num_samples; i++) {
-        movave[i] = (double)ppg[i-(avecap-1)/2-1];
-        for(j = (i-(avecap-1)/2); j < num_samples; j++) {
-            movave[i] = movave[i] + (double)ppg[j];
+
+        // normalize ppg
+        for(i = 0; i < num_samples; i++) {
+            ppg[i] = ppg[i] - (int)movave[i];
+            if(ppg[i] > 1000 || ppg[i] < -1000)
+                ppg[i] = 0;
+        }
+
+        // smoothing curve
+        for(i = 1; i < num_samples; i++) {
+            ppg[i] = ppg[i] + ppg[i-1];
+            //printf("%d   ", ppg[i]);
         }
-        movave[i] = (double)(movave[i]/(num_samples-i+(avecap-1)/2));
-    }
-    for(i = (avecap-1)/2; i < (num_samples-(avecap-1)/2); i++) {
-        movave[i] = (double)ppg[i-(avecap-1)/2-1];
-        for(j = i-(avecap-1)/2; j < i+(avecap-1)/2; j++) {
-            movave[i] = movave[i] + (double)ppg[j];
+
+        // AMPD Algorithm
+        const int kcap = 25;
+        int m[kcap][num_samples];
+        for(int k = 1; k <= kcap; k++) {
+            for(i = 1; i < num_samples; i++) {
+                if(i-1 < 0 || i-k-1 < 0 || i+k-1 >= num_samples)
+                    m[k-1][i] = 0;
+                else if(ppg[i-1] > ppg[i-k-1] && ppg[i-1] > ppg[i+k-1])
+                    m[k-1][i] = 1;
+                else
+                    m[k-1][i] = 0;
+            }
+        }
+        
+        int max_mult[num_samples];
+        for(i = 0; i < num_samples; i++) { // max_mult = m(1,:)';
+            max_mult[i] = m[0][i];
+        }
+
+        for(int k = 1; k < kcap; k++) {
+            for(i = 0; i < num_samples; i++) {
+                max_mult[i] = max_mult[i]*m[k][i];
+            }
+        }
+        
+        int num_max = 0;
+        // extract times that are max
+        for(i = 0; i < num_samples; i++) {
+            num_max = num_max + max_mult[i];
         }
-        movave[i] = (double)(movave[i]/avecap);
-    }
-    
-    // normalize ppg
-    for(i = 0; i < num_samples; i++) {
-        ppg[i] = ppg[i] - (int)movave[i];
-        if(ppg[i] > 1000 || ppg[i] < -1000)
-            ppg[i] = 0;
-    }
-    
-    // smoothing curve
-    for(i = 1; i < num_samples; i++) {
-        ppg[i] = ppg[i] + ppg[i-1];
-        //printf("%d   ", ppg[i]);
-    }
-    
-    // AMPD Algorithm
-    const int kcap = 25;
-    int m[kcap][num_samples];
-    for(int k = 1; k <= kcap; k++) {
-        for(i = 1; i < num_samples; i++) {
-            if(i-1 < 0 || i-k-1 < 0 || i+k-1 >= num_samples)
-                m[k-1][i] = 0;
-            else if(ppg[i-1] > ppg[i-k-1] && ppg[i-1] > ppg[i+k-1])
-                m[k-1][i] = 1;
-            else
-                m[k-1][i] = 0;
+
+//    printf("num_max = %d    ", num_max);
+
+        vector<double> time_of_max;
+        vector<int> index_of_max;
+        vector<int> max_points;
+        for(i = 0; i < num_samples; i++) {
+            if(max_mult[i] == 1) {
+                time_of_max.push_back(TIME[i-1]);
+                index_of_max.push_back(i-1);
+                max_points.push_back(ppg[i-1]);
+            }
+        }
+
+        // calculating HRV
+        vector<double> r;
+        vector<int> index_r;
+        double mean_inter_time = 0;
+        for(i = 0; i < num_max-1; i++) {
+            r.push_back(time_of_max.at(i+1)-time_of_max.at(i));
+            index_r.push_back(index_of_max.at(i+1) - index_of_max.at(i));
+            mean_inter_time = mean_inter_time + r.at(i);
+        }
+
+        //printf("r: %.2f %.2f %.2f   ", r.at(0), r.at(1), r.at(2));
+        mean_inter_time = mean_inter_time/(double)(num_max-1);
+
+        // getting rid of outlier points in r
+        for(i = 0; i < num_max-1; i++) {
+            if(r.at(i) > mean_inter_time + 0.11)
+                r.at(i) = mean_inter_time + 0.11;
+            else if(r.at(i) < mean_inter_time - 0.11)
+                r.at(i) = mean_inter_time - 0.11;
+        }
+
+        // SDNN -- std of normal to normal R-R intervals
+        mean_inter_time = 0;
+        for(i = 0; i < num_max-1; i++) {
+            mean_inter_time = mean_inter_time + r.at(i);
         }
-    }
-    int max_mult[num_samples];
-    for(i = 0; i < num_samples; i++) { // max_mult = m(1,:)';
-        max_mult[i] = m[0][i];
-    }
-    
-    for(int k = 1; k < kcap; k++) {
+        mean_inter_time = double(mean_inter_time/(num_max-1));
+        double SDNN_doub = 0.0;
+        for(i = 0; i < num_max-1; i++) {
+            SDNN_doub = SDNN_doub + (r.at(i)-mean_inter_time)*(r.at(i)-mean_inter_time);
+        }
+        SDNN = (int)(sqrt(SDNN_doub/(num_max-1))*1000);
+        printf("SDNN = %d    ", SDNN);
+
+        // TIME TO CALCULATE HF/LF
+        // FFT: use movave as fftppg 
         for(i = 0; i < num_samples; i++) {
-            max_mult[i] = max_mult[i]*m[k][i];
+            double real = 0;
+            double im = 0;
+            double dum_ppg;
+            for(j = 0; j < num_samples; j++) {
+                dum_ppg = (double) ppg[j];
+                real = real + dum_ppg * cos(2.0*M_PI*(double)(j*i)/(double)num_samples);
+                im = im + dum_ppg * sin(2.0*M_PI*(double)(j*i)/(double)num_samples);
+            }
+            movave[i] = sqrt(real*real + im*im)/(double)num_samples;
+        }
+        
+        
+        // make frequency array: use TIME
+        double sampling_freq = (double)(1/TIME[1]);
+        for(i = 0; i < num_samples; i++) {
+            TIME[i] = sampling_freq*(double)(i)/(double)(num_samples);
         }
+
+        double LF = 0.0;
+        i = 0;
+        while(TIME[i] < 0.15) {
+            LF = LF + movave[i];
+            i++;
+        }
+
+        double HF = 0.0;
+        while(TIME[i] < 0.4) {
+            HF = HF + movave[i];
+            i++;
+        }
+        HF_LF = HF/LF;
+        printf("HF/LF = %.2f    \n", HF_LF);
+        
+        
+        ready = false; // last line
+
     }
-    
-    // extract times that are max
-    int num_max = 0;
-    for(i = 0; i < num_samples; i++) {
-        num_max = num_max + max_mult[i];
-    }
-//    printf("num_max = %d    ", num_max);
-    
-    vector<double> time_of_max;
-    vector<int> index_of_max;
-    vector<int> max_points;
-    for(i = 0; i < num_samples; i++) {
-        if(max_mult[i] == 1) {
-            time_of_max.push_back(TIME[i-1]);
-            index_of_max.push_back(i-1);
-            max_points.push_back(ppg[i-1]);
-        }
-    } 
-    
-    // calculating HRV
-    vector<double> r;
-    vector<int> index_r;
-    double mean_inter_time = 0;
-    for(i = 0; i < num_max-1; i++) {
-        r.push_back(time_of_max.at(i+1)-time_of_max.at(i));
-        index_r.push_back(index_of_max.at(i+1) - index_of_max.at(i));
-        mean_inter_time = mean_inter_time + r.at(i);
-    }
-    
-    //printf("r: %.2f %.2f %.2f   ", r.at(0), r.at(1), r.at(2));
-    mean_inter_time = (double)(mean_inter_time/(num_max-1));
-    
-    // getting rid of outlier points in r
-    for(i = 0; i < num_max-1; i++) {
-        if(r.at(i) > mean_inter_time + 0.11)
-            r.at(i) = mean_inter_time + 0.11;
-        else if(r.at(i) < mean_inter_time - 0.11)
-            r.at(i) = mean_inter_time - 0.11;
-    }
-    
-    // SDNN -- std of normal to normal R-R intervals
-    mean_inter_time = 0;
-    for(i = 0; i < num_max-1; i++) {
-        mean_inter_time = mean_inter_time + r.at(i);
-    }
-    mean_inter_time = double(mean_inter_time/(num_max-1));
-    double SDNN_doub = 0.0;
-    for(i = 0; i < num_max-1; i++) {
-        SDNN_doub = SDNN_doub + (r.at(i)-mean_inter_time)*(r.at(i)-mean_inter_time);
-    }
-    SDNN = (int)(sqrt(SDNN_doub/(num_max-1))*1000);
-//    printf("SDNN = %d\r\n", SDNN);
-    
-    // TIME TO CALCULATE HF/LF
-    // FFT
-    /*
-    double fftppg[num_samples];
-    for(i = 0; i < num_samples; i++) {
-        double real = 0;
-        double im = 0;
-        for(j = 0; j < num_samples; j++) {
-            real = real + ppg[j] * cos((double)(2*M_PI)*(double)((double)(j*i)/num_samples));
-            im = im + ppg[j] * sin((double)(2*M_PI)*(double)((double)(j*i)/num_samples));
-        }
-        fftppg[i] = sqrt
-    }
-    */
-    ready = false; // last line
-    }
-    
+
 }
 
 // main() runs in its own thread in the OS
 int main()
 {
 //    printf("Hello world.\r\n");
-    for(int i = 0; i < num_samples; i++) {
-        TIME[i] = (double)(i*(12.0/num_samples)); // change to 30.0 later
-    }
     
     t.start(callback(&evqueue, &EventQueue::dispatch_forever));
     kw40z_device.attach_buttonLeft(&ButtonLeft);
@@ -413,7 +444,7 @@
     strcpy((char *) text, "Raw PPG:");
     oled.Label((uint8_t *)text,7,0);
 
-    strcpy((char *) text, "SDNN; LF/HF:");
+    strcpy((char *) text, "           LOFI");
     oled.Label((uint8_t *)text,7,40);
     //dynamic text setup
     textProperties.fontColor = COLOR_WHITE;
@@ -421,24 +452,29 @@
     oled.SetTextProperties(&textProperties);
 
     txThread.start(txTask);
+    
+
     while (true) {
 
-        /* Format the time reading */
+        // Format the time reading 
         sprintf(text,"%d",ppg_single_sample);
-
-        /* Display time reading in 35px by 15px textbox at(x=55, y=40) */
         oled.TextBox((uint8_t *)text,55,15,35,15); //Increase textbox for more digits
+        
+/*
+        // Display time reading in 35px by 15px textbox at(x=55, y=40) 
         if (ppg_single_sample > 45) {
             sprintf(text,"%d; %.2f",SDNN, HF_LF);
 
-            /* Display time reading in 35px by 15px textbox at(x=55, y=40) */
+            // Display time reading in 35px by 15px textbox at(x=55, y=40) 
             oled.TextBox((uint8_t *)text,55,55,35,15); //Increase textbox for more digits
         } else {
             sprintf(text,"wait HR");
 
-            /* Display time reading in 35px by 15px textbox at(x=55, y=40) */
+            // Display time reading in 35px by 15px textbox at(x=55, y=40) 
             oled.TextBox((uint8_t *)text,55,55,35,15);
-        }
+        }*/
+
+        
 
         Thread::wait(1000);
     }