Revision 0:805cf4bfa472, committed 2014-12-17

Comitter:

nagasm

Date:

Wed Dec 17 02:02:28 2014 +0000

Child:

1:7bccad3c277f

Commit message:

NucleoF401RE EMG sensor with Artifact Filter, Notch Filter, etc

Changed in this revision

--- /dev/null	Thu Jan 01 00:00:00 1970 +0000
+++ b/FIR_LPF.hpp	Wed Dec 17 02:02:28 2014 +0000
@@ -0,0 +1,72 @@
+// from CQpub0 Mikami / FIR_LPF_Direct
+//  http://developer.mbed.org/users/CQpub0Mikami/code/FIR_LPF_Direct
+//
+//                                   Band1       Band2
+// Lower band edge frequency (kHz) 0.000000    0.600000
+// Upper band edge frequency (kHz) 0.500000    5.000000
+// Gain                            1.000000    0.000000
+// Weight                          1.000000    1.000000
+// Deviation                       0.009734    0.009734
+// Deviation [dB]                  0.084139  -40.234188
+
+const float fc[order+1] = {   
+    -3.566292E-03f,  2.335185E-03f,  1.917338E-03f,  1.681921E-03f,  1.522689E-03f,
+     1.364066E-03f,  1.157961E-03f,  8.803014E-04f,  5.296940E-04f,  1.236180E-04f,
+    -3.047488E-04f, -7.123405E-04f, -1.051992E-03f, -1.279654E-03f, -1.360180E-03f,
+    -1.273725E-03f, -1.018960E-03f, -6.155554E-04f, -1.029986E-04f,  4.621914E-04f,
+     1.013316E-03f,  1.480161E-03f,  1.798055E-03f,  1.916527E-03f,  1.805744E-03f,
+     1.463100E-03f,  9.143371E-04f,  2.130115E-04f, -5.640853E-04f, -1.326518E-03f,
+    -1.978121E-03f, -2.430574E-03f, -2.613837E-03f, -2.486610E-03f, -2.043447E-03f,
+    -1.317959E-03f, -3.802110E-04f,  6.680631E-04f,  1.704594E-03f,  2.600376E-03f,
+     3.235972E-03f,  3.515095E-03f,  3.379365E-03f,  2.817894E-03f,  1.872141E-03f,
+     6.322899E-04f, -7.684087E-04f, -2.167773E-03f, -3.393601E-03f, -4.284042E-03f,
+    -4.705797E-03f, -4.577073E-03f, -3.875274E-03f, -2.648668E-03f, -1.012065E-03f,
+     8.616354E-04f,  2.758435E-03f,  4.447251E-03f,  5.706155E-03f,  6.350451E-03f,
+     6.257734E-03f,  5.386359E-03f,  3.786968E-03f,  1.602121E-03f, -9.435526E-04f,
+    -3.565645E-03f, -5.948320E-03f, -7.781761E-03f, -8.797770E-03f, -8.807067E-03f,
+    -7.725336E-03f, -5.594939E-03f, -2.585050E-03f,  1.010693E-03f,  4.809328E-03f,
+     8.361672E-03f,  1.121444E-02f,  1.295107E-02f,  1.324389E-02f,  1.191081E-02f,
+     8.934171E-03f,  4.488148E-03f, -1.062401E-03f, -7.182842E-03f, -1.321242E-02f,
+    -1.841506E-02f, -2.204769E-02f, -2.343209E-02f, -2.202840E-02f, -1.749386E-02f,
+    -9.730157E-03f,  1.092382E-03f,  1.453105E-02f,  2.989524E-02f,  4.629194E-02f,
+     6.269474E-02f,  7.802411E-02f,  9.123583E-02f,  1.014078E-01f,  1.078185E-01f,
+     1.100078E-01f,  1.078185E-01f,  1.014078E-01f,  9.123583E-02f,  7.802411E-02f,
+     6.269474E-02f,  4.629194E-02f,  2.989524E-02f,  1.453105E-02f,  1.092382E-03f,
+    -9.730157E-03f, -1.749386E-02f, -2.202840E-02f, -2.343209E-02f, -2.204769E-02f,
+    -1.841506E-02f, -1.321242E-02f, -7.182842E-03f, -1.062401E-03f,  4.488148E-03f,
+     8.934171E-03f,  1.191081E-02f,  1.324389E-02f,  1.295107E-02f,  1.121444E-02f,
+     8.361672E-03f,  4.809328E-03f,  1.010693E-03f, -2.585050E-03f, -5.594939E-03f,
+    -7.725336E-03f, -8.807067E-03f, -8.797770E-03f, -7.781761E-03f, -5.948320E-03f,
+    -3.565645E-03f, -9.435526E-04f,  1.602121E-03f,  3.786968E-03f,  5.386359E-03f,
+     6.257734E-03f,  6.350451E-03f,  5.706155E-03f,  4.447251E-03f,  2.758435E-03f,
+     8.616354E-04f, -1.012065E-03f, -2.648668E-03f, -3.875274E-03f, -4.577073E-03f,
+    -4.705797E-03f, -4.284042E-03f, -3.393601E-03f, -2.167773E-03f, -7.684087E-04f,
+     6.322899E-04f,  1.872141E-03f,  2.817894E-03f,  3.379365E-03f,  3.515095E-03f,
+     3.235972E-03f,  2.600376E-03f,  1.704594E-03f,  6.680631E-04f, -3.802110E-04f,
+    -1.317959E-03f, -2.043447E-03f, -2.486610E-03f, -2.613837E-03f, -2.430574E-03f,
+    -1.978121E-03f, -1.326518E-03f, -5.640853E-04f,  2.130115E-04f,  9.143371E-04f,
+     1.463100E-03f,  1.805744E-03f,  1.916527E-03f,  1.798055E-03f,  1.480161E-03f,
+     1.013316E-03f,  4.621914E-04f, -1.029986E-04f, -6.155554E-04f, -1.018960E-03f,
+    -1.273725E-03f, -1.360180E-03f, -1.279654E-03f, -1.051992E-03f, -7.123405E-04f,
+    -3.047488E-04f,  1.236180E-04f,  5.296940E-04f,  8.803014E-04f,  1.157961E-03f,
+     1.364066E-03f,  1.522689E-03f,  1.681921E-03f,  1.917338E-03f,  2.335185E-03f,
+    -3.566292E-03f
+};
+
+float FIR_calc1(float data){
+    int i;
+    float acc = 0;
+    xn1[0] = data;
+    for (i=0; i<=order; i++) acc = acc + fc[i]*xn1[i];
+    for (i=order; i>0; i--) xn1[i] = xn1[i-1];
+    return(acc);
+}
+
+float FIR_calc2(float data){
+    int i;
+    float acc = 0;
+    xn2[0] = data;
+    for (i=0; i<=order; i++) acc = acc + fc[i]*xn2[i];
+    for (i=order; i>0; i--) xn2[i] = xn2[i-1];
+    return(acc);
+}

--- /dev/null	Thu Jan 01 00:00:00 1970 +0000
+++ b/main.cpp	Wed Dec 17 02:02:28 2014 +0000
@@ -0,0 +1,90 @@
+#include "mbed.h"
+#include "sub.hpp"
+#include "FIR_LPF.hpp"
+#include "notch.hpp"
+
+int main(){
+    common_setup();
+    xbee.baud(38400);
+    xbee.attach(&rx_fifoset, xbee.RxIrq);
+    timer_setup.attach_us(&timer_interrupt, 5); // 5usec
+    while(1){
+        tx_fifo_check();
+        if(timer_value[1] > 19){ // 0.1msec sampling
+            timer_value[1] = 0;
+            float data1 = (float)gain * (analog_value2.read() - 0.5f);  // A/D in (3)
+            if(emergence!=0 && data1<0.05 && data1>0.95){
+                data1 = old_1;
+                if(++em_count1 > 10){
+                    em_count1 = 0;
+                    tx_message(0xff0000);   // Emergency 1
+                }
+            }
+            if(data1 < 0) data1 = -data1; // always detection ON
+            if(fir_lpf != 0) data1 = FIR_calc1(data1);      // FIR calc (1) call
+            if(notch_1 != 0) data1 = notch_filter1(data1);
+            old_1 = data1;
+            if(timer_value[2] > 2999){  // 15msec
+                timer_value[2] = 0;
+                if(max_count != 0) data1 = move_mean_calc1(data1);
+                tx_message((uint16_t)((data1 + 1.0f) * 2047)<<4);
+            }
+        }
+        if(timer_value[3] > 19){ // 0.1msec sampling
+            timer_value[3] = 0;
+            float data2 = (float)gain * (analog_value3.read() - 0.5f);  // A/D in (4)
+            if(emergence!=0 && data2<0.05 && data2>0.95){
+                data2 = old_2;
+                if(++em_count2 > 10){
+                    em_count2 = 0;
+                    tx_message(0xff0001);   // Emergency 2
+                }
+            }
+            if(data2 < 0) data2 = -data2; // always detection ON
+            if(fir_lpf != 0) data2 = FIR_calc2(data2);      // FIR calc (1) call
+            if(notch_2 != 0) data2 = notch_filter2(data2);
+            old_2 = data2;
+            if(timer_value[4] > 2999){  // 15msec
+                timer_value[4] = 0;
+                if(max_count != 0) data2 = move_mean_calc2(data2);
+                tx_message(0x400000 + ((uint16_t)((data2 + 1.0f) * 2047)<<4));
+            }
+        }
+        if(timer_value[0] > 199999){ // 1000msec
+            timer_value[0] = 0;
+            myled = !myled;
+        }
+        if(rx_fifo_check() == 1){
+            int sum = 0;
+            for (int i=0; i<6; i++) sum += conv_hex(raw_data[i])<<(4*(5-i));
+            tx_message(sum); // Echo Back 
+            if(sum>>16 == 0x80){
+                switch((sum & 0xff00)>>8){
+                    case 0x00:
+                        fir_lpf = sum & 0x01;
+                        break;
+                    case 0x01:
+                        gain = sum & 0x0f;
+                        break;
+                    case 0x02:
+                        max_count = sum & 0x7f;
+                        if (max_count>100) max_count = 100;
+                        sum_clear();
+                        break;
+                    case 0x03:
+                        notch_1 = sum & 0x01;
+                        break;
+                    case 0x04:
+                        notch_2 = sum & 0x01;
+                        break;
+                    case 0x05:
+                        coef_set(sum & 0x3f);
+                        break;
+                    case 0x06:
+                        emergence = sum & 0x01;
+                        break;
+                }
+            }
+        }        
+    }
+}

--- /dev/null	Thu Jan 01 00:00:00 1970 +0000
+++ b/mbed.bld	Wed Dec 17 02:02:28 2014 +0000
@@ -0,0 +1,1 @@
+http://mbed.org/users/mbed_official/code/mbed/builds/4fc01daae5a5
\ No newline at end of file

--- /dev/null	Thu Jan 01 00:00:00 1970 +0000
+++ b/notch.hpp	Wed Dec 17 02:02:28 2014 +0000
@@ -0,0 +1,47 @@
+/*
+#include <stdio.h>
+#include <stdlib.h>
+#include <math.h>
+int main(){
+    int fs = 10000;
+    double  B0, a1, a2, b1, c0, F0, T;
+    F0 = 60;    // or 50
+    B0 = 100.;
+    T  = 1./fs;
+    a1 = 2. * exp(-M_PI * B0 * T) * cos(2. * M_PI * F0 * T);
+    a2 = -exp(-2. * M_PI * B0 * T);
+    b1 = -2. * cos(2. * M_PI * F0 * T);
+    c0 = (1-a1-a2)/(2+b1);
+    printf("a1 = %f\n", a1);
+    printf("a2 = %f\n", a2);
+    printf("b1 = %f\n", b1);
+    printf("c0 = %f\n", c0);
+    return 0;
+}
+60Hz
+    a1 = 1.936768
+    a2 = -0.939101
+    b1 = -1.998579
+    c0 = 1.642174
+50Hz
+    a1 = 1.937188
+    a2 = -0.939101
+    b1 = -1.999013
+    c0 = 1.938304
+*/
+
+float notch_filter1(float data){
+    y1[0] = data + a1*y1[1] + a2*y1[2];
+    float reault = y1[0] + b1*y1[1] + y1[2];
+    y1[2] = y1[1];
+    y1[1] = y1[0];
+    return(reault);
+}
+
+float notch_filter2(float data){
+    y2[0] = data + a1*y2[1] + a2*y2[2];
+    float reault = y2[0] + b1*y2[1] + y2[2];
+    y2[2] = y2[1];
+    y2[1] = y2[0];
+    return(reault);
+}

--- /dev/null	Thu Jan 01 00:00:00 1970 +0000
+++ b/sub.hpp	Wed Dec 17 02:02:28 2014 +0000
@@ -0,0 +1,129 @@
+unsigned char rxFIFO[256], txFIFO[256], raw_data[6];
+unsigned char rx_top, rx_end, tx_top, tx_end, phase;
+float mean_sum1, mean_sum2, ad_data1[101], ad_data2[101];
+int timer_value[6], max_count, fir_lpf, gain, ad_pointer1, ad_pointer2;
+const int order = 200;
+float xn1[order+1], xn2[order+1];
+float a1, a2, b1, c0, y1[3], y2[3], old_1, old_2;
+int notch_1, notch_2, emergence, em_count1, em_count2;
+const float _60Hz_a1 = 1.936768, _60Hz_a2 = -0.939101, _60Hz_b1 = -1.998579, _60Hz_c0 = 1.642174;
+const float _50Hz_a1 = 1.937188, _50Hz_a2 = -0.939101, _50Hz_b1 = -1.999013, _50Hz_c0 = 1.938304;
+
+RawSerial xbee(PA_2, PA_3);
+Ticker timer_setup;
+AnalogIn analog_value0(A0);
+AnalogIn analog_value1(A1);
+AnalogIn analog_value2(A2);
+AnalogIn analog_value3(A3);
+DigitalOut myled(LED1);
+
+void sum_clear(){
+    int i;
+    for (i=0; i<101; i++) ad_data1[i] = ad_data2[i] = 0;
+    ad_pointer1 = ad_pointer2 = mean_sum1 = mean_sum2 = 0;
+}
+
+void coef_set(int herz){
+    if(herz < 55){
+        a1 = _50Hz_a1;
+        a2 = _50Hz_a2;
+        b1 = _50Hz_b1;
+        c0 = _50Hz_c0;
+    }
+    else{
+        a1 = _60Hz_a1;
+        a2 = _60Hz_a2;
+        b1 = _60Hz_b1;
+        c0 = _60Hz_c0;
+    }
+    for (int i=0; i<3; i++) y1[i] = y2[i] = 0;
+}
+void common_setup(){
+    int i;
+    rx_top = rx_end = tx_top = tx_end = phase = 0;
+    em_count1 = em_count2 = 0;
+    max_count = 100;
+    fir_lpf = emergence = old_1 = old_2 = 0;
+    gain = 5;
+    for (i=0; i<3; i++) timer_value[i] = 0;
+    timer_value[3] = 10;
+    timer_value[4] = 1500;
+    for (i=0; i<=order; i++) xn1[i] = xn2[i] = 0.0;
+    sum_clear();
+    notch_1 = notch_2 = 60;
+    coef_set(60);
+}
+
+void timer_interrupt(){
+    int i;
+    for (i=0; i<6; i++) timer_value[i]++;
+}
+
+void tx_fifo_check(){
+    if(xbee.writeable() == 1){
+        if(tx_top != tx_end){
+            xbee.putc(txFIFO[tx_end]);
+            ++tx_end &= 255;
+        }
+    }
+}
+
+int rx_fifo_check(){
+    unsigned char data;
+    if(rx_top != rx_end){
+        data = rxFIFO[rx_end];
+        ++rx_end &= 255;
+        if (data < 33){
+            phase = 0;
+            return(1);
+        }
+        raw_data[phase] = data;
+        if(++phase > 5) phase = 0;
+        return(0);
+    }
+    return(0);
+}
+
+void rx_fifoset(void){
+    rxFIFO[rx_top] = xbee.getc();
+    ++rx_top &= 255;
+}
+
+void tx_fifoset(unsigned char data){
+    txFIFO[tx_top] = data;
+    ++tx_top &= 255;
+}
+
+unsigned char hex_conv(unsigned char data){
+    data &= 15;
+    if(data < 10) return(data+48);
+    else return(data+55);
+}
+
+unsigned char conv_hex(unsigned char data){
+    if((data > 47) && (data < 58)) return(data-48);
+    else if((data > 64) && (data < 71)) return(data-55);
+    return(0);
+}
+
+void tx_message(int data){
+    int i;
+    for (i=0; i<6; i++) tx_fifoset(hex_conv((data>>(4*(5-i))) & 15));
+    tx_fifoset(13);
+}
+
+float move_mean_calc1(float data){
+    mean_sum1 = mean_sum1 - ad_data1[ad_pointer1] + data;
+    ad_data1[ad_pointer1] = data;
+    ad_pointer1++;
+    if(ad_pointer1 == max_count) ad_pointer1 = 0;
+    return(mean_sum1 / (float)max_count);
+}
+
+float move_mean_calc2(float data){
+    mean_sum2 = mean_sum2 - ad_data2[ad_pointer2] + data;
+    ad_data2[ad_pointer2] = data;
+    ad_pointer2++;
+    if(ad_pointer2 == max_count) ad_pointer2 = 0;
+    return(mean_sum2 / (float)max_count);
+}